Terrain changes from images acquired on opportunistic flights by SFM photogrammetry

Acquiring data to analyze change in topography is often a costly endeavour requiring either extensive, potentially risky, fieldwork and/or expensive equipment or commercial data. Bringing the cost down while keeping the precision and accuracy has been a focus in geoscience in recent years. Structure from motion (SfM) photogrammetric techniques are emerging as powerful tools for surveying, with modern algorithm and large computing power allowing for the production of accurate and detailed data from low-cost, informal surveys. The high spatial and temporal resolution permits the monitoring of geo5 morphological features undergoing relatively rapid change, such as glaciers, moraines, or landslides. We present a method that takes advantage of light-transport flights conducting other missions to opportunistically collect imagery for geomorphological analysis. We test and validate an approach in which we attach a consumer-grade camera and a simple code-based Global Navigation Satellite System (GNSS) receiver to a helicopter to collect data when the flight path covers an area of interest. Our method is based and builds upon (Welty et al., 2013), showing the ability to link GNSS data to images without a complex 10 physical or electronic link, even with imprecise camera clocks and irregular time-lapses. As a proof of concept, we conducted two test surveys, in September 2014 and 2015, over the glacier Midtre Lovénbreen and its forefield, in northwestern Svalbard. We were able to derive elevation change estimates comparable to in-situ mass balance stake measurements. The accuracy and precision of our DEMs allow detection and analysis of a number of processes in the proglacial area, including the presence of thermokarst and the evolution of water channels. 15


Introduction
High spatial and temporal resolution surveys produced by structure from motion (SfM) photogrammetric techniques are emerging as powerful tools for studying geomorphological objects such as glaciers, moraines, and landslides.Modern software and large computing power allow us to produce accurate and detailed data from relatively low-cost surveys that can be conducted with little planning.In this study, we present a method to take advantage of light-transport flights to collect imagery for precise digital elevation model (DEM) and orthoimage generation, to be used in geomorphological analysis.
The typical SfM pipeline (Fig. 1) produces a georeferenced DEM and orthoimage from input imagery and geolocation data (Snavely et al., 2008;Niethammer et al., 2010;Westoby et al., 2012;Kääb et al., 2014;James and Robson, 2014;Ouédraogo et al., 2014;Tonkin et al., 2015;Nolan et al., 2015;Galland et al., 2016;Eltner et al., 2016).Typically, the camera's intrinsic parameters (sensor and optics) are determined first, along with the relative camera locations and pointing angles from the photogrammetry itself.Then, by combining this information with surveyed camera position (from embedded GNSS data) and/or ground control points (GCPs), an accurate georeferenced reconstruction of the terrain is generated through dense multi-view stereo matching.
Our method exploits the malleability of SfM photogrammetry, compared to classical photogrammetry, regarding flight path and aircraft stability.In this study, we test and validate an approach to attach simple cameras and consumergrade GNSS devices to an aircraft on other missions, with Published by Copernicus Publications on behalf of the European Geosciences Union.the aim of collecting data when the flight path is over an area of interest.This paper is motivated and grounded in the concepts introduced in "Cameras as Clocks" (Welty et al., 2013) and herein we aim to present a practical, hands-on application, assessment, and improvement of these concepts within SfM photogrammetry.In addition, we aim to generate accurate DEMs without the use of GCPs, as these may not be available at high enough precisions over regions covered by opportunistic flights.The novelty in the method lies in the way the GNSS information is connected to the imagery, without the physical or electronic link typically found in highprecision purpose-built hardware.The virtual absence of associated costs allows for failed trials and the collection of a surplus of images that could also be used for other projects if properly archived.
We first present our generic method and then our test area and the specific survey (data acquired, specific equipment used, data processing, and additional data sets for comparison).We use this method to generate two DEMs from our low-cost equipment.The accuracy of the generated DEMs is analysed based upon comparisons with professional surveyed DEMs.DEMs are compared in regions considered stable and statistics are used for quantifying their accuracy.In summary, we aim to show the improvements gained by including a few additional steps within the common SfM pipeline.We close with a brief initial geomorphological analysis and interpretation conducted from our products.
The geomorphological research goals are to detect and quantify surface changes of glacial and pro-glacial surfaces.In particular, the measurement of glacier geodetic balance from glacier elevation changes complements traditional in situ measurements of mass balance (e.g.Zemp et al., 2010;Andreassen et al., 2016), for instance by detecting bias in the in situ estimates and providing basis for calibrating or correcting such biases (Zemp et al., 2013).Associated with glacier wastage and climatic change, a number of other geomorphological processes occur in the pro-glacial environment, for example related to buried ice, frozen ground, and hydrological processes.In this study, we present methods for easily collecting DEMs over glaciers and the surrounding terrain, for example, combined with annual mass balance surveys that often require helicopter or fixed wing transport.
2 Generic method

Acquisition hardware and flight plan
The first step is data acquisition, following the general principle that it should allow opportunistic data gathering.The associated hardware list and flight path recommendation presented below is not specific; rather, it includes the key characteristics that should be satisfied.The actual equipment used for our study is described in Sect.4.1.
Cameras: the cameras and lenses used should be easy to attach to any potential aircraft (or come with mounting kits), be weather-resistant, have good battery life, have good image quality (reasonably sharp, good dynamic range, and low electronic noise), and have an integrated time-lapse function.The characteristics of a specific camera might make it adapted for a certain survey and not another.
GNSS: a low-cost consumer-grade code-based GNSS receiver can be used for this method (e.g. a Garmin-style handheld hiking GPS).It must be able to log positions at a relatively high frequency (1 Hz is an acceptable value).
Aircraft: by definition, control over the type of aircraft is limited in our opportunistic method.The aircraft just needs to have space to securely attach the camera in an adequate position (i.e. the camera should look straight down, and have a free view to the survey object), and there should be a place for the GNSS tracker where signal reception is acceptable.
Flight trajectory: for this method to work, the flight trajectory cannot be a straight line at constant speed.Indeed, the process described in Sect.2.3 is based on the fitting of two independent estimate of the trajectory of the flight and fitting a line to a line has a degree of freedom that cannot be solved.If a typical photogrammetric flight path composed of several lines of flight is not an option (because of legal restrictions or time constraints for instance), a zigzag in the trajectory or a couple of gently banked turns should give enough geometry to the acquisition for the method to perform.

Data processing workflow
Once data are acquired, they are processed through a typical SfM photogrammetric pipeline (see Fig. 1).We use the free and open-source integrated photogrammetric library MicMac, developed at the French National Institute of Geographic and Forest Information (IGN) (Pierrot-Deseilligny and Clery, 2011;Pierrot-Deseilligny et al., 2016;Pierrot-Deseilligny, 2016).MicMac is a highly versatile tool in which the user is in control of the processing parameters and can check intermediate results, such as the calibration of the cameras, the relative camera positions, and the residuals of the absolute orientation transformation.Such capability is crucial for our method because the georeferencing information acquired by the (not quite embedded) GNSS system needs to be linked to the images through the non-trivial process discussed below.

Link between the cameras and the GNSS tracker
To estimate absolute camera positions, each picture must be linked to the GNSS track.This should be possible to do by using the time tag in the images' metadata and matching it with the position recorded in the GNSS log for that time.However, camera clocks are neither precise nor accurate (Welty et al., 2013).Typically, time is only recorded with a 1 s precision (some higher-end cameras do embed a tag called SubSecTimeOriginal that provide sub-second time stamps, with varied level of precision).

Existing solution to the problem
In Welty et al. (2013), the issue of absolute accuracy (GNSS Time = Camera Time ) is tackled by assessing the time offset between the GNSS tracker clock and the camera clock (Delay GNSS-Camera ) : where Time i and TimeEXIF i are, respectively, the estimated time of acquisition and the recorded EXIF time tag of image i.The parameter Delay GNSS-Camera can be easily estimated by taking a picture of the screen of the GNSS receiver, by comparing the time displayed on the GNSS receiver image to the time stamp in the image metadata, or by setting the camera time as accurately as possible shortly before the flight, such that Delay GNSS-Camera ≈ 0. Of course, both methods only yield an approximate value that is, for instance, only as accurate and precise as the GPS receiver display.The actual value needs to be refined further.This refinement is performed by minimizing the rootmean-square (RMS) difference between the fitting of the camera's relative positions, as solved from the SfM photogrammetry, to the simulated camera positions along the GPS track assuming different values of Delay GPS-Camera .For each set of such simulated positions, a seven-parameter transformation is computed: where w C rs , w R rs , and λ are the translation vector, rotation matrix, and scaling factor between the relative and absolute coordinate systems, and rs L and w L are the coordinates of a point in the two systems.The parameters computed through least-square adjustment provide the best fit between the two trajectories, and the distances between the simulated positions and the positions resulting from the adjusted transformation provide an estimate of the accuracy of the simulated positions.The fitting error between the position of each point in the time-fitted GPS data and the final estimated camera position is used as the quality indicator and the solution yielding the lowest error provides the value of Delay GPS-Camera .

Improvement of the method
An additional consideration is that even if we could expect a camera time-lapse mode to take pictures regularly, it is hardly the case.We performed in-lab measurements on a static Go-Pro H3+BE pointed at a 60 Hz screen displaying a timer with 0.01 s precision, using a sufficiently fast memory card to prevent the camera's buffer from filling, and found a 0.1 ± 0.1 s variation in the actual lapse time between two pictures, for the constant lapse time programmed in the camera.If the memory card is not fast enough to write a constant flow of images, the buffer of the camera often gets filled, yielding additional irregular delays on the time between images (up to one extra second in our case).Plotting TimeStamp expected vs. TimeStamp EXIF -TimeStamp expected (the deviation from the expected time stamp) shows jumps of 1 s every time the accumulated delay exceeds 1 s (visible in the red line in Fig. 3) as well as irregularly spaced jumps when the buffer gets filled (visible in the blue line).The method presented by Welty et al. (2013) does not solve for sub-second variations or error in the EXIF time stamp, which can lead to local strains on the solution proportional to the flight speed.
Given that the EXIF time stamps can not be trusted, an estimate of the delay between the camera clock and the GPS clock for the first image (Delay GPS-Camera0 ) as well as the actual lapse times are needed to approximate the time stamp Time i at sub-second precision for each picture.We therefore need to find the parameters in where TimeEXIF 0 is the camera clock for the first image and T 0−i the time between the acquisition of the first picture and the ith picture.T 0−i need to be estimated robustly in cases where the lapse times are inconsistent.By performing a piecewise linear regression in the ImageNumber vs. TimeStamp EXIF space, we can estimate the drift of the www.the-cryosphere.net/11/827/2017/The Cryosphere, 11, 827-840, 2017 sub-second camera time between each image.The value of Delay GPS-Camera0 is then estimated in the same way as in Welty et al. (2013).Figure 2 shows the mean residuals of such a fit as a function of different Delay GPS-Camera0 values in our 2015 survey (Sect.4).

Study area
Our study area is in the vicinity of the Ny-Ålesund research base, in northwest Svalbard (78 • 53 N-12 • 03 E; Fig. 4).Ny-Ålesund is characterized by a comparatively mild but nonetheless Arctic climate, with a mean annual temperature of ≈ −6 • C (Førland et al., 2012), and continuous permafrost up to ca. 400 m thick (Liestøl, 1977;Humlum et al., 2003).Glaciers cover ca.60 % of Svalbard's land area and comprise a variety of ice caps, tidewater glaciers, valley glaciers, and cirque glaciers.Svalbard glaciers are typically polythermal (Björnsson et al., 1996), with a temperate accumulation area and cold marginal areas, where the ice is frozen to the bed.This thermal regime gives rise to glacier-marginal land systems in which the glacier's forefield is underlain by permafrost and remnants of glacier ice from previous glacier advances (e.g.Etzelmüller and Hagen, 2005).Icecored moraines predominate and are at present in different stages from stable to areas experiencing downwasting and thermokarst processes (e.g.Etzelmüller, 2000).
Our study site is Midtre Lovénbreen, a well-studied glacier with a wealth of historic data.It is ≈ 5 km 2 and is situated in a north-facing catchment on the Brøggerhalvøya peninsula.Midtre Lovénbreen has one of the longest continuous mass balance records in the High Arctic (Hagen and Liestøl, 1990;Kohler et al., 2007).As part of the long-term mass balance monitoring, the glacier has a system of stakes frozen into the ice.Precise positions of these stakes have been measured using dual-phase GNSS equipment biannually since 1999.The glacier forefield of Midtre Lovénbreen comprises a large moraine system, with some parts potentially underlain by ice (Tonkin et al., 2015).The study area provides an ideal test site for our work since a helicopter is required to access other glaciers nearby and is typically available in late summer and early autumn.Thus we are able to hitchhike on these flights with our simple camera set-up, described in detail below.
4 Data and methods specific to the case study

Data acquisition
Two test surveys (15 and 3 September 2015) were conducted over the glacier and its forefield during routine helicopter flights to other nearby areas.Having surveys on 2 consecutive years allowed us to check the repeatability of the method as well as its ability to detect annual changes.The two surveys were conducted under very different weather conditions, which led to distinctly different glacier surfaces (Figs. 5 and 9): -In 2014, there is a thin layer of fresh snow on the ice and ground.There are long, sharp shadows due to clear skies and the low sun angle of late September at this high latitude.
-In 2015, there is no snow yet, so the images show bare ground and ice.There are no shadows due to overcast conditions on the day of the survey.
The different light conditions led to different camera settings.The surface brightness during the 2014 flight made the camera choose an average exposure time of 1/1500 s, while the darker conditions in 2015 led to an average exposure time of 1/500 s.It is questionable whether this was sufficiently fast to avoid motion blur, and thus the 2015 survey images may contain a slight rolling shutter distortion; that is, images are stretched in the y axis, corresponding to the flight direction.
The aircraft used was an Eurocopter AS350 helicopter.The camera is easily attached to the cargo swing with standard action-camera accessories (green marker in Fig. 6).We acquired data by making a detour on flights going to a neighbouring glacier.The 2015 acquisition was performed in parallel with an IceCam survey (Sect.4.4.2).The surveys took about 10 min each, flying at average altitudes and speeds of 1500 m and 35 m s −1 in 2014, and 1100 m and 50 m s −1 in 2015.
For our camera, we used a GoPro Hero 3+ Black Edition.At a relatively low price, it offers durable construction, decent image quality, and a variety of accessories for attaching to an aircraft's under-body, landing gear, carrying latch, or  other locations.In addition, it offers a wide field of view that, at the cost of fish-eye distortion, ensures good coverage of the whole area of interest.Fish-eye lenses are typically hard to calibrate simultaneously within the SfM process from inflight nadir images; however, one can perform a calibration separately (e.g. in the lab or under more suitable geometric ground conditions) and input it into the orientation step.In a case such as this study where the topography is very rugged and its variation is more than half of the flight height, simultaneous in-flight calibration (also called auto-calibration) is actually possible.It should also be noted that the calibration of a camera is temperature sensitive (especially lower-end models), which pushes the balance in favour of in-flight calibration.Figure 5 shows sample images from the surveys and Table 1 gives a statistic summary of the image parameters.
The camera was set in Protune mode and the images do not display saturation.For both surveys, the camera was set to take a picture every second.
The GNSS system was a hand-held consumer-grade hiking GPS device, the Garmin GPSmap 60CSx.After installing the camera on the helicopter and initializing image acquisitions, we acquired pictures of the Garmin GPS display showing the time for about 1 min.Then, prior to take-off, the GPS receiver was placed in the front of the nose of the helicopter, close to the lower front window (see red marker in Fig. 6).After the flight, we again acquired pictures of the GPS time display in order to calibrate the image acquisition times.

Processing
Once the images and GPS track were obtained, a visual check confirmed that the data were good and ready for processing.Because the helicopter flew higher than expected on both surveys, we were able to crop the image to keep only the 2500 × 2500 pixels around the centre of the GoPro images, where image quality is best and distortion is minimized, while nonetheless maintaining sufficient overlap between images.The first part of the SfM photogrammetric process (automatic tie-point identification, followed by autocalibration and computation of camera relative positions and orientations) was then performed.
Using the methods of Sect.2.3, a piecewise linear regression estimated the T 0−i values measured and the sevenparameter transformation residual minimization method resulted in estimated values for Delay GPS-Camera0 of 4.7 and 2.3 s in 2014 and 2015, respectively.
Calibration errors are common when performing autocalibration using nadir images and may lead to a high-order bias (or doming effect; e.g.James and Robson, 2014) affecting the precision of the DEM.Therefore, we refined the autocalibration through GNSS-guided bundle adjustment using the camera coordinates computed from Delay GPS-Camera0 and all the T 0−i to constrain the x, y, and z axes, with higher weights on the z axis since it is the most sensitive to calibration errors.The z coordinates also vary smoothly along the flight, and a strong temporal coherence in the GNSS data can therefore be expected for that component.The residuals for all images before and after that step are presented in Fig. 7. From this final orientation and calibration, DEMs and orthoimages with 2 m (2014) and 1 m (2015) ground-sampling distance (GSD) were extracted using the MicMac software.

DEMs and orthoimages
The above processing resulted in the creation of DEMs and orthoimages over the glacier and its moraine.However, light snow cover on the surface during the 2014 flight, especially in the upper glacier and combined with the bright sunshine and long shadows, create a high global radiometric dynamic; The Cryosphere, 11, 827-840, 2017 Figure 7. Residuals of the fit between the GNSS-extracted positions and the relative SfM-based orientation of the images, before and after the GNSS-guided bundle adjustment for the 2015 flight.There is a systematic error in the z axis coming from a misestimated focal length (each of the four "bumps" in the blue line representing the residuals in the z axis in the upper panel corresponds to a line of flight; Fig. 9).Righthand panels show the DEM differences to a reference DEM from 2010.The non-compensated DEM (upper-right panel) has a significantly stronger dome-shaped distortion pattern.After correction, the magnitude of the doming is less than the precision of the GPS.this lead to a limited local dynamic range in the images captured by the camera.An absence of coherent features made the correlation phase fail in the upper regions of the glacier (Fig. 8) for the 2014 survey.

Additional data sets
To validate our data and perform surface elevation change analysis, we gathered two additional data sets from different sources.Table 2 compiles the key information and Fig. 9 shows the orthoimages and flight paths of all data sets.

UltracamXP 2010
In 2010, the Norwegian Polar Institute (NPI) commissioned photogrammetric flights over the study area, as part of a Svalbard-wide photogrammetry campaign.The camera system used was a Microsoft Vexcel UltracamXP digital aerial mapping system (Microsoft, 2016) with a 60/20 overlap.This professional photogrammetric system is fully integrated and outputs pan-sharpened images at 0.25 m GSD, with associated GPS/INS information.We processed the images with MicMac to create a DEM and an orthoimage, using the inflight GPS information for georeferencing.

IceCam 2015
The IceCam system (Divine et al., 2016) was designed to retrieve small-scale sea ice topography.It was flown in September 2015 to assess its use for larger-scale land surveys.The hardware component of the system comprises two Canon 5D Mark II cameras, a combined Novatel GPS/INS unit, and a laser altimeter mounted in a single enclosure outside the helicopter.The altimeter has a limited range and could not be used for this survey.The images were also processed with MicMac in a similar way as the UltracamXP 2010 data set, creating a DEM and an orthoimage, using the in-flight GPS information for georeferencing.
5 Results and interpretations

DEM comparisons
First, all DEMs were compared to each other, focusing on areas with assumed stable terrain (i.e.off-glacier).Because of the low absolute accuracy of the low-cost code-based GPS tracker and the remaining errors in the modelled camera capture times, DEM differencing reveals patterns related to translations and rotations between the DEMs.We coregistered the DEMs using a modified version of the algorithm of Nuth and Kääb (2011), which incorporates rotations about the x and y axes, in addition to the 3-D translation: where dh is elevation difference, a and b describe the translation vector amplitude and direction respectively, c is the mean vertical difference, and d and e are the rotation parameters describing a linear plane rotated about the x and y axes.
Translations were typically on the order of a few metres.
To assess the quality of our products and evaluate the improvement that can be attributed to the camera time refinements, we compared two 2015 GoPro products with the 2015 IceCam DEM, which is made from images acquired on the same flight.The first product (Fig. 10, top) is computed without the camera time refinements presented in Sect.2.3 while the second (Fig. 10, bottom) includes the refinements.These differences reveal a standard deviation after co-registration of for the whole scene of 3.87 and 1.93 m, respectively, for the unrefined and refined method, after removal of obvious outliers.These values are however still affected by the remaining doming effect mostly seen at the edges of the survey area, most likely resulting from the imperfect camera calibration (Figs. 7,10,left panels).In the lower glacier area, where correlation is good in both data sets and where the doming effect only minimally affects the data, the standard deviation is 1.19 and 0.47 m, respectively, for the unrefined and refined method.It is quite obvious in the top left panel of Fig. 10 that the strains imposed by the inaccurate EXIF time stamps creates a distortion in the model: the DEM difference shows errors in Z coherent with errors in XY in different directions in different area of the model with a boundary along the centre of the scene (with negative errors in Z moving from southern slopes to northern slopes).
We then co-register and difference the 2014 and 2015 DEMs with each other and with the 2010 DEM.A number of features are visible in our difference maps but the two main geomorphological phenomena visible in the data are the thickness change of the glacier and thermokarst processes at the western lateral glacier moraine.Figure 12 shows profiles The Cryosphere, 11, 827-840, 2017 of the DEMs and DEM difference maps for both phenomena (see positions of the profiles on Fig. 11).

Glacier elevation changes
As part of the Midtre Lovénbreen mass balance program conducted by the NPI, in situ glacier stakes are surveyed and measured biannually.Differential dual-phase GNSS is used to measure elevations of the tops of stakes, which are frozen into the ice surface.Distances to the local ice surface are measured using a tape from the stake top to a pole lying on the uneven ice surface, oriented in two positions, to yield an average elevation of the ice surface in the immediate vicinity of the stake.The difference between two such measurements made in autumn represents the change in the ice elevation over a year, as well as an additional component due to the movement of the stake downslope.We estimate and remove the latter by using the horizontal annual velocity multiplied by the tangent of the slope at the stake, as determined from the 2010 DEM.
Typically, during glacier retreat, greatest ice loss is at the glacier tongue, and decreases upglacier.From 2010 to 2015, the tongue of Midtre Lovénbreen lost up to 14 m of ice thickness.Our 2014-2015 differences reveal significant annual changes of up to −3 m.Comparison of the DEM differences at the position of the stakes with the stake measurements provides a completely independent evaluation of the quality and reliability of our method.Figure 13 shows the position of the stakes and the observed bias in height difference.The bias between the 2010 reference data and the 2015 data that we produced (computed for stakes 1-11) is 0.08 ± 0.19 m yr −1 while the bias between the 2014 and 2015 flight (computed for stakes 1-4 only) is −0.61 ± 0.46 m yr −1 .The residual bias is a combination of the remaining doming of the DEMs from the imperfect camera calibration (±0.46 m) as well as the bias from the glaciological stake measurements, estimated to be 0.28 m yr −1 .The observed trend in bias with elevation in the 2010-2015 difference (Fig. 13) is coherent with the submergence/emergence of ice with about 0.2 m yr −1 on the tongue, 0 m yr −1 at the dynamic ELA, and −0.2 m yr −1 in the accumulation area.
Since the whole glacier is covered by our 2015 survey, we can estimate an average glacier-wide elevation change between 2010 and 2015.The mean elevation change over 5 years at the stake positions is −3.06 ± 0.46 m with our method and −3.51 ± 1.42 m from the stake measurements.The mean elevation change over the entire surface covered by the two DEMs is −2.85 ± 0.46 m.The similarity in pointbased and area-based mean elevation change indicates that the glacier stake measurements are good representatives of the changes occurring over the whole glaciers.

Pro-glacial area
To evaluate the performance of the DEMs produced and make geomorphological interpretations for the pro-glacial area of Midtre Lovénbreen, we focus on the following zones: the left (western) lateral ice-cored moraine, the proglacial dead-ice area, and the main glacier runoff channel (Fig. 14).
The most prominent elevation differences between the 2010 and 2015 DEMs in the proglacial area can be found along the left lateral moraine (Fig. 14b).This moraine is, like many similar ones on Svalbard, ice-cored (Etzelmüller, 2000) and thus potentially subject to thermokarst processes (Kääb and Haeberli, 2001).In fact, strongest elevation losses of up to 2 m yr −1 occurred at thermokarst lakes and ice fronts where substantial ice loss occurred over the 5-year period is observed.Even between the 2014 and the 2015 DEMs ice melt at the margin of thermokarst lakes and retreat of ice fronts can be detected, both on the left lateral ice-cored moraine and at other thermokarst lakes on the proglacial dead-ice area (Fig. 14b).
Also, other sections of the ice-cored moraine without clearly visible thermokarst depressions show surface lower- The Cryosphere, 11, 827-840, 2017 ing of up to 1 m yr −1 and more (Fig. 14b).This points to the fact that the debris cover over the ice deposits is thinner than the local permafrost active layer so that the annual zero-degree isotherm reaches the ice and is able to melt substantial amounts of it (Etzelmüller, 2000;Tonkin et al., 2015).The enhanced lowering rates 2010-2015 on the icecored moraine area give a clear spatial pattern, suggesting that ground-ice contents or debris-cover characteristics could well be mapped and investigated from elevation changes such as derived in this study (cf. Tonkin et al., 2015).For this purpose we also see no difference, regardless of whether the 2015 GoPro or the 2015 IceCam DEM is used.Both reveal similar results because the actual elevation changes clearly exceed the combined uncertainty of the respective DEM differences to the 2010 UltraCam DEM.
The maximum lowering rates seen over intact debris-cover on the ice-cored moraine exceed values found by Tonkin et al. (2015) for 2003-2014 on the lateral moraine/dead-ice area at neighbouring Austre Lovénbreen, hinting to higher ice content or thinner debris cover for the western lateral moraine of Midtre Lovénbreen.Our DEM area also covers the dead-ice/ice-cored moraine area investigated by Tonkin et al. (2015) for 2003-2014, and we find a similar pattern and similar lowering rates as found by these authors (not shown in Fig. 14).
We also assess changes over the (potential) dead-ice areas in front of Midtre Lovénbreen (Fig. 14c and d).Besides some thermokarst lakes growing in the glacier forefield, in which local lowering rates are up to 1 m yr −1 due to expanding lake shores, surface lowering rates elsewhere seem to be close to or within the uncertainty of the DEM differences and can be seen as residual zonal DEM shifts, long-wavelength vertical distortions, or image stripe patterns.We thus choose to only interpret local variations in elevation changes that should not be affected by the above errors with longer wavelengths.
There seem to be not only sections with lowering rates of up to around 0.1 m yr −1 over 2010-2015 but also areas that seem to change little.Most of the local patterns of elevation changes coincide with spatial variations in landforms.For instance, the borders of some (thermokarst?)depressions seem to erode and thus produce elevation loss.Also, where the glacier forefield has a rough topography with depressions and small hills typical for thermokarst areas, a spatially variable pattern of elevation changes is visible.In contrast, where the glacier forefield topography is more smooth, the elevation changes also display less spatial variations.In one area of the forefield there was surface icing present in both 2010 and 2015; DEM differencing showed that the ice accumulation is 1.5-2 m lower in 2015 compared to 2010 (Fig. 14c).In sum, there is a contrasting pattern of mostly stable or lowering surface elevations of up to about −0.1 m yr −1 over the proglacial area that seems to reflect patterns in groundice loss and erosive processes.A detailed quantitative analysis with the DEMs produced and over the time period investigated, however, is difficult as a pattern of higher-order DEM biases of the same order of magnitude overlays the www.the-cryosphere.net/11/827/2017/The Cryosphere, 11, 827-840, 2017 geomorphological signals related to most of the dead-ice processes present in the glacier forefield investigated.Etzelmüller (2000) also investigated elevation changes on icecored moraines and dead-ice areas, but not on Midtre Lovénbreen.Elevation change rates from 1970 to 1990 are fully in agreement with our above findings.As a third landform example we describe elevation changes 2010-2015 over the main glacier runoff stream (Fig. 14d).Over its entire reach in the forefield there are elevation losses exceeding those on the surrounding areas, of up to 6 m over the 5-year period, at places where the stream cut a new path through the moraine material.Lateral channel erosion of up to 5 m over the 5-year period can be seen at the margins of steep stream flanks (Fig. 14d).Also, where the stream area gets wider and more braided in the forefield we find elevation loss of about 0.5 m over the 5-year period, more than on the surrounding, roughly stable terrain, and a loss that is probably due to stream bed erosion during periods of high discharge.Outside of the glacier forefield, imagery shows some sediment fan accumulating with new deposits.Elevation gains at these locations are some tens of centimetres, with a spatially consistent pattern, although the change is not statistically significant with respect to the combined accuracy of our DEM differences.

Conclusions
We present a method to derive high-resolution DEMs and orthoimages using simple, inexpensive off-the-shelf cameras and code-based hand-held GNSS receivers mounted on airborne platforms.The core of the method is synchronizing the inaccurate photograph capture times with the GNSS receiver clock during the SfM processing.By doing so, we refrain from requiring accurate GCPs and can generate DEMs fully automatically without human intervention (i.e.no need to manually define GCPs).We illustrate this approach over an Arctic glacier, Midtre Lovénbreen on Svalbard, using opportunistic helicopter flights performed for other purposes.In our case study, we are easily able to detect glacier thickness changes over a 5-year period, but also over 1 year, although the latter results are more sensitive to the small biases resulting from the imperfection of the photogrammetric solution.
In comparison to the stake-derived in situ elevation differences, the higher-order DEM biases are most significant over 1 year but are reduced dramatically when differencing over many years.
Our results show the potential to derive accurate DEMs of a glacier from opportunistic flights, for example, as combined with in situ mass balance surveys of the glacier in our study site.In this manner, the acquisition of DEMs and subsequent DEM differences for deriving geodetic glacier mass balance is simpler and more easily acquired for those monitoring these glaciers.Our methods can help further the application of photogrammetry to non-experts needing highquality terrain products.As one application scenario, photogrammetric DEMs could be easily acquired as a byproduct whenever helicopter flights are needed to reach a study glacier or any other field sites.
The strongest limitation to detecting subtle elevation changes over our 5-year observational period are higherorder vertical and horizontal biases in the DEMs and thus in the DEM differences derived.These biases restrict the statistical significance of the elevation differences to a range of 0.5-1.0m, although local patterns of elevation differences seem significant to values as low as 0.1 m.
Thermokarst processes and general ice loss on the lateral ice-cored moraine investigated are evident in our results and can be valuable to map potential and degrading ground-ice content or potential thicknesses of its debris cover (Tonkin et al., 2015).Furthermore, other processes such as icing or erosion in and at the margins of stream channels can be quantified.Small elevation changes in other parts of the forefield show reasonable spatial patterns and change rates, but quantification of these changes is complicated by the aforementioned higher-order DEM biases.
Most of the high-order biases visible, such as the doming effect discussed by James and Robson (2014), can be blamed on the GoPro fish-eye lens: its calibration is challenging and the low-quality GNSS data were not sufficient to completely remove the observed biases.However, the GNSS data still helped significantly (Fig. 7).It must also be noted that even if the lens's field-of-view resulted in a very wide swath, it is impossible to use the whole image.This is because the viewing angle to an object on the side of an image is so different to that in the middle of an image that it is recognizable neither in the matching algorithm (which leads to the absence of tie point detection) nor in the correlation step (which fails in its area).Our study uses one of the worst types of camera system for photogrammetry, thereby exemplifying some of the worst-case accuracies that, if already good enough for some analysis, could be improved with better fitted cameras.We therefore recommend using a camera with a more conventional lens with less distortion (i.e. a non-fish-eye lens).Cameras from the Olympus Tough series, Panasonic FT series, or Nikon AW series are good examples of field-worthy cameras that can produce better images (for our purpose) than a GoPro.Even better, higher-end cameras (such as fullframe DSLR or hybrids like the Sony A7S) are of course also available, but their use is subject to the willingness of the owner/user to strap them more or less precariously to an aircraft.
Finally, even with an optically satisfactory camera, there is still the inherent problem of correlating images over the white, contrast-free surface of a landscape covered in fresh snow, as seen in our test survey in the upper glacier sections.This will to some extent always be a problem when using optical imagery, although solutions are being developed, e.g.increased colour depth or cameras without infrared filters (Nolan et al., 2015).
Figure 1.Standard SfM photogrammetry processing workflow.Green boxes are processing steps and grey boxes are products.

Figure 2 .
Figure 2. Mean residuals in x, y, and z and Cartesian distance D of differences between raw GPS positions and SfM-reconstructed relative positions, as a function of different start-time Delay GPS-Camera0 offsets between GPS and camera.Search area (a) from −3.5 to +3.5 s in 0.5 s increments and (b) from 1.5 to 3.0 s in 0.1 s increments.The best value is found at ≈ Delay GPS-0 = 2.3s.Note that the z coordinate is nearly independent of Delay GPS-Camera0 because the flight used here was mostly at a constant altitude.

Figure 3 .
Figure 3. Differences between expected time stamps of pictures for a 30 min period with a lapse time of 1 s and the EXIF time stamp recovered from the first 1800 images of a GoProH3+BE.Results are either for a relatively slow micro-SD card (64 GB SanDisk Ultra microSDXC Class 10 UHS-I) or a fast card (32 GB SanDisk Extrem Plus microSDHC UHS-I/U3).The fast card diverges by only 2 s after 30 min, while the slow card diverges by 134 s at a visibly irregular rate.

Figure 8 .
Figure 8. Correlation scores for 2014 (left) and 2015 (right) GoPro flights.White indicates perfect correlation, shades of grey imperfect correlation, and black no correlation.For both surveys, the accumulation area presents lower correlation scores (correlation even failed for the 2014 survey).

Figure 10 .
Figure 10.DEM differences between simultaneous IceCam and GoPro surveys of 2015: top row without refined camera times and bottom row with refinements.Residual doming is evident even in the refined GoPro survey, as well as areas where either the IceCam or the GoPro surveys failed to produce useful data because of a lack of texture in the images (fresh snow, for instance).Histograms of the differences, after removal of outliers (dH > 10 m), are shown to the right for the whole scene and for the glacier only.

Figure 11 .
Figure11.DEM difference map of the tongue of the glacier between 2010 and 2015.In red are the two lines from which profiles were extracted: "a" for the ice-cored moraine and "b" for the glacier terminus.

Figure 13 .
Figure 13.Top: position of mass balance stakes overlayed on the DEM difference for 2010-2015 (left) and 2014-2015 (right).Bottom: height variation at the stake positions according to the DEMs (x axis) versus according to the stake measurements (y axis).

Figure 14 .
Figure 14.(a) Hillshade of the glacier forefield from the 2015 GoPro-derived DEM.Rectangles (b-d) indicate the location of respective panels of the lower row of the figure.(b) Elevation changes 2010-2015 over a section of the lateral ice-cored moraine with the position of thermokarst lakes indicated.(c) Elevation changes 2010-2015 over a section of the central forefield with icing.(d) Elevation changes 2010-2015 over proglacial channels; glacier terminus to the lower left; elevation changes underlain by the hillshade to better indicate channel positions.

Table 1 .
Statistic summary of the image parameters.

Table 2 .
Characteristic of the products.