Characterizing Arctic sea ice topography using high-resolution IceBridge data

. We present an analysis of Arctic sea ice topography using high-resolution, three-dimensional surface elevation data from the Airborne Topographic Mapper, ﬂown as part of NASA’s Operation IceBridge mission. Surface features in the sea ice cover are detected using a newly developed surface feature picking algorithm. We derive information regarding the height, volume and geometry of surface features from 2009 to 2014 within the Beaufort/Chukchi and Central Arctic regions. The results are delineated by ice type to estimate the topographic variability across ﬁrst-year and multi-year ice regimes. The results demonstrate that Arctic sea ice topography exhibits signiﬁcant spatial variability, mainly driven by the increased surface feature height and volume (per unit area) of the multi-year ice that dominates the Central Arctic region. The multi-year ice topography exhibits greater interan-nual variability compared to the ﬁrst-year ice regimes, which dominates the total ice topography variability across both regions. The ice topography also shows a clear coastal dependency, with the feature height and volume increasing as a function of proximity to the nearest coastline, especially north of Greenland and the Canadian Archipelago. A strong correlation between ice topography and ice thickness (from the IceBridge sea ice product) is found, using a square-root relationship. The results allude to the importance of ice deformation variability in the total sea ice mass balance, and provide crucial information regarding the tail of the ice thickness distribution across the western Arctic. Future research priorities associated with this new data set are presented and discussed, especially in relation to calculations of atmospheric form drag


Introduction
Sea ice is a heterogeneous medium consisting of level and/or deformed ice floes of various spatial scales (hundreds of metres to several kilometres in diameter), separated by cracks and leads.Given sufficient stresses created by the combined forces of atmospheric/oceanic drag and/or ice-ice interaction, an ice floe can break apart (often along a boundary with another ice floe), and the blocks of newly broken ice are redistributed vertically (e.g.Hopkins, 1998;Feltham, 2008).This pattern of deformation is referred to as a pressure ridge, with the upper surface extension commonly known as a sail, and the lower surface extension (into the ocean) known as a keel.Over first-year ice (FYI), distinct pressure ridges are commonly observed against the backdrop of smooth ice.Distinct pressure ridges can often extend laterally for tens/hundreds of metres across an ice floe.Over multi-year ice (MYI), however, networks of sails and rubble fields (at various stages of weathering) dominate the ice surface.Localized regions of deformation are created through convergent stresses within the ice pack (e.g.ice hummocks), while snow redistribution features also distort the ice surface, caused by erosion (sastrugi) and deposition (dunes) (Thomas and Dieckmann, 2009).Snow drift features can build up alongside sails (snow banks), smoothing their slope and extending their areal coverage.Visual imagery of the sea ice surface and a schematic of a typical FYI floe are given in Fig. 1.
In regions where the sail and keel density is high, the resultant obstructions to fluid flow (form drag) are thought to dominate the total drag on the ice cover over frictional (skin drag) effects (Arya, 1973;Leonardi et al., 2003;Tsamados et al., 2014).Ice deformation also impacts the internal strength of the ice pack, further altering the momentum transfer between the atmosphere and ocean (Martin et al., 2014).The sea ice strength is critical for understanding the resultant loads experienced by ice-breaking ships and offshore structures (e.g.Timco and Weeks, 2010).Dynamical ice redistribution also contributes directly to the total thickness of Arctic sea ice (e.g.Thorndike et al., 1975), although this contribution to ice growth (over thermodynamics) has yet to be reliably quantified.In the Arctic, first-order estimates suggest that deformed ice could contribute up to ∼ 50 % of the total sea ice volume (Wadhams, 2000).The ice topography impacts sea ice melt variability through melt pond formation (e.g.Perovich and Polashenski, 2012), where the flatter (variable) topography of FYI (MYI) promotes shallow but extensive (deeper but less extensive) melt ponds to form on the sea ice surface (e.g.Polashenski et al., 2012).Increased understanding of the sea ice topography is also of interest to the satellite (e.g.ICESat and CryoSat-2) and airborne (e.g.Ice-Bridge) altimeter communities, as the interpretation of radar returns over pressure ridges remains challenging (e.g.Newman et al., 2014).
Studies investigating sea ice morphology in detail (i.e.those resolving distinct pressure ridges at the metre scale) have been based predominantly on airborne and underwater measurements (e.g.Tucker et al., 1979;Wadhams, 1980Wadhams, , 1981;;Tucker et al., 1984;Wadhams and Davy, 1986;Haas, 2004;Martin, 2007;Rabenstein et al., 2010).More recently, Doble et al. (2011) used coincident autonomous underwater vehicle (AUV) sonar and airborne laser profiling to perform a high-resolution, three-dimensional analysis of sea ice morphology; however this was limited to one region within the Beaufort Sea.Efforts have been made to compile existing data sets of pressure ridge morphology (Strub-Klein and Sudom, 2012) and airborne surface profiling (Castellani et al., 2014), to increase spatial and temporal coverage.Unfortunately, these data remain sparse (they do not provide annual data on a basin scale), and are predominantly based on linear profiling of surface features.
In this study, we utilize recent, high-resolution data from the Airborne Topographic Mapper (ATM) laser altimeter, flown as part of NASA's Operation IceBridge (OIB) mission (Krabill, 2013), to provide detailed information regarding the sea ice topography over a variety of Arctic sea ice regimes.IceBridge surveys conducted from Fairbanks, Alaska, acquire data over the predominantly FYI cover of the Beaufort and Chukchi seas, while surveys conducted from Thule and Kangerlussuaq, Greenland, sample the thicker, MYI pack of the Central Arctic, north of Greenland and the Canadian Archipelago.IceBridge offers a vast improvement over previous data sets used to investigate ice topography, due to the combination of high spatial coverage and the use of a conical scanner, which allows profiling of the ice surface in three dimensions.The continuous years of data collection (since 2009) also increasingly provide the potential to assess interannual variability within these two regimes.
The paper is organized as follows: Sect. 2 describes the data used in this study; Sect. 3 discusses the surface feature detection methodology; Sect. 4 presents and discusses the Arctic sea ice topography results; and conclusions are given in Sect. 5.

Data
NASA's OIB mission began collecting airborne observations of the polar regions in 2009, bridging the gap between NASA's Ice, Cloud, and land Elevation Satellite (ICESat) mission which retired in 2009, and the future ICESat-2 mission (Abdalati et al., 2010) scheduled for launch in 2017.OIB aircraft carry a suite of instruments designed to measure both land and sea ice, including their overlying snow cover.In this study, we primarily make use of data obtained by the Airborne Topographic Mapper (ATM) which is a conically scanning laser altimeter operating at 532 nm (Krabill et al., 2002).The ATM laser range and aircraft position/orientation are used to assign three-dimensional geographic coordinates to the point where each laser pulse reflects from the surface.The laser elevation data are referenced to the WGS-84 ellipsoid.
The across-track ATM swath width is determined by the maximum off-nadir scan angle, which is normally fixed at 15 • , giving a swath width of ∼ 250 m, assuming a nominal flight altitude of ∼ 460 m.Note that the scan angle was increased to 23 • in 2010, increasing the swath width.Various statistics regarding the IceBridge sea ice flights and ATM data are shown in Table 1.Each elevation measurement has a footprint of ∼ 1 m and a vertical accuracy of 10 cm or better (Krabill, 2013).Martin et al. (2012) showed that for the IceBridge missions specifically, the ATM has an estimated horizontal accuracy of 74 cm, a vertical accuracy of 6.6 cm and a vertical precision of 3 cm.The high vertical resolution of the ATM makes it well suited to the detection of ridges with a characteristic sail height (upper surface extension of the ridge) of around 1-2 m (e.g.Wadhams, 2000).The shotto-shot ATM spacing is variable (due to the conical scan) and depends on the location of the shot within the swath, includ-ing a negligible shot spacing at the edge of the swath, and a variable shot spacing of several metres around the centre of the swath (Krabill, 2013).The shot spacing at the centre of the swath is determined by the off-nadir scan angle, scan frequency and the plane's altitude, pitch, roll and velocity.
The ATM surface elevation data are routinely used in the retrieval of sea ice freeboard, in conjunction with an automated sea ice lead detection algorithm (Onana et al., 2013) based on coincident optical imagery of the surface from the Digital Mapping System (DMS) (Dominguez, 2010) as described in more detail by Kurtz et al. (2013).The DMS provides geolocated, panchromatic or natural colour imagery that features an image resolution (pixel size) of ∼ 10 cm, assuming a nominal flight altitude of ∼ 460 m, and covers the entire width of the ATM scan.An Applanix POS/AV precision orientation system is used to geolocate and orthorectify the images (Brooks et al., 2012).Sea ice thickness is estimated from the sea ice freeboard using snow depth derived from the on-board snow radar system (Kurtz et al., 2013).The sea ice freeboard, thickness and snow depth product, at a 40 m spatial resolution that includes associated uncertainties, is available through the National Snow and Ice Data Centre (NSIDC) (IDCSI4, Kurtz et al., 2015).Since 2012, Ice-Bridge has also provided the community with a quick-look data product, several months in advance of the standard product release (Kurtz et al., 2013).The IDCSI4 (2009IDCSI4 ( -2013IDCSI4 ( ) and quick-look (2014IDCSI4 ( -2015) ) 40 m spatially averaged sea ice data sets are used in the surface feature-ice thickness regression analysis (Sect. 4.3).However, in this study, we mainly utilize the raw ∼ 1 m horizontal resolution ATM elevation data to characterize the surface profile within the entire ATM swath.We primarily analyse the ATM data from 2009 to 2014 in this study, but also make use of the recently released 2015 data in the surface feature-ice thickness regression analysis (Sect. 4.3).The IceBridge sea ice data coverage over the western Arctic from 2009 to 2014 is shown in Fig. 2.
Since 2011, OIB has also operated a "narrow scan" ATM that features a lower across-track swath width of ∼ 45 m, increasing the shot density in the centre of the swath (Krabill, 2014).These narrow scan ATM data will be combined with the regular ("wide scan") ATM data in specific case studies to assess the potential uncertainties in the surface feature detection from the lower mean spatial sampling of the wide scan ATM.
Visual validation of the surface feature detection scheme is carried out using the DMS imagery, while the POS/AV data are used for accurate geolocation of along-track positioning to determine bounds of evenly spaced ATM sections, as discussed later.
In addition to the OIB data, we use the European Organisation for the Exploitation of Meteorological Satellites (EU-METSAT) Ocean and Sea Ice Satellite Application Facilities (OSI-SAF) sea ice type product (Aaboe et al., 2015).This product provides daily sea ice type classification (open water, first-year ice, multi-year ice) based on the analysis www.the-cryosphere.net/10/1161/2016/The Cryosphere, 10, 1161-1179, 2016 of passive microwave and scatterometry data over the entire Arctic Ocean.We also utilize a data set quantifying the distance to the nearest coastline (http://oceancolor.gsfc.nasa.gov/DOCS/DistFromCoast/) to understand sea ice topography/deformation as a function of coastline proximity.Finally, we use the National Snow and Ice Data Center (NSIDC) regional mask of the Arctic Ocean and surrounding regions (http://nsidc.org/data/polar_stereo/tools_masks)(i) to ensure data are over sea ice, and (ii) to exclude regions (e.g. the Canadian Archipelago) from some of our analyses.

Sea ice topography characterization
There has been considerable discussion in the literature regarding pressure ridges and how they should be defined (e.g.Hibler et al., 1974;Wadhams, 1981;Wadhams and Davy, 1986;Martin, 2007).In this study we employ the elevation threshold approach, which has been used extensively in previous studies (e.g.Wadhams, 1980;Dierking, 1995;Martin, 2007;Tan et al., 2012;Castellani et al., 2014).Typically, a ridge (or surface feature) is detected if it has a height above the local level ice/snow surface greater than a chosen elevation threshold.Different elevation thresholds are then used to differentiate different topographic features of the ice cover.Castellani et al. (2014), for example, used 20 and 80 cm thresholds to differentiate "big" sails from "small" sails/snow features.Sastrugi heights were measured during the Sever airborne program (Warren et al., 1999, Fig. 16b).A maximum sastrugi height of 46 cm (north of Greenland) was suggested based on quadratic fits to in situ observations, meaning elevation thresholds higher than this are likely to exclude purely snow drift features.Results based on the lower elevation threshold mean one can not talk solely of deformation features, due to the likely inclusion of snow features.Alter-natively, higher elevation thresholds could result in the exclusion of a significant fraction of the ice topography variability.The choice of cut-off height can provide a significant impact on the sail/feature height distributions (e.g.Tan et al., 2012) and should be considered when analysing the surface feature data derived in this study.
In this study, we choose to focus on a lower elevation threshold of 20 cm, but also provide summary results and discussion of the analysis using a higher 80 cm threshold.Our results are therefore more representative of the total ice and snow topography variability, which is an important factor when considering the potential impact of these results on estimates of atmospheric form drag over sea ice, an expected utility of this data set in the near future.For simplicity, we refer to all measured topographic snow or ice features in this analysis as "features", instead of ridges or sails.Hibler et al. (1972) discussed the concept of a ridge link as the elementary linear segments composing otherwise complex twodimensional deformation features.In fact, our feature detection algorithm (described in the following subsections) selects connected areas around a local maximum in each structure, and our individual features can therefore be thought of as intermediate quantities between an elementary ridge link and the full ridge structure.Visual inspection across several case studies (not shown) demonstrates that for the higher elevation threshold (80 cm), a linear approximation is more valid than for the features detected using a lower (20 cm) threshold.This idea will be explored further in Sect.4.4.
It is worth noting that these features will likely differ from those detected using linear profiling.For one, the Rayleigh criterion (separating peaks by measuring the depth of the crest between them) is not employed in this study, due to the three-dimensional nature of the data.The relatively wide (∼ 200-300 m) swath width also means we are much more The Cryosphere, 10, 1161Cryosphere, 10, -1179Cryosphere, 10, , 2016 www.the-cryosphere.net/10/1161/2016/likely to pick the peaks of the entire surface feature, as opposed to linear profiling studies, which detect the peak of the surface feature along a random (linear) profile.In regions where surface features are sparse, the two-dimensional ATM scan makes it much more likely that we will detect a surface feature (higher than the chosen elevation threshold).These differences in approach, and the impact on the resultant sail heights especially, are discussed in more detail by Lensu (2003).The feature height distributions in this study are thus likely to differ from those presented previously (e.g.Wadhams, 1980).A recent study by Beckers et al. (2015) explored the difference in surface roughness (standard deviation of relative surface elevation) statistics from linear and scanning airborne laser altimetry, for regions north of Svalbard and in the Fram Strait.They found convergence of surface roughness statistics for sampling distances over several kilometres, especially for the drifting ice sampled north of Svalbard.Unfortunately their surface roughness data are different to the surface feature data presented in this study.Future work will attempt to understand, in more detail, the potential differences between the surface feature distributions presented here, and the surface feature distributions generated by linear profiling.

Feature-picking methodology
The following sections detail the surface feature detection scheme that is visually demonstrated in the case study given in Fig. 3. Further case studies are given in the supplementary information, covering a range of ice types (Figs.S1-S3 in the Supplement).Note that these case studies are based on all individual ATM points within the bounds of the DMS image (∼ 350 m along-track in Fig. 3) for visualization purposes.In the processing of all ATM data (all results presented in Sect.4), the size of each ATM section processed is increased to 1 km along-track.This was a balance between having enough data to reliably estimate a level ice surface (discussed in the next subsection), and a small enough region not be influenced by changes in the sea surface height.The Rossby radius, which indicates the length scale of ocean eddies, is 10 km for typical polar latitudes (Chelton et al., 1998), an order of magnitude greater than the 1 km section length chosen.

Level ice surface calculation
To detect features on the ice surface, we first define a level ice surface.Previous approaches include detecting regions where the ice elevation change is less than some threshold over some along-track distance (e.g.Wadhams and Horne, 1980), or detecting the modal ice surface within a given region (e.g.Williams et al., 2015).In this study, we take a similar approach to the recent, three-dimensional Antarctic study of Williams et al. (2015) and detect the most level ice surface within the relevant section.We calculate the cumulative el- evation distribution of all ATM points within a 1 km section and find the percentile bin (using a bin width of 20 %) with the smallest elevation increase.This is equivalent to finding the modal elevation across percentile bins.The level ice surface calculation is demonstrated in Fig. 3c (and other case studies in the supplementary information).In Fig. 3c, the lowest elevation change is found at 15-35 %, meaning the level ice elevation was taken at the 25th percentile of the elevation distribution, corresponding to a level ice elevation of −8.35 m relative to the WGS-84 ellipsoid.Visual inspection using DMS imagery across a variety of case studies showed that a bin width of 20 % proved to be the most reliable.Maps of the calculated level ice elevation percentile from 2009 to 2014 are shown in the supplementary information (Fig. S4 in the Supplement).In the case of a saddle point, where two shallow elevation gradients are separated by a higher elevation gain (see Fig. S1 in the Supplement for an example), the higher of the surfaces is used, as the lower surface is assumed to come from either a lead or a refrozen lead, which could result in an overestimation of the surface features in these sections.

Data interpolation
All the raw, irregularly spaced ATM elevation data (within each 1 km section) are then projected on to a regularly spaced horizontal grid based on the EPSG:3413 polar stereographic projection (https://nsidc.org/data/icebridge/projections_grids.html), using a linear interpolation scheme.
The level ice elevation is subtracted to convert the data to a regularly spaced grid of elevation relative to the level ice surface (Fig. 3d).We note that, due to the on-ice scan pattern of the ATM, grid cell values are informed by a variable number of raw measurements, wherein the effects of spatial sampling and instrument noise will vary across the gridded elevations.Specifically, the higher shot spacing in the middle of the ATM swath poses a potential for over-interpolation, depending on the horizontal grid resolution chosen.To investigate this in more detail, the shot spacing was analysed for several flights across all OIB years, as summarized in Table 1 and demonstrated in the supplementary information (Fig. S5 in the Supplement).We analysed the near-maximum (99th percentile) spacing in each section, as the maximum spacing is often influenced by isolated ATM points caused by adjacent data drop-out.The mean shot spacing is also shown in Table 1.This demonstrates that across all years (2009-2014), most of the data (99 %) have a shot spacing < 4 m, meaning a horizontal grid resolution of 2 m was chosen (over half the near-maximum spacing).Problems can also occur in interpolation around the ATM swath edge within the convex hull (the maximum region bounded by the corners of the ATM section), especially when the plane deviates from a linear trajectory (sections are not analysed if the pitch and/or roll is greater than a set threshold as discussed in Sect.3.1.4).A kd tree algorithm (Maneewongvatana and Mount, 1999) is therefore used to detect the proximity of the projected ATM data to the raw ATM data.If the nearest raw ATM data point is further than a set distance away (5 m), then that data point is discarded.

Identifying unique surface features
All the gridded ATM elevation data below the chosen feature height threshold (20 cm) are then masked.We scan the masked/gridded ATM data for connected data points using a 3×3 structuring element that considers data points to be connected if they touch adjacently or diagonally.Features which occupy an area less than a set threshold (100 m 2 ) are discarded.The information is still retained in the "bulk" ice topography statistics (area fraction/volume of surface features), as discussed later.
Further segmentation is carried out to increase the geometrical characterization of the surface features.We search each of the connected components for local maxima, and a watershed filter (Soille and Ansoult, 1990) is used to find the shallowest contour that separates each local maxima.These local maxima must be separated from each other (horizontally) by at least 10 m, as in previous studies (e.g.Martin, 2007).This segmentation is highlighted by the large feature in Fig. 3 that has been split into several segments, each dominated by a local maxima.This step is especially crucial when using a relatively low elevation threshold (e.g.20 cm as in most of this study) as large features often merge together around their lower elevation bases.
To understand the impact on the surface feature detection from the choice of grid resolution (2 m) used, Fig. 4 shows the feature detection scheme using a 1, 2 and 4 m grid resolution.Figure 4 also shows the feature detection using the default 2 m grid resolution and incorporating the narrow scan ATM data (discussed in the previous section).The results show negligible visual difference in the gridded ATM data, and only small differences in the calculated feature statistics (number and area coverage).The narrow scan data, while successfully filling in some of the lower spot spacing regions in the middle of the swath (shown visually in Fig. 4b), do not appear to provide much additional value, meaning we choose to proceed solely with the wide scan ATM data for all the analysis presented hereafter.Several other case studies were analysed (not shown) and all demonstrate similar The Cryosphere, 10, 1161-1179, 2016 www.the-cryosphere.net/10/1161/2016/results.Note that Fig. 4 demonstrates the feature detection scheme over a typical 1 km section.

Individual feature and bulk topography statistics
Before proceeding with the ATM processing, the POS/AV data are used to assess the pitch, roll and altitude of the plane within the relevant 1 km ATM section.If the mean pitch or roll is greater than 5 • or the mean altitude of the plane is outside the range 300-700 m (based on the nominal sea ice flight altitude of ∼ 460 m), then the ATM section is not processed.
The number of ATM points within the 1 km section is also calculated as low-lying cloud; leads and ATM malfunction can result in significant regions of ATM drop-out.The mean number of ATM points within a 1 km ATM section (summarized in Table 1) varies from ∼ 40 000 points in 2009 to ∼ 20 000 points from 2011 onwards, when the ATM scan angle and frequency were reduced.We therefore use a threshold of 15 000 ATM points (75 % of the minimum) to ensure reasonable data coverage within each ATM section analysed.
We calculate the surface feature height (h f ) by finding the height of all points within each unique surface feature relative to the level ice surface, and define h f as the peak (maximum) value.We calculate the surface feature area (A f ), which is equal to the number of grid points within each feature multiplied by the square of the grid resolution chosen (2 m).We compute the centre of mass of each feature (r c ), which we use as the feature position.Note that we do not weight the surface feature heights based on their areal coverage (A f ).
While not being a major focus of this study, we also compute the covariance matrix (analogous to an inertia matrix) of each feature as where r i is the position of each point within the unique feature, and the integration is performed over the full feature area.C s and C p are the small (secondary) and large (primary) eigenvalues of C respectively, meaning the ratio R = (C p / C s ) 1/2 gives the degree of elongation of the feature (the ratio of long over short axis, assuming an elliptical shape).We present this analysis to highlight further potential applications of this unique data set, and to demonstrate the impact of the elevation threshold on the geometry of the features detected in this study (Sect.4.4).
Several additional bulk properties of the ice topography are calculated directly within the feature detection scheme.For all (1 km) ATM sections, we collect the (i) mean x/y section location, (ii) ATM swath area coverage (used to estimate ice area, assuming minimal open water), (iii) number of features detected, (iv) feature area coverage (all, including small features < 100 m 2 ), (v) feature area coverage (only large features > 100 m 2 ), (vi) mean surface feature height (all, including small features) and (vii) mean surface feature height (only large features).The volume of surface features per unit ice area, V f , is calculated by multiplying the appropriate mean feature height (with or without small features included) by the total feature area coverage within the section, and dividing by the total swath area (units of m).Note that here we use the mean height of all the points included within each feature to calculate the mean feature height and volume (within each section), whereas in the individual feature height analysis, we take the maximum (peak) height of the feature.Using the maximum feature height has the benefit of being independent of the elevation threshold (if the same feature is detected across different thresholds) and the size of the feature detected.The surface feature height, h f , is thought to be more relevant to form drag calculations (discussed in Sect.4.4).The surface feature volume (per unit area), V f , is an integrator of the size (height and areal coverage) and density of the surface features, and is more of an indicator of the total ice topography variability.
Following Richter-Menge and Farrell ( 2013), we analyse the data within the Central Arctic (CA) region, which extends from 210 to 10 • E and 81 to 90 • N, and in the Beaufort/Chukchi (BC) region, which extends from 190 to 240 • E and 69 to 79 • N, as highlighted in Fig. 5.The CA region is dominated by old (Maslanik et al., 2011) and thick (Laxon et al., 2013) MYI, while the BC region contains a variable mix of FYI and MYI (Maslanik et al., 2011).To delineate the results based on the estimated ice type (FYI or MYI), we take the mean of all daily OSI-SAF ice type data within the dates of the relevant OIB yearly sea ice campaign.For ice to be classified as either MYI or FYI, we require over 80 % of the data at a given grid cell to be consistently one ice type (across the daily range).Locations estimated to include a mixture of FYI and MYI (< 80 % of one ice type) are not included in the delineation, but are still kept in the regional analyses.As the OSI-SAF ice type mask excludes some data along the coast, we assume that all locations along the CA (BC) coast are MYI (FYI) in the absence of an OSI-SAF estimate.Note that a polar hole in the OSI-SAF data prevents some ice type discrimination in this region, as shown in the maps of ice type presented in the following section.The ice type mask is projected onto the relevant data set using a nearest neighbour interpolation scheme.The FYI/MYI coverage from the ATM sections used in this study is summarized in Table 1.
www.the-cryosphere.net/10/1161/2016/The Cryosphere, 10, 1161-1179, 2016 Figure 6 shows the probability distributions of surface feature heights within the CA and BC regions for all features, and for the features estimated as either FYI or MYI, using the OSI-SAF ice type mask (discussed in Sect.3).We also exclude data within the Canadian Archipelago and Fram Strait (using the NSIDC Arctic Ocean mask) from this analysis.Statistics from these distributions are summarized in Table 2.Note that a bin width of 10 cm is used in these probability distributions, although the mean and standard deviation are calculated independently.Before interpreting these distributions, it is worth noting that the spatial sampling in 2009 is lower than all other years (Table 1) and is weighted more towards the thick ice directly north of Greenland (Fig. 5).The sampling in the BC region in 2011 is also noticeably sparse.The spatial sampling increases markedly in 2012-2014, allowing for a more reliable discussion of interannual variability within both regions.
The mean feature height in the CA region decreased from 1.46 ± 0.87 m (2009) to 1.24 ± 0.76 m (2013), before increasing to 1.40 ± 0.85 m (2014).The height of the features estimated as MYI showed a similar pattern, decreasing from 1.54 ± 0.90 m (2009) to 1.25 ± 0.77 m (2013), before increasing to 1.42 ± 0.85 m (2014).The mean FYI feature height across all years was 1.03 ± 0.59 m (0.33 m lower than the MYI mean), with no obvious trend or pattern.As shown in Table 1, the ice estimated as FYI is an order of magnitude lower (1-7 %) than the ice estimated as MYI.The FYI feature height in 2014 was anomalously low (0.70 ± 0.39 m), and the distribution was noticeably skewed towards lower feature heights.However, this distribution was influenced by the low sampling of FYI in the CA region that year (MYI/FYI estimated coverage summarized in Table 1).The The Cryosphere, 10, 1161-1179, 2016 www.the-cryosphere.net/10/1161/2016/FYI that was sampled appears to be located to the north-east of Greenland, near to the ice edge.
In the CA region, the number of features classified as MYI is considerably greater than those classified as FYI (2.6 × 10 6 compared to 0.8 × 10 5 ), meaning the changing topography of the MYI is dominating the response of the CA feature height variability over the small changes in MYI coverage.The modal feature height decreased from 0.65 m (2009) to 0.45 m (2010-2014 mean) in both the MYI and all feature distributions.The modal feature height of the FYI and MYI ice is similar (0.45 m mean), meaning the longer tail of the MYI probability distribution is causing the strong difference in the mean surface feature height.Note that a discussion of potential causes of this interannual variability is provided later, in Sect.4.1.3.
To investigate the tail of the distribution in more detail, Fig. 7 shows the distributions on a log-linear scale, clearly highlighting the exponential nature of the surface feature height distributions found in this study.An ordinary exponential distribution of sail heights was proposed by Wadhams (1980), which has been validated (to varying degrees) by further observations of sail/feature height (e.g.Tucker et al., 1979;Dierking, 1995;Martin, 2007;Rabenstein et al., 2010;Tan et al., 2012).As discussed earlier (Sect.3), the feature heights presented here represent the peaks of the unique twodimensional features, so a direct comparison between these earlier studies (that detect the peak of the surface feature along a random (linear) profile) is not appropriate.Figure 7 demonstrates that a higher probability tail is prominent in the CA region in 2009 and 2014, to a lesser extent.
The MYI and FYI surface feature height is lower in the BC region compared to the CA region, suggesting contrasting responses of the ice within these different regimes.The mean feature height in the BC region still shows a similar interannual pattern, decreasing from 1.14 ± 0.74 m (2009) to 0.94 ± 0.57 m (2013), before increasing to 1.03 ± 0.58 m (2014).This appears to be driven, in part, by the decreasing height of the MYI features (1.33 ± 0.85 m in 2009 compared to 1.07 ± 0.65 m in 2013), although the number of features classified as MYI is of a similar order of magnitude to the FYI in the BC.Overall, the number of features detected within the BC region has increased by almost a factor of 3 since 2009 (3.3 × 10 5 in 2014 compared to 1.2 × 10 5 in 2009), consistent with increased IceBridge coverage.The decrease in 2013 appears to be caused by the decreased coverage of MYI in the BC region, as the MYI feature height remained constant, but the relative quantity of features detected as MYI decreased.In 2014, the feature heights in the BC region classified as MYI and FYI were similar (1.03 ± 0.58 and 1.02 ± 0.58 m respectively).Figure 7 shows a similar exponential distribution in the feature height tail, although the probability of "high" features (> 2 m) is consistently lower www.the-cryosphere.net/10/1161/2016/The Cryosphere, 10, 1161-1179, 2016 than the CA region (a steeper gradient in the log-linear trend), as expected.
The feature height probability distributions for all features classified as either FYI or MYI (independent of region) are shown in the supplementary information (Fig. S6 in the Supplement).The distributions again show clear differences across ice types, with the mean feature height higher for the MYI (∼ 1.3 m) compared to the FYI (∼ 1.0 m), although the modal feature height is similar (∼ 0.45 m) across both ice types.
Table 3 provides statistics of the probability distributions of surface feature height, based on the higher (80 cm) elevation threshold processing.The distributions still show differences across regions, with a higher mean feature height in the CA (2.09 ± 0.74 m) compared to the BC region (1.96 ± 0.67 m), although this difference is significantly less than for the 20 cm results.Again, the mean modal feature height is similar (1.65 m for the CA and 1.55 m for the BC).These results further demonstrate the strong impact on the feature height distributions from the choice of cut-off elevation.

Feature volume variability
Figure 8 shows maps of the mean surface feature volume per unit area (V f ) using the surface elevation threshold of 20 cm.Note that while these results include small (< 100 m 2 ) features, V f excluding small features showed similar results, with differences on the order of 0.01 m (not shown).It is worth noting again that V f differs from the individual feature height analysis as it represents the effective thickness of all surface features (total feature volume in the section spread over the entire swath area) within each 1 km ATM section.Figure 8, however, demonstrates a pattern consistent with the surface feature height analysis, including a higher V f ( 0.15 m) in the CA region, and a lower V f ( 0.15 m) in the BC region.V f is greatest along the Greenland coastline (increasing up to ∼ 0.3-0.4m), especially towards northern Greenland (across most years) and along the eastern Greenland coast within the Fram Strait.V f also increases towards the Beaufort Sea coast in 2012.The regional variability in V f appears stronger than the feature height (h f ) variability.Repeating the analysis for both the mean areal coverage and mean height of features (not shown) demonstrates a roughly equal contribution to the regional volume variability from each term (features increasing in area and height concurrently).
To assess the V f variability across regions and ice type, Fig. 9 shows the probability distributions of V f within the CA and BC regions, for all 1 km ATM sections and for the sections estimated as FYI or MYI.Statistics from these distributions are summarized in Table 4.Note that as these data are based on the 1 km ATM sections (as opposed to individual features), the data sampling is significantly reduced.
In the CA region, V f decreased from 0.19 ± 0.11 m (2009) to 0.15 ± 0.15 m (2013), before increasing to 0.19 ± 0.13 m (2014).Similar to the feature height analysis, the number The Cryosphere, 10, 1161Cryosphere, 10, -1179Cryosphere, 10, , 2016 www.the-cryosphere.net/10/1161/2016/  of sections classified as MYI is over an order of magnitude higher across all years than the FYI (4.2 × 10 4 compared to 0.2 × 10 4 ), meaning the changing topography of the MYI is dominating the response of the V f variability in the CA region (over changes in the MYI coverage), as demonstrated by the coincident variability in the MYI V f .The FYI mean V f (0.11 ± 0.07 m) is lower than the MYI mean V f (0.18 ± 0.12 m) and again shows no discernible trend/pattern.The modal V f in the CA experienced a more variable decline from 2009 to 2014 across both FYI and MYI distributions.In 2010 (all sections) and 2010/2014/all (FYI) the modal V f was 0.01 m, highlighting the prevalence of (1 km) ATM sections with a negligible V f (above the 20 cm elevation threshold).Note that this was not demonstrated in the surface feature height analysis, as the size of the features is not taken into account.This highlights how the threedimensional surface feature volume analysis presented is a more useful indicator of the total ice topography variability, compared to linear transects of peak feature heights, as discussed in Sect.3.1.4.
In the BC region, V f demonstrates a similar interannual pattern, decreasing from 0.11 ± 0.08 m (2009) to 0.06 ± 0.07 m (2013), before increasing to 0.09 ± 0.07 m (2014).Similar to the CA, this appears to be driven, in part, by the decreasing MYI V f (0.16 ± 0.08 m in 2009 compared to 0.11 ± 0.07 m in 2013).The decrease in 2013 appears to be driven by a decrease in the FYI V f and, to some degree, by an increased fraction of FYI sections.In the BC region  4.
in 2014, the MYI and FYI V f are similar (0.09 ± 0.07 m to 0.08 ± 0.08 m).The number of sections has increased by a similar ratio (3) than the increase in features detected, suggesting consistency in the density of features detected.

Potential causes of feature height and volume variability
A recent study by Kwok (2015) provided estimates of the relative contribution to sea ice thickness variability from convergence (dynamics) and melt (thermodynamics) north of Greenland and the Canadian Arctic Archipelago.The strong increase in both h f and V f in the CA region between 2013 and 2014 found in this study is consistent with the strong increase in convergence-driven ice growth (within a similar region) during the preceding summer, estimated by Kwok (2015) using sea ice drift and assumptions of mass conservation.Strong ice convergence was also estimated by Kwok (2015) in December 2008, which may be linked to the high features observed in this study along the CA coastline in 2009.The Kwok (2015) study also showed that variability in the convergence-driven ice growth may be higher than thermodynamic (melt-driven) changes, highlighting the important role of ice deformation variability in the Arctic sea ice mass balance.
In the BC region, the feature height and volume variability is driven more by variability in the presence of MYI.The presence of MYI in the Beaufort Sea (e.g. the tongue of MYI extending from the CA to the southern Beaufort Sea in 2012) is the result of a complex interplay between the impact of the Beaufort Gyre on ice drift (e.g.Hutchings and Rigor, 2012;Petty et al., 2016) and the variable melt-out of ice in the Beaufort/Chukchi region (e.g.Hutchings and Rigor, 2012).The strong increase in the BC MYI coverage in 2014 has also coincided with an overall recovery of older ice across the Arctic since 2013 (see Fig. 4.3a in Perovich et al. (2015), based on ice age data from Tschudi et al., 2015).Both Tilling et al. (2015) and Kwok and Cunningham (2015) showed an increase in Arctic sea ice volume in 2014, linked to the retention of MYI.
While our study provides information regarding the surface feature variability, the underside extension of the pressure ridge system, the keel, is thought to be significantly larger in size (e.g.Wadhams, 2000).Strub-Klein and Sudom (2012) recently compiled and analysed several ridge morphology data sets collected over the last few decades.They demonstrated that, on average, the maximum keel depth is around 4 times larger than the maximum sail height, while the keel width is around 6-7 times wider than the sail width.This suggests a keel volume up to ∼ 20-30 times larger than sail volume.The changes in surface feature volume, V f , demonstrated in this study (± 0.05 m) suggest, to a first-order approximation, total deformation variability up to ∼ 1 m, if the keels are taken into account.This simple estimate assumes minimal impact from snow redistribution variability, which will act to reduce the magnitude of this estimate.Unfortunately, detailed information regarding snow variability (spatial and temporal) over Arctic sea ice is lacking.

Sea ice topography as a function of coastline proximity
The maps of surface feature height, h f , (Fig. 5) and mean surface feature volume, V f , (Fig. 8) suggest a strong relationship between surface feature variability and coastline proximity.The convergent, onshore ice drift in the CA region is thought to contribute significantly to increases in ice deformation and thickness across much of this region (e.g.Kwok, 2015).The increased age of the ice along the CA coast (Maslanik et al., 2011) may also provide more time for the ice to thicken through both thermodynamic and dynamic processes.In the BC region, the mean winds (Fig. 2) and ice drift are aligned more parallel to the coast, suggesting less of an impact from convergent coastal boundary stresses.The ice along the BC coastline is also driven, in part, by variability in the import of thicker, older ice from the CA region, linked to variability in the ice circulation around the Beaufort Gyre (e.g.Hutchings and Rigor, 2012;Petty et al., 2016), as mentioned in the Introduction.Tucker et al. (1979) also discussed how relatively thin ice around the landfast ice edge (compared to the thicker, grounded ice, nearshore) can undergo significant ice deformation.Mahoney et al. (2014) used Radarsat Synthetic Aperture Radar (SAR) imagery to show that the landfast ice edge can extend up to 100 km offshore of Alaska, The Cryosphere, 10, 1161-1179, 2016 www.the-cryosphere.net/10/1161/2016/although the results also suggest significant spatial and temporal variability in the width of this BC landfast ice regime.
A detailed analysis of specific IceBridge flight lines (in isolation) is therefore likely needed to detect and estimate the ice topography around the variable landfast ice edge.
In this study, we more broadly analyse the coastal dependency of the surface feature height, h f , and mean surface feature volume, V f , data presented in the previous section.Figure 10 shows h f represented by box and whisker plots, separated into coastline proximity bins (100 km wide) for the BC and CA regions.The coastal proximity data were presented in Sect. 2 and a map of the coastline proximity is given in the Supplement (Fig. S7).Note that less weight should be given to the BC results as there are much fewer data near to the coast (the period of 2012-2014 has the highest coverage of data near to the BC coastline).It is also worth noting that the CA coastal region (northern Greenland and the Canadian Archipelago) is dominated by MYI, whereas the BC coastal region (northern Canada and Alaska) shows greater interannual variability in the dominant ice type, as discussed previously.
Despite the consistent presence of MYI over much of the CA region, Fig. 10 demonstrates a strong increase in h f with increasing coastline proximity (in terms of the 25th, 50th 75th and 95th percentiles) up to 900 km away from the coast.The 0-100 km bin shows a significant fraction (∼ 5 %) of fea-  5), presented using box and whisker plots (5, 25, 50, 75, 95 percentiles).The coastline distance bin width is 100 km.The black boxes (and whiskers) show the results from all features detected in each region, while the colours represent the results from each year.
www.the-cryosphere.net/10/1161/2016/The Cryosphere, 10, 1161-1179, 2016 tures higher than ∼ 3.3 m, compared to the distance bins further from the coast.The results show moderate interannual variability, with 2009 showing higher features (compared to the other years) from 0 to 200 km from the coast, while 2014 shows higher features from 100 to 800 km from the coast, highlighting that the increase in surface feature height in 2014 manifested over much of the CA region, while in 2009, the high surface features were contained mostly along the CA coastline.
The BC region also demonstrates an increase in surface feature height with increasing coastline proximity, although this is mainly observed in the upper percentiles (75th and 95th) of the distributions.The median feature height across the 0-400 km percentile bins shows higher variability than the CA region.The 95th percentile results from 0 to 300 km are lowest in 2013, which may be due, in part, to the thin level ice sampled in the Chukchi Sea north of Point Hope in 2013 (Richter-Menge and Farrell, 2013).The feature heights also tend to increase (across most percentile ranges) at distances greater than 700 km away, which is likely due to the import of MYI from the CA into the northern Beaufort Sea.
The surface feature volume (per unit area), V f , results, shown in Fig. 11, demonstrate a similar and perhaps more obvious coastline relationship.In the CA region, 2009 and 2014 show increases in V f closer to the coastline, similar to the feature height results discussed previously.The median V f across all distance bins shows greater interannual variability compared to h f .In the BC region, the V f increase towards the coast (75th and 95th percentile) is much clearer than the h f results, and the interannual variability is again higher.This suggests that the three-dimensional surface feature volume data are a more useful measure of coastal topographic variability compared to the surface feature heights, especially compared to data compiled from linear transects.Note that reducing the bin width to 10 km and analysing the coastline dependency on this smaller scale did not demonstrate any obvious landfast ice zone (a steep gradient in ice topography) across either region.

Relationship between sea ice thickness and surface feature variability
The relationship between sail height and sea ice thickness has been discussed in several previous studies of sea ice pressure ridging, with varying conclusions drawn.Tucker and Govoni (1981) were perhaps the first to observe the link between sail heights and the thickness of the ice blocks from which they formed, which they assumed to be representative of the parent ice thickness.A square-root relationship was presented, which was validated by additional in situ observations (Tucker et al., 1984) and the two-dimensional particle modelling study of Hopkins (1998).More recently, Martin (2007) found only a weak correlation between sail height and the parent ice thickness using a variety of linear surface profiling data sets and assuming a similar square-root rela- , where b is the calculated regression coefficient, and r is the correlation coefficient for all years of data.σ r is the standard error of the residuals (or root mean squared error) calculated using all years of data.tionship.A stronger, but still only moderate correlation was found when a linear fit was assumed.
Figure 12 shows the correlation between the surface feature height, h f , derived in this study, and the total sea ice thickness, H i , taken from the IceBridge sea ice thickness product (IDCSI4 from 2009 to 2013 and quick-look in 2014, as described in Sect.2).Both data sets are averaged over 10 km (along-track) sections to smooth the data, and the Ice-Bridge thickness data are interpolated onto the mean surface feature height sections using linear interpolation.It is worth noting that the sea ice thickness data used in these regressions are calculated using measurements of sea ice freeboard The Cryosphere, 10, 1161-1179, 2016 www.the-cryosphere.net/10/1161/2016/and assumptions of hydrostatic equilibrium and thus implicitly include the deformed and undeformed ice, meaning h f and H i are not truly independent variables.The regressions are therefore expected to differ from those presented in previous analyses, that correlated the sail heights with the thickness of the ice blocks within the ridge (e.g.Tucker et al., 1984) or the level ice thickness directly (e.g.Martin, 2007).Our likely inclusion of snow drift features, and the expected thermodynamic/dynamic changes over time of surface features (we are not measuring the features as they are formed) will also impact these correlations, and weaken the physical links to pressure ridging constraints.We therefore do not attempt a validation of the square-root relationship found in previous studies, but instead seek to quantify the relationship between the peak surface feature heights found in this study and the local (total) sea ice thickness.Our analysis is therefore more in line with the regressions presented in Beck-ers et al. (2015), between the total ice (plus snow) thickness and their estimated surface roughness (introduced in Sect.3).In that study, a strong (negligible) linear correlation was found over the deformed (drifting) sea ice, although these regressions were limited by their considerably lower spatial/temporal sampling compared to the data presented in this study.
The regressions between the surface feature height, h f , and the total ice thickness, H i , are shown in Fig. 12.We experimented with both linear and square-root relationships using a least-squares fit, and slightly stronger correlations were found with the latter.The square-root relationships also cross the origin and are thus more physically consistent (the linear correlations without a variable y intercept were markedly weaker), so we decided to present and focus on the squareroot regressions in this study, h f = b √ H i , where b is a regression coefficient calculated from the least-squares fit.
www.the-cryosphere.net/10/1161/2016/The Cryosphere, 10, 1161-1179, 2016 The regression using all years of data (2009-2014) demonstrates strong correlation (r = 0.72, b = 0.72).The annual regressions (given in Fig. 12) show that the strongest correlation is observed in 2013 (r = 0.81).Strong correlations are observed across all other years (r = 0.67-0.76),except for 2009, where only a weak (r = 0.35) correlation, and a regression parameter higher than average (b = 0.86) was found.This may be due to the decreased ATM coverage in 2009, although Fig. 12 suggests that the ice thickness results were also skewed low compared to the relationships demonstrated across all other years.Note that changing the averaging length scale (5 and 20 km) resulted in weaker correlations.In general, the consistency of these regressions (similar b value) across different years (2010-2014) suggests consistency in the response of the ice to dynamical forcing.
To demonstrate the potential utility of these findings, Fig. 13 shows the sea ice thickness from the derived Ice-Bridge product, and the sea ice thickness estimated using the surface feature height, h f , and the relationship h f = b √ H i (both using the 10 km mean data).Here we use the recently released 2015 ATM data to calculate the surface feature height (not presented earlier), and the 2015 quick-look sea ice thickness data.A regression parameter of b = 0.72 was used based on the regression analysis across 2009-2014.The maps qualitatively show the close correspondence between the spatial variability in ice thickness across the CA and BC regions from both the IceBridge product and the ice thickness estimated from h f .Differences between observed and predicted ice thickness are up to ± 2 m in some regions, although this is within range of the combined root mean squared error of the regression (1.10 m, given in Fig. 13) and the mean 10 km IceBridge thickness uncertainty of 0.8 m (calculated from the raw IceBridge uncertainty estimates across all years).In general, the results provide a useful means of understanding ice topography and thickness variability in more detail, and demonstrate how the ice thickness estimates could provide a useful proxy for ice thickness, especially in regions where measurements of leads, which are needed to calculate sea ice freeboard, are sparse.

Feature geometry and the potential for additional feature characterizations
As discussed in the introduction, sea ice topography is crucial for estimating atmospheric form drag over Arctic sea ice.Calculations of atmospheric form drag require estimates of the surface feature height (as presented in this study), along with the surface feature density (e.g.Arya, 1973;Tsamados et al., 2014).Linear profiling studies calculating atmospheric form drag (e.g.Castellani et al., 2014) simply measure the spacing between unique surface features along the linear profile, assuming that the features are randomly orientated and sufficiently sampled for this assumption to be valid.Mock et al. (1972) showed that for randomly oriented ridges, the average ridge frequency, µ, and the average ridge density (the ratio of the total length of ridges per unit area), R D , are related via µ = (2/π )R D .In contrast to linear profiling studies, R D can be calculated directly with these data as R D = i L i / A tot = L tot / A tot , where the sum is over all features within the total ice/swath area (given a fully concentrated ice pack).Assuming an elliptically shaped feature, the length of the major axis of a specific feature can be estimated as L i = 2 √ π (RA sf ) 0.5 , where R = (C p / C s ) 0.5 is the degree of elongation of the feature, as mentioned in Sect.3.An average spacing between features can then be estimated from R D as X f = π / (2R D ).
A crucial factor in this calculation is the assumption of linear features in the estimation of ridging density.Figure 14 shows the mean aspect ratio (R) of all features detected across 1 year (2012) using the 20 and 80 cm elevation cut-off thresholds.For the 20 cm elevation cut-off (as used throughout much of this study), the aspect ratio of all features appears to be ∼ 2-2.5 : 1, while for the 80 cm threshold, the estimated aspect ratio is ∼ 3-4 : 1.The assumption of linearity is somewhat arbitrary, but is clearly more questionable in the 20 cm case.We have decided not to present calculations of ridging density and form drag estimates as we believe a more thorough analysis is needed, which is beyond the scope of this current study.Understanding the surface feature geometry variability, and linking this with estimates of feature density relevant to form drag parameterizations and also melt pond formation, will be a crucial next step in the utility of this unique, three-dimensional sea ice topography data set.

Conclusions
We have presented a detailed characterization and analysis of Arctic sea ice topography using high-resolution, threedimensional surface elevation data from the Airborne Topographic Mapper, flown as part of NASA's Operation Ice-Bridge mission.Surface features in the sea ice cover (caused by ice deformation and/or snow redistribution) are detected using a newly developed feature-picking algorithm.We derive information regarding the individual height and volume (per unit area) of surface features from 2009 to 2014 within the Beaufort/Chukchi and Central Arctic regions, across both first-year and multi-year ice regimes.
The results demonstrate that Arctic sea ice topography exhibits significant spatial variability, mainly driven by the increased surface feature height and volume in the multi-year ice.The multi-year ice topography also exhibits greater interannual variability compared to the first-year ice topography.Multi-year ice dominates the Central Arctic region and contributes significantly (but variably) to the Beaufort/Chukchi region.The tail of the surface feature heights (> 2 m) exhibits a clear exponential distribution, further validating previous observational studies.The ice topography also shows a strong coastal dependency, with the feature height increasing as a function of proximity to the nearest coastline, es- The Cryosphere, 10, 1161Cryosphere, 10, -1179Cryosphere, 10, , 2016 www.the-cryosphere.net/10/1161/2016/pecially north of Greenland and the Canadian Archipelago.The coastal proximity results provide useful context regarding interannual variability in the location of surface topography features.A strong correlation between surface feature height and ice thickness (from the IceBridge sea ice product) is found, based on a square-root relationship.The consistency of these regressions across different years (2010)(2011)(2012)(2013)(2014) suggests consistency in the response of the ice to dynamical forcing.Overall, the results allude to the importance of regional and interannual ice deformation variability in the total sea ice mass balance, and provide crucial information regarding the tail of the sea ice thickness distribution across the western Arctic.While this study presents the use of IceBridge data to understand the Arctic sea ice topography, future work will attempt to understand the impact of ice topography on estimates of atmospheric form drag. Another exciting prospect involves the extension of this analysis to Antarctic sea ice, where observations of the sea ice state are extremely lacking.

Figure 1 .
Figure 1.(a) Aerial photograph of the sea ice surface, taken off the coast of Barrow, Alaska.(b) Schematic of a sea ice floe (not to scale) featuring two large pressure ridges, one smaller pressure ridge and a sastruga (snow feature).

Figure 2 .
Figure 2. Top: IceBridge sea ice flight lines and estimated ice type over the western Arctic.The dark grey (light grey) background indicates regions where more than 80 % of the daily data within all IceBridge sea ice campaign dates (across all years) are estimated as MYI (FYI), while the medium grey indicates a mix of FYI and MYI, taken from the EUMETSAT OSI-SAF ice type mask.The coloured stars indicate locations of the various case studies, as highlighted in the relevant figures.Bottom: mean winds from January to March (2009-2014) taken from the ERA-I reanalyses.

Figure 3 .
Figure 3. Example of the surface feature detection algorithm overlaid on a DMS image taken on the 23rd March 2011 as highlighted by the yellow star in Fig. 2. (a) DMS image; (b) raw ATM data overlaid on the DMS image; (c) elevation distribution for all ATM points within the section shown, where the blue line indicates the bounds of the calculated level ice surface and the red line indicates the feature height threshold; (d) gridded (2 m) ATM elevation relative to the level ice surface; (e) unique surface features (> 20 cm) and their elevation relative to the level ice surface; (f) unique surface feature identifier (features larger than 100 m 2 ).

Figure 4 .
Figure 4. Example feature detection algorithm for a 1 km ATM section in 2011.The top row shows the raw ATM data from both the regular wide scan (a) and from the combined wide and narrow scan (b).Panels (c-e) show the features detected using a 1 m (c), 2 m (d) and 4 m (e) gridding of the regular wide scan ATM data, while (f) shows the results from the 2 m gridding of the combined wide and narrow scan ATM data.Panels (c-e) also show the number of surface features (> 20 cm) detected and the total area of these features.

Figure 5 .Figure 5
Figure 5. Surface feature height, h f , from 2009 to 2014, detected using a 20 cm elevation threshold.The dark grey (light grey) background indicates regions where more than 80 % of the daily data within each year's IceBridge sea ice campaign dates are estimated as MYI (FYI), while the medium grey indicates a mix of FYI and MYI.The red (blue) dashed lines represent the Central Arctic (Beaufort/Chukchi) regions used in this study.The data are plotted using hexagonal bins.

Figure 6 .Figure 7 .
Figure6.Probability distributions of the surface feature height, h f , (using a 20 cm elevation threshold) detected within the (a) Central Arctic and (b) Beaufort/Chukchi regions (shown in Fig.5) and for the features estimated as FYI (c and d) or MYI (e and f) using the OSI-SAF ice type mask described in Sect.3. The bin width is 10 cm and the bin values are plotted as lines (joining each value) instead of steps for clarity.The solid (dashed) vertical lines show the mean (mode) of the distributions across each year.The statistics (mean, mode and standard deviation) are summarized in Table2.

Figure 8 .
Figure 8.As in Fig. 5 but for the surface feature volume (per unit area), V f .

Figure 9 .
Figure 9. Probability distributions of the surface feature volume (per unit area), V f , (using a 20 cm elevation threshold) within the (a) Central Arctic and (b) Beaufort/Chukchi regions (shown in Fig. 8) from 2009 to 2014.The bin width is 1 cm.The statistics (mean and mode of each distribution) are summarized in Table4.

Figure 10 .
Figure10.Surface feature height, h f , as a function of distance to the nearest coastline within the (a) Central Arctic and (b) Beaufort/Chukchi regions (given in Fig.5), presented using box and whisker plots(5, 25, 50, 75, 95 percentiles).The coastline distance bin width is 100 km.The black boxes (and whiskers) show the results from all features detected in each region, while the colours represent the results from each year.

Figure 11 .
Figure 11.As in Fig. 10 but for the surface feature volume (per unit area), V f .

Figure 12 .
Figure 12.Correlation between the mean IceBridge sea ice thickness and surface feature height h f , averaged over 10 km along-track sections.The solid lines represent the least-squares fit, assuming a square-root relationship(h f = b √ H i ), where b is the calculated regression coefficient, and r is the correlation coefficient for all years of data.σ r is the standard error of the residuals (or root mean squared error) calculated using all years of data.

Figure 13 .
Figure 13.The 2015 surface feature height, h f , (top left) and the estimated sea ice thickness using the relationship h f = 0.72 √ H i given in Fig. 12 (top right).The bottom left panel shows the quick-look IceBridge ice thickness results, and the bottom right panel shows the difference between the ice thickness estimated in this study and the derived IceBridge thickness (bottom right).

Figure 14 .
Figure 14.Surface feature aspect ratio, R, detected using a 20 cm elevation threshold (left) and 80 cm threshold (right) in 2012.The dark grey (light grey) background indicates regions where more than 80 % of the daily data within each year's IceBridge sea ice campaign dates are estimated as MYI (FYI), while the medium grey indicates a mix of FYI and MYI.The red (blue) dashed lines represent the Central Arctic (Beaufort/Chukchi) regions used in this study.

Table 1 .
IceBridge ATM flight information.Note that all calculated quantities (rows 3-13) are based on the permissible sea ice sections, as described in Sect.3.The ice type classification is also described in more detail in.

Table 2 .
Surface feature height statistics (mean and mode) taken from the probability distributions shown in Fig.6.The value in the brackets (next to the means) equals 1 standard deviation of the relevant distribution.The third column (under each region) shows the number of surface features detected.

Table 3 .
Surface feature height statistics, as in Table2but for the processing using an 80 cm threshold, with the results using all features in each region shown (not delineated by ice type).

Table 4 .
Surface feature volume statistics taken from the probability distributions shown in Fig.9.The value in the brackets (next to the means) equals 1 standard deviation of the relevant distribution.The third column (under each region) shows the number of 1 km ATM sections used in each distribution.