Geodetic Data for USGS Benchmark Glaciers: Orthophotos, Digital Elevation Models, Glacier Boundaries and Surveyed Positions

Photogrammetric methods form the basis of USGS Benchmark Glacier Project geodetic products. Data processing techniques have evolved over the lifespan of the project, but each technique employs photogrammetry to provide spatially extensive glacier surface elevations (DEMs) in a common reference frame. We produce orthorectified mosaics (Sevara, 2013) from the same images used in DEM construction. Here, we describe the three photogrammetric methods used to produce DEMs in historical order from oldest to newest. Additionally, we describe our approach to extract glacier margins from orthomosaics, as well as the format of attribute tables within each glacier margin shapefile. Full citations of all bibliographic references are provided in the document "GeodeticData_README.pdf" included with this data package.

Photogrammetry/Gridded Products: Both airborne and spaceborne (e.g., WorldView, IKONOS) platforms have been used to acquire stereo optical imagery of glaciers in North America. Below we detail the analyses used to produce geodetic data from imagery sources. Techniques and resolution (pixel size) of gridded geodetic products are listed in tables within this data release. The exact process varies by input, as specified in the following process steps.

Analog Recovery: The oldest available DEMs were derived using analog methods. The original topographic maps were created using stereo-plotters guided by plane table benchmarks. High quality scanning of the maps allowed elevation contours to be manually digitized and interpolated using a natural neighbor interpolation routine (Sibson, 1981). For older DEMS at Sperry Glacier, absolute elevation errors, due to the sometimes poorly resolved historic elevation datums, were corrected (i.e., Florentine et al., 2018). Full citations of all bibliographic references are provided in the document "GeodeticData_README.pdf" included with this data package.

Aerial Imagery: Airborne optical stereo imagery has been acquired over glaciers in North America since the mid-1940s. However, most of these images (e.g., Nolan et al., 2017) lack exterior (i.e., camera location and position) and interior (i.e., camera focal length and radial lens distortion) information. Structure from motion (Sfm) algorithms provide means to optimize both exterior and interior camera orientations based on ground control points (GCPs), in lieu of exterior and interior camera information (Verhoeven et al., 2012). Using opportunistic, manually-identified GCPs (e.g., boulders, buildings, intersecting rock joints) visible in both aerial imagery and high-resolution georeferenced images (typically satellite imagery), and available focal lengths if physically printed on aerial imagery, we optimized camera calibrations and locations, tying Sfm point clouds and the resulting DEMs and ortho images to the landscape (e.g., Maurer et al., 2015, Kienholz et al., 2016). Recent increases in the affordability of high-resolution consumer-grade cameras, and access to precision global positioning systems, has permitted aerial imagery to be acquired in bulk, at relatively low cost since 2014. These additional aerial imagery datasets were obtained via a Nikon D810 camera with a Distagon 25 mm f/2.0 ZF.2 Lens, or a Nikon D850 camera with a Zeiss 25 mm installed in a Cessna 180 with a camera port. A Trimble R7 GPS connected to an externally mounted Sensor Systems L1/L2 antenna recorded raw GPS data at 5 Hz for aircraft positioning. The camera remote flash port was used to trigger event markers in the Trimble R7 data recording, thereby precisely timing each shutter actuation. Although the exterior and interior orientations of recently acquired images are known, DEMs and ortho images were still produced using Sfm to maintain consistency with DEMs derived from older aerial imagery and to leverage the automated workflow of commercial Sfm-based commercial software (e.g., Agisoft Photoscan). Dense point clouds derived from Sfm were directly converted to DEMs, with no interpolation/extrapolation applied to unresolved areas of the DEMs. For ortho images, we derive an interpolated mesh from dense point clouds, which was then used to orthorectify and mosaic images. Full citations of all bibliographic references are provided in the document "GeodeticData_README.pdf" included with this data package.

Aerial Light Detection and Range (LiDAR): Since 2017, increased access to consumer-grade LiDAR platforms, coupled to precision global positioning systems has permitted high accuracy elevation data acquisitions. These data were acquired via a Riegl VQ-580 ii laser scanner utilizing a Applanix 610 Global Network Satellite System mounted on the bottom of a Cessna 180. These data were finally processed using commercial software (e.g., Applanix POSPac; Riegl RiProcess; QTModeler) to process a final DEM, with no interpolation/extrapolation applied to unresolved areas.

Satellite Imagery: Commercially sourced satellite imagery (i.e., Digital Globe) was processed using the open source Automated Stereogrammetry Software Ames Stereo Pipeline (ASP) (Shean et al., 2016). Classified satellite imagery (i.e., National Technical Means) was processed using commercially available software (SOCET SET). DEMs derived from classified sources are made publicly available herein, while the ortho images remain classified and hence, not publicly available. Caution should be used in the analysis of DEMs generated by SOCET SET photogrammetry software as they contain interpolated areas where elevations were unresolved. These areas can be easily identified by inspecting hillshade files derived from the DEM. For commercially obtained satellite imagery, processed using ASP, no interpolation was applied in unresolved areas of the DEMs.

Glacier Boundary Shapefiles: Glacier boundaries were either digitized from original topographic maps (Johnson, 1980) or manually delineated along well-defined regions of the glacier margin (e.g., debris bands or between barren ground and ice) using imagery described in the preceding text. This step was sufficient for glaciers that occupied a well-constrained basin. However, glaciers with ice divides that are shared with other glaciers (Taku), required an extra step. For these glaciers, which also tended to be positioned in perennially snow-covered areas, glacier flow velocity fields (Burgess et al., 2013) were used to define glacier boundaries, i.e., glacier outlines were defined along divergent velocity fields (Keinholz et al., 2015). In these locations, we assumed ice divides were stationary and did not move through time. In years with multiple images, we utilized the image closest to the end of the mass balance year. Additional outlines for Glacier National Park have been published previously by USGS in Fagre et al. (2017). The 2005 Sperry Glacier boundary released here is from this previous publication (Fagre et al., 2017), although all other glacier boundary data are novel. For each glacier boundary, we include the year and specific date the boundary represents, the imagery used to produce the glacier boundary, planar metric area, and the length of the glacier measured along a centerline profile. Glacier margin attribute fields include the name of the glacier, the glacier boundary for each year, the specific date of imagery acquisition, glacier area, glacier length, and imagery source. Full citations of all bibliographic references are provided in the document "GeodeticData_README.pdf" included with this data package.

GNSS Data: Global Navigation Satellite System (GNSS) data was collected on mass balance stakes, at and around mass balance index sites, at equipment installations, and at various sites of interest on and around the Benchmark Glaciers. Two receivers were used to collect simultaneous data at a local survey monument (or base station) and at the survey site, or in some cases a single receiver was used if there is an existing CORS station within 10 km. Stop and go kinematic observations were made at and around index sites with a minimum 30-seconds of data at each stop, all other observations are static occupations with a minimum 8 minutes of data. The local survey monument location is fixed using the average from Online Positioning User Service (OPUS) positions of the base station files. Maximum baseline length is less than 10 km. Positions are reported in the local UTM zone for a given glacier, WGS-84, with ellipsoid heights. The base station files associated with local survey monuments were processed using the Online Positioning Service (OPUS) provided by NOAA's National Geodetic Survey. These files are generally 6-12 hrs in length. Multiple files are averaged to determine the local survey monument location, and the RMSE of those positions is used to estimate the uncertainty in local survey monument position. Survey points were then processed against the local base station. We processed baselines using precise ephemerides in Trimble Business Center. We calculated mean positions from the baseline length and direction relative to the local base-station. Observations without coincident local base station data were removed, as were any positions with RMSE greater than 20 cm horizontally and 50 cm vertically. Following this, manual quality control checks were performed. All of the data was plotted and inspected to ensure logical correctness. For example, stake observations were inspected to make sure the name was correct; names were corrected if the correct name could be determined from the time, position, and field notes; the observation was removed if the correct name could not be determined. Repeat surveys of stakes were inspected to make sure the stake appeared to flow down-glacier. Velocities from multiple surveys on the same stake were plotted and inspected for plausibility. Less than 10 observations were removed based on implausible flow directions and velocities, and these were associated with spring (April or May) surveys of older stakes that field notes suggested were badly bent.