Articles | Volume 13, issue 12
Atmos. Meas. Tech., 13, 6675–6689, 2020
Atmos. Meas. Tech., 13, 6675–6689, 2020

Research article 09 Dec 2020

Research article | 09 Dec 2020

A novel lidar gradient cluster analysis method of nocturnal boundary layer detection during air pollution episodes

A novel lidar gradient cluster analysis method of nocturnal boundary layer detection during air pollution episodes
Yinchao Zhang, Su Chen, Siying Chen, He Chen, and Pan Guo Yinchao Zhang et al.
  • School of Optics and Photonics, Beijing Institute of Technology, Beijing 100081, China

Correspondence: Siying Chen ( and Su Chen (


The observation of the nocturnal boundary layer height (NBLH) plays an important role in air pollution and monitoring. Through 39 d of heavy pollution observation experiments in Beijing (China), as well as an exhaustive evaluation of the gradient, wavelet covariance transform, and cubic root gradient methods, a novel algorithm based on the cluster analysis of the gradient method (CA-GM) of lidar signals is developed to capture the multilayer structure and achieve night-time stability. The CA-GM highlights its performance compared with radiosonde data, and the best correlation (0.85), weakest root-mean-square error (203 m), and an improved 25 % correlation coefficient are achieved via the GM. Compared with the 39 d experiments using other algorithms, reasonable parameter selection can help in distinguishing between layers with different properties, such as the cloud layer, elevated aerosol layers, and random noise. Consequently, the CA-GM can automatically address the uncertainty with multiple structures and obtain a stable NBLH with a high temporal resolution, which is expected to contribute to air pollution monitoring and climatology, as well as model verification.

1 Introduction

Air pollution has an important impact on human health, climatic patterns, and the ecological environment (Shi et al., 2019; Su et al., 2020a; Wang et al., 2020). The primary anthropogenic emission source is particulate matter (PM), which is the major source of severe haze in Beijing (Lv et al., 2020; Ma et al., 2019). Many passive and active remote sensing instruments have been combined to observe aerosol optical and microphysical properties (Ji et al., 2018b; Wang et al., 2019), the relationships between PM and meteorology (Li et al., 2019; Zhang et al., 2015), and aerosol–atmospheric boundary layer height (ABLH) interactions (Dong et al., 2017; Su et al., 2020b). With the development in star and moon photometry, continuous day-to-night detection has improved the estimation of column-integrated aerosol properties at night. (Benavent-Oltra et al., 2019; Pérez-Ramírez et al., 2008). Nevertheless, there are a few experiments are observed in the nocturnal boundary layer (NBL). The complexity of weak wind forces, significant stratification, and intermittent turbulence (Stull, 1988; Weil, 2011) results in the continuous accumulation of fine particles near the surface. The turbulent mixing process is accompanied with a strong physiochemical effect, which favours the formation of new particles and worsens the pollution (Hao et al., 2018; Wang et al., 2018). Therefore, accurately acquiring the nocturnal boundary layer height (NBLH) during a polluted episode, especially at night, is of great significance toward combatting air pollution.

Multiple approaches have been developed to determine the ABLH based on various observations, including radio sounding, remote sensing, and parameterization from laboratory experiments (Li et al., 2017b; McGrath-Spangler and Denning, 2012; Nakoudi et al., 2019; Su et al., 2020a). The lidar uses an aerosol as a tracer for mixing processes with high space and temporal resolutions (Kumar, 2006; Leventidou et al., 2013; Yuval et al., 2020). The stable condition shows greater agreement between lidar and radiosonde than the unstable condition because of the complex aerosol structure that complicates NBLH retrieval (Emeis and Schäfer, 2006; Martucci et al., 2007; Sawyer and Li, 2013). At night, the NBLH, as determined by elastic lidar, is either the top of the residual layer or the top of the mechanically driven surface mixing layer (Dang et al., 2019a; Yuval et al., 2020). In the absence of external forces, aerosols in the atmosphere become stratified, resulting in single or multiple layers (elevated or advent aerosols), depending on the location and type of the atmospheric aerosols (Dudeja, 2019; Martucci et al., 2007). A more complex vertical backscatter signal profile can also be formed under specific environmental conditions, such as cloud contamination and local signal noise effect (Dang et al., 2019a; Stull, 1988).

The classical methodologies of lidar-retrieved algorithms are difficult to employ in the identification of multilayer structures in cases of night-time pollution. Gradient-based methods, such as the first-order gradient method (GM) (Hayden et al., 1997), inflexion point method (Menut et al., 1999), logarithm gradient method (Toledo et al., 2017), and cubic-root gradient method (CRGM) (Yang et al., 2017), are sensitive to noisy data unless signal averaging is performed to prevent the loss of some useful instantaneous information. The threshold method is too subjective to set a universal threshold for different weather and terrain (Frioud et al., 2003), while the variance method (Hooper and Eloranta, 1986) is easily affected by lofted aerosol layers and reduces the temporal resolution by calculating the variance profile. The Haar wavelet covariance transform (WCT) (Davis et al., 2000) and the idealized backscatter (Steyn et al., 1999) methods are more robust to noise; however, they can still be affected by low-level clouds and lofted aerosol layers (Caicedo et al., 2017). Some graph theory methods, such as the extended Kalman filter (Banks et al., 2014), Pathfinder and PathfinderTURB (de Bruine et al., 2017; Poltera et al., 2017), k-means clustering (Liu et al., 2018; Toledo et al., 2014), and the STRAT-2D algorithm (Haeffelin et al., 2012) have been proposed to yield promising results via an automated method that reduces the incorrect detection of ABLH. However, these techniques strongly depend on the vertical distribution of particle layers (aerosols and clouds) and are prone to increase the uncertainty under complicated multilayer conditions.

The retrieval of ABLH under cloudy conditions is quite challenging. Some researchers have used the threshold of the attenuated scattering ratio (Campbell et al., 2008; Winker and Vaughan, 1994), the ratio of peaks to the base of the range-corrected signal (RCS) (Wang and Sassen, 2001) to locate cloud tops and bases, while others have employed the objective upper limit of the convective condensation level (CCL) (Li et al., 2017a) and analysed the signal continuity to classify whether the cloud caps within the ABL (Dang et al., 2019b). The height restriction has significant advantages in removing the influence of clouds. Elevated aerosol layers (EALs) are characteristically similar to the aerosol trapped in ABL, using the threshold of lidar backscatter coefficient can distinguish them (Dubovik et al., 2002; Hänel et al., 2012; Peng et al., 2017). More instrument and multi-wavelength lidar systems are combined to obtain more accurate results identify the EALs (Liu et al., 2019; Ortega et al., 2016).

Digressing from these previous efforts to estimate the ABLH, we herein present a new approach – cluster analysis of the gradient method (CA-GM) – to overcome the multilayer structure and remove the noise fluctuation of NBLH with raw data resolution. This study proposes a reasonable parameter to reduce the interference of the cloud layer, EALs, and local noise over the air pollution in megacity regions. The results were evaluated by comparison with the nearby radiosonde site, and they were confirmed through continuous observation via traditional methods in different atmospheric layers.

2 Instruments and datasets

2.1PM2.5 data and radiosonde

Beijing, located in the North China Plain, experienced severe intermittent haze pollution from December 2016 to December 2017. The 39 d lidar and radiosonde data were recorded during that period, and the average concentrations of PM2.5 reached 140 µg m−3. The dataset for lidar, average PM, and air quality index are provided in Sect. S1 in the Supplement. In situ PM2.5 daily measurements in China are primarily obtained from the official website of the China National Environmental Monitoring Centre (CEMC;, last access: 26 November 2020). The radiosonde data are released daily from Nanjiao Station (39.80 N, 116.47 E), which is located southeast of the Beijing Institute of Technology lidar (BIT-lidar) system (39.95 N, 116.32 E). The L-band radiosonde provided a high-resolution profile of temperature, pressure, relative humidity, wind speed, and direction twice a day at 08:00 and 20:00 LST (local standard time) (Guo et al., 2016). The sample temporal resolution is 1.2 s (Zhang et al., 2018), and the vertical resolution is less than 20 m. Previous studies (Hennemuth and Lammert, 2006; Seidel et al., 2012) have adopted the radiosonde as a reference for detecting ABLHs for daily and annual changes in lidar measurements. We resampled the radiosonde data using linear interpolation to achieve the same vertical resolution of lidar and compared it to the 1 h average NBLH centred around the radiosonde launch times. As a result of the complexity of the transition during the morning and early at night, the boundary layer is in a transition between stable and unstable conditions. To determine NBLH from the radiosonde vertical profiles of temperature and humidity, the elevated temperature inversion layer or the height of a significant reduction in moisture is used (Peng et al., 2017). The potential temperature gradient (PTG) should have a good correlation with the relative humidity gradient (RHG), with an allowable error of 100 m (Wang and Wang, 2016). In this study, if the difference between the PTG and RHG is in excess of 100 m, the PTG is considered first, whereas if there is no significant temperature change or the evident changes belong to the cloud or EALs, the result of RHG is referred to as the NBLH.

2.2 BIT-lidar system

A single-wavelength Raman–Mie lidar is operated on the campus of the Beijing Institute of Technology, providing aerosol, cloud, ABLH, and temperature measurements. This lidar system has been continuously enhanced to capture aerosol loading (Chen et al., 2014; Ji et al., 2018a). The standardized RCS is subjected to the correlation lidar factor correction (correction of electronic noise error, background noise error, and overlap factor) and distance correction. The backscatter coefficient can be calculated using the Fernald method (Fernald, 1984), and the assumed lidar ratio is 70 sr (Rosati et al., 2016) owing to the polluted continental aerosol particles. The detailed parameters of the BIT-lidar are listed in Table 1.

Table 1Key parameters of the BIT-lidar.

Download Print Version | Download XLSX

3 Rationale and implementation of the novel algorithm

3.1 Weighted k-means clustering

The NBL shows more complex internal structure at night, the particulate can be used as an important indicator of atmospheric layering because its vertical distribution is strongly affected by the thermal and dynamic structure of the atmosphere (Neff and Coulter, 1986). The assumption of the NBL at which the aerosol concentration and turbulence intensity are significantly higher in the NBL top than in the free atmosphere (FA) (Dang et al., 2019a; Wang et al., 2020). Owing to the influence of the multilayer structure, the minima of RCS gradients are the potential locations of the NBL top. The assembly of these distinguishing peaks with height over time into groups can be considered space- and time-averaged aerosol concentrations. Therefore, it can solve the inadvertent jump between different atmospheric layers. A theoretical schematic of the k-means clustering principle is shown in Fig. 1. To form clusters, the Euclidean distance dis(xi,xj) between two given signal points xi and xj, with coordinates (GMi, hi) and (GMj,hj), is defined according to Eq. (1) as follows:

(1) dis ( x i , x j ) = GM i - GM j 2 + h i - h j 2 1 / 2 ,

where GMi (GMj) is the value of the RCS gradient and hi (hj) is the height of the peak.

Figure 1Theoretical schematic of the k-means clustering. The Euclidean distance dis(xi,xj) between two given signal points xi and xj, with coordinates (GMi, hi) and (GMj, hj), in the cluster Ci and cluster Cj.


Subsequently, we apply the k-means clustering algorithm to classify the datasets with the notable peaks. The cluster number is preset, and the k-means method builds clusters iteratively by moving the centroid until the target function sum of squared errors (SSE) approaches the local minimum (Toledo et al., 2014). The SSE is calculated using Eq. (2) as follows:

(2) SSE = i = 1 k x Ci dis ( c i , x ) 2 ,

where k is the number of clusters and ci and x represent the cluster centroid and all observations in the cluster Ci, respectively. To obtain accurate data for compact and well-separated clusters, the criteria for cluster validation are necessary. The Davies–Bouldin index is employed for the cluster validation analysis and defined as follows:

(3) DB = 1 k i j k max S k ( c i ) + S k ( c j ) S ( c i , c j ) ,

where Sk is the averaged intra-distance between the observations and their cluster centroid and S(cicj) is the distance between cluster centroids ci and cj. The minimum DB index is considered an optimal cluster classification.

Figure 2The theoretical schematic of the weighted k-means clustering. (a) The real profile of a lidar RCS (light grey line) and the hour-averaged RCS (black line). (b) The gradient of RCS (light grey line), the hour-averaged gradient RCS (black line), and the three minima in the profile (yellow points). (c) The distribution of the gradient minima within an hour.(d, e) The results obtained by standard k-means and weighted k-means clustering, where two clusters are differentiated, as shown by hollow and solid red and blue points, respectively.


The standard k-means algorithm must be normalized in cases where the variable is rather different, and data normalization is based on min–max normalization (Virmani et al., 2015). The normalized k-means clustering is “isotropic” in all directions of space, and it tends to capture a spherical shape. Nevertheless, herein we propose putting weight on height and excluding variances that are greater in height. Therefore, the assembling groups of the distribution tend to be separated along variables with greater variances, which is conducive toward setting the upper and lower limiter altitude to classify different atmospheric layers vertically (Fig. 2b). The weighted G is calculated by the difference between the maximum altitude hmax and the minimum altitude, hmin, as shown in Eq. (4), while the weighted height hw is rescaled by the normalized height data hnor, as presented in Eq. (5) as follows:


3.2 Multilayer classification

Because of the presence of the strong gradient signature of the backscatter profile, a dataset of three minima of RCS gradient within an hour works as the dataset of k-means classification (Fig. 2e). The three minima are calculated by a window of 25 m, and selected in orders. The cloud layer (CL) has a larger gradient magnitude of extinction and backscattering coefficient than the aerosol layers (Palm et al., 2012). Additionally, the typical nocturnal clouds are shallow cumulus, stratocumulus, and stratus (Kotthaus and Grimmond, 2018). They have shallow vertical dimensions and are denser than aerosols at the same altitude; hence, they can be distinguished from the aerosol layer (Wang and Sassen, 2001). Meanwhile, the accuracy of the NBLH from GM can be affected by background and electronic noise: it has a non-regular distribution and appears at higher altitudes with lower signal-to-noise ratios. The noise layer lacks a stratified structure but has a GM value similar to that of the lower height. Thus, we calculated the range of vertical extension of different layers, indicating the cluster significance of the noise and other layers. As for the EALs, their presence above the NBL represents a difficulty when retrieving the upper height of the NBL, particularly when the EALs are close to it. Both aerosol layers have a similar characteristic of gradient variance and range of height, which we discover by seeking the empirical threshold value of the EALs in the backscatter coefficient (Hänel et al., 2012). The typical backscatter threshold for a 532 nm wavelength lidar is defined as βth= 1.786 × 10−3km-1sr-1, which is calculated using the Ångstöm parameter as 1.2 under urban-industrial and mixed conditions (Dubovik et al., 2002). The gaps between NBL and EALs in the multilayer structure are determined by Dth= 100 m (Peng et al., 2017).

Figure 3Flowchart of the retrieval method for the CA-GM. hCi(hCj) is the height of centroid of cluster Ci(Cj). Dinter is the inter-cluster distance between minimum hCi and maximum hCj. RhCi is the intra-cluster range from the minimum hCi and maximum hCi.SCi(SCj) is the SD of the GM in the cluster Ci(Cj). Vi(Vj) is the vertical uniformity calculated by RCi/NCi(RCj/NCj), the NCi(NCj) is the amount of peak in the group Ci(Cj). βth is the typical backscatter aerosol layer (1.786 × 10−3km-1sr-1). Dth is the threshold of distance to define a gap between multiple aerosol layers (100 m). Dsig is the empirical threshold, set as 50 m. C is the empirical value, set as 1.5. hnbl is the final location of the nocturnal boundary layer height.


3.3 Implementation of the CA-GM algorithm

The CA-GM method, which is based on the k-means clustering analysis of different types of atmospheric layers, is generally used to retrieve multiple layers in polluted cases. The specific ideas are shown in the flowchart in Fig. 3, and the specific steps are as follows (detailed results are presented as a case study).

The algorithm is divided into three parts: pre-processing, layer attribution classification, and NBL inspection. The CA-GM algorithm is implemented if the data collection exceeds 30 min within an hour period. First, the standardized lidar RCS is applied to a Savitzky–Golay filter for preliminary de-noising. The profile of the backscatter coefficient (β) is calculated, and the reference height (href) is limited by the Fernald method as the theoretical height limiter (Comerón et al., 2017; Ji et al., 2017). Notably, G is a dataset of three gradient minima of the RCS. The cluster is preset as a pair, and k-means clustering is carried out once to seek the minimum DB index as the optimal grouping, Ci and Cj. Second, there is a parameter Dinter that is defined as the minimum inter-cluster distance, which can measure the cluster-stratified significance to classify the cloud and noise mixed in G. If Dinter exceeds the threshold Dsig, it can distinguish the noise from other layers. Dsig is the empirical value to distinguish noise layer for verified stratification. Furthermore, SCi and SCj are a quality control function for noise layer attribution. For the cloud layer, the vertical extension RhCi (RhCj) of the cloud is lower than the aerosol layers; therefore, we define an empirical constant C for this study. In addition, the vertical uniformity parameter Vi(Vj) works as a quality control tool for the features of the cloud and other layers. If a cloud layer or noise exists, the original G is removed from the upper limiter as hcloud or hnoise, respectively. After the elimination of cloud and noise interference, the EALs can be determined from the typical aerosol layer, βth, and the gap distance, Dth. Finally, the new dataset, G, which has been removed from the different attributed layer, goes to the final step with the cluster as Ci and Cj (or Ci and Cj). Owing to the assumption of NBL distribution, the largest deviation of cluster indicates the location of the NBLH (hnbl).

In summary, by k-means clustering analysis of the vertical–temporal gradient of the GM once or twice within an hour, the multilayer NBL structure can be separated according to the physical characteristics of its different layers. The CA-GM method is an objective and robust method for judging the attribution of different layers (NBL, EALs, and CL) and noise.

Table 2Statistical parameters of the lidar-retrieved algorithm compared with radiosonde measurement. Mean bias (MB), correlation coefficient (R), root-mean-square deviation (RMSE), and the percent of relative absolute bias difference (PRD) are shown below.

Download Print Version | Download XLSX

Figure 4Comparison between the radiosonde-determined and lidar-retrieved NBLH measurements via the gradient method (GM, red circle), wavelet covariance transform transition method (WCT, blue triangle), cubic root gradient method (CRGM, orange star), and cluster analysis of gradient method (CA-GM, black circle). The correlation coefficient is represented by R. The solid black line is the 1:1 line.


4 Evaluation and comparative analysis with classical methods

4.1 Evaluation with radiosonde data

The L-band radiosonde provided accurate thermodynamic profiles, and the radiosonde-determined NBLHs were used to evaluate the accuracy of the lidar-retrieved NBLHs. Compared with the two-moment radiosonde with the other three algorithm, it was found that the correlation coefficients (R) ranged from 0.68 to 0.85. The CA-GM had the highest consistency among the classical methods, with the highest correlation coefficient (0.85), the weakest root-mean-square error (RMSE) (203 m), the smaller mean bias (28 m), and the minimum mean relative absolute difference (PRD) (17 %) (Table 2).

The NBLH retrieved by GM and CA-GM (Fig. 4) had a good correlation with the radiosonde approach, and the latter method enhanced the correlation coefficient by 25 %. With the implementation of CA-GM, the data were concentrated, and the RMSE was reduced from 292 to 203 m (Table 2). The means bias of GM is greater than that of the CA-GM, corresponding to the decrease in PRD from GM to CA-GM. Additionally, compared with the WCT and CRGM, the former underestimated the NBLH by approximately 13 m, whereas the latter overestimated the altitude by 186 m. The RMSE of CA-GM is less than that of WCT and CRGM, which is similar to the PRD result. Therefore, the CA-GM showed a good correlation with the radiosonde method, and evinced the least fluctuation and highest consistency in NBLH retrieval.

4.2 Comparison with other classical methods

To enrich our analysis, a comparison of CA-GM with GM, WCT, and CRGM in the 39 d night-time period was applied to compensate for the rare temporal resolution of the radiosonde approach. The results are shown in Fig. 5.

Figure 5Correlation coefficient and RMSE results compared with the CA-GM method under all conditions (see text for details) using the gradient method (GM), wavelet covariance transform transition method (WCT), and cubic root gradient method (CRGM). The sample number is shown at the top of the column, and the condition is represented by the x axis.


Valid CA-GM data were implemented for a total of 256 h, and the data were analysed for comparison with other retrieval algorithms. Under conditions without the interference of multiple layers, the CA-GM had a good correlation of 0.90, 0.70, and 0.82 for GM, WCT, and CRGM, respectively, and the RMSE was the lowest compared to other situations. Consequently, the CA-GM was more similar to the other three methods in cases without a multilayer structure, which proved the feasibility of the CA-GM relative to the classical boundary layer retrieval methods.

Moreover, the extensive results showed that the WCT method was more accurate than the GM during the night (Caicedo et al., 2017), and it was less affected by the low signal-to-noise ratio condition (Brooks, 2003). The dilation and threshold of the WCT method were selected carefully in this study (Mao et al., 2013); thus, the performance of the WCT could ensure the identification of the noise and most of the cloud layers. Notably, compared with the consistency for the WCT to the CA-GM, the improvement of the correlation coefficient from 0.70 to 0.84 in cloud contamination and from 0.70 to 0.72 in noise effect was observed, which proves its ability to remove the attributed layers. Although the fluctuations in noise and cloud layers were relatively large, the CA-GM exhibited an outstanding ability for cloud removal to eliminate noise. As for EALs, because of their ambiguous cluster, as well as the NBL, all the methods had poor correlation coefficients with the CA-GM. Observing EALs is the most challenging part in multilayer structures; hence, more active remote sensing instruments (such as multi-wavelength lidar and polarized lidar), as well as additional methods, are required to determine the accurate layout of EALs.

Table 3 presents the criterion parameters in the CA-GM. The cluster-significant parameter Dinter for noise was 20.75 ± 14.62 m, which was significantly lower compared to other conditions. The typical altitude of NBL, EALs, cloud, and noise in severe haze pollution is 590.49 ± 202.84, 1024.69 ± 166.36, 1252.52 ± 303.28, and 1100.66 ± 253.04 m, respectively. The vertical extension of the cloud layer was shallower than the other layer, with a typical extension of 128.6 ± 82.13 m. The backscatter coefficient of EALs was 1.12 ± 0.76 × 10−3km-1sr-1, which was evidence for choosing a suitable empirical βth value. The cloud had the smallest value in vertical uniformity, which indicated a denser peak distribution than other layers.

Table 3Computed criteria parameters for layer attribution.

Download Print Version | Download XLSX

Figure 6Time–height cross section of range-corrected signal (RCS) with four NBLH retrieved methods on 23 December 2016.


Figure 7Distribution of altitude and normalized gradient method (GM) values at 18:00–21:00 LST. Panels (a–c) indicate hourly intervals.


4.3 Case study with a multilayer structure

4.3.1 Effects of cloud contamination

On 23 December 2016, there was a cloud layer that was 1.3 kma.g.l. (above ground level) between 18:00–23:00 LST (Fig. 6), which is presented as a light blue region in Fig. 6. Below the cloud base, there was a distinct aerosol layer surface and a strong signal negative gradient, indicating the WCT method capture. The cloud significantly influences the GM and CRGM determination and captures the upper edge of the cloud. After 21:00 LST, the cloudiness decreases and the lidar can capture the NBL signal. After defining the minimum in the upper cluster (Ci) as the top limiter altitude, the CA-GM captured slowly increased the NBLH, as shown in Fig. 6. Figure 7 shows the significant two-layer structure distribution hourly for the first k-means clustering distribution. The centroid of the two layers indicated the approximate location at 839 and 1428 m (Fig. 7b), and the cloud located at the upper layer, which had a shallow vertical extension and a relatively dense distribution. The radiosonde measurement had a good correlation with the lidar-retrieved NBLH in Fig. 8. The PTG exhibited the steepest slope at 1.37 km, but it corresponded to the height at the cloud location. Therefore, we selected NBLH, using the RHG method, at 0.78 km, which was less than the CA-GM retrieved height at 20:00 LST.

Figure 8Planetary boundary layer height estimates using radiosonde. Profiles include temperature, potential temperature, relative humidity, potential temperature gradient (PTG), and relative humidity gradient (RHG). Estimated NBLH by PTG (red) and RHG (blue) are shown by dashed horizontal lines.


Figure 9Time–height cross section of RCS with four NBLH retrieval methods on 6 April 2017.


4.3.2 Noise effect

On 6 April 2017, the noise distribution was prone to appear when the low-load aerosol was utilized for the GM. The gradient-based methods were affected by noise and had a wide range of fluctuations (Fig. 9). Conversely, the WCT adequately captured the edge of the aerosol concentration. From the distribution of the GM with height distribution, Fig. 10 shows evident mixing without a stratified layer structure. Therefore, the noise was mixed in the upper layer of the centroid at 1479 and 1452 m at 19:00 and 20:00, which set the upper limiter and recalculated the NBLH in an hourly manner, due to the SD not meeting the requirement of the algorithm (Sci<Scj) at 21:00 UTC. Therefore, the NBL is in the cluster of the upper layer (Sci=0.016, Scj=0.033). The radiosonde data were calculated through the rapid change of the PTG method as 0.79 km (Fig. 11), corresponding to the height retrieval by using CA-GM as 0.74 km.

Figure 10Distribution of altitude and the normalized gradient method (GM) value during 19:00–22:00 LST.


Figure 11Planetary boundary layer height estimates using radiosondes. Profiles include temperature, potential temperature, relative humidity, potential temperature gradient (PTG), and relative humidity gradient (RHG). Estimated NBLH by PTG (red) and RHG (blue) are shown by dashed horizontal lines at 20:00 LST.


4.3.3 Nocturnal aloft aerosol layer

On 2 January 2017, the EALs appeared frequently in the lower troposphere. There was a distinct aerosol layer between 0.7 and 1.2 kma.g.l. between 17:00 and 22:00 LST (Fig. 12). Without any limitation, the GM, CRGM, and WCT captured the height of the EALs when the negative gradient signal at the EALs was stronger than the NBL, corresponding to the lofted aerosol structure from 17:00 to 22:00 LST. As shown in Fig. 13, the two distinct peaks of the cluster were the two aerosol layers; the deviation of the upper layer was larger initially, and both layers gradually exhibited approximately the same gradient magnitude as the time transition. The upward centroid of the upper cluster provides additional evidence for the NBLH with topped EALs. After using the CA-GM method to limit the base of the lofting aerosol layers, the effect of EALs in the polluted cases can be successfully separated. Similarly, in Fig. 14, the first gradient maxima above the surface inversion layer is the NBLH, and both PTG and RHG showed good consistency, while the NBLH was at 0.44 km. The other peak with PTG and RHG corresponded to the height of the EALs.

Figure 12Time–height cross section of RCS with four NBLH retrieval methods. The discontinuity of the RCS at 18:06–18:07 LST is the result of detecting electric noise. The discontinuities of RCS at 20:39–20:58 LST and 23:18–23:39 LST were because of the laser energy adjustment and the signal test.


Figure 13Distribution of altitude and the normalized gradient method (GM) value during 16:00–23:00 LST.


Figure 14Planetary boundary layer height estimates using radiosondes. Profiles include temperature, potential temperature, relative humidity, potential temperature gradient (PTG), and relative humidity gradient (RHG). Estimated NBLH by PTG (red) and RHG (blue) are shown by dashed horizontal lines at 20:00 LST.


5 Discussion and conclusion

Elastic lidars are excellent instruments to determine the NBLH with high space and vertical resolution. Multilayer structures in severely polluted cases impede buoyancy forces and influence pollutant dispersion and dilution. Herein, a novel CA-GM algorithm was developed to capture the multilayer structure and achieve stability at night with raw resolution. A 39 d heavily polluted observation experiment over Beijing (China) thoroughly established the limitations of the current methods employed for boundary layer height determination, and a suitable algorithm for pollution conditions was developed.

Overall, the CA-GM method highlights its high performance relative to the radiosonde approach; the best correlation (0.85), weakest RMSE (203 m), and an improved 25 % correlation coefficient of the GM were established. The possible deviations are due to the different definitions of thermodynamic NBLs from radiosondes and aerosol NBLs. The sound data are also multilayered because of the effect of the aerosol and cloud layers, and the radiosonde-retrieved NBLH combines the PTG and RHG methods to discuss the uncertainty of NBLs in the pollution period. The calculation of the three minimum gradients can be used to determine the potential stratified layer structure, which provides a worst case scenario for estimating the surface concentrations of pollutants released into the NBL. Compared with the 39 d performance of other algorithms, a reasonable parameter selection can distinguish different atmospheric layers, such as cloud layer, elevated (or advected) aerosol, and random noise. The Dinter, RhC, β, and Vi provide a novel idea of classifying multiple layers based on their physical characteristics, which is more objective for automatic clustering under complex conditions.

The correlation coefficient with the CA-GM and WCT had an elevated correlation coefficient from 0.7 to 0.84 and 0.7 to 0.72 in cloud and noise effect, which proved the ability of the CA-GM to ensure the upper edge of the low-level cloud and remove the random noise. The EALs are often located at the top of NBLs, with a similar characteristic to the NBLs. Thus, using the empirical threshold on a single-wavelength elastic lidar is a good way to classify EALs in polluted cases. Consequently, the CA-GM approach can deal with the uncertainty of the multilayered structure and obtain a stable NBLH with a high temporal resolution, which is expected to contribute to the long-term observation of the single-wavelength lidar system and micro-pulse lidar monitoring in air pollution.

The resolution of the CA-GM is at a high resolution compared to previous studies (Martucci et al., 2010; Su and Patrick McCormick, 2019; Tsaknakis et al., 2011). The averaged times for an elastic lidar system are 2, 15, or 30 min, and it will lose the raw time resolution in tracing the aerosol distribution. However, the CA-GM are taking into account the overall set of observations in the effective time (which at least exceed 30 min) and maintain the raw resolution of the data.

The uncertainty of the CA-GM is calculated by the real signal provided in Sect. S2 in the Supplement. Concerning the robustness of the CA-GM approach, the effect of the lidar RCS noise in determining the NBLH has been analysed. Unlike other gradient derivation methods, CA-GM results are slightly affected by lidar signal noise. NBL top height as obtained for “noised” lidar RCS with an additive Gaussian noise coefficient α< 4 % is better than GM results. The intensity of the CLs changes (± 40 %) will not affect with the classification of the CA-GM in the polluted cases; the significantly stratified structure is related to the relative difference on the backscatter signal. As for the EALs, the strict threshold will define the EALs accurately. The limitation of the CA-GM based on the nocturnal boundary layer is stable, hence, we can calculate the distribution of the minima RCS gradients at an hour interval to use weighted k-means clustering for height restriction of the layers. Secondly, based on the limitation of the lidar system, the lower limit of the BIT-lidar is around 300 m. Shallow nocturnal boundary layer heights are not detectable. Thirdly, the method should be used under high signal-to-noise ratio conditions, such as during night-time and where air pollution is present.

Code and data availability

Data can be accessed by contacting Siying Chen (


The supplement related to this article is available online at:

Author contributions

YZ and SC analysed the experimental data and co-wrote the paper. HC and PG supported the experiment and the maintenance of the BIT-lidar. SYC designed and led the study. All co-authors discussed the results and commented on the manuscript.

Competing interests

The authors declare that they have no conflict of interest.


We wish to thank the China Meteorological Administration (CMA) for providing the radiosonde data.

Financial support

This research has been supported by the National Natural Science Foundation of China (grant no. 61505009).

Review statement

This paper was edited by Daniel Perez-Ramirez and reviewed by Yi Yang and two anonymous referees.


Banks, R. F., Tiana-Alsina, J., María Baldasano, J., and Rocadenbosch, F.: Retrieval of boundary layer height from lidar using extended Kalman filter approach, classic methods, and backtrajectory cluster analysis, edited by: Comerón, A., Kassianov, E. I., Schäfer, K., Picard, R. H., Stein, K., and Gonglewski, J. D., Amsterdam, the Netherlands, p. 92420F, 2014. 

Benavent-Oltra, J. A., Román, R., Casquero-Vera, J. A., Pérez-Ramírez, D., Lyamani, H., Ortiz-Amezcua, P., Bedoya-Velásquez, A. E., de Arruda Moreira, G., Barreto, Á., Lopatin, A., Fuertes, D., Herrera, M., Torres, B., Dubovik, O., Guerrero-Rascado, J. L., Goloub, P., Olmo-Reyes, F. J., and Alados-Arboledas, L.: Different strategies to retrieve aerosol properties at night-time with the GRASP algorithm, Atmos. Chem. Phys., 19, 14149–14171,, 2019. 

Brooks, I. M.: Finding boundary layer top: Application of a wavelet covariance transform to lidar backscatter profiles, J. Atmos. Ocean. Tech., 20, 1092–1105, 2003. 

Caicedo, V., Rappenglück, B., Lefer, B., Morris, G., Toledo, D., and Delgado, R.: Comparison of aerosol lidar retrieval methods for boundary layer height detection using ceilometer aerosol backscatter data, Atmos. Meas. Tech., 10, 1609–1622,, 2017. 

Campbell, J. R., Sassen, K., and Welton, E. J.: Elevated Cloud and Aerosol Layer Retrievals from Micropulse Lidar Signal Profiles, J. Atmos. Ocean. Tech., 25, 685–700,, 2008. 

Chen, H., Chen, S., Zhang, Y., Chen, H., Guo, P., and Chen, B.: Experimental determination of Raman lidar geometric form factor combining Raman and elastic return, Opt. Commun., 332, 296–300, 2014. 

Comerón, A., Muñoz-Porcar, C., Rocadenbosch, F., Rodríguez-Gómez, A., and Sicard, M.: Current Research in Lidar Technology Used for the Remote Sensing of Atmospheric Aerosols, Sensors, 17, 1450,, 2017. 

Dang, R., Yang, Y., Hu, X.-M., Wang, Z., and Zhang, S.: A Review of Techniques for Diagnosing the Atmospheric Boundary Layer Height (ABLH) Using Aerosol Lidar Data, Remote Sens.-Basel, 11, 1590,, 2019a. 

Dang, R., Yang, Y., Li, H., Hu, X.-M., Wang, Z., Huang, Z., Zhou, T., and Zhang, T.: Atmosphere Boundary Layer Height (ABLH) Determination under Multiple-Layer Conditions Using Micro-Pulse Lidar, Remote Sens.-Basel, 11, 263,, 2019b. 

Davis, K. J., Gamage, N., Hagelberg, C. R., Kiemle, C., Lenschow, D. H., and Sullivan, P. P.: An objective method for deriving atmospheric structure from airborne lidar observations, J. Atmos. Ocean. Tech., 17, 1455–1468, 2000. 

de Bruine, M., Apituley, A., Donovan, D. P., Klein Baltink, H., and de Haij, M. J.: Pathfinder: applying graph theory to consistent tracking of daytime mixed layer height with backscatter lidar, Atmos. Meas. Tech., 10, 1893–1909,, 2017. 

Dong, Z., Li, Z., Yu, X., Cribb, M., Li, X., and Dai, J.: Opposite long-term trends in aerosols between low and high altitudes: a testimony to the aerosol–PBL feedback, Atmos. Chem. Phys., 17, 7997–8009,, 2017. 

Dubovik, O., Holben, B., Eck, T. F., Smirnov, A., Kaufman, Y. J., King, M. D., Tanré, D., and Slutsker, I.: Variability of Absorption and Optical Properties of Key Aerosol Types Observed in Worldwide Locations, J. Atmos. Sci., 59, 590–608,<0590:VOAAOP>2.0.CO;2, 2002. 

Dudeja, J. P.: Micro-Pulse Lidar for the Determination of Atmospheric Boundary Layer Height, Int. J. Res. Anal. Rev., 6, 810–817, 2019. 

Emeis, S. and Schäfer, K.: Remote Sensing Methods to Investigate Boundary-layer Structures relevant to Air Pollution in Cities, Bound.-Lay. Meteorol., 121, 377–385,, 2006. 

Fernald, F. G.: Analysis of atmospheric lidar observations: some comments, Appl. Optics, 23, 652,, 1984. 

Frioud, M., Mitev, V., Matthey, R., Häberli, C., Richner, H., Werner, R., and Vogt, S.: Elevated aerosol stratification above the Rhine Valley under strong anticyclonic conditions, Atmos. Environ., 37, 1785–1797,, 2003. 

Guo, J., Miao, Y., Zhang, Y., Liu, H., Li, Z., Zhang, W., He, J., Lou, M., Yan, Y., Bian, L., and Zhai, P.: The climatology of planetary boundary layer height in China derived from radiosonde and reanalysis data, Atmos. Chem. Phys., 16, 13309–13319,, 2016. 

Haeffelin, M., Angelini, F., Morille, Y., Martucci, G., Frey, S., Gobbi, G. P., Lolli, S., O'Dowd, C. D., Sauvage, L., Xueref-Rémy, I., Wastine, B., and Feist, D. G.: Evaluation of Mixing-Height Retrievals from Automatic Profiling Lidars and Ceilometers in View of Future Integrated Networks in Europe, Bound.-Lay. Meteorol., 143, 49–75,, 2012. 

Hänel, A., Baars, H., Althausen, D., Ansmann, A., Engelmann, R., and Sun, J. Y.: One-year aerosol profiling with EUCAARI Raman lidar at Shangdianzi GAW station: Beijing plume and seasonal variations: One-year Aerosol Profiling Near Beijing, J. Geophys. Res.-Atmos., 117,, 2012. 

Hao, L., Garmash, O., Ehn, M., Miettinen, P., Massoli, P., Mikkonen, S., Jokinen, T., Roldin, P., Aalto, P., Yli-Juuti, T., Joutsensaari, J., Petäjä, T., Kulmala, M., Lehtinen, K. E. J., Worsnop, D. R., and Virtanen, A.: Combined effects of boundary layer dynamics and atmospheric chemistry on aerosol composition during new particle formation periods, Atmos. Chem. Phys., 18, 17705–17716,, 2018. 

Hayden, K. L., Anlauf, K. G., Hoff, R. M., Strapp, J. W., Bottenheim, J. W., Wiebe, H. A., Froude, F. A., Martin, J. B., Steyn, D. G., and McKendry, I. G.: The vertical chemical and meteorological structure of the boundary layer in the Lower Fraser Valley during Pacific '93, Atmos. Environ., 31, 2089–2105,, 1997. 

Hennemuth, B. and Lammert, A.: Determination of the Atmospheric Boundary Layer Height from Radiosonde and Lidar Backscatter, Bound.-Lay. Meteorol., 120, 181–200,, 2006. 

Hooper, W. P. and Eloranta, E. W.: Lidar measurements of wind in the planetary boundary layer: the method, accuracy and results from joint measurements with radiosonde and kytoon, J. Clim. Appl. Meteorol., 25, 990–1001, 1986. 

Ji, H., Chen, S., Zhang, Y., Chen, H., Guo, P., and Chen, H.: Calibration method for the reference parameter in Fernald and Klett inversion combining Raman and Elastic return, J. Quant. Spectrosc. Ra., 188, 71–78,, 2017. 

Ji, H., Zhang, Y., Chen, S., Chen, H., and Guo, P.: Aerosol characteristics inversion based on the improved lidar ratio profile with the ground-based rotational Raman–Mie lidar, Opt. Commun., 416, 54–60,, 2018a. 

Ji, H., Chen, S., Zhang, Y., Chen, H., Guo, P., and Zhao, P.: Comparison of air quality at different altitudes from multi-platform measurements in Beijing, Atmos. Chem. Phys., 18, 10645–10653,, 2018b. 

Kotthaus, S. and Grimmond, C. S. B.: Atmospheric boundary-layer characteristics from ceilometer measurements. Part 1: A new method to track mixed layer height and classify clouds, Q. J. Roy. Meteor. Soc., 144, 1525–1538,, 2018. 

Kumar, Y. B.: Portable lidar system for atmospheric boundary layer measurements, Opt. Eng., 45, 076201,, 2006. 

Leventidou, E., Zanis, P., Balis, D., Giannakaki, E., Pytharoulis, I., and Amiridis, V.: Factors affecting the comparisons of planetary boundary layer height retrievals from CALIPSO, ECMWF and radiosondes over Thessaloniki, Greece, Atmos. Environ., 74, 360–366,, 2013. 

Li, H., Yang, Y., Hu, X.-M., Huang, Z., Wang, G., and Zhang, B.: Application of Convective Condensation Level Limiter in Convective Boundary Layer Height Retrieval Based on Lidar Data, Atmosphere, 8, 79,, 2017a. 

Li, H., Yang, Y., Hu, X.-M., Huang, Z., Wang, G., Zhang, B., and Zhang, T.: Evaluation of retrieval methods of daytime convective boundary layer height based on lidar data: MEASUREMENT OF BOUNDARY LAYER HEIGHT, J. Geophys. Res.-Atmos., 122, 4578–4593,, 2017b. 

Li, X., Song, H., Zhai, S., Lu, S., Kong, Y., Xia, H., and Zhao, H.: Particulate matter pollution in Chinese cities: Areal-temporal variations and their relationships with meteorological conditions (2015–2017), Environ. Pollut., 246, 11–18,, 2019. 

Liu, B., Ma, Y., Liu, J., Gong, W., Wang, W., and Zhang, M.: Graphics algorithm for deriving atmospheric boundary layer heights from CALIPSO data, Atmos. Meas. Tech., 11, 5075–5085,, 2018. 

Liu, B., Ma, Y., Gong, W., Zhang, M., and Yang, J.: Improved two-wavelength Lidar algorithm for retrieving atmospheric boundary layer height, J. Quant. Spectrosc. Ra., 224, 55–61,, 2019. 

Lv, Z., Wei, W., Cheng, S., Han, X., and Wang, X.: Mixing layer height estimated from AMDAR and its relationship with PMs and meteorological parameters in two cities in North China during 2014–2017, Atmos. Pollut. Res., 11, 443–453,, 2020. 

Ma, X., Wang, C., Han, G., Ma, Y., Li, S., Gong, W., and Chen, J.: Regional Atmospheric Aerosol Pollution Detection Based on LiDAR Remote Sensing, Remote Sens.-Basel, 11, 2339,, 2019. 

Mao, F., Gong, W., Song, S., and Zhu, Z.: Determination of the boundary layer top from lidar backscatter profiles using a Haar wavelet method over Wuhan, China, Opt. Laser Technol., 49, 343–349,, 2013. 

Martucci, G., Matthey, R., Mitev, V., and Richner, H.: Comparison between backscatter lidar and radiosonde measurements of the diurnal and nocturnal stratification in the lower troposphere, J. Atmos. Ocean. Tech., 24, 1231–1244, 2007. 

Martucci, G., Milroy, C., and O'Dowd, C. D.: Detection of Cloud-Base Height Using Jenoptik CHM15K and Vaisala CL31 Ceilometers, J. Atmos. Ocean. Tech., 27, 305–318,, 2010. 

McGrath-Spangler, E. L. and Denning, A. S.: Estimates of North American summertime planetary boundary layer depths derived from space-borne lidar: PBL DEPTH ESTIMATES FROM CALIPSO LIDAR, J. Geophys. Res., 117, D15101,, 2012. 

Menut, L., Flamant, C., Pelon, J., and Flamant, P. H.: Urban boundary-layer height determination from lidar measurements over the Paris area, Appl. Optics, 38, 945,, 1999. 

Nakoudi, K., Giannakaki, E., Dandou, A., Tombrou, M., and Komppula, M.: Planetary boundary layer height by means of lidar and numerical simulations over New Delhi, India, Atmos. Meas. Tech., 12, 2595–2610,, 2019. 

Neff, W. D. and Coulter, R. L.: Acoustic Remote Sensing, in: Probing the Atmospheric Boundary Layer, edited by: Lenschow, D. H., American Meteorological Society, Boston, MA, 201–239, 1986. 

Ortega, I., Berg, L. K., Ferrare, R. A., Hair, J. W., Hostetler, C. A., and Volkamer, R.: Elevated aerosol layers modify the O2O2 absorption measured by ground-based MAX-DOAS, J. Quant. Spectrosc. Ra., 176, 34–49,, 2016. 

Palm, S. P., Hart, W. D., Hlavka, D. L., Welton, E. J., and Spinhirne, J. D.: The Algorithm Theoretical Basis Document for the GLAS Atmospheric Data Products, NASA, NASA/TM-2012-208641/Vol 6, 148 pp., 2012. 

Peng, J., Grimmond, C. S. B., Fu, X., Chang, Y., Zhang, G., Guo, J., Tang, C., Gao, J., Xu, X., and Tan, J.: Ceilometer-Based Analysis of Shanghai's Boundary Layer Height (under Rain- and Fog-Free Conditions), J. Atmos. Ocean. Tech., 34, 749–764,, 2017. 

Pérez-Ramírez, D., Ruiz, B., Aceituno, J., Olmo, F. J., and Alados-Arboledas, L.: Application of Sun/star photometry to derive the aerosol optical depth, Int. J. Remote Sens., 29, 5113–5132,, 2008. 

Poltera, Y., Martucci, G., Collaud Coen, M., Hervo, M., Emmenegger, L., Henne, S., Brunner, D., and Haefele, A.: PathfinderTURB: an automatic boundary layer algorithm. Development, validation and application to study the impact on in situ measurements at the Jungfraujoch, Atmos. Chem. Phys., 17, 10051–10070,, 2017. 

Rosati, B., Herrmann, E., Bucci, S., Fierli, F., Cairo, F., Gysel, M., Tillmann, R., Größ, J., Gobbi, G. P., Di Liberto, L., Di Donfrancesco, G., Wiedensohler, A., Weingartner, E., Virtanen, A., Mentel, T. F., and Baltensperger, U.: Studying the vertical aerosol extinction coefficient by comparing in situ airborne data and elastic backscatter lidar, Atmos. Chem. Phys., 16, 4539–4554,, 2016. 

Sawyer, V. and Li, Z.: Detection, variations and intercomparison of the planetary boundary layer depth from radiosonde, lidar and infrared spectrometer, Atmos. Environ., 79, 518–528,, 2013. 

Seidel, D. J., Zhang, Y., Beljaars, A., Golaz, J.-C., Jacobson, A. R., and Medeiros, B.: Climatology of the planetary boundary layer over the continental United States and Europe: Boundary Layer Climatology: U.S. and Europe, J. Geophys. Res.-Atmos., 117,, 2012. 

Shi, Y., Hu, F., Fan, G., and Zhang, Z.: Multiple technical observations of the atmospheric boundary layer structure of a red-alert haze episode in Beijing, Atmos. Meas. Tech., 12, 4887–4901,, 2019. 

Steyn, D. G., Baldi, M., and Hoff, R. M.: The Detection of Mixed Layer Depth and Entrainment Zone Thickness from Lidar Backscatter Profiles, J. Atmos. Ocean. Tech., 16, 953–959,<0953:TDOMLD>2.0.CO;2, 1999. 

Stull, R. B.: An Introduction to Boundary Layer Meteorology, Kluwer Academic Publishers, Dordrecht, 1988. 

Su, J. and Patrick McCormick, M.: Using multi-wavelength Mie–Raman lidar to measure low-level cloud properties, J. Quant. Spectrosc. Ra., 237, 106610,, 2019. 

Su, T., Li, Z., and Kahn, R.: A new method to retrieve the diurnal variability of planetary boundary layer height from lidar under different thermodynamic stability conditions, Remote Sens. Environ., 237, 111519,, 2020a. 

Su, T., Li, Z., Li, C., Li, J., Han, W., Shen, C., Tan, W., Wei, J., and Guo, J.: The significant impact of aerosol vertical structure on lower atmosphere stability and its critical role in aerosol–planetary boundary layer (PBL) interactions, Atmos. Chem. Phys., 20, 3713–3724,, 2020b. 

Toledo, D., Córdoba-Jabonero, C., and Gil-Ojeda, M.: Cluster Analysis: A New Approach Applied to Lidar Measurements for Atmospheric Boundary Layer Height Estimation, J. Atmos. Ocean. Tech., 31, 422–436,, 2014. 

Toledo, D., Córdoba-Jabonero, C., Adame, J. A., De La Morena, B., and Gil-Ojeda, M.: Estimation of the atmospheric boundary layer height during different atmospheric conditions: a comparison on reliability of several methods applied to lidar measurements, Int. J. Remote Sens., 38, 3203–3218,, 2017. 

Tsaknakis, G., Papayannis, A., Kokkalis, P., Amiridis, V., Kambezidis, H. D., Mamouri, R. E., Georgoussis, G., and Avdikos, G.: Inter-comparison of lidar and ceilometer retrievals for aerosol and Planetary Boundary Layer profiling over Athens, Greece, Atmos. Meas. Tech., 4, 1261–1273,, 2011. 

Virmani, D., Taneja, S., and Malhotra, G.: Normalization based K means Clustering Algorithm, arXiv:1503.00900, 2015. 

Wang, H., Lu, K., Chen, X., Zhu, Q., Wu, Z., Wu, Y., and Sun, K.: Fast particulate nitrate formation via N2O5 uptake aloft in winter in Beijing, Atmos. Chem. Phys., 18, 10483–10495,, 2018. 

Wang, H., Li, Z., Lv, Y., Xu, H., Li, K., Li, D., Hou, W., Zheng, F., Wei, Y., and Ge, B.: Observational study of aerosol-induced impact on planetary boundary layer based on lidar and sunphotometer in Beijing, Environ. Pollut., 252, 897–906,, 2019. 

Wang, H., Li, Z., Lv, Y., Zhang, Y., Xu, H., Guo, J., and Goloub, P.: Determination and climatology of the diurnal cycle of the atmospheric mixing layer height over Beijing 2013–2018: lidar measurements and implications for air pollution, Atmos. Chem. Phys., 20, 8839–8854,, 2020. 

Wang, X. and Wang, K.: Homogenized Variability of Radiosonde-Derived Atmospheric Boundary Layer Height over the Global Land Surface from 1973 to 2014, J. Climate, 29, 6893–6908,, 2016. 

Wang, Z. and Sassen, K.: Cloud Type and Macrophysical Property Retrieval Using Multiple Remote Sensors, J. Appl. Meteorol., 40, 1665–1682,<1665:CTAMPR>2.0.CO;2, 2001.  

Weil, J. C.: Stable Boundary Layer Modeling for Air Quality Applications. Air Pollution Modeling and its Application XXI, NATO Science for Peace and Security Series C: Environmental Security, Springer, Dordrecht, 57–61,, 2011. 

Winker, D. M. and Vaughan, M. A.: Vertical distribution of clouds over Hampton, Virginia observed by lidar under the ECLIPS and FIRE ETO programs, Atmos. Res., 34, 117–133,, 1994. 

Yang, T., Wang, Z., Zhang, W., Gbaguidi, A., Sugimoto, N., Wang, X., Matsui, I., and Sun, Y.: Technical note: Boundary layer height determination from lidar for improving air pollution episode modeling: development of new algorithm and evaluation, Atmos. Chem. Phys., 17, 6215–6225,, 2017. 

Yuval, Levi, Y., Dayan, U., Levy, I., and Broday, D. M.: On the association between characteristics of the atmospheric boundary layer and air pollution concentrations, Atmos. Res., 231, 104675,, 2020. 

Zhang, L., Wang, T., Lv, M., and Zhang, Q.: On the severe haze in Beijing during January 2013: Unraveling the effects of meteorological anomalies with WRF-Chem, Atmos. Environ., 104, 11–21,, 2015. 

Zhang, Y., Zhang, L., Guo, J., Feng, J., Cao, L., Wang, Y., Zhou, Q., Li, L., Li, B., Xu, H., Liu, L., An, N., and Liu, H.: Climatology of cloud-base height from long-term radiosonde measurements in China, Adv. Atmos. Sci., 35, 158–168,, 2018. 

Short summary
Air pollution has an important impact on human health, climatic patterns, and the ecological environment. The complexity of the nocturnal boundary layer (NBL), combined with its strong physio-chemical effect, induces worse polluted episodes. Therefore, we present a new approach named cluster analysis of gradient method (CA-GM) to overcome the multilayer structure and remove the fluctuation of NBL height using raw data resolution.