Advanced characterisation of aerosol size properties from measurements of spectral optical depth using the GRASP algorithm

This study evaluates the potential of using aerosol optical depth (τa) measurements to characterise the microphysical and optical properties of atmospheric aerosols. With this aim, we used the recently developed GRASP (Generalized Retrieval of Aerosol and Surface Properties) code for numerical testing of six different aerosol models with different aerosol loads. The direct numerical simulations (self-consistency tests) indicate that the GRASP-AOD retrieval provides modal aerosol optical depths (fine and coarse) to within 0.01 of the input values. The retrieval of the fine-mode radius, width and volume concentration are stable and precise if the real part of the refractive index is known. The coarse-mode properties are less accurate, but they are significantly improved when additional a priori information is available. The tests with random simulated errors show that the uncertainty in the bimodal log-normal size distribution parameters increases as the aerosol load decreases. Similarly, the reduction in the spectral range diminishes the stability of the retrieved parameters. In addition to these numerical studies, we used optical depth observations at eight AERONET locations to validate our results with the standard AERONET inversion products. We found that bimodal log-normal size distributions serve as useful input assumptions, especially when the measurements have inadequate spectral coverage and/or limited accuracy, such as moon photometry. Comparisons of the mode median radii between GRASP-AOD and AERONET indicate average differences of 0.013 μm for the fine mode and typical values of 0.2–0.3 μm for the coarse mode. The dominant mode (i.e. fine or coarse) indicates a 10 % difference in mode radii between the GRASP-AOD and AERONET inversions, and the average of the difference in volume concentration is around 17 % for both modes. The retrieved values of the fine-mode τa(500) using GRASP-AOD are generally between those values obtained by the standard AERONET inversion and the values obtained by the AERONET spectral deconvolution algorithm (SDA), with differences typically lower than 0.02 between GRASP-AOD and both algorithms. Finally, we present some examples of application of GRASP-AOD inversion using moon photometry and the airborne PLASMA sun photometer during the ChArMEx summer 2013 campaign in the western Mediterranean.


Introduction
Aerosol optical depth (τ a ) measurements that are complemented with angular radiance measurements are routinely used by the scientific community to infer the microphysical and optical properties of atmospheric aerosols (e.g. the AErosol RObotic NETwork or AERONET). The availability of suitable τ a measurements far exceeds the availability of suitable sky radiance scans, mainly because of partly cloudy conditions and scattering angle limitations. Additionally, night-time τ a measurements are becoming increasingly available through moon photometry, while night-time radiance scans are useless. Thus, it is desirable to derive meaningful aerosol information from spectral optical depth measurements when complementary radiance measurements are not available.
The goal of this work is to evaluate and demonstrate the full potential of using spectral optical depth measurement for characterising detailed properties of aerosol. Indeed, the magnitude and spectral dependence of optical depth are known to be related to the number and size distribution of aerosol particles. Therefore, ground-based observations of solar radiation attenuation were one of the first types of measurements analysed in pioneering studies devoted to atmospheric remote sensing algorithms (e.g. see Ångström, 1929, 1961Yamamoto and Tanaka, 1969;Quenzel, 1970;Grassl, 1971;King et al., 1978). However, in the last few decades, the focus of remote sensing retrieval development shifted towards an analysis of more complex observations: angular and polarimetric properties of transmitted and reflected diffuse radiation. In principle, such observations have high sensitivity, which allows complete and accurate characterisation of the aerosol features. However, modelling and interpretation of the diffuse radiation is significantly more complex than the analysis of the direct sunbeam. Therefore, the complexity and efficiency of the retrieval algorithms have significantly improved compared to the those originally developed for the interpretation of aerosol optical depth. At the same time, in a number of practical situations the interpretation of τ a data alone remain interesting for the community. Moreover, due to evolution of With this purpose, we make use of the recently developed GRASP (Generalized Retrieval of Aerosol and Surface Properties; see Dubovik et al., 2014) algorithm and software. The code implements statistically optimised fitting of diverse observations using multi-term least square method (LSM) concept (e.g. see Dubovik, 2004). Correspondingly the retrieval is organised as a solution search in continuous space of solutions without traditional use of pre-calculated look-up tables. The GRASP algorithm is highly versatile and can be applied to a large variety of remote sensing measurements (sun photometer, lidars, satellite imagers, etc). The GRASP concept, which originated from the  algorithm, has been successfully used for 15 years to process observations of AErosol RObotic NETwork (AERONET) ground-based sun/sky-radiometers (Holben et al., 1998), engendering a large number of worldwide climatologies (to cite some Dubovik et al., 2002a;Smirnov et al., 2002a, b;Eck et al., 2005Eck et al., , 2010Giles et al., 2012;Toledano et al., 2012). During this period the algorithm has certainly evolved and several useful modifications were realised Li et al., 2009). The GRASP development inherited of all these retrieval advances.
In addition, the GRASP algorithm has a highly flexible forward model that makes it a convenient tool for sensitivity and tuning studies . This flexibility was one of the motivating factors for performing these studies. Indeed, the lack of scattering information in the τ a observation obliges the retrieval developer to search for an optimum aerosol model and an adequate set of a priori constraints (e.g. parameterisation of size distribution). For example, this study suggests that reasonable results can be obtained by approximating size distributions as bimodal log-normals, which can be described using only six parameters (volume median radius r V i , standard deviation σ V i and volume concentration C V i for fine and coarse mode) instead of binned size distributions (in the case of AERONET standard inversion 22 bins logarithmically equidistant between 0.05 and 15 μm). Moreover, we assume that the complex refractive index and the sphericity parameter are known. The different assumptions and their consequences are discussed in detail in the sensitivity analysis presented in Sect. 3.
A practical motivation for the present study is the large number of "optical depth only" measurements that exist in the ground-based networks. For instance, AERONET includes around forty direct sun measurements per day in its standardised sequence of measurements (in cloudless conditions); only about eight of these sequences are coincident with skyradiance measurements that are suitable as input to the AERONET inversion code 1 . Some sites have new instruments that increase the number of valid clear-sky radiance retrievals to about 16 per day (through the use of Hybrid scans) and the number of τ a measurements to about 200 per day. Nonetheless, there is still a large amount of data that do not contain the sky-radiance scans necessary for full AERONET inversions of size distributions, complex refractive indices and single-scatter albedos. Moreover, many AERONET sites are plagued by several months of partial cloudiness. In these situations, there are no angular The analysis presented here is focused on τ a measurements in the spectral range of the Cimel sun photometers that are used in the AERONET project. This approach allowed us to compare the results obtained by the GRASP-AOD inversion with the standard AERONET inversion. Consequently, the simulation tests proposed in Sect. 3 are done using the wavelengths in the AERONET network, and the main analysis with real measurements in Sect. 4 corresponds to AERONET data. However, we emphasise that the code is not restricted to AERONET measurements, and some example applications of the code to moonphotometer data are also shown in Sect. 4. Finally we present some retrievals from data obtained at different heights during the ChArMEx campaign  with the new airborne sun-tracking photometer PLASMA (Karol et al., 2013) fully designed at the LOA (Laboratoire d'Optique Atmosphérique, University of Lille).

Inversion strategy (GRASP-AOD)
As commented in the introduction, we make use of the recently developed GRASP (Generalized Retrieval of Aerosol and Surface Properties; see Dubovik et al., 2014) algorithm and software (more information and a free version of the code can be obtained in http://www.grasp-open.com/) to assess the potential of spectral aerosol optical depth measurement for characterising aerosol properties. The main concept of GRASP is the multi-term least square method (LSM) (for details of the inversion methodology see King, 2000, Dubovik, 2004 and also and also King and Dubovik, 2013). GRASP is designed to be applied to a broad number of remote sensing measurements (e.g. sun photometer, lidars, satellite images).
The multi-term LSM method solves the following system of equations: f* = f a + Δ f 0* = Ga + Δ g a* = a + Δ a . (1) The first equation in Eq. (1) describes the physical relation between the set of measurements f and the vector of unknowns a. The symbols ∆ f denote the uncertainty in the real measurements f*. Note that, in our case, the set of measurements only contains the spectral aerosol optical depths (defined in Eq. A2). In GRASP code the atmospheric aerosol is modelled as an ensemble of randomly oriented spheroids and the aerosol optical depth is modelled as follows: τ a λ = ∫ lnϵ min lnϵ max ∫ lnr min lnr max C ext λ , n, k, r v r ⋅ dn ϵ dln ϵ dV r dln r dln r dln ϵ , where ν(r) is the volume of particle, C ext is the cross section of extinction, λ is wavelength and n and k are real and imaginary parts of the refractive index. The aerosols are approximated as spheroids (Mishchenko et al., 2002) with ϵ being the axis ratio (ϵ = a/b, a is the axis of spheroid rotational symmetry, b is the axis perpendicular to the axis of spheroid rotational symmetry) and r is the radius of the equivalent sphere. As discussed by Mishchenko et al. (1997), the usage of r and ϵ is convenient for separating the effect of particle shape and size in analyses of aerosol mixture light scattering. Then the functions dV r dln r and dn ϵ dln ϵ denote the volume distribution of the spheroids (for the total column) and the number particle shape (axis ratio) distribution accordingly. Dubovik et al. (2006) have demonstrated that the particle shape distribution dn ϵ dln ϵ for the non-spherical fraction of any tropospheric aerosol can be approximated as constant over the particle size distribution. This assumption simplifies Eq.
(2) and the aerosol extinction is calculated for the retrieval as a mixture of spherical and non-spherical fractions. Moreover, in order to perform fast and accurate calculations, the integrals are replaced by sums of precalculated kernels as follows: τ a λ = τ sph λ + τ nons λ = ∑ i = 1, …, N r C sph K τ sph λ , k, n, r i + 1 − C sph K τ nons λ , k, n, r i dV r i dln r , where N r is the number of bins used to represent the size distribution, C sph is the fraction of the spherical particles and K τ sph and K τ nons are the kernels for spherical and non-spherical particles respectively. The complete information on the forward model and the detailed calculation of the kernels can be gained from Dubovik et al. (2006Dubovik et al. ( , 2011.
The second equation in Eq.
(1) represents a priori smoothness constraints on the retrieved characteristics. They are typically applied to eliminate unrealistic strongly oscillating dependencies in the retrieved characteristics. Specifically, the study of  showed that assuming zeros (0* -zero vector) for the derivatives of retrieved aerosol characteristics allows the elimination of strongly oscillating solutions with high derivatives. The matrix G is composed of coefficients allowing the numerical estimates of derivatives of function y(x) using discrete values a i = y(x i ). These constraints are normally used to smooth the retrieved size distribution and the spectral dependencies of the refractive index. The symbol ∆ g accounts for the uncertainties in the a priori constraints.
However, here in the GRASP-AOD application, the refractive index and the sphericity parameter are assumed as known. The particle size distribution is characterised as bimodal log-normal and is defined by six independent parameters, volume median radius r V i μm , geometric standard deviation σ V i and particle volume concentration C V i μm 3 μm −2 , with i = f, c for the fine and coarse mode respectively. Therefore, the second equation is not used in this particular application of the GRASP code. It should be noted that we initially tested binned size distributions but we rapidly observed that strong smoothness constraints were required in order to assure realistic retrievals. Hence, we decided to model the size distribution in terms of log-normal functions. Note that many physical models use lognormal approximations to represent size distributions (e.g. see Whitby, 1978;Shettle and Fenn, 1979;Koepke et al., 1997;Hess et al., 1998) and this aerosol representation generally agrees well with observations (Tanre et al., 1988;Remer and Kaufman, 1998;Dubovik et al., 2002a). Therefore, applications with a need to represent aerosols using a limited number of parameters naturally choose this concept for modelling particle size distribution. In general, bimodal log-normal functions are fully sufficient to interpret variability of aerosol optical depth (to cite some Eck et al., 1999;O'Neill et al., 2001a, b;Schuster et al., 2006, etc). However, some AERONET retrievals suggest that the particle size distribution is not always perfectly log-normal, as some size distributions are characterised by asymmetrical mode shapes (e.g. Dubovik et al., 2002a;Eck et al., 2005Eck et al., , 2010. Moreover, some size distribution retrievals have a pronounced trimodal structure in the particle size range addressed by the inversion, such as observed in volcanic aerosol plumes (Eck et al., 2010) or aerosol retrievals located near fog or cloud Li et al., 2014). Obviously, our strategy, which is based upon simplified bimodal size distributions would not provide correct retrievals in these specialised situations.
The last equation in Eq.
(1) shows the possibility of using a priori constraints on actual values of any retrieved parameter a i , and a* is the vector of a priori estimates of a i . The symbol ∆ a represents the uncertainty in the a priori constraints of a*. For optimised accounting of those uncertainties, the solution of the system defined in Eq.
(1) is given by the minimisation of the following quadratic form: where the weighting matrices W and the Lagrange parameters γ are defined as follows: where ϵ f 2 and ϵ a 2 are the first diagonal elements of the corresponding covariance matrices C f and C a respectively. During the retrieval, we assume that all input data have log-normal error distribution Dubovik, 2004), which means that the measurements and the retrieved parameters are used in the logarithm space.
Finally, it should be indicated that one of the recent success of GRASP code has been the easy adaptation of the multi-pixel retrieval concept in the methodology proposed by multiterm LSM (Eq. 1). The strategy was developed by Dubovik et al. (2011) in order to overcome the difficulties related to the limited information of the satellite observations over a single pixel. The multi-pixel retrieval regime takes advantage of known limitations on spatial and temporal evolution in both aerosol and surface properties. Specifically, a large group of pixels are inverted simultaneously, using a priori constraints on the temporal and spatial variability of the retrieved parameters. The concept has been expanded in the latest version of the GRASP algorithm, where the measurements can also belong to different remote sensing instruments. The present study focuses on the potential of a single set of aerosol optical depth measurements, so that the multi-pixel inversion will not be used in this work. Nevertheless, in a number of recent promising applications the multi-pixel approach has been used for synergy retrievals when τ a observations were combined with other colocated but not coincident measurements (Lopatin et al., 2013(Lopatin et al., , 2016. Moreover, it appears to be a powerful capability for the interpretation of τ a -only measurements in the future.

Sensitivity studies
The description of any algorithm must answer several inherent questions regarding the retrieval: stability of the inversion, confidence in the retrieved products, dependence up on the a priori assumptions, etc. In our particular case, the main challenges are to (a) identify the reliability of the six parameters which describe the bimodal log-normal size distribution (and of the secondary products derived from them), (b) check the effects of possible errors in the aerosol optical depth measurements, and finally, (c) analyse the consequences of assuming the refractive index and sphericity parameter as "known" during the retrieval process.
Several simulation tests are considered in this section to address these points. First, a selfconsistency analysis, including multiple variations of the initial guess, will be carried out to check the stability of the retrieval. Next, simulated errors in the aerosol optical depth measurements will be introduced to determine the ramifications on the retrieved properties. The last two studies will include biases in the refractive index and sphericity parameter assumptions.

Aerosol models
The numerical tests in this paper are based upon the climatology provided by Dubovik et al. (2002a), which utilises about 10 years of real aerosol retrievals with AERONET's Version 1 processing. AERONET has subsequently implemented a Version 2 aerosol retrieval product, but the single-scatter albedo climatology that is based upon this newer processing scheme is within 0.02 of the Dubovik et al. (2002a) climatology for the same aerosol type . Specifically, we have selected aerosol properties from six different sites with which to carry out the simulations. Three of the sites are dominated by fine-mode aerosols (from less to more absorbing): Goddard Space Flight Center (GSFC, Maryland, USA), which represents urban weakly absorbing aerosol, Mexico City (Mexico) representing urban absorbing aerosol and Mongu (Zambia), which corresponds to a strongly absorbing biomass burning aerosol. Additionally, we have selected Bahrain (Bahrain) and Solar Village (Saudi Arabia) as examples of mixed desert dust and pure desert dust respectively. Finally, we have used the aerosol properties at Lanai (Hawaii, USA) as an example of maritime aerosol. For each example, we have considered at least three aerosol loads: τ a (440) = 0.3, 0.6 and 0.9. The τ a (440) = 0.9 case is omitted for Lanai, however, since this case would be excessively unrealistic (given the typical values observed at that site). On the other hand, we have considered four extra cases which can be used to analyse the retrievals at low-and high-load aerosol conditions: τ a (440) = 0.1 for Lanai and for GSFC and τ a (440) = 1.5 for Mongu and Flight Center), MEXI (Mexico City), ZAMB (Mongu), BAHR (Bahrain), SOLV (Solar Village) and LANA (Lanai). In the same way, we will use different indexes for the different aerosol loads: index 0 for τ a (440) = 0.1, index 1 for τ a (440) = 0.3, index 2 for τ a (440) = 0.6 index 3 for τ a (440) = 0.9 and index 4 for τ a (440) = 1.5.
The aerosol properties of all the examples considered here are represented in Table 1. The first parameter is the reference value of the aerosol optical depth and the rest of the parameters are derived from the expressions in Dubovik et al. (2002a). The parameters in the top part of Table 1 are used to describe the bimodal log-normal size distribution: volume median radius r V i μm , geometric standard deviation σ V i and particle volume concentration C V i μm 3 μm −2 for the fine and coarse modes. It also contains the sphericity parameter; we use a sphericity parameter of 0 (i.e. all particles are non-spherical) for all the Solar Village cases except the one with τ a (440) = 0.3 (SOLV1), which has an  Figure 1 illustrates the size distributions for the examples with the three aerosol loads common to all sites created from the parameters described in Table 1. Each aerosol load for every case is represented by different lines: the cases with τ a (440) = 0.3 by a dashed line, the cases with τ a (440) = 0.6 by a solid line and cases with τ a (440) = 0.9 by a dashed-dotted line.
The three cases with a predominance of fine mode are plotted in the left panel. The two desert dust examples and the case of maritime aerosol are represented in the right subfigure.
In the cases with a predominance of the fine mode, left panel of Fig. 1, the mode radii vary with τ a for both modes. The feature was marked in the climatology study from  and depicted in the summarising formulas, in Table 1 of that study, which are the basis of our aerosol models here. The other three aerosol cases, which are represented in the right panel, do not present this property (as can be gained from the figure or directly from the values in Table 1).
The simulated aerosol optical depth values obtained by running the GRASP forward module with the values from Table 1 are shown in Table 2. We have also included the resulting Ångström exponent, which is computed as a linear regression of τ a (440), τ a (675) and τ a (870). In the last two columns, the values of the fine-mode aerosol optical depth at 500 nm and the effective radius are depicted, since they will be analysed in our different sensitivity studies (as derived secondary products). In order to provide a graphical representation of the tendencies from the τ a simulated values, they are represented for the cases with τ a (440) = 0.3 and 0.6 in Fig. 2, in logarithmic scale. In the figures, we observe a linear relationship with a slight curvature between 340 and 1020 nm for all of the sites. This curvature is negative for aerosol size distributions dominated by fine-mode aerosols and positive for the desert dust cases (especially for SOLV). This is consistent with previous works (Kaufman, 1993;Eck et al., 1999Eck et al., , 2001Reid et al., 1999;Schuster et al., 2006).

Self-consistency analysis
Our strategy for the simulation tests has been adapted from previous work Torres et al., 2014). A general scheme of the procedure is drawn in Fig. 3. The methodology consists basically of inverting the aerosol optical depths from Table 2 (which were obtained by running the forward module using the aerosol properties described in Table 1), but introducing some modifications in order to "test" the sensitivity of the code regarding the different aspects already mentioned: multiple initial guesses, variation of the refractive index and uncertainty in the aerosol optical depth data. Nevertheless, we first present a self-consistency study in which we do not make any modifications during the inversion process.
The first six columns of Table 3 represent the differences between the output of the selfconsistency analysis (obtained from the inversion of the τ a values from Table 2 and considering the values of the refractive indices also from Table 2) and the original values of the size distribution (Table 1). In the last two columns, we show the differences in the selected secondary products: τ f (500) and r eff (theoretical values in Table 2).
The capacity of the GRASP-AOD inversion to discriminate fine-and coarse-mode aerosol optical depths is one of the most important results that can be gained from Table 3: differences between the reference and retrieved values of τ f (500) are always less than 0.004.
This result is even better for those cases with a predominant fine mode, where a maximum difference of 0.001 is obtained.
The retrieved parameters that characterise the fine mode also show a good agreement with the reference values (from Table 1). Thus, the maximum difference observed for r V f is 0.002 μm (LANA1). Once again the comparison is better for the cases with a prevailing fine mode, where no difference are observed up to the third decimal. The differences in σ V f do not exceed 0.01 for the aerosol cases with a predominance of fine mode, while they are a bit higher in the rest of the cases: up to 0.015 for the maritime aerosol and a maximum difference of 0.019 for BAHR1. Finally, the divergences in the concentration C V f are under 0.002 μm 3 μm −2 in all the examples.
The retrieval of the coarse mode is less accurate than the retrieval of the fine mode. We expected this result, since the wavelengths used in this study are less sensitive to the radii in the range of the coarse mode than to those in the range of the fine mode; see, for instance, the discussion about extinction efficiency and size parameter in Lenoble et al. (2013) or in chap. 3 of Coulson (1988). Differences in r V c are up to 0.31 μm for the cases with fine-mode predominance and a bit smaller for desert dust and maritime aerosol, up to 0.21 μm, since those cases have larger coarse-mode aerosol optical depth and therefore more information.
The discrepancies in the concentration C V c are at maximum 0.007 μm 3 μm −2 (MEXI3) for the cases with fine-mode predominance and 0.02 μm 3 μm − 2 (SOLV2) for the cases with a prevailing coarse mode. For both cases these differences represent 6-7 %.
Here we should point out the strong connection between the retrievals of r V c and C V c .
Aerosol examples with the best characterisation of r V c also correspond to those with the best characterisation of C V c : GSFC in the case of fine-mode predominance and SOLV for the cases with a prevailing coarse mode. Those cases with an overestimation in the volume median radius also display an overestimation in the concentration: all cases of MEXI and BAHR1. Note that, for the radius range of the coarse mode, the extinction efficiency diminishes as the radius grows. Therefore an overestimation of r V c needs to be optically compensated increasing the coarse-mode concentration. In the same way, the cases that present an underestimation in the volume median radius show an underestimation in the concentration: all cases of ZAMB and BAHR3.
Finally, the effective radius is also computed and compared with the reference values in Table 3. Note that this parameter is not directly derived from the inversion and is computed from the retrieved values of the bimodal log-normal size distribution. Its accuracy, therefore, will be conditioned by the accuracy of the retrieved parameters. The differences between the effective radii from the reference and those retrieved in the self-consistency are under 0.018 μm, except for the cases of BAHR where differences up to 0.074 μm are obtained.

Initial guess variation
In this section, we compare the sensitivity of the inversion results to the initial guess. More specifically, we test three different initial guesses for each of the six retrieved parameters, which results in a total of 729 inversions (3 6 ) for the 17 aerosol examples. Note that we need not verify that the initial guesses are consistent with the inversion results; rather, we analyse the stability of the code with respect to the initial guess. That is, a perfect retrieval code would give us the same correct results, independently of the initial guess. Table 4 shows the default initial guesses used for the first iteration of the GRASP-AOD code, which depend upon the Ångström exponent and the aerosol optical depth at 440 nm. Note that these initial guesses can be modified if ancillary information for specific applications exists. Table 5 shows the variety of initial guesses chosen for this portion of our study. These values are computed from the expressions in Table 4, with a modal variation of ±25 % for r V i , an absolute variation of ±0.1 for σ V i and ±40 % for C V i . Here, we only show the results obtained for the examples with aerosol load τ a (440)=0.6, since they are similar to the results obtained with the other two aerosol loads.
Before analysing all cases in detail, the 729 inversion retrievals for GSFC2 and BAHR2 are represented in Fig. 4 to get a visual idea of the results obtained in this analysis. Black solid lines are used to represent the reference values, while the grey shadow contains all the retrievals. It can be observed that the fine mode is well characterised, whereas there is a larger uncertainty in the description of the coarse mode. The better performance of the fine mode compared to the coarse mode was expected, since as commented, small particles have higher Mie extinction efficiencies at the selected wavelengths in the analysis.
In Table 6, average values and standard deviations of the retrieved parameters calculated from the 729 inversions are shown for each aerosol type. The last two columns include the secondary products: τ f (500) and the effective radius. We focus on the standard deviations in Table 6 for this second study, since they indicate the stability of the retrieved parameters. If we analyse τ f (500), we see that the standard deviations are between 0.001 and 0.002 when the fine mode dominates, while they are between 0.005 and 0.006 for those cases in which the coarse mode dominates. This result denotes that the separation of the modes is very stable regardless of the initial guess values.
Similarly, we observe that the standard deviations for r V f and C V f are less than 0.01 in all cases (units [μm] and [μm 3 μm −2 ] respectively). The variation is even smaller for the cases with a fine-mode predominance, with a maximum value of 0.002 in the two parameters. In relative terms this variation is less than 3% for both cases. These results indicate that the retrieved r V f and C V f have little sensitivity to the initial guess and the variability is practically negligible in the cases with a dominant fine mode. The coarse modes in Table 6 indicate that the standard deviations are between 10 and 15 % for both r V c and C V c , regardless of whether the fine mode or coarse mode dominates.
The characterisation of the geometric standard deviation of fine mode, σ V f , is more sensitive to the initial guess selection than the other fine-mode parameters: the standard deviation is lower than 0.01 for the cases in which the fine mode dominates and up to 0.06 for the rest. This result denotes a low sensitivity of the retrieval to the width of the fine mode. The results are a bit worse for the standard deviation of σ V c : in cases with a predominant coarse mode the values are between 0.02 and 0.04, while for the cases with a prevailing fine mode the standard deviations are up to 0.08. Therefore, a good selection of the initial guess of σ V f and σ V c is very important for a suited retrieval of both parameters.
The standard deviation of the effective radius is between 5 and 8 % in all cases with a finemode predominance, and between 10 and 14 % in the cases with a prevailing coarse mode. The largest influence of the predominant mode in the calculation of the effective radius explains the result obtained here and, as commented, the fine mode is better characterised than the coarse mode.

Simulation of aerosol optical depth errors
The purpose of this section is to analyse how aerosol optical depth errors affect the inversion and its products. Following the scheme in Fig. 3, we propose to introduce random variations ∆τ(λ) in the aerosol optical depth values (from Table 2) for each channel. The random variations are obtained from a normal distribution which has been differently defined for each wavelength, taking into the account the uncertainty estimation for AERONET field instruments. These uncertainties, ϵ i , are at a maximum of 0.01 in the visible and near infrared and increasing up to 0.02 in the ultraviolet wavelengths (for more information see the Appendix A and the references therein, in particular Holben et al., 1998 andEck et al., 1999). Due to the temperature sensitivity of the detector at 1020 nm and the water vapour absorption at this wavelength, we have also considered ϵ i = 0.02 for this channel. To define the normal noise distribution, we have assigned ϵ i = 3σ which defines −ϵ i ϵ i as a 86 % confident interval.
The number of tests has been set to 1000 for each aerosol case. The multi-term LSM formulation (basis of GRASP-AOD inversion strategy) allows us to account for different uncertainties at each wavelength. Therefore, here and in those applications of GRASP-AOD where the errors have a spectral dependency, the known uncertainties of each wavelength can be introduced in the diagonal elements of the covariance matrix (Eqs. 4 and 5), this information being used in the inversion procedure. This is the case in the analysis that follows. At the same time, we applied a filter to the analysis, so that the retrievals with absolute residual (defined from Eq. 1 as f* − f a /N; with f* = τ) larger than 0.014 or relative residual (the same but dividing by τ(λ) for each wavelength) larger than 15% have been eliminated.
In order to summarise, only the cases with τ(440) = 0.3 and τ(440) = 0.9 for each type are presented in the general analysis. The particular examples with the extreme aerosol loads (τ(440) = 0.1 and τ(440) = 1.5) will be analysed in a different subsection, and finally, the effects of reducing the spectral range to 440-1020 nm will be discussed in the last part of the section. For GSFC, we observe that the red line (all errors +ϵ i ) gives larger values of the size distribution than the reference in the coarse mode, whereas the opposite happens for the blue line. Note that the mode radii are practically the same as those obtained in the self consistency test for both the red and blue lines; the discrepancies are only observed in the concentration values. In the fine mode, however, the analysis is not that simple: for the +ϵ i case, the modal radii are reduced and the fine mode widens so that the concentration increases. All these effects considerably increase the volume size distribution for radii smaller than 0.12 μm. This occurs because ϵ i is larger at the two ultraviolet wavelengths than at the other wavelengths, so the Ångström exponent increases and the retrieved mode radius decreases.

General analysis-So
Similar effects are observed for SOLV3 in the coarse mode. For +ϵ, there is a decrease in the modal radius (0.10 μm) and in the mode concentration (4 %). The opposite occurs for the case in which all the errors equal −ϵ i . In both cases, the variations in the standard deviation are under 0.01. For SOLV1, however, there is a small decrease in the concentration and in the modal radius for both the red and blue lines. The decrease in the radius is more significant for +ϵ i (0.24 μm) than for −ϵ i (0.08 μm), while the reduction in the concentration is similar for both cases: 6 % for +ϵ i and 7 % for −ϵ i . The standard deviation of the mode increases of 0.07 for +ϵ i and decreases 0.06 for −ϵ i .
Looking at the grey area, it seems that the uncertainties in r V f and C V f are significantly reduced for the cases with the highest aerosol load, GSFC3 and SOLV3 (τ a (440) = 0.9), with respect to GSFC1 and SOLV1 (τ a (440) = 0.3). The same effect is observed in the coarse mode, for r V c and C V c , in all cases. Table 7 presents the values and standard deviations of the retrieved parameters calculated for the simulation of random error analysis. The scheme of the table is the same as in the initial guess variation (Table 6), but to the last column we have added the ratio of the inversions that pass the filters out of 1000 simulations for each aerosol case. As in the study of the initial guess variation, standard deviation values will be used to check the stability of the retrieval, though in this case, with respect to the τ a errors.
Average values of τ f (500) obtained in the study do not differ by more than 0.005 from those retrieved in the self-consistency analysis (see Table 2). The standard deviations of τ f (500) are between 0.008 and 0.009. This result is consistent with the normal variations of ϵ i = 3σ introduced in the aerosol optical depth in the study. Figure 6 summarises the differences between the mean values obtained in the error simulation study (Table 7) and the values from the self-consistency test for the parameters representing the bimodal log-normal size distribution. Similarly, the standard deviations obtained in the analysis for the same parameters (also in Table 7) are presented in Fig. 7.
We can observe in Fig. 6 that the largest differences between the mean values of r V f in the error simulation study and the values from the self-consistency test are 0.007 μm and are obtained in the case of LANA2. The standard deviations of r V f (either in Table 7 or Fig. 7) for all of the cases are between 0.013 and 0.018 when τ a (440) = 0.3, and between 0.005 and 0.009 when τ a (440) = 0.9. There is an significant diminution in the standard deviation of r V f for all examples as the aerosol load grows. This effect is more visible in the cases with a prevailing fine mode (GSFC, MEXI and ZAMB), where the variations are around 0.015 for the cases with τ a (440) = 0.3 and diminish to 0.005-0.006 for the cases with τ a (440) = 0.9.
Meanwhile, the cases with a prevailing coarse mode (SOLV, BAHR and LANI) have r V f variability of 0.014-0.018 that diminishes to 0.007-0.009.
The mean values of the fine-mode concentration, C V f , are almost identical to those obtained in the self consistency analysis with maximum differences of 0.002 μm 3 μm −2 or a 4% relative bias (Fig. 7). From Fig. 7, we can observe that the standard deviations are around 12-16% for the cases with τ a (440) = 0.3 and they decrease to 4-6 % for the cases with τ a (440) = 0.9.
For σ V f , the differences between the averages in Table 7 and the self-consistency analysis are typically around 0.01, except for the case of BAHR1 with a maximum difference of 0.025. The standard deviation is around 0.1 for the cases with τ a (440) = 0.3 and they are between 0.03 and 0.05 for the cases with the largest aerosol load in Table 7.
The average values here and those obtained in the self-consistency do not differ by more than 0.15 μm for all cases. Note that if we compare the average values of r V c and the input values, the differences are much larger: around 0.4 μm for the aerosol cases dominated by the coarse mode and around 0.2 μm for the cases with a prevailing fine mode (similar results to those obtained between the input values and the self-consistency analysis in Table 3). The standard deviations of r V c are between 0.03 and 0.1 (see Fig. 7). These variations represent between 2 and 4 % in relative terms. There is, in general, a decrease in the standard deviation as the aerosol load grows but it is not as significant as for r V f .
The analysis of C V c indicates differences of up to 0.002 μm 3 μm −2 between the average values in Table 7 and those obtained in the self-consistency study for the cases with a prevailing fine mode. These differences can represent more than 8 % of the value as can be gained in Fig. 6. For the cases dominated by the coarse mode, the differences are under 5 %. Generally, there is a decrease in the relative differences as the aerosol load grows (except for MEXI and BAHR that do not present any errors when τ a (440) = 0.3). Back to the cases with a fine-mode predominance, Table 7 shows that the standard deviations of C V c are relatively high with values around 0.01 μm 3 μm −2 . Note that, for the cases with the lowest aerosol loads, sometimes these variations represent around 40 % of the total value as we see in Fig.  7. Certainly, variations of ±0.01 for the longest wavelengths should have a large influence on the retrieval of the coarse mode considering the values of the aerosol optical depth at these wavelengths (τ a (1640) ∼ 0.03-0.04; see Table 2 for GSFC1, ZAMB1 and MEXI1). The differences, in relative terms, are considerably smaller for the cases with the highest aerosol loads. For the cases with a prevailing coarse mode, the standard deviations of C V c are around 10 % for the cases with τ a (440) = 0.3 and between 4 and 5 % for the cases with the largest aerosol loads.
The differences between the mean values of the standard deviation of the coarse mode, σ V c , and those retrieved in the self-consistency analysis are under 0.01 in all examples except for MEXI and LANA, which have differences around 0.025-0.03, and GSFC, where differences up to 0.05 are observed. The standard deviations are under 0.025, with the exception of SOLV1 and LANA1, for which a value of 0.05 is registered. For these two cases, there is a significant decrease in the standard deviation when the aerosol load grows that is not noticeable in the rest of the aerosol models.
There is also a good agreement between the average values of the effective radius in Table 7 and those retrieved in the self-consistency test. The differences are under 0.02 μm for all the examples. The standard deviations are between 0.02 and 0.05 for the cases with a predominant fine mode and they are between 0.05 and 0.095 for the cases dominated by the coarse mode. Note that, for all cases, the variations represent between 15 and 20 % in the cases with τ a (440) = 0.3 and between 7 and 10 % for the cases with the largest aerosol load in Table 7.
Almost all the inversions pass the filters for the cases with a prevailing coarse mode. For the examples dominated by the fine mode, the ratio is a little bit lower, around 90 %. The filter that cuts most of the retrievals is the one imposed to the relative differences, and it is due to the relatively low values of τ at infrared wavelengths compared to the values of the random errors.

Low and high aerosol load conditions-The effects of τ a errors on the
examples with the extreme aerosol loads (GSFC0, LANA0, ZAMB4 and SOLV4) are analysed in this subsection. Figure 8 contains the retrievals for two of our particular aerosol test cases: GSFC0 (τ a (440) = 0.1) and SOLV4 (τ a (440) = 1.5). They have been chosen since they were selected to get a visual representation in the previous subsection in Fig. 5. The main output from Fig. 5 was that the uncertainties in the characterisation of the size distribution were significantly reduced for GSFC3 and SOLV3 with respect to GSFC1 and SOLV1. This tendency seems to be confirmed here if we compare SOLV4 to SOLV3, where we observe that the grey area has been significantly reduced. On the contrary, the case of GSFC0 presents larger uncertainties compared to GSFC1. Note that for GSFC0, the two cases in which all the wavelengths have the same error, +ϵ i and −ϵ i , do not pass the filters imposed onto the residuals, and that explains why both lines appear outside the grey area in the coarse mode. We preferred to keep them in the figure to point out that they follow similar patterns to their equivalents in Fig. 5. For instance, in the case of +ϵ i (red line), the modal radius decreases and the concentration increases with respect to the reference. The standard deviation of the mode also rises though not as much as in Fig. 5. On the other hand, red and blue lines in the case of SOLV4 do not differ much from the reference in the right panel of Fig. 8.
The values and standard deviations of the retrieved parameters calculated for the simulation of random errors for our four particular aerosol test cases are shown in Table 8. The structure  of the table is the same as Table 7, including the fraction of the inversions that pass the residual filters for 1000 simulations (last column); note that this ratio is particularly low for the case GSFC0 (0.39). However, ϵ i represents 91 % of the total value of τ a at 1020 nm and 66 % of τ a at 1640 nm for this case with low aerosol loading. The filter imposed to the relative AOD residuals in the inversion process (i.e. relative residual < 15 %) eliminates most of the retrievals.
The differences between the averaged τ f (500) values in Table 8 (with extreme aerosol loading) and those retrieved in the self-consistency analysis are all under 0.001. Moreover, the standard deviations of τ f (500) are between 0.007 and 0.009 as in the general study.
These same results were obtained while analysing the 1000 retrievals without filters in the residuals. Therefore, we can conclude that the retrieval of τ f (500) is consistent with the accuracy of τ a measurements and independent of the aerosol load.
Similar to τ f (500), the differences observed between the averaged values of the bimodal log-normal parameters (r V i , σ V i , and C V i ) in Table 8 and those retrieved in the self- Torres et al. Page 17 consistency analysis are almost negligible. Thus, the main goal is to check the evolution of the standard deviations with respect to the general analysis.
In the cases with τ a (440) = 0.1, GSFC0 and LANA0, we observe that the standard deviations increase their values with respect to the cases with τ a (440) = 0.3. For instance, the standard deviations of r V f were 0.015 and 0.018 (in Table 7) for GSFC1 and LANA1 respectively, representing 10% in both cases. Meanwhile, for GSFC0 and LANA0 the standard deviations are 0.028 and 0.040, which represents 25% of the averaged values of r V f . We find similar increments, from 10 to 25 %, for C V f . For the coarse mode, the increment of the uncertainties is not as strong as in the fine mode though it is significant in the three parameters.
On the contrary, there is a decrease in the standard deviations of all the parameters for the cases SOLV4 and ZAMB4 τ a (440) = 1.5 with respect to the cases SOLV3 and ZAMB3 τ a (440) = 0.9. For example, the standard deviation of r V f for ZAMB3 was 0.005 (around 3.5% of the total value), whereas it is 0.003 for ZAMB4 (2% of the total value). Similar results are obtained in the characterisation for the coarse mode and SOLV. Thus, the standard deviation of r V c was 0.07 for SOLV3 (3% of the total value) and 0.045 for SOLV4 (a bit under 2% of the total value).
Thus the capacity of the GRASP-AOD algorithm to discriminate between fine-and coarsemode extinction is robust, quite stable in absolute terms and has a constant uncertainty linked to the value of the measurement error that does not depend upon the aerosol load. Consequently, the value of τ f can be accurately retrieved irrespective of the value of τ a . On the contrary, the uncertainty in the bimodal log-normal size distribution parameters dramatically increases as the aerosol load decreases. Therefore, we recommend a lower limit of τ a (440) = 0.2 for the bimodal size distribution parameters in order to assure quality retrievals. We emphasise, however, that this lower limit does not apply to the τ f retrieval.

Reducing the spectral range-
The different sensitivity studies presented in this section have been done for the spectral range of 340 to 1640 nm using the AERONET extended sun-photometer channels (for more information see Appendix A). However, some instruments have a reduced spectral range that does not capture the ultra-violet or nearinfrared part of this range (Table A1 in the Appendix). For instance, many sites in the AERONET network were equipped with the Cimel model 318 Polarised photometer (mainly in Africa and Europe), which only includes the 440, 670, 870 and 1020 nm channels. Likewise, measurements in the GAW-PFR network utilise a similar spectral range from 368 to 862 nm. Hence, we repeated the 1000 inversions for each aerosol case using only the wavelengths in the reduced spectral range of 440 to 1020 nm.
The increments of the average values and the standard deviations obtained in the simulation error study for the reduced spectral range (440 to 1020 nm) with respect to the reference study (spectral range 340-1640 nm; see values in Table 7) are represented in Table 9 for the six bimodal log-normal size distribution parameters and τ f (500). In this case, we have not included the ratio of the inversions that pass the filters (now with the maximum absolute residual fixed to 0.01) since it does not change substantially with respect to the one in Table  7.
The differences in the average values in Table 9 are typically quite small and they do not show a clear tendency. Therefore, the discrepancies between the average values obtained in this study and those obtained in the self-consistency analysis (Table 3) are quite similar to those obtained with the full spectral range. However, differences are more notable in some specific examples. For instance in the case of MEXI1, r V f was 0.133 μm in the selfconsistency analysis and an averaged value of 0.131 μm was presented in the analysis with all the wavelengths. Here, in the analysis with the reduced spectral range, the averaged value is 0.145 μm; which represents a difference of 0.012 with respect to the value obtained in the self-consistency analysis and a difference of 0.014 (shown in Table 9) with respect to the analysis with the full spectral range.
If we first observe the standard deviations in Table 9, we see that they are positive or negative depending upon the parameter. The standard deviation decreases for σ V f , σ V c and r V c and it increases for the rest of the parameters. The increment is particularly relevant in the case of τ f (500): the standard deviations increase their value between 0.007 and 0.013 with respect to the values obtained in the tests with all the wavelengths (Table 7). In those tests, the standard deviations were around 0.008-0.009, while here, the standard deviation of τ f (500) exceeds 0.02 for some aerosol cases. In relative terms, the increments of the uncertainty in τ f (500), by reducing the spectral range, increases from 5 to 10 % for the cases with τ a (440) = 0.3 and from 3-4 to 5-7 % for the cases with τ a (440) = 0.9.
For r V f , the standard deviations in the cases with τ a (440) = 0.3 was around 0.015. This value increases on average by 0.01 when reducing the spectral range. The uncertainty of r V f was only 0.005 for the cases with τ a (440) = 0.9 and for the aerosol models dominated by the fine mode (Table 7). With the reduction of the spectral range this uncertainty only rises by 0.003, given a maximum value in the standard deviation of 0.01.

Sensitivity to the refractive index
As it was already discussed in the introduction, the information contained exclusively in the spectral aerosol optical depth measurements is not enough to retrieve the refractive indices.
We have assumed this parameter as known, thus far, and at this point we wish to determine how errors in this assumption affect the characterisation of the size distribution. Thus, we use the same scheme presented in Fig. 3 to answer this question, but this time we modify the values of the refractive indices.
We only analyse the effects on three of the aerosol cases here: GSFC and SOLV represent the cases with fine-or coarse-mode predominance respectively, and ZAMB is also included to review the effects on an absorbing aerosol.
3.5.1 Pre-analysis with the forward code-First, we use the forward code to check the variation in the aerosol optical depth generated by a modification in the refractive index.
Although the retrievals for three aerosol loads were evaluated in the main study, we only evaluate variations for GSFC2, SOLV2 and ZAMB2 (cases with τ a (440) = 0.6) in this portion of the analysis.
Thus, in Table 10, the change in the aerosol optical depth caused by a change of ±2σ in the refractive index (in both parameters n and k) can be determined for the three examples considered. The values of 2σ are also obtained from the climatology study of Dubovik et al. (2002a) (see Table 1 in that work). In the cases GSFC2 and ZAMB2, 2σ n is equal to 0.02, whereas it is 0.06 in the case of SOLV2. For the imaginary part, 2σ k has values of 0.008, 0.006 and 0.002 for ZAMB2, GSFC2 and SOLV2 respectively. In all cases, the values of 2σ are spectrally independent. Note that for SOLV2 and ZAMB2, k −2σ k gives a negative value of the absorption in some wavelengths. In those cases, k is fixed to zero and indicated with an asterisk in Table 10.
An increase in the real part of the refractive index produces an increase in τ a , because more light is scattered, and consequently, the direct beam is reduced (Bohren and Huffman, 1983). This change is symmetric as can be seen in Table 10. Moreover, ∆τ a is linear with the increments of n considered in this study (the same results were obtained in Torres (2012) where a more detailed explanation can be found).
Another interesting result is the larger effect produced by an increment of n for GSFC2 and ZAMB2 compared to SOLV2, although 2σ n is three times larger for SOLV2, the increments of the aerosol optical depth are similar for the three examples (at 440 nm the variations are 0.046 and 0.034 for GSFC and ZAMB, while it is 0.039 for Solar Village). The explanation lies in the fact that most of the increment produced by the real part of the refractive index increments the optical depth of the fine mode without practically changing the optical depth of the coarse mode. Only for SOLV2, is there a slight decrease in the optical depth of the coarse mode at the largest wavelengths. This property also explains the strong spectral dependence in Table 10: the fine-mode optical depth, τ f , is much larger at short wavelengths than for long wavelengths.
Finally, it can be observed in Table 10 that the variation in the imaginary part does not have a big influence in the aerosol optical depth. The maximum variations, always under 0.01, were obtained for the case ZAMB, which is the case with the largest 2σ k . As expected, positive increments of the imaginary refractive index increase the absorption and therefore, the aerosol optical depth. In this case, the variation is not only allocated to the fine mode, but it is distributed between the modes. Nevertheless, the predominant mode gets the largest variations for each particular aerosol example. Fig. 3, GRASP-AOD code has been applied to the optical depth values for the aerosol cases GSFC, ZAMB and SOLV (Table 2), considering variations of ±σ and ±2σ (see Dubovik et al., 2002a) from the original values of their refractive index (in Table 1, only the cases with τ a (440) = 0.3, 0.6 and 0.9) during the inversion procedure.   11 presents the retrieved values of the fine-mode aerosol optical depth at 500 nm for variations of ±2σ in the real and in the complex part of the refractive index. The maximum difference between the retrieved values with and without the variations is observed for the case SOLV3 and is equal to 0.005; however, the differences are typically around 0.001-0.002. A similar result is obtained for the rest of the wavelengths, although it is not presented here. Therefore, the discrimination between the extinction of fine and coarse modes does not depend upon the refractive index assumption. Figure 9 shows the size distributions obtained by varying the real and the imaginary parts of the refractive index by ±σ and ±2σ during the retrieval process, for the three aerosol cases (for each of them with the three different aerosol load τ a (440) = 0.3, 0.6 and 0.9). In the figure, red lines are used to indicate the cases when the refractive index is overestimated and blue lines represent the cases when the parameter is underestimated. Variations of ±σ are depicted by solid lines and the ones of ±2σ by dashed lines. Finally, black solid lines illustrate the reference values, while grey dashed lines represent the retrievals without variations in the refractive index (self-consistency study).

Retrieval analysis-Following the scheme from
The retrieval behaviour of the size distribution presents a similar pattern for the three analysed cases in response to the variations of the refractive index: there is a decrease in the mean radius and the volume concentration of the modes when the real part of the refractive index increases, while for negative variations both parameters increase their values. The same result was obtained in King et al. (1978), where it was pointed out that the shape of the size distribution remains the same but shifts with a varying real part. This is derived from the anomalous diffraction theory by Van de Hulst (Van de Hulst, 1957) and was also discussed in Yamamoto and Tanaka (1969).
Analysing the cases with fine-mode predominance, GSFC and ZAMB, we observe that this tendency is only significant in the fine mode: the increments of r V f are between ∓(4-5) % and between ∓(2-4) % for C V f for variations of ±2σ n . For the coarse mode, the variations in r V c and C V c do not show clear behaviour for GSFC and ZAMB, with values practically negligible (typically under 1 %). These results are expected given the results from the preanalysis with the forward code in the previous section and the stability of τ f observed in Table 11. As mentioned earlier, an increase in the refractive index results in a rise in the aerosol optical depth in the forward model. Since τ f does not change in the retrieval, we expect the increase in the real refractive index to be compensated elsewhere. Figure 9 indicates that increases in the real refractive index are balanced by a decrease in the volume of the particles and also by a reduction of r V f . Note that the extinction efficiency diminishes as the radius decreases in the fine mode, especially for the shortest wavelengths. Figure 9 seems to indicate linear behaviour in the increments of r V f and C V f against the variation of the refractive index. The differences for ±σ n are half as large as for ±2σ n for both GSFC and ZAMB. We have explored whether this behaviour continues for larger variations of the refractive index, in particular, for the case of GSFC. Note that the conditions at the GSFC site are favourable for aerosol hygroscopic growth and the refractive index can suffer large variations with respect to the averaged values. Thus, we introduced differences up to ±5σ n (or 0.05) for the three aerosol loads of GSFC (τ a (440) = 0.3, 0.6 and 0.9). The results of these tests have shown that the linear behaviour continues with these variations of the refractive index and that the differences in r V f and C V f can be approximated as Δ r V f − 0.04 × Δ n, Δ C V f − 0.27 × τ a 440 × Δ n. In relative terms, these differences in r V f represent between ∓(12-13) %, but in C V f they represent between ∓(8-10) % for the maximum variation considered of ±5σ n . The large refractive index variation of ±5σ n also perturbs the coarse mode r V c by ∓(3-5) % and C V c by ∓(2-4) %. The stability of τ f observed in Table 11 is still valid even with a variation of ±5σ n , with a maximum difference of ∆τ f = 0.005 observed for the case GSFC3.
For the case of SOLV, the general tendency is observed for both fine and coarse mode even with a variation of just ±2σ n in the refractive index. In the fine mode, the differences in r V f and C V f are between ∓(8-9) % and ∓(7-8) % respectively. For the coarse mode, the mean radius and the volume concentration decrease between (3-4) % (for +2σ n ) when the real part of the refractive index increases.
The effective radius follows the same tendency as the mean radii of the modes. For the maximum increment of ±2σ n , it varies as ∓0.01 μm in the three aerosol loads of GSFC and ∓0.006μm in those of ZAMB. These changes represent the ∓5 % and the ∓3 % respectively of the total value. In the examples of SOLV, the variations are of ∓3 % for the two lowest aerosol loads and approximatively ∓7 % for SOLV3.
The three subfigures on the right in Fig. 9 depict the retrievals that vary the imaginary part of the refractive index. The size distributions show similar patterns for the three aerosol cases: mean radii of the two modes do not have large variations with difference under 2% in the case of r V f , being a bit higher for r V c around (3-4) %. For both modes, these differences do not present clear behaviour and are randomly positive or negative regardless of the sign in the variation of the refractive index. By analysing the concentration, however, a definite tendency can be observed: they decrease for positive variations of the refractive index and they increase for negative increments. The differences are around 1% for C V f and between (3-4) % for C V c . This result is expected, since when the absorption is increased there is an increase in the aerosol depth that is compensated by a reduction in the number of the particles. The behaviour in the concentration, however, is not so clear in Fig. 9, where it seems that the opposite happens. The explanation is in the standard deviation of the modes, which also diminish for positive increments of the absorption and increases in the opposite case. Thus the maximum values in the size distribution are reached in the cases of +2σ k , though the aerosol concentration are the lowest in these cases. The maximum variations of σ V f and σ V c are between 0.015 and 0.022 (absolute values), independently of the aerosol type. The variations in the effective radius are negligible for GSFC and ZAMB (differences under 0.5 %). In the case of SOLV, differences up to 2 % are observed, though there is no clear tendency in them. This general little sensitivity for the variation of the imaginary part was also obtained in the studies of King et al. (1978) and Yamamoto and Tanaka (1969).

Variations in the sphericity parameter
We analyse the sphericity parameter assumption for our last sensitivity test. This parameter was introduced in Eq.
(3) and defines the percentage of spheres in the aerosol model. The kernels of the extinction are quite similar for spheres and spheroids (unlike the kernels of angular scattering), but some small differences do appear because spheroid cross sections can not be straightforwardly calculated from equivalent sphere radii. It should also be noted that the sensitivity of extinction to particle shape is minor for small particles at the size parameters used in the retrieval, as illustrated in Fig. 23 of Dubovik et al. (2006) (even the P 11 function shows little sensitivity as shown in Fig. 27 of Dubovik et al., 2006). Therefore, the main discrepancies produced by the uncertainty in the sphericity parameter should appear in the coarse mode.
We have selected the aerosol models of SOLV1 and BAHR3 in this analysis, since their values of the sphericity parameter is 0.6 and we can largely decrease and increase it in our tests. At the same time, the use of SOLV1 (τ a (440) = 0.3) and BAHR3 (τ a (440) = 0.9) allows us to see whether the effect is more significant depending on the aerosol load.
We begin with a brief comment on the variations in the aerosol optical depth caused by a modification in the sphericity parameter. We use the forward model to increase the sphericity from 0.6 to 1.0 and decrease it from 0.6 to 0.2 for both SOLV1 and BAHR3. The τ a variations are quite similar in relative terms for both cases -increasing the sphericity parameter from 0.6 to 1 causes the aerosol optical depth to decrease by less than 1 % at wavelengths less than 500 nm, to decrease by (2-3) % at the 675 and 870 nm wavelengths and to decrease by about 5 % at the 1020 and 1640 nm wavelengths. We also analyse the modes separately and find that the coarse-mode optical depth variations are negative and almost constant at around 5 % in both aerosol cases. The fine-mode variations, on the other hand, are positive with values between and (1-2) % and a maximum at the central wavelength. The opposite results but with the same relative increments are obtained when the sphericity parameter is reduced.
We modified the reference value of the sphericity parameter for SOLV1 and BAHR3 during the inversion procedure as a second step. We reduced the sphericity to 0.4 and 0.2 and incremented it to 0.8 and 1.0, and we present the variations in the resulting retrieved parameters in Table 12. As we expected from our tests with the forward mode, the differences for the fine mode are very small -the differences for r V f and C V f are negligible, and σ V f only changes at maximum by ±0.02. Similarly, no changes for the retrieved values of τ f (500) are observed in Table 12.
On the contrary, we observe some significant differences for the parameters that define the coarse mode. These differences are more significant for SOLV1 than for BAHR3. For instance, the maximum difference for r V c is 0.145 μm for SOLV1 (6 % in relative terms) and 0.056 μm for BAHR3 (about 2 %). Differences in C V c are also more significant for SOLV1 than BAHR3, with maximum differences of 0.014 μm 3 μm −2 (about 10 % of total value) for SOLV1 and only 0.004 (1 % in relative terms) for BAHR3. As mentioned earlier, increasing the number of spheroids decreases the sphericity parameter and increases the optical depth. Since the retrieval is constrained by optical depth, the retrieval compensates for the reduction of the sphericity assumption at BAHR3 by increasing r V c (recall that the extinction efficiency diminishes as the radius increases in the coarse mode). For SOLV1 the relation is less clear, since C V c and σ V c also suffer large variations.
Nevertheless, the variations try to compensate for the change in τ a that occurs when the spheroid ratio is altered.

Comparison with AERONET
We present here a comparison between the aerosol products obtained from the new GRASP-AOD inversion and AERONET products in order to validate the results achieved in previous sections. In particular, we compare our results with (a) size distributions obtained from the AERONET inversion code Dubovik et al., 2006), and (b) τ f (500) obtained from the SDA (O'Neill et al., 2003). All the data and products used in the comparison belong to AERONET level 2.0 (quality-assured data; see Smirnov et al., 2000) and can be found in the public AERONET database (http://aeronet.gsfc.nasa.gov).

Data selection-
To homogenise the two different data sets (almucantar and spectral aerosol optical depth measurements), we propose a comparison of daily averages on days with stable aerosol conditions instead of comparing single inversions. This approach allows us to include a large number of τ a inversions which will be useful for checking the stability of the new inversion. We limit our analysis to days that fulfil three requirements:

1.
There are a minimum of 15 τ a measurements per day; each τ a measurement is only eligible if there is a valid almucantar within a maximum of ±30 min delay. This requirement also improves the quality of the τ a data selected, as there are only almucantars at level 2.0 with θ s > 50°. That is, a maximum delay of 30 min ensures that all τ a measurements are obtained when θ s > 40°, which reduces the bias produced by errors in the calibration coefficients (see Eq. A6).

2.
In order to guarantee stable aerosol conditions, the ratio between the standard deviations and the averages of the eligible aerosol optical depths are required to be smaller than 0.1 throughout the day (evaluated for every spectral wavelengths used by the GRASP-AOD inversion).

3.
There should be a minimum of four valid almucantars per day.
Following these conditions, 29 days have been selected from the same six sites used as (climatological) inputs in the simulations and two additional Mediterranean AERONET sites: Rome Tor Vergata (henceforth just Rome) and the ChArMEx super site of Lampedusa. Following Mallet et al. (2013), the dominant aerosol type is background aerosols with frequent dust intrusions for Lampedusa (coarse-mode predominance), while for Rome, the dominant aerosol type is urban-industrialised (fine-mode predominance).
In total, 165 AERONET inversions (almucantar and τ a ) and 744 GRASP-AOD inversions have been compared. Table 13 depicts the information regarding the daily conditions of the analysed cases. The first three columns in Table 13 contain technical information about the day selected (site, photometer number, type and date). Mexico City, Mongu, Lanai, Rome and Lampedusa do not have extended photometers for the measurements, so we have used only standard wavelengths for the analysis at those sites. This is indicated in the second column with the label "STD" (for the standard photometers) or "EXT" (for the extended photometers) together with the AERONET number of the instrument (see Table A1 to get more information about standard and extended photometers). The selected days for each site correspond to the same deployment period (between a common pre-and post-calibration). For Rome and Lampedusa, the chosen interval is coincident with the ChArMEx/ADRIMED (Chemistry-Aerosol Mediterranean Experiment/Aerosol Direct Radiative Forcing on the Mediterranean Climate) summer campaign from 11 June to 5 July 2013 over the western Mediterranean . Unfortunately, at Lampedusa site, the 340 nm channel has not passed the requirements of AERONET level 2.0 during this period, so the data at this wavelength have not been used for the GRASP-AOD inversions at this site.
The daily averaged value and the standard deviation of the aerosol optical depth at 440 nm as well as the Ångström exponent (from 440, 675 and 870 nm) are indicated in columns four and five of Table 13, so as to give an idea of the characteristics of the aerosol analysed: load and predominance of fine or coarse mode. As can be seen from the table, we covered a broad aerosol load with τ a (440) variations by a factor of 2 (LANA) to 3.5 (GSFC) using three or four examples for each site. The last two columns indicate the number of inversions: ALM indicates the number of almucantars used for the AERONET standard inversion, AOD the number of inversions used for both GRASP-AOD inversion and SDA algorithm.
To run the GRASP-AOD inversions, we have assumed climatological values of the refractive index for the different sites in Table 13. Thus, for the first six cases we have taken the values from Table 1, while for Rome and Lampedusa we relied on Table 2 of Mallet et al. (2013) 2 . For the sphericity parameter, we have used values of 0 and 100 for Ångström exponents smaller than 0.6 and larger than 1.1 (respectively), and we have applied a linear interpolation (with respect to the Ångström exponent) to intermediate values. Table 14 represents the daily average of the standard parameters for fine and coarse modes retrieved using AERONET-standard and GRASP-AOD inversion for the chosen days shown in Table 13. While both r V f and r V c are direct outputs of the GRASP inversions, the standard AERONET inversions gives a 22-bin retrieved size distribution. In the latter situation, values of r V f and r V c are estimated using the standard AERONET procedure described in http://aeronet.gsfc.nasa.gov/new_web/Documents/ Inversion_products_V2.pdf. Figure 10 illustrates the absolute differences of the mean volume radii by using GRASP-AOD and AERONET standard inversion calculated from the values in Table 14. In the left panel, the differences for the fine mode are represented, while those of the coarse mode are depicted in the right panel. By analysing the differences in r V f at the four sites with finemode predominance (GSFC, MEXI, ZAMB and ROME), we observe an excellent agreement between the inversions, with differences under 0.015 μm except for the cases of GSFC -A-(maximum difference of 0.028 μm), GSFC -D-(0.018 μm) and ZAMB -A-(0.020 μm). At the four sites with a prevail-ing coarse mode (SOLV, BAHR, LANA and LAMP), the differences are a bit higher, generally ranging from 0.010 to 0.020 μm. The maximum difference of 0.035 μm is observed for the case LAMP -B-. Thus, the averages of the differences of r V f are 0.010 μm at sites dominated by the fine mode and 0.013 μm when we include all the sites.

Results-
The differences in r V c are under 0.3 μm for the cases with τ a (440) > 0.3. The maximum differences are reached in the cases with a prevailing fine mode and τ a (440) < 0.3: MEXI -A-and ZAMB -A-with values of 0.46 μm and 0.48 μm respectively. The analysis in relative terms gives similar results to that in the fine mode, with differences of around 15 % in our particular test cases, although they do not exceed 10 % in the rest of the cases.
In the case of C V f , the largest differences are obtained for the cases with a predominant fine mode, which have a maximum value of 0.016 μm 3 μm −2 found for MEXI -D-(see Table 14). For C V c , the largest differences are found for the desert dust cases, which have a maximum value of 0.043 μm 3 μm −2 (SOLV -D-). The high variability of the concentration for both modes makes the analysis more suitable in terms of relative differences. Figure 11 shows the relative differences (in modulus) of the daily volume concentration retrieved by GRASP-AOD and AERONET standard inversion for the different examples. The average differences in both C V f and C V c are around 17 %. As a general result, we observe that the differences for C V f are the smallest for the cases dominated by the fine mode, and the differences for C V c are the smallest for the cases dominated by the coarse mode. We find that the average differences in volume concentration are around 10 % for the dominant mode, whether it is the fine or coarse mode.
In Fig. 12, the absolute differences for the daily means of σ V f and σ V c , retrieved by GRASP-AOD and AERONET standard inversions (Table 14), are shown. The average difference is 0.071 for σ V f and 0.081 for σ V c , although differences larger than 0.2 are sometimes observed in both cases. These results confirm one of the outcomes from the sensitivity analysis, where we pointed out that the retrieval of the standard deviations (mode widths) is less accurate than the retrieval of the modal radii and concentrations.
Values of τ f (500) are computed using the retrieved parameters from standard AERONET 3 and GRASP-AOD inversions. Both values are presented together with those obtained by the SDA algorithm in Table 15. It can be observed that, in general, the values of τ f (500) retrieved by GRASP-AOD range between those retrieved by the standard AERONET inversion and the SDA algorithm. By analysing the differences for the cases with a prevailing fine mode, we see that a maximum value of 0.020 is found for the case GSFC -Dbetween GRASP-AOD and the standard AERONET inversion. On the other hand, a maximum difference of 0.011 is observed between GRASP-AOD and SDA algorithm for the case GSFC -D-. The analysis of the cases dominated by the coarse mode shows differences of up to 0.034 (SOLV -D-) between GRASP-AOD and the standard AERONET inversion. The maximum difference is a bit larger for the comparison between GRASP-AOD and SDA algorithms with a value of 0.036 (BAHR -A-).
Note that computing τ f (500) also allows us to compare the results obtained by the SDA and Table 15 also contains the values of r eff respectively retrieved by GRASP-AOD and standard AERONET standard inversion. For the cases with a prevailing fine mode, the largest differences are obtained for the cases of MEXI and ROME. Among them, the largest value is 0.06 μm, which is reached in the case MEXI -B-. The differences are higher in the cases with a predominance of the coarse mode. The extreme case, both in absolute and relative terms, is found for LAMP -B-with an absolute difference of 0.196 μm, which is around 30 % in relative terms. The average of the differences is 10 % and it does not depend on the mode predominance.

Night AOD inversion from moon photometers
The recent advances in moon photometry (in terms of instrumentation, calibration and data treatment) have enabled the set-up of night measurements at a few sites. One of the pioneer sites at Lille University has been providing regular night τ a measurements since 2013 from a moon photometer (Cimel model CE318-U). The typical cloudy conditions of the site and the constraints due to the moon cycle (only 14 days out of 28 valuable for night measurements) limit the study cases to a relatively small number of nights. Note that the moon measurements have a spectral range from 440 to 1640 nm (see Table A1), while the sensitivity analyses of Sect. 3 were mainly done for the 340-1640 nm spectral range.
The left panel of Fig. 13 shows the aerosol optical depth measured during one of those nights, specifically on 9 April 2015. In the same figure, we have added the τ a values measured by the sun photometer #741 (AERONET number) during the previous evening and the following morning. Even though there is only 1 h of valid measurements (a total of six measurements), the case represents an interesting example with relatively high τ a values due to the pollution episode that occurred in the north of France in the spring of 2015 on 7-9 April. During those days, the values of τ a (440) typically ranged between 0.4 and 0.6 (with peaks up to 0.9 on 8 April), while the climatological average at the Lille site is < τ a (440)>=0.22.
In the right panel of Fig. 13, size distributions retrieved by GRASP-AOD inversion for the night τ a measurements are represented (six inversions in red solid lines). In the same figure, we have added the GRASP-AOD retrievals from the two τ a measurements of the previous evening (blue lines) and from one of the τ a measurements the next morning (dashed black line). The τ a measurement selected is the one coincident with the almucantar measurement, and together they form the input of the first AERONET-standard inversion of the day and are closest to the night measurement. Note that the refractive indices retrieved from this inversion have been used in the GRASP-AOD retrievals shown in Fig. 13 (day and night). The 22-bin size distribution from this inversion is also illustrated with a black solid line.
Observing the AERONET retrieval, we can see that the episode is characterised by a predominant fine mode C V f /C V c = 4 with relatively high values of r V f = 0.320 μ m and r eff = 0.317 μm compared to the typical urban aerosols in Lille (r eff < 0.2 μm see Mortier, 2013).
The average of the retrieved values during the night shows similar tendencies in the fine mode with r V f = 0.304 ± 0.01 μ m and r eff = 0.296 ± 0.012 μm. The ratio between the concentrations, however, is a bit higher C V f /C V c = 6.94 ± 1.45 than the one found in the AERONET retrieval. This difference is mainly due to the discrepancies in the coarse mode between the retrievals: the C V c retrieved during the night is 0.009 ± 0.002 μm 3 μm −2 , but in the AERONET retrieval it is 0.015 μm 3 μm −2 . The concentration values of the fine mode are quite similar: C V f = 0.057 ± 0.03 μm 3 μm −2 (GRASP-AOD) and C V f = 0.06 μm 3 μm −2 (AERONET). Finally, the τ f (500) = 0.360 ± 0.004 retrieved during the night from GRASP-AOD is in excellent agreement with the value retrieved during the next morning by SDA, τ f (500) = 0.337, and the one derived by the almucantar inversion, τ f (500) = 0.349.

Application of the method to the airborne PLASMA sun photometer
The spectral optical depth measurements obtained by airborne sun photometers represent another interesting application example of the GRASP-AOD inversion. Typically the use of airborne sun photometers is limited to field campaigns developed in singular areas with specific research objectives such as validation of models, characterisation of specific aerosol types or interactions between aerosols and other atmospheric components like clouds or gases. From the τ a measurements taken at different heights, the derivation of the spectral aerosol extinction profile is immediate. The extinction profiles are normally available for every landing/taking off and during pre-scheduled vertical profiles carried out by the airplanes. They give basic information on the characterisation of the aerosol vertical distribution.
In the example presented here, we use the data from the PLASMA airborne sun photometer during the ChArMEx/ADRIMED summer campaign that took place from 11 June to 5 July 2013 over the western Mediterranean . The PLASMA airborne sun photometer was installed in the ATR-42 French research aircraft, which was operated from Sardinia (Italy) and participated in the 18 flights of the campaign. The experimental set-up also involved several ground-based measurement sites on different Mediterranean islands (Corsica, Lampedusa, Minorca and Sicily) and additional measurements from lidar and sun photometers performed on alert during aircraft operations (at Granada and Barcelona). The meteorological conditions observed during the campaign (moderate temperatures and southern flows) were not favourable to producing large concentrations of local polluted smoke particles. However, several moderate dust plumes were observed during the campaign, with the main sources located in the north-west of the Sahara desert. Though peaks in τ a of up to 0.6 (at 440 nm) were registered by ground-based sun photometers during the campaign, the maximum values measured by collocated sun photometers at the time of ATR-42 vertical profiles were around 0.3 (at 440 nm).
The τ a values represented in the left part of Fig. 14 correspond to PLASMA measurements at different heights during one of the aforementioned dust plumes, specifically in the take-off from Minorca airport on the 17 June 2013 between 11:45 and 12:00 UTC (flight 32; Denjean et al., 2016). The average values of the three τ a data measured during this period by the ground-based AERONET sun photometer at Minorca (AERONET site Cap d'en Font; see Chazette et al., 2016) are also plotted as reference. Note that the ground-based sun photometer is standard (see description in Table A1) and does not have the 1640 nm channel. At the same time, the 340 nm channel of PLASMA suffered from instabilities during the campaign and the derived τ a data at this channel are not represented.
One of the first conclusions obtained from Fig. 14 is the vertical homogeneity of the aerosol characteristics in the column up to 4000 m. The same property was observed by analysing the coincident lidar data (Fig. 6 in the study by Chazette et al. (2016) on 17 June 2013). τ a diminishes with height, but the spectral shape of the measurements does not considerably  (870)) varies between 0.46 at 500 m and 0.61 at 4000 m. The Ångström exponent from AERONET sun-photometer reference at ground level is 0.48. This same property is visible in the size distribution obtained by GRASP-AOD inversions and represented in the right part of Fig. 14: neither the mean volume radius nor the standard deviation of the modes change significantly from the retrievals at different heights. For r V f , the values span from 0.120 μm at ground level and 500 m to 0.135 μm at 4000 m, when r V f given by AERONET standard inversion is 0.134 μm. In the case of r V c , the variation goes from 2.06 μm at 4000 m to 2.21 μm at ground level, while r V f given by AERONET is 2.27 μm. The ratio between the

Conclusions
The main goal of the present work was to show the potential of retrieving the total column aerosol size distributions from spectral optical depth measurements without the aid of coincident radiance measurements, and estimating a set of secondary aerosol properties (e.g. effective radius or fine-mode fraction of aerosol optical depth) derived from it. The limited information content in spectral τ a measurements results in the necessity of using a priori constraints. The utilisation of the GRASP public software allowed us to test and evaluate different aerosol model descriptions and a priori constraints. The current analysis indicates that bimodal log-normal size distributions and a priori estimates of refractive indices and sphericity parameter provide a practically efficient retrieval. The validation of the retrieval has been done through (a) a sensitivity analysis using six different aerosol models with three or four aerosol loads for each of them (Sect. 3) and (b) a comparison of the aerosol properties obtained from 744 AERONET observations using GRASP-AOD inversion to those obtained at eight different AERONET sites through 165 almucantar AERONET standard inversion.
The simulated tests have shown that spectral τ a measurements are sufficient for a precise discrimination between the extinction of fine and coarse modes, independent from any assumption, with maximum differences in τ f (500) under 0.01 from the input values found in all the sensitivity tests. Specifically, the characterisations of aerosol fine-mode optical properties are accurate, although they depend on reliable a priori information about the real refractive indices and accurate measurement of aerosol optical depths. The uncertainty observed during the sensitivity tests for the fine mode r V f and C V f is between 5 % for the cases with a fine-mode predominance and 10 % for the cases with a prevailing coarse mode.

NASA Author Manuscript
The characterisation of the optical properties of the coarse mode using τ a measurements is less accurate but can be significantly improved using moderate a priori information about coarse-mode parameters (for example, using as initial guess data retrieved from nearalmucantar inversions or climatological values from the given site). The study showed that good calibration of the long wavelengths (1020 and 1640 nm) is essential due to its influence on the coarse mode and the typical low aerosol optical depth values presented in this spectral range. Nevertheless, in the cases in which the coarse mode is predominant, the uncertainty observed in our sensitivity studies is 10 % for r V c , and around 20 % for C V c . The effective radius has been revealed as a quite stable parameter relative to the calibration errors (uncertainties around 15 %) and bias in the assumed refractive index (differences under 10 % for the analysed cases). The simulated tests including aerosol optical errors showed that the uncertainty of the bimodal log-normal size distribution parameters dramatically increases as the aerosol load decreases. We recommend a lower limit of τ a (440) = 0.2 for the bimodal size distribution parameters in order to assure quality retrievals of aerosol particle bimodal size distribution parameters. However, this lower limit does not apply to the τ f retrieval.
Daily averaged aerosol properties obtained by applying GRASP-AOD inversion to 744 AERONET observations were compared to those retrieved by 165 almucantar AERONET standard inversion at eight selected AERONET sites with varied aerosol characteristics. The retrieved values of τ f (500) using GRASP-AOD are generally between those obtained by AERONET standard inversion and the spectral deconvolution algorithm, with differences typically lower than 0.02 between GRASP-AOD and both algorithms. In some cases with large aerosol load in the coarse mode, maximum discrepancies of up to 0.05 are observed between the three retrieval methods. The study on real cases confirms the retrieval capabilities in the fine mode: maximum differences of 0.028 μm in r V f between GRASP-AOD and AERONET inversions are observed, though they are generally lower than 0.015 μm (10% of the values). In the cases with fine-mode predominance, the differences in C V f are typically under 0.01 μm 3 μm −2 , and in relative terms, the average of these differences is around 10 %. The comparison of coarse-mode parameters shows larger differences, however, especially at the sites with a prevailing fine mode. In these cases differences of 0.2-0.3 μm in r V c are normally obtained, although some extreme differences of up to 0.5 μm are also observed. The comparison improves significantly when limiting the study to those sites characterised by the predominance of desert dust aerosol, with differences typically under 0.2 μm between GRASP-AOD and AERONET. In relative terms, these differences are under 25 % for all the dust cases, which agrees with the uncertainty estimated from the sensitivity tests. The effective radius shows average differences around 10 % and it does not depend on the fine-or coarse-mode predominance.
Finally, we emphasise that the use of GRASP-AOD code is not restricted to the standard AERONET sun-photometer measurements. In the last part of Sect. 4, we also show two practical applications of the GRASP-AOD code: night-time retrievals using moon photometer data and vertically resolved retrievals using the PLASMA airborne sun-tracking photometer.

Data availability.
The AOD data used for the sensitivity studies in Sect. 3 were generated by the aerosol properties depicted in Table 1. The simulated aerosol optical depth values used in the main analyses are shown in Table 2. The AOD values used in Sect. 4 "Real Cases" are publicly available on the AERONET web page (https://aeronet.gsfc.nasa.gov/). GRASP inversion algorithm software used in this work is free and publicly available at http://www.graspopen.com. Moreover, a web application was developed to directly use GRASP-AOD inversion without needing to install all of the GRASP code. The application can be found at www.grasp-open.com/aod-inversion.
Appendix A:: Aerosol optical depth measurements.

A1 General background
Total optical depth of the atmosphere (τ) can be understood as the attenuation of light passing through the atmosphere column containing aerosol particles, molecules and absorbing gases, and it is described by the well-known Beer-Bouguer-Lambert Law: where F(λ) (Wm −2 μm −1 ) and F 0 are the monochromatic direct flux densities at the Earth's surface and at the upper limit of the atmosphere respectively. For atmospheric applications where the sun is the radiation source, F(λ) and F 0 are defined as solar flux densities. We intentionally avoid the adjective "solar" here, since the flux density comes from the moon for the night observations presented in this work. Additionally, optical depth can also be derived from star radiation as well. Finally, the optical air mass (m s ) accounts for the light path in the atmosphere and can be approximated 4 by m s = 1/ cos θ s , where θ s is the zenith angle of the radiation source in the sky.
Under cloud-free conditions, the total optical depth can be separated into the gaseous absorption τ g , the molecular scattering or Rayleigh scattering τ R , and the aerosol scattering and absorption τ a , which is known as the aerosol optical depth. The latter can be derived, therefore, from the following expression: Instruments designed to derive the aerosol optical depth typically utilise wavelengths that lack significant gas absorption in the spectral region from the ultraviolet (UV) to the near infrared (NIR). This spectral region presents the highest sensitivity regarding scattering and extinction by aerosols according to the typical sizes of the natural occurring aerosols (Shaw et al., 1973;Shaw, 1983;Dutton et al., 1994).
The instruments and wavelengths used in the GRASP-AOD inversion tests with real data (Sect. 4) are depicted in Table A1. The main analysis of those tests has been done using data from the two main groups of instruments in AERONET: Cimel-318 standard version (STD) with a spectral range from 340 to 1020 nm, and Cimel-318 extended version (EXT), which includes the 1640 nm channel 5 . The example of application of GRASP-AOD inversion to night measurements has been done with the moon photometer Cimel-318U. Given the low signal of moon radiation in the ultraviolet spectral region, the moon photometer Cimel-318U does not include the 340 and 380 nm channels. The last example is based on the airborne sun-tracking photometer PLASMA (Photomètre Léger Aeroporté pour la Surveillance des Masses d'Air). It has a wider spectral range (340-2250 nm) than AERONET-extended instruments. Nevertheless, the calibration of PLASMA has been done using the same protocols and tools of AERONET, and therefore, only those channels of extended groundbased Cimel photometers were used in GRASP inversions of PLASMA data, with the exception of the 340 nm channel due to instability during the campaign.
To summarise, the sensitivity study presented in Sect. 3 has been mainly done only for the spectral range from 340 to 1640 nm and using the channels of AERONET-extended sun photometer presented in Table A1.
4 Only valid if θ s ≤ 75°. The exact formulation can be found in Kasten and Young (1989).
5 For both instruments, there is an extra spectral band centred at 936 nm which is not used for the description of the aerosol properties. The channel is intentionally selected in an absorption band of water vapour in order to determine the column abundance of this gas.

A2 Calibration and data treatment
Calibration of the instruments is carried out by the Langley plot method (Shaw et al., 1973). The method is based on the Beer-Bouguer-Lambert law and it is directly applicable if the instrument has a linear response. In this case, Eq. (A1) can be transformed to lnV λ = lnV 0 λ − τm s , where V (λ) are the digital counts measured by the instrument and V 0 (λ) is the calibration constant or the extraterrestrial signal of the instrument.
The Langley method is used to derive V 0 (λ) by means of a set of direct sun observations performed over a range of air masses (typically from 2 to 5). The measurements provide a straight line (ln V (λ) vs. m s ), with an intercept (ln V 0 (λ)) from which the calibration constant (V (λ) 0 ) can be extracted.
In moon photometry, calibration presents greater complexity than the common Langley approach for sun photometry. The main difficulty is that the moon is a highly variable source, and the illumination changes continuously with the lunar-viewing geometry. Consequently, a lunar irradiance model needs to be considered. Barreto et al. (2013) developed the Lunar Langley Method (LLM), which is a modification of the usual Langley technique that can be applied to cases with variable illumination conditions. The calibration coefficient V 0 (λ), is variable in time with the LLM and is expressed as a function of the moon's extraterrestrial irradiance (I 0 ) and the instrument calibration constant (κ(λ)): Here, I 0 (λ, t) is taken from the Robotic Lunar Observatory (ROLO) model, developed by Kieffer and Stone (2005). The ROLO model presents a relatively precise I 0 (λ, t), with maximum errors around 1 % at any time and for all of the wavelengths. More information about the calibration process and instrument deployment can be obtained from Barreto et al. (2013) and Barreto et al. (2016). Instrument wavelengths used in the GRASP-AOD inversion tests with real data (Sect. 4).

Instrument type Spectral channels (nm)
Cimel-318 Extended PLASMA 340, 380, 440, 500, 670, 870, 1020, 1640Cimel-318 Standard 340, 380, 440, 500, 670, 870, 1020Cimel-318 Polarized 440, 670, 870, 1020 Cimel-318U (moon photometer) 440,500,670,870,1020,1640 In AERONET, the intercalibration procedure is used for the calibration of field instruments. The method is based on the realisation of simultaneous co-located measurements of a socalled master instrument (calibrated by Langley method) and the field instruments under certain atmospheric conditions. The field instrument can be calibrated simply by a ratio of raw signals of each channel (V field (λ)) with the master raw signal (V master (λ)): The accuracy of the master calibration is about 0.25 % in the visible and NIR and 0.5 % in the UV, whereas for field instruments the calibration uncertainty is 1 % in the visible and NIR and 2 % in the UV, due to uncertainty in the calibration transfer Holben et al., 2006). Deriving Eq. (A3), we obtain the following: This means that, for instance, an uncertainty of 1 % in V 0 represents a maximum error of 0.01 in the total optical depth in the extreme case of m s = 1. This result is taken into the account in the sensitivity tests of Sect. 3.
Finally, all the data used in Sect. 4 are part of the level 2.0 quality-assured AERONET data set  and http://aeronet.gsfc.nasa.gov/new_web/PDF/ AERONETcriteria_final1.pdf). The measurements were done following the standard sequences of AERONET and the final values of τ a obtained through the official data treatment of the network (http://aeronet.gsfc.nasa.gov/new_web/Documents/ version2_table.pdf).    Methodology diagram followed to carry out different sensitivity "tests" on GRASP-AOD code.   Differences between the averaged values obtained in the error simulation study (Table 7) and the retrieved values from the self-consistency analysis (Table 3) for the six parameters that describe the bimodal log-normal size distributions: volume median radius (top subfigures), geometric standard deviation (middle subfigures) and volume concentration (bottom subfigures). The differences for the latter are plotted in relative terms.  Standard deviation (SD) obtained in the error simulation study (Table 7) for the six parameters that describe the bimodal log-normal size distributions: volume median radius (a), geometric standard deviation (b) and volume concentration (c). For the latter, the standard deviations are plotted in relative terms.  The effect of τ a errors on size distributions retrieved with the GRASP-AOD code for the cases of GSFC0 (low aerosol load conditions τ a (440) = 0.1) and SOLV4 (high aerosol load conditions τ a (440) = 1.5). As in Fig. 5   Absolute differences between the daily averages of the median volume radii retrieved by GRASP-AOD and AERONET standard inversions (total values can be found in Table 14) for the examples considered in the analysis as listed in Table 13  Relative differences between the daily averages of the volume concentrations retrieved by GRASP-AOD and AERONET standard inversions (total values can be found in Table 14) for the examples considered in the analysis as listed in Table 13  Absolute differences between the daily averages of the mode standard deviation retrieved by GRASP-AOD and AERONET standard inversions (total values can be found in Table 14) for the examples considered in the analysis as listed in Table 13   Aerosol optical depth values measured by PLASMA airborne photometer at different height levels in the take-off from Minorca airport on the 17 June 2013 (11:45) on the left, and the corresponding particle size distribution retrieved by GRASP-AOD inversion on the right.
The mean values of the three τ a data measured by the sun photometer at AERONET site Cap d'en Font during the PLASMA profile are also represented with its corresponding GRASP-AOD inversion. The size distribution from the closest AERONET standard inversion (09:00) is also represented in the figure as a reference. Description of aerosol properties used for simulating the aerosol optical depth measurements. They are based on the climatology study of Dubovik et al. (2002a). The first part of the table specifies the parameters describing the aerosol load, the size distribution (modelled as a bimodal log-normal function:  Differences obtained from the self-consistency test on GRASP-AOD code obtained by comparing parameters in Tables 1 and 2. The first six columns provide parameters that represent the bimodal log-normal size distributions (volume median radius (r  )). The following columns contain two secondary products derived from the GRASP-AOD inversions and the corresponding reference values: fine-mode aerosol optical depth at 500 nm (τ f (500)) and the effective radius (r eff [μm]). Initial guesses used as defaults for the GRASP-AOD retrieval. The values, as can be seen in the table, are given as a function of Ångström exponent and the aerosol optical depth at 440 nm.  Values used for the initial guess tests. For each aerosol example, there are three possible initial guesses for the six retrieved parameters (volume median radius (r  )).

Model
Model Sensitivity of multiple initial guesses to the GRASP-AOD inversion. The averages and the standard deviations of the retrievals are presented for six primary parameters (volume median radius (r )) in the first six columns (reference values in Table 1). The last two columns contain the averages and the standard deviations of the two secondary products (fine-mode aerosol optical depth at 500 nm (τ f (500)) and the effective radius (r eff [μm])).
Model Results obtained from the simulation of errors in τ a on GRASP-AOD code. The parameters representing the bimodal log-normal size distribution (volume median radius (r  Table 1). The following columns contain two secondary products derived from GRASP-AOD inversion: fine-mode aerosol optical depth at 500 nm (τ f (500)) and the effective radius (r eff [μm]) (reference values in Table 3). The last column R indicates the fraction of the inversions that pass the filters out of the original 1000 simulations.
Model Results obtained from the simulation of errors in τa on GRASP-AOD code for the cases with the lowest and the highest aerosol loads. As in Table 7, the first six columns contain the parameters representing the bimodal log-normal size distribution (volume median radius (r  )). The following columns contain two secondary products derived from the GRASP-AOD inversion: fine-mode aerosol optical depth at 500 nm (τ f (500)) and the effective radius (r eff [μm]) (reference values in Table 3). The last column R indicates the fraction of the inversions that pass the filters out of the original 1000 simulations.  Differences in the mean values and the standard deviations between the studies of the error simulation for the spectral range 440-1020 and 340-1640 nm (these last values represented in Table 7). As in previous tables, the first six columns contain the parameters representing the bimodal log-normal size distribution (volume median radius (r  Change in the aerosol optical depth (references in Table 2) caused by changes of ±2σ in the complex refractive index (n, k) for GSFC2, SOLV2 and ZAMB2 (reference values in Table 1). The value of 2σ n is equal to 0.02 for GSFC and ZAMB, while it is equal to 0.06 for SOLV. For the imaginary part, 2σ k has values of 0.008, 0.006 and 0.002 for ZAMB, GSFC and SOLV respectively. Note that k − 2σ k can be negative at some wavelengths; for those cases, k is fixed as zero and is indicated by an asterisk in the table.  Fine-mode aerosol optical depth at 500 nm (τ f (500)) obtained varying the refractive index by ±2σ during the retrieval process. The value of 2σ n is spectrally independent and equal to 0.02 for GSFC and ZAMB, but is equal to 0.06 for SOLV. For the imaginary part, 2σ k is also spectrally independent and has values of 0.008, 0.006 and 0.002 for ZAMB, GSFC and SOLV respectively. Note that k − 2σ k produces negative values of k for all wavelengths in the aerosol example GSFC and for some in SOLV. In those cases, k is fixed as zero in the retrieval.  Sphericity tests results. Differences between GRASP-AOD code with modified sphericity parameters and GRASP-AOD code with reference sphericity value of 0.6. The first six columns depict the bimodal log-normal size distribution parameters: volume median radius (r ). The last two columns contain the fine-mode aerosol optical depth at 500 nm (τ f (500)) and the effective radius (r eff [μm]).  Summary of the 29 days chosen for comparison between the new GRASP-AOD inversion and the AERONET products. The first columns depict information about the day selected: site, photometer, date, τ a at 440 nm and Ångström exponent. The last two columns contain the number data used for each inversion type: ALM or almucantar used for the AERONET inversion, and τ a used for both the GRASP-AOD inversion and SDA algorithm. STD and EXT refer to the standard and extended version of the Cimel-310 sun photometer respectively (see Table A1). The 340 nm channel of photometer #172, at Lampedusa site, was excluded for quality assurance reasons.