Articles | Volume 14, issue 10
Research article
26 Oct 2021
Research article |  | 26 Oct 2021

Retrieving microphysical properties of concurrent pristine ice and snow using polarimetric radar observations

Nicholas J. Kedzuf, J. Christine Chiu, V. Chandrasekar, Sounak Biswas, Shashank S. Joshil, Yinghui Lu, Peter Jan van Leeuwen, Christopher Westbrook, Yann Blanchard, and Sebastian O'Shea

Ice and mixed-phase clouds play a key role in our climate system because of their strong controls on global precipitation and radiation budget. Their microphysical properties have been characterized commonly by polarimetric radar measurements. However, there remains a lack of robust estimates of microphysical properties of concurrent pristine ice and aggregates because larger snow aggregates often dominate the radar signal and mask contributions of smaller pristine ice crystals. This paper presents a new method that separates the scattering signals of pristine ice embedded in snow aggregates in scanning polarimetric radar observations and retrieves their respective abundances and sizes for the first time. This method, dubbed ENCORE-ice, is built on an iterative stochastic ensemble retrieval framework. It provides the number concentration, ice water content, and effective mean diameter of pristine ice and snow aggregates with uncertainty estimates. Evaluations against synthetic observations show that the overall retrieval biases in the combined total microphysical properties are within 5 % and that the errors with respect to the truth are well within the retrieval uncertainty. The partitioning between pristine ice and snow aggregates also agrees well with the truth. Additional evaluations against in situ cloud probe measurements from a recent campaign for a stratiform cloud system are promising. Our median retrievals have a bias of 98 % in the total ice number concentration and 44 % in the total ice water content. This performance is generally better than the retrieval from empirical relationships. The ability to separate signals of different ice species and to provide their quantitative microphysical properties will open up many research opportunities, such as secondary ice production studies and model evaluations for ice microphysical processes.

1 Introduction

Ice-containing clouds play an important role in Earth's radiation budget and global precipitation (Baran et al., 2009; Field and Heymsfield, 2015; Mulmenstadt et al., 2015; Li et al., 2014). Their formation and evolution involve processes of ice nucleation, ice multiplication, aggregation, and riming, which are closely linked to atmospheric conditions and dynamics (DeMott et al., 2011; Field et al., 2017; Gultepe et al., 2017; Korolev et al., 2020). Such complex interactions make it challenging to complete our understanding of these ice microphysical processes and represent them well in models (Korolev et al., 2017; Morrison et al., 2020).

Polarimetric radar measurements contain information on ice properties and have been proven useful for studying ice microphysical processes (e.g. Kennedy and Rutledge, 2011; Grazioli et al., 2015; Moisseev et al., 2015). Many empirical relationships were developed to provide important bulk properties such as the ice water content (e.g. Ryzhkov et al., 1998; Lu et al., 2015), median volume diameter, and number concentration (e.g. Murphy et al., 2020), but they cannot inform the partitioning between ice species. The ability of the partitioning is of particular importance for studying the aggregation process because it provides information on the size and number concentration of pristine ice and aggregates.

However, separating signals of pristine ice from aggregates in polarimetric radar data is challenging because larger snow aggregates often dominate the radar reflectivity and mask contributions of smaller pristine ice crystals (Hogan et al., 2002; Keat and Westbrook, 2017). As a result, information from horizontal reflectivity (ZH) alone is insufficient to characterize mixtures of ice hydrometeors (Oue et al., 2018), and it is necessary to incorporate other radar observables in retrieval methods. By exploiting distinct fall behaviours between pristine ice and aggregates, Spek et al. (2008) used ZH, differential reflectivity (ZDR), and the Doppler spectrum to retrieve particle size distribution (PSD) parameters of pristine ice and snow aggregates. Without the use of the Doppler spectrum, Schrom et al. (2015) used ZH, ZDR, and specific differential phase shift (KDP) to estimate the PSD of pristine ice in the dendritic growth zone of Colorado winter storms. KDP is a great addition in their approach since it is mainly determined by ice number concentration. Unfortunately, these three radar observables remain insufficient, and their partitioning between pristine ice and aggregates was weakly constrained. To improve the partitioning, Keat and Westbrook (2017) showed that the relative radar signal contributions of pristine ice embedded in snow aggregate populations can be quantified using ZH, ZDR, and the co-polar correlation coefficient (ρhv), but they have not attempted to use their partitioning to provide quantitative retrievals of pristine-ice number concentration, water content, and particle size.

The objective of the paper is to present an ensemble cloud retrieval method (dubbed ENCORE-ice) for simultaneously retrieving the number concentrations, sizes, and ice water contents of concurrent pristine ice and snow aggregates from measurements of ZH, ZDR, KDP, and ρhv. This framework provides full error statistics and characterizes sub-species from radar signals, which is an advance on the existing methods. The polarimetric radar observations and the retrieval method are detailed in Sect. 2. The ancillary datasets for evaluations are introduced in Sect. 3. Section 4 presents evaluation results using synthetic datasets and actual observations from Chilbolton, United Kingdom, in 2018. Finally, Sect. 5 summarizes the key findings and discusses potential applications.

2 Radar observations and ENCORE-ice

2.1 Polarimetric radar data

Our retrieval method uses four polarimetric observables. The first observable is the horizontal reflectivity ZH, which provides information on particle size and concentration, but its dependence on size is much stronger. As such, ZH is dominated by contributions from snow aggregates because their sizes, and thus their backscatter cross-sections, are typically much larger than those of pristine ice crystals. The second observable is the differential reflectivity ZDR, which provides information on particle shape and orientation. A ZDR of 0 dB indicates spherical particles because of equal backscattered power in each polarization. Snow aggregates yield low ZDR values (about 0–0.6 dB; see Hogan et al., 2012) as a result of their sparse and irregular morphology, with the component crystals oriented at a wide range of angles. In contrast, pristine ice particles can yield ZDR values of several decibels because of their aspect ratios and preferential horizontal orientation when falling. Heterogenous regions with concurrent pristine ice and snow aggregates are therefore associated with higher ZDR values than if only snow aggregates were present. The third observable is the co-polar correlation coefficient ρhv, the correlation coefficient between horizontally and vertically backscattered power, which provides information on the diversity of particle shape in a radar sample volume (Kumjian, 2013; Keat et al., 2016). ρhv is unity in homogenous regions but tends towards lower values (e.g.  0.97) in the presence of heterogenous hydrometeor types. Finally, the fourth observable is the specific differential phase shift KDP, which provides information on particle number concentration, shape, and orientation.

Our case study is based on polarimetric radar data (Bennett, 2020) from the Parameterizing Ice Clouds using Airborne obServationS and triple-frequency dOppler radar data (PICASSO) field campaign in Chilbolton, UK, in 2018–2019. During the campaign, the National Centre for Atmospheric Science mobile X-band dual-polarization Doppler weather radar (NXPol; Neely et al., 2018) operated with a 0.98 beam width, 150 m range resolution, and maximum range of 150 km. The radar performed two back-to-back, fixed-azimuth range–height indicator (RHI) scans every 7 min, and each scan was completed in 18 s. Throughout 13 February 2018, RHI scans were performed along the 243 radial and intercepted by the NCAS-managed Facility for Airborne Atmospheric Measurements (FAAM) aircraft on several occasions, providing a unique opportunity for evaluation. Key characteristics of NXPol are summarized in Table 1.

Table 1Characteristics of the NXPol polarimetric radar. Further specifications and details can be found in Neely et al. (2018).

Download Print Version | Download XLSX

2.2 ENCORE-ice

ENCORE is an ensemble-based retrieval method that has previously been used to retrieve three-dimensional cloud microphysical properties (Fielding et al., 2014) and one-dimensional cloud and drizzle properties (Fielding et al., 2015), but several key components are modified here for ice retrieval.

2.2.1 Particle size distribution

We approximate the PSD of pristine ice and aggregates by normalized gamma distributions, given as (Testud et al., 2001)

(1) n ( D ) = N 0 f μ D ; D 0 ,

where N0 is the normalized number concentration and n(D)dD is the number of particles in the range of the maximum particle dimensions (D,D+dD). The choice of the size descriptor in Eq. (1) is because in situ cloud probe data and the ice scattering database are both given based on the maximum particle dimension. The function fμ is defined as

(2) f μ D ; D 0 = 6 3.67 4 3.67 + μ 4 + μ Γ 4 + μ D D 0 μ exp - 3.67 + μ D D 0 ,

where μ is the shape parameter of the PSD and D0 is the diameter used for normalizing D. Following Mason et al. (2018), we assume a constant shape parameter of μ=2. Several studies have shown that the retrieved ice water content is relatively insensitive to the choice of shape parameter (e.g. Delanoë et al., 2005; Spek et al., 2008); we also found that our number concentration retrieval is not sensitive to μ.

Table 2Examples of mass–size relationships (taken from Mason et al., 2018).

Download Print Version | Download XLSX

From the PSD, the total ice number concentration (NI) and the total ice water content (qI) can be computed by

(3) N I = 0 n ( D ) d D = 0 n P ( D ) + n A ( D ) d D = 0 N 0 , P f μ D ; D 0 , P d D + 0 N 0 , A f μ D ; D 0 , A d D = N P + N A


(4) q I = 0 m ( D ) n ( D ) d D = 0 m P ( D ) n P ( D ) d D + 0 m A ( D ) n A ( D ) d D = q P + q A ,

respectively, where n(D) is the combined PSD from nP(D)+nA(D) and the subscripts P and A denote contributions from pristine ice and snow aggregates, respectively. m(D) is the mass at a given maximum particle dimension D. The mass–size relationship can be formulated as

(5) m ( D ) = a D b ,

where a and b are the pre-factor and exponent, respectively. These coefficients depend on ice habit and have been estimated from past aircraft in situ and surface observations as shown in Table 2. From the PSD, we define and calculate the effective mean diameters (Deff) as

(6) D eff = 0 n ( D ) D 4 d D 0 n ( D ) D 3 d D ,

which is the ratio of the fourth to the third moment of the PSD. To compare our retrieval with the empirical estimates (as discussed in Sect. 3), we also calculate an effective mean diameter using the equivalent melted diameter (Dmlt) as the size descriptor, defined as

(7) D eff , mlt = 0 n D mlt D mlt 4 d D mlt 0 n D mlt D mlt 3 d D mlt ,


(8) D mlt = 6 m ( D ) π ρ w 1 3 = 6 a D b π ρ w 1 3

and ρw is water density.

2.2.2 The basis of ENCORE-ice

The state vector (x, i.e. variables to be retrieved) for each ensemble member is defined as

(9) x = log 10 N 0 , P i = 1 G , log 10 D 0 , P i = 1 G , log 10 N 0 , A i = 1 G , log 10 D 0 , A i = 1 G T ,

where the superscript i represents the index of the range gate and the total number of gates to be retrieved is G. Let us use Q members to form an ensemble; i.e.

(10) x 1 , , x Q

such that the mean x represents the best estimate of the state vector and the spread of the ensemble members around the mean represents the uncertainty in the best estimate.

Using the iterative stochastic ensemble Kalman filter approach (Evensen et al., 2019), each ensemble member is updated based on

(11) x k a = x k f + E f w k ,

in which xkf and xka are the prior and posterior ensemble member k, respectively;

(12) E f = x 1 f - x f , , x Q f - x f

is the initial ensemble matrix with the prior mean (xf) subtracted; and wk comprises weight vectors that are calculated from iteratively minimizing the following cost function:

(13) J ( w k ) = 1 2 w k T w k + 1 2 y - h x k f + E f w k - ε k T R - 1 y - h x k f + E f w k - ε k .

In Eq. (13), the observation vector y is defined as gate-by-gate radar observables:

(14) y = Z H ( i = 1 , , G ) , Z DR ( i = 1 , , G ) , - ln K DP i = 1 , , G , - ln ρ hv i = 1 , , G T ,

where ZH and ZDR are in decibels. Radar observations at one range gate will influence the estimation of the state vector at another gate if these two gates are within the pre-defined radius, which will be explained in more detail in Sect. 2.3. h(x) represents the forward model for simulating polarimetric radar observables from the state vector x, and εk is a random perturbation vector drawn from the observation error distribution, which is estimated to be Gaussian with a mean of zero and covariance matrix R (Evensen et al., 2019, with modification from Van Leeuwen, 2020). The covariance matrix R is diagonal with standard deviations given in Table 3.

Table 3Estimated observational errors for X-band observables based on standard and benchmark procedures, adapted from Bringi and Chandrasekar (2001, pp. 359–376) and Wang and Chandrasekar (2009).

a Estimated by the uncertainty of 0.05 km−1 for a typical value of KDP=0.5 km−1. b Estimated by the uncertainty of 0.01 for ρhv=0.95.

Download Print Version | Download XLSX

As detailed later in Sect. 2.3, the prior is assumed Gaussian, and there is no prior correlation between variables N0,P, D0,P, N0,A, and D0,A. But there is correlation in the vertical (i.e. between gates) for each variable in our setup. We have also used a prior with large uncertainty, approximately 1–2 orders of magnitude in the state variables, such that the influence of the prior is minimal. In contrast to the prior, no Gaussian assumption is made in the posterior ensemble members, although the retrieval statistics are largely focused on their means and standard deviations.

2.2.3 Simulating radar observables for h(x)

To model polarimetric radar observables from the assumed PSD, knowledge of the single-scattering properties of ice particles is required. Many scattering databases of realistically shaped ice particles at radar wavelengths are available, and we used Lu et al. (2016) because of the following considerations. Several existing scattering databases assume a total random orientation of the scatterers, e.g. Liu (2008), Hong et al. (2009), Kuo et al. (2016), and Eriksson et al. (2018). Such assumption cannot explain polarimetric radar signals which are produced by non-spherical scatterers with preferred orientations with respect to the zenith direction. The database of Brath et al. (2020) assumes scatterers possess arbitrary fixed orientations relative to the zenith direction but includes only hexagonal plates and aggregates consisting of hexagonal plates. We found that the database described in Lu et al. (2016) fits our needs in the current polarimetric radar study since it contains all necessary polarimetric scattering data in many fixed orientations of a large variety of ice crystal species, including plates, columns, dendrites, and aggregates. The single-scattering properties for each species are available for a range of crystal maximum dimensions, thickness ratios, and types. The pristine habits generally begin at  0.1 mm and do not exceed 6 mm, whereas the aggregates begin at  0.4 mm and extend to approximately 18–45 mm. Multiple morphological realizations per maximum dimension are available for dendrites and aggregates to account for their complexities. Note that for a given size, natural aggregates may have substantially different properties compared to the realizations available in the database.

The scattering calculations were conducted using the generalized multi-particle Mie method (GMM; Xu, 1995) and the discrete dipole approximation (DDA; Yurkin and Hoekstra, 2011). We used properties calculated from GMM because DDA calculations are not available for aggregates. Specifically, we use the amplitude scattering matrix elements in the forward and backward direction for horizontally and vertically polarized radiation, denoted as Shhf,b and Svvf,b, where the superscript and subscript represent the scattering direction (i.e. forward or background) and the polarization status (horizontally or vertically), respectively. From the assumed PSD and the amplitude scattering matrix elements, radar observables for a single sample volume containing multiple ice particle habits can be derived as shown in Appendix A.

2.3 Practical considerations

There are several practical considerations for ENCORE-ice implementation. The first consideration is ice habit. The scattering database provides three habits (plates, dendrites, and columns) for pristine ice. Since the temperature found in PICASSO mostly ranged between 5 and 25 C, all three types of pristine ice can be the preferred habit (see examples in Fig. 1). Currently, we do not predetermine the ice habit. Instead, we ran our retrieval algorithm for all three habits independently and then selected the most appropriate one based on the agreement in the measured and forward-simulated radar observables.

Figure 1Examples of particle images from the Stratton Park Engineering Company Two-Dimension Stereo (2D-S) probe, showing the presence of (a) column, (b) plate, and (c) aggregates of dendrites on 13 February 2018. Each image frame is 1.28 mm high, taken from one of the probe channels only since the other channel was not working properly on this day.


Similarly, the scattering database provides five types of aggregates: two of them were constructed using ice columns (LD-N1e and HD-N1e) and three of them using stellar ice crystals (LD-P1d, LDt-P1d, and HD-P1d). Each aggregate in the database (Lu et al., 2016) was generated by first specifying a reference spheroid with a given horizontal maximum dimension and an aspect ratio of 0.6, defined as the ratio of the lengths of the polar axes to the equatorial axes. Then, small monomers were added to the reference spheroid one at a time; any parts of the monomer that were outside the reference spheroid were removed. This procedure was repeated until the mass of the aggregate reached the desired total mass. As a result, the aspect ratio of the aggregate generated in the database was not necessarily the same as the reference spheroid (0.6).

Figure 2 shows the average aspect ratios for aggregate types available in the database, which were calculated by averaging ratios of the maximum vertical dimension to the maximum horizontal dimension for all realizations within one size bin. Compared to Garrett et al. (2015) and Jiang et al. (2017), who reported an aspect ratio range between 0.3 to 0.6 from observations of falling aggregates at the surface, we found that LDt-P1d and HD-P1d exhibit a similar aspect ratio range. In the mass–size relationship for LDt-P1d we used a=0.000482 and b=1.97 in units of centimetres, grams, and seconds (cgs) as in Table 2, based on aggregates composed of ordinary dendritic crystal (Kajikawa, 1989; Botta et al., 2011), whereas for HD-P1d we used a=0.00145 and b=1.80 in units of cgs, based on aggregates of thin plate (Mitchell and Heymsfield, 2005; Botta et al., 2011). The mass–size relationship of HD-P1d is very close to unrimed aggregates (Erfani and Mitchell, 2017) and more aligned to values in the recent literature listed in Table 2. Hence, we select HD-P1d as the prescribed choice for aggregates.

Figure 2Aspect ratios of various aggregate types available in the scattering database as a function of their maximum dimensions. The aspect ratio is defined as the ratio of the sizes of the minor axes to the major axes. The grey shading between 0.3 and 0.6 represents the typical range of snow aggregate aspect ratios observed in nature (Garrett et al., 2015; Jiang et al., 2017).


The second consideration is the prior used to generate the first guess for ensemble members. Using over 70 h of in situ aircraft observations from a wide range of field campaigns spanning diverse cloud and temperature regimes, Delanoë et al. (2014) characterized the PSDs of ice particles by the normalized gamma distribution. They found that N0 ranged between 1 and 10 000 L−1 mm−1 with a mean of 100 L−1 mm−1, and the median volume diameter (MVD) ranged between 0.2–0.8 mm with a mean of 0.5 mm for the temperature zone of 10 to 20 C. Additionally, Tiira et al. (2016) analysed surface measurements of ice particle number concentration from the Precipitation Imaging Package during the Biogenic Aerosols – Effects on Clouds and Climate field campaign. They found that N0 ranged mainly from 1 to 100 L−1 mm−1 and MVD ranged from 0.5 to 5 mm. As these were surface based, the measured PSDs from Tiira et al. (2016) are more representative of the characteristics of snow aggregates. Note that these values were derived using the equivalent melted diameter as the size descriptor, not the maximum particle dimension. Hence, these values are used only to point out a possible range and serve as a starting point for us to construct the prior.

Table 4The prior and uncertainty used in ENCORE-ice. The means at the lowest radar gate are given in the physical state space, and the rest are in the transformed state space (i.e. log 10). Retrieval is performed using two different sets of the prior; the second set uses values in the parentheses, and the rest remain unchanged. All radar gates above the lowest gates are perturbed by an AR1 red noise process with a vertical correlation of 0.999 and a standard deviation that is half of the standard deviation at the lowest gates.

Download Print Version | Download XLSX

Based on these observational ranges mentioned above, our prior is designed as follows. We started with the lowest radar gate, randomly assigning (N0,P, D0,P, N0,A, D0,A) from normal distributions with the means and standard deviations listed in Table 4. Next, we applied a slope for each ray to provide initial guesses for other radar gates. The slope was randomly selected from a normal distribution described in Table 4. Because the prevalence of active ice nuclei is a function of temperature and thus a function of height as well (DeMott et al., 2010), N0,P likely increases with height, and thus the slopes in the prior are assumed to have a positive mean. In contrast, the dependence of D0,P, N0,A, and D0,A on height is less clear (e.g. Field et al., 2005). For practical reasons, the slopes applied for D0,P, N0,A, and D0,A are assumed to have a slightly negative mean. The slightly negative slope avoids unrealistic priors for radar gates at higher altitudes since we used the logarithm form in the state vector. Finally, red (AR1) noise was added over the vertical with a correlation coefficient of 0.999 and a zero-mean random perturbation with a standard deviation that is half that of the lowest radar gate. Note that without this noise term each ensemble member would be a straight line in the vertical for each variable with a different slope. Since the fundamental idea behind ensemble retrievals is that the true atmospheric profile is drawn from the same distribution as the prior ensemble members and we know the true atmospheric profile is not a straight line, we add random noise with non-zero vertical correlation to each ensemble member profile to make each of them more realistic.

Compared to values reported in Delanoë et al. (2014) and Tiira et al. (2016), we have chosen lower means for N0,P and N0,A to start with. This is because the state vector space is the logarithm of N0. Positive slopes make the changes in N0 much more dramatic in the vertical than those with negative slopes. As a result, large starting values of N0,P and N0,A will lead to unrealistically high concentrations at higher altitudes in the prior. In contrast, we choose larger means for D0,P and D0,A because of the assumed negative slopes in both D0,P and D0,A. In general, the range in our prior is large, approximately 1–2 orders of magnitude across all ensemble members over all the vertical profiles. For such a widely spread prior the solution will be dominated by the observations.

The third consideration is the number of the ensemble members used in ENCORE-ice. Ideally, a large ensemble size is needed to ensure that the sampled prior, as well as the solution, is representative. However, a large ensemble size is computationally expensive. Therefore, we applied a localization scheme to reduce the required number of ensemble members so that we shorten the computational time while achieving the same mean retrieval and associated uncertainty. The localization scheme operates on each gate and takes only observations close to that gate into account to find the solution. This is implemented by multiplying the observation error variance of each observation with an exponential function of the distance between that observation and the gate that is being updated, such that observations far from the gate have less influence. The influence radii vary linearly with height, with one gate at the lower level and about five gates at the upper level. Using our synthetic datasets, we have found that 50 ensemble members with the localization scheme are able to produce similar mean retrievals and associated uncertainty to a non-localized ensemble of size 500. The number of iterations is set to 20, although the solutions often converge at the 10th iteration.

Finally, all radar data underwent the following quality checks and corrections before being used for retrieval:

  • ZH and ZDR were corrected for attenuation due to liquid water, using the method described in Bringi and Chandrasekar (2001, pp. 490–512). The attenuation due to ice at the X-band is negligible and thus ignored here (Vivekanandan et al., 1999).

  • Systematic biases in ZDR were identified using zenith-pointing ZDR observations. As hydrometeors produce a ZDR of 0 dB when viewed at the zenith angle due to their spherical symmetry (e.g. for raindrops) or lack of preferential azimuthal orientation (e.g. for ice particles), any residual ZDR can be treated as bias and removed (Seliga et al., 1981). We have found the ZDR correction factors to be 0.2 dB for the PICASSO cases.

  • KDP was calculated using the method of Wang and Chandrasekar (2009).

  • Once all corrections are applied, measurement noises were removed using a cubic spline approach (Craven and Wahba, 1979).

Additionally, to ensure that gates are associated with sufficient information for our method, we exclude gates that exhibit one or more of the following conditions:

  • gates within 500 m of the 0 C level, avoiding contamination from liquid hydrometeors in the radar sample volume because our state vector is not designed for that;

  • gates where ρhv exceeds 1.0 because these values are unphysical;

  • gates with a signal-to-noise ratio (SNR) of less than 20 dB because of a lack of detectable hydrometeors, a threshold chosen because the precipitating region and its surrounding area typically have SNR values larger than 30–40 dB;

  • radar rays with elevation angles greater than 50 because ZDR tends towards 0 dB at higher elevation angles and polarimetric information becomes ambiguous;

  • gates with ZDR below 0.25 dB, regardless of their elevation angles, because the relative contributions of pristine ice and aggregates become ambiguous at lower values, as indicated in Keat and Westbrook (2017);

  • gates with KDP below 0.1 km−1 to ensure a sufficient number concentration of pristine ice.

Note that negative KDP values indicate the presence of conical graupel (Aydin and Seliga, 1984) or the vertical reorientation of pristine ice crystals in the presence of thunderstorm electric fields (Hubbert et al., 2014). Since our state vector includes only aggregates and horizontally oriented pristine ice, excluding gates with low KDP values is necessary.

3 Independent observations and retrievals for evaluations

3.1 In situ aircraft measurements from PICASSO

During PICASSO, the FAAM aircraft performed multiple transects from Chilbolton to Dorset (50.82 N, 2.56 W) at varied altitudes. Figure 3 depicts the flight path for 13 February 2018, which was a typical pattern during the campaign.

Figure 3Flight paths on 13 February 2018 between 03:26 and 10:21 UTC. The red dot denotes the location of NXPol in Chilbolton, UK, while the path in cyan denotes the path during 06:00–09:00 UTC in which retrievals are evaluated in Sect. 4.2.

To evaluate our cloud retrieval, we use in situ measurements of liquid water content and total water content (i.e. the sum of ice and liquid water contents) from a Nevzorov probe and PSD measurements from a high-volume precipitation spectrometer (HVPS; SPEC Inc, USA). The HVPS is an optical-array particle-imaging probe, which collects images of ice crystals with a pixel resolution of 150 µm. Size distributions of particles between 75 and 19 275 µm were derived from their images and reported here using the maximum particle dimension as the size descriptor. A description of the data processing and quality control can be found in Crosier et al. (2011), and the sources of uncertainties were discussed in O'Shea et al. (2021). All in situ datasets were averaged to 5 s intervals for statistical reliability (Protat et al., 2007). We only use in-cloud samples, defined as having an ice water content (IWC) greater than 0.01 g m−3. Additionally, although our retrieval provides microphysical properties of pristine ice and aggregates separately, we focus on evaluating bulk properties to avoid the ambiguity introduced by applying a threshold to separate these two species in observed PSDs.

Three bulk properties are used for evaluations. Firstly, the total ice number concentration (denoted as NI,HVPS) is calculated by integrating the observed PSD. The associated counting uncertainty is estimated as

(15) 100 % N I , HVPS V HVPS Δ t ,

where Δt is the HVPS sampling time resolution and VHVPS is the sample volume, approximately 310 L s−1. Secondly, IWC, denoted as qI,NEV, was derived by taking the difference between the total and liquid water contents measured by the Nevzorov probe. Similarly to Abel et al. (2014), both the total and liquid water contents were corrected for changes in aircraft altitude and environmental conditions. Finally, effective mean diameters from HVPS PSDs, Deff,HVPS, defined as

(16) D eff , HVPS = 0 n HVPS ( D ) D 4 d D 0 n HVPS ( D ) D 3 d D ,

were calculated, using the same definition as that of Eq. (6).

The evaluations in the total ice number concentration, ice water content, and effective mean diameter all together allow us to indirectly examine whether the partitioning between pristine ice and aggregates is appropriate. A more direct comparison would be ideal but requires classifying each individual particle in image data, which is not trivial and is beyond the scope of this work.

3.2 Bulk ice properties from empirical relationships

As mentioned in Sect. 1, several studies have proposed empirical relationships for estimating IWC, particle size, and ice number concentration. In this study, we compare our retrieval with estimates from Ryzhkov and Zrnic (2019) because of their availability of the ice number concentration estimates. The relationships in Ryzhkov and Zrnic (2019) were based on theoretical calculations, using an assumed exponential size distribution for 12 ice habits. Their method takes advantage of the features that the reflectivity difference between horizontal and vertical polarization (ZDP) is proportional to the third moment of the PSD and that KDP is proportional to the first moment of the PSD. As a result, the ratio of ZDP to KDP is proportional to the second moment of the PSD and can be used to estimate the mean volume diameter of ice particles (Murphy et al., 2020):

(17) D emp = - 0.1 + 2 Z DP K DP λ 1 2 ,

where Demp is the mean volume diameter in millimetres with the subscript emp denoting empirical estimates, ZDP is in units of mm6 m−3, KDP is in  km−1, and λ is the radar wavelength in millimetres. Murphy et al. (2020) also estimated the number concentration and IWC using


where Nemp and qI,emp are in units of L−1 and g m−3, respectively; ZH is in units of decibels relative to Z; and ZDR is unitless. For convenience, we refer to retrievals from these empirical relationships as “Murphy20” hereafter. Based on the evaluation conducted by Murphy et al. (2020) for a stratiform region of a mesoscale convective system over Oklahoma, Nemp scatters significantly with respect to in situ measurements; qI,emp and Demp tend to be systematically biased low but outperformed other empirical relationships. Since, in theory, these derived relationships are not sensitive to ice particle shape and orientation, they remain a good starting point for intercomparisons.

Note that these empirical relationships are designed for radar volumes that include only one species. Hence, if a radar volume is known to include a mixture of different species, caution should be exercised when interpreting their results. Additionally, Eq. (17) was derived using the equivalent volume diameter as the size descriptor in the PSD, and thus Demp cannot be used directly for comparisons to our retrieval that is based on the maximum particle dimension as the size descriptor. Instead, we need to trace back their derivations to find their retrieved PSD, convert the equivalent volume diameter to the equivalent melted diameter (Dmlt), and then calculate the effective mean diameter (Deff,mlt) using Eq. (7) for intercomparisons. The details can be found in Appendix B.

4 Results

4.1 Evaluation using synthetic data

In this section we use synthetic polarimetric radar data to evaluate our retrieval and identify any potential issues. The synthetic dataset was generated as follows. We first generated 501 profiles from the prior used in ENCORE-ice and then randomly selected a profile that has a relatively wide range of ZH and ZDR values for testing. Along with the forward model described in Sect. 2.2.3, this selected profile is used to generate synthetic radar measurements and serves as the “truth” in this evaluation experiment. Because the truth profile and the initial ensemble members were generated from the same prior and used the exact same forward models, any retrieval error found in this experiment is due to the combination of the observed uncertainty and the retrieval method itself only. Hence, the design of this experiment does not allow us to evaluate errors due to the representativeness of forward models or the prior, which likely exist in real-world applications.

Figure 4 shows the synthetic radar measurements over 20 gates with a resolution of 50 m at a given elevation angle of 30, based on the truth profile shown in Fig. 5. The chosen number of gates is arbitrary but represents a frequent scenario in the radar scans collocated with in situ data during PICASSO. In this scenario, the total ice number concentration is dominated by pristine ice and the total ice water content is dominated by aggregates. The truth has an NP range of between 5–20 L−1 and an NA range between 1–3 L−1; D0,P ranges between 0.4–1 mm, while D0,A ranges between  2–5 mm. The combined Deff varies from 1.5 to  5 mm, generally close to D0,A as expected, since it is weighted by size to the third power and mainly controlled by the species with large particle sizes.

Figure 4Profiles of (a) ZH, (b) ZDR, (c) KDP, and (d) ρhv for synthetic observations (black line) calculated from the ice cloud properties given in Fig. 5 and for the mean of forward simulations from the ensemble (red dots). Grey shading denotes the observational uncertainties given in Table 3, while error bars in orange denote retrieval uncertainty calculated as the 1-standard-deviation spread of the ensemble simulations.


Figure 5Profiles of (a) normalized number concentration and (b) normalization diameter for pristine ice and (c) normalized number concentration and (d) normalization diameter for aggregates. Panel (e) represents the total number concentrations, and panel (f) represents ice water contents of pristine ice and aggregates. Truth is denoted by solid black lines in (a)(d) and by dashed lines in (e) and (f). The retrieved ensemble means are denoted by solid blue and red lines with shading that represents the 1-standard-deviation spread of the ensemble members. The habit of pristine ice is plate in this experiment.


Table 5Truth means and the means, root-mean-square error (RMSE), and biases in retrieval in the synthetic dataset experiment.

Download Print Version | Download XLSX

The forward-modelled observables of the ENCORE-ice solution in Fig. 4 agree well within the uncertainty in the synthetic values, providing confidence in retrievals. As shown in Fig. 5, the retrieval captures the vertical trend of the truth; the retrieval uncertainty estimated from the spread of the ensemble members also appears reasonable since the truth falls within the retrieval uncertainty. Because many combinations of N0,A and D0,A could lead to the same ZH, we see some compensatory effects between N0,A and D0,A in aggregates at the lower layer. The errors are compensated for so that the error in the total ice water content is not enhanced, as shown in Table 5. Overall, the retrieval biases in the combined total number concentration, water content, and effective diameter properties are within 5 %, and the root-mean-square errors are small (see Table 5).

To conclude, this evaluation experiment demonstrates that the combination of these four radar observables is appropriate and the current observational uncertainty is sufficient for us to separate signals of pristine ice from aggregates. The errors in retrieved pristine-ice properties are small, and thus further physical interpretation based on the associated vertical profiles can be made to understand the underlying microphysical processes. For aggregates, the errors in retrieved size diameter and water content are small, but the vertical variations in retrieved number concentration may not follow the truth exactly due to the possible compensatory effects between number concentration and particle size.

4.2 Evaluation using PICASSO data

The case of 13 February 2018 from the PICASSO campaign represents a stratiform precipitating cloud system associated with a frontal passage. Using a radar scan at 08:37 UTC as an example, Fig. 6 shows a significant area with reduced ρhv, and enhanced ZDR and KDP at  3 km height, which suggests the presence of enhanced pristine ice embedded in snow aggregates. Based on the temperatures measured by the aircraft (Fig. 7), this area is in a temperature zone approximately between 12 and 18 C, and thus the preferred ice habit is likely to be dendrite and plate for this radar scan. During this radar scan, the FAAM aircraft was too far away to provide meaningful comparison, but cloud images showed that dendrites were present most of time during this period.

Figure 6Height–range plots of observed (a) ZH, (b) ZDR, (c) KDP, and (d) ρhv from the RHI scan at 08:37 UTC on 13 February 2018 during the PICASSO field campaign. The dashed red line denotes the 0 C level, while the dashed black line denotes the approximate flight altitude of the FAAM aircraft during the scan. The black polygon denotes the region that has enhanced ZDR and KDP and reduced ρhv.


Figure 7A temperature profile composited from aircraft in situ data during 06:00–09:00 UTC on 13 February 2018. Data between 06:33:30–06:39:20 UTC were unphysical and thus excluded. The error bars represent 1 standard deviation of sampled observations.


Figure 8Retrieval performance for a radar ray at 08:37:39 UTC. Observed and forward-simulated profiles of (a) ZH, (b) ZDR, (c) KDP, and (d) ρhv. The shading in (a)(d) represents the observational uncertainty. The red dots represent the mean of the ensemble simulations, and the error bars represent 1 standard deviation in forward simulations. Panels (e)(j) represent the retrieved mean normalized pristine-ice number concentration, pristine-ice normalization diameter, normalized aggregate number concentration, aggregate normalization diameter, the total number concentrations, and the total ice water content, respectively. Panels (k) and (l) represent the individual and combined effective mean diameters using the maximum particle dimension and the equivalent melted particle size as the size descriptor, respectively. The shading in (e)(l) represents 1-standard-deviation uncertainty in retrieval. For comparisons, retrievals from Murphy20 are co-plotted in (i), (j), and (l).


Figure 8 shows detailed retrieval performance for a ray taken from the radar scan in Fig. 6. For this case, retrievals using the dendrite habit perform best; the habit suggested by our retrieval is consistent with observed particle images. As shown in Fig. 8a–d, the forward modelled radar observables agree well with the observed vertical profiles. The normalization diameters for pristine ice and aggregates (D0,P,D0,A) are about 2.5 and 5–6 mm, respectively. Retrieved NI is relatively constant at  5 L−1, due to the opposite vertical variations between NP (increasing with height) and NA (decreasing with height). In contrast, qI decreases with height from 1 to 0.2 g m−3 because both NA and D0,A decrease with height.

Compared to Murphy20 retrievals, a few findings stand out. Firstly, retrieved qI profiles from the two methods follow each other closely. This is not surprising because both qI profiles are largely constrained by the same ZH observations. Secondly, the retrieved Murphy20 NI is much larger than that from ENCORE-ice. These results suggest that Murphy et al. (2020) have attributed all radar signals to one species like our pristine ice. Due to the smaller size of pristine ice compared to aggregates, the retrieved Murphy20 NI must be much larger than that of ENCORE-ice to make up for having the same qI value. This also explains why NI and qI in Murphy20 retrievals have similar profile shapes. Considering that the observed ρhv is not close to 1, the attribution to single species is likely inappropriate, leading to a large error in ice number concentration, even though qI may seem reasonable. Finally, Deff,mlt, the effective mean diameter using the equivalent melted diameter as the size descriptor, from Murphy20 retrievals tends to be larger than that from ENCORE-ice. This is partly because Murphy et al. (2020) have used denser ice particles; i.e. both the pre-factor and the exponent in their mass–size relationship are slightly larger than our aggregates. Since Deff,mlt depends on the assumed mass–size relationship, Fig. 8l is used for a qualitative comparison only.

Extending the evaluation from one ray to collocated data, a set of thresholds for matching space and time is needed. Since the enhanced area in Fig. 6 is about 15 km wide and 1 km deep, we use this scale as one of our criteria and consider in situ observations and radar gates collocated if their distance is within 15 km in the horizontal and 1 km in the vertical. The time difference threshold for collocation is set to be within 3.5 min because a pair of back-to-back radar RHI scans were performed every 7 min (see Sect. 2). These spatial and temporal thresholds lead to six clusters of radar scans for intercomparison, which comprise 105 rays with a total of 1675 gates. For a given ray, if the root-mean-square difference between the measured and the forward-simulated radar observable is greater than 0.1 dB in ZDR, 0.1 km−1 in KDP, or 0.01 in ρhv, we consider that the retrieval quality for the entire ray is poor and exclude all the retrievals. After this exclusion, 81 rays with 1237 radar gates remain for the evaluation. Most unsuccessful retrievals are likely due to an inappropriate prior. To make the retrieval method work for those unsuccessful cases, we may need to assume priors with different shapes of vertical profiles. Unfortunately, we do not have good knowledge of those shapes and will need to rely on future campaigns to help gather this information by taking frequent multiple-layer flights around the radar site.

Figure 9Time series of (a) flight altitude; (b) temperature; and observed and retrieved (c) total ice number concentration, (d) ice water content, (e) effective mean diameter using the maximum particle dimension as the size descriptor, and (f) effective mean diameter using the equivalent melted diameter as the size descriptor. Retrievals from ENCORE-ice and Murphy20 empirical relationships are denoted by dots, explained by detailed legends. The dots represent the median of retrievals from all collocated gates, and the vertical bars denote the range between the 25th and 75th percentiles. Note that the counting uncertainty in the total ice number concentration in (c) is plotted but too small to see. All calculations are based on the size range of HVPS observations, except Murphy20 retrievals in (c) and (d). For convenience, we index six retrieval clusters from 1 to 6 as shown in (d).


Figure 9 shows the time series of in situ observations and collocated retrievals. We expect column ice crystals at the beginning and very end of the time series because of the measured temperature higher than 10 C. The flight height was maintained at  2 km from 06:30 to 06:40 UTC, suggesting that the missing temperatures due to a data glitch at  06:40 UTC are likely to be about 5 C. During 07:10–08:45 UTC, the temperatures are between 10 and 20 C and likely favour the presence of both dendrite and plate. These expectations about prevalent ice habits are confirmed by visually checking the in situ cloud particle images (see Fig. 1 for examples).

In our retrievals, 40 % of the collocated radar observables are best fit with plate as the pristine-ice habit, 20 % with dendrite, and 40 % with columns. In general, when the cloud particle images were indeed dominated by columns, we have also found that retrievals with columns as the pristine-ice habit provide the best agreement between the measured and forward-simulated radar observables. In the period between 07:00–08:45 UTC when dendrites appeared much more frequently than plates in cloud particle images, our retrievals suggest the opposite because 40 % of best-fit retrievals are associated with plates and only 20 % of best-fit retrievals are associated with dendrite. Therefore, we consider that there remains a large uncertainty in distinguishing plate and dendrites using our retrievals. Note that even with this habit uncertainty, the choice of plate and dendrite does not lead to significantly different retrievals in NI and qI.

Figure 10Box plots of in situ observations and retrievals from ENCORE-ice and Murphy20 for (a) total ice number concentration, (b) ice water content, (c) effective mean diameter using the maximum particle dimension as the size descriptor, and (d) effective mean diameter using the equivalent melted diameter as the size descriptor. The bottom and top of each box represent the 25 % and 75 % quartiles, and the line inside the box represents the median. The whiskers represent the range of data points within 1.5 times the interquartile distance. The included sample sizes for in situ data and radar gates are 347 and 1237, respectively.


The collocated retrievals in Fig. 9c and d show that NI retrieved from ENCORE-ice is approximately on the same order of magnitude as observations and that the retrieved qI values are close to the Nevzorov probe observations. NI and qI from ENCORE-ice generally perform better than those from Murphy20, but they are both overestimated as indicated by the box plots in Fig. 10. The overestimations in the median of NI and qI are 98 % and 44 %, respectively, for ENCORE-ice, and 445 % and 187 %, respectively, for Murphy20. Note that Murphy20 retrievals in NI and qI are based on empirical relationships derived from size bins between zero and infinity, while HVPS-based and ENOCRE-ice-based estimates are derived using HVPS size bins from between 75 and 19 275 µm. This difference in size ranges is not a concern for comparisons of qI, but it contributes to part of the overestimation in NI in Murhpy20 retrievals. Using our retrieved PSD, we have found that the median NI derived from zero- and infinity-sized bins is  3 % larger than that derived from the HVPS size range. This suggests that the difference in the size range for integration calculations is not the main cause for the 445 % overestimation in Murphy20 NI retrievals. Additionally, similarly to in Fig. 8l, Murphy20 Deff,mlt tends to be larger than our Deff,mlt from aggregates by 0.3 mm in the overall median, as shown in Figs. 9f and 10d.

Recall as per Eq. (6) that the effective mean diameter is weighted by size and thus strongly influenced by large particles. When observations sample both pristine ice and aggregates, the combined mean particle size is expected to be close to the size of aggregates, not pristine ice. The results in Fig. 9e generally match this expectation, showing a reasonable agreement between the observed and the combined mean diameter, except three data clusters during the period between 07:30 and 08:45 UTC. For Cluster 5 around 08:45 UTC, the Deff from HVPS has a large variation, ranging between 1.4 and  5 mm with a median of 2.8 mm. This range is on the same order of the median Deff of pristine ice (1.5 mm) and of aggregates ( 5 mm), though the median combined Deff of 4.7 mm is significantly larger than the observed median. For clusters 3 and 4 between 07:30 and 07:45 UTC, the observed effective mean diameter is closer to the retrieved pristine-ice diameter. This unexpected behaviour might suggest a few scenarios, which are discussed in detail next.

The first scenario is that the radar volume might include pristine ice only or aggregates only. However, as shown in Murphy20 retrievals, assuming single species leads to a large error in NP. The observed ρhv is also too low to support this scenario. The second scenario is that the separation between pristine ice and aggregates in our retrieval is inappropriate. To assess this possibility, we tested various combinations of N0 and D0 and found that the following two conditions must be met for the combined size to be close to the size of pristine ice. The first is that NP needs to be at least 1 order of magnitude larger than NA, and the second is that D0,P cannot be much smaller than D0,A. To meet the first condition, let us assume that our retrieved NA is supposed to be 10 times smaller because our retrieved NP from ENCORE-ice is already overestimated compared to the in situ observations and should not be even higher. Under that assumption, D0,A needs to increase by a factor of 1.5 to maintain the same radar reflectivity observation. An increase in D0,A, however, would make the difference between D0,P and D0,A even larger, which violates the second required condition. The combination of reduced NA and increased D0,A by the factors used above would also reduce qI by a factor of 3. Then, to make up for the reduction in qI, one can increase D0,P to remediate the second required condition. This eventually leads to a scenario where two species are alike, as in the first scenario, which is not supported by the observations.

Figure 11Two-dimensional histograms of occurrences of distances in the vertical and horizontal between in situ measurements and radar gates for Cluster 1–6 in (a)(f), respectively. Note that occurrences are counted for all pairs of in situ data points and radar gates. In calculations of retrieval errors, selected in situ data points and radar gates are only used once with equal weights.


The third scenario is that the discrepancy in Deff is due to a sampling issue. Figure 11 shows two-dimensional histograms of occurrences of the vertical and horizontal distance in the collocated in situ and radar dataset. The distance was calculated with respect to the radar gate; i.e. the positive vertical distance represents the flight altitude being higher than the radar gate of interest. Interestingly, for clusters 1, 2, and 6, in situ samples were taken largely at radar scan heights or below. It is likely that both aircraft and radar have sampled the same regime with notable aggregations, which explains why the observed Deff is close to the retrieved Deff of aggregates. In contrast, in situ samples were taken at higher altitudes over the radar scans for clusters 3–5. In these cases, aircraft may have sampled a pristine-ice growth zone aloft, but the radar gates below sampled the subsequent aggregations, which explains why the observed Deff is closer to the retrieved Deff of pristine ice, rather than of aggregates. Further studies using more datasets and retrievals would be needed to assess the third scenario. Overall, when including clusters 1, 2, and 6, the difference between the observed and the retrieved combined median Deff is about 0.55 mm.

5 Summary

We have introduced a new method for retrieving microphysical properties of concurrent pristine ice and snow aggregates from X-band polarimetric radar observations. The radar observables used here include horizontal reflectivity, differential reflectivity, the co-polar correlation coefficient, and specific differential phase shift. The first observable provides constraints on the combined aggregate and pristine-ice population, while the last three observables provide constraints on the partitioning between aggregates and pristine ice, as well as on the ice number concentration and size of pristine ice. The observations are combined with our prior knowledge via an ensemble retrieval framework to find the best estimates of microphysical properties. Since properties of pristine ice and snow aggregates vary significantly in nature, we apply a widely spread prior, and thus the retrieval is mainly dominated by the observations.

Based on the evaluation using synthetic observations, we have found that the current observational uncertainty is sufficient for quantifying properties of pristine ice and snow aggregates. The retrieval was able to reproduce vertical profiles similar to the truth, and the root-mean-square error with respect to the truth is within the retrieval uncertainty. The biases in the combined total ice number concentration, ice water content, and effective mean diameter are all within 5 %. This exercise demonstrates that our retrieval method works well if the prior and the forward models for simulating radar observables are chosen appropriately and are representative of reality. In general, the appropriateness and representativeness of the prior and forward model can be confirmed by examining the agreement between the observations and the forward simulations.

We have also evaluated our retrieval against in situ cloud probe observations taken from a recent field campaign in Chilbolton, UK, which was coordinated to have collocated X-band radar scans and aircraft flights. We analysed a 3 h long case that had 1237 collocated radar gates. Although the period was not particularly long, the aircraft sampled ice particles in temperature zones from 5 to 35 C, allowing us to assess the retrieval performance for cases that are dominated by column, plate, or dendrite. The collocated in situ data have a median number concentration of 5.6 L−1, ice water content of 0.2 g m−3, and 2.7 mm effective mean diameter. Compared to in situ medians, our retrieved total number concentration and ice water content are overestimated by 98 % and 44 %, respectively. This performance is generally better than that from empirical relationships, which has differences of 445 % and 187 % in the total number concentration and ice water content, respectively, with respect to the in situ medians. For the effective mean diameter, the in situ observations agree with our effective mean diameter combined from pristine ice and aggregate in three data clusters with a difference of 0.55 mm. In other clusters, the observed effective mean diameters agree better with the retrieved size of pristine ice, likely because the aircraft sampled pristine-ice growth zones aloft instead of the aggregation zones that the radar sampled. Since planar crystal growth and subsequent aggregation can lead to zones with distinct ice bulk properties, taking frequent aircraft measurements at multiple vertical layers around the radar location would be particularly helpful to improve collocations and allow us to analyse individual rays in more detail.

Currently, our method is designed to work for conditions with a mixture of pristine ice and aggregates. In the presence of rimed particles, the state vector should be expanded to include additional variables that can accommodate and inform the degree of riming, e.g. the riming factor described in Masson et al. (2018), or to include appropriate rimed species explicitly. When triple-frequency measurements are available and can be used to distinguish particle types effectively (e.g. Kneifel et al., 2015; Barrett et al., 2019), such information on particle types can also be incorporated into our method to provide retrievals for off-zenith radar scans that are more challenging for triple-frequency techniques. It is also possible to expand the observation vector with other radar observables at multiple wavelengths, providing further constraints on retrieval if added information exists.

This work is the first step towards quantifying microphysical properties of concurrent ice species, using a framework that considers our prior knowledge and the observational uncertainties. Since we have focused on radar signals with a reduced co-polar correlation coefficient and enhanced differential reflectivity and specific differential phase shift (i.e. cases with a potentially high ice number concentration), the immediate application will be to study dendritic growth zones commonly found in thick stratiform clouds. In particular, the Atmospheric Radiation Measurement (ARM) user facility has operated X-band polarimetric radars at fixed sites since 2011. Comprehensive X-band measurements from various scan strategies are also available from the Biogenic Aerosols – Effects on Clouds and Climate field campaign in Finland, which took place in 2014. These rich datasets will allow us to study formation of new crystals via either primary nucleation or a secondary ice process, their growth into planar crystals and dendrites, and the subsequent aggregations. The retrieved ice properties can be further compared to model simulations to understand what controls the ice number productions.

Appendix A

Radar equations for a single sample volume containing multiple ice particle habits are given as (Jung et al., 2010)



(A6) K DP = 0.18 λ π i = 1 J 0 C k Re S hh , i f - S vv , i f n ( D ) d D ,

where Zh, Zv, and Zhv are in units of mm6 m−3; D is the maximum particle dimension; λ is the radar wavelength; Kw is the dielectric factor of water and Kw2=0.93; and the amplitude scattering matrix elements (S) are in units of millimetres. The vertical bars represent the magnitude of the terms within them, while Re represents the real part of the complex number, and the asterisk indicates its complex conjugate. The index i represents the species existing in the radar volume, and the index J represents the number of species. Note that the amplitude scattering matrix elements in the database are tabulated for various elevation angles, azimuth angles, and habit realizations. To apply these amplitude elements to Eqs. (A1)–(A6), Shhf,b and Svvf,b are first linearly interpolated with respect to elevation angles for radar rays. The interpolated Shhf,b and Svvf,b are then used to calculate the terms in the parentheses for all azimuth angles and habit variations. Because the azimuthal orientation of hydrometeors relative to the radar is random and unknown and because the exact morphological characteristics of these particles at any given time in nature are also unknown, the terms in the parentheses are averaged over azimuth angles and habit realizations, which are represented by the horizontal bar over the parentheses.

Coefficients A, B, C, and Ck are included to account for the effects of canting on the polarimetric radar moments. Following Jung et al. (2010) and Ryzhkov et al. (2011), the canting angle distributions are assumed to be Gaussian, and their effects can be parameterized using the mean and standard deviation of the distribution. Supposing that all oblate species fall with their major axes preferentially oriented in the horizontal plane, the mean canting angle can be set to zero (Ryzhkov et al., 2011). The width of the canting angle distribution is set to 10 for pristine ice crystals and 60 for snow aggregates, similarly to in Ryzhkov et al. (2011) and Matsui et al. (2019). All detailed equations and coefficients can be found in Jung et al. (2010).

Appendix B

Ryzhkov and Zrnic (2019) used a power-law dependence to describe particle density, given as

(B1) ρ = α D e - 1 ,

where the density ρ is in g cm−3, coefficient α is in g cm−2, and De is the equivalent volume diameter. They also assumed an exponential particle size distribution; i.e.

(B2) N D e = N 0 , s e - D e ,

with an intercept N0,s and the exponent . From their Eqs. (11.32) and (11.33) in Ryzhkov and Zrnic (2019), we can calculate these two parameters by



(B5) N 0 , s = q I , emp 0.0003811 α - 0.2 Z 0.6 ,

where the exponent is in mm−1, N0,s is in m−3 mm−1, α is in g cm−2, Z is the radar reflectivity in mm6 m−3, qI,emp is the retrieved ice water content from Eq. (19) in g m−3, and Demp is the retrieved diameter from Eq. (17) in millimetres.

Once N0,s and are known in Eq. (B2), we further convert the size descriptor De to the equivalent melted diameter (denoted as Dmlt) by

(B6) D mlt = α ρ w D e 2 1 3

and then calculate the effective mean diameter Deff,mlt using Eq. (7).

Data availability

FAAM aircraft observations from PICASSO are available at the Centre for Environmental Data Analysis archive (, FAAM, 2018). NxPol radar observations from PICASSO are publicly available via (Bennett, 2020)​​​​​​​. The ice crystal scattering database used to compute radar moments is available in the Atmospheric Radiation Measurement (ARM) user facility Data Archive ( using Data Product Name = Icepart-mod (Lu et al., 2016). Retrievals are available freely at (, Chiu, 2021).

Author contributions

All authors contributed to the work presented here and to manuscript editing. NJK analysed the radar and aircraft observations, coded the algorithm, performed the retrievals, and provided the initial draft. JCC supervised the project, conceptualized the idea, developed the methodology, analysed the retrievals, and revised the manuscript. SB, SSJ, and VC quality controlled the radar observations and contextualized observed polarimetric features. YL provided guidance in constructing the radar forward model. PJvL introduced the most recent development in the iterative stochastic ensemble Kalman approach and assisted in improving the retrieval algorithm. CW provided and contextualized the PICASSO observations. YB assisted in data analysis and retrieval preparation. SO'S provided and quality controlled the in situ cloud probe data.

Competing interests

The contact author has declared that neither they nor their co-authors have any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


Airborne data were obtained using the BAe-146-301 Atmospheric Research Aircraft (ARA) flown by Directflight Ltd and managed by the Facility for Airborne Atmospheric Measurements (FAAM), which is a joint entity of NERC and the Met Office. This work would not be possible without the efforts of the FAAM group and the PICASSO research team led by Jonathan Crosier. We also thank Ryan Neely III and Lindsay Bennett at NCAS, University of Leeds, for collection of the NXPol radar data.

Financial support

This research has been supported by the Office of Science (BER), DOE (grant no. DE‐SC0018930); the European Research Council via the CUNDA project (grant no. 694509); and the Natural Environment Research Council (NERC) (grant no. NE/P012426/1).

Review statement

This paper was edited by Stefan Kneifel and reviewed by two anonymous referees.


Abel, S. J., Cotton, R. J., Barrett, P. A., and Vance, A. K.: A comparison of ice water content measurement techniques on the FAAM BAe-146 aircraft, Atmos. Meas. Tech., 7, 3007–3022,, 2014. 

Aydin, K. and Seliga,T. A.: Radar Polarimetric Backscattering Properties of Conical Graupel, J. Atmos. Sci., 41, 1887–1892,<1887:rpbpoc>;2, 1984. 

Baran, A. J., Connolly, P., and Lee, C.: Testing an ensemble model of cirrus ice crystals using midlatitude in situ estimates of ice water content, volume extinction coefficient and the total solar optical depth, J. Quant. Spectrosc. Ra., 110, 1579–1598,, 2009. 

Barrett, A. I., Westbrook, C. D., Nicol, J. C., and Stein, T. H. M.: Rapid ice aggregation process revealed through triple-wavelength Doppler spectrum radar analysis, Atmos. Chem. Phys., 19, 5753–5769,, 2019. 

Bennett, L.: NCAS mobile X-band radar scan data from 1st November 2016 to 4th June 2018 deployed on long-term observations at the Chilbolton Facility for Atmospheric and Radio Research (CFARR), Hampshire, UK, Centre for Environmental Data Analysis [data set],, 2020. 

Botta, G., Aydin, K., Verlinde, J., Avramov, A. E., Ackerman, A. S., Fridlind, A. M., McFarquhar, G. M., and Wolde, M.: Millimeter wave scattering from ice crystals and their aggregates: Comparing cloud model simulations with X- and Ka-band radar measurements, J. Geophys. Res., 116, D00T04,, 2011. 

Brath, M., Ekelund, R., Eriksson, P., Lemke, O., and Buehler, S. A.: Microwave and submillimeter wave scattering of oriented ice particles, Atmos. Meas. Tech., 13, 2309–2333,, 2020. 

Bringi, V. N. and Chandrasekar, V.: Polarimetric Doppler weather radar: principles and applications, Cambridge University Press, Cambridge, UK, 2001. 

Chiu, J. C.: cjchiurams/PICASSO-IceSnowProperties: v1.0.0 (v1.0.0), Zenodo [data set],, 2021. 

Craven, P. and Wahba, G.: Smoothing noisy data with spline functions: Estimating the correct degree of smoothing by the method of generalized cross-validation, Numer. Math., 31, 377–403, 1979. 

Crosier, J., Bower, K. N., Choularton, T. W., Westbrook, C. D., Connolly, P. J., Cui, Z. Q., Crawford, I. P., Capes, G. L., Coe, H., Dorsey, J. R., Williams, P. I., Illingworth, A. J., Gallagher, M. W., and Blyth, A. M.: Observations of ice multiplication in a weakly convective cell embedded in supercooled mid-level stratus, Atmos. Chem. Phys., 11, 257–273,, 2011. 

Delanoë, J., Protat, A., Testud, J., Bouniol, D., Heymsfield, A. J., Bansemer, A., Brown, P. R. A., and Forbes, R. M.: Statistical properties of the normalized ice particle size distribution, J. Geophys. Res.-Atmos., 110, D10201,, 2005. 

Delanoë, J. M. E., Heymsfield, A. J., Protat, A., Bansemer, A., and Hogan, R. J.: Normalized particle size distribution for remote sensing application, J. Geophys. Res.-Atmos., 119, 4204–4227,, 2014. 

DeMott, P. J., Prenni, A. J., Liu, X., Kreidenweis, S. M., Petters, M. D., Twohy, C. H., Richardson, M., Eidhammer, T., and Rogers, D.: Predicting global atmospheric ice nuclei distributions and their impacts on climate, Proc. Natl. Acad. Sci. USA, 107, 11217–11222, 2010. 

DeMott, P. J., Möhler, O., Stetzer, O., Vali, G., Levin, Z., Petters, M. D., Murakami, M., Leisner, T., Bundke, U., Klein, H., and Kanji, Z. A.: Resurgence in ice nuclei measurement research, B. Am. Meteorol. Soc., 92, 1623–1635,, 2011. 

Erfani, E. and Mitchell, D. L.: Growth of ice particle mass and projected area during riming, Atmos. Chem. Phys., 17, 1241–1257,, 2017. 

Eriksson, P., Ekelund, R., Mendrok, J., Brath, M., Lemke, O., and Buehler, S. A.: A general database of hydrometeor single scattering properties at microwave and sub-millimetre wavelengths, Earth Syst. Sci. Data, 10, 1301–1326,, 2018. 

Evensen, G., Raanes, P. N., Styrodal, A. S., and Hove, J.: Efficient implementation of an iterative Ensemble Smoother for data assimilation and reservoir history matching, Front. Appl. Math. Stat., 5, 47,, 2019. 

Facility for Airborne Atmospheric Measurements (FAAM), Natural Environment Research Council, and Met Office: FAAM C081 PICASSO flight: Airborne atmospheric measurements from core and non-core instrument suites on board the BAE-146 aircraft, Centre for Environmental Data Analysis, available at: (last access: 25 February 2021), 2018. 

Field, P. R. and Heymsfield, A. J.: Importance of snow to global precipitation, Geophys. Res. Lett., 42, 9512–9520,, 2015. 

Field, P. R., Hogan, R. J., Brown, P. R. A., Illingworth, A. J., Choularton, T. W., and Cotton, R. J.: Parametrization of ice-particle size distributions for mid-latitude stratiform cloud, Q. J. Roy. Meteor. Soc., 131, 1997–2017,, 2005. 

Field, P. R., Lawson, R. P., Brown, P. R. A., Lloyd, G., Westbrook, C., Moisseev, D., Miltenberger, A., Nenes, A., Blyth, A., Choularton, T., Connolly, P., Buehl, J., Crosier, J., Cui, Z., Dearden, C., DeMott, P., Flossman, A., Heymsfield, A., Huang, Y., Kalesse, H., Kanji, Z. A., Korolev, A., Kirchgaessner, A., Lasher-Trapp, S., Leisner, T., McFarquhar, G., Phillips, V., Stith, J., and Sullivan, S.: Secondary ice production: current state of the science and recommendations for the future, Meteor. Mon., 58, 7.1–7.20,, 2017. 

Fielding, M. D., Chiu, J. C., Hogan, R. J., and Feingold, G.: A novel ensemble method for retrieving properties of warm cloud in 3-D using ground-based scanning radar and zenith radiances, J. Geophys. Res.-Atmos., 119, 10912–10930, 2014. 

Fielding, M. D., Chiu, J. C., Hogan, R. J., Feingold, G., Eloranta, E., O'Connor, E. J., and Cadeddu, M. P.: Joint retrievals of cloud and drizzle in marine boundary layer clouds using ground-based radar, lidar and zenith radiances, Atmos. Meas. Tech., 8, 2663–2683,, 2015. 

Garrett, T. J., Yuter, S. E., Fallgatter, C., Shkurko, K., Rhodes, S. R., and Endries, J. L.: Orientations and aspect ratios of falling snow, Geophys. Res. Lett., 42, 4617–4622,, 2015. 

Grazioli, J., Lloyd, G., Panziera, L., Hoyle, C. R., Connolly, P. J., Henneberger, J., and Berne, A.: Polarimetric radar and in situ observations of riming and snowfall microphysics during CLACE 2014, Atmos. Chem. Phys., 15, 13787–13802,, 2015. 

Gultepe, I., Heymsfield, A. J., Field, P. R., and Axisa, D.: Ice-phase precipitation, Meteor. Mon., 58, 6.1–6.36,, 2017. 

Heymsfield, A. J., Schmitt, C., Bansemer, A., and Twohy, C. H.: Improved representation of ice particle masses based on observations in natural clouds, J. Atmos. Scie., 67, 3303–3318,, 2010. 

Hogan, R. J., Field, P. R., Illingworth, A. J., Cotton, R. J., and Choularton, T. W.: Properties of embedded convection in warm-frontal mixed-phase cloud from aircraft and polarimetric radar, Q. J. Roy. Meteor. Soc., 128, 451–476,, 2002. 

Hogan, R. J., Tian, L., Brown, P. R. A., Westbrook, C., Heymsfield, A. J., and Eastment, J. D.: Radar scattering from ice aggregates using the horizontally aligned oblate spheroid approximation. J. Appl. Meteorol. Clim., 51, 655–671,, 2012. 

Hong, G., Yang, P., Baum, B. A., Heymsfield, A. J., Weng, F., Liu, Q., Heygster, G., and Buehler, S. A.: Scattering database in the millimeter and submillimeter wave range of 100–1000 GHz for nonspherical ice particles, J. Geophys. Res., 114, D06201,, 2009. 

Hubbert, J. C., Ellis, S. M., Chang, W.-Y., Rutledge, S., and Dixon, M.: Modeling and interpretation of S-Band ice crystal depolarization signatures from data obtained by simultaneously transmitting horizontally and vertically polarized fields, J. Appl. Meteorol. Clim., 53, 1659–1677,, 2014. 

Jiang, Z., Oue, M., Verlinde, J., Clothiaux, E. E., Aydin, K., Botta, G., and Lu, Y.: What can we conclude about the real aspect ratios of ice particle aggregates from two-dimensional images?, J. Appl. Meteorol. Clim., 56, 725–734,, 2017. 

Jung, Y., Xue, M., and Zhang, G.: Simulations of polarimetric radar signatures of a supercell storm using a two-moment bulk microphysics scheme, J. Appl. Meteorol. Clim., 49, 146–163,, 2010. 

Kajikawa, M.: Observation of the falling motion of early snow-flakes. Part II: On the variation of falling velocity, J. Meteorol. Soc. Jpn., 67, 731–738, 1989. 

Keat, W. J. and Westbrook, C. D.: Revealing layers of pristine oriented crystals embedded within deep ice clouds using differential reflectivity and the copolar correlation coefficient, J. Geophys. Res.-Atmos., 122, 11737–11759,, 2017. 

Keat, W. J., Westbrook, C. D., and Illingworth, A. J.: High-precision measurements of the copolar correlation coefficient: Non-Gaussian errors and retrieval of the dispersion parameter μ in rainfall, J. Appl. Meteorol. Clim., 55, 1615–1632,, 2016. 

Kennedy, P. C. and Rutledge, S. A.: S-Band dual-polarization radar observations of winter storms, J. Appl. Meteorol. Clim., 50, 844–858,, 2011. 

Kneifel, S., vonLerber, A., Tiira, J., Moisseev, D., Kollias, P., and Leinonen, J.: Observed relations between snowfall microphysics and triple-frequency radar measurements, J. Geophys. Res.-Atmos., 120, 6034–6055,, 2015. 

Korolev, A., McFarquhar, G., Field, P. R., Franklin, C., Lawson, P., Wang, Z., Williams, E., Abel, S. J., Axisa, D., Borrmann, S., Crosier, J., Fugal, J., Kramer, M., Lohmann, U., Schlenczek, O., Schnaiter, M., and Wendisch, M.: Mixed-phase clouds: progress and challenges, Meteor. Mon., 58, 5.1–5.50,, 2017. 

Korolev, A., Heckman, I., Wolde, M., Ackerman, A. S., Fridlind, A. M., Ladino, L. A., Lawson, R. P., Milbrandt, J., and Williams, E.: A new look at the environmental conditions favorable to secondary ice production, Atmos. Chem. Phys., 20, 1391–1429,, 2020. 

Kumjian, M.: Principles and applications of dual-polarization weather radar. Part I: Description of the polarimetric radar variables, J. Operational Meteor., 1, 226–242,, 2013. 

Kuo, K.-S., Olson, W. S., Johnson, B. T., Grecu, M., Tian, L., Clune, T. L., Aartsen, B. H. V., Heymsfield, A. J., Liao, L., and Meneghini, R.: The Microwave radiative properties of falling snow derived from nonspherical ice particle models. Part I: An extensive database of simulated pristine crystals and aggregate particles, and their scattering properties, J. Appl. Meteorol. Clim., 55, 691–708,, 2016. 

Li, J.–L. F., Forbes, R. M., Waliser, D. E., Stephens, G., and Lee, S.: Characterizing the radiative impacts of precipitating snow in the ECMWF integrated forecast system global model, J. Geophys. Res.-Atmos., 119, 9626–9637,, 2014. 

Liu, G.: A database of microwave single-scattering properties for nonspherical ice particles, B. Am. Meteorol. Soc., 89, 1563–1570,, 2008. 

Lu, Y., Aydin, K., Clothiaux, E. E., and Verlinde, J.: Retrieving cloud ice water content using millimeter- and centimeter-wavelength radar polarimetric observables, J. Appl. Meteorol. Clim., 54, 596–604,, 2015. 

Lu, Y., Jiang, Z., Aydin, K., Verlinde, J., Clothiaux, E. E., and Botta, G.: A polarimetric scattering database for non-spherical ice particles at microwave wavelengths, Atmos. Meas. Tech., 9, 5119–5134,, 2016 (data available at:, last access: 1 October 2018). 

Mason, S. L., Chiu, C. J., Hogan, R. J., Moisseev, D., and Kneifel, S.: Retrievals of riming and snow density from vertically pointing Doppler radars, J. Geophys. Res.-Atmos., 123, 13807–13834,, 2018. 

Matsui, T., Dolan, B., Rutledge, S. A., Tao, W. K., Iguchi, T., Barnum, J., and Lang, S. E.: POLARRIS: A POLArimetric Radar Retrieval and Instrument Simulator, J. Geophys. Res.-Atmos, 124, 4634–4657,, 2019. 

Mitchell, D. L.: Use of mass- and area-dimensional power laws for determining precipitation particle terminal velocities, J. Atmos. Sci., 53, 1710–1723,<1710:uomaad>;2, 1996. 

Mitchell, D. L. and Heymsfield, A. J.: Refinements in the treatment of ice particle terminal velocities, highlighting aggregates, J. Atmos. Sci., 62, 1637–1644,, 2005. 

Moisseev, D. N., Lautaportti, S., Tyynela, J., and Lim, S.: Dual-polarization radar signatures in snowstorms: Role of snowflake aggregation, J. Geophys. Res.-Atmos., 120, 12644–12655,, 2015. 

Morrison, H., van Lier-Walqui, M., Fridlind, A. M., Grabowski, W. W., Harrington, J. Y., Hoose, C., Korolev, A., Kumjian, M. R., Milbrandt, J. A., Pawlowska, H., Posselt, D. J., Prat, O. P., Reimel, K. J., Shima, S., van Diedenhoven, B., and Xue, L.: Confronting the challenge of modelling cloud and precipitation microphysics, J. Adv. Model. Earth Sy., 12, e2019MS001689,, 2020. 

Mulmenstadt, J., Sourdeval, O., Delanoe, J., and Quass, J.: Frequency of occurrence of rain from liquid-, mixed-, and ice-phase clouds derived from A-Train satellite retrievals, Geophys. Res. Lett., 42, 6502–6509,, 2015. 

Murphy, A. M., Ryzhkov, A., and Zhang, P.: Columnar Vertical Profile (CVP) Methodology for validating polarimetric radar retrievals in ice using in situ aircraft measurements, J. Atmos. Ocean. Tech., 37, 1623–1642,, 2020. 

Neely III, R. R., Bennett, L., Blyth, A., Collier, C., Dufton, D., Groves, J., Walker, D., Walden, C., Bradford, J., Brooks, B., Addison, F. I., Nicol, J., and Pickering, B.: The NCAS mobile dual-polarisation Doppler X-band weather radar (NXPol), Atmos. Meas. Tech., 11, 6481–6494,, 2018. 

O'Shea, S., Crosier, J., Dorsey, J., Gallagher, L., Schledewitz, W., Bower, K., Schlenczek, O., Borrmann, S., Cotton, R., Westbrook, C., and Ulanowski, Z.: Characterising optical array particle imaging probes: implications for small-ice-crystal observations, Atmos. Meas. Tech., 14, 1917–1939,, 2021. 

Oue, M., Kollias, P., Ryzhkov, A., and Luke, E. P.: Toward exploring the synergy between cloud radar polarimetry and Doppler spectral analysis in deep cold precipitating systems in the Arctic, J. Geophys. Res.-Atmos., 123, 2797–2815,, 2018. 

Protat, A., Delanoë, J., Bouniol, D., Heymsfield, A. J., Bansemer, A., and Brown, P.: Evaluation of ice water content retrievals from cloud radar reflectivity and temperature using a large airborne in situ microphysical database, J. Appl. Meteorol. Clim., 46, 557–572,, 2007. 

Ryzhkov, A., Pinsky, M., Pokrovsky, A., and Khain, A.: Polarimetric radar observation operator for a cloud model with spectral microphysics, J. Appl. Meteorol. Clim., 50, 873–894,, 2011. 

Ryzhkov, A. V. and Zrnic, D. S.: Polarimetric microphysical retrievals, in Radar Polarimetry for Weather Observations, Springer, 435–464,, 2019. 

Ryzhkov, A. V., Zrnic, D. S., and Gordon, B. A.: Polarimetric method for ice water content determination, J. Appl. Meteorol. Clim., 37, 125–134,<0125:PMFIWC>2.0.CO;2, 1998. 

Schrom, R. S., Kumjian, M. R., and Lu, Y.: Polarimetric radar signatures of dendritic growth zones within Colorado winter storms, J. Appl. Meteorol. Clim., 54, 2365–2388,, 2015. 

Seliga, T. A., Bringi, V. N., and Al-Khatib, H. H.: A preliminary study of comparative measurements of rainfall rate using the differential reflectivity radar technique and a rain gage Network, J. Appl. Meteorol., 20, 1362–1368, 1981. 

Spek, A. L. J., Unal, C. M. H., Moisseev, D. N., Russchenberg, H. W. J., Chandrasekar, V., and Dufournet, Y.: A new technique to categorize and retrieve the microphysical properties of ice particles above the melting layer using radar dual-polarization spectral analysis, J. Atmos. Ocean. Tech., 25, 482–497,, 2008. 

Szyrmer, W. and Zawadzki, I.: Snow studies. part II: average relationship between mass of snowflakes and their terminal fall velocity, J. Atmos. Sci., 67, 3319–3335,, 2010. 

Testud, J., Oury, S., Black, R. A., Amayenc, P., and Dou, X.: The concept of “normalized” distribution to describe raindrop spectra: A tool for cloud physics and cloud remote sensing, J. Appl. Meteorol., 40, 1118–1140,<1118:tcondt>;2, 2001. 

Tiira, J., Moisseev, D. N., von Lerber, A., Ori, D., Tokay, A., Bliven, L. F., and Petersen, W.: Ensemble mean density and its connection to other microphysical properties of falling snow as observed in Southern Finland, Atmos. Meas. Tech., 9, 4825–4841,, 2016. 

Van Leeuwen, P. J.: A consistent interpretation of the stochastic version of the Ensemble Kalman Filter, Q. J. Royal Met. Soc., 146, 2815–2825,, 2020.  

Vivekanandan, J., Martner, B., Politovich, M., and Zhang, G.: Retrieval of atmospheric liquid and ice characteristics using dual-wavelength radar observations, IEEE T. Geosci. Remote, 37, 2325–2334,, 1999. 

Wang, Y. and Chandrasekar, V.: Algorithm for estimation of the specific differential phase, J. Atmos. Ocean. Tech., 26, 2565–2578,, 2009. 

Xu, Y.-L.: Electromagnetic scattering by an aggregate of spheres, Appl. Optics, 34, 4573,, 1995. 

Yurkin, M. A. and Hoekstra, A. G.: The discrete-dipole-approximation code ADDA: Capabilities and known limitations, J. Quant. Spectrosc. Ra., 112, 2234–2247,, 2011. 

Short summary
Ice clouds play a key role in our climate system due to their strong controls on precipitation and the radiation budget. However, it is difficult to characterize co-existing ice species using radar observations. We present a new method that separates the radar signals of pristine ice embedded in snow aggregates and retrieves their respective abundances and sizes for the first time. The ability to provide their quantitative microphysical properties will open up many research opportunities.