Articles | Volume 14, issue 11
Research article
12 Nov 2020
Research article |  | 12 Nov 2020

Parameterizing anisotropic reflectance of snow surfaces from airborne digital camera observations in Antarctica

Tim Carlsen, Gerit Birnbaum, André Ehrlich, Veit Helm, Evelyn Jäkel, Michael Schäfer, and Manfred Wendisch

The surface reflection of solar radiation comprises an important boundary condition for solar radiative transfer simulations. In polar regions above snow surfaces, the surface reflection is particularly anisotropic due to low Sun elevations and the highly anisotropic scattering phase function of the snow crystals. The characterization of this surface reflection anisotropy is essential for satellite remote sensing over both the Arctic and Antarctica. To quantify the angular snow reflection properties, the hemispherical-directional reflectance factor (HDRF) of snow surfaces was derived from airborne measurements in Antarctica during austral summer in 2013/14. For this purpose, a digital 180 fish-eye camera (green channel, 490–585 nm wavelength band) was used. The HDRF was measured for different surface roughness conditions, optical-equivalent snow grain sizes, and solar zenith angles. The airborne observations covered an area of around 1000 km× 1000 km in the vicinity of Kohnen Station (750 S, 04 E) at the outer part of the East Antarctic Plateau. The observations include regions with higher (coastal areas) and lower (inner Antarctica) precipitation amounts and frequencies. The digital camera provided upward, angular-dependent radiance measurements from the lower hemisphere. The comparison of the measured HDRF derived for smooth and rough snow surfaces (sastrugi) showed significant differences, which are superimposed on the diurnal cycle. By inverting a semi-empirical kernel-driven bidirectional reflectance distribution function (BRDF) model, the measured HDRF of snow surfaces was parameterized as a function of solar zenith angle, surface roughness, and optical-equivalent snow grain size. This allows a direct comparison of the HDRF measurements with the BRDF derived from the Moderate Resolution Imaging Spectroradiometer (MODIS) satellite product MCD43. For the analyzed cases, MODIS observations (545–565 nm wavelength band) generally underestimated the anisotropy of the surface reflection. The largest deviations were found for the volumetric model weight fvol (average underestimation by a factor of 10). These deviations are likely linked to short-term changes in snow properties.

1 Introduction

Surface reflection in the polar regions plays an important role in the Earth's climate system. The snow surface albedo is highly variable on both a temporal and spatial scale. This causes uncertainties in determining the surface solar radiative energy budget in these areas (e.g., Kuipers Munneke et al.2008). As such it is essential to monitor the reflective properties of snow surfaces to accurately predict future climate change in the polar regions. However, their remote and harsh environment makes them difficult to access and requires remote sensing techniques to observe snow surface reflection and its influencing parameters.

The high spatial and temporal variability of the snow surface albedo necessitates continuous observations with global coverage, which are provided by satellite instruments. Polar-orbiting satellites monitor the reflectance (i.e., reflected radiance in units of Wm-2sr-1) at the top of the atmosphere (TOA). However, they are restricted in terms of the number of available observation angles and spectral bands as well as their temporal resolution. The processing of using measurements from polar-orbiting satellites to monitor the broadband surface albedo typically requires three steps (e.g., Qu et al.2015): atmospheric correction, modeling of the angular reflectance, and narrow-to-broadband conversion. During the first step, the TOA reflectance is converted into a surface reflectance by correcting for gaseous and aerosol scattering and absorption applying radiative transfer modeling (e.g., Vermote and Kotchenova2008). For the calculation of the narrowband surface albedo, accurate knowledge of the bidirectional reflectance distribution function (BRDF) of the surface is required due to the limited number of available observation angles. The BRDF describes the scattering of electromagnetic (EM) radiation from one incident direction to another direction of the hemisphere. By using a linear combination of the discrete spectral band measurements with specific weighting, the broadband surface albedo is calculated (e.g., Brest and Goward1987; Klein and Stroeve2002; Liang et al.2003; Pohl et al.2020). The largest uncertainty in this three-step process is introduced by the angular dependence of the surface BRDF, especially when the reflection of the underlying surface is highly anisotropic. This is the case for the reflection of solar radiation at snow surfaces in polar regions due to low Sun elevations and the anisotropic scattering phase function of snow crystals. As the snow surface of the Antarctic ice sheet is also used for the validation and cross-calibration of polar-orbiting satellites (Jaross and Warner2008), a thorough understanding of the anisotropic reflection is needed.

However, as an infinitesimal quantity, the BRDF of a surface cannot directly be measured under natural illumination conditions. From a rigorous physical point of view, most satellite, airborne, and ground-based instruments measure the hemispherical-directional reflectance factor (HDRF; Schaepman-Strub et al.2006) if the reflectance is constant over the full cone angle of the instrument's field of view (FOV). In contrast to the BRDF, the HDRF includes incident irradiance from the entire hemisphere. An effective BRDF can only be derived from HDRF observations if the atmospheric influence is considered and the FOV is sufficiently small (e.g., Gatebe et al.2003; Lyapustin et al.2010).

For example, the retrieval of the aerosol optical depth (AOD) over highly reflecting surfaces is complicated, as it is hard to separate the surface and aerosol contributions to the measured reflectance at the TOA (Tomasi et al.2015), especially for the comparably low values of AOD in polar regions (Mei et al.2013). Due to its high surface albedo, most of the solar radiation incident on a snow surface is reflected, which depends strongly on the BRDF of the snow surface and, thus, on the illumination and viewing geometry. In general, most of the photons are scattered by the snow grains in the forward scattering direction. Yang et al. (2014) compared Moderate Resolution Imaging Spectroradiometer (MODIS) measurements of AOD at 550 nm over eastern China with ground-based data from the Aerosol Robotic Network (AERONET). Using a Lambertian model, 54 % of the regression points fell within the expected error envelope around the 1:1 line. This result is improved to 69 % if a more sophisticated, non-Lambertian model accounting for anisotropic reflection at the land surface is applied.

Yet the BRDF varies with snow grain morphology, the solar zenith and azimuth angles, the liquid water content of the snowpack, the snow impurity type and concentration, the dimension and orientation of surface roughness structures, and the wavelength. To account for these effects, an accurate model of the snow BRDF is needed. A BRDF model can be physical, empirical, or semi-empirical. Physical BRDF models (e.g., Cook and Torrance1982) accurately simulate the scattering of an EM wave at a surface by applying physical laws. The high accuracy comes at the cost of very complex computations. Empirical BRDF models (e.g., Phong1975; Walthall et al.1985) mimic the surface reflection by means of a simple, non-physical formulation. However, the drawback of the rather simple computations is their restricted accuracy. Semi-empirical models (e.g., Martonchik et al.1998; Lucht et al.2000) use simple, direct parameterizations of a more complex physical BRDF with a limited number of independent parameters.

The anisotropic reflection by snow surfaces was investigated by Li (1982) with simulations of the snow BRDF using Mie theory and the doubling and adding method. Applying Mie theory and the discrete-ordinate method, Han (1996) retrieved the surface albedo from satellite measurements in the Arctic. Leroux and Fily (1998) developed a snow BRDF model including the effect of sastrugi by means of regularly spaced, identical, and rectangular protrusions. Leroux et al. (1998) and Leroux et al. (1999) employed the doubling and adding method, Mie theory, and ray tracing to develop a snow BRDF model including polarization. Comparing the simulated values with observations in the principal plane, they found that the BRDF in the near-infrared wavelength range is strongly affected by the snow grain shape, whereby simulations assuming hexagonal particle shapes yield an improved agreement with the observations compared to assuming spheres. Accordingly, Aoki et al. (2000) stressed the importance of the particle shape assumption for the scattering phase function used in their BRDF model (applying both Mie theory and the Henyey–Greenstein phase function) when comparing to observations between 0.52 and 2.21 µm wavelength. Kokhanovsky and Zege (2004) demonstrated the use of an asymptotic analytical equation to model the BRDF of snow. Their approach represents the snow grains for example as fractal particles and, thus, accounts for their non-sphericity. Evaluating this asymptotic model with in situ measurements of the snow reflectance, Kokhanovsky et al. (2005) found generally good agreement but reduced model accuracy in the solar principal plane at large observation angles. In contrast to considering snow grains as independent scatterers of fractal shape, Malinka (2014) provided a framework to calculate the inherent optical properties (extinction coefficient, single-scattering albedo, scattering phase function, and the polarization properties) by applying the concept of stereology, which considers snow as a random mixture of ice and air. Malinka et al. (2016) combined the stereological approach with analytical, asymptotic formulas to calculate the bidirectional reflectance. Their model showed high accuracy when compared to albedo observations of snow-covered sea ice.

The comparison of in situ measurements of the snow reflectance with simulations is essential in terms of model validation. Observations of the HDRF (or effective BRDF in case the FOV is small and the atmospheric influence is considered) are conducted using a variety of different measurement concepts. For ground-based applications, manual or automated gonio-spectrometer systems are employed (e.g., Painter et al.2003; Pegrum et al.2006; Bourgeois et al.2006a). Kuhn (1985) observed a peak in reflectance that is contained within ± 60 to both sides of the solar azimuth angle. This peak becomes more prominent with increasing solar zenith angle and snow grain size. Marks et al. (2015) measured the HDRF of snow surfaces between 400 and 1600 nm wavelength with the Gonio Radiometric Spectrometer System (GRASS; Pegrum et al.2006) at Dome C, Antarctica. Their observations also showed enhanced forward scattering. In addition, they observed a larger anisotropy of the surface reflection at longer wavelengths. Employing a gonio-spectrometer, Bourgeois et al. (2006b) measured strong variations of the HDRF between 0.6 and 13 (in the wavelength range 350–1050 nm) depending on the solar zenith angle and the surface roughness at Summit, Greenland. Measurements with the Automated Spectro-Goniometer (ASG; Painter et al.2003) showed a decrease of snow HDRF at all wavelengths between 400 to 2500 nm when the snow grain size increased from 80 to 280 µm (Painter and Dozier2004b). Further comparisons with the results of a forward discrete-ordinate radiative transfer model (Stamnes et al.1988) revealed larger deviations between the simulations and observations for more complex crystals. The importance of the crystal habit for the anisotropy of reflection at snow surfaces was also emphasized by Dumont et al. (2010) and Stanton et al. (2016). The latter measured increasing anisotropy of the surface reflection during the growth of surface hoar in the laboratory. Several studies observed systematically less anisotropy for a typical snow HDRF than estimated from simulations (Warren et al.1998; Painter and Dozier2004a; Hudson et al.2006; Hudson and Warren2007). In the solar principal plane, the models mainly overestimate the forward scattering and underestimate the backward scattering. Implementing non-spherical grains in the BRDF models (e.g., Mishchenko et al.1999; Xie et al.2006; Jin et al.2008) improves the comparison with observations. The non-spherical model of Jin et al. (2008) agrees within ± 10 % for viewing zenith angles less than 60 with observations in Antarctica performed by Hudson et al. (2006). However, the asymmetry between forward and backward scattering still remains. This highlights the need to further incorporate macroscopic effects such as the roughness of the snow surfaces into the models (Leroux and Fily1998; Hudson and Warren2007).

Ground-based instruments observe the directional reflectance of a characteristic, homogeneous surface, whereas airborne and satellite observations average over a larger measurement area. Thus, the latter are more suitable for studying the influence of macroscopic surface roughness on the surface HDRF. Nolin et al. (2002) employed multi-angular measurements with the Multi-angle Imaging Spectroradiometer (MISR) for the characterization of surface roughness over Greenland and Antarctica. Measurements with the Clouds and the Earth's Radiant Energy System (CERES) showed monthly regional biases between 12 and 7.5 W m−2 in the cloudless TOA solar irradiance. These biases were attributed to the effect of sastrugi, which introduce a significant solar azimuth dependence: Kuchiki et al. (2011) observed a diurnal cycle in MODIS reflectances over the South Pole. In general, sastrugi decrease the forward scattering by casting shadows and increase the backward scattering due to a lower effective incident angle caused by the sastrugi slope, and the snow HDRF loses its azimuthal symmetry (Warren et al.1998). Zhuravleva and Kokhanovsky (2011) simulated a larger effect for a higher density of the sastrugi field.

Airborne measurements can be employed for the validation of BRDF models and the comparison with the large pixel size of satellite observations. Gatebe and King (2016) provided an extensive database of airborne spectral (effective) BRDFs for various surface types, for example, ocean, vegetation, snow, desert, and clouds. The effective BRDFs were acquired by the Cloud Absorption Radiometer (CAR; Gatebe et al.2003) over a 30-year period between 1984 and 2014. The CAR is a scanning radiometer covering 14 spectral channels between 340 and 2324 nm. The effect of surface roughness on the effective BRDF was studied by Lyapustin et al. (2010) with CAR measurements during the Arctic Research of the Composition of the Troposphere from Aircraft and Satellites (ARCTAS) spring campaign in April 2008 (Jacob et al.2010; Matsui et al.2011). Their results showed an agreement between the kernel-based Ross-Thick–Li-Sparse-Reciprocal (RTLSR) BRDF model (Lucht et al.2000) used in the operational MODIS BRDF/Albedo Model Parameters product (MCD43A1) (Schaaf et al.2002) and the CAR measurements to within ±0.05. The RTLSR model performed better in the forward scattering direction, whereas the kernel-based Modified Rahman–Pinty–Verstraete (MRPV) BRDF model (Martonchik et al.1998) used for the processing of data from MISR performed better in the backscattering direction. Jiao et al. (2019) proposed an additional snow kernel within the kernel-driven BRDF model framework to better account for the anisotropic reflection of pure snow surfaces. Incorporating the additional snow kernel yielded a correlation coefficient of 0.9 (compared to 0.65 for the original three-kernel model) and only small biases between the model and different BRDF validation data from ground-based and satellite observations.

Cox and Munk (1954) analyzed radiance-calibrated analog photographs for the parameterization of the ocean HDRF. Nowadays, digital cameras are increasingly applied in vegetation and soil monitoring (e.g., Lebourgeois et al.2008; Koukal and Atzberger2012). The instantaneous measurement of multiple viewing angles facilitates aerial HDRF measurements with digital cameras. The high angular resolution allows the detection of features in the anisotropy of the reflection which could be missed with HDRF measurements allowing only a limited number of observation angles. This was demonstrated by Goyens et al. (2018) for ground-based measurements of the snow HDRF, comparing radiance measurements with a spectroradiometer at discrete viewing angles with hyperangular observations using a fish-eye radiance camera. Ehrlich et al. (2012) used a commercial digital camera equipped with a wide-angle lens with a FOV of 100 for measurements of the HDRF of snow-covered sea ice, ocean, and clouds (viewing zenith angle up to around 60).

This study presents a method to derive the snow HDRF at a visible wavelength (spectral band 490–585 nm) from airborne digital camera measurements in Antarctica during austral summer in 2013/14. Concurrent measurements of the optical-equivalent snow grain size retrieved from spectral surface albedo measurements, and surface roughness determined by means of a laser scanner, allow for investigating their effect on the snow HDRF in separate case studies. Subsequently, the measurements are used to parameterize the snow HDRF applying the RTLSR model (Lucht et al.2000). The presented methodological approach allows for a direct comparison of the airborne HDRF measurements with satellite observations from MODIS.

The definition of reflectance quantities used within this work, the modeling of the BRDF, and the inversion of a semi-empirical kernel-driven BRDF model are presented in Sect. 2. The fieldwork and instrumentation are presented in Sect. 3 together with the detailed calibration involved in the measurement of the snow HDRF. The dependence of the snow HDRF on the solar zenith angle, the surface roughness, and the optical-equivalent snow grain size is studied in Sect. 4 based on three cases. Subsequently, the airborne measurements are compared with satellite-derived BRDF from MODIS. In the concluding Sect. 5, the findings of this work are summarized and perspectives for future studies are given.

2 Methodology

2.1 Definition of reflectance quantities

The definition of the reflectance quantities applied within this work follows Schaepman-Strub et al. (2006). The intrinsic reflectance properties of a surface are given by its BRDF. It quantifies the reflection of the incident radiation at the surface and its scattering from one direction of the hemisphere to another. The spectral BRDF (fBRDF, units of sr−1) provides for each zenith (θi) and azimuth angle (φi) of incident direct irradiance Fi(θi,φi;λ) the reflected radiance Ir(θi,φi;θr,φr;λ) for all directions of reflection (defined by the reflection zenith and azimuth angles θr and φr) by

(1) f BRDF = d I r θ i , φ i ; θ r , φ r ; λ d F i θ i , φ i ; λ .

The bidirectional reflectance factor (BRF, dimensionless) is obtained when the BRDF of a sample surface is divided by the BRDF of a diffuse (Lambertian) standard surface illuminated under the same conditions (identical beam geometry). As the BRDF of an ideal Lambertian surface is (π sr)−1, the BRF is given by

(2) R BRF = d I r θ i , φ i ; θ r , φ r d F i θ i , φ i d F i θ i , φ i d I r ideal θ i , φ i = π sr f BRDF .

In Eq. (2) and in the remainder of this section, the spectral dependence is omitted for reasons of simplicity. In atmospheric conditions, both the BRDF and BRF cannot be measured directly as the global irradiance F0 reaching the surface is composed of a direct (Fi) and diffuse (Fdiff) component. In this case, the measurable quantity is the HDRF. The definition of the HDRF is analogous to that of the BRF but includes irradiance from the entire hemisphere (denoted with 2π):

(3) R HDRF = π sr d I r θ i , φ i , 2 π ; θ r , φ r d F 0 θ i , φ i , 2 π = R BRF ( θ i , φ i ; θ r , φ r ) f dir + R ( 2 π ; θ r , φ r ) 1 - f dir .

fdir denotes the fraction of direct incident radiation (i.e., fdir[0,1]). The second step in Eq. (3) assumes that the incident diffuse radiation is isotropic (for details, see Schaepman-Strub et al.2006). The additional integration over all reflection angles leads to the bihemispherical reflectance, generally called surface albedo α.

2.2 Modeling of the bidirectional reflectance

The shape of the BRDF is described by a weighted sum of trigonometric functions, generally referred to as kernels for volumetric scattering (Kvol), geometric scattering (Kgeo), and isotropic scattering (Kiso).

The kernel-driven semi-empirical Ross–Li model (Lucht et al.2000), which forms the basis of the MODIS 16 d BRDF/albedo product (Schaaf et al.2002), is applied within this study. The BRDF is given as a linear combination of the kernels with corresponding non-negative weighting functions fiso, fvol, and fgeo:

(4) f BRDF θ i , θ r , Δ φ , λ = f iso ( λ ) K iso + f vol ( λ ) K vol θ i , θ r , Δ φ + f geo ( λ ) K geo θ i , θ r , Δ φ ,

with the viewing (reflection) zenith angle θr, incidence zenith angle θi, relative viewing azimuth angle Δφ, and wavelength λ. The kernels,


depend on the scattering angle θ and the function 𝒪:


The functions 𝒞 and 𝒟 depend solely on the viewing and illumination geometry:


Originally, the concept of the model was developed studying typical features in observations of the bidirectional reflectance of various surface types (Roujean et al.1992). First, bare soil surfaces exhibit strong backscattering characteristics and show the effect of geometrical structures on the surface. Secondly, dense leaf canopies typically feature a minimum reflectance close to the nadir direction that increases off nadir for all azimuthal directions. In their approach, Roujean et al. (1992) represented the modeled BRDF as a linear combination of these two observational characteristics – with the former being accounted for by Kgeo, and the latter by Kvol. The volumetric kernel Kvol stems from volume scattering radiative transfer models (Ross1981). Thereby, randomly located facets that absorb and scatter incident radiation are modeled under the single-scattering approximation. The formula given in Eq. (6) corresponds to the Ross-Thick kernel for a dense leaf canopy. The geometric kernel Kgeo is derived from surface scattering and geometric shadow casting theory (Li and Strahler1992) and expresses effects caused by intercrown gaps within vegetation. It represents randomly oriented vertical protrusions on a flat and horizontal surface that isotropically reflect radiation. Equation (7) gives the Li-Sparse geometric kernel in reciprocal form describing the casting of shadows by a sparse ensemble of surface objects.

Hu et al. (1997) validated kernel-based BRDF models with 27 multi-angular data sets from various land cover types. The accuracy of the model was high. The correlation coefficient between model and observations was above 0.7 for all and above 0.9 for more than half of the data sets. Although originally calculated for surfaces covered with vegetation, the Ross-Thick and Li-Sparse-Reciprocal kernels are also applied for snow surfaces for the MODIS BRDF/albedo product. Stroeve et al. (2005) and Stroeve et al. (2013) assessed the accuracy of the 16 d albedo product of 11 years of measurements at 17 automatic weather stations on the Greenland ice sheet. They retrieved physically realistic ice sheet albedo values with an overall mean bias between MODIS and the in situ measurements of 0.022.

2.3 Inversion of semi-empirical kernel-driven model

The main benefit of retrieving the weighting functions fiso, fvol, and fgeo from the HDRF measurements is that they encompass all of the angular radiance information. The retrieval is done by inverting the RTLSR BRDF model. Thus, information about the complete two-dimensional shape of the HDRF is used.

The modeled fBRDFθi,θr,Δφ,λ from Eq. (4) can be written in the form of a sum as

(12) f BRDF , l = k = 1 3 f k K kl .

The spectral dependence is omitted here, and the index l denotes a particular viewing and illumination geometry θi,θr,Δφl. Considering an observation with N directional measurements ρl (l=1,,N), the error function 2 is defined as the difference between the observed and the modeled reflectances such that

(13) E 2 = 1 d l = 1 N ρ l - f BRDF , l 2 w l .

The degree of freedom d is (N−3) and wl denotes weights which are assigned to the respective observations. In general, wl could take the values 1, ρl, or ρl2. The goal of the inversion is the determination of the model weighting functions fk such that 2 is minimized. Strahler et al. (1996), amongst others, presented the analytical solutions for the three fk following this minimization. and Lewis (1995) showed that 2 has a global minimum in fk.

The inversion depends on the choice of the error function 2 and, thus, on the choice of the weights wl. Minimizing the absolute error (using weights equal to unity) leads to smaller errors in angular domains with a large reflectance. In contrast, minimizing the relative error (e.g., wl=ρl,ρl2) performs better in angular regions with lower reflectance. Subsequently, different choices for wl were tested for the retrievals performed in this work, and wl=ρl2 was chosen as it produced the lowest retrieval errors across the entire angular domain.

Figure 1(a) Illustration of research flights with Polar 6 in Dronning Maud Land between 24 December 2013 and 5 January 2014. The five flights on 28 December 2013 (marked in red) and the two flights on 2 January 2014 (green and brown) are used for the case studies. The parameter ranges covered by the case studies are marked with the same colors within the histograms in Fig. 2. (b) Photograph of Polar 6 illustrating the mounting of the digital camera, SMART, and the laser scanner.

3 Measurements and instrumentation

The airborne measurements were performed with the Polar 6 research aircraft (Wesche et al.2016), operated by the Alfred Wegener Institute, Helmholtz Centre for Polar and Marine Research (AWI), Bremerhaven, Germany. Polar 6 was based at Kohnen Station (750 S, 04 E). Figure 1a illustrates the tracks of the research flights with Polar 6 (60 flight hours) covering an area of 1000 km× 1000 km around Kohnen Station in Dronning Maud Land. The flights were performed between 24 December 2013 and 5 January 2014. The observations comprised regions with higher (coastal areas) and lower precipitation amounts and frequencies (inner Antarctica), as well as a variety of surface roughness structures. An overview of the airborne instrumentation is given in Carlsen et al. (2017); the mounting of the digital camera, laser scanner, and the Spectral Modular Airborne Radiation measurement sysTem (SMART) for surface albedo measurements on the underside of the aircraft is illustrated in Fig. 1b. The variability of the solar zenith angle is illustrated in Fig. 2a in the form of a histogram. During the flights, the solar zenith angle varied between 49 and 73.

Figure 2Variability of different parameters during research flights shown as histograms of the probability density function (PDF). The parameter ranges covered by the case studies are marked with the same colors as the corresponding flight tracks in Fig. 1. (a) Variability of solar zenith angle (bin size: 1). (b) Variability of surface roughness (bin size: 4 cm). Note the logarithmic scale. (c) Variability of optical-equivalent snow grain size (bin size: 20 µm).


3.1 Surface roughness measurements using a laser scanner

The snow surface topography was measured using the airborne laser scanner RIEGL VQ-580. In time-of-flight laser ranging, a near-infrared laser beam (1064 nm wavelength) is emitted downward and subsequently reflected upward by the snow surface before the echo is acquired by the sensor. From the time lag between emission and detection, the distance to the ground is calculated with a precision of about 25 mm depending on flight altitude. A fast-rotating polygonal mirror with a FOV of 60 and 10 to 150 scans per second allows for fully linear, unidirectional, and parallel scan lines.

After a correction for the aircraft attitude, a 1 km× 1 km digital elevation model (DEM) is generated from the geotagged laser point cloud with a resolution of 1 m. Subtracting the large-scale topography (smoothed DEM), the residual field contains the roughness information. The standard deviation (SD) of the residual field is interpreted as the surface roughness at the central coordinate of the DEM. Thus, roughness data are given for one data point per 1 km along the flight track of the aircraft. The uncertainty of the absolute height measurements (used for DEM generation) is less than 0.1 m. The relative analysis applied to the measurements yields even higher accuracy. The range in surface roughness is illustrated in Fig. 2b in the form of a histogram. The snow surface was generally smooth with roughness structures mostly below 5–10 cm (note the logarithmic scale; maximum: 2.2 m).

3.2 Optical-equivalent snow grain size as retrieved from spectral surface albedo measurements

Solar spectral radiation measurements were conducted using SMART (Wendisch et al.2001; Ehrlich et al.2008). The upward and downward spectral irradiance [F(λ), F(λ)] is measured within 0.3 to 2.2 µm wavelength applying an active horizontal adjustment system of the optical inlets to compensate for aircraft movement. From the irradiance measurements, the spectral snow surface albedo was obtained (Wendisch et al.2004). The spectral resolution is 2 to 3 nm between 0.3 and 1.0 µm and 15 nm up to 2.2 µm wavelength; the temporal resolution is in the order of 1 s. Combining the errors associated with the signal-to-noise ratio (1.3–3.0 %), the accuracy of the dark correction (0.1 %), the wavelength calibration (1.0 %), the accuracy of the cross-calibration (1.0–4.5 %), the non-ideal cosine response of the optical inlets (3.5 %), and the horizontal stabilization (1.0 %), the surface albedo measurements with SMART have an estimated uncertainty between 4.1 and 8.1 % depending on wavelength (Carlsen et al.2017).

The optical-equivalent snow grain size is defined as the radius Ropt of a collection of spheres with the same volume-to-surface ratio compared to the actual non-spherical snow grains. The retrieval of Ropt from spectral surface albedo measurements is described in Carlsen et al. (2017). In principle, the snow grain size and pollution amount (SGSP) algorithm (Zege et al.2011) was extended to spectral ratios of surface albedo at 1280 and 1100 nm wavelength. Being independent of systematic measurement uncertainties (e.g., cross-calibration of the optical inlets), this approach decreases the uncertainty of the retrieved Ropt compared to the single-wavelength approach. The retrieved Ropt from SMART and analogous ground-based measurements were compared to grain size observations utilizing reflectance measurements with MODIS (Carlsen et al.2017). The grain size as retrieved from SMART measurements (mean value: 105 µm) is slightly higher than the MODIS retrievals (89 µm). The larger data set of ground-based observations (72 µm) agrees well with the MODIS retrievals within the ranges given by the measurement uncertainties (linear correlation coefficient: 0.78; root-mean-square error: 24 µm).

The variability of Ropt during the research flights is illustrated in Fig. 2c in the form of a histogram. The optical-equivalent snow grain size varied between 16 and 480 µm but was mostly below 120 µm depending on precipitation and snow metamorphism processes.

3.3 Directional radiance measurements using a digital camera

3.3.1 Camera specifications

A Canon EOS-1D Mark III digital camera was used for the airborne HDRF measurements of the snow surface. It is a digital single-lens reflex camera with a complementary metal oxide semiconductor (CMOS) image sensor. The CMOS sensor covers 3908 pixel× 2600 pixel on a sensor area of 28.1 mm× 18.7 mm (Advanced Photo System (APS-H) format). During the observations, the camera was configured with a 8 mm F3.5 EX DG Circular Fisheye lens by Sigma. Color information is obtained by means of a color filter array (CFA) in front of the CMOS sensor. The CFA consists of red (R), green (G), and blue (B) color filters, which determine the spectral response of the underlying pixel.

In principle, the measured signal s (in digital numbers, DN) should be saved essentially verbatim in the proprietary RAW format CR2 (Canon RAW version 2). This would ensure the highest flexibility in the data analysis as no interpolation or white balance color correction is applied to the raw data. Full control over the post-processing steps is a prerequisite for radiometric measurements. During the airborne measurements, the camera was set to produce sRAW output format (small RAW) with a resolution of 1944 pixel× 1296 pixel. At some point, internal chrominance subsampling is applied by the camera (Kerr2015).

The sRAW photos are decoded utilizing the open-source tool dcraw (Coffin2017). The darkness level was set to 0 DN (see below) and the saturation level to 16 383 DN as the camera captures images with color depth of 14 bit. In between, the raw data are linearly interpolated. The multipliers for all channels are set to 1, meaning no white balance color correction is applied. The dynamic range of the output file is 16 bit (saturation at 65 536 DN).

To characterize the linearity of the sensor, photos from the radiation emitted by an integrating sphere with varying output intensities were taken with the camera in the laboratory (Carlsen2018). The camera sensor showed a linear sensitivity over a large part of the dynamic range up to 53 000 DN; the coefficient of correlation was exceptionally high at 0.99998. Photographs that measured higher signals (0.4 % of the total number) were excluded from the analysis.

Figure 3Relative spectral response functions for the three camera channels (colored lines: R, G, and B). The vertical dotted lines denote the respective central wavelengths. Dashed black line: RSR function for MODIS spectral band 4.


Compared to the dynamic range, the dark current (5 DN) and the read-out noise (7 DN) can be neglected.

Information on the spectral response of the camera is needed for the radiometric calibration. It determines the extent to which radiation of a certain wavelength passes the fish-eye lens as well as the CFA and is registered by the photodiodes. In this regard, the relative spectral response (RSR) function (in units of nm−1) of the three camera channels is defined as

(14) RSR ( λ ) = T c ( λ ) 0 T c ( λ ) d λ .

Tc(λ) with c=R, G, B denotes the dimensionless spectral transmission coefficients. The RSR function is normalized such that

(15) 0 RSR ( λ ) d λ = 1 .

The RSR function of the camera was measured in the laboratory using an integrating sphere as a radiation source and varying the outgoing radiation between 300 and 750 nm in increments of 5 nm wavelength by means of a grating monochromator. Figure 3 shows the measured RSR functions for the three camera channels. The central wavelengths of the non-Gaussian RSR functions are 602 nm (R), 538 nm (G), and 470 nm (B). Their full width at half maximum (FWHM) varies between 68 nm for the blue channel and 95 nm for the green channel.

The FOV of the camera is 180 due to the optics of the fish-eye lens. However, the camera is installed slightly above the lower aircraft body so that the camera is protected especially during take-off and landing. Therefore, parts of the aircraft frame are constantly in the FOV of the camera, which is why the effective FOV is reduced to approximately 160.

3.3.2 Radiance calibration and image post-processing

The post-processing of each raw image involves (a) the radiometric calibration, (b) the geometric calibration, (c) the aircraft attitude correction, and (d) the calculation of the HDRF. A detailed description of the different steps can be found in Carlsen (2018).

The measured signal was converted into the physical quantity of radiance (units of Wm-2nm-1sr-1) by means of a radiometric calibration in the laboratory. The pixels respond differently to a uniform illumination due to manufacturing tolerances (e.g., irregularities in the silicon used), contamination with dust particles, and optical effects at the edges of the lenses. These deviations lead to the photo response non-uniformity (PRNU) and need to be corrected.

The camera was positioned in front of an integrating sphere that served as a uniform radiation source. However, a distinct decrease in brightness is visible from the center to the edges of the sensor. This vignetting effect is typical for digital cameras (e.g., Lebourgeois et al.2008). From the data sheet, the spectral radiance Isphere(λ) that the integrating sphere emits at the specific optometer current is known. The calibration factor kc is defined at each pixel location on the sensor (row x, column y) as

(16) k c ( x , y ) = I sphere ( λ ) s ( x , y ) t exp

and carries the unit of Wm-2nm-1sr-1(DN/s)-1. During the calibration, the exposure time texp was set to 1∕1000s. Not only does kc correct for the PRNU; it simultaneously performs the absolute calibration, transforming the measured digital signal into the physical quantity of radiance (in units of Wm-2nm-1sr-1).

The geometric calibration of the camera–lens system relates each sensor pixel to a viewing zenith and azimuth angle (θv, φv). Often, the process of geometric calibration of a camera involves calibration equipment or the use of planar targets such as checkerboard patterns (e.g., Tsai1987; Urquhart et al.2016). Within this work, a stellar calibration method is applied (e.g., Schmid1974; Klaus et al.2004; Mori et al.2013; Urquhart et al.2016) utilizing the high precision to which the positions of stellar objects are known. In this regard, 684 different star positions identified in subsequent pictures of the night sky were utilized to calculate θv and φv for each sensor pixel.

Figure 4Exemplary post-processing and retrieval of model parameters for HDRF measurement. (a) Raw image taken on 2 January 2014 at 08:16 UTC (θ0=58.9). (b) Polar plot of observed snow HDRF ρl calculated from raw image for camera channel G. The image is rotated in the azimuthal direction of the Sun such that the Sun is on the left. (c) Averaged HDRF measurement using three images; the footprint is 1670 m. (d) Modeled snow HDRF RHDRF,l after retrieval of model parameters fiso=1.12, fvol=0.17, and fgeo=0.01 from (c). (e)  Relative difference (ρl-RHDRF,l)/RHDRF,l. The accumulated RMSE is 0.04.


3.4 Calculation and uncertainties of measured surface HDRF

The camera took pictures with a temporal resolution of about 8 s. Exemplarily, Fig. 4 demonstrates the derivation of the snow HDRF from a raw image taken by the camera on 2 January 2014 (08:16 UTC; see Fig. 4a). First, the pixels receiving radiation from the direction of the aircraft frame are excluded. For each pixel location (x, y), the radiance I(x,y) (in units of Wm-2nm-1sr-1) is calculated from the measured signal s using the absolute calibration factor kc and the exposure time texp:

(17) I ( x , y ) = s ( x , y ) t exp k c ( x , y ) .

The camera viewing angles are calculated from the geometric calibration for each pixel (x, y). As the camera is fixed to the aircraft frame, a correction for the aircraft attitude was implemented to obtain the reflection angles θr and φr. Thus, beside the geometric calibration, the observed reflection angles are determined by the attitude angles of the aircraft. Utilizing the data from the internal navigation system and the GPS on Polar 6, the viewing angles are corrected depending on the roll, pitch and yaw angles of the aircraft at the time of measurement. In this regard, Euler rotations are applied as described in Ehrlich et al. (2012). Finally, the observed HDRF is calculated by

(18) R HDRF ( θ 0 , φ 0 ; θ r , φ r ) = π sr I ( θ r , φ r ) F θ 0 , φ 0 .

F(θ0,φ0) denotes the downward solar irradiance at the time of measurement (with solar zenith and azimuth angles θ0 and φ0). Based on the resulting reflection angles, a polar plot of the measured RHDRF is created (camera channel G; see Fig. 4b). To achieve comparability, each image is rotated in the azimuthal direction of the Sun.

The downward irradiance measurements from SMART could not be used for the calculation of the HDRF due to calibration issues. Note that this pertains to the radiometric calibration only, which converts the digital numbers registered by the spectrometer into units of irradiance. For the albedo measurements with SMART, a relative calibration of the upper and lower sensors is sufficient and an absolute radiometric calibration is not required. Thus, the albedo measurements and the retrieval of the optical-equivalent snow grain size are unaffected by this calibration issues. For the calculation of the HDRF, the global irradiance was simulated along the flight track with the Library for Radiative Transfer (libRadtran; Mayer and Kylling2005) using the Discrete Ordinate Radiative Transfer (DISORT) solver by Stamnes et al. (1988) instead. Vertical profile data of air temperature and relative humidity measured by radiosondes at Kohnen Station, Neumayer Station III (König-Langlo2014), and Halley Research Station (Durre et al.2016) were used as input for the simulations. The simulated irradiance was integrated over the wavelength range of each camera channel and weighted with the RSR function of the camera. The use of simulations limits the validity of absolute values of the measured HDRF to cloudless conditions. Cloudless cases (approximately 75 % of the camera observations) were identified from visual synoptic observations during the flights as well as using the downward irradiance measured with SMART. Periods with fast fluctuations in the downward irradiance were flagged as cloudy. Although mainly the shape of the HDRF is analyzed within this work (which is independent from the absolute value of F), the analysis is restricted to cloudless conditions only.

The overall relative uncertainty of the HDRF measurements with the digital camera ranges in the order of 4.5 %. The error in the absolute value of RHDRF might be higher depending on the atmospheric conditions due to the usage of simulated values for the global irradiance F in Eq. (18). The uncertainties in the HDRF measurements stem from sensor characteristics (estimated as 0.5 % due to signal-to-noise ratio, dark current, linearity, read-out noise, and chrominance subsampling), the radiometric calibration (4 %), the geometric calibration, and the correction for the aircraft attitude (combined estimate of 1.0 %).

3.5 Averaging

From trigonometric considerations, the footprint of the camera (twice the radius r of the disc on the ground pictured by the camera) depends on the flight altitude z and the FOV. At an altitude of 100 m, the radius is 570 m and the footprint approximately 1 km. The footprint grows with increasing altitude to 5.7 km (z=500m) and 11.3 km (z=1000m). To get independent of local roughness features, averaging is performed over time intervals of around 30 s (including three to four pictures). This interval was chosen trading off the gained independence from local features for a growing footprint and thus reduced comparability with the MODIS BRDF measurements that are provided on a 500 m× 500 m grid. Generally, averaging over 30 s yields footprints between 1 and 2 km. However, individual footprints depend on flight altitude, exact averaging time, and aircraft velocity. The averaged HDRF for the exemplary measurement on 2 January 2014 for a footprint of 1670 m is shown in Fig. 4c.

3.6 Approximation of surface BRF with HDRF measurements

For each averaged HDRF, the weighting functions fiso, fvol, and fgeo are retrieved by inverting the RTLSR BRDF model and setting the HDRF equal to the BRDF in Eq. (12). This is necessary as the BRDF is not measurable under atmospheric conditions. However, the atmospheric influence is assumed to be small due to the high surface elevation and the dry, aerosol-free atmosphere on the Antarctic Plateau: the mean AOD at Kohnen Station was 0.015 between 2001 and 2006 (at 500 nmTomasi et al.2007), and the mean integrated atmospheric water vapor was 1.1 kg m−2 between December 2013 and January 2014 (Carlsen et al.2017). Simulations with libRadtran (Mayer and Kylling2005) of the mean direct fraction of the global irradiance fdir at Kohnen Station between 10 December 2013 and 31 January 2014 indicate rather weak atmospheric scattering effects. For the central wavelengths of the three camera channels, the mean direct fraction is 87 % (for 602 nm wavelength, channel R), 81 % (538 nm, G), and 69 % (470 nm, B). The values for fdir are even higher when contributions from a small scattering cone around the incident direction are considered. Schaepman-Strub et al. (2006) simulated the difference between the HDRF and BRF (RBRF=πsrfBRDF; compare Eq. 2) for snow surfaces for different fractions fdir at 550 nm wavelength. They found that, with an increasing diffuse fraction of the incident irradiance, the shape of the HDRF is smoothed in comparison to the BRF (fdir=1; see Eq. 3). For fdir=0.8 at θ0=30, the shape of the HDRF is still close to that of the BRF. Hence, the measured HDRF shapes can serve as approximations for the intrinsic BRF of the underlying snow surface for the measurements analyzed within this study (fdir>0.8).

3.7 Quality of the inversion

An exemplary inversion is performed on the snow HDRF measurement from 2 January 2014 at 08:16 UTC (θ0=58.9). Figure 4c shows the averaged observations that are used as input for the algorithm. The inversion leads to the retrieved model parameters of 1.12 (fiso), 0.17 (fvol), and 0.01 (fgeo). With the calculated kernels Kkl, the modeled HDRF is obtained (see Fig. 4d). The shapes of modeled and measured HDRF are similar; their relative difference is mostly within the range of ± 5 % (see Fig. 4e). The largest deviations occur in the forward scattering direction (up to 3 %) because the location of maximum is not mimicked perfectly by the modeled HDRF. However, up to viewing zenith angles of around 75, the relative difference is below 1 %.

Figure 5Statistics of the quality of the inversion. RMSE of all airborne retrievals (PDF; bin size: 0.01). Vertical dashed line shows RMSE threshold of 0.1 for MODIS full inversions.


The root-mean-square error (RMSE),

(19) σ RMSE = E = 1 N - 3 l = 1 N ( ρ l - R HDRF , l ) 2 w l ,

serves as a criterion for the quality of the inversion. Only RMSEs below 0.1 are considered full inversions within the MODIS retrieval (e.g., Stroeve et al.2005). For the example shown in Fig. 4, it is 4.0 %. σRMSE strongly depends on the anisotropy of the observed HDRF. Figure 5 shows the RMSE for all airborne retrievals as a histogram, with 96 % of the retrievals showing an RMSE below 0.1. Hence, the upper bound for the uncertainty of the inversion is in the range of 10 %.

The influence of the angular sampling on the retrieval of the model parameters can be investigated using the weights of determination (WoD; Lucht and Lewis2000), which depend both on the angular sampling and the number of directional measurements N. WoD (also: noise amplification factors) below 1 indicate a retrieval uncertainty that is smaller than the RMSE of the inversions. Due to the large N and consistent sampling of the viewing hemisphere, the WoD for the model parameters fk within this study are mostly below 0.0002. For comparison, the sparse angular sampling of the reflectance measurements from MODIS observations can lead to WoD above 1, and only model fits with WoD below 2.5 are considered full inversions.

To study the stability of the inversion as a result of the uncertainties in the reflectance measurements, the stability of the retrieval of fk was studied for a 5 min flight leg over a homogeneous snow surface on 28 December 2013: the surface roughness varied between 2 and 4 cm, Ropt varied between 80 and 90 µm, and the solar zenith angle was around 55.3. The inversion was performed for each individual camera picture (no averaging). The mean and SD of the retrieved model parameters along the flight leg were 1.104 ± 0.030 (fiso), 0.296 ± 0.027 (fvol), and 0.031 ± 0.003 (fgeo). The same area was sampled five times at different solar zenith angles on 28 December 2013; the influence of θ0 on the SDs of the fk is negligible. The mean SDs of 0.031 (fiso), 0.020 (fvol), and 0.003 (fgeo) are used as error bars for the fk within this study.

3.8 Comparison with BRDF derived from MODIS satellite observations

The use of the RTLSR model in the retrieval of the model parameters fk from airborne observations provides the framework to compare with satellite observations. The MODIS BRDF/Albedo Model Parameters product (MCD43A1) provides the fk at 500 m resolution utilizing multi-date, atmospherically corrected, and cloud-cleared input data over a period of 16 d (Schaaf et al.2002). Version 6 of the MCD43A1 product uses both Terra and Aqua data and is temporally weighted to the ninth day of a 16 d retrieval window, thereby putting greater emphasis on the actual day of interest than previous versions. The combination of MODIS on Terra and Aqua provides a higher number of high-quality reflectance measurements, resulting in better temporal and angular sampling. As the retrieval uncertainty increases with larger solar zenith angle, a thorough quality assessment of the satellite data is needed especially in polar regions. Stroeve et al. (2005) compared the MODIS albedo product with in situ observations on the Greenland ice sheet and recommended only using full-inversion data (RMSE < 0.1, WoD < 2.5). Hence, this study restricts the satellite data to full inversions based on the quality flag information provided in the MCD43A2 product (quality flag values: 0 or 1). Lower quality-magnitude inversions would facilitate archetypal BRDF parameters to improve the retrieval in case of less than seven high-quality reflectance observations or poor angular sampling. The BRDF data from MODIS spectral band 4 (0.545–0.565 µm) are used as it coincides best with the green camera channel (see Fig. 3).

The MODIS data were resampled on the flight track using the nearest neighbor with respect to the great circle distance. For the comparison, the aircraft measurements are filtered analogously: only clear-sky observations with RMSE < 0.1 (96 % of total measurements), WoD < 2.5 (valid for all observations), and a solar zenith angle lower than 70 were used. Stroeve et al. (2005) restricted the analysis to θ0< 75; however, the number of comparable observations between MODIS and the aircraft is not reduced by the stricter condition. In fact, only 434 (21 %) MODIS observations can be compared to the 2078 airborne observations that in principle fulfill the quality requirements. For the remaining observations, no full-inversion retrieval could be performed from the MODIS observations.

4 Results and discussion

To study the influence of the solar zenith angle, the surface roughness, and the optical-equivalent snow grain size on the snow HDRF, the effects of the individual parameters need to be separated. For this purpose, several case studies were selected: (a) measurements at different θ0 (marked with red circles in Fig. 2a, five flights on 28 December 2013), (b) measurements at different surface roughness (marked with green circles in Fig. 2b, morning flight on 2 January 2014), and (c) measurements at different Ropt (marked with brown circles in Fig. 2c, afternoon flight on 2 January 2014). The remaining two parameters were kept constant throughout the individual cases. After discussing the individual cases, the airborne observations are compared with the BRDF derived from MODIS.

4.1 Influence of solar zenith angle

The research flights with Polar 6 on 28 December 2013 were exploited to study the influence of θ0 on the snow HDRF. Five consecutive flights were conducted in an area northeast of Kohnen Station (see red flight track in Fig. 1a).

Figure 6(a) Averaged snow HDRF (490–585 nm wavelength) over 30 s segment during easternmost flight leg for the five consecutive research flights on 28 December 2013. (b–d) Dependence of the retrieved fk on the solar zenith angle θ0 for the five flights. The dashed lines represent linear regression fits of the parameters with respect to θ0. (b) The isotropic weight fiso. (c) The volumetric weight fvol. (d) The geometric weight fgeo.


For each flight, the snow HDRF was averaged over a 30 s segment during the easternmost flight leg; the results are shown in the top panel of Fig. 6. The solar zenith angle varied between 51.8 and 71.6. As Polar 6 took identical routes throughout the flights, the respective flight leg always covered the same area. The central points of the five consecutive segments are separated by less than 1 km. Therefore, the surface roughness (retrieved from the laser scanner) of approximately 6 to 7 cm remained constant. Within the retrieval uncertainties, the same holds true for the optical-equivalent snow grain size. Ropt varied between 70 and 85 µm. This minimizes a possible influence of Ropt or lrough on the snow HDRF that would be superimposed on the effect of θ0. Cloudless conditions prevailed during the flights.

With increasing θ0, the maximum of the HDRF in the forward scattering direction becomes more pronounced (see Fig. 6). Correspondingly, the anisotropy gets larger. Figure 6 shows the dependence of the fk on θ0 for the green camera channel (490–585 nm wavelength). fiso shows no clear dependence on θ0 and varies between 1.06 and 1.10. fgeo weakly increases with θ0 from 0.03 at θ0=51.8 to 0.05 at 71.6. During the first three research flights, θ0 varied only slightly. Therefore, the strongest trends are visible between flights 3 and 5. In particular, fvol increases from 0.25 to 0.36.

The probability of photons entering the snowpack to leave it after just a few scattering events increases for lower Sun elevation. With increasing θ0, the reflection properties of the snow layer converge to the single-scattering properties of ice crystals. In addition, longer shadows are cast by surface roughness structures at larger θ0. The expected increase in the anisotropy of the snow HDRF for increasing θ0 is obvious in the changing model parameters fk and has been observed in earlier studies (e.g., Dirmhirn and Eaton1975; Kuhn1985; Warren et al.1998; Hudson et al.2006).

4.2 Influence of surface roughness

On 2 January 2014, a research flight was conducted from Kohnen Station in the direction of the coastline to Neumayer Station III. The flight track is shown in Fig. 1a (green). Between 08:50 and 09:00 UTC, Polar 6 crossed a sastrugi field at an elevation between 100 and 250 ma.s.l.. This fostered large variations in the surface roughness between 42 cm and 2.2 m which became distinguishable by visual observation from the aircraft. At the same time, the solar zenith angle was constant at around 55.5. For these values, only little variation of the snow HDRF with θ0 was found (compare Fig. 6). The optical-equivalent snow grain size varied between 240 and 320 µm; the flight altitude remained at around 600 m above the surface, leading to a footprint of approximately 3.5 to 4.0 km. During this flight leg, 14 measurements were used to study the influence of the surface roughness on the snow HDRF. Two examples of the snow HDRF are shown in the top panel of Fig. 7.

Figure 7(a) Two examples of snow HDRF (490–585 nm wavelength) at different times during an overflight over a sastrugi field on 2 January 2014 (morning flight). (b–d) Dependence of the retrieved fk on the surface roughness. The dashed lines represent linear regression fits of the parameters with respect to lrough. (b) The isotropic weight fiso. (c) The volumetric weight fvol. (d) The geometric weight fgeo.


For the smooth surfaces, the maximum in forward scattering direction is more pronounced. With increasing surface roughness, the anisotropy gets lower. Roughness structures enhance the backscatter by changing the effective angle of incidence. In addition, they cast shadows that reduce the forward scatter. The evolution of the retrieved model parameters fk as illustrated in Fig. 7 demonstrates the decreasing anisotropy: fiso and fgeo show a slight decreasing trend, whereas, most prominently, fvol decreases from 0.33 to 0.24 with increasing surface roughness. This is in accordance with Warren et al. (1998), who observed a reduction of the forward reflection peak due to sastrugi during tower measurements at South Pole Station. The effect of surface roughness on the BRDF has been previously investigated in observational (e.g., Grenfell et al.1994; Hudson and Warren2007) and modeling studies (e.g., Leroux and Fily1998; Zhuravleva and Kokhanovsky2011) and is facilitated by the multi-angular reflectance data from MISR to map surface roughness over ice sheets and sea ice from space (Nolin et al.2002).

Figure 8(a)  Two examples of snow HDRF (490–585 nm wavelength) at different times on 2 January 2014 (afternoon flight). (b–d) Dependence of the retrieved fk on the optical-equivalent snow grain size. The dashed lines represent linear regression fits of the parameters with respect to Ropt. (b) The isotropic weight fiso. (c) The volumetric weight fvol. (d) The geometric weight fgeo.


4.3 Influence of optical-equivalent snow grain size

The separation of the two parameters θ0 and lrough showed that they influence the snow HDRF in opposite ways. Ice absorption becomes dominant in the near-infrared part of the solar EM spectrum. Hence, even though a wide variety of optical-equivalent snow grain sizes was covered during the flights (16 to 480 µm), no large effect on the snow HDRF was expected as the camera is only sensitive in the visible wavelength range (Warren et al.1998, compare Fig. 3). However, 32 observations between 14:45 and 15:15 UTC on 2 January 2014 (see brown flight track in Fig. 1a), depicted in Fig. 8, show a decrease of fvol from 0.34 at 184 µm to 0.17 at 325 µm snow grain size. This is consistent with the decrease in anisotropy as visible in the averaged HDRF measurements (top panel in Fig. 8).

Figure 9Comparison between airborne (490–585 nm wavelength) and MODIS (545–565 nm) retrievals of fk for the research flight on 29 December 2013. (a) Map of the area covered by the research flights; the colored contours correspond to the MODIS retrieval of fiso. Black: flight track on 29 December 2013. Colored circles on black flight track: airborne retrieval of fiso. (b) As in (a) but for fvol. (c) As in (a) but for fgeo. (d) Corresponding combined scatter plot and histograms between airborne and satellite (resampled on flight track) retrieval for fiso. Error bars correspond to the mean SD of the airborne retrieval. Dashed line in scatter plot: 1:1 line. Dashed line in histograms: center of axis. (e) As in (d) but for fvol. (f) As in (d) but for fgeo.

However, this stands in contrast to findings by Steffen (1987), who measured the bidirectional reflectance on a glacier in northwestern China for powder snow (grain size: 0.15 mm), new snow (0.5–1.0 mm), and old snow (1–3 mm) and found higher anisotropy for larger snow grains. It should be noted that the range of snow grain sizes covered within this study is significantly smaller and the footprint of the measurements is larger compared to the local study by Steffen (1987). Even though other parameters were kept constant during this case study (θ0: 51.3–53.4; lrough: 4–7 cm; footprint: 1.3–1.7 km), effects of other influencing parameters cannot be ruled out completely. For example, Steffen (1987) reported changing ice crystal shapes between the measurements and suggested a possible influence of dust particles within the snow. This could partly explain the opposite conclusions drawn within this study.

4.4 Parameterization of snow HDRF

The most significant effects of θ0, lrough, and Ropt were observed on the model parameter fvol, for which the parameterizations in terms of linear regression fits (dashed black lines in Figs. 6c–8c) are given by


The linear correlation coefficients (R2) are 0.98 (for θ0), 0.75 (lrough), and 0.68 (Ropt). Note that these parameterizations are only valid for the wavelength range of about 490–585 nm (green camera channel) and the ranges of all three parameters as given within the descriptions of the case studies.

4.5 Comparison with satellite observations

Figure 9 shows the comparison between retrieved fk from airborne and satellite measurements for the research flight on 29 December 2013, which showed the largest amount of full-inversion retrievals for MODIS (235). The flight covered a triangular pattern southeast of Kohnen Station. The top panel shows contour plots with color coding for the three parameters as retrieved from MODIS on that day: (a) fiso, (b) fvol, and (c) fgeo. The colored circles along the flight track (black) correspond to the airborne-retrieved fk. The corresponding scatter plots between the MODIS and aircraft retrievals are shown in Fig. 9d–f. The isotropic and geometric model weights fiso and fgeo from airborne and MODIS retrievals are in the same order of magnitude and partly agree within the measurement uncertainties. fiso values from aircraft are slightly higher than the corresponding MODIS retrievals. This is evident from the difference in the probability density functions (PDFs) between aircraft and MODIS retrievals as illustrated in Fig. 9d–f. The comparison of fgeo values shows no robust difference between airborne and MODIS retrievals. However, the lowest correlation is found for the volumetric weights fvol. fvol values retrieved from aircraft range mostly between 0.1 and 0.3 and are much larger than the MODIS retrievals. The map in Fig. 9b shows that the flight covered an area where MODIS retrieved fvol values of exactly or close to 0. In particular, this holds true for a latitudinal band around 76 S. Even though full inversions were reported for this area, this might still be a measurement artifact. But even for the cases with nonzero values, the aircraft retrievals lead to much higher fvol values and, thus, larger anisotropy of the surface reflection than shown in the MODIS retrieval.

The combined scatter plots and histograms for the three model weights including all 434 observations from seven different research flights are shown in Fig. 10 and indicate that the conclusions drawn from the flight on 29 December 2013 are valid for all research flights: fvol shows the strongest discrepancies between airborne and MODIS retrievals, whereas fiso values are slightly higher when retrieved from the aircraft measurements. The correlation between airborne and satellite measurements is low for all three model parameters; R2 is 0.18 (for fiso), 0.004 (fvol), and 0.055 (fgeo).

Various factors can contribute to the observed differences: (1) data quality, (2) short-term changes in snow properties, (3) footprint size differences, and (4) inherent challenges comparing the camera and satellite data.

First, a thorough quality assessment for both data sets formed the prerequisite for the comparison. Only airborne data from cloudless conditions (at solar zenith angles below 70) with a RMSE lower than 0.1 were compared with full-inversion MODIS retrievals (quality flag values: 0 or 1; RMSE < 0.1). However, possible artifacts such as the latitudinal band around 76 S with values of exactly 0 cannot be ruled out entirely.

Figure 10(a) Scatter plot between airborne (490–585 nm wavelength) and MODIS (resampled on flight track, 545–565 nm) retrieval for fiso including all seven research flights used for the comparison. Error bars correspond to the mean SD of the airborne retrieval. Dashed line: 1:1 line. In addition, the histograms for both airborne and satellite retrievals are given. Dashed line: center of axis. (b) As in (a) but for fvol. (c) As in (a) but for fgeo.


Secondly, the case studies demonstrated the strong influence of the solar zenith angle, the surface roughness, and the optical-equivalent snow grain size on the anisotropy of the reflection. Precipitation events or drifting snow altering the optical-equivalent snow grain size and the surface roughness can change the surface properties on short timescales. These short-term changes in snow properties are not captured well by the 16 d BRDF product from MODIS. Even though the 434 observations stem from seven different research flights, most of the data originate from the two flights on 29 December 2013 and 2 January 2014 (morning flight). For most of the remaining flights, no full-inversion retrieval was possible for MODIS. Being restricted to mainly two different days, short-term changes in snow properties can indeed explain the low correlation between the airborne and MODIS retrievals. Also, the solar zenith angle at the time of the aircraft measurements varies from the different zenith angles during the satellite overpasses used for the MODIS retrieval.

Thirdly, the airborne measurements still included varying footprint sizes (mostly between 1 and 3 km) compared to the 500 m resolution of the MODIS product. This was the result of trading off comparability with MODIS for data coverage. Also, no correlation between the fk retrieved from the camera measurements and the footprint size was found (not shown here).

Lastly, comparing the airborne camera data with the MODIS product poses challenges that are inherent to the measurement setup. The observation times of the MODIS product do not coincide with the times of the research flights. MODIS spectral band 4 was chosen as its RSR function agrees best with the green camera channel. As the HDRF is wavelength dependent, differences naturally occur due to different spectral coverage but should be small.

Thus, despite the thorough quality assessment analogously applied to both airborne and satellite data, a quantitative comparison between the two retrievals seems challenging. Nevertheless, the comparison highlights the difficulties to accurately measure the snow BRDF from MODIS as well as how quickly snow properties can change, leading, in this study, to an underestimation of the anisotropic reflection. For the analyzed cases, MODIS underestimates the volumetric model weight fvol on average by at least a factor of 10 compared to the airborne measurements (neglecting MODIS values of exactly 0). Many satellite retrievals of surface albedo and atmospheric parameters over snow surfaces rely on a correct angular modeling of the surface reflectance. Thus, inaccuracies in the surface BRDF can propagate and influence other retrieval results.

5 Conclusions

The HDRF of snow surfaces in the wavelength band 490–585 nm was derived from airborne measurements in Antarctica during austral summer in 2013/14. A digital 180 fish-eye camera installed on a research aircraft was used for this purpose. The surface roughness was retrieved using an airborne laser scanner, while the optical-equivalent snow grain size was calculated from spectral surface albedo measurements. While based at Kohnen Station (750 S, 04 E) at the outer part of the East Antarctic Plateau, the airborne measurements with the Polar 6 research aircraft covered an area of around 1000 km× 1000 km in Dronning Maud Land. Thus, the snow HDRF was obtained for a variety of conditions with different solar zenith angle (θ0), surface roughness (lrough), and optical-equivalent snow grain size (Ropt).

The digital camera provides airborne radiance measurements with high angular resolution. The technical characterization of the digital camera in the laboratory revealed excellent linearity as well as negligible dark current and read-out noise. From the measured angular radiances, the HDRF was calculated applying simulations of the global irradiance with libRadtran using the DISORT solver. The relative uncertainty of the HDRF measurements is estimated as 4.5 %. The footprint of the snow HDRF measurements analyzed within this study varies between 1 and 4.5 km.

Three case studies were investigated to separate the effects of θ0, lrough, and Ropt on the snow HDRF. With increasing solar zenith angle, the HDRF maximum in the forward scattering direction becomes more pronounced. Conversely, the backscatter is enhanced over rougher surfaces by changing the effective angle of incidence and by casting shadows. This was quantified by inverting the semi-empirical kernel-driven RTLSR model and parameterizing the snow HDRF with respect to θ0, lrough, and Ropt (Eqs. 2022). The uncertainty of the inversion is estimated as 10 %. The increased anisotropy (with increasing θ0) is mainly shown by the increase in the model parameter fvol. Vice versa, fvol decreases with increasing lrough. The snow grain size reveals an effect similar to surface roughness structures in terms of a decrease of the anisotropy with increasing Ropt, which differs from earlier findings (e.g., Steffen1987). Possible reasons are discussed (footprint size, ice crystal shape, and contamination with dust particles). However, to yield stronger dependence on the optical-equivalent snow grain size, a camera sensitive to radiation in the near-infrared part of the EM spectrum needs to be employed.

Continuous satellite observations of the snow surface albedo with global coverage critically depend on the angular modeling of the surface BRDF, which may introduce large errors especially in polar regions due to low Sun elevations and the anisotropic scattering phase function of snow crystals. The RTLSR model was chosen for the inversion as it forms the basis of the MODIS 16 d BRDF/albedo product. Thus, the retrieved model parameters fk from the HDRF measurements allow for direct comparison with BRDF products from satellite remote sensing. Despite the similar quality assessment applied to airborne and satellite data, large deviations were found especially for the volumetric model weight fvol. The airborne values for fvol are larger than the corresponding MODIS retrievals (by at least a factor of 10). Short-term changes in snow properties (precipitation and drifting snow), which are not captured by the 16 d MODIS retrieval, and different solar zenith angles at the time of measurements are discussed as likely reasons for the observed differences. Although the effects are presumably small, influences of a varying footprint size for the airborne observations and the different RSR functions of the green camera channel and MODIS spectral band 4 can further contribute to this discrepancies.

Generally, MODIS observations underestimated the anisotropy of the reflection by the snow surfaces. This has possible implications for satellite retrievals of surface albedo and atmospheric parameters over snow surfaces (e.g., AOD and cloud properties), which strongly depend on a correct angular modeling of the surface reflectance (e.g.,  Qu et al.2015).

The analysis of the snow HDRF measurements can be readily applied to other airborne campaigns. For example, the same digital camera with a similar setup of the Polar 6 research aircraft was used (a) during the Arctic CLoud Observations Using airborne measurements during polar Day (ACLOUD) between 22 May and 28 June 2017 based at Longyearbyen, Svalbard (Wendisch et al.2019; Ehrlich et al.2019; Jäkel and Ehrlich2019), and (b) during the Polar Airborne Measurements and Arctic Regional Climate Model Simulation Project (PAMARCMiP) between 12 March and 6 April 2018 at Station North, Greenland (Herber et al.2012).

Data availability

The tabulated data including the retrieved fiso, fvol, and fgeo as well as θ0, lrough, Ropt, and auxiliary information (time and position of the measurements) are available for all research flights in the Publishing Network for Geoscientific & Environmental Data (PANGAEA): (Carlsen et al.2020). The MODIS BRDF/Albedo Model Parameters product (MCD43A1, Version 6) and the quality data set (MCD43A2) were downloaded from NASA's Land Processes Distributed Active Archive Center (LP DAAC) website (, Schaaf and Wang2015a and, Schaaf and Wang2015b).

Author contributions

TC performed the post-processing and data analysis of the camera and surface albedo measurements, produced the figures, and wrote the paper with help from all authors. GB and MS conducted the measurements at Kohnen Station and with the Polar 6 research aircraft and helped write the paper. VH analyzed the laser scanner data for the surface roughness measurements and helped write the paper. EJ gave substantial input with regard to the data analysis and helped write the paper. AE and MW planned the measurements, designed the study, and helped write the paper.

Competing interests

The authors declare that they have no conflict of interest.


This work was supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) in the framework of the priority program “Antarctic Research with comparative investigations in Arctic ice areas” (SPP 1158) through grants WE1900/29-1 and BI 816/4-1. We gratefully acknowledge the funding by DFG – Project-ID 268020496 – TRR 172, within the Transregional Collaborative Research Center “ArctiC Amplification: Climate Relevant Atmospheric and SurfaCe Processes, and Feedback Mechanisms (AC)3”. This work was partly supported by the European Research Council (ERC) through grant StG 758005. We are grateful to the Alfred Wegener Institute, Helmholtz Centre for Polar and Marine Research, Bremerhaven, Germany, for supporting the campaign with logistics, the aircraft, and man power in Antarctica. In addition, we would like to thank Kenn Borek Air Ltd., Calgary, Canada, for the great pilots who made the complicated measurements possible. For his insights into the sRAW format of the digital camera, we would like to express our gratitude to Douglas Kerr. The map in Fig. 1a was generated using the Generic Mapping Tools (Wessel et al.2019).

Financial support

This research has been supported by the Deutsche Forschungsgemeinschaft (grant nos. WE1900/29-1, BI 816/4-1, and 268020496) and the European Research Council (grant StG 778005).

Review statement

This paper was edited by Marie Dumont and reviewed by two anonymous referees.


Aoki, T., Aoki, T., Fukabori, M., Hachikubo, A., Tachibana, Y., and Nishio, F.: Effects of snow physical parameters on spectral albedo and bidirectional reflectance of snow surface, J. Geophys. Res., 105, 10.219–10.236,, 2000. a

Bourgeois, C., Ohmura, A., Schroff, K., Frei, H.-J., and Calanca, P.: IAC ETH goniospectrometer: A tool for hyperspectral HDRF measurements, J. Atmos. Ocean. Tech., 23, 573–584,, 2006a. a

Bourgeois, C. S., Calanca, P., and Ohmura, A.: A field study of the hemispherical directional reflectance factor and spectral albedo of dry snow, J. Geophys. Res., 111, D20108,, 2006b. a

Brest, C. and Goward, S.: Deriving surface albedo measurements from narrow-band satellite data, Int. J. Remote Sens., 8, 351–367,, 1987. a

Carlsen, T.: Influence of snow properties on directional surface reflectance in Antarctica, PhD thesis, Faculty of Physics and Earth Sciences, Leipzig University, available at: (last access: 8 November 2020), 2018. a, b

Carlsen, T., Birnbaum, G., Ehrlich, A., Freitag, J., Heygster, G., Istomina, L., Kipfstuhl, S., Orsi, A., Schäfer, M., and Wendisch, M.: Comparison of different methods to retrieve optical-equivalent snow grain size in central Antarctica, The Cryosphere, 11, 2727–2741,, 2017. a, b, c, d, e

Carlsen, T., Birnbaum, G., Ehrlich, A., Helm, V., Jäkel, E., Schäfer, M., and Wendisch, M.: Hemispherical-directional reflectance factor (HDRF) of snow surfaces retrieved from airborne camera measurements in Antarctica during austral summer in 2013/14, PANGAEA,, 2020. a

Coffin, D.: Decoding raw digital photos in Linux, available at: (last access: 8 November 2020), 2017. a

Cook, R. and Torrance, K.: A Reflectance Model for Computer Graphics, ACM T. Graphic., 1, 7–24,, 1982. a

Cox, C. and Munk, W.: Measurement of the roughness of the sea surface from photographs of the sun's glitter, J. Opt. Soc. Am. A, 44, 838–850, 1954. a

Dirmhirn, I. and Eaton, F.: Some characteristics of the albedo of snow, J. Appl. Meteorol, 14, 375–379,<0375:SCOTAO>2.0.CO;2, 1975. a

Dumont, M., Brissaud, O., Picard, G., Schmitt, B., Gallet, J.-C., and Arnaud, Y.: High-accuracy measurements of snow Bidirectional Reflectance Distribution Function at visible and NIR wavelengths – comparison with modelling results, Atmos. Chem. Phys., 10, 2507–2520,, 2010. a

Durre, I., Xungang, Y., Vose, R. S., Applequist, S., and Arnfield, J.: Integrated Global Radiosonde Archive (IGRA), Version 2. Radiosonde data from Halley research station (AYM00089022), NOAA National Centers for Environmental Information, (last access: 24 October 2018), 2016. a

Ehrlich, A., Bierwirth, E., Wendisch, M., Gayet, J.-F., Mioche, G., Lampert, A., and Heintzenberg, J.: Cloud phase identification of Arctic boundary-layer clouds from airborne spectral reflection measurements: test of three approaches, Atmos. Chem. Phys., 8, 7493–7505,, 2008. a

Ehrlich, A., Bierwirth, E., Wendisch, M., Herber, A., and Gayet, J.-F.: Airborne hyperspectral observations of surface and cloud directional reflectivity using a commercial digital camera, Atmos. Chem. Phys., 12, 3493–3510,, 2012. a, b

Ehrlich, A., Wendisch, M., Lüpkes, C., Buschmann, M., Bozem, H., Chechin, D., Clemen, H.-C., Dupuy, R., Eppers, O., Hartmann, J., Herber, A., Jäkel, E., Järvinen, E., Jourdan, O., Kästner, U., Kliesch, L.-L., Köllner, F., Mech, M., Mertes, S., Neuber, R., Ruiz-Donoso, E., Schnaiter, M., Schneider, J., Stapf, J., and Zanatta, M.: A comprehensive in situ and remote sensing data set from the Arctic CLoud Observations Using airborne measurements during polar Day (ACLOUD) campaign, Earth Syst. Sci. Data, 11, 1853–1881,, 2019. a

Gatebe, C., King, M., Platnick, S., Arnold, G., Vermote, E., and Schmid, B.: Airborne spectral measurements of surface-atmosphere anisotropy for several surfaces and ecosystems over southern Africa, J. Geophys. Res., 108, 8489,, 2003. a, b

Gatebe, C. K. and King, M. D.: Airborne spectral BRDF of various surface types (ocean, vegetation, snow, desert, wetlands, cloud decks, smoke layers) for remote sensing applications, Remote Sens. Environ., 179, 131–148,, 2016. a

Goyens, C., Marty, S., Leymarie, E., Antoine, D., Babin, M., and Belanger, S.: High Angular Resolution Measurements of the Anisotropy of Reflectance of Sea Ice and Snow, Earth and Space Science, 5, 30–47,, 2018. a

Grenfell, T. C., Warren, S. G., and Mullen, P. C.: Reflection of solar radiation by the Antarctic snow surface at ultraviolet, visible, and near-infrared wavelengths, J. Geophys. Res., 99, 18669–18684,, 1994. a

Han, W.: Remote sensing of surface albedo and cloud properties in the Arctic from AVHRR measurements, PhD thesis, University of Alaska, Fairbanks, available at: (last access: 8 November 2020), 1996. a

Herber, A., Haas, C., Stone, R., Bottenheim, J., Liu, P., Li, S.-M., Staebler, R., Strapp, J., and Dethloff, K.: Regular airborne surveys of Arctic sea ice and atmosphere, Eos T. Am. Geophys. Un., 93, 41–42,, 2012. a

Hu, B., Lucht, W., Li, X., and Strahler, A.: Validation of kernel-driven semiempirical models for the surface bidirectional reflectance distribution function of land surfaces, Remote Sens. Environ., 62, 201–214,, 1997. a

Hudson, S. R. and Warren, S. G.: An explanation for the effect of clouds over snow on the top-of-atmosphere bidirectional reflectance, J. Geophys. Res., 112, D19202,, 2007. a, b, c

Hudson, S. R., Warren, S. G., Brandt, R. E., Grenfell, T. C., and Six, D.: Spectral bidirectional reflectance of Antarctic snow: Measurements and parameterization, J. Geophys. Res., 111, D18106,, 2006. a, b, c

Jacob, D. J., Crawford, J. H., Maring, H., Clarke, A. D., Dibb, J. E., Emmons, L. K., Ferrare, R. A., Hostetler, C. A., Russell, P. B., Singh, H. B., Thompson, A. M., Shaw, G. E., McCauley, E., Pederson, J. R., and Fisher, J. A.: The Arctic Research of the Composition of the Troposphere from Aircraft and Satellites (ARCTAS) mission: design, execution, and first results, Atmos. Chem. Phys., 10, 5191–5212,, 2010. a

Jäkel, E. and Ehrlich, A.: Radiance fields of clouds and the Arctic surface measured by a digital camera during ACLOUD 2017, Leipzig Institute for Meteorology, University of Leipzig, PANGAEA,, 2019. a

Jaross, G. and Warner, J.: Use of Antarctica for validating reflected solar radiation measured by satellite sensors, J. Geophys. Res., 113, D16S34,, 2008. a

Jiao, Z., Ding, A., Kokhanovsky, A., Schaaf, C., Breon, F.-M., Dong, Y., Wang, Z., Liu, Y., Zhang, X., Yin, S., Cui, L., Mei, L., and Chang, Y.: Development of a snow kernel to better model the anisotropic reflectance of pure snow in a kernel-driven BRDF model framework, Remote Sens. Environ., 221, 198–209,, 2019. a

Jin, Z., Charlock, T. P., Yang, P., Xie, Y., and Miller, W.: Snow optical properties for different particle shapes with application to snow grain size retrieval and MODIS/CERES radiance comparison over Antarctica, Remote Sens. Environ., 112, 3563–3581,, 2008. a, b

Kerr, D.: The Canon sRaw and mRaw output formats, Issue 3, avialable at: (last access: 8 November 2020), 2015. a

Klaus, A., Bauer, J., Karner, K., Elbischger, P., Perko, R., and Bischof, H.: Camera calibration from a single night sky image, in: Proceedings of the 2004 IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR 2004), IEEE Computer Society, 151–157,, 2004. a

Klein, A. and Stroeve, J.: Development and validation of a snow albedo algorithm for the MODIS instrument, Ann. Glaciol., 34, 45–52,, 2002. a

Kokhanovsky, A. A. and Zege, E. P.: Scattering optics of snow, Appl. Optics, 43, 1589–1602,, 2004. a

Kokhanovsky, A. A., Aoki, T., Hachikubo, A., Hori, M., and Zege, E. P.: Reflective properties of natural snow: approximate asymptotic theory versus in situ measurements, IEEE T. Geosci. Remote, 43, 1529–1535,, 2005. a

König-Langlo, G.: Radiosonde measurements from Neumayer Station (2014-01), Alfred Wegener Institute, Helmholtz Center for Polar and Marine Research, Bremerhaven, PANGAEA,, 2014. a

Koukal, T. and Atzberger, C.: Potential of multi-angular data derived from a digital aerial frame camera for forest classification, IEEE J. Sel. Top. Appl., 5, 30–43,, 2012. a

Kuchiki, K., Aoki, T., Niwano, M., Motoyoshi, H., and Iwabuchi, H.: Effect of sastrugi on snow bidirectional reflectance and its application to MODIS data, J. Geophys. Res., 116, D18110,, 2011. a

Kuhn, M.: Bidirectional reflectance of polar and alpine snow surfaces, Ann. Glaciol., 6, 164–167,, 1985. a, b

Kuipers Munneke, P., Reijmer, C. H., van den Broeke, M. R., König-Langlo, G., Stammes, P., and Knap, W. H.: Analysis of clear-sky Antarctic snow albedo using observations and radiative transfer modeling, J. Geophys. Res., 113, D17118,, 2008. a

Lebourgeois, V., Bégué, A., Labbé, S., Mallavan, B., Prévot, L., and Roux, B.: Can commercial digital cameras be used as multispectral sensors? A crop monitoring test, Sensors, 8, 7300–7322,, 2008. a, b

Leroux, C. and Fily, M.: Modeling the effect of sastrugi on snow reflectance, J. Geophys. Res., 103, 25779–25788,, 1998. a, b, c

Leroux, C., Deuze, J.-L., Goloub, P., Sergent, C., and Fily, M.: Ground measurements of the polarized bidirectional reflectance of snow in the near-infrared spectral domain: Comparisons with model results, J. Geophys. Res., 103, 19721–19731,, 1998. a

Leroux, C., Lenoble, J., Brogniez, G., Hovenier, J., and De Haan, J.: A model for the bidirectional polarized reflectance of snow, J. Quant. Spectrosc. Ra., 61, 273–285,, 1999. a

Lewis, P.: The utility of kernel-driven BRDF models in global BRDF and albedo studies, Proc. Int. Geosci. Remote Sens. Symp.  1186–1187,, 1995. a

Li, S.: A Model for the Anisotropic Reflectance of Pure Snow, M.Sc. thesis, Department of Geography, University of California, Santa Barbara, 1982. a

Li, X. and Strahler, A.: Geometric-optical bidirectional reflectance modeling of the discrete crown vegetation canopy: effect of crown shape and mutual shadowing, IEEE T. Geosci. Remote, 30, 276–292,, 1992. a

Liang, S., Shuey, C., Russ, A., Fang, H., Chen, M., Walthall, C., Daughtry, C., and Hunt, R.: Narrowband to broadband conversions of land surface albedo: II. Validation, Remote Sens. Environ., 84, 25–41,, 2003. a

Lucht, W. and Lewis, P.: Theoretical noise sensitivity of BRDF and albedo retrieval from the EOS-MODIS and MISR sensors with respect to angular sampling, Int. J. Remote Sens., 21, 81–98,, 2000. a

Lucht, W., Schaaf, C. B., and Strahler, A. H.: An algorithm for the retrieval of albedo from space using semiempirical BRDF models, IEEE T. Geosci. Remote, 38, 977–998,, 2000. a, b, c, d

Lyapustin, A., Gatebe, C. K., Kahn, R., Brandt, R., Redemann, J., Russell, P., King, M. D., Pedersen, C. A., Gerland, S., Poudyal, R., Marshak, A., Wang, Y., Schaaf, C., Hall, D., and Kokhanovsky, A.: Analysis of snow bidirectional reflectance from ARCTAS Spring-2008 Campaign, Atmos. Chem. Phys., 10, 4359–4375,, 2010. a, b

Malinka, A.: Light scattering in porous materials: Geometrical optics and stereological approach, J. Quant. Spectrosc. Ra., 141, 14–23,, 2014. a

Malinka, A., Zege, E., Heygster, G., and Istomina, L.: Reflective properties of white sea ice and snow, The Cryosphere, 10, 2541–2557,, 2016. a

Marks, A., Fragiacomo, C., MacArthur, A., Zibordi, G., Fox, N., and King, M. D.: Characterisation of the HDRF (as a proxy for BRDF) of snow surfaces at Dome C, Antarctica, for the inter-calibration and inter-comparison of satellite optical data, Remote Sens. Environ., 158, 407–416,, 2015. a

Martonchik, J. V., Diner, D. J., Kahn, R. A., Ackerman, T. P., Verstraete, M. E., Pinty, B., and Gordon, H. R.: Techniques for the retrieval of aerosol properties over land and ocean using multiangle imaging, IEEE T. Geosci. Remote, 36, 1212–1227,, 1998. a, b

Matsui, H., Kondo, Y., Moteki, N., Takegawa, N., Sahu, L. K., Zhao, Y., Fuelberg, H. E., Sessions, W. R., Diskin, G., Blake, D. R., Wisthaler, A., and Koike, M.: Seasonal variation of the transport of black carbon aerosol from the Asian continent to the Arctic during the ARCTAS aircraft campaign, J. Geophys. Res., 116, D05202,, 2011. a

Mayer, B. and Kylling, A.: Technical note: The libRadtran software package for radiative transfer calculations - description and examples of use, Atmos. Chem. Phys., 5, 1855–1877,, 2005. a, b

Mei, L., Xue, Y., de Leeuw, G., von Hoyningen-Huene, W., Kokhanovsky, A., Istomina, L., Guang, J., and Burrows, J.: Aerosol optical depth retrieval in the Arctic region using MODIS data over snow, Remote Sens. Environ., 128, 234–245,, 2013. a

Mishchenko, M., Dlugach, J., Yanovitskij, E., and Zakharova, N.: Bidirectional reflectance of flat, optically thick particulate layers: an efficient radiative transfer solution and applications to snow and soil surfaces, J. Quant. Spectrosc. Ra., 63, 409–432,, 1999. a

Mori, Y., Yamashita, A., Tanaka, M., Kataoka, R., Miyoshi, Y., Kaneko, T., Okutomi, M., and Asama, H.: Calibration of fish-eye stereo camera for aurora observation, in: Proc. of the International Workshop on Advanced Image Technology (IWAIT2013), Nagoya, Japan, 729–734, 2013. a

Nolin, A., Fetterer, F., and Scambos, T.: Surface roughness characterizations of sea ice and ice sheets: Case studies with MISR data, IEEE T. Geosci. Remote, 40, 1605–1615,, 2002. a, b

Painter, T. and Dozier, J.: Measurements of the hemispherical-directional reflectance of snow at fine spectral and angular resolution, J. Geophys. Res., 109,, 2004a. a

Painter, T. and Dozier, J.: The effect of anisotropic reflectance on imaging spectroscopy of snow properties, Remote Sens. Environ., 89, 409–422,, 2004b. a

Painter, T., Paden, B., and Dozier, J.: Automated spectro-goniometer: A spherical robot for the field measurement of the directional reflectance of snow, Rev. Sci. Instrum., 74, 5179–5188,, 2003. a, b

Pegrum, H., Fox, N., Chapman, M., and Milton, E.: Design and Testing a New Instrument to Measure the Angular Reflectance of Terrestrial Surfaces, in: IEEE International Geoscience & Remote Sensing Symposium, IGARSS 2006, 31 July–4 August 2006, Denver, Colorado, USA, Proceedings, 1119–1122,, 2006. a, b

Phong, B.: Illumination for Computer Generated Pictures, Commun. ACM, 18, 311–317,, 1975. a

Pohl, C., Istomina, L., Tietsche, S., Jäkel, E., Stapf, J., Spreen, G., and Heygster, G.: Broadband albedo of Arctic sea ice from MERIS optical data, The Cryosphere, 14, 165–182,, 2020. a

Qu, Y., Liang, S., Liu, Q., He, T., Liu, S., and Li, X.: Mapping Surface Broadband Albedo from Satellite Observations: A Review of Literatures on Algorithms and Products, Remote Sens. Basel, 7, 990–1020,, 2015. a, b

Ross, J.: The Radiation Regime and Architecture of Plant Stands, Tasks for Vegetation Sciences 3, Dr. W. Junk Publishers, The Hague, the Netherlands, 1981. a

Roujean, J.-L., Leroy, M., and Deschamps, P.-Y.: A bidirectional reflectance model of the Earth's surface for the correction of remote sensing data, J. Geophys. Res., 97, 20455–20468,, 1992. a, b

Schaaf, C. and Wang, Z.: MCD43A1 MODIS/Terra+Aqua BRDF/Albedo Model Parameters Daily L3 Global – 500m V006 [Data set], NASA EOSDIS Land Processes DAAC,, 2015a. a

Schaaf, C. and Wang, Z.: MCD43A2 MODIS/Terra+Aqua BRDF/Albedo Quality Daily L3 Global – 500m V006 [Data set], NASA EOSDIS Land Processes DAAC,, 2015b. a

Schaaf, C. B., Gao, F., Strahler, A. H., Lucht, W., Li, X. W., Tsang, T., Strugnell, N. C., Zhang, X. Y., Jin, Y. F., Muller, J. P., Lewis, P., Barnsley, M., Hobson, P., Disney, M., Roberts, G., Dunderdale, M., Doll, C., d'Entremont, R. P., Hu, B. X., Liang, S. L., Privette, J. L., and Roy, D.: First operational BRDF, albedo nadir reflectance products from MODIS, Remote Sens. Environ., 83, 135–148,, 2002. a, b, c

Schaepman-Strub, G., Schaepman, M. E., Painter, T. H., Dangel, S., and Martonchik, J. V.: Reflectance quantities in optical remote sensing-definitions and case studies, Remote Sens. Environ., 103, 27–42,, 2006. a, b, c, d

Schmid, H.: Stellar calibration of the orbigon lens, Photogramm. Eng., 40, 101–115, 1974. a

Stamnes, K., Tsay, S., Wiscombe, W., and Jayaweera, K.: A numerically stable algorithm for discrete-ordinate-method radiative transfer in multiple scattering and emitting layered media, Appl. Optics, 27, 2502–2509,, 1988. a, b

Stanton, B., Miller, D., Adams, E., and Shaw, J. A.: Bidirectional-reflectance measurements for various snow crystal morphologies, Cold Reg. Sci. Technol., 124, 110–117,, 2016. a

Steffen, K.: Bidirectional reflectance of snow at 500–600 nm, Large Scale Effects of Seasonal Snow Cover, IAHS Publ., 166, 415–425, 1987. a, b, c, d

Strahler, A., Wanner, W., Schaaf, C., Li, X., Hu, B., Muller, J., Lewis, P., and Barnsley, M.: MODIS BRDF/Albedo product: Algorithm Theoretical Basis Document Version 4.0, 1996. a

Stroeve, J., Box, J., Gao, F., Liang, S., Nolin, A., and Schaaf, C.: Accuracy assessment of the MODIS 16-day albedo product for snow: comparisons with Greenland in situ measurements, Remote Sens. Environ., 94, 46–60,, 2005. a, b, c, d

Stroeve, J., Box, J., Wang, Z., Schaaf, C., and Barrett, A.: Re-evaluation of MODIS MCD43 Greenland albedo accuracy and trends, Remote Sens. Environ., 138, 199–214,, 2013. a

Tomasi, C., Vitale, V., Lupi, A., Di Carmine, C., Campanelli, M., Herber, A., Treffeisen, R., Stone, R. S., Andrews, E., Sharma, S., Radionov, V., von Hoyningen-Huene, W., Stebel, K., Hansen, G. H., Myhre, C. L., Wehrli, C., Aaltonen, V., Lihavainen, H., Virkkula, A., Hillamo, R., Strom, J., Toledano, C., Cachorro, V. E., Ortiz, P., de Frutos, A. M., Blindheim, S., Frioud, M., Gausa, M., Zielinski, T., Petelski, T., and Yamanouchi, T.: Aerosols in polar regions: A historical overview based on optical depth and in situ observations, J. Geophys. Res., 112, D16205,, 2007. a

Tomasi, C., Kokhanovsky, A., Lupi, A., Ritter, C., Smirnov, A., O'Neill, N., Stone, R., Holben, B., Nyeki, S., Wehrli, C., Stohl, A., Mazzola, M., Lanconelli, C., Vitale, V., Stebel, K., Aaltonen, V., de Leeuw, G., Rodriguez, E., Herber, A., Radionov, V., Zielinski, T., Petelski, T., Sakerin, S., Kabanov, D., Xue, Y., Mei, L., Istomina, L., Wagener, R., McArthur, B., Sobolewski, P., Kivi, R., Courcoux, Y., Larouche, P., Broccardo, S., and Piketh, S.: Aerosol remote sensing in polar regions, Earth-Sci. Rev., 140, 108–157,, 2015. a

Tsai, R.: A versatile camera calibration technique for high-accuracy 3D machine vision metrology using off-the-shelf TV cameras and lenses, IEEE J. Robotic. Autom., 3, 323–344,, 1987. a

Urquhart, B., Kurtz, B., and Kleissl, J.: Sky camera geometric calibration using solar observations, Atmos. Meas. Tech., 9, 4279–4294,, 2016. a, b

Vermote, E. and Kotchenova, S.: Atmospheric correction for the monitoring of land surfaces, J. Geophys. Res., 113, D23S90,, 2008. a

Walthall, C., Norman, J., Welles, J., Campbell, G., and Blad, B.: Simple equation to approximate the bidirectional reflectance from vegetative canopies and bare soil surfaces, Appl. Optics, 24, 383–387,, 1985. a

Warren, S., Brandt, R., and O'Rawe Hinton, P.: Effect of surface roughness on bidirectional reflectance of Antarctic snow, J. Geophys. Res., 103, 25789–25807,, 1998. a, b, c, d, e

Wendisch, M., Müller, D., Schell, D., and Heintzenberg, J.: An airborne spectral albedometer with active horizontal stabilization, J. Atmos. Ocean. Tech., 18, 1856–1866,<1856:AASAWA>2.0.CO;2, 2001. a

Wendisch, M., Pilewskie, P., Jäkel, E., Schmidt, S., Pommier, J., Howard, S., Jonsson, H. H., Guan, H., Schröder, M., and Mayer, B.: Airborne measurements of areal spectral surface albedo over different sea and land surfaces, J. Geophys. Res., 109, D08203,, 2004.  a

Wessel, P., Luis, J. F., Uieda, L., Scharroo, R., Wobbe, F., Smith, W. H. F., and Tian, D.: The Generic Mapping Tools version 6, Geochem. Geophy. Geosy., 20, 5556–5564,, 2019. a

Wendisch, M., Macke, A., Ehrlich, A., Lüpkes, C., Mech, M., Chechin, D., Dethloff, K., Velasco, C. B., Bozem, H., Brückner, M., Clemen, H.-C., Crewell, S., Donth, T., Dupuy, R., Ebell, K., Egerer, U., Engelmann, R., Engler, C., Eppers, O., Gehrmann, M., Gong, X., Gottschalk, M., Gourbeyre, C., Griesche, H., Hartmann, J., Hartmann, M., Heinold, B., Herber, A., Herrmann, H., Heygster, G., Hoor, P., Jafariserajehlou, S., Jäkel, E., Järvinen, E., Jourdan, O., Kästner, U., Kecorius, S., Knudsen, E. M., Köllner, F., Kretzschmar, J., Lelli, L., Leroy, D., Maturilli, M., Mei, L., Mertes, S., Mioche, G., Neuber, R., Nicolaus, M., Nomokonova, T., Notholt, J., Palm, M., van Pinxteren, M., Quaas, J., Richter, P., Ruiz-Donoso, E., Schäfer, M., Schmieder, K., Schnaiter, M., Schneider, J., Schwarzenböck, A., Seifert, P., Shupe, M. D., Siebert, H., Spreen, G., Stapf, J., Stratmann, F., Vogl, T., Welti, A., Wex, H., Wiedensohler, A., Zanatta, M., and Zeppenfeld, S.: The Arctic Cloud Puzzle: Using ACLOUD/PASCAL Multiplatform Observations to Unravel the Role of Clouds and Aerosol Particles in Arctic Amplification, B. Am. Meteorol. Soc., 100, 841–871,, 2019. a

Wesche, C., Steinhage, D., and Nixdorf, U.: Polar aircraft Polar5 and Polar6 operated by the Alfred Wegener Institute, J. Large-Scale Res. Facilities, 2, A87,, 2016. a

Xie, Y., Yang, P., Gao, B. C., Kattawar, G. W., and Mishchenko, M. I.: Effect of ice crystal shape and effective size on snow bidirectional reflectance, J. Quant. Spectrosc. Ra., 100, 457–469,, 2006. a

Yang, L., Xue, Y., Guang, J., Kazemian, H., Zhang, J., and Li, C.: Improved Aerosol Optical Depth and Angstrom Exponent Retrieval Over Land From MODIS Based on the Non-Lambertian Forward Model, IEEE Geosci. Remote S., 11, 1629–1633,, 2014. a

Zege, E. P., Katsev, I. L., Malinka, A. V., Prikhach, A. S., Heygster, G., and Wiebe, H.: Algorithm for retrieval of the effective snow grain size and pollution amount from satellite measurements, Remote Sens. Environ., 115, 2674–2685,, 2011. a

Zhuravleva, T. B. and Kokhanovsky, A. A.: Influence of surface roughness on the reflective properties of snow, J. Quant. Spectrosc. Ra., 112, 1353–1368,, 2011. a, b

Short summary
The angular reflection of solar radiation by snow surfaces is particularly anisotropic and highly variable. We measured the angular reflection from an aircraft using a digital camera in Antarctica in 2013/14 and studied its variability: the anisotropy increases with a lower Sun but decreases for rougher surfaces and larger snow grains. The applied methodology allows for a direct comparison with satellite observations, which generally underestimated the anisotropy measured within this study.