SNOw-COvERED SURFACE vARIAbILITY AND DEM GENERATION USING AERIAL PHOTOGRAMMETRY IN MOUNT ODIN, CANADA

. Seasonal snow-covered surface has a critical role in global water resource supplement especially provi-ding fresh water for humankind and flora’s consumptions as well as local underground water storages. The in situ measurements of seasonal snow-covered variability are extensively prodigal and costly particularly in existence of severe climate conditions such as high latitude regions and polar areas. It is therefore necessary to apply remote sensing techniques and observations to estimate accurately the snowpack melting and accumulation for different seasons. In this paper, we estimate snow-covered surface variability for four different seasons of year in Mount Odin, Canada using aerial photos. In order to do this, firstly Digital Elevation Model (DEM) with respect to Earth Gravitational Model 1996 (EGM96) for each flight mission of A, B, C and D from these aerial photos by applying Bundle Adjustment (BA) triangulation is being generated precisely. Moreover, the displacement of each two DEMs is computing in order to determine snow-covered surface variability between each two flight missions. The results demonstrate that flight mission C has the highest elevation topographically compare to the missions A, B and D while mission C was planned in February 2011 in existence of vast snow throughout Mount Odin area as well as mission C’s DEM which has higher elevation values than the others. The proposed methodology and problem solution and the case study information with the details of each flight mission are discussed in expatiation. Reference to this paper should be made as follows: Najibi, N.; Arabsheibani, R. 2013. Snow-covered surface variability and DEM generation using aerial photogrammetry in Mount Odin, Canada, Geodesy and Cartography 39(3): 113–120.


Introduction
Seasonal snow is a vital resource to supply aquatic environment life like agriculture farms as well as underground water storages, continental water and water sources for different industries and factories, hydroelectricity generation and drinking urban water. The estimation of seasonal snow-covered surface variability in terms of the snow melting volume and accumulation plays a critical role in different levels of management especially natural hazard disaster warning system, management and planning when the huge amount of snowpack is going to be melted (Lee et al. 2008). Therefore, it is necessary to estimate precisely the snow-covered surface variations spatially and temporally.
Additionally, the variability in seasonal snow-covered surface as well as ice and permafrost can affect on the climate conditions locally and then globally especially the snow and ice surfaces which have higher albedo rather than the other natural land surfaces and thus are capable of reflecting a greater percentage of the incident solar radiation in spite of the decreasing coverage of snow and ice will result in an increase in the amount of solar radiation absorbed by the earth.
Since from last three decades there have been difficulties and unsuccessful missions for determination desirable accuracy in height in snow-covered terrain image matching and Digital Elevation Model (DEM) creations domain (Al-Rousan et al. 1997;Höhle 2009); nowadays the current enhancements applied to the digital camera systems as well as considering high accurate flight planning encouraging developments for high accuracy measurements and estimations (Wolf, Dewitt 2000).
Many studies have been done to determine snowcovered surface variability by using different techniques; for instance Bacher et al. (1999) used aerial photography to generate DEMs of glacial regions; they concluded their project with accuracy up to 2 m. Moreover, Ledwith and Lunden (2001) estimated DEMs of snow-covered areas in Norway through using aerial photos at a scale of 1:30000 and reached to an average difference of 2.8 m between derived DEMs and Global Positioning System (GPS) field measurements. Though Sun et al. (2003) used Shuttle Radar Topography Mission (SRTM) products including surface height to compare with the surface height from Shuttle Laser Altimeter (SLA) but they considered SRTM and airborne mission in order to minimize the random errors of creating snow-covered Siberia landscape through comparisons between the measurements from SRTM and Light Detection And Ranging (LiDAR) which provided them with the height error caused by the vegetation and snow-covered mountains; since firstly they assessed SLA data by field observations and statistical examinations of functional waveforms. Besides, assimilation effects from the Moderate Resolution Imaging Spectroradiometer (MODIS) of snow cover fraction products as well as snow cover fraction and snow depth estimation from the Air Force Weather Agency (AFWA)-National Aeronautics and Space Administration (NASA) Snow Algorithm (ANSA) product are analyzed and discussed in details by Liu et al. (2013) for Alaska region.
Furthermore, redistribution analysis of snow transport and integration measurements of snow sublimation and snow-covered surface melt rates have been done recently by Bernhardt et al. (2012) in Berchtesgaden National Park of Germany as the case study area. They used several complex models (e.g. SnowModel, SnowTran-3D, SnowSlide) in order to compare analytically the model results and the relevant sorted satellite images to validate the general proposed snow model which followed by this conclusion that the snow transport processes inclusion into the model environment is improving the accuracy of the predicted snow distribution.
To monitor snow cover variability, Paudel and Andersen (2011) firstly analyzed the spatial and temporal variability of snow cover in the Trans Himalayan region of Nepal during last decade and then proposed five practical steps in order to remove cloud obscuration from MODIS daily snow products including combining data from the Terra and Aqua satellites, adjacent temporal deduction, spatial filtering based on orthogonal neighbouring pixels, spatial filtering based on a zonal snowline approach and finally temporal filtering based on zonal snow cycle.
Although nowadays the DEM generation derived by means of at least two images which its method and results are described in details in many studies, for example (Gülch 1994) and (Heipke 1995); the general steps are as the following: flight planning, aerial photos acquisition, georeferencing, photos matching, applying necessary parameters and triangulation solution, DEM generation, necessary editing and finally the quality control. The overall steps used during this study in order to measure snow-covered surface variability using aerial photogrammetry technique are being shown in Fig. 1.
In this study, we concentrate on Mount Odin located at Northern Canada by generating 4 DEMs of snow-covered terrain where there is a vast and sizable snow on the ground. The aerial photos captured during 4 photogrammetry missions defined here as mission A, B, C and D and the related photos have been taken in August 17, 2010, November 10, 2010, February 19, 2011 and May 10, 2011 flight missions respectively. This paper represents firstly the proposed methodology in Section 1. In Section 2, the data and observations are placed. The results and discussions are discussed in Section 3 which eventually followed by conclusion remarks.

DEM generation using Bundle adjustment of aerial photogrammetry
In method of Bundle adjustment (BA) simultaneously the three dimension position (X o , Y o , Z o ) and rotations (ω, φ, κ) of all the block's bundles are adjusting to get the desired possible intersection of conjugate rays of each point on the ground (X, Y, Z) and fit to the given control points as the following (Wolf, Dewitt 2000): (1) where c is the calibrated focal length; x, y are the coordinates of each point in the certain photo, λ p is the scale factor and R is the well-known rotation matrix consists of ω, φ and κ. It is assumed that each light ray is to be a straight line passing through entire of the projection center, image point and terrain point straightforwardly. Therefore, the observations equations are given as: , (2) where x ij and y ij are the photo's coordinates of each point in the photo of j th , v is the possible random error residual and c is the calibrated focal length, the coordinates of projection centre are , , , the elements of point i th are defined as , … and the terrain coordinates of point i th considered as X i , Y i , Z i .
In order to solve Equation (2), there are many methods; the convenient method is to linearize it according to Taylor expansion (Wolf, Dewitt 2000); and thus the derivatives of Equation (2) will be concluded with respect to BA's orientation parameters (X 0 , Y o , Z o , ω, φ and κ) and three dimension coordinates of the considered tie points.
Obviously, this procedure is iterative due to applying the linearization to the observation equations, but we applied Powell's dog leg non-linear least squares technique in BA solution instead of tradition method of Levenberg-Marquardt optimization algorithm in order to avoid storing and operating the zero entries. Practically the number of iterations depends critically on how good we know beforehand the approximate value of the unknowns (the approximate value of each X, Y, Z, X 0 , Y 0 , Z 0 , ω, φ and κ).
The general process of BA is illustrated in Fig. 2.

From DEM to snow-covered surface variability estimation
The proposed methodology to determine snow-covered surface variability requires at least two DEMs with time span differentially. These two DEMs are provided in two different time spans which it is decided to be measured the changes upon targeted snow-covered surface. Clearly, the DEM difference between each two of these DEMs can provide us with snow-covered surface variability between this two time span in terms of average snow volume estimation and the seasonal snow covered region boundary. In this study, using BA approach leads us to generate the DEMs in 4 different time spans in one year. The related DEMs have been generated based on EGM96 geopotential reference frame.

Case study, data and observations
The aerial photogrammetry has been planned in 4 flight missions as mission A, B, C and D. These aerial photos have been captured by digital mapping camera (DMC) with focal length 120 mm over Mount Odin in Baffin Island, Nunavut, Canada (Coordinates in WGS84: 63.54611°N and 65.43028°W). Nowadays Mount Odin is the highest mountain within the Baffin mountains and also one of the highest mountain in the Arctic Cordillera (Baffin mountains are a mountain range located at the northeastern coast of Baffin Island and Bylot Island as a part of the Arctic Cordillera) which holding a topographic prominence of 2147 m and thus Odin is the third highest mountain in Nunavut, Northern Canada by the mentioned topographic prominence. Fig. 3 illustrates the location of Mount Odin.
In Table 1, the number of captured aerial photos in each flight mission, numbers of tie and control points, mission's date and mean flight elevation have been represented respectively.
This area is holding a dynamic behaviour with respect to the seasonal snow-covered surface during different months of year.
Additionally, Fig. 4 (a, b, c and d) represents the in situ captured aerial photos as a schematic for separated flight line during each mission respectively.  The minimum 60 percent aerial photos overlapping for each mission is been followed completely (Gogineni 2012). The number excess than 60 percent will lead to extra cost for this project as the less overlapping value can cause uncertainty and additional problems in image matching as well as the stability of the proposed mathematical solution.
However, natural brightness conditions and illumination in the moments of capturing photos or handling them should be followed as an identical situation for each pair images during BA solution as much as it could be feasible. Since Toutinm (2004) demonstrated that the extracted DEM from Ikonos images based on terrain slope and azimuth resulted in error reduction as linear evolution with respect to the slope; moreover those elevations in sun-facing slopes conditions were more precise about 1 m rather than those elevations in against direction of the sun approximately.
Although the direct sunlight due to causing very high differences in the light intensity especially creating hard shadows may practically contribute to very large systematic errors through confining the darker areas in the stereo images (Heinz 2002;Orun, Alkis 1996), but some large noises can be eliminated by using high quality photography lenses as well as by covering the excavation area with a width roof inside of an absolute dark room.
In this study we assume that these errors sources are in minimum values and cannot cause a significant deviation in BA solution; although the all flight missions were taken place in the similar day timing conditions.

Flight runs mission simulation and overlapping
As it is illustrated in Table 1  According to Fig. 5, flight run for mission C was planning in higher elevation especially when it is getting close to Mount Odin area (see Table 1). It is clear that the corresponding aerial photos for this mission will be in larger scale as well.
The reason that additional aerial photos have been captured for mission C is describing as the following; since the exact boundary of snow-covered surface in Mount Odin was not known precisely, extra photos Similarly to the other aerial photogrammetry missions, the captured photo sensor board (aircraft) might be affected inevitably by its sudden movements due to the wind drift and small movements of the aircraft which led into further undesirable errors; Fig. 6 presents this inhomogeneity during aircraft flight time period. However the idealistic situation is happening when the both values of pitch (rotation about the Y-axis or North direction) and yaw (heading) (rotation about the X-axis or East direction) are equal to zero approximately, but if they are exactly zero, it is a so-called nadir photo (Moffit, Mikhail 1980). In fact during practical photogrammetry missions, this will never occur due to the mentioned unavoidable obstacles during the mission execution. Consequently, the result of this irregularity can contribute to a challenging procedure in DEM creations and then DEM comparisons to understand the snow-covered surface variability.

DEM generation with respect to EGM96
DEM generation is performed using BA approach for mission A, B, C and D with respect to EGM96. The BA procedures to generate DEM are discussed completely in Section 2. In this study, we transformed the DEM relative values to EGM96 due to simplicity and reality. The generated DEM corresponding to each mission A, B, C and D are given in Fig. 7 respectively.
Since mission A was performed in August 2010; in this time of the year there were not a dense and deep snow on the ground in Mount Odin region as it can be understood from aerial photos for mission A which is illustrated in Fig. 4 as well. Mission A's DEM in Fig. 7 clearly depicts an average elevation for Mount Odin topography which can be concluded as a region without snow accumulation or with a little snow accumulation placed in this region caused mainly by snow melting during the spring for Baffin Island.
There were snow falling in November 2010 in Mount Odin which resulted in snow deposition over this area as flight mission B is shown in Fig. 4 and mission B's DEM in Fig. 7.
Moreover, snow falling is reaching to its maximum rate in February in Mount Odin as the snow melting is in its minimum one. Since mission C is done in February 2011 (see also Fig. 4); it is decided to take more covered ground surface due to importance of this time span in the year and lack of field information about the exact boundary between bare land and snow surface. Additionally, the mission C's DEM demonstrates a high elevation topography level as there was a sizable snowpack on throughout Mount Odin (Fig. 7).
Mission D is planned in May 2011. Commonly this time span is a commencement to snow melting while the spring is started. The aerial photos captured during mission D (Fig. 4) is presenting that these less deep snow-covered layers are melted.
Moreover, it is given in Fig. 7 as Mission D's corresponding DEM that demonstrates a sudden change in the topography elevation levels in Mount Odin especially compare to mission C's DEM.

DEM comparisons and snow-covered surface variability monitoring
In order to determine the snow-covered surface variability caused by snow falling and accumulation or snow melting over Mount Odin, the comparison of generated DEMs is proposed. It is assumed that with a good approximation each two DEMs belong to the same ground footprint of Mount Odin.
The comparison of each two DEMs will result in the average displacement of topography which means the snow-covered variations during this time span of two DEMs generation by using the related aerial photos. Obviously, the positive average difference (D) Mission D value presents the increase in ground surface elevation which mostly caused by snow accumulation. In the other words, the snow melting will lead to a sudden decrease in surface elevation of the region which can demonstrate by average negative value in difference of two DEMs; for example the difference of mission A's DEM from C's DEM gives us a negative value in elevation direction which means the snow-covered surface were diminished in mission A's calendar time (August 2010) rather than mission C which was in February when there were a vast and deep snow-covered over the ground in Mount Odin. According to these comparisons, Table 2 illustrates the average snow-covered surface variability during one calendar year for these 4 missions of A, B, C and D respectively.
Since it is presented in Table 2, mission C's DEM has positive values completely which demonstrates its highest elevation compare to the other DEMs and thus snow-covered surface is in its maximum elevation during February 2011. For instance, the mean difference between mission C's DEM and D's one is 1.4 m; this says that there was 1.4 m snow in mission C (February 2011) more than mission D (May 2011) averagely.
However, we should consider the general weather condition especially the average air temperature average during the flight missions. The nearest meteorological center to Mount Odin is named QIKIQTAR-jUAQ CLIMATE. As Fig. 8 presents, during the whole four flight missions the air temperature as well as snow depth had their seasonal trends and there was no abnormal event like sudden increase in air temperature or even unusual decrease in air temperature time series which might be resulted in the unusual changes on snow-covered surface and thus led to incorrect DEM generating and eventually followed by unreasonable conclusions based on it.

Conclusions
This paper presents a practical study on snow-covered surface variability measurement using remote sensing technique especially DEM generation from aerial photogrammetry over snowy ground. As it is being discussed, different DEMs generated for different time point during a calendar year can result in a reasonable estimation for snow-covered variability for each time span. The accuracy of generated DEM here is in order of 0.1 m which compare to the snow depth over Mount Odin may be a precise estimation. Although satellite remote sensing mission's products as well as a planned field measurement can provide us with a high advantage in terms of accuracy temporally and spatially but due to special condition of Mount Odin the aerial photogrammetry has been preferred.
The results represented that mission C has the maximum elevation topographically rather than the other missions; additionally the air temperature and snow depth values derived from the nearest meteorological center to Mount Odin (see Fig. 8) are in consistent to this conclusion. Moreover, the generated DEM for mission C has higher elevation values as well as it is illustrated in the mission's C aerial photos which there are a vast snowpack throughout Mount Odin. However the proposed methodology has been done successfully but it is necessary to apply and use additional missions and different remote sensing sensors and techniques especially GPS-Reflectometry technique (Najibi, jin 2013) and Synthetic Aperture Radar (Interferometric) (SAR or InSAR) over different regions in order to better understand snow-covered surface variability during different certain time span. It is important to  consider the derived accuracy, total costs and mission applicability of each methodology in the whole mentioned approaches in order to make a good decision upon snow-covered surface variability measurements for the water related issues management.