Retrieving ice-nucleating particle concentration and ice multiplication factors using active remote sensing validated by in situ observations

Understanding the evolution of the ice phase within mixed-phase clouds (MPCs) is necessary to reduce uncertainties related to the cloud radiative feedback in climate projections and precipitation initiation. Both primary ice formation via ice-nucleating particles (INPs) and secondary ice production (SIP) within MPCs are unconstrained, not least because of the lack of atmospheric observations. In the past decades, advanced remote sensing methods have emerged which provide high-resolution data of aerosol and cloud properties and could be key in understanding microphysical processes on a global scale. In this study, we retrieved INP concentrations and ice multiplication factors (IMFs) in wintertime orographic clouds using active remote sensing and in situ observations obtained during the RACLETS campaign in the Swiss Alps. INP concentrations in air masses dominated by Saharan dust and continental aerosol were retrieved from a polarization Raman lidar and validated with aerosol and INP in situ observations on a mountaintop. A calibration factor of 0.0204 for the global INP parameterization by DeMott et al. (2010) is derived by comparing in situ aerosol and INP measurements, improving the INP concentration retrieval for continental aerosols. Based on combined lidar and radar measurements, the ice crystal number concentration and ice water content were retrieved and validated with balloon-borne in situ observations, which agreed with the balloon-borne in situ observations within an order of magnitude. For seven cloud cases the ice multiplication factors (IMFs), defined as the quotient of the ice crystal number concentration to the INP concentration, were calculated. The median IMF was around 80, and SIP was active (defined as IMFs > 1) nearly 85 % of the time. SIP was found to be active at all observed temperatures (−30 to −5 ^∘C), with the highest IMFs between −20 and −5 ^∘C. The introduced methodology could be extended to larger datasets to better understand the impact of SIP not only over the Alps but also at other locations and for other cloud types.

1 Introduction and background

The increase in the Earth's global mean temperature in recent years is unequivocal, yet the extent of a cloud cooling effect remains most uncertain (IPCC, 2021). The radiative feedback of a cloud is a strong function of the hydrometeor phase (Sun and Shine, 1994). Mixed-phase clouds (MPCs) consisting of both water phases (liquid and solid) contribute strongly to uncertainties in the radiative feedback. They are thermodynamically unstable (see e.g., Korolev et al., 2017) because of the lower vapor pressure with respect to ice than with respect to liquid water. This causes the ice crystals to grow at the expense of the evaporation of cloud droplets, which is referred to as the Wegener–Bergeron–Findeisen process (Wegener, 1911; Bergeron, 1935; Findeisen, 1938). The complexity of phase partitioning adds to difficulties simulating MPCs in models (e.g., McCoy et al., 2016; Matus and L’Ecuyer, 2017) where experimental observations could reduce uncertainties (see e.g., Baumgardner et al., 2011; Mahrt et al., 2019). In the evolution of an MPC, the ice phase plays an important role as it controls precipitation initiation and consequently cloud lifetime (see e.g., Field and Heymsfield, 2015; Heymsfield et al., 2020). Thus, the ice phase determines not only how a cloud impacts the radiative budget but also for how long. The first ice crystals in a cloud either formed within the cloud by homogeneous or heterogeneous ice nucleation or can be externally introduced by, for example, sedimentation from a higher cloud (seeder–feeder process, e.g., Proske et al., 2021; Ramelli et al., 2021a, and references therein) or levitated from the ground (blowing snow, especially in mountainous terrain, e.g., Beck et al., 2018; Walter et al., 2020, and references therein). Approximately below −38 ^∘C, supercooled droplets can freeze homogeneously or heterogeneously. Above that temperature, heterogeneous ice nucleation is favored on sparsely abundant aerosols called ice-nucleating particles (INPs; see e.g.,Wegener, 1911; Vali, 1971; Pruppacher and Klett, 2010; Murray et al., 2012). Until today, a variety of aerosol particles acting as INPs in the atmosphere have been identified such as desert and soil dust, organics from biomass burning, marine or terrestrial biogenic particles, atmospherically aged soot, bacteria, and others (Kanji et al., 2017; Huang et al., 2021). However, due to the geospatial variability of INP sources, atmospheric INP concentrations feature a high spatiotemporal variability, complicating their quantitative assessment (see e.g., DeMott et al., 2010; Kanji et al., 2017; Murray et al., 2021, and references therein). Understanding heterogeneous ice formation and subsequent ice crystal growth is key to understanding the link between aerosols and precipitation formation and is therefore an important step towards constraining weather and climate predictions (Ansmann et al., 2019a; Bühl et al., 2019). After first ice crystals are found in a cloud, secondary ice production (SIP) can enhance the ice crystal number concentration (ICNC) by various processes, e.g., fragmentation during ice–ice collisions (e.g., Vardiman, 1978; Takahashi et al., 1995) or splintering during riming (Hallett–Mossop process, active between −8 and −3 ^∘C, Hallett and Mossop, 1974), given favorable environmental conditions (Korolev and Leisner, 2020). The Hallett–Mossop process requires the presence of supercooled liquid droplets. Thus, the Hallett–Mossop process can play an important role, especially over orographic terrain, given that the orographic-forcing-induced updrafts could sustain the availability of liquid water droplets (Lohmann et al., 2016). Various field studies have observed ICNC exceeding the ambient INP concentration by several orders of magnitude (see e.g., Koenig, 1963; Auer et al., 1969; Hobbs and Rangno, 1985, 1990, 1998; Gayet et al., 2009; Crosier et al., 2011; Stith et al., 2011; Crawford et al., 2012; Heymsfield and Willis, 2014; Lawson et al., 2015; Lasher-Trapp et al., 2016; Ladino et al., 2017; Mignani et al., 2019; Lauber et al., 2021; Pasquier et al., 2022). However, the underlying SIP processes are weakly constrained (Korolev and Leisner, 2020). Therefore, the prevalence of SIP in the atmosphere as well as the environmental conditions for SIP to be active remains uncertain.

Uncertainties in predicting atmospheric INP concentrations and ice multiplication are (partly) related to spatiotemporally limited field observations, because of the large effort needed to obtain field data and their point-like characteristic. Remote sensing techniques can be a suitable solution to overcome the issue by providing continuous data in time and at least one spatial dimension. Starting in the 1930s, remote sensing techniques have evolved to very advanced and reliable instruments today (Wandinger, 2005). Using lidar (light detection and ranging) and radar (radio detection and ranging) instruments, a broad variety of aerosol and cloud properties can be retrieved (e.g., Ansmann et al., 2019a). Lidar systems are key instruments for the investigation of aerosol optical properties and their vertical layering in the atmosphere. The retrieved backscatter and extinction coefficients allow for the determination of physical properties of aerosol particles in the atmosphere, such as size and particle number concentration (see e.g., Müller et al., 1999; Veselovskii et al., 2005; Mamouri and Ansmann, 2016). Advanced lidar systems also deliver other products such as extinction-to-backscatter ratio (lidar ratio), linear depolarization ratio (LDR), or Ångström exponent from multiwavelength observations of the backscatter and extinction coefficients that allow the determination of more detailed aerosol information, e.g., aerosol type (see e.g., Burton et al., 2012; Groß et al., 2013; Baars et al., 2017). The accuracy of lidar retrievals has been validated in many studies, with a focus on different properties comparing lidar measurements with in situ observations using aircraft or unmanned aerial vehicles (UAVs) (Ferrare et al., 1998; Wandinger et al., 2002; Sakai et al., 2003; Müller et al., 2012; Sawamura et al., 2017; Schrod et al., 2017; Marinou et al., 2019; Haarig et al., 2019; Düsing et al., 2021). Furthermore, the LDR was utilized to determine the contribution of dust to the observed aerosol (Shimizu et al., 2004; Tesche et al., 2009; Mamouri and Ansmann, 2016; Haarig et al., 2017). Precise knowledge about the aerosol type is essential to use the lidar retrieval to estimate INP concentrations, owed to the fact that the INP concentration is deduced only from a physical property such as particle number concentration and surface area concentration. However, the parameterizations are often proposed for specific aerosol types such as dust (e.g., Niemand et al., 2012; DeMott et al., 2015; Ullrich et al., 2017; Harrison et al., 2019), marine aerosol (e.g., McCluskey et al., 2018), soot (e.g., Ullrich et al., 2017), and global aerosol (e.g., DeMott et al., 2010). The feasibility of predicting INP concentration from lidar profiles has been shown by Mamouri and Ansmann (2016). In lidar INP studies, typically three aerosol type categories are used: mineral dust (Ansmann et al., 2003; Mamouri and Ansmann, 2016; Schrod et al., 2017; Haarig et al., 2019; Marinou et al., 2019), marine aerosol (Mamouri and Ansmann, 2016; Haarig et al., 2019), and continental aerosol (Mamouri and Ansmann, 2016; Schrod et al., 2017; Düsing et al., 2018; Marinou et al., 2019; Düsing et al., 2021). If necessary, these categories can be extended by wildfire smoke (Ansmann et al., 2021; Engelmann et al., 2021) or volcanic ash (Ansmann et al., 2011). Whereas the aerosol sources within the first two categories can be comparably narrowed down by source, continental aerosol sources are much more diverse, featuring (among others) biological, organic, or lifted soil particles and also anthropogenic emitted particles from biomass burning or combustion for example (Kanji et al., 2017; Huang et al., 2021, and references therein), complicating the retrieval of continental INP concentrations. Thus, finding the most ideal INP parameterization for a given aerosol type is key for minimizing uncertainties in predicted INP concentration. Schrod et al. (2017) proposed calibration factors to optimize the INP parameterization by DeMott et al. (2015) for the retrieval of INP concentration in condensation and deposition mode from dust-dominated air masses. Marinou et al. (2019) showed that separating the dust-carrying air masses into a dust and a continental component and applying the parameterizations of DeMott et al. (2015, 2010), respectively, yields the best agreement between lidar retrieval and in situ observations (basing on Schrod et al., 2017). However, Marinou et al. (2019) stated further that “additional measurements are required in order to define the optimum INP parameterizations for non-dust atmospheric conditions (e.g., continental, marine, smoke)”, suggesting that an optimal INP parameterization for the retrieval of INP concentration from continental air masses is yet to be proposed.

Remote sensing instrumentation has recently increasingly been used for cloud observations and is becoming a promising method for retrieving ice-crystal-related parameters (see e.g., Seifert et al., 2010; Bühl et al., 2013; Zhang et al., 2014; Ansmann et al., 2019a). It allows for a continuous monitoring of the vertical cloud structure and therefore permits the observation of spatiotemporal cloud evolution (Seifert et al., 2010). This is a crucial addition to in situ measurements, which only allow for measurements in a constrained time and height frame (Seifert et al., 2010). The ice water content (IWC) can be retrieved solely based on radar measurements (Hogan et al., 2006), and ICNC can be estimated by linking the IWC to an assumed particle shape and size distribution (see e.g., Bühl et al., 2019). Complementing radar retrievals with size information from a collocated lidar improves the accuracy of cloud property retrieval further (Delanoë et al., 2013). In previous studies, remote sensing was used to study (primary) ice nucleation in the atmosphere (see e.g., Sakai et al., 2003; Ansmann et al., 2008; Eidhammer et al., 2010; Seifert et al., 2010, 2011; Ansmann et al., 2019a; Engelmann et al., 2021) and to estimate SIP (see e.g., Auer et al., 1969; Luke et al., 2021; Sotiropoulou et al., 2021). The identification of a specific SIP process from atmospheric in situ and remote sensing observations is a challenge. Lauber et al. (2021) recently showed that recirculation of melted ice crystals can enhance the concentration of small ice crystals in the proximity above the melting layer. Whereas deepening the understanding of individual processes seems best achievable in controlled conditions during laboratory studies (Korolev and Leisner, 2020), the contribution and magnitude of SIP still needs to be assessed in the (complex interacting) atmosphere. From field observations, the ice multiplication factor (IMF) defined as the ICNC divided by the INP concentration can be utilized to quantify SIP. Whereas the IMF allows quantification of the excess of ice crystals, it cannot be used to directly identify a specific SIP process due to the potential (temporal and spatial) difference between the origin of the first ice crystals and SIP to occur.

In this study we combined a suite of in situ and remote sensing instruments to understand cloud formation and evolution in orographic terrain. In February and March 2019, we performed an intensive field campaign in the Swiss Alps in the region of Davos. In a high valley, a combined lidar–radar system was employed along with ground-based in situ aerosol (including INP) observation (Wieder et al., 2022b) and balloon-borne in situ cloud observations (Ramelli et al., 2021a, b). A second in situ aerosol site was located on a nearby mountaintop (height difference 1.1 km, Mignani et al., 2021; Wieder et al., 2022b). The near collocation of the lidar beam and the mountaintop site (horizontal displacement 3.65 km) is ideally suited for aerosol and thus INP closure as also previously shown by Bedoya-Velásquez et al. (2018), who studied the hygroscopicity of aerosol particles. Based on aerosol data collected over 8 weeks and seven observed cloud events, we address the following. First, we validate the lidar retrieval of aerosol number concentration and surface area concentration, which is the basis for the INP concentration retrieval (Sect. 3.1). Second, based on INP concentrations measured during a Saharan dust event, we assess the accuracy of different dust parameterizations and evaluate their performance (Sect. 3.2.1). Third, we evaluate two INP parameterizations (DeMott et al., 2010; Ullrich et al., 2017) used in the lidar community for the retrieval of INP concentration from continental air masses (Sect. 3.2.2). Based on our in situ measurements, we propose a calibration factor to optimize one INP parameterization and validate the tuning with the lidar observations (Sect. 3.2.3). Fourth, we validate radar-retrieved IWC and ICNC with balloon-borne in situ observations (Sect. 3.3). Lastly, using the tuned parameterization we present a methodology to obtain estimates of IMFs in MPCs (Sect. 3.4).

Figure 1 Overview of the measurement locations: (a) Davos-Wolfgang (black star) in the east of Switzerland (map source: Federal Office of Topography). (b) Measurement locations with the local topography around. The OCEANET container and the lower aerosol monitoring site were located in a high valley at Wolfgangpass (blue dot, WOP, 1631 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). A second aerosol measurement site was located on the mountaintop Weissfluhjoch (orange dot, WFJ, 2693 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$), which is located approximately 1.1 km higher than WOP and horizontally displaced by 3.65 km. Lidar product retrieval heights for the comparison with in situ measurements at WFJ and WOP are indicated in orange and blue, respectively, right of the lidar beam (green). The topography was extracted from the digital height model DHM200 from the Federal Office of Topography swisstopo. (c) View from WFJ in direction of WOP. The approximate location of the lidar beam is indicated in green. (d) Placement of the aerosol trailer right of the OCEANET container's lidar and radar at WOP. The tethered balloon used for in situ cloud observations is seen in the front right.

2 Measurement setup and methodology

The RACLETS (Role of Aerosols and CLouds Enhanced by Topography on Snow) campaign took place in the region of Davos, Switzerland (Fig. 1a), in February and March 2019, where an extensive set of aerosol, cloud, precipitation, and snow measurements were conducted (Envidat, 2019; Walter et al., 2020; Mignani et al., 2021; Ramelli et al., 2021a, b; Lauber et al., 2021; Georgakaki et al., 2021; Wieder et al., 2022b). The two main measurement sites were located on a saddle at the entrance of a high valley (Wolfgangpass, 1631 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$, hereafter referred to as WOP) and the other on a mountaintop (Weissfluhjoch, 2693 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$, hereafter referred to as WFJ) as shown in Fig. 1b. At both locations, two similarly equipped aerosol measurement sites were set up. In addition, the remote sensing instrumentation was installed at WOP.

2.1 In situ aerosol and INP measurements

The aerosol measurement sites at WFJ and WOP were previously presented (Wieder et al., 2022b; Georgakaki et al., 2021; Mignani et al., 2021). At both sites, ambient air was sampled through a 46 ^∘C heated inlet. The heating was a preventive measure to avoid icing of the outside inlet parts, to evaporate activated cloud droplets, and to sublimate ice crystals. The evaporation of volatile compounds of the aerosol cannot fully be excluded. However, the effect is expected to be minor given the high flow rate through the inlet (300 L min⁻¹), such that the temperature of sampled air was likely below 46 ^∘C. Furthermore, the degradation of relevant INPs (mostly biological) should only occur at temperatures above 46 ^∘C (Kanji et al., 2017; Huang et al., 2021) and is hence regarded as unlikely (see Wieder et al., 2022b, for further details). Downstream, an aerodynamic particle sizer spectrometer (APS; model 3321, TSI Inc., US) and a scanning mobility particle sizer spectrometer (SMPS; model 3938, TSI Inc., US) recorded aerosol size distributions between approximately 10 nm (electric mobility) and 20 µm (aerodynamic diameter). In this study, electric mobility (SMPS) and aerodynamic diameter (APS) were converted to physical diameter assuming a shape factor χ=1.2 and assuming a particle density ρ=2 g cm⁻³ (Thomas and Charvet, 2017). After conversion, the observed size range of the SMPS and the APS covered particles with physical diameters between 10 and 400 nm and between 400 and 15 µm, respectively. The total surface area concentration (hereafter referred to as s) of all aerosols was calculated from the entire size distribution utilizing both SMPS and APS data. The number concentrations of particles with radii ≥ 250 nm (hereafter referred to as n₂₅₀) were obtained by integrating the APS distributions for physical diameters ≥ 500 nm. The error in the obtained concentrations is given by ±10 %. As described above, we cannot fully exclude the loss of volatile compounds which could shift the size distributions to smaller sizes, thus underestimating in situ predicted INP concentrations when using size-dependent aerosol number concentrations as an input parameter (Sect. 3.2.3). However, as we discussed above the effect is likely negligible.

Ambient aerosol was collected over a time span of 20 min into pure water (15 mL, W4502-1L, Sigma-Aldrich, US) for offline INP analysis using high-flow-rate impingers (Coriolis^® μ, Bertin Instruments, France, 300 L min⁻¹) attached to the end of the inlets (Wieder et al., 2022b; Mignani et al., 2021). The (immersion mode) INP analysis was done on site utilizing the drop-freezing apparatuses LINDA (Stopelli et al., 2014) at WFJ and DRINCZ (David et al., 2019) at WOP. In short, the liquid sample is split into small aliquots of the same volume (V_a) and cooled down in a cryostat. During cooling a camera mounted above the cryostat takes pictures of the droplets in which frozen droplets appear more opaque than unfrozen ones. After all droplets froze, the fraction of frozen droplets as a function of the cryostat temperature (FF(T)) is retrieved in post-processing from the pictures. Consequently, the INP concentration is calculated according to Vali (1971) as

$$\begin{array}{}\text{(1)}& n_{\mathrm{INP}}\left(T\right)=-{\displaystyle \frac{\ln [\mathrm{1}-\mathrm{FF}(T\left)\right]}{V_{\text{a}}\times C}},\end{array}$$

with C as the conversion factor from INP concentration found in the sample liquid (in mL⁻¹) to the concentration in ambient air (in StdL⁻¹). The error in the obtained INP concentration is temperature dependent and varies per sample, but extends on average to ±50 %. Further details about the aerosol setups as well as the INP data processing can be found in Wieder et al. (2022b).

2.2 Lidar and radar measurements

During RACLETS, a 35 GHz (Ka-band) cloud radar and the multi-wavelength polarization Raman lidar Polly^XT-OCEANET (Engelmann et al., 2016, hereafter referred to as Polly^XT) of the Leibniz Institute for Tropospheric Research (TROPOS) were deployed at the WOP site (Fig. 1d) as part of the mobile observation platform OCEANET-Atmosphere (Griesche et al., 2020). The instruments provided a continuous overview on the temporal evolution of the vertical distribution of aerosol particles and clouds from 8 February to 16 March 2019 and are used in the following for remote-sensing-based retrievals of aerosol and cloud microphysical properties.

Polly^XT measured vertical profiles of the particle backscatter coefficient (at 355, 532, and 1064 nm), the particle extinction coefficient with the Raman method (nighttime, 355 and 532 nm), and the particle linear depolarization ratio (LDR, 355 and 532 nm). Based on the Polly^XT observations, INP number concentration was calculated following the POLIPHON method as described in Mamouri and Ansmann (2016). Here, only a short description of POLIPHON will be provided. In a first step, the dust and non-dust contributions to the backscatter coefficient are separated using the particle LDR at 532 nm. The backscatter contributions are transferred to extinction contributions using an extinction-to-backscatter ratio (lidar ratio) of 55 sr for dust and 50 sr for continental particles (rural background, pollution). These extinction coefficients σ_i at height z are converted for each aerosol type i (dust or continental aerosol) to n_250,i and s_i following the relation described in Mamouri and Ansmann (2016).

$$\begin{array}{}\text{(2)}& {\displaystyle }{n}_{\mathrm{250},i}\left(z\right)& {\displaystyle }=c_{\mathrm{250},i}\times {\mathit{\sigma }}_i\left(z\right)\times C\left(z\right)\text{(3)}& {\displaystyle }{s}_i\left(z\right)& {\displaystyle }=c_{\mathrm{s},i}\times {\mathit{\sigma }}_i\left(z\right)\times C\left(z\right)\end{array}$$

The additional factor C(z) defined as $p_{\mathrm{0}}/p\left(z\right)\times T\left(z\right)/T_{\mathrm{0}}$ converts the retrieved concentrations to standard conditions (Stdcm⁻³ for n₂₅₀ and µm² Stdcm⁻³ for s) at standard pressure p₀ and standard temperature T₀ according to the ambient pressure p(z) and temperature T(z) at height z. The meteorological data were taken from the COSMO-1 reanalysis data (see Sect. 2.4). The needed conversion factors ($c_{\mathrm{250},i},c_{\mathrm{s},i}$) are aerosol type dependent and are obtained from long-term sun-photometer (Aerosol Robotic Network, AERONET,  Holben et al., 1998) datasets from the Sahara (dust, Ansmann et al., 2019b) and the Alpine site of Davos (continental aerosol, calculated for the present study) following the method described in Mamouri and Ansmann (2016). The applied extinction-to-number-concentration (particle radii ≥ 250 nm) conversion factors are $c_{\mathrm{250},\mathrm{dust}}=\mathrm{0.19}$ Mm cm⁻³ for Saharan dust and $c_{\mathrm{250},\mathrm{cont}.}=\mathrm{0.0828}$ Mm cm⁻³ for continental aerosol (Alps), and the extinction-to-surface-area conversion factors are $c_{\mathrm{s},\mathrm{dust}}=\mathrm{2.4}\times {\mathrm{10}}^{-\mathrm{12}}$ Mm m² cm⁻³ for Saharan dust and $c_{\mathrm{s},\mathrm{cont}.}=\mathrm{2.48}\times {\mathrm{10}}^{-\mathrm{12}}$ Mm m² cm⁻³ for continental aerosol (Alps). The uncertainties are 30 % for n₂₅₀ and 30 %–50 % for s. Further assessments of the uncertainties are provided in Mamouri and Ansmann (2016) and Haarig et al. (2019). With this method, the vertical profiles of the basic input parameters in various INP parametrizations are derived (see Sect. 3.2). The uncertainty in the obtained INP concentrations is about a factor of 3 (Haarig et al., 2019).

The 35 GHz cloud radar was of the type Mira-35 (Görsdorf et al., 2015). During RACLETS, the radar was operated in vertical-stare mode. Pulses with a length of 208 ns were emitted at a repetition frequency of 6000 Hz, resulting in a vertical resolution of 31.17 m and a maximum unambiguous velocity range of 25.6 m s⁻¹, which spans from −12.8 to 12.8 m s⁻¹. The return signals of the emitted linearly polarized pulses were detected separately in the co- and cross-polarized planes. For both channels, Doppler spectra are derived from Fourier transformations of the return signals from a series of 512 consecutive pulses, corresponding to a Doppler-velocity resolution of 0.05 m s⁻¹. The final temporal resolution of the acquired cloud radar dataset of 10 s is obtained from incoherent averaging of 100 consecutive Doppler spectra.

From the cloud radar's reflectivity the IWC was derived according to Hogan et al. (2006) as

$$\begin{array}{}\text{(4)}& \mathrm{IWC}={\mathrm{10}}^{\mathrm{0.000242}ZT+\mathrm{0.0699}Z-\mathrm{0.0186}T-\mathrm{1.63}}\end{array}$$

with the radar reflectivity (Z) and the ambient temperature (T; see Sect. 2.4). For this approach, the expected error in the retrieval is temperature dependent up to a factor of 2. Furthermore, ICNCs were retrieved from the cloud radar observations with the method described in Bühl et al. (2019). The procedure to derive ICNC from the RACLETS cloud radar observations was described by Ramelli et al. (2021a) and is only explained briefly in the following. Observations of LDR, radar reflectivity, Doppler velocity, and Doppler spectral width are compared with a lookup table of forward modeled cloud radar spectra holding microphysical and observable quantities. In this way, ICNC and the distribution of maximum particle diameter are estimated. The uncertainty of the retrieved values is estimated by varying the measurement values within their measurement errors and thereby retrieving many retrieval results at once. The most common value in the resulting distribution of retrieved values is considered the result and the width of the distribution represents its uncertainty. Within this study, the uncertainty in ICNC is a factor of around 4 (see Ramelli et al., 2021a).

2.3 In situ cloud measurements

During the campaign, cloud properties were also measured in situ using the tethered balloon system HoloBalloon (Ramelli et al., 2021a, b). On HoloBalloon, a HOLographic Imager for Microscopic Objects (HOLIMO) is installed which images an ensemble of cloud particles in the size range from 6 µm and 2 mm in a well-defined sample volume to obtain phase-resolved cloud properties (Ramelli et al., 2020). The captured particles are classified into cloud droplets, ice crystals, and artifacts using a convolutional neural network (Touloupas et al., 2020) for particles larger than 25 µm and a decision tree for particles smaller than 25 µm. The differentiation between cloud droplets and ice crystals is only done for particles larger than 25 µm based on the particle shape (circular versus non-circular), whereby all particles smaller than 25 µm are classified as cloud droplets. All ice crystals predicted by the neural network are manually confirmed or reclassified to ensure a high classification accuracy of ice crystals. Following this approach, the ICNC can be estimated with an uncertainty of ±5 % for ice crystals larger than 100 µm and ±15 % for ice crystals smaller than 100 µm (Beck, 2017). The IWC (in kg m⁻³) is calculated using the mass–diameter relationship given in Cotton et al. (2013):

$$\begin{array}{}\text{(5)}& \mathrm{IWC}=\left\{\begin{array}{ll}\sum _i\frac{\mathit{\pi }}{\mathrm{6}}{\mathit{\rho }}_{\mathrm{ice}}{D}_{\mathrm{i},\max }^{\mathrm{3}}\times V_{\mathrm{cloud}}^{-\mathrm{1}}& \text{for }{D}_{\max }\le \mathrm{70}\mathrm{\mu }\text{m}\\ \sum _i\mathrm{0.026}{D}_{\mathrm{i},\max }^{\mathrm{2}}\times V_{\mathrm{cloud}}^{-\mathrm{1}}& \text{for }{D}_{\max }>\mathrm{70}\mathrm{\mu }\text{m},\end{array}\right.\end{array}$$

where D_max is the maximum dimension of the ice crystals (in m), ρ_ice is the effective ice density of small ice (700 kg m⁻³), and V_cloud is the cloud sampling volume (in m³). Due to the underlying assumptions involved in the calculation of the IWC (e.g., mass–diameter relationship, uncertainties in the effective ice crystal density), the uncertainty in the IWC can be up to a factor of 2 (Beck, 2017). In a comparison to radar observations (Sect. 3.3), we use IWCs and ICNCs measured with HoloBalloon on 8 March 2019, for which further details on the measurements and the synoptic situation are provided in Ramelli et al. (2021a).

2.4 Cloud temperature data

Temperatures at cloud top and within the clouds were retrieved from the COSMO-1 analysis data at corresponding heights above WOP. The COSMO-1 data were provided by MeteoSchweiz for the processing of the remote sensing data.

2.5 Investigated INP parameterizations

DeMott et al. (2015)Ullrich et al. (2017)Harrison et al. (2019)DeMott et al. (2010)Ullrich et al. (2017)

Table 1 Overview of INP parameterizations for immersion freezing mode considered in this study along with the aerosol type they are applicable to, the aerosol property they are based on (n₂₅₀: aerosol number concentration of particles with radii ≥ 250 nm; s: aerosol surface area concentration), the applicable temperature range, and the explicit form equation number in this paper. This overview was adapted from the more extensive summary given by Marinou et al. (2019, Table 1).

^∗ Note that D10 was not primarily developed for predicting continental INP concentrations and includes samples from different field observations in North America, the Amazon, and the Pacific also featuring dust-carrying air masses. It was shown in recent publications (Mamouri and Ansmann, 2016; Marinou et al., 2019) that it is suitable to deduce continental INP concentration from lidar observations.

INP concentrations can be obtained from remote sensing data by applying a suitable INP parameterization to retrieved aerosol properties. By now, a wide range of INP parameterizations specific to a certain aerosol type (e.g., dust, continental, marine, soot) and a freezing mode (e.g., immersion or deposition) exist. Here, we compare lidar-estimated INP concentrations with in situ observations in immersion mode at temperatures $\ge -\mathrm{20}$ ^∘C during times of a Saharan dust event and otherwise continental aerosol over the region of Davos. In Table 1, we present the investigated existing parameterizations for immersion freezing which are applicable to air masses during our observations. The overview is based on the summary of Marinou et al. (2019), which gives detailed information about the parameterizations (Marinou et al., 2019, Sect. 2). In the presence of Saharan-dust-carrying air masses, we investigate the predictions of INP concentration (given in StdL⁻¹) from the parameterizations of DeMott et al. (2015) (D15) defined as

$$\begin{array}{}\text{(6)}& n_{\mathrm{INP},\mathrm{D}\mathrm{15}}(T,n_{\mathrm{250}})=\mathrm{cf}\times (n_{\mathrm{250}}{)}^{\mathrm{1.25}}\times e^{\mathrm{0.46}\times (\mathrm{273.16}-T)-\mathrm{11.6}},\end{array}$$

with cf being a calibration factor (unity by default); Ullrich et al. (2017) (dust parameterization, immersion freezing mode, U17d) defined as

$$\begin{array}{}\text{(7)}& n_{\mathrm{INP},\mathrm{U}\mathrm{17}\mathrm{d}}(T,s)={\mathrm{10}}^{-\mathrm{9}}\times s\times e^{-\mathrm{0.517}\times T+\mathrm{150.577}};\end{array}$$

and Harrison et al. (2019) (H19) defined as

$$\begin{array}{}\text{(8)}& n_{\mathrm{INP},\mathrm{H}\mathrm{19}}(T,s)=f_{\mathrm{K}}\times {\mathrm{10}}^{-\mathrm{5}}\times s\times {\mathrm{10}}^{{\sum }_{i=\mathrm{0}}^{\mathrm{5}}{c}_i\times (T-\mathrm{273.15}{)}^i}\end{array}$$

where f_K is the fraction of K-feldspar on the particle surface, $c_{\mathrm{0}}=-\mathrm{3.25}$, $c_{\mathrm{1}}=-\mathrm{0.793}$, $c_{\mathrm{2}}=-\mathrm{6.91}\times {\mathrm{10}}^{-\mathrm{2}}$, $c_{\mathrm{3}}=-\mathrm{4.17}\times {\mathrm{10}}^{-\mathrm{3}}$, $c_{\mathrm{4}}=-\mathrm{1.05}\times {\mathrm{10}}^{-\mathrm{4}}$, and $c_{\mathrm{5}}=-\mathrm{9.08}\times {\mathrm{10}}^{-\mathrm{7}}$. The latter parameterization is a more recent dust parameterization scaling with the relative contribution of K-feldspar, which was found to correlate most strongly with the INP concentration in their study. Note that temperature T is given in kelvin, and both n₂₅₀ (given in Stdcm⁻³) and s (given in µm² Stdcm⁻³) refer to dry particle dimensions. In the other cases, which are considered air masses carrying continental aerosol, we investigate the INP concentration prediction based on DeMott et al. (2010) (global parameterization, D10) defined as

$$\begin{array}{}\text{(9)}& \begin{array}{rl}& n_{\mathrm{INP},\mathrm{D}\mathrm{10}}(T,n_{\mathrm{250}})=\mathrm{5.94}\times {\mathrm{10}}^{-\mathrm{5}}\times (\mathrm{273.16}-T{)}^{\mathrm{3.33}}\\ & \times (n_{\mathrm{250}}{)}^{\mathrm{0.0264}\times (\mathrm{273.16}-T)+\mathrm{0.0033}}\end{array}\end{array}$$

and Ullrich et al. (2017) (soot parameterization, immersion freezing mode, U17s) defined as

$$\begin{array}{}\text{(10)}& \begin{array}{rl}& n_{\mathrm{INP},\mathrm{U}\mathrm{17}\mathrm{s}}(T,s)=\mathrm{7.463}\times {\mathrm{10}}^{-\mathrm{9}}\times s\\ & \times e^{\mathrm{0.7667}-\mathrm{0.8525}\times (T-\mathrm{273.15})-\mathrm{0.0101}\times (T-\mathrm{273.15}{)}^{\mathrm{2}}}.\end{array}\end{array}$$

The latter was developed on soot aerosol, thus making it suitable to specifically capture anthropogenic contributions to continental aerosol. Notably, the application ranges of all parameterizations presented in Table 1 are applicable mainly at temperatures $\le -\mathrm{15}$ ^∘C. Therefore, we extrapolate the parameterizations to −5 ^∘C in accordance with Marinou et al. (2019).

3 Results and discussion

Atmospheric INP concentrations and IMF are estimated from the lidar and radar measurements. Three main uncertainties are involved in the estimation: (i) uncertainties linked to the measurement of the aerosol extinction coefficient and its conversion to number or surface area concentration, (ii) uncertainties linked to the INP parameterization itself, and (iii) using a parameterization not suitable for the dominant aerosol constituent. In the following, we first validate the aerosol properties as input and consequently determine the most suitable INP parameterization for dust and continental aerosol during RACLETS.

3.1 Aerosol concentration comparison

Figure 2 Comparison of aerosol properties observed in situ and with lidar. (a) Time series of aerosol number concentration with radii ≥250 nm (n₂₅₀) measured in situ (black) and with lidar (orange) at WFJ. (b) Same as (a) but measured at WOP (lidar data in blue). (c) Time series of surface area concentration (s) measured in situ (black) and with lidar (orange) at WFJ. (d) Same as (c) but measured at WOP (lidar data in blue). Right of each time series panel, a scatter plot compares the in situ to the lidar observations of the data presented on the left. Uncertainties are indicated by the shading for the time series and at representative data points of the scatter plots in black. Note the lidar retrieval height at WOP was 2052 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ (see Fig. 1).

In the following, we investigate the accuracy of the lidar retrieval of aerosol properties with in situ observations. For the comparison of lidar observations to in situ observations at WFJ (mountaintop site) and WOP (high valley site), lidar retrievals were taken from the closest height bins at 2695 and 2052 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$, respectively (Fig. 1b). In the case of WOP, the lowest bin with complete overlap of the lidar (see e.g., Wandinger and Ansmann, 2002) was taken for the comparison, which is still around 400 m above the in situ site. The incomplete overlap of the emitted and received beam is a general issue in the comparison of lidar measurements to ground-based in situ observations. The comparison to WFJ, more than 1000 m above the lidar site where complete overlap is surely reached, is a great advantage of the present study.

Figure 2a, b show n₂₅₀ from in situ observations and lidar retrievals for WFJ and WOP, respectively. Note that the accuracy will be discussed only qualitatively. A quantitative assessment and improvement of the lidar aerosol retrieval is beyond the scope of this publication. For both locations in situ and lidar observations agree qualitatively well, in particular for higher concentrations (dusty conditions). At relatively low aerosol concentrations and in the presence of continental aerosol, a plateauing of the lidar-retrieved aerosol concentrations was observed. The clear atmosphere over the Alps with very low values of the extinction coefficient (<10 Mm⁻¹) could be responsible for deviations from the assumed linear relationship of extinction to n₂₅₀ and s, respectively, (Mamouri and Ansmann, 2016). The larger diurnal variability between the in situ observation and lidar retrieval at WOP (high valley site) compared to WFJ (mountaintop site) can be explained by the diurnal changes of aerosol concentration near the ground (in situ observations) not affecting the air masses on the lidar retrieval height (height difference approx. 400 m; see Fig. 1b). This difference in height and therefore air mass commonly limits a quantitative conclusion between ground-based in situ observations and remote sensing instruments as the well-mixed boundary layer could at times not extend up to the lowest retrieval height. For the retrieval of s (Fig. 2c, d) the aforementioned observations hold equally true, which is not surprising as the surface area relates to the square of particle radius. However, comparing the retrieval accuracy at WFJ (Fig. 2a, c), a stronger bias of the lidar-retrieved surface area concentrations is apparent, which has also been previously reported by Haarig et al. (2019), based on similar observations carried out at Barbados.

Despite the fact that the comparison in this section was purely qualitative, two important insights are gained: (i) due to the better agreement of in situ observations and lidar retrieval for n₂₅₀ compared to s, INP parameterizations based on n₂₅₀ should be preferred. (ii) Comparing ground-based in situ observations (collocated with the lidar instrument) does not allow for a comparison even to the lowest lidar retrieval height. Consequently, we restrict the comparison to in situ INP concentrations in the following to observations made at WFJ (mountaintop site), where the lidar retrieval overlaps the in situ observations.

Figure 3 Overview of the Saharan dust event over the Davos region starting on 19 February 2019 and the continental aerosol conditions before. (a) Lidar attenuated backscatter coefficient (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$) measured above WOP (high valley site). Note that the gap in data on 20 February is due to technical maintenance of the lidar system. Retrieval height of WFJ (mountaintop site) is indicated in dashed white. (b) Lidar attenuated backscatter coefficient (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$) retrieved from the lidar height bin closest to WFJ (2695 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$, Fig. 1b). (c) In situ aerosol number concentrations for particles with radii ≥ 250 nm (n₂₅₀, left axis) at WFJ (dashed orange) and WOP (dashed blue) and INP concentrations at −13 ^∘C ($n_{\mathrm{INP},-\mathrm{13}{}^{\circ }\mathrm{C}}$, right axis) measured in situ at WFJ. INP samples during continental background are indicated with crosses, and samples during Saharan dust are indicated with hexagons. Uncertainties in INP concentration are indicated for representative samples in black.

3.2 INP concentration comparison

Throughout the campaign, the region around Davos was mostly exposed to continental aerosol (e.g., 14–18 February 2019 as seen in Fig. 3). The situation changed 19–22 February 2019, when the synoptic wind situation promoted transport of Saharan dust from North Africa towards Davos. The arrival of the dust plume is clearly visible in the attenuated backscatter coefficient (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$) of the lidar (Fig. 3a, b). Additionally, an increase by more than an order of magnitude in n₂₅₀ occurred at both sites on 19 February (in situ observations, Fig. 3c). Before 19 February, the higher variability of n₂₅₀ at WOP (high valley site) compared to WFJ (mountaintop site) can be explained by local sources and accumulation of aerosol beneath a nighttime inversion along with valley and mountain breezes (Wieder et al., 2022b). Over the course of 19 February, n₂₅₀ at WOP followed the same trend and magnitude as at WFJ and continued doing so in the following days. This drastic change in trend underlines that the region was affected and dominated by the long-range-transported Saharan dust. INP concentrations at −13 ^∘C measured in situ at WFJ also showed an increase of approximately 1 order of magnitude already on 18 February (a day before the strong Saharan dust signal in Fig. 3a, b). INP concentrations could have already been higher on 17 February; however, no in situ measurements are available on that day. It is conceivable that air masses carrying only a smaller fraction of Saharan dust had already influenced the region of Davos. Analysis of the air mass properties and history prompted us to distinguish continental aerosol and Saharan dust for 19 February 11:30 UTC. The FLEXPART dispersion model (Pisso et al., 2019) indicated that air masses before that time originated from eastern Europe and Italy (not shown). Between 11:00 and 12:00 UTC (time corresponding to measurements at Davos) the flow over southern Italy was affected by intrusion of air masses from North Africa. To identify when exactly air masses from North Africa carrying Saharan dust entered the flow, we utilized the relative humidity measured at WFJ (Fig. A1a in the Appendix). A local minimum in relative humidity was identified on 19 February at 11:30 UTC. Therefore, we classify samples before that time as still continental. Nonetheless, the Jungfraujoch Research Station (approximately 140 km WSW of Weissfluhjoch) released a Saharan dust event warning only at 13:37 UTC. Brunner et al. (2021) reported that the determination of a Saharan dust event varies based on the underlying parameter used as proxy. The three samples with higher INP concentrations (orange crosses on 18 and 19 February in Fig. 3c) coincided with local maxima in relative humidity and attenuated backscatter (indicated by 1 and 2 in Fig. A1a, b, respectively). Therefore, the higher INP concentrations observed at WFJ at the end of the continental aerosol period could be a result of a stronger exchange with air masses at lower height resulting in uptake of biogenic and soil particles from fertile lands over southern Italy and eastern Europe (Conen et al., 2015).

Table 2 Overview of the number of samples (N) available for INP closure in dust-dominated and continental air mass cases. Additionally, the number of days over which the samples are taken is given.

^∗ Note that due to technical maintenance of the lidar instrument only 9 of the 11 total samples during the Saharan dust event shown in Fig. 3 could be used.

In the following we will investigate the performance of INP parameterizations for dust-carrying and continental aerosol air masses. For the dust comparison (Sect. 3.2.1) we use 9 of the 11 in situ samples indicated by hexagons in Fig. 3c. For the two samples in the morning of 20 February, no lidar data are available for a comparison. For the comparison with continental aerosol (Sect. 3.2.2) we use the 14 samples presented with crosses in Fig. 3c and additional samples that were collected individually during clear-sky conditions before and after cloud events between 23 February and 12 March 2019, resulting in a total of 23 samples collected on 14 d (see Table 2). The grouping into Saharan dust and continental samples is consistent with Mignani et al. (2021) (see Fig. S1 in the Supplement of Mignani et al., 2021).

In the in situ and lidar observations, four sources of uncertainty are present.

• The first source is the measurement uncertainty of the in situ observation (drop freezing technique), which depends on various parameters such as freezing temperature and number of droplets. Vali (1971) states the uncertainty to be between a factor of 2 and 4 (based on an assay of 150 droplets).

• The second source is the retrieval uncertainty of the lidar that are passed on to the INP concentration through the applied parameterization. The retrieval uncertainty of the INP concentration from the lidar measurements is given by a factor of 3–10 (Haarig et al., 2019; Mamouri and Ansmann, 2016).

• The third source is air mass differences due to the horizontal distance between in situ measurements at WFJ (mountaintop site) and the lidar beam, which was 3.65 km (Fig. 1b). However, being located at mountaintop height, it is conceivable to assume no drastic difference between the two locations. Thus, the uncertainties caused by the spatial distance are negligible compared to the other sources of uncertainties.

• Ultimately, the natural variability of INP is typically referred to as 1 order of magnitude (Kanji et al., 2017).

The latter being the largest source of uncertainty, we consequently view data points of in situ measurements and lidar-based observations within 1 order of magnitude as acceptable.

Figure 4 Comparison of in situ INP concentrations ($n_{\mathrm{INP},\mathrm{obs}.}$) measured at WFJ (mountaintop site, Fig. 1) to lidar-retrieved INP concentrations (n_INP,lidar) during the Saharan dust event (19–21 February 2019) for three dust parameterizations: (a) D15 (DeMott et al., 2015), (b) D15 with a calibration factor (cf=0.086) applied as proposed by Schrod et al. (2017), (c) H19 (Harrison et al., 2019), and (d) U17d (Ullrich et al., 2017). Notably, we compare INP concentrations at temperatures between −20 ^∘C (dark orange) and −5 ^∘C (light orange). The 1:1 line is shown in dashed gray. The dashed red lines show a linear regression fit through the logarithmically transformed data points, including the 95 % confidence interval in dotted cyan. Shown in every panel are the mean absolute error (MAE), the fraction of data points falling within a factor of 2 and 5 with respect to the 1:1 line (f₂ and f₅), and the slope of the fit (including the standard error). The gray shaded area indicates a deviation of 1 order of magnitude from the 1:1 line, which is the upper limit uncertainty of the lidar-based INP retrieval (Mamouri and Ansmann, 2016). For D15 (a), the data fit is also presented, applying the calibration factor from (b) for comparison (black dashed). For H19 (c), the data fit is also presented, assuming a K-feldspar fraction of 10 % (black dashed) and 1 % (black dotted), respectively. Uncertainties in INP concentration are indicated for representative samples in black.

3.2.1 Retrieval during times of Saharan dust

In situ INP observations are compared against lidar-retrieved INP concentrations using different INP parameterizations to evaluate their performance during Saharan dust presence. In Fig. 4, we present comparisons to the parameterizations of D15, H19, and U17d (Table 1) along with performance measures such as the mean absolute error (MAE), the fraction of data points falling within a factor of 2 and 5 with respect to the 1:1 line (f₂ and f₅), and the slope of a linear regression performed on the logarithmically transformed data. Among all investigated parameterizations, D15 (Fig. 4a) results in a slope of the linear regression fit that is closest to the desired unity. However, D15 overestimates the INP concentration, which could partially be attributed to the slight overestimation in n₂₅₀ retrieval (Fig. 2a). Schrod et al. (2017) proposed a calibration factor of cf=0.086 based on a lidar comparison with Saharan dust samples collected with UAVs over Cyprus. Applying this calibration factor to the data pushes the data beyond the 1:1 line, reducing the MAE and increasing f₂ and f₅ (Fig. 4b). Previously, Mignani et al. (2021) compared INP observations at −15 ^∘C of the same samples to predictions of D15 using in situ aerosol data and applying the calibration factor of Schrod et al. (2017). In this study, we extend the analysis to the whole temperature spectrum (−20 to −5 ^∘C) and relied on aerosol retrievals from the lidar. For concentrations between 10⁻¹ and 10¹ StdL⁻¹, which is the concentration range Schrod et al. (2017) used to determine the calibration factor, the fit cuts the 1:1 line (black dashed extension of the red dashed data fit). Note that Schrod et al. (2017) also proposed a slope correction; however, we did not correct the slope because all data points already fall within the uncertainty range (gray shaded area). H19 (Fig. 4c) results in a higher f₂ than D15 (using cf=0.086) but with a more shallow slope leading to a lower f₅. For the retrieval we assumed a K-feldspar fraction of 20 % in accordance with Kandler et al. (2011), who investigated the relative K-feldspar fractions on Cabo Verde west of central Africa. H19 scales log-linearly with the percentage of K-feldspar. During measurements at Barbados Kandler et al. (2018) found that Saharan dust featured only 1 % K-feldspar. However, using 1 % K-feldspar results in an underestimation of INP concentration (Fig. 4c). A fraction of approximately 10 % K-feldspar would match our observations best, yet the slope of the relation would be too small (0.51). U17d (Fig. 4d) overestimates the INP concentration the most strongly, leading to the highest MAE (1.93) among the investigated parameterizations.

Figure 5 Same as for Fig. 4 but during continental background (see Fig. 3) for the two parameterizations: (a) D10 (DeMott et al., 2010), and (b) U17s (Ullrich et al., 2017).

Our comparison suggests D15 in combination with the calibration factor (cf=0.086) proposed by Schrod et al. (2017) as the preferable choice for INP concentration retrieval from dust air masses (slope closest to unity, lowest MAE). Moreover, D15 should be preferred due to its dependence on n₂₅₀ instead of s, which is associated with lower lidar retrieval uncertainty (see Sect. 3.1).

3.2.2 Retrieval during times of background aerosol

Aside from a Saharan dust event, the measurement region was not affected by a long-range-transported dominating aerosol species. Thus, the samples collected outside the Saharan dust episode represent a general mix of continental aerosols. In Fig. 5, we compare the in situ observed INP concentration of these samples with lidar retrievals using the parameterizations (a) D10 for global aerosol and (b) U17s for soot, which are the best parameterizations available to mimic continental aerosol (Table 1). D10 results in a large MAE (1.84) and with no data points agreeing within a factor of 2 or 5 ($f_{\mathrm{2}}=f_{\mathrm{5}}=\mathrm{0}$). Yet, the slope is close to unity (1.05). INP concentrations based on U17s result in smaller MAE (0.66) and higher f₂ (0.27) and f₅ (0.53), but the slope is more shallow (0.47). Purely based on the error assessment of the parameterizations in their original form, U17s would be the parameterization of choice for continental INP concentration lidar retrieval. The slope of the D10 comparison being close to unity indicates that D10 captures the composition of a continental aerosol mix. Together with the large MAE, the comparison suggests that the expected fraction of INPs among the aerosol mix is overestimated. In the next section, we propose a tuning of D10 to improve the retrieval of continental INP concentration from lidar measurements.

3.2.3 Tuning of D10 based on in situ INP and aerosol measurements

Figure 6 Same as for Fig. 5a but comparing in situ INP concentrations ($n_{\mathrm{INP},\mathrm{obs}.}$) measured at WFJ (mountaintop site, Fig. 1) to predicted INP concentrations (n_INP,APS) based on D10 now applied to in situ aerosol concentrations (from APS) during continental background (see Fig. 3). (a) D10 applied to in situ aerosol concentrations (from APS) in its standard form. (b) D10 shifted towards the 1:1 line, maximizing the fraction of data points within 1 order of magnitude around the 1:1 line (see Eq. 11).

The slope in the comparison between in situ observations and lidar-retrieved INP concentrations suggests that D10 captures the feature of a continental INP mix (Fig. 5a). The large overprediction by nearly 2 orders of magnitude on average (MAE = 1.84) indicates an overestimation in INP active aerosol fraction. D10 was developed based on nine datasets from different locations (e.g., Alaska, Pacific, Amazon Basin) and thus contains various aerosol types (DeMott et al., 2010). Although a large number of non-dust samples contributed to the parameterization, dust samples, e.g., from the Pacific Dust Experiment (PACDEX, Stith et al., 2009), were included in the D10 parameterization. Mineral dust is more ice active compared to (non-biogenic) continental aerosol. Therefore, it is not surprising that D10 overestimates the INP concentration when used for non-dust continental aerosol. Using the collocated in situ observations at WFJ (mountaintop site) of n₂₅₀ (from APS) and INP concentration (Sect. 2.1), we could derive a tuned D10 parameterization optimized for continental aerosol. The optimization is based purely on the in situ observations in order to be independent of the lidar data (which are the test data) as we cannot fully rule out a difference in air mass. Analogously to the lidar retrieval, INP concentration is predicted by n₂₅₀ (now taken from the APS) and D10 applied to it. The comparison of the obtained INP concentrations (Fig. 6a) exhibits very similar performance measures as in the lidar comparison (Fig. 5a), underlining that the difference in air mass was minuscule. Furthermore, as discussed in Sect. 2.1, the use of a heated inlet could potentially lead to the loss of volatile compounds. This could subsequently decrease the obtained n₂₅₀, causing an underestimation of in situ predicted INP concentrations. We note that this effect does not apply to lidar-retrieved aerosol number concentrations and therefore lidar-retrieved INP concentrations. Therefore, the similarity of the performance measures also indicates the validity of the assumption that evaporation of volatile compounds is inconsequential. We optimize D10 with a multiplicative calibration factor (cf) as a free parameter (similarly to DeMott et al., 2015; Schrod et al., 2017) for maximizing the number of data points within an order of magnitude around the 1:1 line (gray shaded area in Fig. 6b). Based on the performance parameters it is apparent how the prediction has improved while the slope remained unchanged. To obtain the tuned INP concentrations ($n_{\mathrm{INP},\mathrm{D}\mathrm{10},\mathrm{tuned}}$) from INP concentrations using D10 (n_INP,D10, as defined in Eq. 9), we calculated

$$\begin{array}{}\text{(11)}& n_{\mathrm{INP},\mathrm{D}\mathrm{10},\mathrm{tuned}}=\mathrm{cf}\times n_{\mathrm{INP},\mathrm{D}\mathrm{10}},\end{array}$$

with the calibration factor cf=0.0204 and the INP concentrations given in StdL⁻¹. The calibration factor indicates that under clean continental aerosol conditions just 2 % of the INPs predicted with D10 are present. The improvement in prediction that was found using the in situ aerosol data is also found for the lidar retrieval when applying the calibration factor (Fig. 7). The now tuned D10 parameterization outperforms U17s (compare Figs. 5b and 7b). Thus, we propose the tuned D10 parameterization (Eq. 11) to retrieve atmospheric INP concentrations from wintertime continental air masses.

Figure 7 Same as for Fig. 4 but comparing in situ INP concentrations ($n_{\mathrm{INP},\mathrm{obs}.}$) measured at WFJ (mountaintop site, Fig. 1) to lidar-retrieved INP concentrations (n_INP,lidar) using (a) the standard D10 parameterization and (b) the tuned D10 parameterization (see Eq. 11) for continental aerosol. Note that (a) is the same as Fig. 5a.

Figure 8 Comparison of the remotely retrieved and in situ measured (a) ice water content (IWC) and (b) ice crystal number concentration (ICNC, n_ice as variable) on 8 March 2019. The IWC data were retrieved from the lowest size bin (1783 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). Due to measurement gaps, ICNC between 1900 and 2500 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ was averaged (median). Left: time series of available remotely retrieved data (blue) and in situ observed data (black). Right: scatter plot of data during times the balloon was at an elevation higher than 1750 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ (gray shading in the time series plots, maximal elevation: 1910 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). The gray shading in the scatter plots indicates an error of 1 order of magnitude. The values f₁₀ and f₃ indicate the fraction of data points within a factor of 10 and 3, respectively. Uncertainties are indicated by the shading for the time series and at representative data points of the scatter plots in black.

3.3 Validation of remotely retrieved cloud properties

In this section, remotely retrieved ice water contents (IWCs) and ice crystal number concentrations (ICNCs) are compared with balloon-borne in situ observations to validate the accuracy of the retrievals. On 8 March 2019, the balloon performed vertical profiles between 1650 and 1900 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ next to the remote sensing instruments (see overview plot (Fig. 6b) in Ramelli et al., 2021a). Since the lowest data bin of the radar and lidar was at 1783 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$, we considered in situ data only for times when the balloon was at an elevation of 1750 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ and higher (gray shading in time series plots of Fig. 8, maximum elevation 1910 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). IWC is deduced solely from radar observations and could continuously be taken from the lowest height bin (1783 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). The retrieval of ICNCs needs parallel measurements of the radar and lidar. As the lidar featured data gaps, it was necessary to average (median) the ICNC data between 1900 and 2500 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ The uncertainty of remotely retrieved cloud properties is typically given by a factor of 3 (Bühl et al., 2019). Within a factor of 3, 68 % and 46 % of the data points agreed for IWC (Fig. 8a) and ICNC (Fig. 8b), respectively. For a factor of 10, the agreement was 92 % and 88 % for IWC (Fig. 8a) and ICNC (Fig. 8b), respectively. Given the inhomogeneous distribution of hydrometeors within an MPC (Korolev et al., 2017) and that the balloon could have been horizontally displaced up to 200 m from the remote sensing beams, we assess the comparability to be satisfactory. In the following analysis we will consider remotely retrieved ICNC associated with an uncertainty of 1 order of magnitude.

3.4 IMF estimation in orographic MPCs

Here we propose a method to assess secondary ice production (SIP) on a single cloud basis by calculating ice multiplication factor (IMF) histograms using combined lidar and radar retrievals. First, the procedure is explained for one observed cloud case (Sect. 3.4.1). Second, the effect of sublimation below cloud base on ICNC and the temporal evolution of ICNC profiles is presented. Observed IMF histograms with height and temperature are presented and discussed with regard to the ICNC input data (Sect. 3.4.2). Third, the limitations and caveats of the method are discussed (Sect. 3.4.3). Finally, the results of this analysis are compared to observations of previous field measurements (Sect. 3.4.4).

Figure 9 Exemplary illustration of data processing workflow to obtain in-cloud ice multiplication factors (IMFs) for the cloud observed on 12 March 2019. (a) Reflectivity of the cloud radar is given in shades of gray. Ice crystal number concentration (n_ice) is overlaid in color. The limit of used data for the IMF calculation with height is indicated in red (height span h_d, time span t_d). Additionally, the ambient temperature at the heights of the upper and lower bounds of selected data is given in red (see T_dt and T_db in Table D1). The INP concentration (n_INP,ct) was taken from the last lidar retrieval (solid green shading) at cloud top (dashed blue line). (b) Attenuated backscatter coefficient at 1064 nm (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$) from the lidar in color. Cloud base estimates determined from the attenuated backscatter coefficient at 532 nm (${\mathit{\beta }}_{\mathrm{532}\mathrm{nm}}^{\mathrm{att}.}$) are indicated by red dots. Regions where the Cloudnet algorithm detected liquid water are indicated in bright green. The height of cloud top is indicated by the dashed blue line. The region within the yellow rectangle was classified as below cloud base and the region within the cyan rectangle as above cloud base. Overviews of all individual cloud events can be found in Figs. B1 and C1.

3.4.1 Procedure description

Locally in a cloud, the magnitude of ice crystal enhancement can be estimated with the IMF defined as ICNC divided by INP concentration. Note that the effect of blowing snow is neglected, and sedimenting ice crystals from a cloud layer aloft are ruled out based on manual inspection of cloud radar observations. For this analysis, only ICNCs in ice saturated regions were considered, as otherwise ICNCs can decrease due to sublimation. The INP concentration is indicative of the number of primary ice crystals formed in the cloud that could subsequently evoke SIP. Taking a too low INP concentration results in an overestimation of the IMF. The (cumulative) INP concentration increases with decreasing temperature, resulting in the highest INP concentrations at cloud top (assuming that the entire cloud formed by the same ascending air mass). ICNCs at cloud top are generally created due to primary ice formation (see e.g., Crosier et al., 2014). At sufficient size, the primary ice crystals will no longer be levitated by the updrafts and start to sediment to lower heights. Thus, using INP concentrations at cloud top is an upper estimate for the primary ice crystals in the cloud and avoids an overestimation of the IMF. In the next section, we determine IMFs for seven individual cloud events based on the following procedure.

1. Based on the aforementioned considerations, we retrieve the INP concentration at cloud top height from the lidar retrieval closest to the cloud event (Fig. 9a), using the tuned D10 parameterization (Eq. 11), which is the best estimate using ground-based lidar. ICNCs were taken over a height h_d (Fig. 9a), excluding heights near the radar echo boundaries to avoid retrieval errors and to avoid zones of cloud top entrainment. Additionally, data were trimmed to a time t_d (Fig. 9a) to account for echo-free regions at the beginning and end of the clouds. Echo-free regions in the selected data frame (e.g., at a height of approximately 3.6 $\mathrm{km}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ at around 13:40 UTC in Fig. 9a) are thought to be negligible for the consequent statistical analysis.

2. Ice crystals sedimenting out of the cloud will start to sublimate upon entering the unsaturated environment below cloud base. Thus, in a second step, lidar variables are used to determine the cloud base height to discard data from regions where sublimation could occur. In our study, the classification was done manually by inspecting (i) the attenuated backscatter coefficient at 1064 nm (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$), (ii) a cloud base estimate derived from the attenuated backscatter coefficient at 532 nm (${\mathit{\beta }}_{\mathrm{532}\mathrm{nm}}^{\mathrm{att}.}$), and (iii) the presence of liquid water flag by the Cloudnet algorithm (Fig. 9b). In most cases criteria (ii) and (iii) were used to determine cloud base because these estimates are deemed to be more reliable than criterion (i). However, in cases where these estimates provided frequent jumps in cloud base, for example, due to snowfall blurring the backscattered signal (see e.g., 15:00–15:30 UTC in Fig. 9b), we used (i) as a proxy for the cloud base. As the cloud base height in general changes with cloud evolution, the data were divided into time slices of approximately 1 h for which individually averaged cloud base heights were determined. The individual sizes of the time slices were set to match distinctive regime changes (e.g., the lowering of the cloud base in Fig. 9b). Consequently, each time slice was vertically divided and classified into below cloud base and above cloud base. Strong precipitation can hamper the lidar signal close to the instrument, disallowing cloud base determination above. Time slices during strong precipitation or without lidar information are therefore classified as unknown.

3. The IMF at a given height h and time t is consequently calculated as

   $$\begin{array}{}\text{(12)}& \mathrm{IMF}(h,t)={\displaystyle \frac{n_{\mathrm{ice}}(h,t)}{n_{\mathrm{INP},\mathrm{ct}}}},\end{array}$$

   with n_ice(h,t) being the ICNC at height h and time t and n_INP,ct being the INP concentration at cloud top (as illustrated in Fig. 9a).

Figure 10 Median ice crystal number concentrations (ICNCs) as a function of altitude for each cloud and for time slices of approximately 1 h. (a) Temporal evolution of the ICNC (n_ice) profiles: the faintest colors correspond to the first slice of data, and the strongest colors correspond to the last slice of data. ((b) Cloud regime classification for each time slice: ICNCs (n_ice) taken above cloud base (solid), below cloud base (dashed), and for an unknown cloud base (faint color). For cloud 3, the cloud base could not be determined.

3.4.2 IMF observations of seven orographic MPCs

IMFs were calculated for seven cloud events according to the procedure described above. An overview of the ICNC data and the structure of each cloud is presented in Appendix B. The time slice division and classification are provided in Appendix C. Our methodology uses a fixed INP concentration for all IMFs (Eq. 12). Thus, changes in the IMFs are driven by the evolution of the ICNCs. In a generalized case of cloud evolution, one would expect the ICNCs to increase with time due to SIP as long as conditions for SIP prevail. At cloud top primary ice crystals will nucleate, start growing, and eventually sediment. Due to the sedimentation, ICNC increases downward. If the environmental conditions permit, SIP processes start to enhance the ICNC. With the presence of more ice crystals and warmer temperatures, aggregation becomes more likely, reducing the ICNC towards cloud base and at higher temperatures. With ice crystals sedimenting below cloud base, not only aggregation but also sublimation of entire ice crystals could reduce the ICNC. In this simple picture, the recirculation of ice crystals back into the cloud due to updrafts is not considered. In reality it is not possible a priori to know where in the cloud ICNCs increase or decrease as this strongly depends on the environmental conditions and hence varies from cloud to cloud and throughout a cloud's lifetime.

Figure 10 gives an overview of the median ICNC profiles for each of the seven observed clouds in this study as time slices of approximately 1 h. A general trend of increasing ICNC with time manifests itself across different clouds (e.g., clouds 2, 3, 5, 6, and 7). For clouds 5 and 6 ICNCs increase from cloud top downwards before they start to decrease towards the ground. In contrast, the ICNC profiles of the other clouds (clouds 1, 2, 3, 4, and 7) generally decrease with decreasing height. Note that the top of the profile is below cloud top (see discussion in Sect. 3.4.1). Above cloud base, aggregation can reduce the ICNC. Below cloud base, sublimation can reduce the ICNC in addition to aggregation (Fig. 10b). Despite no cloud base determination, ICNC profiles of clouds 5 and 6 during strong precipitation increase with time near the surface (especially for cloud 6), suggesting that the relative humidity was close to 100 % with respect to ice everywhere below cloud. Note that for cloud 5 in the last hour, there could be ice crystals originating from homogeneous freezing sedimenting from higher up (see Fig. B1d). Cloud base for cloud 3 could not be determined due to a liquid layer in a lower cloud fully attenuating the lidar signal (see Fig. A1b). It is worth mentioning that the ICNCs do not differ strongly between precipitating and non-precipitating clouds.

Figure 11 Probability density functions (PDFs) of obtained ice multiplication factors (IMFs, Eq. 12) for ice crystal number concentration data above cloud base (solid black), below cloud base (dashed black), and of unknown cloud base (solid gray). The solid red line indicates unity. The median IMF and the percentage of IMF larger unity for each PDF are given in the legend.

Following the definition of the IMF (Eq. 12), those below unity are linked to regions in the cloud where the INP concentration at cloud top exceeds the local ICNC. Potential regions of IMF close to and lower than unity are (i) in the beginning stage of cloud formation and (ii) near cloud top, where primary ice crystals have just formed and SIP processes are not active yet. Furthermore, ICNC and thus IMF could be reduced (iii) towards cloud base due to aggregation and (iv) below cloud base due to sublimation of ice crystals. Note that sublimation could also act as a SIP process as ice crystals could fragment during sublimation (Korolev et al., 2020). The obtained IMF distributions for the three categories (above cloud base, below cloud base, and unknown cloud base) for all clouds combined are presented in Fig. 11. Above cloud base, IMFs >1 (and thus active SIP) were observed nearly 84 % of the time, with a median IMF of around 80 and an interquartile range of 2 orders of magnitude (see Table D1). As expected, IMFs below cloud base are generally found at lower values along with a wider spread than above cloud base. In contrast, the distribution of IMF where no cloud base could be determined is shifted to even higher values than the IMF distribution above cloud base (median IMF higher by approximately a factor of 3). Additionally, the shape resembles the above-cloud-base distribution. These observations suggest that ICNCs in the unknown category were potentially not strongly affected by sublimation and could be classified as above cloud base. Table D1 summarizes the IMF calculations per cloud as well as individual key numbers of the clouds (e.g., concentration of INP and aerosol at cloud top and temperatures at different heights). The individual IMF distributions per cloud and all calculated IMFs combined are presented in Fig. D1.

Figure 12 Frequency histogram of all observed ice multiplication factors (IMFs) above cloud base combined with altitude. Frequencies are normalized to the number of observations (right axis) per height bin. Gray background indicates where no ice crystal number concentration data were observed. The dashed white line indicates an IMF of unity. The median IMF of all data is indicated in solid white.

In Fig. 11, it is conspicuous that the IMF distribution (above cloud base) features a heavy tail below unity signifying higher INP concentration at cloud top than local ICNCs. These IMFs can be explained when looking at the IMFs as a function of height (Fig. 12). IMFs below unity are found in two height ranges (one between 3.5 and 4 $\mathrm{km}\mathrm{a}.\mathrm{s}.\mathrm{l}.$ and one between 4.5 and 6 $\mathrm{km}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). The two regions belong to the early stages of clouds 5 and 6 (see Fig. 10a) where the cloud was still developing. This suggests that regions of IMFs smaller than unity could be a proxy for regions where primary ice generation dominates. After the first hour, median ICNCs of clouds 5 and 6 increased by over 2 orders of magnitude and 1 order of magnitude, respectively. Accordingly, the IMFs also increased at the corresponding heights above unity, potentially resembling a change from a primary ice-dominated regime towards the onset of SIP.

Figure 13 As for Fig. 12, but with the altitude converted to temperature using data from the COSMO-1 analysis.

Despite all uncertainties related to SIP, the effectiveness of the underlying processes varies with temperature (e.g., the Hallett–Mossop process is thought to be active only between −8 and −3 ^∘C). To elucidate the obtained IMF with regard to ambient temperatures, retrieval height was converted to ambient temperature using the COSMO-1 analysis. The temperatures at the top and base of the data frame were averaged over the time of available ICNC data (as listed in Table D1 and indicated in Fig. B1) and interpolated linearly with height. The resulting histogram of IMFs from all cloud events is presented in Fig. 13. Individual histograms per cloud can be found in Fig. E1. SIP was active in all observed clouds (IMFs>1) and at all observed temperatures, highlighting that SIP processes are not only restricted to warm subfreezing temperatures, where the Hallett–Mossop process prevails. High frequencies of IMF above the median IMF were observed between −10 and −5 ^∘C. Hanna et al. (2008) found that for cloud tops at temperatures between −20 and −10 ^∘C frequently precipitation is initiated. Strong SIP leading to high ICNC could promote precipitation initiation given a favorable environment. However, no substantial difference in IMF distributions between precipitating and non-precipitating clouds was found (highest IMFs were even found for a non-precipitating cloud, cloud 7 in Fig. E1). Such a link of the observed IMFs to the ambient environmental conditions is not made since our measurements are restricted to only seven cloud events, with the number of observations biased towards temperatures warmer than −15 ^∘C.

3.4.3 Uncertainties and caveats of the proposed methodology

As stated in the introduction, using the IMF (as defined in Eq. 12) itself does not allow the inference of the underlying SIP process due to the potential (temporal and spatial) displacement of first ice crystal formation and occurrence of SIP. For the same reason, it is also not possible to quantitatively describe the uncertainty of the calculated IMFs. Three main factors add uncertainty to this method:

• (i) uncertainties from the ICNC retrieval, which were shown to be within an order of magnitude (Fig. 8),

• (ii) uncertainties from the INP concentration retrieval at cloud top, which are also within an order of magnitude (Fig. 7b), and

• (iii) uncertainties from assuming a constant INP concentration with time throughout the cloud.

Based on the retrieved INP concentrations for the seven clouds (Table D1) ranging from approximately 10⁻³ to 10⁻¹ StdL⁻¹, we estimate this uncertainty to also be of 1 order of magnitude. The combined influence of the uncertainties on the calculated IMF is not possible to determine due to the missing spatial (and thus logical) link between cloud top INP concentration and in-cloud ICNC following the considerations above. Our method has additional caveats.

• All ICNC data above cloud base are assumed to not be affected by sublimation. However, it is not necessarily true that above the detected cloud base no sublimation occurred as dry air layers could have been entrained.

• In most SIP processes, small ice splinters are created that could then start growing in favorable cloud conditions. This small ice would create a small but sharp peak at small velocities in the Doppler velocity spectrum (see e.g., Fig. 1c in Li et al., 2021). In the ICNC retrieval, the entire Doppler velocity spectrum is fit to a unimodal distribution for an assumed ice crystal shape (Bühl et al., 2019). Consequently, the secondary ice would only considerably affect the retrieved ICNC when the ice splinters have grown to a size when their distribution approaches and ultimately merges with the background ice distribution.

• Excluding ICNC found in sublimation conditions could exclude contributions of the SIP process of sublimation fragmentation (Korolev et al., 2020). However, the produced secondary ice crystals will not fully sublimate if they re-enter the cloud again where they would contribute to the ICNC above cloud base.

• In our approach, a constant cloud top height and thus a constant INP concentration at cloud top were assumed for an entire cloud event. We tested that the effect of individual cloud heights did not substantially affect the obtained IMFs.

• Calculating the ICNC from the remote sensing observations requires a sufficiently strong detectable lidar echo, limiting our method to optically thinner clouds (where the signal is not attenuated), which could induce a bias. Thus, all quantitative findings above need to be treated with caution. The statement that can be conveyed with higher certainty is the observation of SIP being active at all temperatures.

The data of each cloud were selected manually (see red data frames in Fig. B1). A robustness test of the data selection is provided in Appendix F.

3.4.4 Comparison to previous field observations

Figure 14 Comparison of obtained ice multiplication factors (IMFs) in this study (gray shading indicating the 10th and 90th percentiles of observations with the median indicated by vertical lines per temperature bin) to previous field observations. The polygons envelop the observed IMFs per study. Circles indicate observations of one IMF at one temperature. Note that for the study of Crawford et al. (2012), no lower limit in IMF could be determined. Note that the retrieved IMFs from Ladino et al. (2017) tend towards 10⁹ for temperatures towards 0 ^∘C. The dashed black line and the solid black line indicate an IMF of unity and the observed median IMF in this study, respectively.

Atmospheric ICNC and INP concentrations have been studied for more than 50 years, most often utilizing aircraft as a measurement platform. In the beginning of the millennium it was brought to the cloud physics community's attention that previous measurements could have recorded artificially enhanced ice crystal concentration due to shattering of ice crystals on the tips of the cloud probes (see e.g., Field et al., 2006; Korolev et al., 2011). As a response, correction algorithms were developed (e.g., Field et al., 2006; Lawson, 2011), and the probe design was adapted (e.g., Korolev et al., 2011, 2013). Thus, we compare our findings only to field observations of selected studies from the past decade in which the effect of shattering was mitigated. Details about the individual studies and the procedure of obtaining IMFs from the ice crystal and INP data are described in Appendix G.

Figure 14 comprises the range of observed IMFs in this study along with derived IMFs from previous studies which observed maritime convective systems (Stith et al., 2011), continental clouds (Crawford et al., 2012; Taylor et al., 2016), tropical clouds (Lasher-Trapp et al., 2016; Ladino et al., 2017), orographic MPCs (Mignani et al., 2019; Lauber et al., 2021), and Arctic MPCs (Pasquier et al., 2022). The IMFs derived from the previous studies not only exhibit great variability among each other but also span over a considerable value range in each individual study. The ranges of IMF from previous studies generally coincide with the observations from our study, strengthening the applicability and validity of our method. Most studies observed clouds at temperatures warmer than −15 ^∘C, and some acknowledge the Hallett–Mossop process as the most likely mechanism producing secondary ice concentrations around −5 ^∘C (Crawford et al., 2012; Taylor et al., 2016; Ladino et al., 2017). Independent of temperature we observed IMFs between 10¹ and 10³ (see Fig. 11). At temperatures warmer than −15 ^∘C, previous studies consistently confirm this observation and suggest that even IMFs of 10⁴ are frequently observed in the atmosphere. The observations of Stith et al. (2011) and Pasquier et al. (2022) suggest SIP also occurring at temperatures much lower than −15 ^∘C, supporting our findings. Especially at these colder temperatures, the contribution of SIP processes remains uncertain due to the lack of available measurements. At these low temperatures, which are often encountered higher up in the troposphere, our method is a promising tool to investigate IMFs for many clouds and associate ranges of IMFs with regimes of environmental conditions.

4 Conclusions and outlook

In this study we retrieved atmospheric ice-nucleating particle (INP) concentrations in dust-dominated and continental air masses and ice multiplication factors (IMFs) in wintertime orographic mixed-phase clouds (MPCs) using active remote sensing and in situ observations obtained during the RACLETS field campaign in the Swiss Alps in February and March 2019. Using in situ aerosol and INP observations on a mountaintop in proximity to a lidar located in a nearby high valley, we investigated the remote sensing retrieval performance of aerosol properties and different INP parameterizations for dust and continental-aerosol-loaded air masses. In addition, we validated the retrievals of ice crystal number concentration (ICNC) and ice water content (IWC) derived from the combination of radar and lidar data with in situ cloud observations. We proposed a methodology to estimate IMFs from collocated lidar and radar data using a tuned INP parameterization for the continental aerosol. We conclude with the following key findings.

1. The aerosol retrieval of the lidar could be validated. We found that the retrieval of particle number concentration n₂₅₀ (number concentration of particles with physical radii ≥ 250 nm) is less prone to bias than the retrieval of surface area concentration.

2. The dust parameterization by DeMott et al. (2015) in combination with the calibration factor derived by Schrod et al. (2017) was found to predict the INP concentrations from Saharan-dust-carrying air masses best (within a factor of 10). Although the calibration factor was derived for condensation and deposition mode INP concentrations at temperatures $\le -\mathrm{20}$ ^∘C and below, we could extend its applicability for immersion mode INP concentrations at temperatures $\ge -\mathrm{20}$ ^∘C.

3. We found that the INP parameterization by DeMott et al. (2010) (D10) in its original form overpredicts continental INP concentration in the Swiss Alps. We propose a multiplicative calibration factor (cf=0.0204) to significantly improve its performance, making it applicable to wintertime continental aerosol.

4. Retrievals of ICNC and IWC from collocated lidar and radar measurements were found to agree with in situ cloud observations within an order of magnitude.

5. We presented a methodology to estimate IMF from collocated lidar and radar measurements. Secondary ice production (SIP) was active nearly 84 % of the time and at all observed temperatures roughly between −30 and −5 ^∘C. Overall, the median IMF found in wintertime MPCs over the Swiss Alps was around 80, with an interquartile range of 2 orders of magnitude. High frequencies of IMF above 80 were observed between −10 and −5 ^∘C based on the restricted number of observed clouds.

6. The IMFs observed with the proposed method showed substantial overlap to IMFs retrieved from previous field observations, especially at temperatures warmer than −15 ^∘C, where IMFs between 10¹ and 10⁴ were frequently observed.

Our study profited from an ideal setting for the closure of continental INP concentration retrieval from the lidar observations combined with in situ observations. The proposed calibration factor should be validated during future field campaigns or on large datasets from past field campaigns with suitable setup. Future development of INP parameterizations intended for lidar applications should be based on aerosol number concentration (n₂₅₀) due to the reduced bias in the lidar retrieval. Different assumptions are included in the proposed methodology to determine IMF from collocated lidar and radar observations. However, the applicability of the methodology is supported by the consistency to the results from previous field observations. The tuning of the D10 parameterization was essential to reduce errors induced by the INP concentration retrieval. As pointed out by Field et al. (2017), reliable quantification of primary ice formation is crucial to constrain SIP. Without the tuning, the INP concentration would have been overestimated by nearly 2 orders of magnitude, and IMFs would have consequently been underestimated by 2 orders of magnitude. The retrieval of ICNC could be further improved by varying the assumed particle shapes. Ultimately, the efficiency and occurrence of proposed processes that cause SIP vary based on the environmental conditions (Field et al., 2017; Korolev and Leisner, 2020). It is thus important to link the observed IMFs back to the environmental conditions the cloud formed in. The application of our methodology to large datasets of ground-based or space-borne aerosol and cloud observations could help to assess IMFs in different cloud types and constrain the environmental conditions for different SIP regimes. Machine learning tools such as random forests or gradient-boosted decision trees could be valuable tools in constraining environmental regimes impacting SIP. This could not only accelerate our understanding of ice multiplication in the atmosphere but also be of immense help to motivate and plan laboratory experiments that investigate individual ice multiplication processes.

Appendix A: Categorization of continental background and Saharan dust samples

Shortly before the Saharan dust event on 19 February 2019, the Davos region was influenced by air masses already featuring high INP concentrations and a higher aerosol content suggested by an increased signal in the attenuated backscatter (see events 1 and 2 in Fig. A1b). The fact that these air masses do not necessarily belong to the following plumes of Saharan dust is further supported by the peaks observed in relative humidity before the phase of Saharan dust (Fig. A1a).

Figure A1 (a) Relative humidity observed during the comparison period (14–22 February 2019) on 19 February 2019 at WFJ (mountaintop site). The two peaks in relative humidity during the transition time are indicated with 1 and 2. (b) Observed lidar attenuated backscatter coefficient ${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$ at the height of WFJ (2687 $\mathrm{m}\mathrm{a}.\mathrm{s}.\mathrm{l}.$). The two characteristic time periods for continental background aerosol and Saharan dust are indicated. Two perturbations in ${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$ during the transition time are indicated with 1 and 2, respectively.

Appendix B: Overview of individual clouds

The exemplary procedure for obtaining IMFs was described for cloud 7 observed on 12 March 2019 in Fig. 9. Figure B1 shows the overview for all observed clouds. For cloud 4 (Fig. B1c) the last lidar retrieval was more than 1 d before the ICNC data were available. After the ICNC retrieval time the next lidar profile was also more than a day later. However, the profile used (as indicated in Fig. B1c) and the next lidar profile only differed by a factor of 2; thus a strong deviation during the available ICNC data times seems unlikely. For cloud 5 the last lidar profile prior to the event was only available 2 d before. Thus, the closest profile in the morning of the next day was used. Small parts of cloud 5 reach temperatures below −38 ^∘C where potentially homogeneous freezing could occur and additional ice crystals could sediment into higher temperature levels.

Figure B1 As in Fig. 9a for each panel except for the parameter labeling of h_d, t_d, and n_INP,ct. The green shading indicates the period during which the INP concentration at cloud top (dashed blue line) was retrieved. Calculated and retrieved parameters for each cloud are presented in Table D1. The blue dotted line in (d) indicates the height above which homogeneous freezing is possible ($T<-\mathrm{38}$ ^∘C).

Appendix C: Cloud regime classification

The exemplary procedure for obtaining IMFs was described for cloud 7 observed on 12 March 2019 in Fig. 9. Figure C1 shows the classification for all observed clouds. Cloud 3 (Fig. C1b) was entirely classified as unknown due to the liquid water layer at around 2500 m.

Figure C1 Time slicing and classification of the ice crystal number concentration data indicated in Fig. B1 (red rectangle encompasses rectangles presented in this plot). White background represents regions where the lidar either was fully attenuated or underwent automatic maintenance. The attenuated backscatter coefficient at 1064 nm (${\mathit{\beta }}_{\mathrm{1064}\mathrm{nm}}^{\mathrm{att}.}$) is presented in the background. A (liquid) cloud base estimate derived from the attenuated backscatter coefficient at 532 nm (${\mathit{\beta }}_{\mathrm{532}\mathrm{nm}}^{\mathrm{att}.}$) is indicated by red circles. The automated liquid water classification by the Cloudnet algorithm is presented by lime tiles. Empty rectangles in cyan, yellow, and magenta indicate regions classified as above cloud base, below cloud base, and unknown cloud base, respectively. Cloud top derived from radar observations (Fig. B1) is indicated by the dashed blue line. The dotted light blue line in (d) indicates the height above which homogeneous freezing is possible ($T<-\mathrm{38}$ ^∘C).

Appendix D: IMF summary and additional cloud information

Table D1 provides a summary of obtained IMF for all ICNC data (red rectangles in Fig. B1) and using only ICNC data above cloud base (cyan rectangles in Fig. C1) to demonstrate the effect of the selection. Generally, the selection does not change the median IMF (cloud 1) or increases it (clouds 2, 4, and 7). In the cases of clouds 5 and 6 the median IMF decreased due to the selection. This is explained by the exclusion of high ICNCs that were present during strong precipitation when the lidar was attenuated close to the ground and no cloud base could be determined. Figure D1 shows IMF histograms per cloud and for all clouds combined using all ICNC data.

Table D1 Overview of the secondary ice production (SIP) analysis per investigated cloud and environmental information for each cloud. For each cloud the total number of ice crystal number concentration (ICNC) observations (N_obs.), and the median ICNC (${\stackrel{\mathrm{̃}}{n}}_{\mathrm{ice}}$), the times when SIP was active (SIP active, i.e., percentage of IMFs>1), and the median IMF ($\stackrel{\mathrm{̃}}{\mathrm{IMF}}$) and its interquartile range (IQR_IMF, in orders of magnitude) are presented for all ICNC data available and ICNC data above a detectable cloud base. Further, the concentration of particles with radii ≥ 250 nm and INPs at cloud top before the event (n_250,ct, n_INP,ct), the cloud height (h_ct), the height (h_d) and time (t_d) over which ICNC was obtained, and the temperatures at cloud top (T_ct) and top (T_dt) and bottom (T_db) of available ICNC data are given. In the last column, the total number of observations (N_obs.) and averages (weighted by number of observations) of the variables are given combined for all clouds.

Figure D1 Probability density functions (PDFs) of obtained ice multiplication factors (IMFs, Eq. 12) for each of the seven observed clouds and all IMFs combined (thick solid red) using all ice crystal number concentration data (not only above cloud base). The dashed black line indicates unity. The median IMF and the percentage of IMF larger than unity for the combined histogram are given in the legend.

Appendix E: IMF temperature histogram per cloud

IMF histograms with temperature for each cloud for ICNC data above cloud base are presented in Fig. E1. Note that for cloud 3 no cloud base could be determined (see Fig. C1b).

Figure E1 Frequency histogram of observed ice multiplication factors (IMFs) above cloud base with temperature per cloud. Frequencies are normalized to the number of observations (right axis) per temperature bin. Radar height was converted to temperature using data from the COSMO-1 analysis. Gray background indicates where no data were observed. The dashed white line indicates unity. The median IMF of all data is indicated in solid white.

Appendix F: Sensitivity of obtained IMF based on data selection

To assess the robustness of our methodology we performed a sensitivity study by changing the upper and lower heights of the selected ICNC by 200 m and the starting and ending times by 20 % of the entire time span towards the center (see h_d and t_d in Fig. 9a). Note that this test was done on the entire ICNC data available, i.e., the data frames indicated in red in Fig. B1. As seen in Table F1 the obtained SIP parameters did not drastically change upon changes in the data selection. The fraction of time when SIP was active did not change by more than five percentage points (highest upon changing the starting time by 20 %, t_s+20 %), and the IMF-related parameters did not change by more than a factor of 1.5 (highest upon changing the starting time by 20 %, t_s+20 %).

Table F1 Benchmarks of sensitivity studies of obtained secondary ice production (SIP) and ice multiplication factor (IMF) based on data selection for both estimates. The data selection limits in height and time were changed by ±200 m and by ±20 % of the time span, respectively.

Appendix G: Extraction of IMFs from previous publications

Airborne ice crystal measurements could have overestimated the ICNC due to shattering of ice crystals on the cloud probes (see e.g., Field et al., 2006; Korolev et al., 2011). The retrieval of ICNC from radar observations has also become more reliable with the further technical development of the instruments. Hence, we considered publications from the last decade to compare our obtained IMF to them (Table G1). However, most studies provide information of ambient temperature, INP concentration, and ICNC separately. Additionally, INP concentrations are typically given at the same location or temperature as the observed ICNC – in contrast to our approach. In the following we briefly describe how IMFs were retrieved from each study. From Stith et al. (2011) Fig. 8, maximum and minimum IMF were obtained for the two legs 05:22 and 05:42 UTC. The ambient temperature ranges during the legs were taken from Stith et al. (2011) Fig. 5. Crawford et al. (2012) estimate an upper limit of 0.01 L⁻¹ in INP concentration based on source considerations. From Crawford et al. (2012) Fig. 3, maximal ICNC with temperature were taken. A minimal ICNC could not be distinguished from zero such that IMFs smaller than 1 are possible, but not further determinable. Lasher-Trapp et al. (2016) Fig. 5 presents the excess of ICNC over INP concentration combined as a function of temperature. We considered ICNC data that were not prone to possible side contamination as the authors write. Taylor et al. (2016) Fig. 13 provides ICNC with estimates on INP concentration with temperature. Note that we consider Run 11.1 as a point measurement. Ladino et al. (2017) Fig. 4 encompasses an extensive overview of their ICNC measurements and predicted INP concentration with temperature. IMFs were calculated using the envelope of all ICNC observations and the INP concentration at corresponding altitude. Mignani et al. (2019) provide IMF between 4 and 25 (Fig. 3) at temperatures between −17 and −12 ^∘C (Fig. 2). Lauber et al. (2021) Fig. 6 provides a median ICNC concentration close to the melting layer. INP concentrations below the detection at this temperature are estimated to be smaller than the lower detection limit of the drop freezing instrument. For measurements close to the melting layer, a temperature of −0.5 ^∘C was used for plotting. IMFs from Pasquier et al. (2022) were provided by personal communication. Li et al. (2021) do not provide correlated data of ICNC, INP concentration, and temperature. They state that “ice number concentrations tend to be 1–3 orders of magnitude higher than expected INP concentrations for clouds at temperatures of  ^∘C or warmer”.

Stith et al. (2011)Crawford et al. (2012)Lasher-Trapp et al. (2016)Taylor et al. (2016)Ladino et al. (2017)Mignani et al. (2019)Lauber et al. (2021)Li et al. (2021)Pasquier et al. (2022)

Table G1 Overview of the studies considered in the ice multiplication factor comparison and the corresponding cloud type sampled in the study.