Contributors: J. Wuite (ENVEO IT GmbH), T. Nagler (ENVEO IT GmbH)

Issued by: ENVEO / Jan Wuite

Date: 10/05/2023

Ref: C3S2_312a_Lot4.WP2-FDDP-IS-v1_202212_IV_ATBD-v4_i1.1

Official reference number service contract: 2021/C3S2_312a_Lot4_EODC/SC1

History of modifications

List of datasets covered by this document

Related documents

Acronyms

General definitions

Single Look Complex (SLC): Level-1 Single Look Complex (SLC) products are Synthetic Aperture Radar (SAR) images in the slant range by azimuth imaging plane, in the image plane of satellite data acquisition. Each image pixel is represented by a complex magnitude value and therefore contains both amplitude and phase information. The imagery is geo-referenced using orbit and attitude data from the satellite. SLC images are produced in a zero Doppler geometry.
(Adapted from: https://sentinels.copernicus.eu/)

Interferometric Wide (IW): The Interferometric Wide (IW) swath mode is the main acquisition mode of Sentinel-1 over land, including ice sheets. It acquires data with a 250 km swath at 5 m by 20 m spatial resolution. IW mode captures three sub-swaths using Terrain Observation with Progressive Scans SAR (TOPS). (Adapted from: https://sentinels.copernicus.eu/)

Terrain Observation with Progressive Scans (TOPS): With the TOPS(AR) technique, in addition to steering the radar beam in range, the beam is also electronically steered from backward to forward in the azimuth direction for each burst, avoiding scalloping and resulting in homogeneous image quality throughout the swath (De Zan and Guarnieri, 2006)

Offset Tracking (OT): OT refers to several related methods that include speckle tracking, coherence tracking and amplitude tracking or feature tracking. Feature tracking uses cross-correlation of image patches to find the displacement of surface features such as crevasses or rifts and edges, that move with the same speed as the ice and are identifiable on two co-registered amplitude images, to derive ice flow velocity. In coherence tracking, the offset which maximises the interferometric coherence within a certain window size is determined and used to derive the ice velocity. Speckle tracking uses the cross-correlation function of radar speckle patterns, rather than visible features, to derive ice flow velocity.

Scope of the document

This document is the Algorithm Theoretical Basis Document (ATBD) for the Ice Velocity (IV) product, v1.4 produced as part of the Copernicus Ice Sheets and Ice Shelves service. It describes the satellite-mounted instruments, auxiliary datasets, auxiliary models and basic algorithms used to create the data products.

Executive summary

The v1.4 Ice Velocity product is a gridded product representing the mean annual ice surface velocity (IV) of the Greenland Ice Sheet. The main input data consist of Sentinel-1A and Sentinel-1B synthetic aperture radar (SAR) data. In Chapter 1 we introduce the satellite source for IV retrieval (Sentinel-1) including the mission's acquisition strategy, technical details of the SAR sensor and characteristics of the Interferometric Wide (IW) swath mode. Chapter 2 describes the input and auxiliary data needed for the IV retrieval algorithm and quality assessment. The primary input data are Level-1 Single Look Complex (SLC) images acquired in IW swath mode. Auxiliary data, used in the processing chain for IV production and quality assessment, include precise orbit files, a digital elevation model (DEM), ocean/ice land masks and validation data.

The retrieval algorithm and processing line for IV and the uncertainty characterization is specified and explained in Chapter 3. The algorithm is based on offset tracking (OT), utilizing long stripes of Sentinel-1 data. Offset tracking refers to several related methods that include amplitude tracking or feature tracking, coherence tracking and speckle tracking. OT requires less operator interaction than interferometric techniques (InSAR), and can be more easily and reliably automated. Another advantage of the OT technique over InSAR is the provision of two components of the velocity vector (along track and line-of-sight, LOS) from data of a single track. The generation of a regional ice velocity map requires the combination of OT results, from several tracks, each providing slant range and azimuth displacements at local incidence and heading angle. For various reasons the algorithm sometimes fails to find matching features or erroneous matches are returned leading to gaps and/or errors in the velocity fields. These are to some extent dealt with using a filtering and gap filling approach. Since the calving front location (CFL) is in fact a highly dynamical environment, dynamic ice/ocean masking for outlet glaciers based on annually/periodically updated CFL's is added as a technical enhancement for future Climate Data Records (CDR's).

Details and examples of the output Greenland Ice Sheet IV product are provided in Chapter 4. We conclude with the reference list in Chapter 5.

1. Instruments

The primary satellite sources for generation of the ice velocity (IV) fields are Copernicus Sentinel-1A and Sentinel-1B. The Sentinel-1 (S1) satellites carry a C-band synthetic aperture radar (SAR) instrument providing high resolution SAR images in different acquisition modes. The Interferometric Wide-Swath mode (IW) is the main acquisition mode over land areas, including Greenland. The IW mode combines a large swath width (250 km) with a moderate ground resolution (5 m x 20 m in range and azimuth, respectively).

Sentinel-1A provides continuous coverage for the entire Greenland Ice Sheet margin since October 2014 with a 12-day repeat pass period. Sentinel-1B was launched in April 2016 and, in combination with Sentinel-1A, reduced the repeat pass period from 12 to 6 days. Besides the continuous coverage of the margin, the ice sheet is covered completely during dedicated annual (winter) campaigns, with 4 to 6 repeat observations per track. In 2019, further expansion of the continuous coverage in Greenland commenced including also the interior ice sheet. This provides an opportunity to produce Greenland wide velocity maps at sub-annual, and even monthly, intervals. The S1 dedicated acquisition plan for Greenland makes the constellation the primary source for monitoring of IV, providing year-round, day/night, all-weather observing capability. The main instrument parameters and IW mode are described in Table 1 and Table 2.

Note: Since 23rd December 2021, there has been a technical problem with Sentinel-1B related to the power supply of the SAR system and unfortunately no data has been acquired since. Recovery attempts have been unsuccessful and European Space Agency (ESA) has announced the end of the mission. Sentinel-1A remains fully operational and Sentinel-1C, which will replace Sentinel-1B, is scheduled for launch in April 2023.

Table 1: Sentinel-1 sensor specifications

Table 2: Characteristics of Sentinel-1 IW mode

2. Input and auxiliary data

2.1. SAR data

Primary input data for IV retrieval are repeat pass Sentinel-1A and Sentinel-1B Level-1 Single Look Complex (SLC) images acquired in the Interferometric Wide (IW) swath mode. These products consist of focused SAR data geo-referenced using orbit and altitude data from the satellite and provided in zero-Doppler slant-range geometry. IW implements the Terrain Observation by Progressive Scans (TOPS) mode, consisting of 3 or 5 sub-swaths, each consisting of a series of bursts. Each burst has been processed as a separate SLC image. Sentinel-1 data is provided free of charge through the Copernicus Open Access Hub¹(SciHub) and various mirror sites.

2.2. Auxiliary data

Auxiliary data, used in the processing chain for IV production and quality assessment, include precise orbit files, a digital elevation model, ocean/ice land masks and validation data. Details are provided in the following:

2.2.1. Orbits

Orbit files (Precise Orbit Ephemerides, POE) are required for IV processing to enable accurate co-registration. Precise orbit files for Sentinel-1 are provided in ASCII XML format by the Copernicus Precise Orbit Determination (CPOD) service available from the Copernicus Sentinels POD Data Hub² approximately 3 weeks after acquisition. POE files cover periods of approximately 26 hours (one complete day in Global Positioning System (GPS) time, plus overlap of two hours between consecutive files) and contain orbit state vectors at 10-second intervals.

2.2.2. DEM

A contemporary high-resolution digital elevation model (DEM) is required in the processing chain for co-registration, calculation of the vertical component of velocity (v_z) and final geocoding of the product. The DEM used for the annual Greenland Ice Sheet velocity map is compiled from TanDEM-X 90m DEM tiles provided by the German Aerospace Center DLR³. This is a product variant of the original 12m (0.4 arcsec) DEM product with reduced pixel spacing (Wessel et al., 2018a). The extent and grid spacing of the DEM is equal to the IV product. Further information can be found in the TanDEM-X Ground Segment DEM Products Specification Document (Wessel, 2018b).

2.2.3. Land/Ice/Ocean Mask

This is an ESRI-shapefile of the ocean, land and ice boundaries, based on recent optical imagery and used for masking the retrievals in the ocean (usually coming from tracked sea ice). The ice sheet and glacier boundaries are based on the Randolph Glacier Inventory (RGI 6.0, RGI Consortium, 2017) with updated glacier fronts. The inventory has been compiled from more than 70 Landsat scenes (mostly acquired between 1999 and 2002) using semi-automated glacier mapping techniques (Rastner et al., 2012). For C3S2 dynamic ice/ocean masking for outlet glaciers based on annually updated calving front locations (CFLs) is added as technical development (See Chapter 3.4).

2.2.4. Validation data

The quality assessment for ice velocity includes detailed validation with contemporaneous in-situ GPS data at various sites across the ice sheet (see: RD.1 and RD.4). The GPS instruments are attached to Automatic Weather Stations (AWS) operated by the Geological Survey of Denmark and Greenland (GEUS) in collaboration with the National Space Institute at the Technical University of Denmark (DTU Space) and ASIAQ Greenland Survey as part of the Danish Programme for Monitoring of the Greenland Ice Sheet (PROMICE; Fausto and Van As, 2019). The product is also evaluated against publicly available ice velocity products covering the same area and time span. These IV maps were produced as part of the NASA 'Making Earth System Data Records for Use in Research Environments' (MEaSUREs) program (Joughin et al., 2021; Joughin, 2021).

3. Algorithms

3.1. Coherent and Incoherent Offset Tracking using Sentinel-1

Offset tracking (OT) is the principal method for generating velocity fields on glacier surfaces from Sentinel-1 operationally. OT refers to several related methods that include amplitude tracking or feature tracking, coherence tracking and speckle tracking. Feature tracking uses cross-correlation of image patches to find the displacement of surface features such as crevasses or rifts and edges, that move with the same speed as the ice and are identifiable on two co-registered amplitude images, and subsequently ice flow velocity (Figure 3.1). In coherence tracking, the offset which maximises the interferometric coherence within a certain window size is determined and used to derive the ice velocity. Speckle tracking uses the cross-correlation function of radar speckle patterns, rather than visible features, to derive ice flow velocity.

Figure 3.1: Concept of Feature tracking. Two correlation plots derived from comparing a reference image template with a larger search window extracted from two different SAR scenes. The peak of the surface corresponds to the reported match. The shape of the correlation function is an important indicator of the measurement accuracy. Sharp pronounced peaks (a) have a higher accuracy then broader peaks (b). (Figure adapted from Wuite, 2006)

OT requires less operator interaction than interferometric techniques (such as Interferometric SAR (InSAR)), and can be more easily and reliably automated. Another advantage of the OT technique over InSAR is the provision of two components of the velocity vector (along track and line-of-sight, LOS) from data of a single track. Also, the Greenland Ice Sheet is subject to highly variable meteorological conditions, resulting in rapid temporal decorrelation of the radar signal phase impairing the comprehensive application of InSAR. OT is less sensitive or insensitive to temporal decorrelation of the radar phase signal, albeit at a lower accuracy of velocity. In areas with distinct and stable surface features, as common on mountain glaciers and outlet glaciers of ice sheets, coherence of the repeat-pass SAR data is not required because the cross-correlation of image-templates is based on the amplitude signal. In the level interior sections of ice sheets, without distinct amplitude features, temporal stability of the target phase (coherence) is required for speckle tracking. Loss of stability of the target phase (temporal decorrelation) is less critical than for the InSAR method because a continuous path is not required for deriving velocities of an isolated patch in an image, as is the case with InSAR.

The ice velocity algorithm is developed as part of the ESA Climate Change Initiative (CCI) and Austrian Space Applications (ASAP) Programmes. The IV system's core module performs coherent and incoherent OT, utilizing long stripes of Sentinel-1 data acquired in interferometric wide swath (IW) mode. For Sentinel-1 Terrain Observation by Progressive Scans (TOPS) mode acquisitions, OT algorithms need to support burst handling and stitching. Also, the azimuth-varying phase introduced by the TOPS beam steering must be accounted for, both for complex offset tracking and for incoherent offset-tracking, where a factor of two oversampling of the SLCs is required prior to detection.

OT accuracies depend on various factors, including on the level of coherence, on the correlation window size (Bamler and Eineder, 2005) as well as on the type of features that are being tracked. In addition, ionospheric scintillations can have a large impact on offset tracking by causing azimuth shifts that are introduced by the fluctuating electron density along the sensor path. Misregistration is not a problem, but the large azimuth shifts are interpreted as ice motion and observed offsets can exceed the azimuth pixel size. DEM errors have less impact on the algorithm performance compared to InSAR, since there is no phase unwrapping. However, since offset-tracking can be applied with longer temporal baselines, and even though larger spatial baselines can be used for offset-tracking, in the typical scenario the sensitivity to DEM errors is much smaller.

For further details the reader is referred to the Algorithm Theoretical Basis Documents of the Greenland and Antarctic Ice Sheet CCI projects (Wuite et al. 2021, Wuite et al. 2020) and Nagler et al. (2015).

3.2. Generation of Regional Ice Velocity Maps from Multiple Tracks using OT

The generation of a regional ice velocity map requires the combination of offset tracking results, from several tracks, each providing slant range and azimuth displacements at local incidence and heading angle. Figure 3.2 illustrates the high-level processing line for IV production. The system includes 3 main modules: 1) the IV Module; 2) the Merge Module; and 3) the Validation Module. These modules are detailed in the System Quality Assurance Document, but briefly explained here.

Figure 3.2: High-level flow chart of the IV processing system. Green – input data, Blue – processing modules, Red - product and intermediate products, Yellow – product data base.

3.2.1. IV Module

Within the IV module SAR data and orbit data are imported into the system and velocity maps are generated for pairs of repeat pass data of the same track. The long tracks are geocoded into the common map projection of the output grid, providing observations of slant range displacement and azimuth for each track separately. The local incidence and heading angle are calculated using the annotated image and orbit parameters and a DEM as input. The IV Module includes an outlier removal/filter and gap filling interpolation scheme (See Chapter 3.3).

3.2.2. Merge Module

The Merge module combines all IV products from all tracks and image pairs over a specified time span (i.e. 1 year), applying a weighted average approach. The combination of displacement observations of multiple tracks applies a weighted least squares model to fit the multiple observed displacements according to:

where x is the horizontal velocity vector (in local Easting/Northing coordinates), y is a vector with the observed velocities (in slant range/azimuth), A is the matrix projecting the horizontal velocity to slant range/azimuth geometry, and ε is a noise vector. The matrices are defined as:

where f() is the master, t() is the template image (sub-image you are searching for in matrix form) at the pixel position (u,v). The peak in the correlation matrix corresponds to the best matching shift of the slave window with respect to the master window.
The weighted least squares estimate of the horizontal velocity vector is then:

An advantage of the approach above is that if crossing tracks (ascending/descending) are used, the displacement measurement in the direction with higher resolution dominates in the velocity retrieval.

3.2.3. Validation Module

Lastly, the Validation Module facilitates quality assessment of the IV products, by automating standard validation tests, consisting of internal consistency checks and intercomparisons with independent data sets (e.g. in-situ, GPS, or other published velocity maps). The output is statistical information (e.g. scatter plots, bias, root-mean-square error) on the intercomparisons compiled in the Product Quality and Assessment Report (PQAR, RD.1).

3.3. Filtering and interpolation

For various reasons the OT algorithm sometimes fails to find matching features or erroneous matches are returned leading to gaps and/or errors in the velocity fields. Most frequently this is caused by a lack of traceable surface features or low coherence between image pairs due to for example summer surface melt. In order to smooth and remove outliers a 3x3 simple median filter is applied, after which a 9x9 first order plane fit filter is applied to fill small data gaps. Further filtering/gap filling is left to the user if required.

3.4. Dynamic ice/ocean masking

The existing Greenland IV products applied a fixed ice-ocean mask derived from Landsat images acquired in summer 2016. Since the calving front location (CFL) is in fact a highly dynamic environment, dynamic ice/ocean masking for outlet glaciers based on annually/periodically updated CFL's is added as a technical enhancement for future CDR's. For CDR v4.0 we updated the existing ice/ocean mask which is applied to the ice velocity map as a final step. For this we generated a 10 m geocoded and cloudless Copernicus Sentinel-2 mosaic (band 8 Near Infrared (NIR)) covering the ice sheet margin. All Sentinel-2 images were acquired in the period August-September 2022. This is at the end of the summer when most of the fjords⁴ are ice free and hence ice fronts are easily distinguishable from open water. For CDR v4.0, the main focus was on updating the calving fronts for 28 large glaciers. These glaciers have been identified in the ESA Greenland Ice Sheet CCI project as key outlet glaciers and are all major active glacier systems or ice streams (Figure 3.3). The delineation of the calving front is performed manually in QGIS (open-source GIS software) by an experienced operator. Figure 3.4 shows an example of the updated ice mask in the vicinity of Steenstrup Gletsjer (NW-Greenland).

Figure 3.3: Cloudless 10 m Copernicus Sentinel-2 mosaic of the Greenland Ice Sheet margin. The images are acquired in the period August-September 2022. Shown are the selection of glaciers for updating the calving front/ice edge. All of these glaciers are marine terminating, except Isunnguata Sermia.

Figure 3.4 Left: Ice velocity map with old land/ice/ocean mask in the vicinity of Steenstrup Gletsjer. Right: same with updated masking. Inset: location in Greenland. Background: Sentinel-2 image mosaic from August-Sept 2022.

3.5. Uncertainty characterisation

For incoherent offset tracking (the general case), the uncertainty is characterised for each pixel by calculating a local offset-map standard deviation in a 5x5 neighbourhood. A plane fit to the offset map in the 5x5 neighbourhood is subtracted prior to calculating the standard deviation, so that an actual velocity gradient is not interpreted as a noise signal. The standard deviation estimate is corrected for any averaging carried out, as well as correlation between neighbouring samples (i.e. if the radar data are oversampled). The Greenland Ice Sheet velocity map is accompanied by its associated error standard deviation. The latter is also a map, in the same geometry as the associated measurement, providing a measure of uncertainty on a per-pixel basis.

4. Output data

The ice velocity product covers the Greenland Ice Sheet and is provided as a NetCDF file with 7 separate layers/maps representing different variables: the 3 velocity components: v_x, v_y, v_z as well as v_v (the magnitude of the horizontal components), and gridded maps showing the valid pixel count and the uncertainty/standard deviation (easting & northing) (Table 3, Figure 4.1). The ice velocity map is annually averaged and provided at 250m grid spacing in North Polar Stereographic projection (EPSG: 3413). A new NetCDF file of annually averaged IV is produced each year. The horizontal velocity is provided in true meters per day, towards easting (v_x) and northing (v_y) direction of the grid, and the vertical displacement (v_z), is derived from a digital elevation model. More details are given in Section 1.3 of the C3S Product User Guide for Ice Velocity [RD.2].

Table 3: The variables provided in the IV output product for the Greenland Ice Sheet. 

Figure 4.1: Example IV product covering the Greenland Ice Sheet (2019-2020), depicted are the easting component (v_x), the northing component (v_y), the magnitude of velocity (v_v), the standard deviation of the easting velocity component (v_x STD) and the valid measurement count (Count).

References

Bamler, R. and Eineder, M. (2005): Accuracy of differential shift estimation by correlation and split-bandwidth interferometry for wideband and delta-k SAR systems, IEEE Geosc. Rem Sens. Lett. 2, 151-155.

De Zan, F., and Guarnieri, A. M. (2006): TOPSAR: Terrain Observation by Progressive Scans. Geoscience and Remote Sensing, IEEE Transactions on, 44(9), 2352-2360. doi:10.1109/TGRS.2006.873853

Fausto, R.S. and van As, D., (2019): Programme for monitoring of the Greenland ice sheet (PROMICE): Automatic weather station data. Version: v03, Dataset published via Geological Survey of Denmark and Greenland. DOI: https://doi.org/10.22008/promice/data/aws [Date Accessed: 15th December 2022].

Joughin, I. (2021): MEaSUREs Greenland Annual Ice Sheet Velocity Mosaics from SAR and Landsat, Version 3. Boulder, Colorado USA. NASA National Snow and Ice Data Center Distributed Active Archive Center. doi: https://doi.org/10.5067/C2GFA20CXUI4 [Date Accessed: 15th December 2022].

Joughin, I., Howat, I., Smith, B. and Scambos, T. (2021): MEaSUREs Greenland Ice Velocity: Selected Glacier Site Velocity Maps from InSAR, Version 4. Boulder, Colorado USA. NASA National Snow and Ice Data Center Distributed Active Archive Center. doi: https://doi.org/10.5067/GQZQY2M5507Z [Date Accessed: 15th December 2022].

Nagler, T, Rott, H., Hetzenecker, M., Wuite, J. and Potin, P., (2015): The Sentinel-1 Mission: New Opportunities for Ice Sheet Observations. Remote Sens., 7, 9371-9389.

Rastner, P., Bolch, T., Mölg, N., Machguth, H., Le Bris, R., and Paul, F. (2012): The first complete inventory of the local glaciers and ice caps on Greenland, The Cryosphere, 6, 1483-1495, https://doi.org/10.5194/tc-6-1483-2012 [Date Accessed: 15th December 2022], 2012.

RGI Consortium (2017): Randolph Glacier Inventory – A Dataset of Global Glacier Outlines: Version 6.0: Technical Report, Global Land Ice Measurements from Space, Colorado, USA. Digital Media. DOI: https://doi.org/10.7265/N5-RGI-60 [Date Accessed: 15th December 2022]

Wessel, B., Huber, M., Wohlfart, C., Marschalk, U., Kosmann, D., & Roth, A. (2018a): Accuracy assessment of the global TanDEM-X Digital Elevation Model with GPS data. ISPRS Journal of Photogrammetry and Remote Sensing, 139, 171–182. doi:10.1016/j.isprsjprs.2018.02.01

Wessel, B. (2018b): TanDEM-X Ground Segment – DEM Products Specification Document”, EOC, DLR, Oberpfaffenhofen, Germany, Public Document TD-GS-PS-0021, Issue 3.2. Available: https://geoservice.dlr.de/web/dataguide/tdm90/pdfs/TD-GS-PS-0021_DEM-Product-Specification.pdf [Date Accessed: 15th December 2022]

Wuite, J. et al., (2020): Algorithm Theoretical Baseline Document (ATBD) for the Antarctic Ice Sheet CCI+ Phase 1, version 1.0. Available from https://climate.esa.int/en/projects/ice-sheets-antarctic/ [Date Accessed: 15th December 2022]

Wuite, J. et al., (2021): Algorithm Theoretical Baseline Document (ATBD) for the Greenland Ice Sheet CCI+ Phase 1, version 1.4. Available from https://climate.esa.int/en/projects/ice-sheets-greenland/ [Date Accessed: 15th December 2022]

Related articles

• Page:
  ECV - Glaciers, Ice Sheets, Lakes, Soil Moisture
• Page:
  IV Version 1.2; SEC Version 2.0; GMB Version 2.0: Product User Guide Specification (PUGS)
• Page:
  IV Version 1.2; SEC Version 2.0; GMB Version 2.0: System Quality Assurance Document (SQAD)
• Page:
  IV Version 1.2; SEC Version 2.0; GMB Version 2.0: Product Quality Assurance Document (PQAD)
• Page:
  IV Version 1.2; SEC Version 2.0; GMB Version 2.0: Product Quality Assessment Report (PQAR)