Retrieval of microphysical dust particle properties from SALTRACE lidar observations: Case studies

Saharan dust is a major natural atmospheric aerosol component with significant impact on the Earth radiation budget. In this work we determine the microphysical properties of dust particles after a long-range transport over the Atlantic Ocean, using input from three depolarization channels of a multi-wavelength polarization Raman lidar. The measurements were performed at Barbados in the framework of the Saharan Aerosol Long-Range Transport and Aerosol-Cloud-Interaction Experiment (SALTRACE) in the summers of 2013 and 2014. The microphysical retrievals are performed with the software tool 5 SphInX (Spheroidal Inversion Experiments) which uses regularization for the inversion process and a new two-dimensional (2-D) extension of the Mie model approximating dust with spheroids. The method allows us to simultaneously retrieve shapeand size-dependent particle distributions. Because dust particles are mostly non-spherical this software tool fills the gap in estimating the non-spherical particle fraction. Two cases measured on 10 July 2013 and 20 June 2014 are discussed. 2-D radius-bimodal shape-size distribution are retrieved. The ratio of spherical-to-non-spherical contributions to the particle num10 ber concentration was found to be about 3/7. A volume-weighted effective aspect ratio of 1.1 was obtained, indicating slightly prolate particles. The total effective radius for the two cases in the preselected radius range from 0.01-2.2 μm was found to be, on average, 0.75 μm. The stronger dust event (10 July 2013) showed about 24% higher values for the total surface-area and volume concentration. Finally, we compare our results with the ones from the polarization lidar-photometer networking (POLIPHON) method and ground-based photometers as well as with airborne in situ particle counters. Considering all differ15 ences in these independent approaches, we find a qualitatively good agreement between the different results and a consistent description of the dust cases. Such an extensive comparison is a novel and fruitful exercise and corroborates that the mathematical retrieval based on Raman lidar data of particle backscattering, extinction, and depolarization is a powerful tool even in the case of dust particles.

1 Introduction 20 In recent years, mineral dust has received growing interest in atmospheric and climate research because of its significant impact on cloud formation and the radiation budget of the Earth. It contributes to about 50% of the global annual particle emission by optical properties to derive the dust particle size distribution and microphysical properties of mineral dust after long-range transport from the Saharan desert across the Atlantic Ocean to the Caribbean island of Barbados. Part of the retrieval is the separation of fine and coarse mode dust. To solve the so-called ill-posed problem we use the inversion technique with regularization. Such methods are a widely used and permit an efficient retrieval of the size distribution for spherical particles (Müller et al., 1999;Böckmann , 2001;Veselovskii et al. , 2004;Müller and Böckmann et al. , 2016). In Böckmann and Wauer 25 (2001), algorithms for spherically multilayered and absorbing particles are presented. Simulations show that the influence of particle shape properties and chemical composition is not negligible and need to be carefully considered when inversion algorithms are used to determine the particle microphysical properties. Lopatin et al. (2013) present the GARRLiC (Generalized Aerosol Retrieval from Radiometer and Lidar Combined data) algorithm which assumes height-independent microphysical particle properties and permits retrievals of the vertical profiles 30 of both fine and coarse mode aerosol concentrations as well as the size distribution and the complex refractive index for each mode.
In most attempts, the Mie scattering model is used. This model offers the computation of the optical properties of fine mode aerosol (Böckmann et al. , 2005) such as biomass burning smoke and urban haze with good accuracy (Osterloh et al. , 2013;Samaras et al. , 2015a;Ortiz-Amezcua et al. , 2017;Ritter et al. , 2018). However, the Mie model no longer provides a viable 35 2 https://doi.org /10.5194/acp-2020-459 Preprint. Discussion started: 30 June 2020 c Author(s) 2020. CC BY 4.0 License. description of optical properties in the case of mineral dust particles since dust particles are clearly non-spherical in shape. The consideration of non-sphericity, even in its simplest form by assuming a spheroidal shape for the dust particles, adds multiple layers of complexity to the retrieval procedure. Numerical stability becomes an important issue as well.
Nonetheless, particle shape is known to have substantial impact on the scattering characteristics, especially with respect to backward scattering (Mishchenko et al. , 1996). The assumption of spheroidal particle geometry seems to reproduce the 5 optical properties of dust particles significantly better than the use of spherical shape model (Mishchenko et al. , 1996;Kahnert and Kylling , 2004). Ensemble averages based on the T-matrix method for spheroids agree better with measurements than computations based on Mie theory (Bohren and Singham , 1991). Simulations with a mixture of simple spheroids are found to agree well with measured optical properties of non-spherical aerosol particles. The exact knowledge of the shape characteristics in the modeling efforts is obviously not needed in many cases of measurements with contributions from individual particles 10 with different orientations, shapes, and composition (Wiscombe and Mugnai , 1986;Bohren and Singham , 1991;Mishchenko et al. , 1997;Kahnert and Kylling , 2004;Veselovskii et al. , 2010). Several non-spherical-particle approximations are present in the literature which focuses on optical parameter investigations. In contrast, non-sphericity is usually considered in the microphysic retrieval branch by using the spheroid shape model (Dubovik et al. , 2002;Veselovskii et al. , 2010;Osterloh , 2011;Böckmann and Osterloh , 2014;Samaras et al. , 2015b;Samaras , 2017). 15 The microphysical parameter database hosted by the Aerosol Robotic Network (AERONET) is such an example. Aerosol particle size distributions are assumed to consist of an ensemble of polydisperse, homogeneous spheres and a mixture of polydisperse, randomly-oriented homogeneous spheroids with size-independent aspect-ratio distribution (Giles et al. , 2019).
A customary assumption in non-spherical microphysical retrieval codes, adopted also by AERONET, is that shape and size are independent parameters with the consequence that particle shape is given by a separate aspect ratio distribution. This approach 20 was used by Dubovik et al. (2006) to invert AERONET sun-photometer data into microphysical particle properties. Dubovik et al. (2006) show that mixtures of spheroids can well explain measured spectral and angular dependencies of scattering intensity and polarization. Furthermore, the fraction on non-sphercial particles for aerosol mixtures with a dominant coarse mode (particle radii ≥∼ 1µm) can be determined. The retrieval results indicate that nonspherical particles with aspect ratios ∼ 1.5 and higher dominate in desert dust plumes. In Müller et al. (2013) this concept was then also used in lidar inversion 25 procedures. Bi et al. (2018) present the invariant imbedding T-matrix method to assess the backscattering linear depolarization ratio (LDR) of non-spherical particles in a super-ellipsoidal shape space. Super-ellipsoids have inherent flexibility regarding modelling of particle aspect ratio, roundness, and concavity, which are the salient characteristics of, e.g., sea salt and dust aerosols.
The results provide comprehensive insight in the relationship between LDR of atmospheric aerosols and particle shape and 30 refractive index characteristics. This can be used in the interpretation of lidar observations. Gasteiger et al. (2011) modelled lidar-relevant optical properties of mineral dust aerosol and compared the simulated results with optical properties derived from lidar data during the Saharan Mineral Dust Experiment (SAMUM) and found that irregularly shaped dust particles with typical refractive index characteristics cause a higher linear depolarization ratio than corresponding spheroids. Simulations with a variety of irregular-shape models improved the agreement with lidar observations. 35 In this work we are using the new concept of a two-dimensional (2-D) spheroidal particle model implicit in the software tool SphInX (Spheroidal Inversion eXperiments) developed by Samaras (2017). In contrast to Dubovik et al. (2006) and Müller et al. (2013), no restriction regarding shape-size separability is given. It is left to the retrieval to reveal the shapesize relationship. In fact, it appears that there is a link between the aspect ratio distributions and particle size as it was found in independent studies during the SAMUM campaign Koepke et al. , 2015). The SphInX tool was 5 already successfully used in Soupiona et al. (2019) to investigate aerosol optical and microphysical properties during selected Saharan dust events over Athens (Greece; NE Mediterranean) and Granada (Spain; NW Mediterranean) focusing on short and long-range dust processes.
SphInX is able to include depolarization ratio information at any wavelength. Depolarization-ratio measurements with advanced polarization lidars at four wavelengths (355, 532, 710 and 1064 nm)  were performed during 10 SAMUM. The mean linear particle depolarization ratio at 532 nm was 0.31, in the range between 0.27 and 0.35 with a slight decrease of the linear particle depolarization ratio with wavelength between 532 and 1064 nm from 0.31 ± 0.03 to 0.27 ± 0.04.
Other statistical properties of Saharan dust are provided in Tesche et al. (2009) during SAMUM. The 500 nm dust optical depth ranged from 0.2-0.8 at the field site south of the High Atlas mountains, Ångström exponents derived from photometer and lidar data were between 0-0.4. On average, extinction-to-backscatter ratios of 53-55 sr (±7-13 sr) were obtained at 355, 15 532 and 1064 nm. The question arises which depolarization channel is the most valuable for the microphysical retrieval. It is not surprising that the answer depends on the aerosol type. Gasteiger and Freudenthaler (2014) found for transported Eyjafjallajökull volcano ash, e.g., that additional 1064 nm depolarization measurements significantly reduce the uncertainty in the retrieved mass concentration and effective particle size. This significant improvement in accuracy is the result of the increased sensitivity of the measured optical data to the impact of larger particles. 20 Here we investigate two measurement cases from SALTRACE Barbados campaign. The measurements were conducted on 10 July 2013 (case 2) between 23:15-0:45 UTC and 20 June 2014 (case 1) between 23:10-02:10 UTC. The used lidar of the Leibniz Institute for Tropospheric Research (TROPOS), shortly called BERTHA (Backscatter Extinction lidar Ratio Temperature Humidity profiling Apparatus), provides depolarization ratio profiles at 355, 532 and 1064 nm (Haarig et al. , 2017a) in addition to three backscatter and two extinction coefficient profiles, the so-called 3β + 2α data set. These measurements were 25 basically performed to characterize optical and microphysical properties of aged dust plumes after long-range transport across the Atlantic Ocean. Both cases were already investigated by Mamouri and Ansmann (2017); Haarig et al. (2017a). Furthermore, case 1 was also investigated in Tesche et al. (2019) with focus on the retrieval of microphysical properties and the guiding question to what extent uncertainties can be reduced or avoided by including specific depolarization-ratio information in the data analysis. In contrast to our approach, Tesche et al. (2019) used the model of Dubovik et al. (2006) with 30 a fixed axis-ratio distribution. Tesche et al. (2019) found that the consideration of light-depolarization information strongly improves the quality of the retrieval products.
Case 2 provides a unique opportunity to compare our lidar inversion results (profiles) with airborne in situ particle counter measurements. In addition, we compare our results with respective ones obtained by means of the polarization lidar-photometer networking (POLIPHON) method (Mamouri and Ansmann , 2017). Additional comparisons are made with inversion products from AERONET sun-photometer observations. The structure of this paper will be as follows. In Section 2, we provide an overview of the 2-D model methodology and the used instruments. Section 3 presents the results of the two case studies. The final Section 4 summarizes our findings and contains concluding remarks.

2-D Model and Algorithm
The extinction-(α) and backscatter coefficients (β), both denoted with Z(λ) with wavelength λ are commonly defined as the collective interaction (scattering/absorption) probability over a particle number size distribution (PSD) n(r). Denoting with C α/β (r, λ; m) the associated (interaction) cross sections, we can express this relation through a Fredholm integral of the first where m is the complex refractive index (RI). A Fredholm integral equation of first kind is a compact operator and therefore the inversion of it is an ill-posed problem. We refer to Eq. 1 shortly as the Mie model, since the optical efficiencies are calculated by the (Lorenz-) Mie theory. 15 We extend this model considering a spheroid-particle approximation, which needs only one additional parameter, i.e. the aspect ratio η, to be fully defined. It is formalized as follows: where [r min , r max ] and [η min , η max ] are sensible and sufficiently wide radius-and aspect ratio ranges determined experimentally as part of the initial parameters of a microphysical retrieval as discussed in Samaras et al. (2015b). Thus, the integral over 20 size, i.e. the radius r, is integrated again over a parameter representing the different aspect ratios η of spheroids, i.e. oblate-(η < 1), sphericalη = 1, and prolate (η > 1) shapes. This means that we are now looking for 2-D particle number shape-size distributions. Note that there is no particular shape-size separation relation for the shape-size distribution, as for instance in Veselovskii et al. (2010). The radius here makes sense as the one of a sphere with equal volume to a spheroid. The latter is found by V = 4πa 2 b/3, so that r = 3 √ a 2 b, where a and b are the semi-minor and semi-major axis of the spheroid, respectively.

25
Regarding the PSD, we see that the case of spheroids defines a two-dimensional (2-D) generalization of the spherical model, but the optical parameters Z still remain only wavelength-dependent, which is why we shall consider this, from a mathematical point of view, a quasi-2-D model. Particle orientation comes about naturally when addressing non-sphericity but is generally unknown. For this, we further assume randomly oriented spheroidal particles for the calculations of the scattering cross sections in Eq. 2 providing further simplifications to the model. Replacing the scattering cross sections with dimensionless efficiencies Q α/β and the number-with volume distribution we can reformulate the generalized model Eq. 2 to the one we use in practice for our application where A is the particle surface area, where we used the fact that in a convex particle ensemble the average area per particle is 5 equal to A/4 (van de Hulst , 1958).
A precalculated database will be used, created by the software tool Mieschka (Rother and Kahnert , 2014). Mieschka is able to perform scattering-related T-matrix calculations for spherical particles and rotationally symmetric non-spherical particles with a convergence strategy based on Barber and Hill (1990); Wiscombe and Mugnai (1986). Additionally, it provides an extensive database of scattering quantities for spheroidal geometries. Mieschka's look-up tables include scattering efficiencies 10 for a 6 × 7, ( (m) × (m)) refractive index grid, a total of 42 RI values, 7 different aspect ratios and a size parameter range [0.02, 40 µm] with a resolution of 0.2 µm, see Table 1. While the maximum size parameter is reasonably large, its potential cannot be fully exploited here in terms of the radius extent since the lidar wavelength at 355 nm restricts the maximum radius (r max ) to about 2.2µm. All working formulas are with respect to r. The resolution gap in the aspect ratio needed for the integration is handled by interpolation to the nearest neighbor. 15 The microphysical parameters can be redefined in 2-D analogously. Similarly to the 1-D case we set the number shape-size distribution n(r, η) = 3 4πr 3 v(r, η), where v(r, η) is the volume shape-size distribution. We define the number concentration the total volume concentration 20 and the surface-area concentration The effective radius is calculated by r eff = 3 vt at , using Eq. 5 and 6. The function G(r, η) denotes here the spheroidal geometrical cross section of the particle, which can be explicitly computed as follows 2π a 2 + ab e sinh −1 (e) , where e = 1 − a 2 /b 2 if η > 1.
(7) 25 Furthermore some new parameters are introduced (Samaras , 2017) in order to study the shape of the size distribution in more detail. We define the volume-weighted-effective aspect ratio and the aspect ratio width The parameters η eff and η var give us an estimation of a central tendency of the aspect ratio and the spread of the values from this type of the mean.
The following parameter specializes in identifying individual (non-)spherical contributions (%) in the volume concentration.
We define the non-spherical volume concentration fraction by 10 where χ is a positive small, yet not negligible number called the shape-transition parameter. Similarly we call 1−γ the spherical volume concentration fraction (svcf).
Since there is no equivalent in the literature of a two-dimensional particle distribution, we introduce the reduced volume size distribution v η (r), defined as the volume shape-size distribution v(r, η) integrated over the aspect-ratio domain, i.e. 15 where r is the volume-equivalent particle radius. The function v η (r) is able to provide the collective trend of all contributing particle geometries in the particle distribution. However, we note that this limited (in terms of information as compared to v(r, η)) particle distribution is not directly comparable with the usual size distribution used in literature, but it can be used in order to have a general sense qualitatively, see Fig. 1 for an overview.
Furthermore, in our analysis we use data derived by inversions of sun-photometer measurements provided by the AERONET 20 database (Holben et al. , 1998). AERONET offers a series of inversion products including the effective radius, the volume concentration, the complex refractive index, the single scattering albedo and the aerosol optical depth, which can be available in two quality levels, namely 1.5 and 2.0, for cloud screened and quality assured data respectively. Sun-photometers are passive remote sensing instruments with different operation principles from lidars and use different inversion techniques based on the theory of optimal estimation (Dubovik and King , 2000;Dubovik et al. , 2000). Although we will not expand on the latter, 25 we will point out some incompatibilities with our approach, which are essential to recognize for our subsequent microphysical analysis. The first and most important difference is that AERONET retrievals relate to the whole atmospheric column while lidar data target specific layers and therefore the comparison cannot be quantitative. Moreover, the volume concentration size distribution of the particle volume u(r) retrieved by AERONET is defined as the derivative du/d ln r with the associated total volume concentration u t , both given in µm 3 µm −2 . Following the LM-model, the volume size distribution (for lidars) is found by v(r) = dv t /dr and measured in µm 3 µm −1 cm −3 = µm 3 µm −2 m −2 .
In order to make sense of these two different measures and have some kind of comparison we turn to the quantityṽ(ln r) = dv t /d ln r = rv(r), which indicates that the difference in units with du/d ln r lies within a multiple of the meter. Here we apply 5 the same concept also for our generalized model just by replacing v(r, η) by v η (r). In practice, we often multiplyṽ(ln r) with the aerosol layer thickness, usually several kilometers, (Osterloh et al. , 2013;Samaras et al. , 2015a). This is the motivation behind the units given in the explaining flowchart, see Fig. 1, for the so-called volume concentration size distribution (VCSD), following the terminology from Samaras et al. (2015a). Converting a reduced size distribution v η (r) to a VCSD is even more involved, but we use these functions as a bridge (in terms of units) between lidar and sun-photometer. AERONET's VCSD is 10 retrieved in 22 logarithmically equidistant bins between 0.05 and 15µm. For the shape-size distribution we use 30 × 30 (r × a) grid points with the radius range [0.01, 2.2] (µm) and the aspect ratio range [0.67, 1.5]. For clarity we note that whenever we show here an one-dimensional size distribution (or VCSD) associated with the lidar, it is always implied that there was originally a shape-size distribution v(r, η) which was first converted to it.
Second, many of AERONET's inversion products are given as a function of wavelength at least at the spectral bands of 442, 15 675, 872, 1020 nm where the four almucantar scans take place, in this category fall also the parameters aerosol optical depth and the refractive index (RI). Especially for the RI we will consider a spectral average of these values in order to compare with the respective parameters from our retrievals. The ranges of the real (RRI) and imaginary part (IRI) of the refractive refractive index, 1.33 ≤ RRI ≤ 1.6 and 0.0005 ≤ IRI ≤ 0.5 form the predefined grid necessary for AERONET's inversion.
Following AERONET's mode distinction in the inversion products we calculate in addition to the total effective radius, the 20 one for the fine and coarse mode. This is done in SphInX by looking for the minimum of the shape-size distribution between r = 0.5 and 1 µm and setting it as the higher and the lower integration boundary for the fine and coarse mode, respectively. In contrast, here, we used a fixed separation limit r = 0.5µm since later we are going to compare our results with the POLIPHON method and in situ particle counter devises where essentially r = 0.5µm is used.
The inversion procedure does not offer a unique solution. Regularization reduces the solution space, keeping only a small 25 number of solutions which minimize the errors. The retrieval products are compared in the following error-related quantities (ERQ), all of which are calculated as a percentage (%) using the standard deviation.
1. Variability of the solution space (Var).

Randomness uncertainty (Unc).
Definitions and explanations are given in the Appendix. 30 We extend the traditional spherical lidar setup 3β + 2α which consists of 3 backscatter coefficients at 355, 532 and 1064 nm and 2 extinction coefficients at 355 and 532 nm, to 3β + 2α + (n)β ⊥ , where β and β ⊥ denote the horizontally and vertically polarized backscatter coefficients, respectively, and n is either 1, 2 or 3 depending on the availability of Raman channels of the lidar system in 355, 532 and 1064 nm. The two cases in this paper pertain to a setup 3β + 2α + 3β ⊥ for retrieving a shape-and size-dependent particle distribution for the first time to the best of our knowledge. Alternatively, instead of β ⊥ , the particle linear depolarization ratio δ = β ⊥ /β may appear in the setup notation to imply the same thing, and further α, β and δ may be omitted (3 + 2 + n). We note that in the case of spherical approximation, the backscatter coefficient (β ) in the aforementioned 6-, 7-or 8-point datasets coincides with the total backscatter coefficient β.
After discretizing of the model Eq. 3 we solve the resulting linear system with regularization, which is the first step to 5 counteract the ill-posedness of this inverse problem. SphInX is able to use different regularization techniques and parameter choice rules. The following combinations are available: -Truncated singular value decomposition (TSVD) with the discrepancy principle, -Tikhonov regularization with the discrepancy principle, with the L-curve method or with the generalized cross validation, -Padé iteration with the discrepancy principle or with the L-curve method. 10 Details on the widely used methods TSVD and Tikhonov and the parameter choice rules can be found in most books about regularization, e.g., Hansen (2010). Padé iteration is used by the software in the context of the so-called generalized Runge-Kutta regularization methods (Böckmann and Kirsche , 2006;Böckmann and Osterloh , 2014). For our two case studies in the next section we use Padé iteration since preliminary numerical tests revealed an overall superior behavior of this iterative regularization as compared to the other built-in methods, a fact also confirmed by further massive theoretical simulations. 15 Extensive details on the latter as well as full functional demonstration of the software tool SphInX (a graphical-user-interface program), designed by the first author, can be found in Samaras (2017).

Leipzig triple-wavelength polarization Raman lidar and AERONET Sun-Photometer
We give a brief description of the triple-wavelength lidar system which provided the optical parameters. The ground-based  (Haarig et al. , 2017b(Haarig et al. , , 2019. The BERTHA lidar system is a container-based, 25 multi-wavelength polarization Raman lidar. It has been upgraded since the SAMUM campaigns (Althausen et al. , 2000;Tesche et al. , 2011) to enable the measurement of the depolarization ratio at three wavelengths (355, 532 and 1064 nm) simultaneously, which is crucial to characterize the dust after long-range transport. A more detailed description of the 13-channel lidar system and the polarization characteristics can be found in Haarig et al. (2017a).
Currently it operates as a 3+2+3 lidar system (3 backscatter coefficients, 2 extinction coefficients and 3 depolarization ratios) 30 with an additional water vapor channel (407 nm) and a high-spectral-resolution channel at 532 nm. In recent times, it has been used in a 3+3+2 configuration to determine the extinction at 1064 nm (Haarig et al. , 2016b). The signals are detected with a range resolution of 7.5 m and a time resolution of 10 s. Moreover, the lidar data are smoothed with a sliding window (742.5 m) for the backscatter-and depolarization-profile and as well with 742.5 m for the calculation of the extinction-profile (range for the linear regression). Note that this suppresses the noise in the lidar data, but smooths out also the fine structure of the layers.
Additionally, an AErosol RObotic NETwork (AERONET) sun-photometer (see AERONET web page https://aeronet.gsfc. nasa.gov/, Barbados_SALTRACE site), and a Vaisala radiosonde station (RS92 for profiling of pressure, temperature, RH, and 5 the vector of the horizontal wind component) were operated at the field site.

Airborne Measurements
During the first SALTRACE period in June-July 2013, the ground-based measurements were complemented by airborne in situ measurements. The platform for these in situ measurements was the Falcon research aircraft of the German Aerospace Center (DLR), which was for this purpose equipped with an extensive aerosol instrumentation. An overview of the complete 10 instrumentation and related data products is given in Weinzierl et al. (2017) (Supplement). Information on size-resolved particle number concentrations between about 10 nm and 50 µm are acquired by a combination of wing-mounted optical particle spectrometers (OPS) and condensation particle counters (CPC) and OPS operated inside the aircraft cabin behind an isokinetic aerosol inlet. Lognormal number size distributions (in terms of volume-equivalent particle diameter) are fitted to these data using a Bayesian inversion method . Within this method, uncertainties in each instrument's 15 response and in the particle optical properties are propagated to gain realistic size distribution uncertainties. For comparison with the results of the retrieval described in Section 2.1, the particle number size distributions derived for the inside of the dust layer are first converted into particle volume size distributions, using πd 3 p /6 for the volume of each (volume-equivalent) particle diameter d p . The volume size distributions are then integrated for the sub-and supermicron diameter range to obtain the total "fine" and "coarse" particle volume concentrations. In doing so, the initial number size distribution uncertainties are 20 fully propagated.

SALTRACE campaign case study 20 June 2014
The optical properties are thoroughly investigated in Haarig et al. (2017a). Here we give a brief overview, necessary for our analysis. Fig. 2(a) shows the color plot of the range-corrected signal at 1064 nm (cross-polarized channel) for 20 June 2014 25 (case 1). We can distinguish a well defined intense aerosol layer up to 4.3 km. Later on, we will concentrate on two specific layers, namely a mixed layer containing dust and marine aerosol between 1.5 and 2.8 km and a pure dust layer between 3.3 and 4.0 km. The full overlap between the laser beam and the receiver field of view is reached at 1 km. This is the reason that several lidar profiles are only available at heights above about 1 km (more details ar given by Haarig et al. (2017a)).
Looking at the optical profiles in Fig. 2(b)-(e), the particle linear depolarization ratio at all available wavelengths is well 30 above 20% for heights >1.5 km and reaches up to 32% (δ532) at about 3.7 km, a clear indication for the presence of non-spherical particles. Furthermore, the particle linear depolarization ratio at all available wavelengths in Fig. 2(b) shows a steady increase in the region 0.5-1.5 km from 8% to 20%. This is a clear signature of turbulence and vertical mixing of marine and dust aerosol. These meteorological processes separate the mixing layer from the dust layer according to the studies of Groß et al. (2015). We will mainly focus on the dust layer above 1.5 km height.
A reasonable layer selection is a trade-off between slowly varying lidar ratio and Ångström exponent (AE) and intense 5 backscattering and extinction by dust particles. Focusing on the extinction-to-backscatter ratio profiles at 355 and 532 nm (not shown here, see Haarig et al. (2017a)) in the altitude range 1.5-4.3 km we have values of about 50±9 sr, typical of dust presence (Ansmann et al. , 2003;Groß et al. , 2015;Heese et al. , 2009). The values of the extinction-related Ångström exponent (AE) AE α 355/532 and of the backscatter-related Ångström exponent AE β 532/1064, see Table 2, are also characteristic of Saharan dust particles Müller et al. (2007). Moreover, the depolarization ratio values and the lidar ratios at 355 and 532 nm fall within the ranges of aged Saharan dust found during SALTRACE campaign in 2013 (Groß et al. , 2015;Haarig et al. , 2017a). In the following, we will focus on the layers from 1.5-2.8 km and 3.3-4.0 km height. The average values and standard deviations of the 3β + 2α + 3β ⊥ optical and intensive properties for both layers are tabulated in Table 2.

Shape-Size Distribution and Microphysical Properties
The inversion of the optical into microphysical properties was performed with the software tool SphInX. We used Padé reg-20 ularization with 100 iterations, 9-14 spline points combined with a spline degree within the range 2-5; more details on initial parameter determination can be found in Samaras et al. (2015a). Less spline points, e.g., six were ruled out because they did not behave well with the strong tendency towards radius-bimodality, which was indicated as most probable by preliminary tests. The refractive index grid (RIG) was fixed to RRI × IRI = [1.4, 1.5, 1.6, 1.7, 1.8] × [0, 0.001, 0.005, 0.01, 0.05]. The upper integration boundary of the Fredholm equation was set to the maximum available r max = 2.2 µm as mentioned before. 25 The retrieved volume shape-size distributions are shown in 3-D in Fig. 3 and in 2-D in Fig. 4. Focusing on the altitude range 1.5-2.8 km in Fig. 4(a), we see two very well separated modes, namely a narrow fine mode and a much broader coarse mode, with maxima of 58 and 52 µm 3 µm −1 cm −3 at a radius of about 0.43 and 1.45 µm, respectively. There are contributions of all kinds of spheroidal particle geometries including spherical ones, and for higher aspect-ratios (η ≥ 1.20) the peak difference fades, so that the maxima equalize at about 30 µm 3 µm −1 cm −3 . The 2-D plot in Fig. 4(a) shows a shift to higher aspect ratios 30 for the coarse mode, i.e., to the prolate shape. Turning now to the 2-D plot for the layer 3.3-4.0 km in Fig. 4(b), we see that while the shape of the first mode is similar to the one in the lower layer, the second mode differs significantly. Large prolate particles with high aspect ratios disappear and the volume magnitude drops.
The retrieved microphysical parameters a t , v t , r eff (total, fine, coarse) η eff , η var , svcf, and RI for both layers with a step of 100 m are shown in Table 3. In this table, the variability (Var %) pertains to 5 best solutions for a single dataset corresponding to a specific altitude range, and the mean variability when more than one data sets are involved. The uncertainty (Unc %) of the mean parameter value, found for every altitude range, is also given, corresponding to different altitude ranges. Table 3 shows that our method provides very good stability. The uncertainty of all parameters is very low and thus we could consider any 5 of these height ranges (100 m) to arrive to quite the same result. This could be also the result of good data quality and the smoothing procedure mentioned before. It should be noted that the inversion is further aided by using the maximum number of depolarization data available to date (dataset 3β + 2α + 3β ⊥ ), the benefit of which was also demonstrated in Samaras (2017) as part of the investigations for an ideal lidar setup.
Our algorithm predicts large particles. This is expected for an intense dust event. An effective radius of 0.73 µm (layer: 10 1.5-2.8 km) for the entire size distribution is obtained. For the fully separated modes (see Figs. , 4(a), and 5(b)), we derived an effective radius of r eff =0.33 µm (fine mode) and r eff =1.45 µm (coarse mode). Turning now to the layer from 3.3-4 km, the effective radii (total/fine/coarse) of 0.71/0.35/1.47 µm are similar to the values found for the lower dust layer (see Table 3, second panel). Another notable characteristic is that the effective radius in (3.3-4 km) is still well-retrieved despite the significant attenuation of signal in these heights. Note, however, that the uncertainties (Unc) for the upper layer are larger, see Table 3. The 15 lower particle concentration is already visible above 2.8 km in the lower values for the extinction and backscatter coefficients (Fig. 2). From the point of view of the microphysical retrieval it translates into a diminished surface-area-and total volume concentration ( Fig. 6(d)-(i)). Fig. 6 shows the column resolved microphysical properties from 1.5 to 4.3 km. There is an overall decreasing tendency with height for the microphysical properties v t and a t , but in the lower layer we see a rather steady pattern compared to the one in the upper layer. The total surface area and volume concentration in the lower layer is larger than in the

25
The refractive index is found to be 1.4 + 0.0497i for all layers between 1.5 and 4 km. Experimental findings have shown that such a high absorption can mainly be found near the sources of dust events. For instance, a soil sample from Burkina Faso (rich in hematite and kaolinite) was found in Wagner et al. (2012) using also a spheroid-particle-based inversion scheme to have IRI = 0.0495 ± 0.0206 at 305 nm. Moreover, Kandler et al. (2007) at Izana (Tenerife) found that the IRI increases a lot (up to 0.03) with decreasing particle size (< 0.5 µm), again due to the predominant hematite / soot component in smaller particles.

30
In addition, Kandler et al. (2011) found at Praia (Cape Verde) between 100 and 250 nm particle diameter a mode with high absorption of around 0.1i, which is produced by soot-sulphate mixture particles. It is also noted that the imaginary part varies depending on the source region for these dust periods. In particular, high absorbing small particles came from Mauritania, as seen also in the backward trajectories, whereas less absorbing ones came from Mali/Niger. AERONET's retrieval suggests, on the contrary, a refractive index with RRI = 1.52 (Var: 1.02%, Unc: 1.54%) and IRI = 0.0017 (Var: 21.13%, Unc: 64.24%), i.e., weakly absorbing particles but the uncertainty is high. The given values of the complex refractive index retrieved by AERONET represent spectral mean values in the four almucantar wavelengths. AERONET observes the whole column including the predominant marine particles in the altitudes below 1.5 km which are only weakly absorbing. Var (%) corresponds to a mean spectral variability, and the given Unc (%) corresponds to the uncertainty of the 5 mean retrieved parameters (spectrally) between the consecutive measurements. The retrieved IRI is lower than the usual values considered for Saharan dust particles . The latter study also finds a strong increase in IRI (up to five times higher at 637 nm) when the dust concentration is lower so that a soot-type absorber prevails. We should note that the limited resolution of the refractive index grid (Mieschka database) is an apparent shortcoming of our approach. The smoothing effect caused by high absorption raises the degree of ill-posedness (Samaras , 2017) which at the end can lead to (under-) 10 overestimation of the (RRI) IRI.
In this respect, Tesche et al. (2019) found that the use of depolarization input at any wavelength, i.e. 355 or 532 or 1064 nm, generally increases the retrieved values of the 532 nm SSA compared to the one obtained with the 3+2 input. The use of depolarization information leads to lower values of the imaginary part compared to the inversion in which the traditional 3+2 data set is used. 15 As another contribution to the discussion, Bi et al. (2018) found that a detailed shape change process from spherical particles to non-spherical particles is critical to the backscattering LDRs. To achieve high LDRs for nearly-spherical particles, the real part of the refractive index should be in the region of 1.3-1.7, and the imaginary part of the refractive index should be less than 0.01. The maximum size parameter of the high depolarization ratio depends on the imaginary part of refractive index. 20 The fine and coarse mode fraction of the volume concentration derived with our method was compared with the products obtained with the POLIPHON separation method Mamouri and Ansmann (2017) and with respective AERONET results. The comparison is shown in Fig. 7. With respect to the retrieval of the particle size distribution the spheroidal particle approximation used in our inversion approach is now well attested as the superior method compared to the Mie-scattering-based method.

Comparison with the POLIPHON Method and AERONET retrievals
The comparison with the POLIPHON shows reasonable agreement. This is a remarkable finding when keeping in mind that 25 spheroids by no means capture the true morphological nuances of dust and that the inversion is an ill-posed problem, i.e., small deviations on input may produce high deviations on output. The graphs in Fig. 7 show the same qualitative behaviour through the whole range for both modes. The average absolute difference in volume concentration (POLIPHON/AERONET vs inversion method) with respect to fine mode fraction is about 5.5 µm 3 cm −3 only. Note that for the coarse mode fraction the inversion makes a cutoff at r = 2.2 µm, which already explains the found differences to the POLIPHON and AERONET 30 results.

SALTRACE campaign case study 10 July 2013
A strong and long-lasting Saharan dust outbreak reached Barbados between 9 and 13 July 2013, substantially stronger than in June 2014 (case 1). On 10 July 2013, the BERTHA lidar system and the Falcon aircraft measured the dust plume which reached almost 5 km in height (Fig. 8). More specifically, this intense Saharan dust layer extended from 2.0-4.5 km height and showed a 532 nm particle depolarization ratio of 0.27±0.015 (Fig. 8, Table 4), and a lidar ratio of 51±3 sr. Downward mixing 5 of mineral dust to the marine aerosol layer was observed below 2 km height. The evening AOD was again high with values from 0.35 (1640 nm) to 0.50 (340 nm) according to the AERONET observations.
In the following, the lidar measurements from 10-11 July 2013 between 23:15 and 00:45 UTC are investigated. Local sunset was at 22:29 UTC. We focus on the height range from 2-3 km. The dust optical properties are summarized in Table 4.

Microphysical Properties and Comparison with Case 1
It is not surprising that the shape-size distributions in Fig. 9 and Fig. 10 are similar to the ones for case 1. The HYSPLIT backward trajectories showed similar long-range transport features (Haarig et al. , 2017a). SphInX produces again stable retrievals in the same sense as presented earlier, i.e, running twice for the same layer considering the range from 2-3 km (i) as a whole (single run) and (ii) breaking it into ranges of 100 m (multiple runs). The results are summarized in Table 5. There is virtually no difference between the retrieved values for whole range and the mean retrieved values for the split range and 15 with very low uncertainties. For the total volume concentration we get 57.23 and 56.10 µm 3 cm −3 , respectively. Fig. 12 shows the microphysical properties for the whole aerosol column from 1.5-5.0 km. Within the considered layer 2-3 km the volume concentration shows its largest value (Fig. 12(d)-(f)). The partly strong variability in the retrievals results and peaking structures in the profiles have to be interpreted with caution and are attributed to the sensitive impact of small uncertainties in the optical properties on the inversion products. 20 Fig. 11 shows the VCSD for three different layers. The decrease of the coarse mode particle concentration in the layer from 3-4.5 km (blue line) when compared to the respective coarse mode number concentration in the lower layer from 2-3 km (orange line) is in line with the reduced backscatter and extinction strength (Fig. 8). Focusing on the whole layer from 1.5-5 km height (pink lines in Fig. 5(b) for case 1, Fig. 11 for case 2), the maximum values for fine and coarse mode are 8 and 12 µm 3 cm −3 (case 1) and 29 and 44 µm 3 cm −3 (case 2). The values for case 2 are again larger since the dust event was more 25 intense. All this is in agreement with the results shown in Fig. 2 and 8 and presented in the Table 2 and 4.
When comparing the retrieval parameters in Tables 3 and 5, we find that the total effective radii are very close to each other in all three layers, namely, 0.73, 0.71 and 0.78 µm and that the parameters svcf, η eff and η var are equal. This points to the fact that the aged dust plumes originated from similar dust sources. On the other hand, the values for total surface-area (174.09, 70.32 and 215.61 µm 2 cm −3 ) and volume concentration (42.46, 16.51 and 56.10 µm 3 cm −3 ) differ significantly. They are larger for case 2, 10 July 2013, reflecting the stronger dust outbreak measured in July 2013.

Comparison with POLIPHON method
The retrieved fine and coarse mode volume concentration were compared with the respective POLIPHON results. The graphs in Fig. 13 show almost the same qualitative behaviour (monotonicity) throughout the entire height range for both modes and agree quite well. We should mention again that the coarse mode fraction derived by the POLIPHON method (thick green line) is not limited by the cutoff (r = 2.2 µm) and is therefore larger. The absolute difference in volume concentration between 5 POLIPHON method and inversion with respect to fine mode fraction is, on average, 7.81 µm 3 cm −3 .
Here we compare the inversion retrieval results with airborne measurements of particle volume concentrations. The FAL-CON observations are included in Fig. 13. By comparing the coarse mode volume concentrations, we found good agreement. 15 The mean relative difference is 39 %, and 14 % if we compare the inversion results with the particle counter 84-percentile value. The inversion algorithm overestimates the fine mode volume concentration. However, taking into account all differences, e.g., time and horizontal distances between the in situ and lidar measurements and measurement errors in the lidar and particle counter data the agreement is good. 20 We compare the volume concentration size distributions found by lidar-based-and AERONET retrievals. The flowchart in Fig. 1 provides an overview of the different retrieval procedures and products. Fig. 5(b) and Fig. 11(b) display the lidarbased volume concentration size distribution for three layers and Fig. 5(a) and Fig. 11(a) show the ones retrieved by the AERONET data analysis method. The latter figures show merely a very prominent coarse mode predicting very large particles.

Comparison with AERONET Sun-photometer
The complete absence of a fine mode is often a mathematical artifact when one of the expected modes is much more dominant 25 than the other. The smaller one is either suppressed or smoothed out. Although the lidar-based-and AERONET VCSDs cannot be directly compared, we can see that the volume concentration values from the two derivations (lidar and AERONET) are within the same order of magnitude for the layers 1.5-4.3 km (case 1, Fig. 5a) and 1.5-5.0 (case 2, Fig. 11a). Thus, considering the thickness of the whole dust layer of about 4 and 4.5 km, respectively, the estimated maximum lidar-based column volume concentration is about 0.12 µm 3 µm −2 and 0.20 µm 3 µm −2 which is comparable with the respective AERONET maximum 30 values of about 0.16 µm 3 µm −2 and 0.19 µm 3 µm −2 . Because the dust volume concentration was not height-independent within the dust layer from 1.5 to 4-5 km height and there were additional contributions from the layer below 1.5 km in both cases, the comparison of profile and column-integrated observations remains difficult. Nevertheless, the absolute differences are very small in particular in case 2, and the maximum peaks in Fig. 5 (a) and (b) and Fig. 11(a) and (b) occur at similar radius values, namely at 1.71 and 1.75 µm.
Finally, we included the column-integrated volume concentrations from the sun-photometer, separated into fine and coarse mode, in Figs. 7 and 13. Evidently, the qualitative comparison between the sun-photometer retrieval and the POLIPHON 5 method shows good agreement for the fine modes for both cases and a larger deviation for the coarse modes, especially for case 2.

Conclusions
In this study, we derived the microphysical properties for two dust events of different intensity by using input data from a Raman lidar system with 3 depolarization channels (355, 532, 1064 nm). The measurements were performed at Barbados 10 during the SALTRACE campaign in 2013 and 2014. Furthermore we conducted a series of comparisons with independent retrievals with the POLIPHON method, AERONET and in situ airborne measurements.
We followed a natural generalization of the Mie model in two dimensions by considering spheroidal optical efficiencies which refines the traditional size distribution to a shape-size distribution (radius and aspect-ratio dependent) and introduces new shape-related parameters. The inherently unstable nature of the retrieval calls for a careful selection of layers so that the 15 aerosol type (characterized by lidar ratio, AE) is reasonably confined, potentially turbulent parts are left out, and at the same time, the most intense parts of the optical profiles are included. This gave the opportunity to focus on the stability of the algorithm and its response to particle size differences and signal attenuation. By performing retrievals to subsequent sub-layers (100 m) within the entire dust layer we showed that the algorithm is able to sense the variation of particle size for case 1 as predicted previously by looking at AE profiles. We also saw that as long as the intensive properties are relatively constant the 20 microphysical properties calculated either within a layer as a whole or within smaller sub-layers have tenuous differences. A good quality of the optical data and adequate smoothing procedure obviously play a role as well.
Equally large particles were found for both cases (r mean eff = 0.75 µm), but case 2 was a more intensive dust event with 24% larger surface-area and volume concentrations. The shape parameters svcf, η eff and η var were identical which is an indication for similar origin of the aged dust particles. The shape-size distribution revealed two distinct modes with spherical and non-25 spherical contributions for both aged dust cases. The coarse mode was clearly prevailing as expected for dust-like particles. The retrieved IRI (0.05) indicate high absorption, which may occur under particular circumstances we briefly explored. In marked contrast, AERONET produces a possible underestimated IRI (0.0017).
A rough comparison for the retrieved volume concentrations between AERONET and our inversion could be done by summing up all shape contributions of a shape-size distribution to produce an analogue of a (1-D) size distribution and then scaling 30 it to the corresponding columnar range. We found that for both cases the orders of magnitude match, and notably there was even some quantitative resemblance in the volume concentration values (low absolute differences) and the radii of the maximum points of the coarse modes, all of which, again, are merely of qualitative interest. A comparison between the POLIPHON method and our inversion revealed remarkable similarities in the trend of the volume concentration with respect to height, with an average absolute difference of 6.6 µm 3 cm −3 (mean for both cases). In the same coarse-grained sense and taking into account the essential temporal and procedural differences, the in airborne situ measurements, available only for case 2, were consistent with the retrieved volume concentration, showing an average (84% percentile) deviation of 14 %.
In conclusion, the microphysical retrieval results derived with our approach provide an adequate description of the aged 5 dust cases, and the multiple comparisons with different and independently obtained data analysis and measurement approaches corroborate the strength of the inversion algorithm. As stressed throughout this paper, the substantial raise in complexity by considering the spheroidal approximation and the difficulty in retaining numerical stability forces the use of offline calculations (kernel database and discretizations). This in turn restricts significantly the solution framework in particle shape (aspect ratio), size (radius) and chemical composition (RI). A well-grounded and argued extension of the software tool for the future concerns 10 the considered size range. The current radius cutoff seems to underestimate the retrieved parameters and also force the shapesize distribution to have a steeper tail. Although all known approaches have several restrictions on chemical composition via a predefined refractive index grid and particle shape (e.g. use of experimental aspect ratio distributions), it is well-known that this practice has severe effects on the retrieval outcome and therefore these ranges have to be extended as well. For instance, the oscillatory trend of backscattering efficiencies fades away (smoothing effect) when particles become ever more non-spherical or 15 absorbing (Samaras , 2017) and therefore a greater resolution and diversity in RI and aspect ratios may facilitate the inversion.
Finally, it also remains open for future investigations whether it is necessary to extend the retrieval database with respect to a wavelength-dependent refractive index.

Appendix A
Here, we give definitions and explanations of ERQ reported in Section 2.1 (Samaras , 2017): 20 1. Variability of the solution space (Var). The algorithm is completed by ordering the solutions with respect to increasing error level and choosing a few of the first ones (calculate the mean solution out of a few least-residual solutions). This uncertainty percentage is the standard deviation of a sought parameter, which is derived by the chosen best solutions, divided by the mean value of the parameter. It describes how much a mean value of a parameter varies from all best solutions in the same solution space i.e. for a specifier error level . In this regard there are two interpretations of 25 such a statistical measure, both of which have to do with potential solution clusters. We rely on these clusters to find physically meaningful solutions, since previous experience with the sphere-particle approximation, showed that most of the mathematically acceptable solutions do not qualify physically. Therefore, on the one hand, the solution space should be variable enough in its full extent (all solutions), i.e. to produce clusters of more physically probable solutions. On the other hand, there should be a relative homogeneity for a small sample of "best" solutions (small residual error), reflecting 30 the ability of the examined method to recognize such clusters. More details on clusters and patterns in solutions spaces for real-data inversions can be found in Samaras et al. (2015a). For more than one data-sets, Var represents the mean variability of all produced solution spaces.
2. Randomness uncertainty (Unc). In simulation studies this is related to the stability of the examined method with respect to several repetitions of a numerical experiment of the same simulated atmospheric scenario but with different (random) instances of the same error level. By extension, this obviously characterizes the capacity of the method to reproduce well a possibly accurate result. The value of Unc is derived by first calculating the mean value of a parameter for every data-set of different data error, and then divide the standard deviation of these values by their mean. In other words, Unc 5 is a form of Var with respect to the different-error data-sets. For measurement cases, these data-sets could consist of optical data values related to consecutive smaller "sub-layers" of a layer which is partitioned in order to keep intensive parameters (e.g. AE, LR) relatively constant, and therefore Unc still makes sense as an additional measure of variability among the retrieved solutions.  The lidar observations (level 0 data and measured signals) and photometer data as well as the analysis products are available at TROPOS upon request (info@tropos.de). In situ data are available on request to bernadett.weinzierl@univie.ac.at.  volume concentration between the lidar-based retrieval (thick blue single line) and the AERONET sun-photometer retrieval (thick orange single line): A two-dimensional particle distribution is retrieved as a function of the radius and the aspect ratio using lidar data from a specified layer. By integrating over all aspect-ratio contributions we end up with an one-dimensional size distribution (reduced VCSD). The difference in units of the latter and the AERONET-based size distribution motivates the final step, i.e. taking the product of the reduced VCSD with the particular column height. properties from left to right: (b) total particle backscatter coefficient at 355 (blue), 532 (green) and 1064 nm (red line), particle linear depolarization ratio (PLDR) at 355, 532 and 1064 nm (dashed lines in the same colours), (c) backscatter coefficiens compited from the cross and parallel-polarized signal profiles, (d) extinction coefficient at 355 and 532 nm as well as (e) three Ångström coefficients. A vertical smoothing window of 99 bins (742m) was applied to the profiles in (b)-(e).