Sea ice drift from Sentinel-1 SAR imagery using open source feature tracking

A computational ::::::::::::: computationally : efficient, open source feature tracking algorithm, called ORB, is adopted and tuned for sea ice drift retrieval from Sentinel-1 SAR images. The best suitable setting and parameter values have been found using four representative Sentinel-1 image pairs . A new quality measure for feature tracking algorithms is introduced utilising the distribution of the resulting vector field ::::::::::: representative :: of :::: sea ::: ice ::::::::: conditions :::::::: between ::::::::: Greenland :::: and ::::::::: Severnaya 5 :::::: Zemlya :::::: during :::::::::::: winter/spring. The performance of the algorithm is compared with :: to : two other feature tracking algorithms (SIFT and SURF). Applied on a test image pair acquired over Fram Strait, the tuned ORB algorithm produces the highest number of vectors (6920, SIFT: 1585 and SURF: 518) while being computational ::::::::::::: computationally most efficient (66 s, SIFT: 182 s and SURF: 99 s using a 2,7 :: 2.7 GHz processor with 8 GB memory). For validation purpose, 350 ::: 314 manually 10 drawn vectors have been compared with the closest calculated vectors : , : and the resulting root mean square distance is 609.9 :::: error :: of ::: ice :::: drift :: is :::: 563 m(equivalent to 7.5pixel). All test image pairs show a significant :::::::::: significantly : better performance of the HV channel ::: due :: to :::::: higher ::::::::::::: informativeness. On average, around 4 ::: four : times more vectors have been found using HV polarisation. All software requirements necessary for applying the presented feature tracking algorithm are open source to ensure 15 a free and easy implementation.


Introduction
Sea ice motion is an essential variable to observe from earth observation data, because it strongly influences the distribution of sea ice on different spatial and temporal scales.Ice drift causes advection of ice from one region to another and export of ice from the Arctic Ocean to the sub-Arctic seas (IPCC, 2013).Furthermore, ice drift generates convergence and divergence zones that cause formation of ridges and leads.Presently, sea ice drift data do not provide sufficient resolution to estimate convergence and divergence fields on a spatial scaling of a few kilometres.
The main region of interest for this work are the ice covered oceans between Greenland and Severnaya Zemlya, i.e.Greenland Sea, Barents Sea, Kara Sea and the ad-Introduction

Conclusions References
Tables Figures

Back Close
Full jacent part of the Arctic Ocean.These oceans are characterised by a strong seasonal cycle of sea ice cover, a large variation of different ice classes (Multi Year Ice, First Year Ice, Marginal Ice Zone etc.) and a wide range of drift speeds (e.g.strong ice drift in Fram Strait).
With systematic acquisition of space-borne Synthetic Aperture Radar (SAR) data over sea ice areas, Kwok (2010) and Kwok and Sulsky (2010) have demonstrated that high resolution ice drift fields can be derived from SAR data.SAR is an active microwave radar which can acquire data independent of solar illumination and weather condition.Sea ice motion fields of the Arctic Ocean with a grid spacing of 5 km have been produced on a weekly basis between 1997-2012 using Radarsat and ENVISAT (Environmental Satellite) SAR data and the geophysical processor system introduced by Kwok (1998).Thomas et al. (2008) have used pattern recognition to calculate sea ice drift between successive ERS-1 (European Remote-sensing Satellite) SAR images with a resolution of 400 m.This work has been continued by Hollands and Dierking (2011) using ASAR data from ENVISAT.
With the successful launch of Sentinel-1A in April 2014 and the planned launch of Sentinel-1B in early 2016, high resolution SAR data will be delivered for the first time with open and free access for all users and a never before reached repeat frequency of less than 1 day in the Arctic (ESA, 2012).This introduces a new area in SAR earth observation, but no sea ice drift algorithm using Sentinel-1 data has been published so far.The objective of this paper is to identify and develop the most efficient open source method for high resolution sea ice drift retrieval from Sentinel-1 data.
Our goal is to exploit recent improvements and developments in computer vision by adopting a state of the art feature tracking algorithm to derive sea ice drift.The advantage of feature tracking to algorithms based on pattern recognition is that each drift vector is independent of the surrounding motion, which leads to better performance e.g.along shear zones.For application on large data sets and for operational use, a computational efficient algorithm, called ORB (Oriented FAST and Rotated BRIEF) Introduction

Conclusions References
Tables Figures

Back Close
Full ( Rublee et al., 2011), has been considered, tuned and compared with other available feature tracking algorithms.The software requirements necessary for deriving ice drift fields from Sentinel-1 data (python with openCV and the python toolbox Nansat) are all open source to ensure a free, user friendly and easy implementation.

Data
The Sentinel-1 mission, an initiative of the European Space Agency (ESA), is composed of a constellation of two identical satellites sharing the same near-polar, sunsynchronous orbit: Sentinel-1A, launched in April 2014, and Sentinel-1B, planned launch in early 2016.Sentinel-1 carries a single C-band Synthetic Aperture Radar (SAR) instrument operating at a centre frequency of 5.405 GHz and supporting dual polarisation (HH+HV, VV+VH).With both satellites operating, the constellation will have a repeat frequency of less than 1 day in the Arctic.Radar data are delivered to Copernicus services within an hour of acquisition with open and free access for all users (ESA, 2012).
The Sentinel-1 product used in this paper is called "Extra Wide Swath Mode Ground Range Detected with Medium Resolution".These images cover an area of 400 × 400 km with a resolution of 40 m and provide both HH and HV polarisation.Four image pairs (Table 1) representative for our region of interest have been chosen, covering a range of different sea ice conditions (pack ice, fast ice, leads, ridges, marginal ice zone, ice edge etc.) and time spans between the acquisitions.
Sentinel-1 datasets were opened and processed with the open source software Nansat (https://github.com/nansencenter/nansat)(Korosov et al., 2015).Nansat is a scientist friendly Python toolbox for processing 2-D satellite earth observation data.It is based on the Geospatial Data Abstraction Library (GDAL) and provides easy access to geospatial data, simple and generic interface to common operations including reading, geographic transformation and export.Nansat proves to be efficient both for Figures

Back Close
Full development and testing of scientific algorithms and for fast operational processing.To extend the functionality of GDAL, Nansat reads metadata from XML files accompanying Sentinel-1 data and supplements the GDAL data model with georeference information stored as ground control points (GCPs).Originally GCPs are pairs of latitude/longitude and corresponding pixel/line coordinates.In order to increase the accuracy of the geographic transformation, the projection of GCPs is changed from cylindrical to stereographic centered at the center of the scene.The reprojected GCPs are then utilised by GDAL to calculate geographic coordinates of any pixel in the raster using spline interpolation.Reprojection of GCPs does not require much additional computational effort, but improves the result significantly, particularly in high latitudes.

Method
Our algorithm for sea ice drift detection includes three major steps: (a) resampling of raw data to lower resolution, (b) detection and matching of features and (c) a new introduced quality measure.
(a) To decrease the influence of speckle noise and increase the computational efficiency, the resolution is reduced before applying the ice drift algorithm.Various resolution resampling algorithms have been tested: Gaussian, nearest neighbour, bilinear, cubic, cubic spline and Lanczos.Best performance and computational efficiency was achieved by using simple averaging from 40 m to 80 m resolution.
(b) For detection and tracking of features on large data sets and for operational use, a computational efficient algorithm, called ORB (Oriented FAST and Rotated BRIEF), has been utilised (Rublee et al., 2011).In our numerical experiments we tuned the parameters of ORB for optimal SAR sea ice drift application.The best suitable parameter set (including spatial resolution of SAR image, patch size of FAST descriptor, number of pyramid levels, scale factor etc.) has been evaluated.(c) A new quality measure using the amount and deviation of vectors in a grid cell is introduced for feature tracking algorithms.The introduced ORB setup is tested against other available OpenCV feature tracking algorithms for comparison and manually drawn vectors for validation.

ORB algorithm
ORB (Oriented FAST and Rotated BRIEF) is a feature tracking algorithm introduced by Rublee et al. ( 2011) as a computationally-efficient replacement to SIFT (Lowe, 2004) with similar matching performance and less affected by image noise.ORB builds on the FAST keypoint detector (Rosten and Drummond, 2006) and the binary BRIEF descriptor (Calonder et al., 2010) with many modifications to enhance the performance.
It uses FAST to find multiscale-features on several pyramid levels and applies Harris corner measure (Harris and Stephens, 1988) to find the best ones among them.To achieve rotation invariance, the orientation of the feature is calculated by using the intensity weighted centroid of a patch with the located keypoint at the centre.The ORB descriptor performs as well as SIFT and better than SURF (Bay et al., 2006), while being almost two orders of magnitude faster (Rublee et al., 2011).To match features, we use a Brute-Force matcher and Hamming-distance.An additional benefit of ORB is that it is free from licensing restrictions, unlike SIFT and SURF.Before the feature tracking algorithm can be applied on a satellite image, the SAR backscatter values σ 0 have to be transformed into the intensity i range [0,255] used in openCV.This transformation is done by using Eq. ( 1) and setting all intensity values below and above the range to 0 and 255.
Converting the linear backscatter values into decibel units before the transformation has been tested, but decreased the algorithm performance for both HH and HV chan- After the transformation into intensity values, keypoints are detected on both SAR scenes using the FAST-9 keypoint detector (Rosten and Drummond, 2006).FAST-9 compares the intensity of a centre pixel to the intensities of the surrounding circle with perimeter of 16 pixels.If 9 contiguous pixel in the circle have an intensity difference greater (and with the same sign) than a certain threshold, the centre pixel is recoginsed as keypoint.The threshold is set low enough to get more than the predefined amount of retained keypoints.
Each keypoint is addressed a score R using the intensity variation around the keypoint (Harris corner measure Harris and Stephens, 1988).A high intensity variation in both dimensions returns a high R value.The predefined amount of keypoints with the highest R values are utilised.
To detect features of different scales, the keypoint search is performed on several pyramid levels.The number of pyramid levels in combination with the scale factor defines the range and increment of the keypoint detection scaling.
ORB adds an orientation θ to each keypoint, derived from connecting the keypoint and the intensity weighted centroid of a patch p with the keypoint at the centre (Rosin, 1999).
To describe keypoints, ORB applies a modified version of the binary keypoint descriptor BRIEF (Calonder et al., 2010).A binary test τ on a patch p is defined by: τ(p; x, y) := 1 if p(x) < p(y) where p(x) is the intensity value at point x.A feature f can be described by a vector of n binary tests:

Conclusions References
Tables Figures

Back Close
Full A set of n binary tests with sampling pair location (x i , y i ) can be written in a 2 x n matrix: To be invariant to in-plane rotation, Rublee et al. ( 2011) steers S according to the orientation θ using the corresponding rotation matrix R θ : A good set of sampling pairs needs to be uncorrelated, so that each pair adds new information to the descriptor, and have high variance, to make features more discriminative.Rublee et al. ( 2011) applied a greedy search on a training dataset to obtain a set of 256 relatively uncorrelated tests with high variance.
After detection and description, each feature of the first image is compared to all features in the second image (Brute Force matching) using the number of positions in which the feature vectors have a different value (Hamming distance).Our setting returns the best two matches and applies the ratio test from Lowe (2004) to decide whether the best match is accepted or rejected.

ORB setting and parameter tuning
Achieving the best possible performance of ORB for sea ice drift from Sentinel-1 images, requires a good setting and tuning of the parameters shown in Table 2.As a compromise between performance and computational efficiency, the resolution of the Sentinel-1 image is reduced from 40 to 80 m (resize factor = 0.5) and the amount of maximum retained keypoints is set to 100 000.
It is not recommended to re-project one image onto the projection of the second image before applying the ORB algorithm, since this is computational very expensive.Instead, geographic coordinates of the matched start and end point shall be calculated independently using the georeference information from GCPs of the first and second image.To reject less reliable matches, we use the ratio test explained in Lowe ( 2004).Manual interpretation of ice drift results (using the training data from Table 1), suggest a good compromise between amount of vectors and correct results with a ratio test threshold of 0.75, meaning that the Hamming-distance of the best match has to be less than 0.75 × Hamming-distance of the second best match.Tested on the image pairs from Table 1, the ratio test showed clear better performance and is computational less expensive than the alternative cross-check, where keypoints are matched in both directions (first image to second image and vice versa) and rejected, if the matches are not correlated.
Unreasonable high velocities above 10 km + 1 km h −1 are removed in a postprocessing step of the drift field.
To tune the remaining parameters in Table 2, we assume that the amount of wrong matches relative to the amount of correct matches does not increase with increasing number of matches.Using this assumption, it can be concluded that more matches equals better algorithm performance.ORB is computational very efficient, making it possible to test the remaining parameters over a wide range with high resolution using both HH and HV polarisation.This has been done to find the best suitable values for patch size, HH and HV brightness boundaries, pyramid levels and scale factor.

Quality measure
The resulting sea ice drift vector field is not evenly distributed, but according to the recognition performance of the respective area.Regions with few vectors represent low reliability, whereas regions with many vectors suggest high reliability.By using a grid and calculating the amount (N) and root mean square distance (RMSD) of all vectors appearing in one grid cell, the distribution of the vector field can be used as a quality measure.Introduction

Conclusions References
Tables Figures

Back Close
Full

Results
Before testing the individual parameters, the remaining parameters were set to the following values: pyramid levels = 8, scale factor = 1.2, HH limits = [0,0.12],HV limits = [0,0.012],ratio test = 0.8.An additional low speed filter with 2.5 km is applied during the testing to reject matches over land.This filter does not influence the number of correct matches, since the sea ice velocities in all considered test images are above 2.5 km.Once the best suitable value for a certain parameter is found, it is applied for further testing.

Patch size
Figure 1 shows that patch size values between 10 and 60 pixel can vary the resulting amount of vectors by an order of magnitude.To resolve velocity gradients with high resolution, the patch size shall be as small as possible.Taking that into account and the performance represented by the amount of matches, the best suitable patch size was chosen to be 34 pixel (2.72 km).This yields around 1000 and 6000 vectors for HH and HV respectively for the training dataset from Table 1.

Brightness boundaries
The performance of the algorithm (represented by the amount of matches) for different backscatter limits σ max and σ min for HH and HV polarisation is shown in Fig. 2. Within the chosen backscatter range, the amount of vectors can vary by an order of magnitude.As a compromise between the different results of the four image pairs, the maximum backscatter σ max is suggested to be set to 0.08 and 0.013 for HH and HV.The chosen minimum σ min is 0 for both HH and HV, since the number of matches is decreasing towards higher values for most training images.Applying this setting on the training dataset yields on average around 1500 and 6000 vectors for HH and HV.

TCD Figures Back Close
Full

Pyramid levels and scale factor
Figure 3 displays the number of matches using 1-16 pyramid levels and the scale factors 1.1, 1.2, 1.3 and 1.4.As a compromise between performance, i.e. number of matches, and computational efficiency (linked to the number of pyramid levels), a scale factor of 1.2 with 7 pyramid levels was chosen.

Recommended parameter set
Table 2 shows the recommended parameter set for ORB Sentinel-1 sea ice drift application.Using these parameters yielded the best compromise between performance and computational efficiency for the 4 representative image pairs from Table 1.

HH and HV comparison
Figure 1, 2 and 3 display the HH and HV results with solid and dashed lines respectively.All image pairs show a significant better performance of the HV channel.On average, around 4 times more vectors have been found using HV.Even the image pair "Franz Josef Land" (Table 1), which has the best HH and the worst HV performance, shows more than two times more vectors using HV channel.However, due to the different appearance of sea ice in the HH and HV image, the spatial distribution of the resulting drift vectors is also slightly different, supporting the usefulness of a combination of both results.

Comparison with SIFT and SURF
To compare the introduced ORB setup with other available OpenCV feature tracking algorithms, SIFT (Lowe, 2004) and SURF (Bay et al., 2006) are considered.The performance of the three algorithms is tested on the image pair "Fram Strait" (Table 1) using the introduced quality indexes number of vectors (N) and root mean square distance (RMSD) on a grid with cell size 1 Full

Computational efficiency
The OpenCV feature tracking algorithms ORB, SIFT and SURF in combination with the python-toolbox 'Nansat' are sufficient computational efficient to compute sea ice drift fields from datasets with large temporal and spatial extent.The processing time indications shown in Table 3 are based on testing the algorithms on a MacBook Pro from early 2013 with a 2.7 GHz Intel Core i7 processor and 8 GB 1600 MHz DDR3 memory.Applying the introduced ORB algorithm needs 36 and 67 % of the processing time necessary to compute drift fields with SIFT and SURF respectively.

Validation
In order to validate the calculated ORB drift field shown in Fig. 4 a, 350 features have been identified by a sea ice expert in both images and manually connected using Ar-cGIS.Figure 5 shows the manually drawn vectors (green) and the respective nearest neighbour vectors from ORB (red).The resulting RMSD between manual and calculated vectors is 609.9 m (equivalent to 7.5 pixel).The RMSD represents a combination of the manually produced error and the displacement variation between the manual and calculated vector.Using the tuned ORB algorithm, a total of 6920 vectors have been calculated within 3 minutes, whereas identifying 350 sea ice drift vectors manually takes several hours.The presented ORB algorithm outperforms SIFT and SURF not only in processing time, but also in quantity and quality, measured by the two introduced indexes N and RMSD.This shall proof that ORB is the best option for feature tracking of sea ice on Sentinel-1 SAR imagery.In addition, the ORB results compare very well with manually drawn drift vectors, proofing good reliability of the algorithm.Current algorithms for calculating drift vectors from consecutive image pairs are based either on feature tracking or pattern recognition.Feature tracking provides vectors, which are independent from each other, whereas pattern recognition includes the surrounding drift information.Therefore pattern recognition is more prone to errors in areas with high velocity gradients.The resulting drift fields from feature tracking are generally not evenly distributed, but according to the feature recognition performance of the respective area.This can be utilised for quality measure (as shown in this paper), but is a disadvantage when it comes to estimation of divergence, shear and deformation.Pattern recognition algorithms however, deliver evenly distributed vector fields.
Hence, combining the two different drift calculation approaches and making use of the respective advantages is recommended to improve high resolution ice drift estimation.Full Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Figure 4 a shows the combined vector fields of the HH and HV channel for ORB, SIFT and SURF respectively.Figure4b and c display N and RMSD on the considered grid.A total of 6920, 1585 and 518 vectors are found using ORB, SIFT and SURF respectively.ORB covers the largest area with many (here > 50) vectors per grid cell and corresponding low RMSD values.
Discussion Paper | Discussion Paper | Discussion Paper |

Figure 1 .
Figure 1.Patch size of descriptor versus average number of matches of the 4 test image pairs from Table1.Solid and dashed lines represent results for HH and HV polarisation respectively.Vertical grey line at 34 pixel (2.72 km) represents chosen parameter.

Table 3 .
Processing time of steps during computation of sea ice drift field from Sentinel-1 imagery, NB: per channel (apart from creating Nansat object).