High-spectral-resolution lidar (HSRL) is a powerful tool for atmospheric aerosol remote sensing. A ground-based high-spectral-resolution lidar (HSRL), operated at 532 nm wavelength, has been developed at Zhejiang University (ZJU) for aerosols and clouds studies. This lidar provides vertical profiles of aerosol scattering ratio together with lidar ratio and particle depolarization ratio at 532 nm. Determination of overlap function is a key step in the calibration of a high-spectral-resolution lidar (HSRL) and important guarantee of data retrieval, an iterative-based general determination (IGD) method for overlap function in HSRL is proposed. The standard method to retrieve the extinction coefficient from HSRL signals depends heavily on the signal-to-noise ratio (SNR). An iterative image reconstruction (IIR) method is proposed for the retrieval of the aerosol extinction coefficient based on HSRL data under low SNR condition. With the optical properties, a state-of-the-art method for feature detection and classification is proposed to automatically identify the features attributed to dust/polluted dust, urban/smoke, maritime aerosols, as well as ice and liquid water cloud during day and night.
Aerosol particles play an important role in Earth’s climate system by absorbing and scattering solar radiation and by modifying cloud properties . Currently, the aerosol impact is the largest source of uncertainty in estimation of the radiation flux forcing . For reducing these uncertainties, the continuous and global monitoring of atmospheric aerosol is needed . Lidar has been proved to be a promising tool for vertical profiling the cloud and aerosol properties . Numerous types of lidars nowadays are used for aerosol observations, including NASA-CNES Cloud-Aerosol Lidar and Infrared Pathfinder Satelite Observation (CALIPSO) , the airborne NASA Langley high-spectral-resolution lidars (HSRL) instrument , the ground-based Raman lidars (RL) from European Aerosol Research Lidar Network (EARLINET) , the Micro-Pulsed Lidar (MPL) network , the Asian Dust and aerosol lidar observation Network (AD-Net) , etc. Elastic backscatter lidars are widely used for aerosol studies, however for the evaluation of the aerosol extinction a priori assumption of the ratio of aerosol extinction to backscatter (i.e., the lidar ratio) is needed . Meanwhile, HSRLs and RLs can independently retrieve aerosol backscattering and extinction coefficients , providing the basis for aerosol classification.
However, lidar measurements cannot be taken from the ground level due to the incomplete geometrical overlap of the receiver field-of-view (FOV) and the laser beam. This overlap is usually described by the Overlap Function, varying from 0 to 1 . It may not be a problem for spaceborne lidar or middle and upper atmospheric research but is indeed a problem for ground-based lidars that study the dynamics of the planetary boundary layer (PBL) or air polluted. To reduce the impact of the overlap function in ground-based lidars, improvements in optical design are also advanced . Moreover, in order to achieve a higher detection range and better signal-to-noise ratio, receivers with larger aperture have been widely employed, leading to a higher incomplete overlap range. Thus, accurate measurement of overlap function is of great importance for better lidar detection of aerosols.
An iterative-based general determination (IGD) method for retrieving the overlap function in HSRL is proposed to refer to the method used in Raman lidar . In the iterative process, the signals of the combined and molecular channels are employed to produce backscatter coefficients with the HSRL method  and Fernald algorithm . A correction to the overlap function, whose initial value can be 1, is then produced by the difference between the two retrieved backscatter coefficients. The new overlap function is then used to modify the old and produce a new signal of the combined channel. New backscatter coefficients will be retrieved again, and the steps followed continue. The iterative process stops until the backscatter coefficients got by the Fernald algorithm is close enough to that from the HSRL method. No crucial requirement of atmospheric condition nor additional auxiliary optical system is needed in this method. Moreover, compared with Raman lidar, the frequency shift of HSRL backscatter signal is negligible, so the application is simpler. Validation experiment is performed with the HSRL system at Zhejiang University (ZJU HSRL). The obtained overlap function meets well with that from the dual-FOV method, which shows that the IGD method can be used for monitoring of HSRL overlap function with real-time lidar signals.
The intensity of lidar signals attenuate inversely as the square of the distance. The signal-to-noise (SNR) of the lidar signals is affected by the solar background light and the noise of the acquisition circuit . The SNR decreases dramatically with the increase of the detection distance, and to retrieve the aerosol extinction coefficient accurately, the standard method in HSRL requires sufficient SNR . The molecular channel in HSRL, which eliminates Mie scattering by the spectral discrimination filter, has a lower signal intensity. In addition, the standard method based on the numerical derivative makes the retrieval process extremely sensitive to the SNR, which greatly affects the accuracy of the aerosol extinction coefficient in HSRL .
A novel method based on an iterative image reconstruction (IIR) technique is proposed to retrieve the aerosol extinction coefficient in HSRL with improving accuracy under low SNR conditions. For a photomultiplier tube (PMT) operating in the analog detection regime, a previous study had proved that the noise of the lidar signals is the compounded Poisson noise . As the atmospheric conditions are relatively stable in a short period, lidar signals reflecting similar atmospheric conditions are related. Taking advantage of these properties of lidar signals, the IIR method reconstructs the “image” of the “feature” of the lidar signals on the basis of reasonable estimates of the lidar ratio by minimizing the difference between the reconstructed signals with the observations under an appropriate assumption of signal noises. The “feature” is identified with a scattering ratio profile greater than one standard deviation from the expected clear-sky . Using an iterative regularization method to carry out the best optimization estimates of the lidar ratio, the IIR method transforms the ill-posed problem of the numerical derivative into a well-posed problem of the image reconstruction. The IIR method has the potential to retrieve aerosol optical properties reasonably and automatically as well as significantly improve the detection capability of the HSRL under the condition of low SNR.
The applications of lidar systems are still limited by the lack of capabilities to classify types of aerosols and clouds automatically. Firstly, considering that measurements from lidars contain a large amount of invalid data, the raw data should be analyzed to distinguish valid data associated with aerosols and clouds (i.e., feature) efficiently. In other words, the analysis must locate all layer boundaries, which is usually referred to as the feature detection. Afterwards, the automatic classification of features should be performed, also called feature classification.
Observations over eastern China from a ground-based 532 nm HSRL system will be presented by Zhejiang University (ZJU) lidar group for the very first time. Hereafter, we will focus on the feature detection and classification based on measurements of the ZJU HSRL, as a detailed reference of the data processing for other HSRLs at single wavelength. Generally, the principle of feature detection from previous works is to discriminate layers from the pre-estimation of aerosol free height range based on the change in the raw lidar signals. Unlike previous methods adopted by elastic lidars, the direct measurements of the ratio of total backscatter to molecular backscatter (i.e., scattering ratio) from HSRLs are exploited to determine vertical extent of feature in this work. It only relies on the fact that the gain ratios of detectors are stable enough to retrieve backscatter accurately, which have been proved in the practical HSRL system . As for feature classification, this study demonstrates that different aerosol types can be distinguished using the measurements combined with lidar ratio and aerosol depolarization ratio analogous to the method proposed by Groβ et al. , modified to fit ZJU HSRL’s needs and capabilities. Our results are validated with the sun photometer and supported by comparable results of typical aerosols from established literatures. In addition, our measurements were performed over eastern China, which also extend the database of aerosols especially with respect to urban and marine and dust near the Yangtze River Delta region of China.
The basic idea of the IGD method is to calculate the overlap function from the difference of backscatter coefficients from HSRL and those from the Fernald algorithm. The IGD method, as an iterative algorithm, is summarized in Figure 1a. This iterative process can be considered as a reverse optimization procedure that can be used in many aspects . According to Figure 1a, the IGD method takes combined signals and molecular signals as inputs and calculates the particulate backscatter coefficients based on the HSRL retrieval and Fernald algorithm, respectively. The particulate backscatter coefficients based on the HSRL retrieval are not affected by the overlap function; thus, it can be regarded as the exact values. While the particulate backscatter coefficients based on the Fernald algorithm are affected by the overlap function. The particulate backscatter coefficients' difference between these two methods can be converted into the correction of the overlap function, thereby correcting combined signals. Through iterative optimizations, the particulate backscatter coefficients calculated by the Fernald algorithm approach the exact values. Therefore, the signals in the last iteration can be considered as the combined signals corrected the overlap function.
Figure 1. Flow diagram for the overlap function retrieval in the IGD method, (a) the general flow of the iterative process; (b) the flow chart to calculate the correction of the overlap function based on the difference of the backscatter coefficients.
In ground-based HSRL systems, Eloranta et al. proposed a dual-FOV lidar to measure the overlap function . In the dual-FOV method, the difference between the two overlap functions corresponding to two receiving channels is used to determine the overlap function of the wide-FOV receiving channel. To verify the IGD method, a wide-FOV receiving channel is introduced into the HSRL system to obtain the overlap function. The accuracy of the retrieval overlap function from the IGD method is verified by comparing the results of two different methods.
In general, the lidar signal detected by the wide-FOV channel can be expressed as
In Eq. (1), is the energy received by the wide-FOV receiving channel and is the overlap function corresponding to the wide-FOV receiving channel. Through special design, such as increasing the FOV of the wide-FOV channel or reducing the distance between the emission axis and the receiving axis, the transition zone of the narrow-FOV channel already falls into the overlapping zone of the wide-FOV channel. Therefore, , the overlap function can be derived from the signals of the dual-FOV receiving channels,
Where and represent the energy received by the narrow-FOV and wide-FOV receiving channel, respectively.
Experimental demonstration is performed by our dual-FOV ZJU HSRL system  locates at Yuquan Campus of Zhejiang University (30.26 N, 120.12 E，Hangzhou, China) on October 28th, 2019, as an example to verify the IGD method. Here, one profile is taken as an example to illustrate the overlap function retrieval processes by two methods. Figure 2a corresponds to the dual-FOV method. Signals collected by the auxiliary receiving system () and the main receiving system () are used to obtain the overlap function of the narrow-FOV channel. The overlap function and backscatter coefficients profiles obtained from the IGD method are shown in Figure 2(b). According to the observations based on the dual-FOV HSRL system in Hangzhou, China, the lidar ratio is assumed to be 40 sr in the IGD method. With the increasing of iteration steps, the backscatter coefficient profile generated in each iteration step converges to the results obtained from the HSRL system.
Figure 2. (a) The retrieval processes of the dual-FOV method. The blue dash line is the signal measured by the narrow-FOV channel, and the light blue dash line is the signal measured by the wide-FOV channel. The solid black line is the retrieved overlap function. (b) The retrieval processes of the iterative method. The red dash line is the backscatter coefficient of the HSRL. The grey dash lines indicate some of the backscatter coefficient profiles during the iteration. The light blue dashed line and the green dashed line indicate the backscatter coefficients generated during the 1st and 20th iterations, respectively. The solid blue line represents the backscatter coefficient profile generated from the 318th iteration, which meets well with the HSRL inversion results. The solid black line is the retrieved overlap function.
The 100 profiles (about 3.3 hours) obtained from the experiment are used to calculate the error of the overlap function obtained by the IGD method and the dual-FOV method, as shown in Figure 3 and Table 1. The iterative retrievals show good consistency with the results of the dual-FOV method. The black line in Figure 3 is the mean overlap function profile obtained by the iterative method, and the three standard deviations are labeled out with error band. The solid red line is the mean retrieval result of the dual-FOV method. As shown in Figure 3, the retrieval results show good accordance with each other, and the mean MAD value is 4.56%. In addition, the mean SD of the overlap function calculated by the iterative method is 0.0177, close to the MC simulation value of 0.0091, which shows that the method and system both have strong stability in practical use. The dual-FOV method has no extra assumed atmospheric parameters in most stable atmospheric conditions. Therefore it can be considered as an accurate real-time overlap function. By comparing real-time overlap functions obtained from two methods, the effects of changes in the optical axis and atmospheric environment can be excluded. However, due to the deviation of the signal-to-noise ratio and the lidar ratio assumption in the experiment, the accuracy of experimental results might be slightly worse than the simulation results.
Figure 3. The retrieval results of the overlap function with two different methods. The red line represents the result of the dual-FOV method. The black line represents the result of the iterative method, and the error band denotes three standard deviations.
Table 1. RMSE, MAD, and SD results of the dual-FOV method and the IGD method.
It should be noticed that the aerosol backscatter coefficient could be used to constrain the retrieval of the aerosol extinction coefficient, since there is a linear relation (i.e., lidar ratio) between the aerosol backscatter coefficient and extinction coefficient. The range of the lidar ratio is about 0-100 typically. Once the correct estimate of the lidar ratio is obtained, the aerosol extinction coefficient could be algebraically computed by the product of the aerosol backscatter coefficient and the lidar ratio. In this work, the IIR method is proposed to reconstruct the “image” of the “feature” of the lidar signals on the basis of a reasonable estimates of the lidar ratio under appropriate assumption of signals noise. In this study, the Gaussian PMF is utilized when the observations are from the analog mode lidar measurements. The characteristics of the proposed method can be summarized into two points: 1) considering the relation of the lidar signals and utilizing the feature detection method to obtain the “feature” of the “image”, 2) transforming the ill-posed problem of the numerical derivative into a well-posed one of the estimates of the lidar ratio. This “feature” is identified with the Feature Detection Algorithm from a previous work to represent the data associated with aerosols and clouds . The region where the scattering ratio profile is greater than one standard deviation from the expected clear-sky is determined as the “feature.” The reconstruction of the “image” of the “feature” must fit with the observations in HSRL to ensure the accuracy of the optical property retrieval. This process would be achieved by using an iterative optimization method on the basis of an adaption of the sparse Poisson image reconstruction algorithm (SPIRAL) .
Figure 4 shows the flow diagram for the retrieval of the optical properties in the IIR method, where is the retrieval window of profile and is the half width of the window. The different region of clear-sky area we utilized may change the results of the “feature” and the retrieval based on the feature detection is easy to reproduce once the “feature” is determined. Therefore, combining the feature detection algorithm with the SPIRAL of Gaussian model, the IIR method has the potential to obtain the optical properties of the “feature” reasonably and automatically. The lidar ratio of each profile could be retrieved through the corresponding retrieval window. Hereafter, the final result is obtained by averaging the lidar ratio in every . More details about SPIRAL can be seen in .
Figure 4.Flow diagram for the aerosol extinction coefficient and the lidar ratio retrieval in the iterative image reconstruction (IIR) method. The details of the adaption of sparse Poisson image reconstruction algorithm (SPIRAL) can be seen in .
In order to verify the feasibility of the IIR method, MC simulations are employed, which allow us to compare the results obtained by different methods with the true values input into the simulations. For this section, the synthetic dataset was created on the basis of real observations of cirrus clouds between 8.5–10 km with 109 profiles and used to estimate the performance of different methods. The measured backscatter coefficient of the cirrus clouds is employed and the extinction coefficient of the cirrus clouds is created with the simulated lidar ratio. The lidar ratio of the cirrus clouds is modeled in Gaussian with a value of 26.28 ± 3.63 sr, while it is set to sr for the clear-sky region . The modeled inputs are presented in Figure 5.
Figure 5. (a) The backscatter coefficient, (b)the extinction coefficient, and (c) the lidar ratio used in the simulation experiment. The lidar ratio inside the clouds is modeled in Gaussian with a value of 26.28 ± 3.63 sr, while it is for the background region.
The lidar signals are generated from forward MC simulations based on the HSRL model mentioned in Section 3.1. Then, a Gaussian distributed random noise with standard deviation equal to the square root of the signals is incorporated into the lidar signals. For the simulation, the scene represents a real dataset with low SNR where the resolution is 7.5 m by 10 s to compare the performance of two methods. The system specifications employed by MC simulations are listed in Table 2, the is a constant for simplicity, and the is obtained by the measured length and temperature of the iodine molecular absorption filter.
Table 2. System specifications employed by Monte-Carlo (MC) simulations.
In this subsection, the performance of the standard and the IIR method is compared against each other. The root mean square error (RMSE) is used to indicate the performance of different methods, along with the relative bias. The RMSE can be defined as:
where ; the superscript “true” and “retrieval” represent the true value and the retrieval value, respectively.
Hereafter, we only focus on the extinction coefficient retrieval of aerosols and clouds (i.e., features) identified by a feature detection method since the accurate measurements of the clear-sky region is beyond the capabilities of the HSRL. The SG filter is used as the low-pass filter in the standard method for both temporal and spatial axes. The SG filters of both temporal axis and spatial axis are first-order polynomials, whose window sizes are nine profiles and nine altitude bins, respectively. The moving SG window for optical thickness is 71 altitude bins. These SG filter parameters are selected based on the minimization of the extinction coefficient RMSE. The window of the IIR method is nine profiles as a balance of efficiency and accuracy.
Figure 6 shows the results obtained by the standard method (left panel) and the IIR method (right panel). The true values are shown in Figure 5. Figure 6a, c, and e are the results obtained by the standard method, which represent the extinction coefficient, the absolute errors of the extinction coefficient of the true values and the retrieval values, and the lidar ratio, respectively. Figure 6b, d, and f are obtained by the IIR method. The dashed lines in f are selected for further analysis, whose profile numbers are 19, 48, 78, and 99, respectively. The extinction coefficient obtained by the standard method is smoother than that obtained by the IIR method. However, the absolute errors of the extinction coefficient from the IIR method (Figure 6d) are smaller than those obtained by the standard method (Figure 6c). Comparing to the input values in Figure 5b, it could be noticed that the structure of the extinction coefficient from the IIR method is closer to the input values than that from the standard method. Furthermore, the comparisons of the true values in Figure 5c and the results obtained by two methods indicate that the IIR method is able to achieve better results of the lidar ratio by the image reconstruction technique. The standard method would require more accumulated profiles and altitude bins (i.e., a higher SNR but a lower spatial and temporal resolution) to reduce the variance of the noise to a sufficient level until the slope of the optical thickness are non-negative.
Figure 6.Comparisons of the results obtained by the standard method (left panel) and the IIR method (right panel) based on the simulation experiments. (a), (c), and (e) are obtained by the standard method, which represent the extinction coefficient, the absolute error of the extinction coefficient of the true values and the retrieval values, and the lidar ratio, respectively. (b), (d), and (f) are obtained by the IIR method. The dashed lines in (f) are selected for further analysis, whose profile numbers are 19, 48, 78, and 99, respectively.
The reason why the IIR method does not suffer from the same problem (negative slope of the optical depth) as the standard method is because it does not try to retrieve the aerosol extinction coefficient algebraically. The IIR method transforms the ill-posed problem of the numerical derivative into a well-posed problem of the image reconstruction. The lidar ratio is always limited to the bounds of 0–100, as it should be physically. Figure 7a and c show four profiles of the extinction coefficient and the lidar ratio obtained by different methods, respectively. The locations of four profiles are shown in Figure 6d. The lines in black are the true values, while the blue and red lines represent the results of the “feature” obtained by the standard and the IIR method, respectively. What we concerned most is the region containing a “feature,” i.e., dense aerosols and clouds. Therefore, a feature detection method is utilized to exclude the clear-sky region . It is obvious that the structures of the extinction coefficient and the lidar ratio obtained by the IIR method are closer to the true values than those obtained by the standard method. The relative biases of the retrieval of the extinction coefficient and the lidar ratio are presented in Figure 7b and d, respectively. It could be found that the extinction coefficient and the lidar ratio obtained by the IIR method have smaller biases than those from the standard method. The relative bias of the extinction obtained by the IIR method is mostly within 20%, while that obtained by the standard method is much bigger. From the Table 3, the mean relative bias of the extinction coefficient obtained by the IIR method is 17.0%, while it is up to 28.5% with the standard method. The maximum relative bias of the lidar ratio obtained by the standard method is up to 100%. It indicates that the standard method is not suitable for the further aerosol analysis under the condition of low SNR, which is sensitive to the extinction coefficient errors, such as the feature classification, the retrieval of microphysical properties, etc. . The relative biases of profile 48 (the second profile) from the IIR method at the start of “feature” may be caused by the system constant and the initial value of the lidar ratio. The mean relative bias of the lidar ratio obtained by the IIR method is 8.5%, while it is up to 22.6% with the standard method. The green lines in Figure 7 represent the results obtained by the IIR method based on the Poisson model. The comparisons between the red lines and the green lines indicate that the results of the Gaussian model are similar to those of the Poisson model and the relative bias of the Gaussian model is smaller than that of the Poisson model.
Figure 7.Comparisons of the extinction coefficient and the lidar ratio obtained by different methods. (a) The profiles of the extinction coefficient, whose profile numbers are 19, 48, 78, and 99, respectively. The lines in black are the true values, while the blue, red, and green lines represent the results obtained by the standard method, the IIR method, and the IIR method based on the Poisson model, respectively. The area within the dotted pink line is the valid area of clouds. (b) The relative biases of the extinction coefficient. The blue, red, and green lines represent the relative biases obtained by the standard method, the IIR method, and the IIR method based on the Poisson model, respectively. The black dotted line represents a relative bias of 10 percent, while the yellow dotted line represents a relative bias of 20 percent. The area within the dotted pink line is the valid area of the clouds. (c) and (d) represent the profiles of lidar ratio and the relative bias of the lidar ratio, respectively.
Table 3. Error analysis of optical properties. RMSE：root mean square error.
Here we follow the CALIPSO terminology and refer to aerosol and cloud layer as “feature” . Generally, a feature in the lidar profile is characterized by enhanced atmospheric backscattering and by boundaries in the form of a more or less rapid change. Therefore, the method of feature detection could be normally based on threshold or derivative method. Benefited from the direct measurements of scattering ratio from HSRL, the expected feature is quite straightforward to identify without pre-estimation of aerosol free height range. After retrieval of optical properties, we expect a scattering ratio R of 1 where is aerosol free height range.
However, as the scattering ratios are derived from the noisy lidar data, the errors which came from insufficient SNR can easily lead to wrong estimation and fluctuation of scattering ratio. Considering the effect of lidar signals’ SNR, the threshold to determine whether features or aerosol free height range should be dynamically varied with SNR to prevent the misidentification of aerosol plume. Based on the theory of general error propagation , it is mainly accomplished by the evaluation of uncertainty of scattering ratio in our method.
A threshold profile is then defined as
where is the expected uncertainty of aerosol backscatter.
That is, the scattering ratio profiles greater than one standard deviation from the expected aerosol free height range are considered to potentially contain a feature :
An example of the threshold profile determined by lidar signals is given in Figure 8. The mask of feature that is denoted by black indicates the clear-sky region obviously where the scattering ratio is smaller than 1. However, the mask of feature that has been determined to be false detection is denoted by pink, which means that the scattering ratio of the bin is greater than 1 but lower than the threshold profile. And true features are denoted by blue. It can be noticed that the only two features are identified correctly on the basis of our method.
Figure 8. An example of feature detection from the scattering ratio and threshold profile. The theoretical threshold profile (dashed yellow) and the practical threshold profile (dashed red) and the scattering ratio profile (solid blue) are given. The masks of the features are noted by the color in the bars on the right-hand side. Features that are denoted by black indicates the aerosol free height range where the scattering ratio is smaller than 1. Features that have been determined to be false detections are denoted by pink. And true features are denoted by blue.
After pre-processing and retrieval of optical properties and feature detection which have been mentioned in previous sections, the next step in lidar data evaluation with respect to optical data is the identification of aerosol or cloud layers. Cloud-aerosol discrimination (CAD) depends on the information available from the lidar measurements. Only backscatter information at a single wavelength (532 nm) is considered here. Therefore, the different signature of clouds and aerosol must be used. Clouds typically show stronger backscattering, lower lidar ratio, higher attenuation, steeper boundaries than aerosols. Generally, our method relies on the fact that backscattering from clouds is much larger than aerosol backscattering. The classification between clouds and aerosol here is based on a set of rules and empirical thresholds determined by both the expected scattering properties . With one exception, ice clouds could have small scattering ratio and high depolarization ratio which are similar to dust and polluted dust aerosol. But the lidar ratio of ice cloud is smaller than dust and polluted dust significantly, and it can be separated from dust and polluted dust based on this fact. To illustrate the choices for the empirical thresholds, the corresponding distributions of the retrieved feature properties themselves are given in Figure 9a from the measurements of ZJU HSRL. The classification starts by making an initial partitioning of clouds and aerosol using the particulate backscatter and depolarization ratio threshold depicted in Figure 9b.
Figure 9. (a) The histograms of depolarization ratios and scattering ratio from selective ZJU HSRL measurements. (b) The scheme of CAD by aerosol depolarization ratio and scattering ratio. Each of data located in the green region or the gray region is classified as aerosol or cloud, respectively. The red region is considered as the mixture of aerosols and ice clouds where the lidar ratios of aerosols (a threshold of 40 sr) is significantly higher than the ice clouds.
After classifying Level 1 data as cloud or aerosol, the final goal is to identify the type of the Level 1 data further. Cloudy data are further classified into liquid clouds or ice clouds using the depolarization ratio threshold in Figure 10a. All clouds that occur above a temperature of 0 ℃ are considered to be liquid or below a temperature of -40 ℃ are considered to be ice oppositely. When the temperature is below 0 ℃ and above -40 ℃, both liquid and ice could exist. In this situation, depolarization is especially useful as only non-spherical particles like ice induce a strong depolarization . Generally, higher aerosol depolarization ratio usually indicates ice cloud. The particle depolarization ratio of water cloud is lower than 0.05 typically based on our measurements.
Figure 10. Flow diagrams for the feature classification of (a) cloud and (b) aerosol, respectively.
For the classification of aerosols, the categories of our method are similar to the previous works from CALIPSO  and Burton et al . The types of aerosols could be divided into three main types in this study followed previous works: dust/polluted dust, marine/polluted marine, and urban/smoke. On the basis of polarization sensitive HSRL measurements, two intensive parameters, that is lidar ratio Sa and aerosol depolarization ratio δa, are derived. It has been generally noticed that lidar ratio and aerosol depolarization ratio only vary with aerosol type and do not depend on concentration. Past analyses of these quantities reveal that characteristic value can be attributed to the different aerosol types [39. A summary of δa and Sa found from lidar measurements at 532 nm for different aerosol types has been listed in Table 4.
Table 4. Summary of key finding concerning the lidar ratios and aerosol depolarization ratios at 532 nm found in established literatures and our finding.
It is obvious that δa and Sa show the large differences for the different aerosol types. The spread ranges from low δa values for marine aerosols of 2-10% to high δa values of 15-35% for dust and polluted dust. Dust and polluted dust are clearly distinguishable by their higher aerosol depolarization ratios. Other three aerosol types show quite similar range of δa values, but differ in other optical parameters. While Sa for marine aerosols is low with a mean value of ～20 sr, urban aerosols show a significantly higher mean value of ～50 sr. The similar results of lidar ratio were found between smoke and urban aerosols. Measurements of aerosol depolarization and lidar ratio at 532 nm are not sufficient for the clear separation between urban and smoke. Additional information such as the spectral depolarization ratio may be used to separate both types of aerosols , however, such information is beyond the current measurement capabilities of our instruments. In this work, we treat the two types of aerosol (i.e., urban and smoke) as one combined type for convenience.
Based on the discussion above, to classify aerosols via the characteristics of different aerosols, Figure 10b demonstrates the classification scheme of aerosols. Firstly, dust and polluted dust could be separated from other aerosols using threshold segmentation of 0.15 aerosol depolarization ratio. Also, it could be noticed in Table 4 that marine and polluted marine show significant lower lidar ratios than other types of aerosols. The rest of aerosols, the lidar ratios of which are lower than 30, should be classified as marine or polluted marine. The remaining aerosols are classified as the urban/smoke type in the end.
The raw results of feature classification could be finalized by applying a filter. For each bin of the data, the classification of aerosols or clouds is replaced by the most common type surrounded the bin. This filtering is done both to reduce noise in the classification and to correctly identify the edges of ice clouds, which are commonly initially misidentified as dust and polluted dust with low backscatter and high aerosol depolarization ratio .
To validate the efficiency of the method we proposed in practical lidar signals, Figure 11 gives the inputs into and the example of feature classification on 6 May 2019 at Zhoushan, Zhejiang, China. All feature categories identified by our method are presented in this example. As can be seen from the observation in Figure 11d, it is presented by a three-layer structure roughly: clouds, dust/polluted dust, and marine/polluted marine. The scattering ratio of the two layers at the bottom are relatively low which indicate that these layers are aerosol layers that shown in Figure 11a. The top layer, which has a much higher scattering ratio and drastic change on a short time scale, should be classified as clouds obviously. First, we take a look at the clouds. The supercooled liquid clouds are observed at the altitude of ~5 km. It should be noticed that the depolarization ratio below the clouds is significantly higher, pointing to the presence of ice particles (Figure 11b). Such coexistence of supercooled liquid and ice layer was reported by RL before . Cirrus clouds appear momently at 17:10 (UTC+8) at a height of 9 km, soon after the mixed-phase clouds were observed again.
Figure 11. Feature classification for 6 May 2019 in Zhoushan. Shown are the inputs used by the feature detection and classification scheme: (a) scattering ratio, (b) aerosol depolarization ratio, (c) lidar ratio. (d) The results of feature classification with each bin classified as either water cloud (green), ice cloud (red), marine/polluted marine (purple), or dust/polluted dust (yellow).
Next, an elevated aerosol layer appears in the 3000 m - 4000 m range and stays well separated from the bottom aerosol layer. The depolarization ratio of the dust/polluted dust aerosols observed in this case are ~20%, while the depolarization ratios of other aerosols at 532 nm normally are below 15%. Hence depolarization measurements provide a convenient way to separate the dust and polluted dust aerosol layers from other aerosols. The bottom marine/polluted marine aerosol layer is often observed under the boundary layer in Zhoushan. This is mainly due to the fact that the observation site is very close to the coast and there is little intensive anthropogenic activity. The lidar ratios at 532 nm for the marine/polluted marine aerosols are about 15-30 sr and are intermediate between the pure marine and other aerosols that shown in Figure 11c, which is consistent with past observations of polluted marine .
Validating all of our retrievals of optical properties is difficult since no instrument at our site could make comparable measurements (such as RL). We do have compared our measurements of backscatter and depolarization ratio with MPL (provided by Hangzhou Meteorological Bureau) and CALIPSO, and a good agreement is found. But it is still not convincing since the retrievals based on measurements from Mie lidar rely on the assumption of lidar ratio. However, for aerosol optical depth, the Hangzhou Meteorological Bureau is equipped with Cimel CE-318 sun photometer allowing the comparison of daytime aerosol optical depth with HSRL shown in Figure 12. The coincident sun photometer optical depth is calculated at 532 nm by applying the Ångström exponent to the 550 nm channel optical depth. The AOD of the clouds in the HSRL measurements is excluded by the CAD (see in the section 4.1.2). The results are then compared with the level 1.5 cloud-screened sun photometer data. The data of optical depth provided by Hangzhou Meteorological Bureau are available from 18 Dec 2019 through 08 Feb 2020 since the sun photometer was sent to Beijing for recalibration for a period of time.
Figure 12. Comparison of AOD derived from ZJU HSRL and the Cimel CE-318 sun photometer. Summary statistics are given in the top-left corner: then sample data points (N), the slope of the linear regression line, the bias error, and the correlation coefficient (R2).
The comparison in Figure 12 depends on several portions of our method, including the correct optical properties retrieval, feature detection and classifying the feature as aerosol or cloud. In this comparison, a good agreement is found with a bias error about -5.9%. We consider the bias error presented in Figure 12 to be reasonably good considering that the sun photometer was measuring optical depth along a different path through the atmosphere and the lack of direct extinction measurement from HSRL at lower heights (<200 m).
The investigation of optical properties of aerosols is a major topic in atmospheric research. Therefore, a number of established publications also provide the opportunity to validate whether our measurements are reasonable. Figure 13 shows the results of feature classification from measurements of ZJU HSRL displayed as a two-dimensional graph of lidar ratio versus aerosol depolarization ratio. Those points, which are mainly from selective measurements in coastal city Zhoushan classified into marine and polluted marine, are colored blue. The dust/polluted dust and urban/smoke aerosols are colored yellow and purple, respectively. Figure 13 also shows the aerosol intensive properties at 532 nm from some other measurements from established literatures . There are general qualitative agreements showing, for example, that dust/polluted dust aerosols have higher depolarization ratio and lidar ratio, while marine/polluted marine aerosols typically have lower lidar ratio. Figure 13 also clearly shows that there can be considerable spread in these observations for the same type of aerosols observed in different locations. Detailed statistical results including our findings and established literatures of different aerosols have been listed in Table 4.
Figure 13. The results of the classification from the measurements of ZJU HSRL. Also indicated in this figure are the aerosol types identified by Burton et al. (2012), Groβ et al. (2011 and 2013), Ansmann et al. (2003), Xie et al. (2008), and Stephanie et al. (2018).