A variational technique to estimate snowfall rate from coincident radar, snowflake, and fallspeed observations

Estimates of snowfall rate as derived from radar reflectivities alone are non-unique. Different combinations of 10 snowflake microphysical properties and particle fallspeeds can conspire to produce nearly identical snowfall rates for given radar reflectivity signatures. Such ambiguities can result in retrieval uncertainties on the order of 100-200% for individual events. Here, we use observations of particle size distribution (PSD), fallspeed, and snowflake habit from the Multi-Angle Snow Camera (MASC) to constrain estimates of snowfall derived from Ka-band Zenith Radar (KAZR) measurements at the ARM NSA Barrow Climate Facility site. MASC measurements of microphysical properties with uncertainties are introduced 15 into a modified form of the optimal-estimation CloudSat snowfall algorithm (2C-SNOW-PROFILE) via the a priori guess and variance terms. Use of MASC fallspeed, MASC PSD, and CloudSat snow particle model as base assumptions resulted in retrieved total accumulations with a -18% difference relative to nearby National Weather Service observations over five snow events. Use of different but reasonable combinations of retrieval assumptions resulted in estimated snowfall accumulations with differences ranging from -64% to + 94% for the same storm events. Retrieved snowfall rates were 20 particularly sensitive to assumed fallspeed and habit, suggesting that in-situ measurements can help to constrain key snowfall retrieval uncertainties. More accurate knowledge of these properties dependent upon location and meteorological conditions should help refine and improve ground and space-based radar estimates of snowfall.


Introduction
The high-latitude regions play a critical role in shaping climate response to anthropogenic forcing.Model predictions suggest that it is these areas that are most susceptible to change and will experience the most dramatic temperature increase in response to the release of greenhouse gases into the atmosphere (Manabe and Stouffer, 1980;Holland and Bitz, 2003;Serreze and Francis, 2006).Observed changes in the Arctic over the late 20th century and early 21st century have been dramatic and have included increased surface temperature and decreased sea ice, permafrost, glacial ice sheet, and spring Arctic snow cover extents (Serreze et al., 2000;Frauenfeld et al., 2004;Dyurgerov and Meier, 2005;Stroeve et al., 2008Stroeve et al., , 2014;;Brown et al., 2010;Perlwitz et al., 2015, among numerous others).
Snowfall can dramatically change surface conditions at high-latitude regions, acting to increase shortwave surface albedo while impacting sensible and latent heat fluxes and longwave emission (Cohen and Rind, 1991).Changes in snow cover can feed back on snow cover, sea ice, and permafrost distributions (Brown, 2000;Ramonovsky et al., 2002;Holland et al., 2006;Vavrus, 2007), where the effects on permafrost in turn may affect high-latitude carbon storage.Snowfall changes at higher north latitudes may also impact freshening of the North Atlantic Ocean and the strength of the Atlantic Meridional Overturning Circulation.Peterson et al. (2006) attributed more than half of the cumulative freshwater input anomaly in the Arctic and high-latitude North Atlantic oceans over the previous 50 years to increases in net precipitation over land and oceans, exceeding the estimated contributions from glacial melt and sea ice reduction.Snowfall is also important as it provides mass influx for vast glacial ice sheets such as those found in Antarctica and Greenland (Lenaerts et al., 2013;Palerme et al., 2014).The distribution of snowfall helps define ice sheet dynamics, and the relative difference between snowfall and melt will impact the long-term survival of these glaciers (Van Tricht et al., 2016).Significant loss of ice sheets may have deleterious effects on human society through a corresponding rise of global ocean levels (Gardner et al., 2013;Jacob et al., 2012).
The quantitative estimation of snowfall at the global scale from spaceborne measurements has occurred only recently.Initial retrieval approaches were based on passive microwave measurements (Skofronick-Jackson et al., 2004;Noh et al., 2006) with a shift in emphasis to radar observations with the launch of the CloudSat Cloud Profiling Radar (CPR) in 2006.Matrosov et al. (2007) and Liu (2008a) demonstrated the first-order capability of the CloudSat CPR to retrieve vertical profiles of dry snowfall.Kulie and Bennartz (2009), in turn, used CloudSat CPR reflectivities to estimate global dry snowfall rate for a year of CloudSat data.They found that snowfall estimates depended critically upon assumed relationships between radar reflectivity and snowfall particle size distributions (PSDs) and shapes.
These studies suggest that estimates of snowfall as derived from radar reflectivities alone are non-unique.Numerous different combinations of snowflake microphysical properties and snow particle fall speeds may yield nearly identical surface snowfall rates for a given reflectivity profile.As such, use of traditional Z e -S relationships to quantify snowfall cannot be expected to produce accurate results for all snowfall events.
From an operational retrieval perspective, clever selection of representative snowflake microphysical properties may produce estimates of snowfall amounts that agree well with reported values for climate or regional applications.For example, Wood et al. (2015; but see also Wood, 2011) developed snowflake models for the CloudSat snowfall retrieval algorithm based upon field campaigns focused on coldseason clouds and precipitation, principally the Canadian CloudSat-CALIPSO Validation Project (C3VP, Hudak et al., 2006).They used this a priori knowledge of snowfall microphysics to refine expected snowfall-radar reflectivity relationships for the optimal-estimation-based (Rodgers, 2000) CloudSat snowfall retrieval scheme used for the CloudSat 2C-SNOW-PROFILE (Wood et al., 2013) product.The product provides estimates of snow size distribution and bulk properties (water content and snowfall rate) at the surface and aloft over land, ice, and ocean surfaces.These estimates have proven valuable for assessing snowfall budgets in remote regions (Palerme et al., 2014;Kulie et al., 2016), providing data for testing global climate models (Palerme et al., 2016;Christensen et al., 2016), and evaluating the performance of ground-based radar measurements of snowfall (Norin et al., 2015).The Global Precipitation Mission con-tinues the global monitoring of snowfall with dual-frequency Ka-band and Ku-band radar and extensive ground measurement activities (Skofronick-Jackson et al., 2015).
Despite such efforts, validation activities suggest that uncertainties in retrieved snowfall rates can still be on the order of 200 % for individual snow events (Wood, 2011).In our retrieval scheme, instead of using a priori guesses of snowflake microphysical properties from field campaigns such as C3VP, we used coincident observations of snowflake microphysical properties from the Multi-Angle Snowflake Camera (MASC; Garrett et al., 2012) to constrain the radarbased retrieval approach.The MASC takes multiple images of snowflakes in free fall while simultaneously measuring fall speed.From these pictures, estimates of snowflake maximum particle dimension, habit, and other properties that could be used to refine a radar retrieval scheme are estimated.It should be noted, however, that any instrument that measures snow particle properties or fall speed, e.g., the Precipitation Imaging Package (Newman et al., 2009), could replace the role of the MASC in the variational approach presented in this work.
We describe here a retrieval technique for snowfall rate and its application to five snow events as observed at the Atmospheric Radiation Measurement (ARM) North Slope Alaska (NSA) Climate Research Facility site at Barrow.Our scheme uses ground-based Ka-band ARM zenith radar (KAZR) measurements for reflectivity profiles and was directly modified from the W-band CloudSat 2C-SNOW-PROFILE algorithm.The flexible optimal-estimation framework was used to incorporate coincident MASC observations into the radar-based scheme through the retrieval a priori terms.A primary objective was to identify those combinations of retrieval assumptions that allowed the best match with snowfall observations at the Barrow site for five snow events during April and May 2014.Another objective was to quantify the sensitivity of snowfall rate retrieval results to the key microphysical assumptions of particle size distribution, habit, and fall speed from the MASC.These results were contrasted with snowfall rates found using alternate assumptions such as the Locatelli and Hobbs (1974, hereafter LH74) particle dimension-fall speed parameterizations, observed KAZR Doppler velocities, and previous field campaign estimates of snow particle size distributions.In Sect.2, we discuss the methodology including the CloudSat snowfall retrieval scheme, the MASC observations, and the combined radar-MASC retrieval approach.In Sect.3, we explore the application of the combined scheme to five snowfall events observed at the Barrow site.In Sect.4, we discuss the implications of our results for the utility of ground-based in situ measurements to refine and improve radar retrievals of snowfall.

CloudSat snowfall retrieval scheme
The CloudSat snowfall retrieval algorithm uses profiles of Wband 94 GHz radar reflectivity to estimate vertical profiles of snow properties.The flexible optimal-estimation approach is used to combine measurements and a priori microphysical information into a common retrieval cost function.Specifically, the scheme assumes an exponential form for snow particle size distributions for each radar reflectivity bin as in Eq. (1), where λ is PSD slope, N 0 is its intercept, and D is particle maximum dimension (Heymsfield et al., 2008).For our implementation, the PSD slope parameter is allowed to vary with height but the number density is held constant given the limited number of independent observations from the radar.(Note that the uniform N 0 used in this work represents a divergence from the current CloudSat algorithm in which number density is allowed to vary with height.)A priori assumptions of particle mass-diameter and diameter-fall speed relationships for D in Eq. ( 1) allow the determination of snow water content and snowfall rate for each radar range bin.
For the CloudSat algorithm, the required a priori microphysical and scattering properties were determined from analyses of snow observations from field campaigns focused on cold-season clouds and precipitation, principally the Canadian CloudSat-CALIPSO Validation Project (Hudak et al., 2006).Snowflake particle models were constructed based upon observed mass and horizontally projected area as a function of particle size (Wood et al., 2015).These efforts were somewhat unique in that they explicitly use ground-based observations to refine forward-model assumptions for the retrieval, an idea we will return to for our combined radar-MASC scheme as described below.Scattering properties for the snowflake particle models were developed through use of the discrete dipole approximation (DDA; Draine and Flatau, 1994) method.DDA replaces the true par- ticle shape with an approximate shape constructed as a threedimensional array of small, cubic ice dipoles.
Figure 1 shows the particle models developed for the CloudSat algorithm.These particle models, intended to simulate the coarse features of snow particles, consist of solidice dipoles intermixed with empty (i.e., air-occupied) dipole locations to meet observed mass and horizontally projected area (Ap) constraints.Analyses of disdrometer, X-band radar, 2-D video disdrometer, and Snowflake Video Imager observations for four significant C3VP snow events provided best-estimate mass-dimension and Ap-dimension relationships (Wood et al., 2015).Scattering calculations for a range of shapes conforming to these constraints suggested that backscatter cross sections can vary by 2 orders of magnitude between these models for a given particle size typical of snowfall.
For the C3VP snow events, these particle models were evaluated by testing against coincident, observed W-band radar reflectivities (Wood et al., 2015).DDA scattering properties were used in combination with video disdrometer observations of snow particle size distributions (the Snowflake Video Imager, Newman et al., 2009) to calculate synthetic W-band reflectivities.These synthetic reflectivities were then plotted versus observed reflectivities as in Fig. 2. Minimization of these differences suggested that the B8pr-30 model best matched observations for the C3VP snow events.
For our study, we will use this B8pr-30 particle model as a base assumption for our Barrow storm events.Given the nonuniqueness of the snowfall retrieval problem, however, it is possible that sources of error other than particle model may have compensated to allow observed and synthetic reflectivities to match as in Fig. 2. We will examine the sensitivity of the retrieved snowfall rate to MASC observations of shape, particle dimension, and fall speed to better understand this issue as discussed in Sect.2.3 below.

Multi-Angle Snowflake Camera
The Multi-Angle Snowflake Camera (MASC) takes highresolution multi-angle photographs of snowflakes as they settle near the surface.It simultaneously measures snowflake fall speed.The MASC consists of three cameras each pointing at an identical focal point approximately 10 cm away in a ring opening.This ring houses a system of near-infrared emitter-detector pairs, arranged in two arrays that are separated vertically by about 3.2 cm.If a hydrometeor passes through both upper and lower arrays, the MASC will trigger each of the three cameras and flash a bank of lights aimed at the center of the camera depth of field.Fall speed is calculated from the time it takes to traverse the distance between the upper and lower triggering array.
Given images from the MASC, it is possible to derive estimates of properties such as particle shape, aspect ratio, maximum dimension, complexity, and orientation (Garrett et al., 2012(Garrett et al., , 2015;;Garrett and Yuter, 2014).Here, specifically, we use observations of maximum dimension and habit to constrain uncertainties inherent to radar-based retrievals of snow rate.For example, Fig. 3 shows typical snowflakes observed at the NSA Barrow Climate Research Facility site during the Spring 2014 MASC deployment.Observed habits include graupels, columns, plates, and aggregate combinations of each.Such observations of shape can be used to select the most appropriate particle model for a given precipitation scene or to identify those scenes in which none of the currently available particle models will be expected to produce a good fit.
Observations of fall speed from the MASC, in turn, can be used to translate retrieved cloud snow water content to a snowfall rate.In this work, we contrast retrieved snowfall rates found using MASC observations, near-surface Doppler measurements, and LH74 particle dimension-fall speed relationships.

Combined radar-MASC retrieval
The MASC deployment to the ARM NSA Barrow Climate Research Facility site in Spring 2014 provided an ideal opportunity to employ our combined radar-MASC snowfall retrieval scheme.In this section, we present sample observations from a Barrow snowstorm to illustrate how we merge MASC information into our radar-based retrieval scheme.We focus specifically on MASC observations of particle maximum dimension, particle model, and fall speed.Given the temporary nature of the deployment at Barrow, the MASC was left unshielded during the storm events presented in this work.
Given the lack of W-band radar measurements during this time period, we use Ka-band ARM zenith radar general mode observations for this study.The KAZR is a vertically pointing radar that measures reflectivity, vertical velocity, and spectral width at 34.8 GHz with a resolution of 30 m from near surface up to about 20 km (Bharadway et al., 2013).The radar had a 2 m antennae with 0.31 • beam width.KAZR systems have a minimum detectable signal in general mode near −20 to −25 dBZ dependent upon target range (Feng et al., 2014;Chandra et al., 2015).Such performance is sufficient for these studies given we need reflectivities near −10 dBZ during the entire day to generate measurable snowfall (0.01 inches) as reported by the National Weather Service (NWS). Figure 4 shows KAZR reflectivities for an all-day 23 April snow event at Barrow.The MASC images shown below the reflectivity plot indicate typical observed snowflakes selected as a function of time.For example, these images suggest rimed graupel-like particles around 04:00 UTC when the KAZR reflectivities suggest the most cloud vertical development.The snowflakes transitioned to sector-plate-type crystals and aggregates as the cloud tops lowered and became more homogenous.
To understand precisely how we incorporate desired MASC information into the modified KAZR snowfall scheme, it is useful to consider the optimal-estimation ap-proach used in the CloudSat retrieval scheme.Letting x denote the vector of snowfall properties to be retrieved, the optimal-estimation technique consists of minimizing a combination of the variance between the set of observations, y, and a corresponding set of simulated measurements, F (x), and between x and a suitable a priori guess, a.Assuming Gaussian statistics, this is accomplished by minimizing the scalar cost function, with respect to x. F denotes the physical model relating snowfall parameters to the radar observations and is called the "forward model", S a is the a priori error covariance matrix, and S y is the measurement error covariance matrix.S y represents not only random instrument noise but also the impact of uncertainties in forward-model assumptions on simulated measurements, F (x).We use an S y standard deviation value of 2 dBZ for this study for the diagonal matrix elements based upon Hammonds et al. (2014), who quantified the uncertainty in forward-modeled radar reflectivity due to assumptions such as mass-dimensional relationships at Kaband frequency.The optimal-estimation approach weights the magnitude of the error covariances to determine the relative impact of both a priori guesses and observations on the final retrieval estimates.The observation vector, y, is the vertical profile of KAZR reflectivities where the retrieval vector, x, would be the vertically varying snowfall particle PSD slope parameter and vertically uniform number density.The values of x that minimize Eq. ( 2) are found through Newtonian iteration.
We introduce MASC PSD information into the retrieval scheme through use of the a priori estimate, a, and a priori covariance matrix, S a .For the PSD slope parameter a priori guess, the MASC images were processed to quantify maximum particle dimension for each snowflake according to the techniques developed in Garrett et al. (2012), Garrett andYuter (2014), andGarrett et al. (2015).We then fit a slope parameter for an assumed exponential particle size distribution to the tail of the size distribution > 1 mm in size.We derived estimates of the slope parameter either for the entire storm event or for subsections of the event based upon differences in storm morphology.For example, for the 23 April event, we found separate estimates of the slope parameter corresponding to the vertically developed part of the storm before 08:00 UTC (1.11 mm −1 ) and for the stratiform part of the storm after 08:00 UTC (0.74 mm −1 ).The use of such lengthy time periods was necessary given the infrequent sampling of snowflakes by the MASC at the Barrow site, typically from 1 to 10 snowflakes per minute for these events.
Values for uncertainties in the slope parameter as defined in the a priori error covariance matrix were calculated for the individual storm subsections.These values depended upon Table 1.MASC-observed slope parameters for exponential particle size distributions for the Barrow snow events.The C3VP observations were derived from field observations with similar snow rates as the Barrow storms (Wood et al., 2014) both the measure of fit of the exponential slope parameter to MASC observations and expected uncertainties in MASCderived estimates of particle size from snowflake images.Uncertainties in the slope fit, σ FIT , were defined as regression 95 % confidence values.Uncertainties in MASC estimates of particle size, σ SIZE , were assumed to be 15 % of the median maximum particle dimension based upon the preliminary work of Kleinkort et al. (2016).That study used 3D printed synthetic snowflakes with known geometry to evaluate MASC estimates of particle dimensions using both a three-camera and five-camera system.For our study, the variance of the slope parameter guess, σ 2 a , is then defined as Table 1 lists estimated slope parameters with uncertainties for the five storm events presented in Sect.3. The relatively small uncertainties in the MASC slope parameter observations dictate a solution close to the a priori guess.As a consequence, the uncertainty in the number density was then allowed to vary over several orders of magnitude to ensure retrieval convergence, i.e., allow the particle number to vary so that the forward-model-simulated reflectivities can match observed radar reflectivities regardless of particle size.The number density term for this scheme is assumed constant with height.Although this approach may be a problem for space-based radars that must see through the entire storm to estimate surface snowfall rate, it is more applicable for use with the upward-looking KAZR.Retrieved snowfall rates are based only on the lowest non-noise radar range bins from 160-310 m above the surface.Given a retrieved slope parameter and number density, snow water content was estimated for each radar range bin assuming a knowledge of snow particle mass-dimension relationships.As in Fig. 4, MASC images are used to select the most appropriate particle model with associated scattering properties for a given scene.The development of more particle models, e.g., lightly rimed aggregates or graupels to better represent the natural variability of observed snowflakes, is a current line of research.In terms of the optimal-estimation technique, this will reduce a major source of uncertainty in the forward model (S y term in Eq. 2).Observations of particle fall speeds from the MASC were used to translate retrieved snow water contents into surface snowfall rates.Figure 5 shows MASC-observed snowflake fall speeds plotted as hourly averages for the 23 April snow event.This plot shows that fall speeds from the MASC observations were almost always lower than those measured from the KAZR Doppler and those predicted from LH74 parameterization schemes (graupel-like snow and aggregates of dendrites).Discrepancies in fall speed between the MASC and KAZR Doppler observations were expected.Both of these measurements depend upon still air-particle fall speed relationships, but they also depend upon atmospheric vertical motion which is expected to vary between the near-surface MASC and aloft Doppler measurements.The MASC observations would also be influenced by near-surface turbulence and disruptions to airflow as the snowflakes pass through the sampling ring and infrared sensors.Indeed, the hourly MASC and KAZR Doppler fall speeds had a Pearson correlation coefficient of −0.34.We note that these discrepancies could have been exacerbated, in part, from the unshielded status of the MASC.
The LH74 fall speeds (m s −1 ) were estimated through Eqs. ( 4) and ( 5) as below: V = 1.1D 0.28 for graupel, (4) where the snowflake maximum dimensions (mm) used for the LH74 parameterization schemes were calculated using MASC observations.The discrepancies in Fig. 5 are in general agreement with Garrett and Yuter (2014), who noticed differences between MASC observations of fall speed at the Alta Ski Resort in Utah and LH74 parameterizations.The LH74 schemes predicted fall speeds poorly for lowtemperature and highly turbulent environments that might be expected during high-latitude Arctic snow events.However, the hourly MASC observations still had a higher Pearson correlation coefficient with the LH74 schemes (0.66) than with the Doppler fall speeds (−0.34).
3 Snowfall retrieval results

Snow events and retrieval assumptions
The KAZR-MASC retrieval approach was applied to five snow events as observed at the Barrow site during Spring 2014.KAZR reflectivities for these storms are shown in Figs. 4 and 6.These events were selected as they produced measurable snowfall at the nearby Barrow National Weather Service site and triggered coincident MASC snowflake images.Retrieved estimates of total snowfall accumulation were compared with NWS snowfall observations to evaluate retrieval performance.
We examined the impact of snowflake habit, slope parameter (λ), and fall-speed assumptions on retrieved snowfall liquid water equivalent (LWE) for these events.In terms of habit, we selected the CloudSat B8pr-30 particle model, sector plates, and hexagonal columns.The CloudSat particle model not only performed well in CloudSat snowfall validation studies but also was visually consistent with MASC images of the general shape of snowflakes seen at Barrow.Sector plates and hexagonal columns were chosen as they too were observed during these snow events.Scattering properties for these pristine habits were derived following the constrained discrete dipole modeling method described by Wood et al. (2015) with particle dimension relationships from Auer and Veal (1970) and particle mass constraints from Mitchell (1996).
For PSD assumptions, we used MASC estimates of the slope parameter with uncertainties as listed in Table 1 as introduced through the a priori guess and covariance terms in Eq. ( 2).We also employed an a priori PSD λ assumption (2.8 mm −1 ) as derived from snowflake observations during the C3VP snow measurement field campaign (Wood et al., 2014) for snow events with similar snowfall rates as the Barrow events.For fall-speed assumptions, we used MASC and We performed retrievals with permutations of habit, PSD, and fall-speed assumptions to determine which combinations produce the best match with snowfall rates observed at the Barrow NWS site.Our best guess or base assumptions were the CloudSat B8pr-30 particle model, MASC fall speed, and MASC PSD λ.Given the non-unique nature of snowfall retrievals and the difference in location between radar and the NWS observations, we do not pretend that agreement between retrieval results and NWS snowfall observations validates our retrieval assumptions.Furthermore, the accurate measurement of snowfall amount from precipitation gauges is a challenging research topic in its own right (Goodison et al., 1998;Rasmussen et al., 2012).Numerous studies over the past century have sought to quantify uncertainties in these measurements and then provide a correction factor for environmental conditions (e.g., Black, 1954;Larson and Peck, 1974;Yang et al., 1995Yang et al., , 2005)).During windy and turbulent conditions, for example, the gauges will likely underestimate snowfall amount as the snowflakes cannot settle into the gauge opening.Overall, collection efficiency is a complex function of instrument design, shielding, wind, turbulence, temperature, and topography among other factors (Folland, 1988;Goodison et al., 1998;Rasmussen et al., 2012;The-riault et al., 2015) and can result in snowfall uncertainties on the order of 25-50 % over the course of a winter season (Wolff et al., 2015).As such, we understand that the reported NWS snowfall values cannot be taken at absolute face value.However, the numerical experiments presented in this work do quantify the sensitivity of snowfall retrieval results to snow microphysical parameters that can be observed by in situ instrumentation.A better understanding of the impact of these parameters, in turn, should provide a platform in which we can then examine other sources of error in the snowfall problem.

23 April snow event
Radar reflectivities and MASC images for a 23 April snow event that produced a total snow accumulation of 4.57 mm (0.18 inches) liquid equivalent are shown in Fig. 4. Retrieval results assuming different particle model, slope parameter, and fall-speed combinations as described in Sect.3.1 are presented in Table 2.The percentage difference term in the table is defined through Eq. ( 6): Error values ranged from −51 to +132 % for this snow event.We found the best agreement with the NWS observations using our "base" assumptions of MASC fall speed and PSD λ, CloudSat B8pr-30 particle model with a retrieved accumulated snowfall of 4.88 mm liquid equivalent, which is a difference of +7 %.This accuracy was closely matched with the use of the combination of hexagonal columns, KAZR Doppler fall speeds, and the PSD λ from the MASC.For these two scenarios, the highly reflective per unit mass hexagonal columns and high-fall-speed Doppler observations offset the lower reflective per unit mass CloudSat particle model and low-fall-speed MASC observations to produce nearly equivalent snow accumulations.Since hexagonal columns were not observed during the snowstorm, however, that solution is not valid.The worst agreement (+132 %) arose from the use of the MASC PSD λ, the CloudSat B8pr-30 particle model, and the LH74 graupel fall speed.This large discrepancy was driven by the high fall speed of the LH74 graupel parameterization as seen in Fig. 5.These calculations allowed us to quantify the sensitivity of retrieval results to assumptions of PSD λ, fall speed, and crystal habit for a given storm event.The assumed PSD slope parameter had the least impact on variability in estimated snowfall accumulation.Snowfall totals of 4.88 and 6.07 mm were found when using the MASC PSD λ and the C3VP-field-campaign-observed PSD λ, respectively, holding fall speed and habit fixed.In terms of the optimal-estimation scheme, changes in the PSD λ are offset by corresponding changes in the PSD number density.Therefore, larger particles (smaller slope parameter) require fewer particles to match radar reflectivity, and vice versa.Such a relationship modulates the impact of changes in particle size distribution on estimated snowfall rates.These results suggest that use of the C3VP PSD λ provides a reasonable alternative when lacking coincident MASC measurements.
By contrast, differing assumptions of fall speed led to a factor of 2 differences in retrieved snowfall rates.Accumulations varied from 4.88 mm with MASC fall-speed observations to 10.62 mm with LH74 graupel-like snow parameterizations, holding base habit and PSD fixed.These results were driven by the variability in fall speeds as plotted in Fig. 5 and the linear impact of fall speed on retrieved snowfall rate.Interestingly, retrieved snowfall rates using low level Doppler fall speeds were twice those found using MASC fall-speed observations.Such results imply the need for further research to determine which fall-speed observation is appropriate for a given scene.
Differing assumptions of particle model also led to a factor of 2 differences in retrieved snowfall rates.Estimated snowfall accumulations varied from 4.88 mm for the CloudSat B8pr-30 particle model to 2.23 mm for hexagonal columns, given fixed fall-speed and PSD assumptions.NWS "truth" (4.57mm) fell between results assuming the branched CloudSat particle (4.88 inches) and sector plates (3.33 mm) where both of these habits were seen during the storm event as in Fig. 4. Given the presence of rimed particles throughout the storm, however, we would not expect perfect agreement given use of non-rimed particle models for the retrieval.
In addition to total accumulation, we also considered the temporal distribution of the snowfall during the event.Figure 7 plots NWS hourly snowfall accumulations for the 23 April storm against retrieved values arising from different assumption permutations.Estimates from the CloudSat B8pr-30 particle model (with MASC fall speed and MASC PSD) agreed reasonably well in trend and magnitude with the NWS hourly observations with a Pearson correlation coefficient of 0.65 and a +7 % overestimate of snowfall over the entire event.Retrieval results had a higher correlation with the hourly observations for the first part of the storm (0.92) than the second (0.12) even though differences in subsection snowfall totals were greater for the first section (+11 %) than for the second (+1 %).However, we do not necessarily expect perfect correlation between the NWS observations and the retrieved values given the differences in location between the snow gauge and KAZR.For example, the NWS snow gauges observed significant snowfall (0.76 mm) over the last 2 h of the event even though low KAZR reflectivities suggested perhaps the weakest part of the storm.
The sector plate model (with MASC fall speed and MASC PSD) produced slightly too little snowfall throughout the course of the event.The hexagonal column model (with MASC fall speed and MASC PSD) produced even less.The three curves for the retrieval scenarios using the MASC fall speeds in Fig. 7 were highly correlated with each other (coefficients of 0.99), given common reflectivity and fall-speed assumptions.These curves had a weaker correlation (coefficients near 0.8) with the Doppler scenario given marked discrepancies in fall speed as seen in Fig. 5.

15 May snow event
Table 3 lists retrieved accumulated snowfall amounts for a 15 May snow event as depicted in the upper left panel of Fig. 6.This storm produced approximately 7.62 mm (0.30 inches) liquid equivalent of snowfall as observed by the NWS before transitioning to light rainfall.Although the use of our base CloudSat particle model, MASC fall speed, and MASC λ assumptions produced snowfall estimates close to NWS observations for the 23 April event, these same assumptions performed poorly for this storm.We found a −43 % difference relative to NWS observations.For all retrieval permutations, differences between retrieved snowfall and NWS observations varied from −74 to +56 %.The best agreement with NWS values came with the permutations of sector plates, MASC λ, and KAZR Doppler fall-speed assumptions (+6 %) and CloudSat B8pr-30 particle model, MASC λ, and LH74 dendritic aggregate fallspeed assumptions (−10 %).Again, it is the compensating nature of fall speed and particle model for these two different retrieval permutations that allows them to produce results consistent with NWS observations.Agreement was particularly poor for those cases that assumed the highly reflective Table 4.Total retrieved snowfall amounts over five Barrow snow events for designated retrieval assumption combinations.Nearby NWS observations suggested 16 mm (0.63 inches) of snowfall liquid equivalent.The MASC-derived PSD slope parameters used for these events are listed in Table 1.The " % Difference" value refers to the difference between the NWS and retrieved accumulations over the five snow events.The "Avg % difference" value refers to the average difference for the individual events.hexagonal column particle model, e.g., a difference of −74 % was found with columns when used with MASC-observed fall speeds and PSD λ.The use of these columns (which were not seen during the event) could not generate enough retrieved snowfall to match observations regardless of fallspeed assumptions.
A possible reason for these observed discrepancies in this case is the presence of heavy riming in observed snowflakes.MASC images suggest lumpy graupel-like snowflakes and aggregates throughout the entire snow event as shown in Fig. 8.Such particles would be expected to be denser and thus have more mass for given maximum dimension than the un-rimed particle models employed in this study.Such conditions would cause the snowfall retrievals to underestimate snow rate, with all other variables fixed.Similarly, the particle-backscatter relationship would change given a coating of frozen or possibly liquid water on the particles, further biasing results.Current work by the authors focuses on scattering calculations for rimed particle models assuming different densities.Ideally, these properties could be used with a hydrometeor classification scheme that estimates degree of riming and melt water such as Praz et al. (2017) or similar for better retrieval performance.

Totals for five snow events
Table 4 lists total retrieved liquid water equivalent over the five snow events at Barrow for different retrieval assumption permutations.The approximate accumulated snowfall was 16.0 mm (0.63 inches) liquid equivalent for these five events as measured by the nearby NWS site.Use of our base assumptions (MASC fall speed, MASC PSD λ, and CloudSat particle model) led to total accumulated snowfall of 13.11 mm, which represents a difference of −18 % relative to NWS observations.Use of the C3VP a priori PSD λ led to a slightly better match with NWS observations (−3 %) than the use of MASC PSD λ (−18 %) when used with base particle model and fallspeed assumptions.This relatively close agreement in results is not surprising given the limited sensitivity of retrieval results to changes in PSD λ.Additionally, the a priori λ value of 2.8 mm −1 falls within the range of observed MASC PSD λ values from 0.74 to 3.42 mm −1 .Therefore, it would be expected that retrieved values using the a priori assumption would be greater than those found using MASC values for some cases but less than for others, producing similar results when totaled over the five events.
Other assumption permutations found retrieved liquid water accumulations with differences ranging from −64 to +122 % relative to NWS observations.These results again highlight the likelihood of compensating errors when inverting radar observations to estimate snowfall.The best agreement (3 %) came for two retrieval scenarios with different assumptions for each habit, PSD λ, and fall speed.The Cloud-Sat particle model, C3VP a priori PSD λ, and MASC fall speed combination resulted in 3 % less retrieved snowfall than NWS observations.Sector plates, MASC PSD λ, and a 1 m s −1 fall speed resulted in 3 % more retrieved snowfall than NWS observations.To further illustrate this idea of compensating uncertainties, consider the use of assumptions that are all demonstrably wrong (hexagonal columns, 1 m s −1 fall speed, and C3VP a priori PSD λ) for the vast majority of these observed snow events.This retrieval permutation generates a value for accumulated snowfall of 11.94 mm that is very close to the 13.11 mm found using our base assumptions.
The use of total accumulation over multiple snow events provides a practical metric to evaluate retrieval performance over time.However, it is also important to evaluate the average error for the individual storm events since a retrieval scheme could produce perfect results over multiple storm events while being markedly wrong for each event.Table 4 lists the accumulation-weighted average absolute differences for the individual snow events for different retrieval permutations.We chose this accumulation-weighted approach since two of the storm events had snowfall accumulation less than 1.05 mm liquid equivalent.If we had equally weighted the events, small differences in retrieved snowfall values on the order of 0.5 mm for these events could heavily skew overall average differences.
Table 4 shows accumulation-weighted average absolute differences ranging from 21 to 122 % relative to NWS observations for the different retrieval permutations.For our base assumptions (MASC fall speed, MASC PSD λ, and Cloud-Sat particle model), we found an accumulation-weighted average absolute difference of 36 %.This result is similar to the 35 % value found with use of the C3VP a priori PSD λ, again demonstrating the limited sensitivity of retrieval results to changes in PSD λ for these events.Two of the smallest average difference values (21 and 27 %) were found using hexagonal columns.Again, since columns were very rarely seen during these events, it is likely that such agreement arises from compensating errors from other retrieval assumptions.Similar average difference values (23, 27, and 27 %) were found for a variety of retrieval combinations based upon a sector plate particle model.Such consistent results identify the sector plate model as a likely habit candidate for a Barrow-based snowfall retrieval scheme.
In theory, analyses of retrieval performance could be broken down as a function of meteorological conditions.For example, we might expect that the retrievals would match snow gauge observations poorly under windy conditions.Snowfall gauge measurements can suffer significant uncertainties that are a strong function of wind speed (Black, 1954;Larson and Peck, 1974;Goodison et al., 1998;Yang et al., 2005, Wolff et al., 2015;among numerous others).Likewise, sampling artifacts for MASC particle size and fall-speed observations are likely as increased wind speed and turbulence disrupt the settling of snowflakes as they pass through the MASC sampling volume.In practice, we had only a handful of snow events during the Barrow deployment, making a thorough exploration of this topic difficult.For example, retrievals with our base assumptions matched the snow gauge observations better for the slightly less windy 23 April event (2.85 m s −1 ) than the 15 May (4.40 m s −1 ) event.For the 23 April event, however, the hourly snow gauge and retrieval estimates were more highly correlated during the windier first part of the storm (3.96 m s −1 ) than the calmer second part (2.17 m s −1 ).Both events were also complicated by habit uncertainties, making it difficult to identify a clear wind signature without more observations of snow storms.
Perhaps of more interest, the 21 May event had the strongest average winds at 9.41 m s −1 (21.0 mph), the greatest measured slope parameter, and thus the smallest particles.These results would be consistent with blowing snow composed of shattered ice crystals.If so, we would expect poor retrieval performance as the PSD and habit measurements near the surface would not be representative of the snowfall conditions aloft.The total snowfall accumulation for this high wind event was only 1 mm (LWE), so again it would be unwise to make conclusions on retrieval performance until we get more data.
Finally, we contrast snowfall estimates found with our variational approach with those found from two different 35 GHz reflectivity-snowfall relationships.Our use of in situ observations from the Barrow event, in theory, should produce more accurate results than those derived from these parameterizations.Kulie and Bennartz (2009, hereafter KB09) estimated snowfall as in Eq. ( 7), Z e = 24.04S 1.51 , where they assumed a three-bullet rosette particle model (Liu, 2008b) with scattering particles derived from DDA. S is snowfall in mm h −1 .Matrosov (2007, hereafter M07) assumed aggregate snowflakes with scattering properties modeled from spheres and T-matrix theory to parameterize snowfall through Eq. ( 8).
Z e = 56.00S 1.2 (8) These equations were derived for vertically pointing radar and dry snowfall conditions that lack significant attenuation at 35 GHz.Table 5 lists total estimated snowfall and accumulationweighted average absolute differences for the KB09 and M07 schemes for the five Barrow storm events.KB09 found a total accumulation with a +37 % difference relative to NWS observations with an average error of 50 % for the individual events.M07 found a total accumulation with a −48 % difference relative to NWS observations with an average error of 49 % for the individual events.These values fall well within the range of results presented in Table 4, although they are slightly larger than those found in the better case optimal-estimation retrieval scenarios.Likewise, the correlation coefficient between hourly Z e -S estimates and NWS observations (0.41) was slightly worse than that found between the base optimal-estimation scheme and NWS observations (0.65).Again, it is hoped the continued refinement of our combined radar, in situ approach will produce results consistently more accurate than such Z e -S relationships.
Table 5. Retrieved snowfall amounts as derived from two different Z e -S relationships for the 23 April, 15 May, and all snow events seen at Barrow.These results are contrasted with those found using our base optimal-estimation scheme and NWS observations.Difference values are defined as in Table 4 4 Discussion and conclusions In this work, we present an optimal-estimation retrieval scheme to calculate surface snowfall rates using coincident radar and in situ snow particle observations.The scheme was modified from the W-band CloudSat snowfall algorithm and applied to measurements from the ground-based Ka-band ARM zenith radar located at the ARM NSA Barrow Climate Research Facility site.Multi-Angle Snowflake Camera estimates of particle size distribution and fall speed were used to constrain the inverse calculations based upon KAZR reflectivities.Images of snowflakes from the MASC were used as a guide to select the most appropriate particle model for a given storm event, e.g., branched aggregates, sector plates, or columns.Retrieved snowfall accumulation were compared with snowfall measurements at the nearby NWS Barrow office for a first-order evaluation of our results.Retrieval snowfall values found using MASC observations as assumptions were contrasted with those found using alternate assumptions such as the Locatelli and Hobbs fallspeed parameterizations, Doppler fall-speed observations, and field campaign observations of snow particle size distributions.Use of these different permutations of retrieval assumptions (habit, PSD λ, and fall speed) allowed us to determine which combination of assumptions best matched nearby NWS snowfall observations for five different snow events.Differences between these approaches also quantified the sensitivity of estimated snowfall amounts to PSD λ, particle fall speed, and snowflake particle model.Although the number of events and snowfall totals with coincident MASC and KAZR observations is limited, they do provide an initial data set to demonstrate the technique and to quantify retrieval performance across multiple snow events.
Use of the base assumptions (CloudSat particle model, MASC fall speed, and MASC PSD λ) resulted in estimated snowfall totals over the five events with a −18 % difference relative to nearby NWS snowfall observations.This agreement, of course, could result in part from compensating errors in individual storm totals.We also calculated the accumulation-weighted average absolute difference for the individual events and found a value of 36 %.Such results demonstrate that modification of the CloudSat particle model to Ka-band frequency and use of MASC observations can produce reasonable snowfall values for spring conditions at Barrow.This agrees in spirit with the validation studies from the C3VP program in which simulated reflectivities using the CloudSat particle model matched well with field observations as demonstrated in Fig. 2. The use of the sector plate model with different fall-speed and PSD assumptions often produced average error values lower than 30 %, suggesting the sector plate as a candidate for a Barrow-based snowfall scheme.These results are not surprising since both aggregates and sector plates were seen during these snow events.These better case retrieval scenarios matched NWS observations more closely than two 35 GHz Z e -S relationships that had average errors near 50 % for the Barrow events.
Other combinations of retrieval assumptions found differences in accumulated snowfall ranging from −64 to +122 % relative to NWS observations as listed in Table 4.The accumulation-weighted average absolute difference values for these permutations ranged from 21 to 122 %.The nonunique nature of the snowfall retrieval problem, however, makes it difficult to determine if we get good results for the "right reasons".For example, the use of a trio of assumptions (hexagonal column particle model, C3VP a priori PSD λ, and a 1 m s −1 fall speed) that are all demonstrably wrong yielded overall differences (−25 %) that are similar to those found with our base assumptions (−18 %).Further, this combination had an average difference (28 %) for the individual storms that was indeed better than the value found for our base assumptions (36 %).For this "wrong" scenario, compensating errors in fall speed (bias results high) and particle model (bias results low) likely offset to produce snowfall results that are reasonable across multiple snow events.The use of MASC or other habit observations would eliminate these column scenarios as valid solutions in our retrieval approach.
Use of alternate fall-speed and particle model assumptions led to factor of 2 differences in retrieved snowfall rates when averaged over the five snow events.The assumed PSD slope parameter had less impact on variability in estimated snowfall accumulation.Snowfall totals of 13.11 and 15.47 mm were found when using the MASC PSD λ and the C3VPfield-campaign-observed PSD λ, respectively, holding fall speed and habit fixed.In terms of the optimal-estimation approach, changes in the PSD λ are offset by corresponding changes in the particle number density.Larger particles require fewer particles to match radar reflectivity, and vice versa, thus limiting the impact on estimated snow water content with changing particle size.
Significant sensitivities to multiple variables imply difficulties for the design and evaluation of radar-based snowfall retrieval schemes.Those studies that focus on one or even a few variables may likely lead to misleading conclusions.Again, we hope the use of coincident in situ observations of snowfall microphysical properties will reduce the non-uniqueness of the retrieval problem.
Future work will focus on expanding the coincident radar and ground-based in situ instrumentation approach presented here.We will modify existing particle models and DDA scattering properties to match the graupels and rimed particles often seen in the MASC images as in Figs. 4 and 8. Differences in the scattering properties and mass-dimensional relationships for rimed and un-rimed particles are anticipated to explain the large discrepancies in retrieved and observed snowfall rates observed for the rimed 15 May case.Use of a hydrometeor classification scheme that estimates degree of riming from MASC images such as that developed by Praz et al. (2017) would provide a quantitative manner to link snowflake images with particle scattering properties.Likewise, we will seek to determine the most appropriate fall-speed metric given available instrumentation.Large discrepancies were observed between near-surface MASC fallspeed measurements, lowest radar range bin Doppler velocities, and particle fall speed-dimension parameterizations.Meaningful results from such studies will require the quantification of sampling artifacts for the MASC or other in situ microphysical instrumentation through efforts similar to Kleinkort et al. (2016).They will also require use of state-ofthe-art snow gauge measurements to the extent possible.
Given the limited data available from the NSA Barrow site, we stress the need for more data sets with coincident radar, snowfall microphysical, and snow gauge observations.The general technique presented here for the KAZR and MASC could be adapted to any set of coincident radar and PSD instrumentation.Along these lines, the authors are deploying a Micro Rain Radar, MASC, and Precipitation Imaging Package to two snowfall measurement sites in Scandinavia over the next two winters.These sites, one run by the Norwegian Meteorological Institute near Haukeliseter Fjellstue and one by the Swedish Institute of Space Physics near Kiruna, experience numerous snowfall events representing a diverse range of synoptic and mesoscale conditions.Such a comprehensive set of observations should allow us to refine our retrieval approach and to gain insights into the state-dependence of snowfall microphysics.In turn, these improved snowfall estimates could be used to explain observed differences between ground-based radar, satellite-based radar, and snow gauge estimates of snowfall (Cao et al., 2014;Smalley et al., 2014;Norin et al., 2015;Saltikoff et al., 2015;Speirs et al., 2017) or as input for weather and climate studies.
Data availability.All the data used and produced for the purposes of this paper are freely available and can be requested from the corresponding author.
Competing interests.The authors declare that they have no conflict of interest.

Figure 2 .
Figure 2. Comparisons of observed W-band reflectivities with synthetic reflectivities derived from the DDA scattering properties as taken from a modified figure from Wood et al. (2015).

Figure 3 .
Figure 3. MASC-captured images of snowflakes as observed at the NSA Barrow Climate Research Facility site.

Figure 4 .
Figure 4.The top scene shows KAZR reflectivities (dBZ) for a 23 April snow event at the NSA Barrow site.MASC images of snowflakes suggest a graupel-like structure around 04:00 UTC which transitioned to more pristine structures and aggregates by 13:00 UTC.

Figure 6 .
Figure 6.KAZR reflectivities (dBZ) for snow event systems as observed at the NSA Barrow Climate Research Facility site.Areas of orange and red indicate likely periods with snowfall.

Figure 8 .
Figure 8.Typical rimed and graupel-like particles observed by the MASC on 15 May. .
Figure 5. MASC-observed fall speeds for the 23 April snowfall event at Barrow are plotted in black.Results are compared to calculations from LH74 parameterization schemes for crystal types similar to those observed by MASC during the storm event.Green represents graupel-like snow, while red represents dendritic aggregates.Near-surface Doppler fall speeds from the KAZR are plotted in blue.The MASC fall-speed observations had correlation coefficients of −0.34 and 0.66 with the KAZR and LH74 fall speeds, respectively.

Table 2 .
Retrieved snowfall amounts for a 23 April snow event for designated retrieval assumption combinations.Nearby NWS observations suggested 4.57 mm (0.18 inches) of snowfall liquid equivalent.The MASC-derived PSD slope parameters (λ) used for this event are listed in Table1.

Table 3 .
Retrieved snowfall amounts for a 15 May snow event for retrieval assumption combinations.Nearby NWS observations suggested 7.62 mm (0.30 inches) of snowfall liquid equivalent.The MASC-derived PSD slope parameters (λ) used for this event are listed in Table1.
Figure 7. Hourly NWS snowfall accumulations (mm) plotted against retrieved values for the 23 April event for different retrieval assumption permutations.MASC PSD slope parameter (λ) is assumed for all curves. .