Relating optical and microwave grain metrics of snow: The relevance of grain shape

The first sentence could probably be improved. I would start with something like: � “At first order, specific surface area (or optical grain size) is the primary parameter used to simulate snow optical and microwave properties. However, the latter also depend on grain shape....” I would also suggest being more explicit in the results: � l.5-8 : e


Introduction
Use either correlation function or two-point correlation function, but try not to alternate.Reply: We agree.We systematically use the two-point correlation function.P2 l.25 : Picard et al. (2009) do not really mention B and g.They use Monte Carlo ray-tracing on different collections of geometrical shapes instead.Hence for sake of clarity, the reference should be put after "attributed to shape", rather than at the end of the sentence.It might be useful to add The Kokhanovsky and Zege (2004) reference after introducing B and g, because this is where they originate from (at least B).Reply: We agree.The suggested references are placed accordingly.P2 l.27 : "the question remains which..." Reply: We agree.The sentence is changed accordingly.p5 l.17 : in terms of the ?". Reply: We agree.The sentence is changed accordingly.p7 l.22 : g is the asymmetry factor.In fact I would say "(phase function, single scattering albedo, etc)" because g is just computed from the phase function.Reply: We agree.The sentence is changed accordingly.Fig. 2 caption : parameter Reply: Yes.The typo is corrected.p18 l.8 : are calculated ?Reply: We agree, sentence is changed accordingly.P20 l.18 : in the obtaining ?Reply: 'the' is removed to make the sentence correct.P20 l.22 : I do not understand the argumentation.Why should this ratio be constant?And constant for all samples?On the contrary, I would expect λ2 to depend on the samples, because this is a shape parameter.To me, λ2 should be resolution insensitive.I would have expected you used different resolutions for the same sample and check that the retrieved λ2 does not change.Maybe I just did not understand well your point, but it might be useful to rephrase this part of the paragraph.
Reply: The confusion might come from the fact that λ 2 is still a length scale, and not a dimensionless shape parameter.So it is naturally affected by resolution.In general the fraction of any physical length-scale and resolution should be high to sufficiently resolve that length scales by a discrete representation.By showing that the fraction of λ 2 /voxelsize is approximately the same over the data set (and sufficiently large), we can exclude a resolution bias.We will explicitly stress in the discussion that our shape metrics are still length scales in contrast to a perception of a dimensionless parameter.Sentence in 5.3.1 last paragraph is added: "Note that within this definition grain shape is not a dimensionless parameter.With this perception of shape we now…" P21 l 1-2 : is there a problem with the syntax of the sentence ?Yes, sentence restructured to "The present dataset was previously used to study the anisotropic properties of snow ( Löwe et al. 2013).Therefore it is necessary to elaborate on the impact … " P21 l.3 : remove "is", remove "the" We agree.Sentence is adjusted to: "Therefore it is necessary to elaborate on the impact of anisotropy in the present analysis that exclusively involves isotropic two-point correlation functions."p21 l.12 : a corollary is whether anisotropic media can be satisfactorily represented by "equivalent" isotropic media, for microwave and optical properties.This is probably beyond the scope of this paper, but one sentence at the end of this section 5.1.3might be relevant if you have an opinion on this.
Reply: In general it will depend on the quantity of interest if anisotropy would play a role.If the quantity of interest "performs" a similar type of volume or directional averaging the 'isotropic' approach is fine.However, commenting on these type of issues is beyond the scope of the paper and very speculative.P25 l.2 : reference in parentheses We agree.Changed accordingly.p25 l.18 : remove parentheses for reference We agree.Changed accordingly.p25 l.22 : statistically We agree.Changed accordingly.p25 l.24 : an exponential We agree.Changed accordingly.p25 l.29 : correlation We agree.Changed accordingly.p27 l.24 : in Libois et al. (2013) We agree.Changed accordingly.p27 l.31 : although the range of B obtained experimentally is larger than that resulting from Malinka (2014), because the latter implies a shape independent B at weakly absorbing wavelengths, it is worth noting that the actual values are very similar, which suggests that the random two-phase medium is not inconsistent with laboratory and field measurements.Reply: The sentence is changed to include this consistency to "The predicted values for B (Fig. 7) are very similar to the values obtained by experiments (Libois et al., 2014) but show a smaller variation.It should be noted that, as the authors discuss, the correlation between the experimentally obtained B and shape, as defined by Fierz et al. (2009), is statistically not significant and variations might be attributed to shadowing effects relevant at higher densities."P28 l.25 : "length scales () of snow samples, which" Reply: We agree.Changed accordingly.p30 l.11 : involves Reply: We agree.Changed accordingly.p31 l.10 : moments Reply: We agree.Changed accordingly.

Introduction
Linking physical properties of snow to the microstructure always requires to identify appropriate metrics of grain size.In this regard , the two-point correlation function has become a key quantity for the prediction of various properties such as thermal conductivity, permeability and electromagnetic properties of snow (Wiesmann and Mätzler, 1999;Löwe et al., 2013;Calonne et al., 2014b;Löwe and Picard, 2015).The :::::::: two-point correlation function carries, in essence, information about a distribution of relevant sizes in the microstructure.For microwave applications, the analysis of :::::::: two-point correlation functions was already used in the era before micro-computed tomography (µCT), where thin section data and stereology were employed to obtain the required information (Vallese and Kong, 1981;Zurk et al., 1997;Mätzler and Wiesmann, 1999).The recently gained interest in :::::::: two-point correlation functions is mainly driven by available data from µCT, from which the ::::::: two-point : correlation function can be conveniently estimated.The relevance of the two-point correlation function for microwave modeling originates from the connection between its Fourier transform and the scattering phase function in the Born approximation for small scatterers (Mätzler, 1998;Ding et al., 2010;Löwe and Picard, 2015), or the connection to the effective dielectric tensor via depolarization factors (Leinss et al., 2015).
A common practical way to characterize the ::::::: two-point : correlation function is a fit to an exponential, such that the fit parameter, the so called exponential correlation length ⇠, can be used to model the decay of microstructural correlations in snow by a single size parameter.This approach dates back to Debye et al. (1957) in the context of small angle scattering of heterogeneous materials.However the characterization of snow by a single length ⇠ is only an approximation since the occurrence of multiple length scales (Löwe et al., 2011) are known to play a role, in particular to characterize anisotropy (Löwe et al., 2013;Calonne et al., 2014b).Despite this caveat, ⇠ still constitutes the main microstructural parameter for microwave modeling of snow (Proksch et al., 2015a;Pan et al., 2016) when the Microwave Emission Model of layered snowpacks (Wiesmann et al., 1998) is used.
The exponential correlation length is often inferred from measurements of the optical equivalent diameter d opt or, equivalently, from the specific surface area (SSA).This link was established statistically (Mätzler, 2002) leading to the empirical relation where is the ice volume fraction.This relation facilitates using the measured optical diameter as the primary input for microwave modeling (Durand et al., 2008;Proksch et al., 2015b;Tan et al., 2015).However, this link between ⇠ and d opt can only serve as a first approximation.The numerical prefactor in Eq. ( 1) seems to depend on snow type (Mätzler, 2002) which causes a significant scatter in estimating :: the :::::::::: exponential : correlation length from optical diameter.This poses the question which additional size metric captures variations in grain shape and explains the scatter.
The two examples from microwave or optical modeling above reflect the known fact that the optical diameter as a single metric of grain size is not sufficient to characterize the microstructure for many physical properties.It is thus necessary to account for additional grain size metrics which implement the idea of grain shape.A key requirement for potential, new shape metrics is a well-defined geometrical meaning.Present snowpack models (Vionnet et al., 2012;Lehning et al., 2002) contain empirical shape descriptors such as sphericity (Brun et al., 1992).An objective definition of these quantities for arbitrary twophase materials is, however, not possible.New shape metrics should thus ideally seek to replace empirical parameters by an objective, measurable and geometrically comprehensible metrics ::::: metric.
One appealing route to define shape is via curvatures of the ice-air interface because curvatures i) have already been used to comprehend snow metamorphism via mean and Gaussian curvatures (Brzoska et al., 2008;Schleef et al., 2014;Calonne et al., 2014a) ii) are natural quantities to assess shape via deviations from a sphere, very close to the definition of sphericity in Lesaffre et al. (1998) and iii) naturally emerge as higher order terms in the expansion of the ::::::: two-point : correlation function (Torquato, 2002).The latter fact can be used in turn to assess variations of the microwave parameter (⇠) from µCT images which links back to the aforementioned microwave modeling problem.
Another appealing route to define shape is via chord length distributions because they i) naturally implement the idea of size dispersity and ii) have been recently :::::: recently ::::: been : put forward by Malinka (2014) to derive closed-form expressions for the averaged optical properties of a porous medium.Again, the latter fact can in turn be used to assess variations in the optical parameters (g G , B) from µCT images which links back to the aforementioned optical modeling problem.
The motivation of the present paper is to investigate and interconnect these two routes of objectively ::::::::::: (objectively) defining grain shape.First, we will assess the curvature-length in the expansion of the :::::::: two-point : correlation function.We will be guided by the question if and how the well-known statistical relation Eq. ( 1) between the exponential correlation length and the optical diameter can be improved by incorporating curvatures.Second, we will characterize the microstructure in terms of chord length distributions in order to make contact to aspects of shape in snow optics.An interconnection between the two routes can be established by an approximate relation between the :::::::: two-point correlation function and the chord length distribution that was originally suggested in the context of small angle scattering (Méring and Tchoubar, 1968).By means of this approximate relation we establish various statistical links between all involved size metrics, the moments of the chord length distributions, optical diameter, surface areas, curvatures and the exponential correlation length.The established links imply a microstructural connection between geometrical optics and microwave scattering via size dispersity, which constitutes one aspect of grain shape.
The paper is organized as follows.In Section :::::: section 2 we present the theoretical background for the :::::::: two-point correlation function, the chord length distribution, the connection between both quantities and the governing length scales.In Section :::::: section 3 we provide a summary of the µCT image analysis methods.To provide confidence of the interpretation of the curvature metrics derived from the :::::::: two-point correlation function, we present an independent validation of these quantities via the triangulation of the ice-air interface.The results of the statistical models are presented in Section 4 and discussed in Section :::::: section 5.
2 Theoretical background

Two-point correlation function and microwave metrics
The interaction of microwaves with snow are :: is commonly interpreted as scattering at permittivity fluctuations in the microstructure which can be described by the two-point correlation function (Vallese and Kong, 1981;Mätzler, 1998;Ding et al., 2010;Löwe and Picard, 2015).The :::::::: two-point : correlation function can be derived from ::: the spatial distribution of ice and air that is characterized by the ice phase indicator function I(x), which is equal to 1 for a point x in ice and 0 for x in air.From that, a covariance function can be defined which is often referred to as the ::::::: two-point : correlation function (2) In the following we disregard anisotropy by stating that C(r) only depends on the magnitude of r = |r|.To interpret snow with this approach, an average over different coordinate directions must be carried out.
For small arguments r, also rigorous results for the decay of the correlation can be inferred since the expansion of A(r) can be interpreted in terms of geometrical properties of the interface.According to Torquato (2002), the expansion for an isotropic medium reads in terms of the length scales 1 , 2 .The first order term is the slope of the :::::::: two-point : correlation function at the origin and can be expressed in terms of the interfacial area per unit volume s (Debye et al., 1957).The size metric 1 is one of the most fundamental lengths scales for a two-phase medium and commonly referred to as the Porod length in small angle scattering, or correlation length in Mätzler (2002).We will adhere to Porod length here to clearly distinguish 1 from the exponential correlation length ⇠.The metric 1 can be also related to the SSA, defined as the surface area per ice mass (m 2 kg 1 ), or in turn to the equivalent optical diameter d opt of snow via with ⇢ i representing the density of ice.The last equality is obtained when the definition of d opt = 6/⇢ i SSA is inserted (see Mätzler (2002)).
For a two-phase material with a smooth interface, the second order term ⇠ r 2 is missing in the expansion Eq. ( 5) and the next non-zero term is the cubic one with a prefactor 1/ 1 2 2 .Here the length scale 2 has a geometric interpretation in terms of interfacial curvatures and is therfore referred to as the curvature length hereafter.As originally shown by Frisch and Stillinger (1963), the following identity holds in terms :: of ::: the average squared mean curvature H 2 and the averaged Gaussian curvature K.The quantity 2 2 is proportional to the orientationally averaged normal curvature of an interface (Tomita, 1986). (a) Illustration of the chord lengths obtained from an ice sample.The mean chord length is defined as the average length of the green line lengths.A stereological approach (Underwood, 1969) to calculate s is to count the number of blue dots per unit length.The estimation for s mf is given by the red contour.b) Illustration of the ::::::: two-point : correlation function A(r) and the method obtaining an estimate for the Porod length 1 to get s cf by fitting the slope at the origin, and the exponential correlation length ⇠ by fitting A(r) to exp ( r/⇠) over a larger span.

Chord length distributions and optical metrics
In snow optics the microstructural characterization within radiative transfer theory (Kokhanovsky and Zege, 2004) commonly involves a single metric, the optical diameter.An interesting approach for geometrical optics in arbitrary two-phase media was recently put forward by Malinka (2014).Thereby, the microstructure is taken into account by the chord length distribution of : a : medium which can be unambiguously defined for arbitrary two-phase random media (Torquato, 2002).Chord lengths in illustrated in the schematic in Fig. 1a.The chord length distribution p(`) of the ice phase denotes the probability (density) for finding a chord of length `.
In contrast to the Born approximation for microwaves, where the microstructure enters as the Fourier transform of the :::::::: two-point correlation function, the theoretical approach Malinka (2014) ::::::::::::: (Malinka, 2014) relates the key optical quantities (absorption, phase function, asymmetry-factor) to the Laplace transform of the chord length distribution p(`) which is denoted by with Laplace variable z.For small z, the Laplace transform can be approximated by the expansion where µ i denotes the i th moment of the chord length distribution, viz Hence : , within the approach from (Malinka, 2014) ::::::::::::: Malinka (2014) , the optical response of snow can be systematically improved by successively including higher moments of the chord length distribution.According to the theory Malinka (2014), the Laplace transform has to be evaluated at z = ↵, with the absorption coefficient ↵ = 4⇡/ .Here is the wavelength and  the imaginary part of the refractive index of ice.It is generally sufficient (Malinka, 2014) to retain only a few terms in Eq. ( 10).
It is straightforward to show (Underwood, 1969) that the first moment, i.e, the mean chord length µ 1 is given by and thus related to the surface area per unit volume s from Eq. ( 6), or the optical diameter d opt via Eq.( 7).Therefore, in lowest order, the Laplace transform Eq. ( 9) only contains the Porod length or specific surface area of snow.The next order correction involves the second moment µ 2 for which no geometric interpretation has been hitherto given for arbitrary two-phase random media.

Connection between chord lengths and the Porod length and the curvature-length
Following the previous two sections, a link between optical and microwave metrics of snow thus requires to establish a link between :::::::: two-point correlation functions and chord length distributions.To this end we employ a relation between the :::::::: two-point correlation function and chord length distribution that was put forward in the early stages of small angle scattering (Méring and Tchoubar, 1968) to interpret the scattering curve in terms of particle properties.In the present notation the relation can be written as which was also used by Gille (2000).
Although Eq. ( 13) is only valid under certain assumptions which will be discussed in sec.:::::: section : 5, it has already some non-trivial implications that can be exploited for the subsequent analysis.As a first consistency check of the approximation Eq. ( 13), we can compute the first moment of the chord length distribution from Eq. ( 11) for n = 1, by inserting Eq. ( 13) and integrating by parts.This yields µ 1 = µ 1 A(0) which is correct by virtue of Eq. ( 3).As a next step, we aim at an expression for the second moment of the chord length distribution in terms of interfacial curvatures by using Eq. ( 11) for n = 2. Again, inserting Eq. ( 13) and integrating by parts yields Though f is an unknown function here, this link shows that the chord length metric µ 2 must be somehow related to the :::::::: two-point correlation function metrics 1 and 2 .In section :: 4 we will statistically investigate the dependence of f on its arguments.

Data
For the following analysis we used an existing µCT dataset of 3D microstructure images described and used in Löwe et al. (2013) for a thermal conductivity analysis and Löwe and Picard (2015) for a comparison of microwave scattering coefficients.
All samples were classified according to Fierz et al. (2009) as described in the supplement of Löwe et al. (2013).The data set comprises 167 different samples including two time series of isothermal experiments, four time series of temperature gradient metamorphism experiments and a set of 37 individual samples.In total, the set includes 62 samples of depth hoar (DH), 54 of rounded grains (RG), 33 of faceted crystals (FC) 10 of decomposing and fragmented precipitation particles (DF), 5 of melt forms (MF) and 3 of precipitation particles (PP).
3.2 Geometry from :::::::: two-point : correlation functions Obtaining the normalized :::::::: two-point correlation function A(r) from a µCT image can be conveniently done by using the Fast Fourier Transform (FFT) as e.g.described in Newman and Barkema (1999).The FFT is typically used for performance issues to evaluate the convolution integral Eq. ( 2) since direct methods can be very slow.The spatial resolution of the :::::::: two-point correlation function depends on the voxel size of the µCT image which ranges from 10 to 50 µm.
Since the snow samples in the data set are anisotropic (Löwe et al., 2013), the normalized ::::::: two-point : correlation function is first obtained in the x, y and z direction and then averaged arithmetically over the three directions i.e, A(r) = (A x (r) + A y (r) + A z (r)) /3.From the normalized :::::::: two-point correlation function two types of parameter fittings are performed.First, the exponential correlation length ⇠ is obtained by fitting the µCT data to the exponential form Eq. ( 4).Technically, we estimated the inverse parameter k by least-squares optimization of the model A(r) = exp( kr) to the data in a fixed range of 0 < r < 50 .An illustration of this method is shown in Fig. 1b.In the following we denote by ⇠ the inverse of the optimal fit parameter ⇠ := 1/k.
Second ::::::: Secondly, we estimated the expansion parameters 1 and 2 of the :::::::: two-point : correlation function by a least-squares regression to the expansion Eq. ( 5).Technically, we fitted ) in the fixed range of 0 < r < 3 which determines the derivatives at the origin.We denote by cf 1 and cf 2 the inverse of the optimal fit parameters cf 1 := 1/k 1 and The superscript is added to discern these :::::::: two-point correlation function based estimates from those presented in the next section for a validation.The influence of resolution and anisotropy to the estimates of 1 and 2 is discussed in section 5.

Geometry from triangulations
To confirm the geometrical interpretation of cf 1 and cf 2 we use an alternative and independent method to estimate these parameters by measuring the surface area and the local curvatures with a VTK-based image analysis as described in Krol and Löwe (2016).In short, a triangulated ice-air interface is obtained by applying the VTKContour filter.After this step, the interface still resembles the underlying voxel structure.Therefore, in a second step the triangulated interface is smoothed by applying the VTKSmoothing filter which involves a smoothing parameter S which is the number of iterations a Laplacian smoothing on a mesh is repeated.For further details we refer to Krol and Löwe (2016).

Accuracy of surface area and curvatures estimates
The measured total surface area is obtained by integrating (summing) the surface area of the triangles over the surface and the estimate vtk 1 which naturally depends on the smoothing parameter.A comparison of the triangulation and the :::::::: two-point correlation function based length scale is shown in Fig. 2 (middle row).A higher value of the smoothing parameter implies a lower surface area s (caused by shrinking of the enclosed volume upon smoothing) and in turn higher estimates for vtk 1 .Using higher smoothing also results in a higher variance in the data.This is likely due to filtering of small perturbations in the surface causing the individual samples to react differently.
It is illustrative to note that even without smoothing for S = 0 the obtained triangulated surface is still different from the voxel surface s mf , which is obtained by the union of ice-air transition faces in the voxel based image (as illustrated by the red contour in Fig. 1a).The quantity s mf is one of the four Minkowski functionals and can be computed by standard counting algorithms (Michielsen and Raedt, 2001).For isotropic systems, and statistically representative samples, the relation between the surface obtained from the :::::::: two-point correlation function s cf = 4 (1 )/ cf 1 and the Minkowski functionals is known to be s cf = 2s mf /3 as discussed in Torquato (2002, p. 290).correlated : with the size of the structures and the resolution.In the end, we chose a smoothing parameter S = 200 that is, on average, acceptable for all involved samples.
Overall, the comparison provides reasonable confidence that the geometrical interpretation of the :::::::: two-point correlation function parameters is correct, though uncertainties inherent to the smoothing operations must be acknowledged.In the following we solely use the quantities derived from the :::::::: two-point correlation function, viz. 1 = cf 1 and 2 = cf 2 where the superscripts are omitted for brevity.

Chord length distribution
To compute the ice chord length distribution from the binary images, all linear lines through the sample in all three Cartesian directions = x, y, z are considered and all ice chords were measured and binned to obtain direction dependent counting densities n (`).Here n x (`) denotes the total number of chords in x direction which have length `.For a binary CT image, can take integer values 0 < `< L x which are restricted by the sample size L x = N x and the voxel size of the image.The mean chord length and other moments µ i are then computed from `in (`). (15)

Statistical models
The main part of the following analysis comprises statistical relations between the length scales derived from the chord length distribution and the ::::::: two-point : correlation function in section 2. In total, we will consider a few statistical models that first relate the exponential correlation length ⇠ and µ 2 to the geometrical length scales 1 and 2 and second, relate ⇠ to µ 1 and µ 2 .We will start with a one-parameter statistical model and compare the results to the two parameter models.We will assess and compare the quality of the fits with the adjusted correlation coefficient R 2 .

Relating exponential correlation length to the Porod length and curvature-length
As a starting point for the statistical analysis we revisit the empirical relation which is equivalent to Eq. ( 1) by virtue of Eq. ( 7), as suggested by Mätzler (2002).To this end we fitted ⇠ and 1 and obtained an average slope of 0.79 with a correlation coefficient of R 2 = 0.733, shown by the green dashed line in Fig. 3a.In the next step we fitted the same data to include an intercept parameter Here the adjusted correlation coefficient, accounting for the inclusion of extra parameters, is R 2 = 0.731 and the parameters are given by a 0 = 5.93 ⇥ 10 2 mm, a 1 = 0.794, with very low p-values (p < 5 ⇥ 10 4 ) for the intercept and the slope ensuring the significance of the parameters of the fit.The order of magnitude of the intercept a 0 is negligible.To understand the remaining scatter we have plotted the residuals ⇠ (a 0 + a 1 1 ) versus the curvature-length 2 as shown in Fig. 3b.The correlation coefficient is given by R 2 = 0.644 and suggest ::::::: suggests : that including the curvature lengths can improve Eq. ( 17).For an overview, this and all other statistical models will be listed in Table 1.
In the next step we include the curvature-length 2 where we fitted the exponential correlation length ⇠ to the model The results are shown in Fig. 3c.Here we find an improvement compared to Eq. ( 17).The correlation coefficient is R 2 = 0.922 and the fit parameters are given by b  9).Qualitatively, the characteristic form (i.e, single maximum), the location of the maximum, and the width of the distribution are correctly predicted by Eq. ( 13).On the other hand, there are obvious shortcomings, such as the oscillatory tail for the RG example when the chord length distribution is derived via Eq.( 15).We will revisit these characteristics in the discussion.

Relating the second moment of the chord length distribution to the Porod length and the curvature-length
Using the previous results we can derive an approximate relation between the second moment of the chord length distribution and the interfacial curvatures.To motivate a statistical model, we start from Eq. ( 14), We investigate the dependency of the function f on parameters 1 , 2 and of this expression by successively including 1 , 2 and in a statistical model.In a first step we approximate f by a statistical model including only The optimal parameters for model Eq. ( 20) are l 0 = 2.40 ⇥ 10 2 mm and l 1 = 1.25, with negligible p values and a correlation coefficient of R 2 = 0.898.The results are shown in Fig. 5a.In view of the inclusion of the curvature-length 2 , we analyzed the residuals of the previous statistical model and plotted them as a function of 2 (Fig. 5b).The correlation coefficient (R 2 = 0.295) is small but including 2 in the analysis further improves the fit.The respective statistical model yields optimal parameters n 0 = 3.95 ⇥ 10 3 mm, n 1 = 1.50 and n 2 = 2.46 ⇥ 10 1 with a correlation coefficient R 2 = 0.949.The p-value for the intercept n 0 is 0.36.For n 1 and n 2 the p-values are again very low.We have heuristically found a possibility of improving Eq. ( 21) even further.This was achieved by including a factor (1 ) on the left-hand side.More precisely, we tried as a statistical model.Here the optimal parameters are q 0 = 1.23 ⇥ 10 2 mm, q 1 = 1.03, :::::::::::::::::::: q 0 = 1.23 ⇥ 10 2 mm, ::::::::: q 1 = 1.03, and q 2 = 1.98⇥10 1 .The p-values for all coefficients are negligible and the correlation coefficient is R 2 = 0.980.The results are shown in Fig. 5c.

Relating microwave metrics and optical metrics
In the previous sections we found a statistical relation between the exponential correlation length ⇠ and the geometrical lengths    20) versus the curvature-length scale parameter 2, c) the statistical model predicting (1 )µ2/2µ1 (see Eq. ( 22)) depending on the Porod length 1 and the curvature-length 2. and 2 on the other hand.Both findings can be recast into a direct connection between the moments of the chord lengths µ 1 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35  23) that depends on the first and second moment of the chord length distribution, µ1 and µ2.
The summary of all models is given in Table 1.To ensure that the inclusion of an additional parameter : , : e.g. by going from model Eq. ( 17) to model Eq. ( 18), is indeed an improvement, we have employed the Akaike information criterion (AIC) (Akaike, 1998).The AIC measure allows to discern if the improvement of the correlation coefficient is trivially caused by an increasing number of fit parameters or an actual improvement on the likelihood of the fit due to the relevance of the added parameters.Absolute AIC-measures have no direct meaning, however a decrease of at least 2k between two models, where k is the number of extra parameters, implies a statistical improvement.For our case k = 1 the difference in the AIC-measure between Eq. ( 17) and Eq. ( 18) is 177 confirming the statistical relevance significance of 2 .

Shape factors g G and B
As an application of the values obtained for the moments of the chord length distribution we can now compute the "shape diagram" of the optical parameters (g G , B) suggested in Libois et al. (2013) derived from (Malinka, 2014, Eq. 60), and Eq. ( A4).
The results depend on the value of the Laplace transform at the absorption coefficient ↵, and thus on wavelengths.For most wavelengths in the visible and near infrared regime ↵µ 1 ⌧ 1 is small and therefore the Laplace transform Eq. ( 9) can be ap-  Warren and Brandt, 2008), the first order, the fraction of the first and second order of Eq. ( 10), and the obtained estimates for B and g G averaged over all snowsamples, including the standard deviation .proximated by a few terms in the expansion Eq. ( 10).Taking typical values for ↵ allows us to estimate the relative importance ↵µ 2 /2µ 1 of the second-order term compared to the first-order term in the expansion Eq. ( 10).These values are obtained by using the values for  provided by Warren and Brandt (2008).The first order ↵µ 1 and ratio ↵µ 2 /2µ 1 is :: are : calculated for typical wavelengths and shown in Table 2.The values and standard deviations denote averages taken over all samples.Wavelengths are selected to match common optical methods, namely 0.9 µm (Matzl and Schneebeli, 2006), 1.31 µm (Arnaud et al., 2011), and the SWIR wavelengths 1.63 µm, 1.74 µm and 2.26 µm used by Domine et al. (2006).We added the wavelength 2.00 µm, which is not used by any instrument, but has the highest value for ↵ in this range.Note that for this wavelength ↵µ 1 is not small and the expansion of the Laplace transform, Eq. 10, likely not a good approximation.The standard deviations are high as a result of the variations due to grain type.The lowest values of µ 2 /2µ 1 ::::::: are found for fresh snow (PP) and highest for depth hoar (DH) and melt forms (MF).

10
The values in Fig. 7 for g G and B are computed for wavelength 1.3µm and shown as a scatter plot of B versus 1 g G similar to Libois et al. (2013) ) 2 [0.2, 0.5] obtained by ray-tracing simulations for different geometrical shapes (Libois et al., 2013).The variations of the values for different snow types is however very small.To complete the analysis we have computed g G and B for higher absorbing wavelengths for which the shape signature might be higher, but the expansion of Eq. ( 10), less reliable.The results are averaged over all snow samples and included in Table 2.

Methodology
Before turning to the discussion of physical implications of the results, we first address methodological details.Retrieving parameters from µCT images must be taken with care.In addition to the uncertainties related to filtering and segmentation pointed out by Hagenmuller et al. (2016), the present method also requires to discuss the interface-smoothing for the validation of 1 and 2 , the image resolution, and the anisotropy of the samples.

Geometrical interpretation
The present analysis and cross-validation of the curvature metric imposes requirements on the smoothness of the interface.
The subtle influence of the smoothing parameter on the surface area s and averaged mean and Gaussian curvatures H and K is apparent from Fig. 2. Naturally, H 2 is most sensitive to smoothing.We found a competing performance of 1 and 2 with the smoothing parameter when comparing the triangulation based estimates with the ::::::: two-point : correlation function based values.The agreement for the surface area seems to be best with smoothing parameter S = 50.In contrast, more smoothing is required to obtain an agreement for the curvature-length.This higher sensitivity on the smoothing parameter is reasonable, since curvatures are defined by surface gradients which are more sensitive to a smooth mesh representation than the surface area.The competing behavior is caused by the smoothing filter, which neither preserves the volume nor the surface area of the enclosed ice upon smoothing iterations.This causes the drop in agreement for 1 in Fig. 2 (left, middle) with increased smoothing.As a remedy, more sophisticated smoothing filters could be used which, for example, ensure the conservation of the enclosed volume (Kuprat et al., 2001).Such problems could be partly avoided by computing normal vector fields and curvatures directly from voxel-based distance maps (Flin et al., 2005).A detailed comparison of all these different methods however, is beyond the scope of this paper.In contrast to 1 and 2 , the interpretation of first and second moments of the chord length distribution, µ 1 and µ 2 , is rather straightforward, where µ 1 is directly related to the optical diameter d opt , and µ 2 is a measure of the variations of this size metric.

Resolution
Resolution plays an important role in the obtaining estimates for 1 and 2 .For a µCT measurement the resolution is commonly chosen appropriately depending on snow type.While fresh Snow (PP) is typically reconstructed with 10µm voxel size, melt forms (MF) and larger particles have larger voxel sizes of 35µm or 54µm.Since we have obtained 1 and 2 with two independent methods that agree reasonably well we conclude that the resolution is generally sufficient to estimate the involved length scales.To further confirm that that there is no remaining bias with resolution we assessed the ratio 2 /voxelsize.Ideally this would be constant for all samples, implying that 2 is equally well resolved for all snow samples.For our data, this this ratio is 9.8 with a standard deviation of 2.6.The correlation coefficient with the voxel size is R 2 = .2,which implies that there is a slight dependence on resolution.A systematic assessment is however difficult since snow types and grain sizes are not equally distributed over the resolution.
The image resolution plays another important role in the interpretation of the expansion of the :::::::: two-point correlation function.
As pointed out by Torquato (2002), a missing r 2 term is generally equivalent to a smooth interface while discontinuities, like sharp edges, would lead to a second order term.Fresh snow and depth hoar crystals are known to have these discontinuities, at least visually.But it remains questionable if these features can be detected objectively at the micrometer scale from image analysis.In an image, discontinuities are always smeared out, virtually contributing to the third order term.

Anisotropy
The present data set of snow samples embodies a large number of anisotropic samples, which was specifically the subject of Löwe et al. (2013) the data is based on.It is thus ::::: dataset :::: was ::::::::: previously :::: used :: to ::::: study ::: the :::::::::: anisotropic ::::::::: properties :: of ::::: snow :::::::::::::::: (Löwe et al., 2013) .::::::::: Therefore :: it :: is : necessary to elaborate :: on : the impact of anisotropy on : in : the present analysis which is exclusively involves isotropic ::::::: two-point : correlation functions.It is important to note that the our analysis does not assume isotropy, but it rather includes the orientational averaging in the three Cartesian directions as a part of the method.Such a procedure is principally valid for arbitrary samples.Moreover, also the geometrical interpretation of the quantities remains valid.This was rigorously shown for 1 Berryman (1998) which relates the slope of the ::::::: two-point : correlation function at the origin for arbitrary anisotropic structures after orientational averaging to the surface area per unit volume s.Though we did not find a mathematical proof for the corresponding statement for 2 , the agreement of cf 2 (obtained from the :::::::: two-point correlation function, orientationally averaged) with vtk 2 (obtained from direct computation of the interfacial curvatures) strongly suggests its validity.In addition, we assessed that the residuals between vtk 2 (where anisotropy does not play a role) and cf 2 are not correlated with anisotropy (R 2 = .026).Overall, we are confident that the method can be applied to arbitrary anisotropic samples to provide orientationally averaged length scales with the correct geometric interpretation with acceptable uncertainties due to image resolution.

Linking size metrics in snow
Accepting the methodological uncertainties, we shall now discuss our findings of the statistical analysis and their relevance for the interpretation of snow microstructure.

Including size dispersity to estimate the exponential correlation length
By construction, the exponential correlation length ⇠ must be understood as a proxy to characterize the entire :::::::: two-point : correlation function with a single length scale.This single length scale contains signatures of both, ; : properties that dominate the behavior of the :::::::: two-point : correlation function for small arguments ( 1 and 2 ) and other properties that dominate the tail-behavior of the :::::::: two-point correlation function for large arguments.
To discuss the statistical relations we found we will start with recovering Mäzler :::::: Mätzler's model (Mätzler, 2002).This statistical model covers a relation between the exponential correlation length and the optical grain size, or in their nomenclature: the correlation length.Mäzler :::::: Mätzler's model predicts the slope to be a 1 = 0.75, which is an average of a 1 = 0.8 for depth hoar and a 1 = 0.6 for other snow types.This is consistent with our finding a 1 = 0.79 since we have many depth hoar samples in the data set, suggesting that grain shape has a direct influence on the statistical relation.This influence was made quantitative by including the curvature-length to the statistical analysis, resulting in the statistical model Eq. ( 18) (Fig. 3c).The quantitative improvement on the statistical model Eq. ( 16) by using Eq. ( 18) is given by the increase in the correlation coefficient from R 2 = 0.733 to R 2 = 0.922.
In addition we established a new statistical relation Eq. ( 23) between ⇠ and the moments of the chord length distribution, µ 1 and µ 2 .This model performs even better when the correlation coefficient R 2 = 0.985 is taken as a quality measure.We confirmed that the inclusion of an additional parameter in Eq. ( 18) and Eq. ( 23) indeed improves on eq. ( 16), by employing the Akaike information criterion (AIC) measure (Akaike, 1998).
All statistical models showing improvements of ::::::: proposed ::::::::: statistical :::::: models ::::: show ::: an ::::::::::: improvement ::: to ::: Eq. : (1) indicate :::::::: indicating : that at least two different length scales 1 and 2 or µ 1 and µ 2 are required to obtain a reasonable prediction of the exponential correlation length.While 1 and µ 1 are both trivially related to the optical radius via Eq.( 7) and Eq. ( 12), the two other size metrics µ 2 or 2 are the origin of performance increase.This seems surprising at first sight.Why should local aspects of the interface ( 1 and 2 ) determine the non-local decay of structural correlations (⇠)?To illustrate our explanation for this finding, we resort to a particle picture and consider a dense, random packing of monodisperse hard spheres.For such a packing, the particle "shape" is trivial and fully determined by the sphere diameter d, which determines the slope of the :::::::: two-point correlation function at the origin.However, also particle positions and thus the decay of correlations is fixed by d.This becomes obvious from the representation C(r) = nv int (r) + n 2 v int (r) ⇤ h(r) for the :::::::: two-point correlation function for such a system at number density n (Löwe and Picard, 2015).In this representation, the spherical intersection volume v int and the statistics of particle positions h(r) both depend on d.Now imagine that each sphere is deformed by a hypothetical, volume-conserving re-shape operation to an irregular, non-convex particle, which is still located at the center of the original sphere.Due to re-shaping, the parameter H 2 would increase.After the re-shape, neighboring particles would overlap (on average), since their maximum extension must have been increased compared to the sphere diameter.To recover a non-overlapping configuration, all particle positions must be dilated.The latter, however, also affects the tail of the :::::::: two-point correlation function.This is exactly what we observe: the "shape of structural units" in snow, as exemplified by H 2 is always correlated with the "position of the structural units" in space.We note that this particle analogy has clear limitations and only serves here to illustrate the rather abstract statistical relations between different length scales.Snow remains a bi-continuous material where individual particles cannot be distinguished.
Overall, we conclude that both, 2 or µ 2 can be used to significantly improve estimates of ⇠ when compared to optical diameter based estimates.

Linking moments of the chord length distributions to Porod and curvature-length
Hitherto no geometrical interpretation for the second moment µ 2 of the chord length distribution was known.Our results suggest an empirical relation, Eq. ( 22), that involves the two geometrical length scales 1 and 2 .In the following we provide supporting arguments for the link between µ 2 and 1 and 2 by discussing the relation Eq. ( 13) between the chord length distribution and the :::::::: two-point correlation function.
The relation Eq. ( 13) was originally raised in the context of small angle scattering long time ago (Méring and Tchoubar, 1968) and later revisited e.g. by Levitz and Tchoubar (1992), revealing two different approximation steps.A first simplification comes from the assumption that consecutive chords on the random ray in Fig. 1 are statistically independent.This issue has been discussed in detail also by Roberts and Torquato (1999), who established an exact relation between the Laplace transforms of the ::::::: two-point : correlation function, the chord length distribution, and a surface-void correlation function based on this assumption.Their results however show that for level-cut Gaussian random fields, where this assumption is violated, the prediction of the chord length distribution can be still very accurate.This indicates that assuming independent chords is per se not a serious limitation.Secondly, Eq. ( 13) is actually an approximation for dilute systems which is generally not valid for snow.
To test the range of validity of the relation ( 13) for snow, we have taken three samples and computed the chord length distribution directly to compare them to the prediction of Eq. ( 13) as shown in Fig. 4.An obvious drawback of Eq. ( 13) can be seen for the rounded grains (RG) sample.Due to the quasi-oscillations in the :::::::: two-point correlation function (cf.Löwe et al. (2011)), A(`) and its second derivative assume negative values, which would imply negative values for p(r) ::: p(`) : via Eq.( 13).

Grain shape for microwave modeling
Thus far, the exponential correlation length ⇠ as a key parameter for MEMLS based microwave modeling (MEMLS) was mainly predicted from the optical diameter.Our conclusions from section 5.2.1 could now be restated: The inclusion of a grain shape parameter, 2 or µ 2 improves the prediction of the exponential correlation length significantly.Or, according to the conclusion from the previous section, one may alternatively restate that size dispersity has an influence on microwave properties.This is known from other models than MEMLS, where an influence of polydispersity on the effective grain scaling parameter within DMRT-ML microwave modeling was found Roy et al. (2013) :::::::::::::: (Roy et al., 2013) .
This equivalence of shape and size dispersity at the level of :::::::: two-point correlation functions can be further illustrated by an interesting example.Consider a microstructure of polydisperse spherical particles.The definition of grain shape from the classification (Fierz et al., 2009) would assign a spherical shape to this microstructure, while the averaged squared mean curvature H 2 would instead vary depending on the variance of particle radii.As pointed out by Tomita (1986), for low density, such a system of polydisperse spherical particles can always be mapped uniquely onto an assembly of monodisperse but irregularly shaped particles by solving an integral equation, if only the :::::::: two-point correlation function is considered.Shape can be equivalent to polydispersity, and snow types which are visually very different might still have very similar physical properties.This example also explains why the objective size dispersity parameters 2 or µ 2 cannot be mapped onto the classical definition of grain type from Fierz et al. (2009).

Grain shape in geometrical optics
Finally, we turn to the implications of size dispersity or grain shape on geometrical optics within the scope of (Malinka, 2014) :::::::::::: Malinka (2014 on chord length distributions. As pointed out by (Malinka, 2014) ::::::::::::: Malinka (2014) , if consecutive chords were statistical ::::::::: statistically : independent i.e. a Markovian process, then the obtained distribution would be :: an exponential, and all optical properties solely determined by the optical diameter (or µ 1 ).To quantify the deviation from an exponential chord length distributions we calculated the fraction µ 2 /2µ 2 1 which is unity for a exponential chord length distribution.This fraction is on average 0.75 for rounded grains (RG), 0.76 for melt forms (MF), 0.77 for precipitation particles (PP) and defragmented particles (DF), 0.79 for faceted crystals (FC) and the closest value to unity is 0.876 for depth hoar (DH).This implies that the chord length distribution for depth hoar is closest to an exponential, which can be visually confirmed by Fig. 4. We reach a similar conclusion for the :::::::: two-point correlation function where 1 is already a fairly good predictor for the exponential corrrelation :::::::: correlation : length when depth hoar is considered (see Fig. 3)a).But due to the deviations from an exponential ::::::::: distribution, an influence of shape via µ 2 on the optical properties would be expected from :::::::: according :: to : Malinka (2014).
Using the chord length distributions we were able to calculate the shape factors B and g G from Malinka (2014) and  2013) ) in Fig. 7 was obtained for wavelength 1.3 µm where the Laplace transform Eq. ( 10) can be approximated by the first and second order.The variations of the absolute values for B, g G shown in Fig. 7 predominantly stem from corrections which are linear in µ 1 (by virtue of (A5)), while the small, scattered deviations from a perfect straight line are caused by µ 2 .If B and g G were evaluated for wavelength 0.9 µm, the influence of µ 2 would be even smaller.Our results show that the values for B and g G are exactly within the range that is suggested by ray-tracing simulations for various geometrical shapes for a wavelength of 0.9µm Libois et al. (2013) :::::::::::::::: (Libois et al., 2013) , but show a much smaller variation over the entire set of snow samples.Comparing our results to ray-tracing of geometrical shapes is however not straightforward, since the 3D microstructures cannot be mapped on an ensemble of regular geometrical objects.
Overall, our analysis indicates a smaller variation of optical properties with shape via µ 2 according to Malinka (2014) ::::::::::::: Malinka (2014) wh compared other methods.We can only hypothesize potential origins which are connected to the present analysis.A crucial assumption made in the geometrical optics framework (Malinka, 2014) is the statistical independence of the chord length and the consecutive ice-air incidence angle for a ray which passes through a grain.Such an assumption might be progressively violated for lower absorption where a higher number of internal reflections in fact probes this assumption more often.Hence the true effect of shape on B and g G might be still more pronounced as captured :::: more :::::::::: pronounced :::: than :::::::: predicted : by size dispersity via µ 2 within (Malinka, 2014).Further details on the discrepancies between measurements, simulations and theory remain to be elucidated by combining tomography imaging and shape analysis together with optical measurements and ray-tracing simulations in the future.

Conclusions
We have analyzed different microstructural length scales ( 1 , 2 and µ 1 , µ 2 ) :: of :::: snow ::::::: samples : which were derived from the :::::::: two-point correlation function and chord length distribution, respectively.All length scales have a well-defined geometrical meaning.While the first order quantities (µ 1 , 1 ) are both related to the mean size (optical equivalent diameter), their higher order counterparts ( 2 , µ 2 ) are objective measures of size dispersity present in the snow microstructure.For the ::::::: two-point : correlation function, the length scale 2 is essentially determined by the second moment of the mean curvature distribution.For the chord lengths, µ 2 is the second moment of the chord length distribution.Both quantities naturally extend the concept of mean grain size as covered by the optical equivalent diameter.The statistical relation established between ( 1 , 2 , µ 1 , µ 2 ) indicates that practically :: in ::::::: practice the two measures of size dispersity can be used interchangeably.
32 and b 2 = 3.85 ⇥ 10 1 .The p-values are very small for all coefficients b i .The order of magnitude of the improvement can already be roughly estimated from the ratio of the prefactors b 1 and b 2 .4.2 Connection between chord length distributions and :::::::: two-point correlation functions To relate the chord length metrics to the Porod length and the curvature-length, we first assess :::::: assessed : the relation between the chord length distribution p(`) and the ::::::: two-point : correlation function A(`) as suggested by Eq. (13).To this end we compared the chord length distribution obtained directly from the µCT image (cf.section 3.5) with the prediction of Eq. (13) via the :::::::: two-point correlation function for a few examples of different snow types.The results are shown in Fig. 4. The selected snow samples are the same as those used in Löwe and Picard (2015, Fig. 8 and Fig.

Figure 3 .
Figure 3. Scatter plots of a) the exponential correlation length ⇠ versus the Porod length 1.A linear fit is plotted in green.Additionally the prediction of Eq. (16) (MM) is plotted in red.b) The residuals of ⇠ and the statistical model Eq.(17), versus the curvature-length 2. c) The statistical model Eq.(18) predicting ⇠ depending on the Porod length 1 and the curvature-length 2.

Figure 4 .
Figure 4. Comparison of the chord length distributions computed from Eq. (13) (symbols) and by direct analysis of the µCT data (solid-line) for three examples of snow types (PP, RG and DH).
1 and 2 on one hand and a relation between the first and second moment of the chord length distribution (µ 1 and µ 2 ) and 1

Figure 6 .
Figure 6.Scatterplot of the exponential correlation length ⇠ versus the statistical model Eq.(23) that depends on the first and second moment

Figure 7 .
Figure7.Scatterplot of the asymmetry factor g G and the optical shape factor B evaluated for :: the : refractive index at wavelength = 1.3 µm.

Table 1 .
Summary Statistical Models

Table 2 .
Determination of the absorption coefficient ↵ (