**Aerosol Optical Radiation Properties in Kunming (the Low–Latitude Plateau of China) and Their Relationship to the Monsoon Circulation Index**

**Haoyue Wang 1,2, Chunyang Zhang 1,3, Ke Yu 4, Xiao Tang 5, Huizheng Che 6,\*, Jianchun Bian 7, Shanshan Wang 2, Bin Zhou 2, Rui Liu 8, Xiaoguang Deng 1, Xunhao Ma 1, Zhe Yang 1, Xiaohang Cao 1, Yuehua Lu 1, Yuzhu Wang <sup>1</sup> and Weiguo Wang <sup>1</sup>**


Received: 10 October 2019; Accepted: 2 December 2019; Published: 5 December 2019

**Abstract:** Based on the Langley method and the EuroSkyRad (ESR) pack retrieval scheme, we carried out the retrieval of the aerosol properties for the CE–318 sunphotometer observation data from March 2012 to February 2014 in Kunming, China, and we explored the possible mechanisms of the seasonal variations. The seasonal variation of the aerosol optical depth (AOD) was unimodal and reached a maximum in summer. The retrieval analysis of the Angstrom exponent (α) showed the aerosol types were continental, biomass burning (BB), and urban/industrial (UI); the content of the desert dust (DD) was low, and it may have contained a sea–salt (SS) aerosol due to the influence of the summer monsoon. All the aerosol particle spectra in different seasons showed a bimodal structure. The maximum and submaximal values were located near 0.2 μm and 4 μm, respectively, and the concentration of the aerosol volume was the highest in summer. In summer, aerosol particles have a strong scattering power but a weak absorption power; this pattern is the opposite in winter. The synergistic effect of the East Asian monsoon and the South Asian monsoon seasonal oscillations can have an important impact on the variation of the aerosol properties. The oscillation variation characteristic of the total vertical columnar water vapor (CWV) and the monsoon index was completely consistent. The aerosol types and sources in the Yunnan–Kweichow Plateau and the optical radiation properties were closely related to the monsoon circulation activities during different seasons and were different from other regions in China.

**Keywords:** Yunnan–Kweichow Plateau; low–latitude plateau monsoon climate; aerosol type and source; aerosol properties; monsoon index; seasonal variation

#### **1. Introduction**

The complexity of aerosol sources (primary and aged aerosol) determines the diversity of the composition types and particle sizes, and it has important effects on the weather, climate, air quality, and human health. The direct and indirect effects of aerosols are some of the most poorly–understood variables affecting atmospheric and water circulation. It is necessary to combine aerosol measurements and radiation techniques with model simulations to accurately determine the aerosol effect [1–4]. Model simulation uses atmospheric radiation transfer models, such as the Santa Barbara DISORT(Discrete Ordinates Radiative Transfer) Atmospheric Radiative Transfer (SBDART), Global Atmospheric Model (GAME), Moderate Resolution Atmospheric Transmission (MODTRAN), and RSTAR. To estimate the accuracy of the models, the modeled global, diffuse, and Sun's direct radiation at ground level have been compared with experimental measurements. The Sun–sky radiometric measurements at ground level were also needed to characterize the aerosol properties. The CE–318 sunphotometer is used internationally. Systems such as Prede POM can be used to retrieve the optical properties of aerosols from solar and sky radiance measurements [5]. Aerosols with different particle sizes have different effects on the environment, transport, and human health [6,7]. However, due to the complex composition of aerosols and the uncertainty of their spatiotemporal distribution, the physicochemical and optical properties of aerosols vary [8], which makes environmental research and assessing the radiation effects of climate change challenging [9–13]. Therefore, quantitative remote sensing analysis of aerosol properties in different geographical locations has important theoretical and practical implications.

At present, a number of international and regional aerosol ground–based observation networks have been established. For example, the Aerosol Robotic Network (AERONET), established by the National Aeronautics and Space Administration (NASA), provides the most extensive aerosol database in the world [14]. The European Skynet Radiometer Network (ESR) is a new type of network in partnership with SKY–Radiometer NETwork–SKYNET (Skynet–Asia) in Japan [15,16]. There is also the Global Atmosphere Watch Precision Filter Radiometer Network (GAWPFR NET) [17]. The Chinese Sun Hazemeter Network (CSHNET), the China Aerosol Remote Sensing Network (CARSNET) [18], and the Campaign for the Atmospheric Aerosol Research Network of China (CARE–China) were established in China [19,20]. ground–based network observation can directly provide basic data for research; it can also provide a reference for the satellite detection data and numerical simulation results and provide an observational basis for the impact of environmental, weather and climate change [14,21–25].

The sunphotometer is one of the most widely used instruments for ground–based passive telemetry to accurately characterize aerosol properties, and its retrieval algorithm is an important component. Dubovik et al. evaluated the physical quantities, errors, and information from different sources. They gave different weights to the data in the retrieval process and applied advanced numerical optimization techniques to obtain the final statistical optimal solution. Through this algorithm, they obtained parameters such as the aerosol particle volume size distribution, complex refractive index, and single scattering albedo by using ground–based observation data [26]. Olmo et al. applied the Nakajima algorithm to incorporate the randomly distributed ellipsoid approximation to retrieve the aerosol parameters [27]. The ESR.pack retrieval scheme proposed by Estellés et al. [28] was based on the Nakajima algorithm and Skyrad.pack to improve and compile the applications of the CE–318, POM, and various sunphotometers. He et al. [29] compared the retrieval results of the aerosol properties generated by Nakajima et al. [30] and Dubovik et al. [26]. Huang et al. used the sun direct radiation data to retrieve the aerosol optical depth (AOD) and Angstrom exponent as well as the single scattering albedo (SSA) and scattering phase functions from the sky scattering data [31]. At present,

the internationally–used, comprehensive aerosol retrieval results come from AERONET (for business retrieval) [26] and the ESR.pack retrieval scheme that covers ESR and Skynet–Asia [28]. Estellés et al. performed a comparison of the aerosol properties derived by the ESR.pack with AERONET and obtained good results [32].

The Yunnan–Kweichow Plateau (100–111◦E, 22–30◦N) is located in the southeastern part of the Qinghai–Tibet Plateau in southwestern China and is a typical low–latitude plateau monsoon climate zone. Its geographical location (it is adjacent to Southeast Asia and South Asia, and it is an important source of aerosols in southwest China due to the strong influence of the monsoon circulation) [33,34] and climate (the annual temperature difference is small and the daily range is large, the dry season and wet season are distinct, and the sunshine and ultraviolet radiation are strong) are unique. It is also an important path for water vapor transport of convection and advection in China. Research shows that aerosols have important effects on global and regional water vapor variation [35,36] and influence the regional climate and environment. However, there are few ground–based stations in the Yunnan–Kweichow Plateau, and there is a lack of systematic research on aerosols. Therefore, this paper explores the relationship between the seasonal and interannual variation characteristics of the East Asian monsoon and South Asian monsoon and the variation of the aerosol optical radiation properties, which were based on the combination of the aerosol properties retrieved from the CE–318 observation data from the Kunming Atmospheric Ozone Monitoring Station, No. 209 of the Global Ozone Observing System (GO3OS), and the monsoon circulation index. It is important to know and understand the aerosol variation and monsoon activities in the Yunnan–Kweichow Plateau and their impact on the environment and climate in specific areas.

#### **2. Instruments and Data**

#### *2.1. Instruments*

The Kunming atmospheric ozone monitoring station (25.03◦N, 102.68◦E; 1917 m above sea level (a.s.l.)) is located in the center of the city, and it is equipped with a CE–318 sunphotometer (CE318NTS8, France). It has eight channels with central wavelengths λ at 340, 380, 440, 500, 670, 870, 1020 and 1640 nm. The instrument can track the Sun for direct radiation observation, and it scans the sky for the scattered radiation of the Almucantar–azimuth angle, principal plane standard–scattering angle, and principal plane polarization–zenith angle.

The CE–318 sunphotometer can be used for atmospheric environmental monitoring and the radiometric calibration of remote sensing satellite sensors. Calibration is carried out every 6 months using the calibration facilities at the Chinese Academy of Meteorological Sciences [18,37]. We conducted sun direct radiation calibration by comparing the instrument with the master sunphotometers in Beijing. The master sunphotometers were calibrated using the Langley method at either Izaña (Spain, 28.31◦N, 16.50◦W; 2391.0 m a.s.l.) or Mauna Loa (HI, USA, 19.54◦N, 55.58◦W; 3397.0 m a.s.l.) [10]. The sky scattered radiation channel is calibrated by integrating the sphere radiation source. Tao et al. described the sphere calibration methods and protocols for CARSNET. The CARSNET sphere calibration results were compared with the original values provided by Cimel, the manufacturer; the linear interpolation method can be used to obtain the calibration coefficient for each period [37,38].

#### *2.2. Data*

We used the CE–318 observation data of the Kunming atmospheric ozone monitoring station from March 2012 to February 2014 (including eight bands of Sun direct radiation data and four bands of equal zenith angle scanning scattered radiation data; the minimum observation time point is in minutes, and there are missing measurements when it is cloudy) and the total ozone column data observed in real time by the Dobson ozone spectrophotometer (No. D003). The absorption coefficients of ozone and the water vapor at different wavelengths and temperatures were derived from the ESR.

We used ERA–Interim daily four times (00,06,12,18 UTC) over the same period of the wind field and relative humidity reanalysis data provided by the European Center for Medium–range Weather Forecasts (ECMWF). The data have a total vertical resolution of 37 layers, with a horizontal resolution of 0.25◦ × 0.25◦. The daily surface meteorological data provided by the National Meteorological Base Station in Kunming included the hourly ground pressure, temperature, and relative humidity. The East Asian monsoon index (EAMI) and South Asian monsoon index (SAMI) were provided by the China National Climate Center.

NASA provided the water vapor data observed by the Atmospheric Infrared Detector (AIRS) carried by the Aqua satellite, which was compared with the vertical columnar water vapor (CWV) retrieved by the CE–318 observation data. The MODIS Fire Points Product Data (MCD14ML) were provided by the University of Maryland website (ftp://fuoco.geog.umd.edu), and the description and validation of the data are shown in the literature [39].

#### **3. Retrieval Algorithms**

#### *3.1. Sun Direct Radiation*

The direct solar irradiance of the ground–based measurement of a specific wavelength (λ) is based on the Beer–Bouguer–Lambert rule, which is expressed as the output voltage *V* on the CE–318 [40]:

$$V = V\_0 R^{-2} T\_\mathcal{S} \exp(-m\tau) \tag{1}$$

where *V*<sup>0</sup> is the calibration voltage constant and can be obtained by extrapolating from a series of observations to *m* = 0. *R* = *rt*·*rm*−<sup>1</sup> is the Sun–Earth distance factor at the time of measurement, *rt* is the distance between the Sun and the Earth at the time of measurement, and *rm* is the mean distance between the Sun and the Earth. *Tg* is the gas absorption transmission rate (mainly considered ozone and water vapor). In the CE–318 channel, only the water vapor absorption band at 936 nm cannot be ignored, so *Tg* of the other channels is 1. *m* = (cosθ) <sup>−</sup><sup>1</sup> is the air mass factor, and θ is the solar zenith angle. τ is the total optical depth of the vertical atmosphere. Taking the logarithm of Equation (1) on both sides, the reciprocal of the slope of the line drawn by ln*V* + 2ln*R* and *m* is τ, which is the Langley method [41,42]:

$$\tau = -\frac{1}{m} [\ln(V/V\_0) + 2\ln R] \tag{2}$$

In the formula, τ = τ*<sup>a</sup>* + τ*<sup>g</sup>* + τ*<sup>r</sup>* is composed of aerosol scattering τ*a*, gas absorption τ*<sup>g</sup>* (such as ozone and water vapor), and Rayleigh scattering τ*r*. Except for the 936 nm channel, τ*<sup>g</sup>* of the other channels can be ignored, and τ*<sup>r</sup>* can be calculated from the actual measured value of the ground pressure. Then the AOD is τ*<sup>a</sup>* = τ − τ*r*.

Assuming that the aerosol particles follow the Junge volume size distribution [43], the calculation is:

$$m(r) = \frac{dN(r)}{dr} = c(z)r^{-(v+2)}\tag{3}$$

where *r* is the radius of the spherical particle, *N*(*r*) is the total number of aerosol particles per unit area, *v* is the Junge parameter, and *c*(*z*) is a function of height *z*. If λ is independent of the Junge spectrum type and the complex refractive index, Ångström generalizes the relationship between τ*<sup>a</sup>* and λ [44]:

$$
\pi\_{\mathfrak{a}}(\lambda) = \beta \lambda^{-\alpha} \tag{4}
$$

Here, α is the Angstrom exponent. Usually, 0 < α < 4 reflects the scale characteristic of the aerosol particle size and is inversely proportional to the particle size. When the coarse particles dominate, α tends to 0, and when the particle size is on the molecular scale, α is close to 4 [45]. β is the Angstrom

turbidity coefficient, and the AOD at λ = 1 μm. When β ≥ 0.20, the atmosphere is turbid; if β ≤ 0.10, the atmosphere is relatively clean [46]. From Equation (4), the following can be obtained:

$$\alpha = -\frac{\ln\left[\pi\_a(\lambda\_1) / \pi\_a(\lambda\_2)\right]}{\ln(\lambda\_1/\lambda\_2)}\tag{5}$$

$$\beta = \exp[\ln \tau\_{\mathfrak{a}}(\lambda) + a \ln \lambda] = \tau\_{\mathfrak{a}}(\lambda)\lambda^{\mathfrak{a}} \tag{6}$$

In the formula, if τ*a*(λ1) and τ*a*(λ2) of two wavelengths λ<sup>1</sup> and λ<sup>2</sup> are known, α and β are obtained. Thus, τ*a*(λ) of arbitrary λ under the same conditions is calculated.

The formula of the water vapor transmission rate *T*w on the channel is [47]

$$T\_w = \exp(-aw^b) \tag{7}$$

Here, *w* is the total amount of water vapor in the slant path, and *a* and *b* are constants. The Sun direct radiation response of CE–318 at the 936 nm water vapor absorption band is:

$$V = V\_0 R^{-2} T\_w \exp(-m\tau\_{ar}) \tag{8}$$

where τ*ar* = τ*<sup>a</sup>* + τ*r*, and τ*<sup>a</sup>* is obtained by interpolating two channels at 870 and 1020 nm. At the same time, *w* = *m*·*W*C, and *W*<sup>C</sup> is the total amount of the vertical water vapor column (CWV). By combining this with Equations (7) and (8), we obtain:

$$
\ln V + m\tau\_{ar} = \ln(V\_0 R^{-2}) - am^b \mathcal{W}\_c^b \tag{9}
$$

*W*<sup>C</sup> (CWV) can be obtained, which is an improvement of the Langley method because it includes the influence of water vapor, making it more accurate and reliable [42].

#### *3.2. Equal Zenith Angle Scattered Radiation*

In the single channel of the radiation transmission model, the surface sun direct radiation is defined as *E*, and the surface scattered radiation is *F* [30]:

$$\begin{aligned} E &= E\_0 \exp(-m\_0 \tau) \\ F(\theta, \phi) &\equiv F(\Theta) = E m\_0 \Delta \Omega [\omega \tau P(\Theta) + q(\Theta)] \end{aligned} \tag{10}$$

where *E*<sup>0</sup> is the sun direct radiation of the channel at the upper boundary of the atmosphere (unit: <sup>W</sup>·m−2·μm<sup>−</sup>1), and *<sup>m</sup>*<sup>0</sup> is the atmospheric optical mass. For <sup>θ</sup> <sup>≤</sup> <sup>75</sup>◦, *<sup>m</sup>*<sup>0</sup> <sup>=</sup> (cosθ) <sup>−</sup>1. ΔΩ is the stereoscopic observation angle of the sunphotometer, calculated at 1.2◦, ω is the SSA of the entire atmosphere [48], Θ is the scattering angle, *P*(Θ) is the total phase function when the scatter angle is Θ [49,50], and *q*(Θ) represents the contribution of multiple scattering (MS). The relationship between Θ and θ is observed in the actual equal zenith angle scan observation:

$$\cos(\Theta) = \cos^2 \theta + \sin^2 \theta \cos(\phi - \phi\_0) \tag{11}$$

At this time, ϕ<sup>0</sup> is the solar azimuth angle, and 0 ≤ Θ ≤ 2θ0. We introduce *G*(Θ) and define it as:

$$G(\Theta) \equiv \frac{F(\Theta)}{Em\_0 \Delta \Omega} = \omega \tau P(\Theta) + q(\Theta) \equiv \beta(\Theta) + q(\Theta) \tag{12}$$

where β(Θ) = ωτ*P*(Θ) is a single scattering equal to the total differential scattering coefficient, including Rayleigh scattering and Mie scattering. The AOD is defined as:

$$\pi\_d(\lambda) = \int\_{r\_m}^{r\_M} \pi r^2 Q\_{ext}(\mathbf{x}, \widetilde{m}) n(r) dr \tag{13}$$

Here, *Qext* is the extinction effective factor of spherical particles, *x* = (2π/λ)·*r* is the Mie size parameter (*x* 1 is Rayleigh scattering, *x* ≥ 1 is Mie scattering), *rm* and *rM* are the minimum and maximum values of the aerosol particle radius, respectively, and *<sup>m</sup>* <sup>=</sup> *mR* <sup>+</sup> *imI* represents the complex refractive index of the aerosol. If *Qext* is only a scattering effective factor, the corresponding calculation result is the scattering optical depth τ*as*, from which the formula of the SSA can be obtained:

$$
\omega\_a = \frac{\tau\_{as}}{\tau\_a}.\tag{14}
$$

The differential scattering coefficient β*<sup>a</sup>* of the entire layer of the atmospheric aerosol is defined as:

$$\beta\_a(\Theta) = \frac{\lambda^2}{2\pi} \int\_{r\_m}^{r\_M} [i\_1(\Theta, \mathbf{x}, \widetilde{m}) + i\_2(\Theta, \mathbf{x}, \widetilde{m})] n(r) dr \tag{15}$$

where *i*<sup>1</sup> and *i*<sup>2</sup> are the Mie scattering intensity functions, thereby defining the aerosol phase function as:

$$P\_{\mathfrak{a}}(\Theta) = \frac{\beta\_{\mathfrak{a}}(\Theta)}{\omega\_{\mathfrak{a}}\tau\_{\mathfrak{a}}} \tag{16}$$

The asymmetry factor *g* is the first moment of the phase function and is used to describe the relative intensity of the forward scattering:

$$\mathcal{g} = \frac{1}{2} \int\_{-1}^{+1} P(\Theta) \cos \Theta d \cos \Theta \tag{17}$$

Here, *g* varies from −1 to 1. Usually, the Mie scattering has a peak *g* > 0 in the forward direction, and Rayleigh scattering has the same property in all directions *g* = 0. g is also an important parameter for discussing radiation transmission problems and aerosol properties.

The aerosol particle size spectrum *<sup>n</sup>*(*r*) is the number of particles (cm−2·cm<sup>−</sup>1) in a unit cross-section of the gas column and a unit radius interval [51]. The volume spectrum distribution *V*(*r*) of the particles is defined as the volume of the aerosol (cm3·cm<sup>−</sup>2) in a unit cross–section gas column and a unit logarithmic radius interval:

$$V(r) = \frac{4\pi}{3}r^4 n(r) \tag{18}$$

In the ESR retrieval scheme, τ*<sup>a</sup>* and β*<sup>a</sup>* can be summarized as:

$$\begin{array}{l} \pi\_{\mathfrak{a}}(\lambda) = \frac{2\pi}{\lambda} \int\_{r}^{r\_{M}} K\_{\text{ext}}(\mathbf{x}, \widetilde{m}) v(r) d(ln r) \\ \beta\_{\mathfrak{a}}(\Theta) = \frac{2\pi}{\lambda} \int\_{r\_{\mathfrak{m}}}^{r\_{M}} K(\Theta, \mathbf{x}, \widetilde{m}) v(r) d(ln r) \end{array} \tag{19}$$

where *Kext* and *K* are core functions:

$$K\_{\rm ext}(\mathbf{x}) = \frac{3}{4} \frac{Q\_{\rm ext}(\mathbf{x})}{\mathbf{x}} K(\Theta, \mathbf{x}, \widetilde{m}) = \frac{3}{2} \frac{i\_1 + i\_2}{\mathbf{x}^3} \tag{20}$$

The values of the above two equations determine the radius of the aerosol particles, which in turn have a greater impact on the physical and optical properties of the aerosol [52,53].

#### *3.3. Retrieval Scheme*

The Sun direct radiation data are filtered by clouds, and the instrument retrieves the AOD of 8 channels (340, 380, 440, 500, 670, 87, 1020 and 1640 nm), α, β, CWV, and other properties. Smironv et al.'s [54] cloud filtering algorithm was applied to remove artifacts introduced by clouds. When retrieving *<sup>V</sup>*(*r*), *<sup>P</sup>*(Θ), *<sup>g</sup>*, SSA, *<sup>m</sup>* and the other parameters using the sky scattering data observed by equal zenith angle scanning, it was necessary to eliminate cloud interference by examining the data symmetry. Holben et al.'s [14] quality control and cloud–removal scheme was applied to obtain high

precision retrieval information: (1) The scattered radiation data used for retrieval had to satisfy the symmetry angle of more than 21 angles (the maximum quantity was 28); (2) the fitting error of the retrieval result was less than 5%; and (3) it needed to meet the condition of θ > 50◦ to perform the retrieval calculations.

#### **4. The Characteristic Analysis of the Seasonal Variations**

#### *4.1. Optical Depth*

Figure 1 shows the seasonal mean variation of the eight λ channel AOD year by year and the variation characteristics of AOD440 from March 2012 to February 2014. In Figure 1a, the AOD decreases with the increase of λ, and both have the same characteristics. The value for the AOD from June to August (JJA) in summer is greater than that of March to May (MAM) in spring, which is greater than that of autumn or September to November (SON), which is greater than that for winter or December to February (DJF). The seasonal mean is less than 0.80, but there are some differences in the seasonal mean variation year by year. This seasonal variation is related to the monsoon circulation activities (as shown in Figure A1) and environmental pollution emissions.

**Figure 1.** (**a**) Seasonal mean variation of the AOD in the 8 bands year by year. (**b**) Seasonal mean columnar distribution of AOD440 and comparison of the seasonal variations year by year. (**c**) Distribution of scatter data between AOD440 and the Angstrom Exponent (α).

The seasonal mean value of AOD440 from Figure 1b is between 0.2 and 0.6 and is less than 0.50, except for the summer of 2012. After analyzing the AOD of the eight channels, we found that they reach their maximum value in the summer of 2012; there are two reasons for this. First, combined with the seasonal variation of the 700 hPa level atmospheric circulation and relative humidity (RH) field in Figure A1, the summer monsoon circulation in 2012 brought more water vapor than that in 2013, and the summer RH around Kunming was 90% in 2012 and 80% in 2013. Second, it may be related to the increase of aerosol particles in the atmosphere caused by the large–scale municipal construction in the urban area of Kunming in the summer of 2012. However, we cannot rule out that the increase of aerosol particles was caused by other factors.

In Figure 1c, the value of AOD440 is concentrated between 0.10 and 1.0, and the Angstrom Exponent (α) is uniformly distributed in the range of 0.2–1.6. α increases with the increase of AOD440, indicating that it is mainly affected by fine particle aerosols generated by human activities, and the particle types and sources are different in different seasons (Table A1).

#### *4.2. Angstrom Exponent and Turbidity Coe*ffi*cient*

Figure 2 shows a seasonal variation of the Angstrom Exponent (α) and frequency distributions at different intervals. The statistical results show the α value is distributed between 0.2 and 1.7. The highest frequency of 0.6–1.0 is 39.90%, and 1.4–1.8 only occurs 2.54% of the time. Figures 5 and 11 show radius *r* of the main control aerosol particles is more than 0.5 μm.

*Remote Sens.* **2019**, *11*, 2911

**Figure 2.** Retrieval results of the Angstrom Exponent (α). (**a**) Seasonal variation of α values. (**b**) The frequency distribution of α values at different intervals. (**c**) Frequency distribution of α in different seasons.

In the summer of 2012, the α value was slightly higher than that in the spring, while the autumn and winter values decreased in turn, reflecting that the mean aerosol particle size *r* in the spring and summer (autumn and winter) was relatively small (large). In 2013, the α value gradually decreased from spring to winter, indicating that there is a certain difference in the mean aerosol particle size *r* in different seasons. The grain size *r* is slightly larger (small) in autumn and winter (spring and summer). On the whole, the seasonal mean of the Angstrom Exponent (α) is between 0.6 and 0.9, which is slightly greater in spring and summer than in autumn and winter. It is noted that the mean value of α in the autumn and winter of 2012 decreased significantly, and the annual mean α value was slightly lower than that in 2013, mainly because the water vapor transport in 2012 was greater than that in 2013, so the hygroscopic growth effect was more significant.

For frequency statistics in different seasons, the distribution frequency of 0.6 ≤ α < 1.0 varies with the season and is the most stable and highest in autumn. In spring and summer, 1.0 ≤ α < 1.4 (0.2 ≤ α < 0.6) increases (decreases), while in autumn and winter, it decreases (increases), which shows the opposite frequency seasonal variation. The distribution frequency of 1.4 ≤ α < 1.8 is the lowest and gradually decreases from spring to winter.

A further analysis of the data in Figure 2, Figure A2, and Figure A3 showed that the mean value of α is 0.88 in spring, and the frequency in the range of 0.6–1.4 is more than 75%. This result is mainly due to the frequent biomass burning around the Yunnan–Kweichow Plateau (Southeast Asia–South Asia) (Table A2), which produces a large amount of fine particle aerosols. It is also related to local anthropogenic emissions (industrial pollution, coal burning, motor vehicles, and other human activities). The mean value of α in summer is 0.87, and the frequency in the range of 0.6–1.4 is more than 80%; compared with the data from spring, the frequency of the small value range 0.6–1.0 increased, while that of the large value range 1.0–1.4 decreased slightly. These results are mainly because the summer monsoon circulation in East Asia and South Asia brings adequate water vapor to make fine particle aerosols grow hygroscopically, which leads to the large AOD in summer (Figure A1). The mean value of α in autumn is 0.76, and the frequency in the range of 0.6–1.4 is about 65%; the frequency in the small value range increased significantly. The mean value of α in winter is 0.70. The frequency between 0.6 and 1.4 is less than 60%, while the frequency between 0.2 and 0.6 is more than 40%, and the frequency between 1.4 and 1.8 is almost 0. In winter, Southeast Asia–South Asia also experienced more frequent biomass burning, resulting in the transport of aerosol particles to the Yunnan–Kweichow Plateau.

By comparing the values of the different intervals' α in Table A1, we found that there are main aerosols, such as UI (urban/industrial), BB (biomass burning), continental, and DD (desert dust), with some SS (sea–salt) aerosols being imported in the wet season. However, the most dominant content is that of the continental aerosol, followed by the BB and UI aerosols; the DD aerosol content is relatively low. From spring to winter, the dominant particles are coarse particles with a high frequency and fine particles with a low frequency.

Figure 3 shows the seasonal mean variation of the β coefficients for the 440 nm and 870 nm bands. The β values at the 2 channels differ little and vary substantially. The comparison shows that the β value is consistent with the seasonal variation of the AOD in Figure 1. The seasonal mean is between 0.10 and 0.30, which is the largest in summer, followed by the values for spring, autumn, and winter. The mean value of β in the summer of 2012 is the largest, indicating that the degree of atmospheric turbidity was the highest.

**Figure 3.** Seasonal mean variations in the Angstrom turbidity coefficient β at 440 nm and 870 nm and seasonal mean variations year by year.

Table A3 shows the division of different β values and atmospheric turbidity. Figure 3 and Table A3 illustrate that there is little atmospheric turbidity. The degree of atmospheric turbidity in spring and summer is significantly higher than that in autumn and winter; the β value in winter is slightly higher than 0.1, and the atmosphere is the cleanest.

#### *4.3. Total Column Water Vapor*

Figure 4 shows a comparison of the retrieval results of the CWV and observations with AIRS. The seasonal variation of the CWV is obvious, with an annual mean of about 1.0 g·cm−<sup>2</sup> that reaches the maximum of more than 3.0 g·cm−<sup>2</sup> in summer. Spring and autumn are similar, and the minimum is in winter. This is consistent with the climatic characteristics of the outbreak and end of summer monsoon over the Yunnan–Kweichow Plateau, and the retrieval results are consistent with the CWV from the AIRS detection.

**Figure 4.** Retrieval results of column water vapor (CWV). (**a**) The comparison of the CWV inter-monthly variation in CE–318 and AIRS retrieval. (**b**) The comparison of CWV seasonal variation in CE–318 and AIRS retrieval.

It is not difficult to see that the seasonal variation of the AOD and CWV is similar, but slightly different. Overall, it reaches the maximum in summer and the minimum in winter, but the difference is

that the AOD is significantly higher in spring than that in autumn, while the CWV has no significant size difference in spring and autumn.

#### *4.4. Particle Spectrum Distribution and Complex Refractive Index*

Figure 5 shows the mean annual aerosol particle spectrum distribution in different seasons. The distribution of the particle spectra varies in different seasons. However, the variation situation is basically the same, and it shows a bimodal shape. The maximum values of the fine mode and the coarse mode are located near 0.2 μm and 4 μm, respectively. The particle spectral structure is similar to the ratio of the continental and UI aerosol models in the Standard Chinese Radiation Atmosphere [55]. The seasonal variation analysis shows the volume concentration in summer, and the distribution of fine particles and coarse particles reaches the maximum. The maximum value of the fine particles is 0.07, and the maximum value of the coarse particles is 0.06. The main reason is the same as that for the seasonal variation of the AOD. In winter it is the smallest, with a maximum particle value of 0.05 and a maximum coarse particle value of 0.02. In winter, due to the influence of the dry inland winter monsoon (Figure A1), the aerosol content is the lowest, and the particle volume concentration is also the smallest.

**Figure 5.** The mean annual aerosol particle spectrum distribution in different seasons.

In contrast, the retrieval results of the aerosol volume spectrum (the fine mode of 0.1–0.2 μm has a maximum value) in the Beijing area [56] are slightly smaller, reflecting the relatively light pollution in Kunming. Fine particles with a radius *r* below 1 μm occupy the main body and are mainly particles in the Aitken nuclei mold and the accumulation mode; the components of these particles are UI and continental aerosols.

Figure 6 shows the relationship between the annual mean of the FR and FI in the complex refractive index and λ, and only the seasonal mean variation of the 440 nm complex refractive index is extracted. The variation of the FR and FI is the opposite when they vary with λ, and its absolute value increases with the increase of λ. The variation of the FR is between 1.40 and 1.50, which is more obvious than the variation of FI, which is between 0.005 and 0.015. The variations in the FR and FI are not significant in the relatively short 440–670 nm band, while the FR and FI are significantly different in the relatively longer bands. The difference in the FR for different aerosol particles is not large, and the FI values can differ by several orders of magnitude.

**Figure 6.** (**a**) Relationship between *<sup>m</sup>***<sup>R</sup>** and *<sup>m</sup>***<sup>I</sup>** of the annual average aerosol complex refractive index and wavelength (λ). (**b**) Seasonal variation of *<sup>m</sup>***<sup>R</sup>** and *<sup>m</sup>***<sup>I</sup>** of the aerosol complex refractive index at 440 nm.

The FR and FI characterize the scattering and absorption properties, respectively, of an aerosol; they offer a comprehensive reflection of the aerosol absorption properties of different components. Black carbon is the most absorbent component in aerosols, and minerals and dust are also important absorbent components [57]. Table A4 shows the statistical mean results of the values of the aerosol FR and FI for different components. Figure 6 illustrates that the content of water vapor and ammonium sulfate in the aerosol is relatively high; it also contains a very low amount of dust and the black carbon aerosol.

Studies [58] have shown that the value of the FR in urban/industrial areas is between 1.4 and 1.47, and if the area affected by the ocean is large, the value of the FR will be low. From the seasonal variation of the FR in Figure 6, it can be determined that Kunming is an urban/industrial area. The seasonal mean variation of the FI is below 0.01; it is the closest to 0.01 in winter, which is relatively large, while it is relatively small in summer, and thus the absorption aerosol content is low. However, there is also some uncertainty. Because the aerosol content in winter is low, the FI value is relatively large, so a more in–depth discussion of the retrieval results is needed to explain the variation of the FI. It may be worthwhile to use the refractive index to study the seasonal cycles of the aerosol types.

#### *4.5. Single Scattering Albedo and Asymmetry Factor*

Figure 7 shows the seasonal mean variation of the SSA in the 440 nm band and the *g* in the 4 bands. In Figure 7a, the seasonal variation of the SSA is a unimodal type; in summer (winter), it reaches a maximum (minimum) of about 0.96 (0.90). The water vapor content is high in summer, and the particle size and volume increase after the hygroscopic aerosol absorbs water, which enhances the scattering. The summer monsoon brings abundant water vapor and a little SS aerosol, and the non–absorbing sulfate aerosol enhances the scattering ability to some extent. The amount of wind and sand near the ground in spring is relatively large, which increases the amount of the flying dust and floating dust particles; the single scattering albedo (SSA) is slightly higher than that in autumn. The SSA is smaller in winter; the main reason is that the strong absorption of the aerosol scattering ability is relatively weak, which is consistent with the analysis of the complex refractive index in Figure 6.

**Figure 7.** (**a**) Seasonal variation of the aerosol single scattering albedo (SSA) at 440 nm (the deviation line is the standard deviation). (**b**) Seasonal variation of the asymmetry factor (*g*) in 4 bands.

Figure 7b reflects the scattering rule of the aerosol particle *g* > 0, and it exhibits the variation that *g* decreases with the increase of λ. The shorter λ the stronger is aerosol forward scattering. The *g* value of the same λ in different seasons is not much different and is not affected by the seasonal variation. The result is reflected in that *g* has no obvious seasonal variation rule, but the inverse relationship between *g* and λ is particularly obvious.

#### **5. Case Analysis**

#### *5.1. Retrieval of Direct Radiation*

Figure 8 shows the daily variation of the AOD, Angstrom Exponent (α), and CWV from 9:00 to 18:00 (Beijing Time) on January 9, 2014. The AOD varies with λ and shows the typical Mie scattering characteristics. The Mie scattering process exists when the aerosol scale parameter *RmM* = 2π*r*·λ−<sup>1</sup> is 0.1–50 [59]. Particles with *r* in the range of 0.3–0.7 μm have the greatest influence on the extinction of visible light. This result shows that the shorter λ the larger AOD, and the stronger is extinction effect of the particles.

**Figure 8.** Daily variation of the AOD, Angstrom Exponent (α), and CWV from 9:00 to 18:00 (Beijing Time) on January 9, 2014. (**a**) The 8 wavelengths (λ) correspond to the hourly mean AOD. (**b**) The AOD at 440, 670, 870, and 1020 nm at 11, 13, 15, and 18 o'clock, respectively. (**c**) The hourly mean Angstrom Exponent (α) (the deviation line is the standard deviation) and the CWV (unit: g·cm<sup>−</sup>2).

The daily variation of the hourly mean of the Angstrom Exponent (α) ranges from 0.6 to 1.2, which is significantly higher in the morning than in the afternoon; it decreases significantly in the afternoon and increases after 16:00. In the morning (afternoon), *r* of the aerosol's dominant particle is smaller (larger). By combining the AOD data and Table A5, we see that the AOD reaches the maximum at around 11:00, when the extinction effect of the aerosol is the strongest. At this time, the Angstrom Exponent (α) is 0.967. It is speculated that at the moment when the extinction effect is the strongest throughout the day, the dominant particle *r* should be around 1 μm. The CWV hourly mean is lower and varies from 1.0 g·cm−<sup>2</sup> to 1.3 g·cm−2, which is consistent with the climatic characteristics of the winter monsoon (dry and little rain).

#### *5.2. Retrieval of Scattered Radiation*

In order to comprehensively analyze the equal zenith angle observation data, the retrieval calculations were performed on the observation data at 10:00, 11:00, 12:00, and 13:00 (Beijing Time) on 9 January, 2014, and combined with the retrieval results of the direct radiation for analysis.

#### 5.2.1. The Detection of Sky Radiation

Figure 9 shows the clear sky data after the equal zenith angle scan and the cloud detection; the abscissa indicates the difference angle with the solar azimuth (the positive and negative signs indicate the right–handed and left–handed scans, respectively, when the CE–318 is observing), and the ordinate indicates the response value corresponding to the amount of radiation received by the instrument when scanning through the azimuth. The data points scanned in 2 directions at different time points have good symmetry, and the data points in the curve satisfy the condition of (*A*<sup>l</sup> − *A*r)·(*A*<sup>l</sup> + *A*r) <sup>−</sup>1·0.5 < 10%, where *A*l and *A*r represent the response values of the radiation received during the left–handed and right–handed scans of the instrument, respectively.

**Figure 9.** Data from the equal zenith angle scan observations at 4 times (Beijing Time) on 9 January, 2014.

#### 5.2.2. SSA and Phase Function

Table 1 shows the variation of the SSA (ω*a*) and *g* with λ at the 4 points. At 12:00, ω*<sup>a</sup>* reaches a maximum value, and AOD440 is greater than 0.40. Table 1 also shows that ω*<sup>a</sup>* is greater than 0.87 and shows a decrease with the increase of λ. Scattering plays a dominant role in the extinction effect of the aerosol on the radiation; λ is shorter, and the proportion of scattering is greater. By combining the CWV and Angstrom Exponent (*a*) data in Figure 1c, we find that the aerosol grows hygroscopically at 12:00 and 13:00, resulting in the enhanced scattering effect; therefore, the SSA reaches the maximum at 12:00. *g* decreases slightly with the increase of λ, but the difference at different times is small.

**Table 1.** Variation of the SSA (ω*a*) and *g* with λ at 10:00, 11:00, 12:00, and 13:00 (Beijing Time) on 9 January 2014.


Figure 10 shows the variation of *P*(Θ) with Θ for different λ. *P*(Θ) with different λ has a good consistency at each Θ angle and reaches the minimum near Θ = 120◦. The exponential curve of Θ and *P*(Θ) accords with the Junge model of the aerosol phase function.

**Figure 10.** The variation of the scattering phase function *P*(Θ) with Θ at 4 times (Beijing Time) on 9 January, 2014.

#### 5.2.3. Particle Spectrum and Turbidity

The range of particle *r* was determined to be 0.05–15 μm, and it contains 22 intervals in the ESR retrieval scheme with reference to the Skyrad scheme. The initial *<sup>m</sup>* of the particle is given during the retrieval process, and theoretical studies have shown that the effect on the retrieval results is small, so the empirical value of *<sup>m</sup>* is 1.500−0.005i. Figure 11 shows the relationship between the particle volume size distributions and the particle radius *r*. The measurement results at different times are from the bimodal spectrum and basically keep the variation in sync. The first peak mode *r* is located at 0.1–0.5 μm, and the second peak area *r* is located at 3–8 μm. For the 2 peak areas, the particle volume concentration is more concentrated on the fine particles with a smaller particle size *r*. The main reason is that anthropogenic aerosol particles are mostly fine particles such as those from the combustion processes, which produce a large number of submicron particles [60]. At 15 μm, the number of particles at 10:00 is greater than the number of particles at 11:00, 12:00, and 13:00; the radius of those particles can be more than 15 μm, which indicates that coarse aerosol particles may exist in the atmosphere. In addition, if the CWV content does not reach the saturation growth of the hygroscopic aerosol, it may also show an increase in the aerosol concentration of fine particles, resulting in a decrease in the effective radius of the particles.

**Figure 11.** Aerosol particle spectrum at 4 times (Beijing Time) on 9 January, 2014.

Figure 12 shows a daily variation of the turbidity coefficient. The β440 and β870 coefficients calculated using the AOD of 440 nm and 870 nm, respectively, are very close. Since the β coefficient represents the AOD at the wavelength of 1 μm (Equation (4)), it is proved that the calculated β coefficient is highly reliable. It was found that β and the AOD maintain almost simultaneous characteristics. It is also verified by the physical meaning of β that the AOD in Figure 1a decreases with the increase of λ. Obviously, the β value can reflect both the degree of turbidity in the atmosphere and the degree of extinction of the atmosphere. Table A3 shows the division of different β values and atmospheric turbidity. In Figure 12, the daily variation range of the β value is 0.10 ≤ β ≤ 0.20, the air quality is neither turbid nor clean, and the degree of turbidity at noon is relatively high.

**Figure 12.** The Angstrom turbidity coefficient β and the daily variation of the AOD corresponding to λ at 440 nm and 870 nm at 9:00 to 18:00 on 9 January, 2014.

#### *5.3. Atmospheric Circulation and Specific Humidity Field*

Through the analysis of the weather map of the 600–700 hPa levels on 9 January, 2014 (figure omitted), it can be seen that the dry west wind and southwest wind from the Indian Peninsula-Bangladesh–Myanmar prevailed throughout the day over the Yunnan–Kweichow Plateau. Not only is it less humid, but the air is very dry. The 600 hPa in Kunming is the westerly circulation, and the atmosphere is very dry and has a specific humidity of about 2, the speed of southwest wind at 700 hPa gradually increases, and there is a significant increase at 14:00 compared with the value at 8:00; the specific humidity is significantly higher than 600 hPa but only 4–6. This shows that the dry weather with little rain in winter is also conducive to the long distance transport of aerosols from Southeast Asia–South Asia (aerosol pollutants in the atmosphere) to the southwestern area of China. However, the contribution and influence of the aerosol load in the Kunming area need to be further observed and verified by numerical simulation.

#### **6. Analysis of the Causes of the Seasonal Variation**

Kunming is located in the central part of the Yunnan–Kweichow Plateau and is affected by the East Asian monsoon and the South Asian monsoon. The annual sunshine is about 2200 h, and the ultraviolet radiation is strong. The monthly mean temperature is between 9.1 ◦C and 20.7 ◦C, and the monthly precipitation is between 11.3 mm and 204.0 mm. The precipitation in the rainy season from May to October accounts for 85% of the whole year and for more than 60% in summer. The southeast's warm and humid airflow from the western Pacific and the southwest's warm and humid airflow from the Indian Ocean meet over the Yunnan–Kweichow Plateau (which contains Yunnan), and a variety of aerosol species combine with water vapor to generate aerosol particles, which are transported over long distances (Southeast Asia–South Asia and the Iranian plateau to North Africa) [33,61], or the aerosol particles are hygroscopically grown. Under the conditions of a high temperature, high humidity, and strong ultraviolet radiation, chemical and photochemical reactions are favored to generate new aged aerosol particles. Therefore, due to the increase in the type and quantity of aerosol particles in

summer (which makes the AOD reach its maximum), there is an increase in the volume concentration of the aerosol particles and the concentration of the coarse and fine particles. Diversification of the particle size (the Angstrom Exponent (α) is higher than the annual mean value) affects the degree of the atmospheric turbidity (the β coefficient reaches the largest value). Because the physical, chemical, and optical properties (absorption and scattering cross section) of the aerosol changed greatly, the extinction effect was enhanced. Moreover, the increase of the water vapor (the CWV reaches the maximum value) enhances the extinction of the solar radiation absorption and scattering, which leads to the enhancement of the secondary or multiple scattering of the aerosol particles and promotes the variation of the aerosol properties.

Precipitation in the dry season from November to April of the following year accounts for only 15% of the whole year. There is little rain and plenty of sunshine in winter. It is mainly influenced by the west wind circulation and cold air from the higher latitude continent. The source and the physicochemical properties of the aerosol are different from those in the summer, which leads to the characteristic parameters of the AOD, CWV, and β coefficient change with the seasons and the AOD, CWV, and β coefficient reach the minimum values in winter.

In the spring, there are few clouds and little rain, the air is dry, the evapotranspiration is strong, the solar ultraviolet radiation is enhanced, and the diurnal range of the temperature is large. The temperature drops rapidly in the autumn, which is about 2 ◦C lower than that in the spring, the precipitation is reduced and less than 30% of summer (but more than that in spring and winter), and the air is dry. The values of the aerosol properties in spring and autumn are between those of summer and winter, but they are slightly different in the two seasons. Comparing the spring and autumn, the aerosol properties are higher (larger) to different degrees. The main reason is that the wind speed in the near–surface layer in spring is much larger than that in the autumn, so that the flying dust and floating dust generated by the exposed surface can enter the atmosphere, which makes the amount of aerosol particles increase slightly. The increase of the desert dust aerosol in spring makes the average AOD higher than that in autumn, and the particle radius changes obviously make the turbidity stronger. Due to the influence of the summer monsoon circulation, the CWV values in spring and autumn are quite different (there is no obvious difference and regularity). The wet season starts at the end of spring and ends in autumn, but the variations of the CWV in the two seasons are slightly different. A more specific difference analysis is needed to combine the characteristics of the actual monsoon circulation evolution with different years, characteristics of water vapor transport, and variation in meteorological factors.

Figure 13 shows the characteristics of the seasonal variation and inter–annual activities of the EAMI and SAMI and aerosol properties normalized time series from March 2012 to February 2014. EAMI and SAMI are consistent with the variation in the aerosol properties, and the CWV and CWV (AIRS) are exactly the same as the seasonal variation of the monsoon index. The positive (negative) phase of the summer (winter) oscillation variation is significant, while those of spring and autumn are relatively weak, reflecting the difference in the aerosol type and source and the optical–radiation properties in different seasons. Except for the influence of other factors, it is mainly affected by the variation in the monsoon circulation activities (Figures A1 and A3). The interannual variation difference of the monsoon circulation can also be reflected in the variation of the aerosol properties. The positive phase of the EAMI and SAMI variations in the summer of 2012 is much longer than that in 2013, resulting in the difference of the variation of the aerosol properties in the summer during those two years.

**Figure 13.** The characteristics of the seasonal variation and inter–annual activities of EAMI and SAMI and aerosol properties normalized time series from March 2012 to February 2014. The aerosol properties in sub-graphs are (**a**) AOD1640, (**b**) AOD1020, (**c**) AOD870, (**d**) AOD670, (**e**) AOD500, (**f**) AOD440, (**g**) AOD380, (**h**) AOD340, (**i**) Angstrom Exponent(α), (**j**) SSA(440), (**k**) CWV, (**l**) AIRS CWV.

The linear fitting regression analysis of the aerosol properties and the monsoon index normalized sequence shows that *y* (the response of the aerosol properties to the monsoon circulation) = *kx* (the variation of the monsoon circulation). Only the influence of the monsoon circulation variation is considered here. The physical meaning of the slope *k* is the sensitivity of the aerosol properties to the monsoon circulation variation, while *kx*×% characterizes the relative influence rate of the monsoon circulation on the variation of the aerosol properties.

Figures 14 and 15 show the correlation analysis of the normalized sequence between EAMI and SAMI and aerosol properties containing CWV (AIRS CWV). The CWV and CWV (AIRS) have the best linear relationship with the monsoon index. The relationship not only reflects the transport of the water vapor from the atmospheric circulation but also verifies the reliability of the aerosol properties retrieval. The linear relationship between the AOD at 8 bands, SSA at 440 nm, CWV (AIRS CWV), EAMI, and SAMI all pass the significance test of more than 99%, but the linear relationship between the Angstrom Exponent (α) and the monsoon index is relatively insignificant.

**Figure 14.** The correlation analysis between the standardized sequence of the EAMI and aerosol properties: (**a**) AOD1640, (**b**) AOD1020, (**c**) AOD870, (**d**) AOD670, (**e**) AOD500, (**f**) AOD440, (**g**) AOD380, (**h**) AOD340, (**i**) Angstrom Exponent(α), (**j**) SSA(440), (**k**) CWV, (**l**) AIRS CWV. Except for the Angstrom Exponent (α) (which passed 95% for the significance test), all passed the significance test with more than 99%.

**Figure 15.** Correlation analysis between the normalized sequence of the SAMI and aerosol properties: (**a**) AOD1640, (**b**) AOD1020, (**c**) AOD870, (**d**) AOD670, (**e**) AOD500, (**f**) AOD440, (**g**) AOD380, (**h**) AOD340, (**i**) Angstrom Exponent(α), (**j**) SSA(440), (**k**) CWV, (**l**) AIRS CWV. The except for the Angstrom Exponent (α) (which only passed 69% of the significance test), all passed the significance test with more than 99%.

Table 2 shows the sensitivity *k*<sup>E</sup> of the aerosol properties to the EAMI variation and the influence rate (%) of the variation of the East Asian monsoon circulation on the relative variation of the aerosol properties.


**Table 2.** Sensitivity *k*<sup>E</sup> of the aerosol properties to the EAMI variation and the influence rate (%) of the East Asian monsoon circulation variation on the relative variation of the aerosol properties.

During the dry season, from November to March (or April) of the following year, the aerosol characteristic parameter value decreases with the weakening of the East Asian monsoon circulation (negative phase) and reaches the minimum value in January. During the spring/summer transition from March (or April) to May, the positive phase of the East Asian monsoon circulation variation gradually increases. At this time, the variation of the aerosol characteristic parameter value under the influence of the monsoon circulation starts to change from decreasing to increasing. During the rainy season from May to October, the values of the aerosol properties increase rapidly in the period of the East Asian monsoon circulation enhancement (positive phase). It should be noted that the aerosol properties reached a maximum value in August 2012. However, the maximum value not only was delayed by two months but was relatively small. It shows that the interannual activity of the East Asian monsoon circulation is significantly different. The negative phase of the East Asian monsoon circulation variation gradually increases during the autumn–winter transition from October to November. At this time, the variation of the aerosol characteristic parameter value begins to change from increasing to decreasing due to the variation of the monsoon circulation. The transformation of the positive and negative (negative and positive) phases of the monsoon circulation affects the seasonal variation characteristics of the aerosol properties. However, the activity intensity of the monsoon circulation in different years and the difference in the transition period have different effects on different aerosol properties. Except for the CWV (AIRS CWV), the effect of the SSA on the 440 nm band is the most significant, and the effect on the AOD at different wavelengths is more consistent.

Table 3 shows the sensitivity *k*<sup>S</sup> of the aerosol properties to the SAMI variation and influence rate (%) of the South Asian monsoon circulation on the relative variation of the aerosol properties. Compared with the East Asian monsoon circulation, the South Asian monsoon circulation weakens 1 month earlier. From October to April of the following year, the South Asian monsoon circulation weakens, and the aerosol characteristic parameter values decrease accordingly, also reaching a minimum value in January. During the transition period from April to May, the monsoon circulation has little effect on the variation of the aerosol properties due to the aerosol characteristic parameter values' change from decreasing to increasing.


**Table 3.** Sensitivity *k*<sup>S</sup> of the aerosol properties to the SAMI variation and the influence rate (%) of the South Asian monsoon circulation on the relative variation of the aerosol properties.

From May to September, the South Asian monsoon circulation increases, and the aerosol characteristic parameter values continued to increase, reaching a maximum value in August. During the transition period from September to October, the South Asian monsoon circulation begins to weaken, which changes the aerosol characteristic parameter value from increasing to decreasing and the impact of the variation of the monsoon circulation on the aerosol characteristic parameter values is also small.

By combining Tables 2 and 3, it is clear that the effects of the East Asian monsoon and the South Asian monsoon on the variation in aerosol properties are similar. Although there is a seasonal difference, the synergy between the East Asian monsoon and the South Asian monsoon is obvious, and the variation in the aerosol properties is consistent with their seasonal variations. The sensitivity of the CWV and CWV (AIRS) to the EAMI and SAMI variations is 0.83–0.90 and 0.93–0.95, respectively. The CWV (AIRS CWV) is more significantly affected by the South Asian monsoon than the East Asian monsoon, which indicates the contribution of the water vapor from the atmospheric circulation. The sensitivity of the SSA in the 440 nm band to the EAMI and SAMI is 0.72 and 0.77, while the sensitivity of the AOD for different λ to the EAMI and SAMI ranges from 0.58 to 0.66 (mean value is 0.62) and 0.48 to 0.69 (mean value is 0.56), which indicates that the extinction of the aerosol particles is deeply affected by the variation of the monsoon circulation activity. It is also noted that the sensitivity of the AOD to the EAMI is greater than that of the SAMI, and the influence of the inter–monthly variation of the EAMI on the aerosol properties is more drastic than that of the SAMI. This may be due to the complexity of the East Asian monsoon and the South Asian monsoon in the process of converging over the Yunnan–Kweichow Plateau, which needs to be studied in detail.

In summary, aerosol properties that vary with the seasons are affected and controlled by the variation in the monsoon circulation activity. For example, the AOD is closely related to the seasonal variation of the CWV. The increase of the water vapor content can affect the variation of the aerosol particle size and quantity, and the physical and chemical reactions of aerosol particles, which promotes the increase of the AOD caused by the extinction of the aerosol. In 2012 and 2013, the variation difference of the monsoon circulation can affect and control the difference in the aerosol sources and types, and the variation in the intensity of the seasonal and interannual activity may also change the number of aerosol particles and the formation difference of the aged aerosol. However, since only two years of data were used, it is necessary to verify whether it is credible through the evaluation of longer–term data. In addition, it should be noted that except for the influence of the monsoon circulation, there is an impact from the increased local emissions of fugitive dust on the aerosol properties. For example, the differences between the interannual and monthly variation in the aerosol properties in Kunming in 2012 and 2013 are not only related to the interannual differences in monsoon circulation but also related to the increase in the local emissions.

#### **7. Conclusions**

The relationship between the influence of the East Asian and South Asian monsoon index on the seasonal variation of the aerosol properties is discussed below based on the retrieval of the aerosol optical and radiation properties from the observational data of the CE–318 sunphotometer at the Kunming atmospheric ozone monitoring station. The conclusions are as follows.

The AOD decreases with the increase of λ and is consistent with the typical Mie scattering properties. The seasonal variation of the AOD, Angstrom turbidity coefficient β*,* and CWV is unimodal, the summer (winter) is the largest (smallest), and the spring is greater than the autumn. The CWV values in spring and autumn showed no significant difference, and the seasonal variation of the Angstrom Exponent (α) was not obvious. The turbidity in the atmosphere was relatively low. This indicates that the aerosol properties of the low–latitude plateaus in southwestern China are different from those in other parts of China. Although they are affected by anthropogenic emissions from the Southeast Asia–South Asia region and biomass burning aerosols, the correlation between the intensity of the seasonal activities of the monsoon circulation and the local anthropogenic emissions and meteorological factors will affect the aerosol properties, so further in–depth simulation studies are needed.

The seasonal and annual variations of the aerosol particle size distributions are bimodal. The complex refractive index and single scattering albedo showed that the summer aerosol particles had a stronger scattering power and a relatively weak absorption power; in winter, they exhibited the opposite. In summer, the water vapor content was high, and after the hygroscopic aerosol particles absorbed water, the particle size and volume increased, resulting in enhanced scattering. In order to explain the enhancement of the aerosol absorption capacity in winter, the retrieval result for the FI needs to be further discussed.

The variation of the aerosol properties with the season was significantly affected by the monsoon circulation activity, and the synergistic effect of the East Asian monsoon and South Asian monsoon made the variation of the aerosol properties have the same seasonal oscillation characteristics as that of the monsoon. In the transition period of the season (especially during the transition periods of the dry and wet seasons in Kunming), the monsoon circulation had little effect on the variation of the aerosol properties. The AOD was closely related to the seasonal variation of the CWV. The increase or decrease of the water vapor content did affect the variation of the aerosol particle size and quantity, thus affecting the physical and chemical reactions of the aerosol particles. Furthermore, the extinction of the aerosol was promoted to cause an increase in the AOD. Only in 2012 and 2013 did the variation difference of the monsoon circulation affect the source and type of aerosol. The variation in the intensity of the seasonal and interannual activity may have also changed the number of aerosol particles and the formation difference of the aged aerosol. It should also be noted that, in addition to the influence of the monsoon circulation, there was an impact from the increased local emissions of fugitive dust on the aerosol properties.

#### **8. Discussion**

This study only dealt with the retrieval and discussion of the aerosol properties in Kunming. Follow–up work needs to compare the aerosol properties of the retrieval between the satellite data and CE–318 observation data. We need to simulate the seasonal variation for the monsoon circulation activity variation and aerosol type and source in greater depth, as well as the direct radiative forcing

effects of aerosols. This study investigated a limited time span of the observational data from the Kunming atmospheric ozone monitoring station, so the rule of seasonal change needs data from a longer period of time, and the long–term evolution trend of the interannual and interdecadal variations needs further research. Furthermore, a more in–depth statistical analysis of the relationship between the aerosol properties and meteorological elements is needed.

**Author Contributions:** Conceptualization, H.C., B.Z., H.W. and W.W.; Data curation, C.Z., X.T., J.B., S.W. and R.L.; Formal analysis, K.Y., Z.Y., X.C. and Y.L.; Methodology, Y.W.; Software, X.D. and X.M.; Original draft, H.W. and C.Z.

**Funding:** This work was supported by the "Second Tibetan Plateau Scientific Expedition and Research Program (STEP)" (grant no. 2019QZKK0604), the National Science Fund for Distinguished Young Scholars and the National Natural Science Foundation of China (grant nos. 41825011, 41807308, 21777026, 41641044, and 21477021), and the Open Project of Shanghai Key Laboratory of Atmospheric Particle Pollution and Prevention (LAP3) (grant no. FDLAP19009).

**Acknowledgments:** We want to thank ESR for the ESR.pack retrieval program, ECMWF for the ERA–Interim reanalysis data, the University of Maryland for the MODIS fire point product data, NASA for the AIRS water vapor observation data, the Kunming National Meteorological Base Station for the ground meteorological elements, and the China National Climate Center for the EAMI and SAMI data.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Appendix A**

**Figure A1.** Seasonal variation of the 700 hPa wind field and relative humidity field from 2012 to 2013; the red point is the location of Kunming, the colored area is the Relative Humidity (%), and the reference wind vector is 5 m·s<sup>−</sup>1.

**Figure A2.** The emissions of the seasonal biomass burning (forest fires, straw burning, and burning) in South Asia–Southeast Asia from 2012 to 2013; the red part is the number of fire points for biomass burning. Refer to the color bar on the right side of the figure for more information.

**Figure A3.** The 72-h backward trajectory of the Hybrid Single–Particle Lagrangian Integrated Trajectory (HYSPLIT) mode of Kunming at 750 hPa from March 2012 to February 2014.


**Table A2.** The frequency of the biomass burning fire points in different seasons in Southeast Asia and South Asia (0–43◦N, 59–140◦E) from 2012 to 2013.


**Table A3.** Characterization of Angstrom's Atmospheric Turbidity Coefficient β [63].


**Table A4.** Mean statistical results of the complex refractive index of aerosol in different compositions\* [58,64–66].


**\*** Combined with the spectral properties of the complex refractive index of the aerosol, it is assumed that the imaginary parts at 870 nm and 1020 nm are approximately equal to those at 670 nm.

**Table A5.** Statistical relationship between the Angstrom exponent (α) and mean radius *r* of the aerosol particles [59].


#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **New Global View of Above-Cloud Absorbing Aerosol Distribution Based on CALIPSO Measurements**

**Wenzhong Zhang 1,2, Shumei Deng 3, Tao Luo 1,\*, Yang Wu 1,4, Nana Liu 1,4, Xuebin Li 1, Yinbo Huang <sup>1</sup> and Wenyue Zhu <sup>1</sup>**


Received: 3 September 2019; Accepted: 12 October 2019; Published: 16 October 2019

**Abstract:** Above-low-level-cloud aerosols (ACAs) have gradually gained more interest in recent years; however, the combined aerosol–cloud radiation effects are not well understood. The uncertainty about the radiative effects of aerosols above cloud mainly stems from the lack of comprehensive and accurate retrieval of aerosols and clouds for ACA scenes. In this study, an improved ACA identification and retrieval methodology was developed to provide a new global view of the ACA distribution by combining three-channel CALIOP (The Cloud–Aerosol Lidar with Orthogonal Polarization) observations. The new method can reliably identify and retrieve both thin and dense ACA layers, providing consistent results between the day- and night-time retrieval of ACAs. Then, new four-year (2007 to 2010) global ACA datasets were built, and new seasonal mean views of global ACA occurrence, optical depth, and geometrical thickness were presented and analyzed. Further discussion on the relative position of ACAs to low clouds showed that the mean distance between the ACA layer and the low cloud deck over the tropical Atlantic region is less than 0.2 km. This indicates that the ACAs over this region are more likely to be mixed with low-level clouds, thereby possibly influencing the cloud microphysics over this region, contrary to findings reported from previous studies. The results not only help us better understand global aerosol transportation and aerosol–cloud interactions but also provide useful information for model evaluation and improvements.

**Keywords:** above-cloud aerosol; low-level cloud; CALIPSO

#### **1. Introduction**

The long-range transport of aerosols plays an important role in several regions of the world, having a potential impact on aerosol–cloud interactions, atmospheric chemistry, and air quality [1–8]. In particular, aerosols often overlay lower level clouds [9], for example, biomass burning aerosols and wind-blown dust overlay low-level cloud deck over the Atlantic. The above-low-level cloud aerosols (ACAs) occupy about 25% of the mean aerosol optical depth (fine mode) at a global scale [9], and this fraction could be much higher regionally and seasonally [10]. Current models experience significant inter-model discrepancies in aerosol forcing assessments, especially over the aerosol–cloud overlap regions [11], that result from inter-model differences in both aerosol and cloud properties [11,12]. A recent evaluation showed that most models cannot reproduce the observed large aerosol load episodes [13]. Therefore, improved observations are needed to better understand ACA–cloud interactions and to constrain the aerosol–cloud radiative processes in models.

The ACA has gained increasing attention in recent years because of its important, but not well-understood, radiative and microphysical effects on clouds [14–19]. Different from clear-sky aerosols, which are usually associated with a negative (cooling) direct radiative effect [15,19–23], the absorption effect of the ACA is significantly amplified due to the strong reflection light of low clouds, resulting in a less negative or even positive (warming) direct radiative effect. By warming the free troposphere and cooling the surface below, the ACA can increase the low cloud cover by enhancing atmospheric stability [24]. However, if the ACA mixes directly with the cloud layer, warming of the ACA could reduce the relative humidity, dissipating the cloud [24]. Therefore, the radiative and microphysical effects of ACAs on clouds depend on both their loadings and their positions relative to the clouds [25].

Despite its importance, quantifying the effect of the ACA on the clouds from satellite observations is still challenging. For passive remote sensing, the ACA has usually been neglected, and only clear-sky aerosol concentrations can be derived from passive satellites [26,27]. This is because passive satellites employ the reflection of natural sunlight to retrieve aerosol and cloud properties separately. The neglect of the ACA in passive remote sensing can result in uncertainties in the retrieval of cloud micro- and macro-properties, such as the liquid water path, cloud optical thickness, and effective radius of cloud droplets [26,27]. Until recently, based on different absorption effects of the ACA in the visible and near-infrared channels, a "color ratio" (CR) method was used to retrieve the ACA optical depth (ACAOD) and aerosol-corrected cloud optical depth (COD) [26] and this was further improved to a multi-channel method to simultaneously retrieve ACAOD, COD, and cloud effective radius data [27]. The CR between a pair of wavelengths is a function of both the aerosol and cloud optical thicknesses, and the measured reflectance can be related to pairs of aerosol and cloud optical thicknesses. The results indicate that the mean liquid cloud optical thickness can be increased by roughly 6%, and the mean liquid effective radius can be increased by roughly 2.6% after correcting for the effect of the ACA [22,28–32]. An inter-satellite comparison of the ACAOD retrieved from NASA's A-train sensors revealed a good level of agreement between the passive sensors over homogeneous cloud fields [33]. However, passive satellites can only provide daytime ACAOD and are unable to determine the spatial position of aerosols with respect to the clouds.

Other than passive sensors, the active lidar Cloud–Aerosol Lidar with Orthogonal Polarization (CALIOP) onboard the Cloud–Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) can provide valuable vertically resolved information about aerosols and clouds. CALIPSO Level-2 data have already been employed in studies of ACA climate effects [19,22,23,34–36]. It was shown that the vertical structure of ACA and clouds was critical in determining the aerosol–cloud radiative forcing. However, inter-satellite comparisons of ACAOD [33,37,38] showed that the CALIPSO-derived ACAOD data was consistently lower than those derived from other satellite sensors in the A-Train, by a factor of four to six, as compared with passive retrieval. The main reason for this is that the current CALIPSO Level-2 aerosol retrieval algorithm employs the 532 nm channel, which cannot detect the true aerosol layer base because of the strong attenuation by the ACA in this channel [33,38]. Based on the CALIPSO Level-2 product, most of the ACA resides above the cloud deck at a distance of about 1 km over the tropical Atlantic region where both the ACA and low cloud frequently occur, indicating that weak cloud microphysical effects occur due to the aerosols not mixing with the cloud [24]. However, this could be misleading because the detected ACA layer base in the CALIPSO Level-2 product has been suggested to be higher than the real position [33,38]. Some other issues also limit the use of CALIPSO Level-2 data in ACA studies, such as the poor consistency between day- and night-time retrieval due to the poor signal-to-noise-ratio (SNR) performance of CALIPSO daytime observations.

This paper aims to provide a new global view of the ACA distribution based on an improved ACA detection method by combining the CALIPSO 532 and 1064 nm channels. Furthermore, efforts to address the low SNR issue in CALIPSO daytime observations were made to minimize the difference between day- and night-time retrieval. The data sources used in this study are introduced in Section 2. The new ACA identification and retrieval method is detailed in Section 3. Section 4 presents the results and discussions, and the conclusions are summarized in Section 5.

#### **2. Data**

CALIPSO, CloudSat and operational meteorology datasets were used to perform the study. CloudSat has not been able to operate at night-time and has operated in daylight-only mode since a spacecraft battery anomaly occurred in 2011. Therefore, the 2007–2010 observations were employed in this study to build both day- and night-time ACA datasets.

The CALIPSO Level-1B Version 4 product provides 532 nm total attenuated backscatter (β 532), a perpendicular polarization component (β 532*p*), and 1064 nm total attenuated backscatter (β 1064) [39], with a 30 m vertical resolution below 8.2 km and a 1/3 km resolution in the horizontal direction along the ground track. The photomultiplier tube exhibited a nonideal recovery at 532 nm after encountering a strong backscattering objective, which was corrected following Li et al. [40]. The CALIPSO Level-2 Vertical Feature Mask (VFM) and Aerosol Profile Products (APro) files were used for comparison in this study.

The cloud type, cloud top height, and cloud bottom height were provided by the 2B-CLDCLASS-Lidar product [41]. The 2B-CLDCLASS-Lidar product used combines the CloudSat Cloud Profiling Radar (CPR), which is a 94 GHz microwave radar, and CALIPSO lidar measurements and thus can provide more accurate and complete cloud mask information than radar-only or lidar-only measurements [42].

The meteorological reanalysis data MERRA2 (The Modern Era Retrospective-analysis for Research and Application, Version 2 [43]) assimilated the meteorological data using a modern satellite database, which was released by the Global Modelling and Assimilation Systems of Goddard Space Flight Centre of NASA. It was used to provide temperature and pressure profiles for molecular backscattering estimation in this study.

#### **3. Methodology**

Previous studies have shown that the CALIPSO operational aerosol product tends to miss the bottom of dense ACA layers and underestimates the 532 nm aerosol optical depth because of the relatively strong aerosol attenuation in the 532 nm channel [33,37,44]. In contrast, both the aerosol attenuation and molecular backscattering in the 1064 nm channel are smaller than those in the 532 nm channel. Therefore, CALIPSO 532 and 1064 nm channel lidar observations were combined to develop an ACA retrieval methodology in this study, as detailed below. A flow chart of the methodology is shown in Figure 1.

**Figure 1.** Flow chart of the new above-low-level cloud aerosol (ACA) identification and retrieval methodology.

Firstly, daytime CALIPSO observations experience lower SNRs than night-time observations, resulting in a strong apparent day- and night-time difference in the ACA global distribution (refer to Figure 3 in Alfaro-Contreras et al. [45]). To overcome the poor SNR issue in daytime CALIOP signals, different smoothing scales were selected for day- and night-time observations to provide consistent retrieval. For the β <sup>532</sup> and β <sup>1064</sup> night-time observations, 20 km horizontal averaging was adopted. For the β <sup>532</sup> and β <sup>1064</sup> daytime observations, multiscale averaging in different signal-strength bins was selected according to the night-time noise level. As shown in Figure 2, daytime noise was calculated with horizontal averaging of 20, 40, 60, 80, and 100 km and compared to the 20 km night-time results. According to Figure 2, β <sup>532</sup> daytime signals of higher than 2.3 <sup>×</sup> 10−<sup>2</sup> km−1sr−1, between 2.3 <sup>×</sup> 10−<sup>2</sup> and 2.8 <sup>×</sup> 10−<sup>3</sup> km−1sr−1, and lower than 2.8 <sup>×</sup> 10−<sup>3</sup> km−1sr−<sup>1</sup> were exhibited at 20, 40, and 60 km, respectively. The maximum averaging scale was selected to be 60 km, after which the SNR performance improved a little. To avoid involving cloud signals in the moving-smoothed aerosol signals, only non-cloud data were smoothed using the cloud mask from the CloudSat 2B-CLDCLASS-Lidar product.

**Figure 2.** Comparison of noise between the 20 km moving–smoothing in the night-time and multiscale moving–smoothing in the daytime.

Secondly, after multiscale smoothing, ACA profiles were screened out using the CloudSat 2B-CLDCLASS-Lidar product and the estimated two-way transmittance at 1064 nm. Low cloud was defined as a cloud top lower than 3 km above ground level (AGL). Only single-layer low cloud or multi-layer thin cirrus cloud above the low cloud cases were considered in this study. Low-level cloud sample quantity is presented for each 2.5◦ × 2.5◦ grid box and averaged for winter (DJF, i.e., December, January, and February), spring (MAM, i.e., March, April, and May), summer (JJA, i.e., June, July, and August), and autumn (SON, i.e., September, October and November) in Figure 3. For each low cloud profile, the two-way transmittance at 1064 nm within 6 km above the cloud top was estimated by ignoring the molecular backscattering, as

$$\mathbf{x} = \mathbf{e}^{-2\zeta} = 1 - 2\mathbf{S} \int \beta \nu\_{1064} \mathbf{d} \mathbf{r} \tag{1}$$

Here, <sup>x</sup> is the two-way transmittance, and β dr is the integrated attenuated backscatter at 1064 nm within 6 km above the cloud top, and S is the lidar ratio at 1064 nm. For cirrus cloud, S was chosen to be 19 sr. For aerosol, the S at 1064 nm has a smaller range than that in the 532 nm channel [46] and was assumed to be 40 sr in this step. By considering the noise, e−2∗0.015 = 0.97 was chosen as the threshold, and values smaller than that were identified as possible ACA cases.

**Figure 3.** Low cloud samples quantity in (**a**) winter (DJF, i.e., December, January, and February), (**b**) spring (MAM, i.e., March, April, and May), (**c**) summer (JJA, i.e., June, July, and August), and (**d**) autumn (SON, i.e., September, October, and November).

Thirdly, following the screening-out of possible ACA cases, the first-guess ACA extinction at 1064 nm was retrieved with the forward integration scheme [47] to build the ACA mask. In this step, the 1064 nm S was chosen to be 40 sr. The first-guess ACA mask was set to unity when the first-estimate ACA extinctions were larger than four times the measured noise, indicating the possible presence of aerosols. Then, the first-guess ACA top and bottom were identified as the highest and lowest points above the cloud top where the ACA mask equals 1.

Finally, using the first-guess ACA top and bottom from the last step, the ACA extinctions from the 532 nm channel were retrieved, with the lidar ratio determined according to the aerosol type classification, identical to the method used by Omar et al. [46]. The ACA mask was thus further refined according to the 532 nm extinctions using a threshold of four times of the measured noise.

Figure 4 presents comparisons of an ACA case derived from the new method and the CALIPSO level-2 V4 aerosol and cloud product. This is an outflowed smoke case immediately above the top of low marine clouds, as can be seen from the smoothed β<sup>532</sup> and data in Figures 4b and 4c, respectively. It should be noted that the same case was demonstrated in Jethva et al. [33]. The retrieval from the new method (Figure 4e) was close to that from passive methods and the CALIOSP "DR" and "CR" methods (refer to Figure 2 in Jethva et al. [33]). In contrast, the ACAOD from the CALIPSO Level-2 product was much smaller than that from the new method, with a strong underestimation by a factor of up to about 2. This is because the standard CALIOP Level-2 aerosol and cloud products use the 532 nm signal to detect respective layers, which, in the presence of thick smoke layers, can be incorrectly assigned due to strong attenuation at 532 nm. Figure 4h shows that only the upper portion of the layer with dense aerosols was detected, while most of the lower portion of the ACA was missed. We tested the projection of a new extinction (Figure 4f) to the CALIOP Level-2 aerosol mask in Figure 4h, and this produced a similar ACAOD to L2. Benefitting from the aid of the 1064 nm extinction coefficient, the new method was shown to provide a complete aerosol mask. Therefore, the new method can provide more accurate ACA detection and retrieval results than the CALIPSO Level-2 product.

**Figure 4.** ACA case on 2 August 2007, 13:12:11 UTC: (**a**) geo-location; (**b**) β532; (**c**) β1064; (**d**) ACA mask from the new method in this study; (**e**) aerosol optical depth at 532 nm determined using the new method (red) and CALIPSO Level-2 V4 APro files (blue); (**f**) extinction coefficient of ACA at 532 nm; (**g**) extinction coefficient of ACA at 1064 nm; (**h**) ACA and cloud layer from the CALIPSO Level-2 VFM product.

#### **4. Results**

The new methodology was applied to the 2007–2010 observations to build a new global ACA dataset. Statistical results are also presented for each 2.5◦ × 2.5◦ grid box and averaged for DJF, MAM, JJA, and SON. Section 4.1 presents a comparison of day- and night-time cloudy-sky ACA occurrences (ACAOC) and ACAOD, and Section 4.2 analyses the global distribution of the seasonal-mean ACAOC and ACAOD with different aerosol types. The ACAOC is defined as the number of ACA profiles divided by the number of low cloud profiles in each 2.5◦ grid box.

#### *4.1. Day- and Night-Time Comparison of the New ACA Dataset*

Figure 5 presents a comparison between night- and day-time global annual mean distributions of low-level cloud fractions (left column), ACAOC (middle column), and ACAOD (right column) derived from the new method using new four-year global ACA datasets. Table 1 summarizes the night- and day-time statistics of the annual mean ACAOC and ACAOD in different seasons globally and over the dust region (5–30◦ N, 60–16◦ W), smoke region (22◦S–5◦ N, 18◦ W–15◦ E), and Eastern Asia (19–40◦ N, 100–140◦ W), respectively. These places have also been documented in previous studies [28,48–51].

As shown in Figure 5 and Table 1, the ACAOC and ACAOD show very similar patterns and values between day- and night-time retrieval. The global mean day- and night-time ACAOC are 0.125 and 0.108 respectively, and the global mean ACAOD values are 0.146 and 0.143. The correlation coefficient between day- and night-time is 0.938 for ACAOC and 0.796 for ACAOD. This shows that the new method developed in this study can produce reasonable day- and night- time results, better than the prior results derived from the CALIPSO level-2 product (refer to Figure 3 in Alfaro-Contreras et al. [45]). The night-time results from the CALIPSO level-2 product are quite similar to the results produced in this study. However, in the day-time, the CALIPSO level-2 products missed lots of weak aerosol layers due to their poor SNRs and thus gave quite different global distributions of ACAOC and ACAOD from those of night-time CALIPSO level-2 product and our method. The new method described in

this study can detect more ACA cases over both source regions and long-range transport regions. Furthermore, the ACAOD values derived by the new method are larger than those from the CALIPSO Level-2 product. The CALIPSO Level-2 product retrieves ACAOD values that are too low over the tropical Atlantic region, a region where other passive methods suggest that high above-cloud aerosol loading could exist (i.e., refer to Figure 2 in Devasthale et al. [52]).

**Figure 5.** Day- and night- time comparison of the global annual mean distribution of low-level cloud fractions (**a**(1)-(2)), cloudy-sky ACA occurrences (ACAOC, **b**(1)-(2)), ACA optical depth (ACA optical depth (ACAOD), **c**(1)-(2)), and (**d**) is the global annual mean distribution of OMI daytime ACAOD.


**Table 1.** ACAOC and ACAOD data from different regions in the night- and day-time.

A preliminary comparison between our ACA dataset and passive satellite product was also done. Figure 5d shows the global distribution of annual-mean daytime 500 nm ACAOD derived from OMACA product version 3, which retrieves ACAOD from OMI's two near-UV observations (354 and 388 nm) [17]. The data processing of Figure 5d follows the same way as Figure 8 in Jethva et al. [10], except that only 2007–2010 data were used in this paper. Comparing to the daytime ACAOD result of

this paper (Figure 5c(1)), it shows that both products have quite similar pattern of global distribution to each other, while with regional differences in values. Our ACA dataset retrieves larger ACAOD than OMACA product over some major aerosol source regions, such the main dust band region including Sahara and Middle-Eastern dust regions and related long-range transport regions (tropic Atlantic, Arabian Sea, and Southern Asia), central Africa wildfire region and eastern Asia region except southernmost part of China. In contrast, OMACA product retrieves larger ACAOD than our ACA dataset over some weak aerosol source regions, such as Siberia, Alaska and Southern Oceans. Further detailed evaluation of the ACA dataset developed in this paper with passive satellite products and its possible implements in improving passive ACA and cloud retrievals will be our study interest in the near future.

#### *4.2. Global Distribution of the Seasonal Mean ACA Properties*

Figures 6 and 7 illustrate the global distribution of the seasonal mean ACAOC and ACAOD data for all aerosol types (left column), dust aerosols (pure and pollute dust, middle column), and smoke and polluted continental aerosols (right column). Smoke and polluted continental aerosols were analyzed together in these figures, because these two types of aerosol have similar optical properties and are quite difficult to discriminate each other for the ACA case [46]. Furthermore, the existing classification algorithm is not suitable for the detection of above-cloud marine aerosols. As can be seen in those figures, above-cloud marine aerosol occurs over the western tropical Pacific region but is misclassified as smoke or polluted continental aerosols. The global mean ACAOC in each season was found to be 13% (DJF), 20% (MAM), 18% (JJA), and 15% (SON), and the global mean ACAOD was 0.12 (DJF), 0.14 (MAM), 0.15 (JJA), and 0.14 (SON), respectively.

As shown in Figures 6 and 7, ACA frequently occurs near the source regions, such as over the Sahara and Middle-Eastern dust regions (main dust band region), the Africa smoke region, and the Eastern Asia region. Those aerosols are mobilized far from the sources and travel along transoceanic pathways, such as dust transport over the tropical Atlantic, region and smoke transport over the southeast Atlantic region, resulting high ACAOc values there. For convenience, the ACA and its long-distance transport are discussed with regard to its main source as in the following text.


over the Pacific Ocean. The above-cloud smoke and polluted continental aerosols over Eastern Asia are also most active in MAM, representing ~32% of occurrences and having an AOD of 0.31 above low clouds. These aerosols can also be transported over the Pacific and contribute to the ACA presence there. The occurrence of above-cloud smoke and polluted continental aerosols is ~25% and the ACAOD value is 0.23 over the Pacific. In other seasons, the ACA occurrence is about 7–25%, and this is associated with very weak long-range transport.

**Figure 6.** Global distribution of the seasonal mean ACAOC for all aerosol types (left column), pure dust and polluted dust aerosols (middle column) and smoke and polluted continental aerosols (right column).

**Figure 7.** Global distribution of the seasonal mean ACAOD for all aerosol types (left column), pure dust and polluted dust aerosols (middle column) and smoke and polluted continental aerosols (right column).

Other than these major aerosol source regions, our results also indicate some weak aerosol source regions, such as South America, Central America, and Timor Sea, where ACA frequently occurs in some seasons, but this has rarely been studied in previous studies. South America and the Timor Sea can be considered to be anthropogenic source regions due to the biomass burning emissions [53,54]. During SON, the above-cloud smoke aerosol occurrence over South America is about 24% with an AOD value of about 0.16. Trade winds will carry the smoke aerosol produced from Northern Australian savannah fires and agricultural/peat burning fires in Indonesia over the Timor Sea [54], resulting in an above-cloud smoke aerosol occurrence of 57% with an ACAOD value of 0.20. The ACA over Central America mostly occurs in MAM as a result of both anthropogenic sources due to human activity and biomass burning activities [55]. The occurrence of above-cloud polluted continental and above-cloud smoke aerosol in Central America is about 35% with an AOD value of about 0.21 in MAM. Our results also indicate that smoke from southeastern Africa transports in an easterly direction across the Indian Ocean in SON, and may reach the west coast of Australia.

#### **5. Discussion**

The relative positions of the ACA and low clouds determine the response of low clouds to the ACA [24]. Figure 8 shows the 4-year zonal mean ACAOC (Figure 8a), ACAOD (Figure 8b), vertical distance between the ACA and cloud (VDAC, Figure 8c), and ACA geometrical thickness (AGT, Figure 8d). The ACAOC and ACAOD values have quite large inter-season variations at low latitudes, which are related to the dust and smoke activity. At low latitudes (between –40◦ and 40◦N), the ACAOC can reach about 40%, and the ACAOD can reach 0.2–0.4 during high ACA activity seasons (i.e., MAM in the northern tropical region and SON in the southern tropical region). Additionally, the VDAC is smaller than about 0.5 km and the AGT can reach 1–2 km, because of the ACA being close to the main source region (Saharan dust and African smoke). When transported away from the source to high latitudes (<–40◦ and >40◦), the ACAOC decreases as the latitude increases. However, the VDAC increases from 0.5 to about 1–1.5 km during long-range transport, associated with the decrease in AGT, except in MAM. In MAM, the AGT in the northern hemisphere is generally about 2 km, as a result of the long-range transport of Asian dust at mid-latitudes [56]. The different behaviors that occur at low and high latitudes indicate different cloud responses to the ACA. The stabilization effect of the ACA dominates the cloud response at the mid-to-high-latitude regions, because the ACA generally resides away from the low clouds [24]. In contrast, at low latitudes, the ACA is more likely to directly mix with low clouds and influence the cloud microphysics.

**Figure 8.** (**a**) ACAOC line chart with latitude. (**b**) ACAOD line chart with latitude. (**c**) Vertical distance between the ACA and cloud (VDAC) line chart with latitude. (**d**) ACA geometrical thickness (AGT) line chart with latitude.

A more insightful view of the relative positions of the ACA and low clouds at low altitudes is shown in Figure 9 by the 4-year mean longitude–altitude distribution of cloud (gray contour line) and ACA (color shade) occurrences over the dust (left column) and smoke regions (right column). The latitude ranges of the dust (tropical Atlantic) and smoke (southeast Atlantic) regions are the same as those stated in Section 4.1. Similar results by the CALIPSO level-2 product can be found in Figure 1 in Zuidema et al. [24], in which most of the ACA layer was shown to reside above the cloud deck at a mean distance of about 1 km [24]. According to the CALIPSO level-2 product, a semi-direct effect (atmospheric stabilization, reference) due to absorption in the aerosol layers above the cloud deck could dominate the ACA–cloud interaction over these regions. However, our results show that the bottom of the high ACA occurrence layer resides immediately above the most frequent cloud height (Figure 9). Our findings are different to those derived from the CALIPSO level-2 product. The main reason for this is the misdetection of the ACA layer bottom in the CALIPSO level-2 product (as stated in Section 3), which was also reported by Rajapakshe et al. (2017) using two seasons of night-time CATS 1064 nm observations over the southeast Atlantic region [44]. Therefore, as well as the semi-effect, smoke and dust layers can mix into the low cloud over the transport region of the Atlantic and influence the cloud microphysics (indirect effect) due to its proximity to the cloud layer [34].

**Figure 9.** Positional relationship between the aerosol and cloud layers in the dust and smoke area with longitude.

#### **6. Conclusions**

This paper presented an improved ACA identification and retrieval methodology that combines CALIOP 532 and 1064 nm observations. The ACA was mainly detected with the 1064 nm channel by taking advantage of the weaker aerosol attenuation compared with that of the 532 nm channel. The selection of the 1064 nm signal for ACA identification allowed the detection of the full column of aerosols above the clouds, which the 532 nm signal misses due to its stronger attenuation. Effort was made to address the low SNR issue in the CALIOP day-time observation to provide consistent results between day- and night-time retrieval. Another feature of this methodology is the reliable cloud–aerosol distinction that occurs with the combination of CALIPSO and CloudSat observations.

Then, new four-year (2007 to 2010) global ACA datasets were built, and new global seasonal-mean views of ACA properties were presented and analyzed, including the occurrence of ACA and its optical properties. The results indicate that the new method can not only capture the main ACA occurrence regions, such as the tropical Atlantic region which is associated with outflows of African dust and smoke, but also weak ACA occurrence regions, such as Eastern Asia, and long-range transport regions, such as the Pacific, using both day- and night-time observations. Compared with the CALIPSO L evel-2 product, the newly derived day- and night-time global distributions of ACA properties show good agreement, benefitted from multiple moving-smoothing processes.

The relative vertical positions of the ACA and low clouds were presented and discussed. The results indicate that at high latitudes (<–40◦ and >~40◦), the ACA mainly resides above the low clouds at a distance greater than 0.5 km, indicating the semi-direct influence of the ACA on the cloud by stabilizing the atmospheric temperature profile. At low latitudes, the distance between the ACA and low clouds is smaller than 0.5 km. Especially over the Atlantic near the Saharan dust and African smoke source regions, the ACA mainly resides immediately above the clouds; this contrasts with the results derived from the CALIPSO operational product in previous studies. The ACA over the tropical Atlantic (dust region) and southeast Atlantic (smoke region) can directly mix with the cloud, influencing the lower cloud deck through both semi-effect and indirective effects. This highlights the complexity of the aerosol–cloud interactions over these regions, where the model still has difficulty reproducing the vertical aerosol structure and aerosol–cloud radiative effect [11–13].

To qualify the roles of these two effects, improve the passive cloud property retrieval, and constrain the relative model process, a height-resolved ACA database like the one developed in this study is needed. The new method developed in this study can provide a more complete and accurate global view of the ACA. Together with other satellite measurements such as MODIS [57] and CERES [58] and cloud property retrieval methods such as that described by Luo et al. [59], the radiative effects of the ACA and its influences on the macro- and micro-physics of low clouds will be further studied in future work. These results not only help us better understand global aerosol transportation and aerosol–cloud interactions but also provide useful information for model evaluation and improvements.

**Author Contributions:** Conceptualization, S.D.; Data curation, Y.W.; Formal analysis, W.Z. (Wenzhong Zhang); Funding acquisition, T.L.; Investigation, N.L.; Methodology, W.Z. (Wenzhong Zhang); Project administration, X.L.; Software, W.Z. (Wenzhong Zhang); Supervision, Y.H. and W.Z. (Wenyue Zhu); Validation, W.Z. (Wenzhong Zhang); Visualization, W.Z. (Wenzhong Zhang); Writing—original draft, W.Z. (Wenzhong Zhang); Writing—review and editing, T.L.

**Funding:** This research was funded by the Strategic Priority Research Program of Chinese Academy of Sciences, grant number XDA17010104. Tao Luo's work was also supported by the National Key Technology Support Program of China, grant number 2017YFC0209801 and Natural Science Foundation of China, grant number 41875041.

**Acknowledgments:** The authors thank the CALIPSO Team for providing data products, which were obtained from the NASA Langley Research Center Atmospheric Science Data Center. The authors are also grateful to the CloudSat team for making the data available; these data were obtained from the CloudSat Data Processing Center (http: //www.cloudsat.cira.colostate.edu). The authors would also like to thank the MERRA team for providing the data obtained from the Global Modeling and Assimilation Systems. The authors thank the OMACA team for providing data products, which were obtained from https://avdc.gsfc.nasa.gov/pub/data/satellite/Aura/OMI/V03/L2/OMACA. The authors acknowledge the anonymous reviewers and the editor for their valuable comments and suggestions.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


Solomon, S., Qin, D., Manning, M., Marquis, M., Averyt, K.B., Tignor, M., Miller, H.L., Chen, Z., Eds.; Cambridge University Press: Cambridge, UK; New York, NY, USA, 2007; pp. 19–91.


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Vertical Structures of Dust Aerosols over East Asia Based on CALIPSO Retrievals**

**Di Liu 1, Tianliang Zhao 1,\*, Richard Boiyo 1,2, Siyu Chen 3,\*, Zhengqi Lu 1, Yan Wu <sup>4</sup> and Yang Zhao 5,6**


Received: 1 February 2019; Accepted: 21 March 2019; Published: 23 March 2019

**Abstract:** The spatiotemporal and especially the vertical distributions of dust aerosols play crucial roles in the climatic effect of dust aerosol. In the present study, the spatial-temporal distribution of dust aerosols over East Asia was investigated using Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) retrievals (01/2007–12/2011) from the perspective of the frequency of dust occurrence (FDO), dust top layer height (TH) and profile of aerosol subtypes. The results showed that a typical dust belt was generated from the dust source regions (the Taklimakan and Gobi Deserts), in the latitude range of 25◦N~45◦N and reaching eastern China, Japan and Korea and, eventually, the Pacific Ocean. High dust frequencies were found over the dust source regions, with a seasonal sequence from high to low as follows: spring, summer, autumn and winter. Vertically, FDOs peaked at about 2 km over the dust source regions. In contrast, FDOs decreased with altitude over the downwind regions. On the dust belt from dust source regions to downwind regions, the dust top height (TH) was getting higher and higher. The dust TH varied in the range of 1.9–3.1 km above surface elevation (a.s.e.), with high values over the dust source regions and low values in the downwind areas, and a seasonally descending sequence of summer, spring, autumn and winter in accord with the seasonal variation of the boundary layer height. The annual AOD (Aerosol Optical Depth) was generally characterized by two high and two low AOD centers over East Asia. The percent contribution of the Dust Aerosol Optical Depth to the total AOD showed a seasonal variation from high to low as follows: spring, winter, autumn and summer. The vertical profile of the extinction coefficient revealed the predominance of pure dust particles in the dust source regions and a mixture of dust particles and pollutants in the downwind regions. The dust extinction coefficients over the Taklimakan Desert had a seasonal pattern from high to low as follows: spring, winter, summer and autumn. The results of the present study offered an understanding of the horizontal and vertical structures of dust aerosols over East Asia and can be used to evaluate the performance aerosol transport models.

**Keywords:** CALIPSO; dust top height; frequency of dust occurrence; pure dust; polluted dust; extinction coefficient

#### **1. Introduction**

Dust aerosol, one of the most important aerosol species, modifies the energy balance and the hydrologic cycle directly by absorbing and scattering solar radiation and indirectly by altering cloud microphysical properties [1–3]. East Asia, the second largest contributor of dust aerosols in the world, emits nearly 600 tons of dust particles into the atmosphere annually through the erosion of soil mainly from natural conditions and partly from human activities [4,5]. Approximately 30% of the dust particles redeposit back into emission regions, 20% are transported to eastern China, and the rest are transported to the Pacific Ocean and beyond by westerly jets which can significantly affect the Asian monsoon system [6,7]. The deposition of dust particles serves as a key mineral supplement for the marine biopshere and for remote rainforests, thus altering the global carbon cycle [8].

It is well known that the vertical structure of dust aerosols plays a crucial role in the atmospheric thermal structure and the aerosol radiative forcing [9,10]. The vertical distribution of dust aerosols largely determines their residence and transport time [11]. Dust aerosols can also uplift to the upper troposphere and travel around the world via westerlies [12]. However, previous work using aircraft, upper-air sounding balloons and ground-based platforms are limited in space and time, offering little information on the vertical structure of dust aerosols. These limitations have been partially addressed with the launch of the Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) satellite [13,14]. Huang et al. provided an analysis of the distribution of aerosols of different types on a global scale on the basis of five-year CALIPSO data [15]. Guo et al. used CALIPSO, MODIS (Moderate Resolution Imaging Spectroradiometer) and OMI (Ozone Monitoring Instrument) measurements to build three-dimensional (3D) patterns of the occurrence frequency of aerosol, dust and smoke in China [16]. Proestakis et al. described the 3D distribution of dust aerosols over southeastern Asia [17]. Nan et al. investigated the vertical distribution of dust particles for the Taklimakan Desert and for the downwind regions using CALIPSO measurement [18].

In this study, we used five years (January 2007 to December 2011) of CALIPSO Level 3 data to gain further insights into the generation, emission, transport, distribution and speciation of aerosols over East Asia. East Asia was divided into five regions to achieve a stronger understanding of the spatial distribution of dust aerosols in relation to emission sources and transport processes. For each region, the analysis was made of the frequency of dust occurrence (FDO), dust Top Height (TH), Aerosol Optical Depth (AOD), percent contribution of dusts to the total AOD (D\_AOD), light extinction coefficient, speciation (pure dust versus pollution aerosols) and horizontal dust transport flux.

A previous study conducted over our study domain concentrated mainly on individual dust storm events involving pure dust [19]. However, human health concerns and policy on pollution controls require detailed knowledge of the speciation of dust (pure dust versus dust mixed with biomass burning smoke) and the mean behaviors of total dust (sum of pure and polluted dust). Another related study documented the seasonal variations of the aerosol extinction profile and occurrence frequency using five-year Level 3 CALIPSO products, focusing on global patterns [15]. In the present study, we also used the Level 3 products, but focused on regional patterns and specifically on the contrast between pure dusts generated from dust source regions in western China and those mixed with smoke from pollution source regions in eastern China. By analyzing spatial patterns of AOD, D\_AOD and profiles of extinction coefficients of various aerosol types, we further investigated into long-distance dust transport, which is known to impact the public health of populated centers in eastern China, Korea and Japan.

Section 2 outlines the study domain, data and methodology related to the present work. The results and discussion are given in Section 3. Section 4 summarizes the main findings of the present study.

#### **2. Data and Methods**

CALIPSO was launched into a sun-synchronous orbit in April, 2006 with a repeating cycle of 16 days. The main instruments aboard CALIPSO include a Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP), a Wide Field Camera (WFC) and an Infrared Imager Radiometer (IIR). These instruments can not only retrieve the vertical profiles of clouds and aerosols, but also provide information concerning cloud-aerosol interactions and surface radiative budgets [20]. CALIPSO continuously retrieves profiles of attenuated backscatter at 532 and 1064 nm and of polarized backscatter at 532 nm in the latitude range of 82◦N–82◦S with high horizontal and vertical resolutions of 333 m and 30–60 m, respectively [21].

Three-level products are included in CALIPSO data. Level 1 (L1) products offer a large volume of raw signals of high spatial resolutions. They are classified with a cloud-aerosol discrimination (CAD) algorithm into clouds and aerosols products, with negative and positive values representing aerosols and clouds, respectively. Level 2 (L2) products consist of six aerosol subtypes classified using a scene classification algorithm [22]. The six aerosol subtypes are pure dust, polluted dust, smoke, clean continental, polluted continental and clean marine aerosols with the corresponding Lidar ratios for the retrieval of aerosol light extinction coefficient. Level 3 (L3) products provide specific monthly variables, such as monthly mean AOD and extinction coefficients, from the L2 products. Relative to data products from other space-borne sensors, the CALIPSO data offer three main advantages for the investigation of the vertical distribution and transport of dust particles. First, CALIPSO can retrieve profiles of backscatter in different atmospheric layers to estimate the vertical profiles of aerosols and clouds, can distinguish between dust and other types of aerosols based on the observed depolarization ratios, and can minimize the interference associated with prevailing surface conditions, and hence, are able to retrieve aerosol profiles where other space-borne sensors may not be able to. Second, being an active remote sensor, it can retrieve information on aerosols and clouds in both daylight and night hours [23]. Third, CALIPSO can retrieve the vertical distributions of different aerosol subtypes; such information is crucial for understanding the long-range transport of aerosols at the regional and the global scales.

L3 is a globally monthly product with a 2◦ × 5◦ latitude-longitude grid, and a vertical resolution of 60 m from the ground to 12 km with a total of 208 layers. Uncertainties exist when CALIPSO retrieves optical properties [24]. Data on extinction below the height of 180 m were excluded as a precautionary measure to avoid biases resulting from the physical interpretation of near-ground aerosols [25]. Following [15], we also only considered extinction coefficients of >0.001 km−<sup>1</sup> to ensure high levels of accuracy.

To better understand the distribution of dust aerosols over different surface types, five homogenous regions over East Asia were analyzed in the present study (Figure 1). They include the Taklimakan Desert (TD; dust source region), the Gobi Desert (GD; dust source region), northern China (NC; a region affected by anthropogenic aerosols), southern China (SC; also under the influence of anthropogenic aerosols) and Korea-Japan (KJ; under the influence of marine aerosols).

The FDO is defined as the ratio of the number of CALIPSO overpasses with dust observations to the total number of CALIPSO overpasses (including conditions of clear air and aerosols) [17]:

$$\text{FDO} = \frac{P\_{\text{dust}}}{P\_{all}} \tag{1}$$

where *Pdust* represents the number of CALIPSO dust overpasses and *Pall* denotes the total number of CALIPSO overpasses.

The seasonally averaged CALIPSO L3 dust profile Top Height (TH) is by definition the height at which dust AOD contribution (D\_AOD) aggregates to 98% [17].

The dust column concentration (Mdu) and D\_AOD (τdu) are calculated with Equations (2) and (3) [26];

$$\mathbf{M}\_{\rm du} = \left(\frac{\mathfrak{g}4\pi}{3}\right) \int r^3 n(r) dr\tag{2}$$

$$
\pi\_{\rm du} \pi\_{\rm du} = \pi \int Q(r) r^2 n(r) dr \tag{3}
$$

where *n*(*r*) represents particle size distribution of dust aerosol, ρ is dust aerosol density and *Q* is extinction index. Following other researchers [26], we obtain:

$$\mathbf{M}\_{\rm du} = 2.7 \times \tau\_{\rm du} \,\mathrm{g} \,\mathrm{m}^{-2} \tag{4}$$

according to Kaufman et al. (2005), the horizontal dust transport flux (F) is calculated with Equation (5) [26]:

$$\mathbf{F} = \mathbf{M}\_{\text{du}} \times \mathbf{W} \times \mathbf{L} \text{ g s}^{-1} \tag{5}$$

where W represents the monthly mean west wind speed (m s−1), Mdu is the monthly mean column dust concentration (g m<sup>−</sup>2), and L is longitudinal length (m). Following other researchers [26–28], we adopted wind fields of 500 hPa over the TD and GD and of 850 hPa over the NC and SC.

For the purpose of validation, the CALIPSO AOD and pure dust D\_AOD were compared with those from the Modern-Era Retrospective analysis for Research and Applications, version 2 (MERRA-2), over the same time period. MERRA-2 calculation of these parameters is driven by anthropogenic emission inventories and natural dust emission parameterization and is bias-corrected by assimilation of MODIS and AVHRR satellite AOD products [29].

**Figure 1.** Topographical map of the study domain. Homogenous regions highlighted in black boxes include: (1) the Taklimakan Desert (TD), (2) the Gobi Desert (GD), (3) Northern China (NC), (4) Southern China (SC) and (5) Korea-Japan (KJ).

#### **3. Results and Discussion**

#### *3.1. Frequency of Dust Occurrence (FDO)*

The TD is the largest desert in China and the second largest shifting desert in the world. It emits large volumes of dust particles into the atmosphere annually. Most of the emitted particles are deposited back to the ground primarily as the result of weak tropospheric winds and complex topography: the center of the Tarim Basin is surrounded by mountains toward the south, north and west [30]. Dust episodes are less frequent but more severe over the GD. Most dust particles emitted by the GD are easily uplifted to the troposphere and transported downwind primarily due to the plateau topography (910–1520 m) of the GD. Additionally, as the main desert region of East Asia, the GD includes stationary deserts, such as the Tengger Desert and Badain Jaran Desert, and shifting deserts, such as the Mu Us Sandland. Over the NC and SC, abundant anthropogenic aerosols are mixed with pure dust through long-range transport, forming polluted dusts [31]. The KJ region suffers not only from polluted dust transported by westerlies from upwind regions (e.g., the NC and SC), but also from locally-formed marine aerosols. To understand the generation and transport of dust aerosols over East Asia, we divided the study period into four seasons, March–April–May (MAM), June–July–August (JJA), September–October–November (SON) and December–January–February (DJF), based on prevailing climatological conditions observed over the study domain.

FDO values show significant levels of spatiotemporal heterogeneity, in order from high to low as follows: MAM, JJA, SON and DJF (Figure 2). The highest FDOs of 18.4% were observed over the TD in MAM and decreased in JJA with the highest values of 17.1%. Distributions of FDOs observed during SON corresponded to those observed in JJA but with a lower value of 14.6% over the TD and higher values observed over the NC. The FDOs observed during DJF were the lowest relative to those of other seasons with the highest value of 11% found over the NC. These results were consistent with those reported previously for East Asia [17]. The spatial pattern of seasonally averaged FDOs was characterized by two high FDO centers over dust source regions (the TD and GD). There existed a dust belt from the dust source regions (the TD and GD) to eastern China, Japan, Korea and beyond, in the latitudinal range of 25◦N to 45◦N. In MAM, FDOs over KJ reached as high as 11 %, which was close to the value for the GD dust source region, indicating efficient long-distance transport in this season. That two distinct FDO patterns exited over the TD and GD implied that the diffusion and long-range transport of dust aerosols from the TD were limited, or to put it differently, most dust aerosols emitted from the TD could not be transported to downwind regions since the TD was mainly surrounded by rugged mountains [3]. Since the whole study domain is under the control of westerlies and obstacles from the Himalayan Mountains, dust aerosols from South Asia and India enter continental China from southwestern China, affecting the NC, SC and KJ in MAM and DJF [32–34]. In the NC, the highest value of FDOs of 16.1% was observed during MAM.

**Figure 2.** The seasonal variations in the frequency of dust occurrence over East Asia based on Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) data.

The vertical distributions of the FDO over the dust source regions (TD and GD) follow a comparable pattern, with significant seasonality (Figure 3). The profile peaked at a height of about 2 km above the surface; beyond this height, the FDO decreased with increasing height. In other words, large volumes of dust aerosols were uplifted into the lower troposphere from the dust source regions, where they could undergo mixing and long-range transport (to eastern China and even to the western Pacific Ocean). The peak FDO values were 52.0% and 40.0% in MAM and 43.5% and 28.5% in JJA over the TD and GD, respectively.

**Figure 3.** Vertical distribution of the frequency of dust occurrence over: (1) the Taklimakan Desert (red line), (2) the Gobi Desert (magenta line), (3) Northern China (green line), (4) Southern China (light blue line) and (5) Korea-Japan (dark blue) based on CALIPSO data.

Over the regions downwind of these sources (the NC, SC and KJ), FDO decreased with increasing height from the ground to the upper troposphere. The closer to the dust source regions, the higher the FDO: of these three regions, the overall FDO profile showed the highest values over the NC and lowest values over KJ, with the profile over SC generally falling in between. At the NC, the near surface FDO was greater than 60% in MAM, SON and DJF with the highest FDO of 67.1% recorded during DJF. This high value could be attributed to the proximity to the dust regions and to the dominance of aerosols generated by use of fossil fuels [35]. FDOs over the SC ranged from 27.0% to 29.9%, with peak values of 29.9% observed during MAM. The highest FDO observed over the KJ (40.2%) occurred during MAM, a pattern that resembled that of the SC. These profiles confirmed the presence of a dust belt in the latitude range of 25◦N to 45◦N that transported significant volumes of dust from the dust source regions.

#### *3.2. Seasonal Distribution of Dust Top Height*

The information on the dust TH can help elucidate mechanisms of long-range transport of the dust aerosols. The dust TH, which is defined as the height above surface elevation (a.s.e.), shows significant seasonal variations over East Asia (Figure 4). Over the TD and GD, the dust TH was the largest in MAM with an average value of 3.1 km (a.s.e.) and the lowest in DJF with an average value of 1.9 km (a.s.e.). Note that the variations here resembled those of the boundary layer height [36]. In comparison, over the downwind regions (the NC, SC and KJ), the dust TH values were still the highest in MAM with a range of 3.9~5.0 km, and the lowest in SON with a range of 2.7~3.2 km. The dust belts were evident in all four seasons. During MAM, the dust TH increased along the dust belt from the dust source regions to the downwind regions, implying progressive vertical expansion of the dust layer as the airmass absorbed dust particles from the ground along its transport pathway. Although not shown in Figure 4, other studies have demonstrated that, aided by the ascending movement of the East Asian and the North American troughs, the vertical expansion of the dust layer can continue all the way to North America [37]. Thus, these patterns of dust storm explain why the highest dust TH values were recorded during MAM in the downwind regions.

**Figure 4.** The seasonal distribution of dust Top Height (km) over East Asia based on CALIPSO data.

#### *3.3. Seasonal Distribution of AOD and the Percentage of D\_AOD to the Total AOD*

To develop a stronger understanding of the contributions of dust aerosols to atmospheric aerosols, it is necessary to further analyze dust aerosols generated from anthropogenic activities by investigating the variations of AOD and the percent contribution of dust aerosols to the total AOD (D\_AOD). The seasonally averaged AOD is generally characterized by two high AOD centers and two low AOD centers over East Asia (Figure 5). Economically developed and industrialized areas of eastern China and areas under the influence of natural dust sources constituted the high AOD centers. In contrast, less developed areas with smaller populations across the Tibetan Plateau and eastern Inner Mongolia constituted the low AOD centers. In addition to natural dust aerosols, aerosols emitted by anthropogenic activities, such as smoke particles from burning of agricultural biomass, sulphate and black and organic carbon aerosols from industrial activities could be responsible for enhanced levels of aerosol loadings during DJF over the SC with a high AOD value of 0.88 [38–40].

The data on D\_AOD reveal once again a distinct dust belt, especially in MAM and DJF. In these two seasons, high percentages (>70%) were found along the belt that extends from the dust source regions to the coastal line of eastern China. The large percentage of greater than 95% observed over the dust source regions indicated that natural dust aerosols served as the most important component of atmospheric aerosols. Meanwhile, obvious seasonal variations were observed. The highest D\_AOD during MAM was closely related to frequent and intensive dust storm events in the dust source regions. The distribution of D\_AOD in SON was similar to that observed in MAM but with lower percentage values [41]. In JJA, the dust belt was much smaller in extent, limited mostly to the dust source regions (the TD and GD) but with abnormally high percentages (>85%). The small spatial extent of the dust belt in JJA can be partially explained by the prevailing climate: JJA is the rainy season in Eastern China, Japan and Korea, and most of the dust aerosols are removed by wet deposition.

**Figure 5.** Seasonal distribution of Aerosol Optical Depths and percentages of D\_AOD to the total AOD over East Asia based on CALIPSO data.

#### *3.4. Extinction Coefficient of Pure and Polluted Dust*

The vertical profiles of extinction coefficients for the three main aerosol types (total aerosol, pure dust and polluted dust) are presented in Figure 6. Here, polluted dust is defined as the mixture of pure dust and smoke particles generated from biomass burning, and total aerosol includes clean marine, pure dust, polluted dust, polluted continental, clean continental, polluted dust and smoke components [42]. Except for a small proportion of polluted dust in DJF, pure dust was the only component of the total aerosols over the dust source regions (TD and GD), showing that pure dust dominated throughout the year. In these two regions, high dust extinction coefficients (TD: 0.17 km<sup>−</sup>1; GD: 0.06 km<sup>−</sup>1) were observed at an altitude of 1~2 km during MAM.

As evidenced in Figure 6s,t, extinction coefficients of polluted dust (<0.04 km−1) were observed over the dust source regions (the TD and GD) in the winter, largely as the result of fossil fuel burning in the winter heating period. Over the NC, the region downwind of the deserts, extinction coefficients of pure dust were much lower than those of polluted dust, indicating the dominate role of anthropogenic activities in this region across all four seasons. Another notable feature about NC is that extinction coefficients of pure dust were higher than those of polluted dust above the height of 1.1 km in MAM (Figure 6c), again confirming that emission and transport of dusts from the deserts were large contributors to air quality problems during this time of the year.

Over the SC and KJ, aerosol extinction profiles did not show considerable seasonal variations, and major aerosol subtypes could be listed in descending order as follows: total aerosols, polluted dust and pure dust. The extinction coefficients of pure dues were less than 0.04 km−<sup>1</sup> year-round while polluted dust accounted for 40~60% of total aerosols, illustrating that local pollution sources dominated over the role of dust long-distance transport.

To further understand the dust transport pattern, we present the vertical distribution of the dust extinction coefficient along the W-E transact at 40◦N (Figure 7). The highest dust extinction coefficients were observed over the TD (80–97.5◦E). If we used a threshold value of 0.1 km<sup>−</sup>1, the dust layer in this longitudinal range extended to a height of 4 km in MAM and JJA, 3 km in SON, and roughly 2.5 km in DJF. Fewer dust storm events occur in JJA than in MAM. However, it appears that such dust devils were able to generate large volumes of dust aerosols through localized disturbances, and because of strong convections, the emitted dust particles were uplifted to up to 4 km beyond the boundary layer height (BLH), where strong wind caused some dust aerosols to be transported to the downwind regions [36,37,43]. On the other hand, due to prevailing topographic conditions, only some particles were successfully transported to the downwind regions on account of topographic obstacles, especially in SON and DJF [3].

**Figure 6.** Profiles of the seasonal average aerosol extinction coefficient (km<sup>−</sup>1) derived from CALIPSO observations over 1. The Taklimakan Desert (**a**,**f**,**k**,**p**), 2. The Gobi Desert (**b**,**g**,**l**,**q**), 3. Northern China (**c**,**h**,**m**,**r**), 4. Southern China (**d**,**i**,**n**,**s**) and 5. Korea-Japan (**e**,**j**,**o**,**t**).

**Figure 7.** Seasonal cross sections of total dust extinction coefficients (km−1) along a latitude of 40◦N and longitudinal range of 73◦E to 105◦ E based on CALIPSO data.

#### *3.5. Horizontal Dust Transport Flux*

Figure 8 shows the horizontal dust transport flux for the five homogenous regions of East Asia. Unsurprisingly, horizontal dust transport flux was the highest over the dust source regions: the TD followed by the GD. Seasonally, the horizontal total dust transport flux over all the five regions was, in descending order: DJF, MAM, SON and JJA. The two dust source regions emitted 95.2 Tg, 35.8 Tg, 46.9 Tg and 82 Tg of dust aerosols into the atmosphere during MAM, JJA, SON and DJF, respectively. The NC and SC sources contributed 35.4%, 15.6%, 31.7% and 54.8% to the total in MAM, JJA, SON and DJF, respectively, and the corresponding flux fractions for the KJ amounted to 9.9%, 6.7%, 6.1% and 12.4%. Notably, since a large percent of the dust over the KJ region originated from the TD, the seasonal variation of the horizontal dust transport flux over the KJ resembled that observed over the TD, with higher levels observed in MAM and DJF and lower levels recorded in JJA and SON.

**Figure 8.** Seasonal distribution of horizontal dust transport flux (Tg) over East Asia (dark blue box: Taklimakan Desert, magenta box: Gobi Desert, green box: Northern China, light blue box: Southern China and yellow box: Korea-Japan).

#### *3.6. Comparison with MERRA-2*

Supplementary Figure S1 presents the MERRA-2 AOD and the percent contribution of pure dust aerosol to AOD for our study region. To allow easy comparison, in Figure S2 we present the same quantities from the CALIPSO retrieval. In Figure 5, the percentages of D\_AOD include both pure dust and polluted dust. In Figure S2, the percentage values are for pure mineral dust originated from the ground, and so a direct comparison can be made of the MERRA-2 D\_AOD percentage values. MERRA-2 can be used as independent evaluation of the CALISPO products because MERRA-2 aerosol variables have been extensively validated against Aerosol Robotic Network (AERONET) measurements [44] and because MERRA-2 does not use CALIPSO to do bias correction. A previous study has shown that the vertical profile of the CALIOP attenuated backscatter coefficient agrees well with that derived by MERRA-2 for Continental USA, South America, Northern and Southern Africa [29]. Here we find that the CALIPSO captured the broad spatial patterns of both AOD and D\_AOD for East Asia as calculated by MERRA-2. In the case of AOD, both show two pollution centers, with a heavily polluted region in eastern China and another less pollution region in western China. In the case of D\_AOD of pure dust, the highest values were located in the dust source regions (Taklimakan Desert and Gobi Desert) according to both data products.

#### **4. Conclusions**

The present study analyzed the seasonal variation of dust aerosols and transport over East Asia using CALIPSO retrievals. To develop deeper insight into the spatiotemporal distribution of dust aerosols, the study domain was divided into five homogenous regions including the Taklimakan Desert (TD), the Gobi Desert (GD), northern China (NC), southern China (SC) and Korea-Japan (KJ) with distinct emissions sources. The frequency of dust occurrence, dust Top Height and extinction coefficients were assessed for each study region.

Our results confirmed that large amounts of dust particles were emitted from the dust regions (TD and GD). The emitted dusts formed a dust belt within a latitudinal range of 25◦N to 45◦N and extended to eastern China, Japan, Korea and the Pacific Ocean. The dust belt was strongest in the spring outside the rainy season. High frequencies of dust occurrence were found over the dust source regions sequenced in descending order as follows: spring, summer, autumn and winter. FDOs showed to isolated centers over the Taklimakan Desert and Gobi Desert, implying limitations of the diffusion and long-range transport of the dust aerosols over the Taklimakan Desert, meaning that most of uplifted dust aerosols were blocked by the mountains around the Tarim basin and eventually settled back to the ground.

The dust top height over East Asia presents significant seasonal variations with values of roughly 3.5~5.1 km (a.s.e.) over the dust source regions and recorded in descending order as follows: summer, spring, autumn and winter. Over the dust belt, the dust Top Height increased from dust source regions to downwind regions; furthermore, frequent dust storms occurring over East Asia during the spring showed why the dust Top Height was the highest in the spring across the downwind regions.

The AOD was generally characterized by two high AOD centers and two low AOD centers over East Asia. The dust belt was the shortest and had the lowest AOD in the summer with higher dust to AOD percentages (>85%) observed only over the Taklimakan and Gobi Deserts, showing that dust aerosols can hardly travel long distances in the rainy season due to wet deposition.

Pure dust was predominant across all four seasons. Over northern China and in areas direct downwind of the dust source regions, polluted dust generated from anthropogenic activities also played an important role. The dust extinction coefficients over the Taklimakan Desert were sequenced in descending order as follows: spring, winter, summer and autumn. The west-to-east transect of the extinction coefficient further confirmed that only some of the emitted dust particles successfully reached the downwind regions due to topographical obstacles and this blocking effect was especially strong in the autumn and winter.

The main advantages of the CALIPSO system lie in its ability to obtain vertical profiles of aerosols, which are crucial for understanding the distribution and transport of dust. However, the satellite completes a cycle over 16 days, producing data gaps that may miss detection of some dust events. A preliminary comparison with MERRA-2 data product, which is continuous in time, showed that when averaged over multiple years, the horizontal distributions of AOD and D\_AOD from the CALIPSO captured the broad patterns calculated by MERRA-2. Future work should use MODIS, MERRA-2 and CALIPSO data to build a 3D structure of the dust trajectory over East Asia to gain STRONGER insight into long-range transport of East Asian dust to other regions of the globe.

**Supplementary Materials:** The following are available online at http://www.mdpi.com/2072-4292/11/6/701/s1, Figure S1: Seasonal distribution of Aerosol Optical Depth (left panels) and percentage contributions of pure dust to the total AOD (right panels) over East Asia based on MERRA-2 data, Figure S2: Seasonal distribution of Aerosol Optical Depth (left panels) and percentage contributions of pure dust to the total AOD (right panels) over East Asia based on CALIPSO data.

**Author Contributions:** D.L., T.Z., and S.C. designed the study. D.L. carried out THE research and wrote the manuscript. T.Z., R.B., and S.C. contributed to the preparation of the manuscript through review, editing, and comments. Z.L. wrote the results section. Y.W. and Y.Z. helped polish the manuscript. And all authors were involved in modifying the paper, and the literature review.

**Funding:** This research was supported by the Foundation for National Natural Science Foundation of China (No. 41830965; No. 91837103).

**Acknowledgments:** We thank TC Chakraborty, Yale University, for providing the MERRA-2 data graphs. **Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
