Measurements of ice crystal fluxes from the surface at a mountain top site

New observations of anomalously high cloud ice crystal concentrations at the Jungfraujoch research station (Switzerland, 3.5 km a.s.l.) are presented. High-resolution measurements of these ice crystals using a high-speed 2D imaging cloud particle spectrometer confirm that the concentrations far exceed those expected from any known primary ice production mechanisms and are at temperatures well below those for known secondary ice production processes to contribute. A detailed 10 analysis of the cloud ice crystal sizes and habits with respect to wind speed and growth regimes confirms, as hypothesised in previous studies, that their origin is not due to blowing snow or related wind speed influences. As first hypothesised by Vali et al. (Vali et al., 2012) using remote sensing and by Lloyd et al., (2015) with in situ measurements, the most likely explanation is due to a strong surface source generated by the interaction of turbulent deposition of supercooled droplets to fragile icecovered snow surfaces. This process enhances the detachment of crystal fragments wherein the smaller size mode is turbulently 15 re-suspended even at low wind speeds below expected blowing snow thresholds. These then continue to grow, adding significantly to the ice crystal number concentrations whose size and habit is determined by the transport time between the ice crystal source and measurement location and liquid water profile within the cloud. We confirm, using eddy covariance measurements of ice particle number fluxes, that the likely source is significantly far upwind to preclude flow distortion effects such that the source plume has homogenised by the time they are measured at the mountain top summit. 20


Introduction
It has been previously reported that very high ice crystal number concentrations can be observed during cloud events in alpine regions (Lloyd et al., 2015). These events cannot be explained by known number concentrations of ice nucleating particles (INPs) and known primary ice mechanisms (e.g. Conen et al. 2015). In some cases the responsible process has been identified 25 as a consequence of so-called secondary ice processes (SIP), such as the Hallet-Mossop process (Hallett et al., 1974) and aeolian generated fragmentation and resuspension of large surface ice particles (blowing snow), Vionnet et al. (Vionnet et al., 2013) and Geerts et al. (2015). In other cases, it was not possible to link high ice number concentrations to these processes mentioned. In this chapter, we attempt use of the well-documented Eddy Covariance (EC) method to differentiate between cloud ice events according to the dominant vertical ice number fluxes to identify the source footprint of the high ice number 30 concentrations that are measured in alpine regions when the clouds are near or in contact with the surface. Although the EC method is well developed technical constraints make it difficult to apply in complex terrain and inside a mixed phase cloud with relatively low concentrations of large ice particles to determine cloud particle number fluxes. Nevertheless, we present results from several cases using this approach that strongly suggests the upwind surface was the likely source of these particles, at least in the selected cases. In the first part of the chapter, we describe the instrumentation and methodology used at the 35 Jungfraujoch station where it was applied. The next part explains the Eddy Covariance method used for calculating the fluxes of ice particles. We then provide a critical assessment of the results and their current limitations with a view to highlighting the requirements for more detailed experiments.  Table 1) and their position inside the experimental setup. The Jungfraujoch ridge is located at 3580 m a.s.l. in the south of Switzerland between the two mountains Jungfrau and Mönch. 45 The research station was selected because of its accessibility, high cloud frequency (37 % per year) and the extensive characterisation of its local aerosol and cloud microphysics prevalent at the site (Baltensperger et al., 1998;Hinz et al., 2005;Lloyd et al., 2015;Hoyle et al., 2016;Kupiszewski et al., 2016) and aerosol composition, CCN and ice nuclei . On a terrace situated on the station roof (see Figure 2), a tower system used for the measurements has been set up as shown in Figure 1. The use of an aero-wing attached to the top of the tower allows for the automated 3D positioning 50 of the attached instruments to orient them into the mean 3D wind vector on a continuous basis. This arrangement has been already successfully used in previous experiments, e.g. the CLoud Aerosol Characterisation Experiments (CLACE) campaigns (Lloyd et al., 2015). For this field campaign the instruments, described in Table 1, have been used and form an EC system to determine cloud particle number fluxes. The 3-View Cloud Particle Imager (3V-CPI) is a hybrid instrument comprising a twodimensional stereoscopic (2D-S, Lawson et al., 2006) shadow imaging spectrometer and a Cloud Particle Imager (CPI, Stratton 55 Park Engineering Company Inc, Boulder Colorado (SPEC)), a charge-coupled device (CCD) imaging spectrometer (3V-CPI Model 2, SPEC Inc.). The use of a 3D sonic anemometer (Metek USA-1, Meteorologische Messtechnik) fixed to the rotating aero-wing next to the inlet of the 3V-CPI is essential for the flux measurement that connects the fluctuations of vertical wind speed with the change in the number concentration of the measured ice particles over a specific time period. A second Metek sonic anemometer mounted on a separate tower at the same height provided measurements for the alignment of the aero-wing 60 with the mean wind vector. It provides the direction of cloud events. The Cloud Droplet Probe (CDP-100, Droplet Measurement Technologies, DMT, Lance et al., 2010), as well as a Particle Volume Monitor (PVM, Gerber Scientific, Inc.) for bulk liquid water content measurement only (Gerber, 1991), provides information about possible contributions from the liquid phase during the measurements. The CDP, PVM and similar cloud spectrometers have been used previously to determine size dependent surface fluxes of cloud droplet number as well as total liquid water to surfaces in complex terrain, (e.g. Beswick 65 et al. 1991, Vong et al. 1995, Kowalski et al. 1997, Klemm et al. 2005, Pryor et al. 2008. Both the CDP and PVM sensors are also known to be sensitive to ice particles, Gerber & Demott (2014), and as such fluxes from these were not considered here due to uncertainties associated with their lack of ice-liquid discrimination, but rather to quality control total particle number and volume concentrations from the other spectrometers. Vaisala and Rotronics sensors provided measurements of the temperature and the humidity during cloud events. 70 are responsible for creating a strong temperature gradient in the upper part of the snowfield. This creates a warm point just below the surface (around 1 cm). The water vapour originating from this point forms an upward flux that eventually reaches the colder surface and is responsible for a thin region of supersaturation above the surface. This zone of supersaturation is caused by the strong non-linear relationship between the temperature and saturation vapour pressure. In this region, the water 80 vapour freezes and produces part of the crystal-clusters. This way, the sublimation of snow does not lead to a mass loss but to a formation of sublimation crystals. This effect is similar to the production of frost flowers (Style and Worster, 2009) over young sea ice. In contrast to the studies over ice where frost flowers only appear in spots over the ice that provide nucleation sites, the surface of snow provides an endless amount of these nucleation sites. This way, the crystal covered areas can extend over very large areas. The whole process of sublimation crystal production is very sensitive to the density of the snowpack as 85 it alters the thermal profile and hence the release of water vapour. In addition to the mentioned processes, the latent heat shows a diurnal cycle that leads to sublimation of snow during the day and condensation of water vapour from the atmosphere during the night. This condensation process generates surface hoar crystals and further contributes to the formation of crystal clusters.
Both processes produce similar crystals that can persist over several weeks before being removed by strong winds.

Eddy covariance method
The eddy covariance (EC) technique provides a direct method of measuring the vertical exchange between the surface and its surrounding (Lee et al., 2005). It has been used to measure the deposition of cloud droplets to hill sites, moorland and forests in sizes up to 15.5 µm (Gallagher et al., 1988(Gallagher et al., , 1992Beswick et al., 1991;Kowalski et al., 1997;Klemm et al., 2005;Pryor et al., 2008), of fog (Holwerda et al., 2006), and aerosol particles (Schmidt and Klemm, 2008). The original theories 95 have been established under the conditions of stationarity in time and homogeneity in space. These conditions are hardly met in practice and therefore, data quality checking, instrument corrections and coordinate rotation of the wind field must be evaluated to obtain meaningful fluxes. The basis of the EC approach is to define a control volume over a representative, ideally homogeneous, surface footprint, recording the accumulation within it, measure the exchange across all the aerial faces of the volume and infer the surface exchange by the difference. If the turbulent mixing acts as a physical averaging operator, the flow 100 field is effectively one-dimensional, and the flux can be written as where c(t) is a generic scalar, w(t) is the vertical component of the velocity vector, the overbar denotes the averaging operation, ̅ is the surface flux, z the vertical or surface normal coordinate, and ( ) is the Dirac delta function. If on the one hand there is no accumulation of c over time, the first term on the left-hand side of formula 1 becomes zero (condition of stationarity). 105 Integrating the formula from the ground z=0, to the sensor height h, taking into account any zero-plane displacement height for momentum transfer, the function reduces to the total covariance of w(t) and c(t):

( ) ######## = 5
(2) In the end, the total covariance needs to be replaced with the measured eddy flux. The Reynolds decomposition is applied to separate the expectation value of a quantity from its fluctuations with the decomposition written as: 110 The basic formula for the flux calculations is therefore: The first term on the right-hand side of equation (4) often gets substituted by 0. The argument in the literature is that 5 → 0 for the conditions of stationarity in time and homogeneity in space. As the measurement location in this experiment is atop a ridge with aspect ratios that differ depending on wind direction, equation 4 is the final version and the remaining term is used to test if the corrections made were valid or not. In this experiment the calculations presented are performed on particle optical size ranges (> 50 µm) discriminated as far as possible by particle shape as a proxy for ice versus water phase. As such there 120 may be significant uncertainties associated with counting statistics based on sample volume limitations and corrections to particle concentrations which need to be carefully addressed.
The next question that has to be addressed is the averaging period. The eddy fluxes need to be calculated over a sufficient period such that any significant contributing eddies can be sampled. On the other hand, the flux can be influenced by meteorological changes rather than the fast turbulent changes, if the averaging time is too long. To address this problem, the 125 ogive presentation of the co-spectrum of the vertical wind component and another property (e.g., sonic temperature) can be used. An ogive shows the cumulative contribution of different eddy scales of increasing importance to the total transport as a point on the ogive plot which is the integral under the spectral density curve between the highest frequency recorded and the frequency of interest. If the curve of the ogive reaches its asymptote, the corresponding time to this frequency represents an adequate flux averaging time. Another important task to address is definition of an appropriate wind coordinate system. 130 Recently, the planar fit coordinate system has been favoured in eddy flux measurements (Wilczak et al., 2001). Unfortunately, this coordinate system cannot be used if the sonic anemometer has been moved frequently as is required to minimise possible instrument inlet sampling issues due to misalignment with the mean flow as in this experiment where the slope flow has a significant non-vertical component. In this experiment the rotating wing, used on the Jungfraujoch to minimise artefacts in particle size distribution measurements arising from instrument inlet alignment, moves every minute into the mean of the wind 135 vector. Therefore, we are able to adopt the natural wind coordinate system forcing the mean lateral and cross wind components to zero without the need for a third rotation as has been used in low wind environments where the vertical component may be significant compared to the horizontal wind components (see Finnigan, 2004). In addition, each periods of coordinate rotation have to be carefully evaluated for (unnatural) over-rotations. A separate fixed 3D sonic anemometer on a parallel tower at the same height is therefore used for quality control of the rotations. 140 Table 2 Summary of days with significant upward (+) and downward (-) fluxes for particles with a diameter of more than 50 µm.

Date
Time ( time of typically 8 minutes for flux calculations. The wind direction has been corrected for two minutes only, as the 3D wing rotates into the mean wind field every minute to ensure correct, unperturbed measurements of the particle number concentration and size distributions. For the concentration, only particles larger 50 µm are considered in the calculations of the flux due to uncertainty in particle shape analysis. Shape and sample volume corrections to the optical array probes (OAP) instrument data 150 products (2D-S) are discussed in detail by Crosier et al. (2011) and Korolev (2007) whilst the CCD imaging probe used to quality control the shape analysis interpretation is discussed by Connolly et al. . A more recent discussion of the limitations and comparisons between various imaging probes and optical scattering probes, e.g. CDP, are discussed in further detail by O'Shea et al. (2019). However it should be emphasised that many of these discussions are for aircraft based measurements which can be problematic due to much higher sampling speeds leading to increased artefacts such as particle 155 shattering on instrument surfaces e.g. Korolev et al. (Korolev et al., 2013 b) which has led to such instruments being modified to minimise these effects, e.g. Lawson et al. (2011) .

(circles) represent data that extend beyond the whiskers.
To focus the discussion, here we concentrate on periods where significant upward particle fluxes, over 800, or downward fluxes, of less than -800 m -2 s -1 are observed (20 % of all periods with fluxes larger or smaller than 0 m -2 s -1 ). The box plots for all such periods above this threshold (see Figure 3 from the south. Downward cases tend to have higher absolute values for the calculated flux > -10,000 m -2 s -1 than upward cases (+6,000 m -2 s -1 ). Outlying cases with fluxes up to +-30,000 m -2 s -1 are often reached during cloud events with an overall enhancement of ice particle concentrations but not out of cloud. The larger liquid water content (LWC, based on CDP data) in 170 downward cases in comparison to upward cases is also seen in the analysis of the CPI images where prevailing water droplets that are larger than 50 µm can be seen. These large droplets do not occur in cases with winds from the south and hence for periods with temperatures lower than in northern cases. This is also reflected in the mean wind speeds that are typically greater for southern cases than northern cases.
In conclusion, the differences between northern and southern flux cases may be related to the orography of the Jungfraujoch 175 ridge which has a steep and narrow flank in the north down to the base of the mountain with just a few snow covered fetches by comparison with a gently sloped southern fetch dominated mainly by the Aletsch glacier. Only the last 100 m presented to the station from the southern sector is extremely steep . This difference might lead to the general increase in speed up factor from southern slope. In addition air masses from this fetch experience a longer cooling period than the more rapid, short vertical ascent experience on the northern side. 180

Figure 4 Scatter plot between concentration and flux vs. wind speed (A-F), concentration vs. LWC (G-H), and concentration vs. temperature (I-J) for particles larger 50 µm (left panels) and smaller 50 µm (right panels). Each point represents an 8-minute mean (blue: from the south, red: from the north).
Looking into the summary above in more detail, out of the 28 days of the campaign, 14 days show relevant flux periods 185 (more/less than +-800 m -2 s -1 , see Table 2). Scatter plots of the data of these 14 days show more interesting periods (see Figure   4).  Figure 4 show the mean concentration of the 8-minute periods used to calculate the flux from. They reveal that there are only 7 out of 213 periods from the north with mean wind speeds larger than 8 ms -1 , with maximum concentrations around 10 3 l -1 for particles larger 50 µm. In contradiction, flux periods from the south cover the whole range of wind speeds 190 up to 20 ms -1 . The same observation can be done for particles smaller 50 µm.

Panel A and B of
Comparing flux calculations versus wind speed, it appears that there may be a trend with larger upward fluxes with increasing wind speed up to 8 ms -1 for both, northern and southern wind directions (Figure 4, C). A similar relationship does not exist for upward fluxes and downward fluxes (see Table 3 while fluxes of small particles reach higher values by two orders of magnitude. Looking at the LWC and concentration (Figure 4, G), a linear correlation seems to exist up to 0.2 gm -3 for fluxes coming from 200 the south. Almost no fluxes from the south exist for larger LWC (3 out of 196). Two out of 3 are negative fluxes from the 2 nd and 9 th of February 2017. The size distribution of the downward flux cases looks similar and differ from the positive flux case (see Figure 5). The negative fluxes have a smaller concentration for particles larger 100 µm. For smaller particles, it seems like there is a peak for the positive flux period in comparison to the negative flux periods. There is also no cut off for larger particles after this maximum and the maximum size for the CDP. Instead, there are particles up to 10 3 l -1 µm -1 . This can also 205 be seen in the CPI pictures of the 3V-CPI (see Figure 6). Very large particles, such as intact dendrites or fragments of them, plates and particles with no specific structures around 50 µm coexist with very small particles, presumably water droplets. As we used 50 µm to distinguish water droplets from ice particles to calculate the fluxes from the larger than 50 µm ice particles. This is a limitation of the approach used here. Due to technical limitations of the imaging probes and the uncertainties associated with the scattering spectrometers it becomes technically challenging to discriminate water droplets from non-210 spherical ice particles for particles smaller 50 µm. This is especially true for the 2D-S part of the 3V-CPI. The resolution of the 2D-S, whose data was used to calculate the flux from, is significantly smaller than that of the CPI (10 vs 2.3 µm), which is a general compromise associated with all cloud instruments between sample volume and field of view and hence counting statistics, e.g. the maximum value of the sample volume of the CPI is only 4 % of the 2D-S sample volume (Lawson et al. 2006). These are the reasons why the data of the 2D-S have been used. Comparing the particle images of the 10 th of February 215 with the 2 nd and 9 th of February (with negative fluxes) shows also larger ice particles but no large dendrite like ones and also plates in case of the 2 nd (even double plates) and rimed particles on the 9 th of February. Liquid water droplets on the 9 th exceeded the 50 µm threshold. This way, they may appear as ice particles using the 50 µm threshold. This could be the reason for the negative ice particle flux. In reality the situation could be reversed if the liquid water droplets didn't appear in the number concentration of the ice particles. The 10 th of February shows rimed ice particles measured at lower temperatures than the 220 negative flux periods on the 2 nd and 9 th of February. The size distribution (see Figure 5) shows larger concentrations, around one to two orders of magnitude larger at the second peak for larger particles and also a wider spread size distribution. This distribution is not as defined as the one for the negative fluxes. Here, water droplets also existed but smaller than the cut-off threshold of 50 µm and therefore do not affect the flux measurement for the larger ice particles. Examining the LWC versus the concentration panels (see Figure 4, G-H), there is enough evidence to conclude there is a significant linear relationship between the LWC and the measured concentration from the South (at least for some cases) as 230 well as for the North (see Table 4). This seems to correspond to the riming effect that not have been ruled out by Rogers & Vali (1987). More liquid water droplets interact with the surface (snow, stone, and trees) and produce small ice crystals that get lofted into the cloud. This correlation can be seen even more for small particles. The concentration versus temperature panels (see Figure 4,     • Fluxes can be obtained from measured ice particle number concentrations. In cloud free cases, the most reasonable source is ice particles grown on the surface particularly when the wind speed threshold does not reach values necessary for initiation of blowing snow. CPI images recorded during suspected blowing snow events confirm that 255 the particles are similar to those observed previously during blowing snow events, namely aged particles with rounded edges due to the transport through unsaturated air with respect to ice. Still, the low wind speeds as well as a missing correlation between wind speed and ice number concentration allow us to discard blowing snow as a mechanism for the measured high ice number concentrations. We, therefore conclude that frost crystals growing on the snow surface are thye likely source of the ice crystals following Lloyd et al 260 • Also, fluxes (with relevant values) do not always occur. It seems that it takes time for the surface source to form and reach a number density sufficient for the wind to begin generating a measurable upward flux. These upward fluxes also have to interact with an orographic cloud near enough to the surface to influence the ice number concentrations inside the cloud. The observed fluxes and conditions in which thy occur are not generally consistent with blowing snow 265 • In cloud cases indicate strong positive and negative fluxes during an event are present most of the time. CPI images show similar ice particle habits present during each event. This leads to the conclusion that the Jungfraujoch research station remote from the source of the ice particles fluxes so that the air mass sampled at the summit has become well-mixed.
• Also, the high ice number concentration during cloud events masks the fluxes. On the one hand, downward fluxes 270 of ice particles subdue local upward fluxes from the surface. On the other hand, events from the south produce high liquid water content with larger droplets that have downward motion. Because they are larger than 50 µm they are included into the calculation of the large particle fluxes. Hence, they influence the upward fluxes.
• There are cases, with a mixture of up-and downward fluxes but also cases that only present positive or purely negative fluxes. The reason for these cases with one dominant direction from or to the surface has yet to be found. 275 • The use of a 3D rotating wing is suited for the measurement of ice particles using a closed path shadow imager that normally is located underneath a wing of an airplane that has to be faced parallel to the mean wind field. For future work, we suggest improvements to such setups can include use of open-path cloud instruments as closed path instruments are more prone to shattering, although such artefacts can be removed using inter-particle arrival time analyses. Other potential particle sources, e.g., from above the Jungfraujoch, could be examined using 280 contemporaneous ceilometer measurements along with local cloud base location. Particle concentrations upwind of the station would also confirm the results presented here.
Calculating the fluxes for particles larger 50 µm with reasonable high values similar to the values that Farrington et al. (2016) modelled (around 1e3 m -2 s -1 ) strengthens and reaffirm the conclusion from previous publications for the existence of the surface fluxes mechanism. Hereby, it is still dubious if the origin of the upward flux is the growth of surface hoar and 285 sublimation crystals. The riming effect might be an additional possible mechanism as the correlation between LWC and concentration of large but also for small particles suggests as it also happens in the close proximity of the surface. More experiments closer to the origin or in the lab are therefore necessary to reveal more of the mechanism in detail.

Author contribution 290
WS, GL, KB, MF made measurements at the mountain top site. WS analysed the data from the campaign and prepared the manuscript. All authors were involved in interpretation of the findings from the campaign.