Determination and climatology of the diurnal cycle of the atmospheric mixing layer height over Beijing 2013–2018: lidar measurements and implications for air pollution

The atmospheric mixing layer height (MLH) determines the space in which pollutants diffuse and is thus conducive to the estimation of the pollutant concentration near the surface. The study evaluates the capability of lidar to describe the evolution of the atmospheric mixing layer and then presents a long-term observed climatology of the MLH diurnal cycle. Detection of the mixing layer heights (MLHL and MLH′L) using the wavelet method based on lidar observations was conducted from January 2013 to December 2018 in the Beijing urban area. The two dataset results are compared with radiosonde as case studies and statistical forms. MLHL shows good performance in calculating the convective layer height in the daytime and the residual layer height at night. While MLH′L has the potential to describe the stable layer height at night, the performance is limited due to the high range gate of lidar. A nearly 6year climatology for the diurnal cycle of the MLH is calculated for convective and stable conditions using the dataset of MLHL from lidar. The daily maximum MLHL characteristics of seasonal change in Beijing indicate that it is low in winter (1.404±0.751 km) and autumn (1.445±0.837 km) and high in spring (1.647±0.754 km) and summer (1.526±0.581 km). A significant phenomenon is found from 2014 to 2018: the magnitude of the diurnal cycle of MLHL increases year by year, with peak values of 1.291±0.646 km, 1.435±0.755 km, 1.577± 0.739 km, 1.597± 0.701 km and 1.629± 0.751 km, respectively. It may partly benefit from the improvement of air quality. As to converting the column optical depth to surface pollution, the calculated PM2.5 using MLHL data from lidar shows better accuracy than that from radiosonde compared with observational PM2.5. Additionally, the accuracy of calculated PM2.5 using MLHL shows a diurnal cycle in the daytime, with the peak at 14:00 LST. The study provides a significant dataset of MLHL based on measurements and could be an effective reference for atmospheric models of surface air pollution calculation and analysis.

Abstract. The atmospheric mixing layer height (MLH) determines the space in which pollutants diffuse and is thus conducive to the estimation of the pollutant concentration near the surface. The study evaluates the capability of lidar to describe the evolution of the atmospheric mixing layer and then presents a long-term observed climatology of the MLH diurnal cycle. Detection of the mixing layer heights (MLH L and MLH L ) using the wavelet method based on lidar observations was conducted from January 2013 to December 2018 in the Beijing urban area. The two dataset results are compared with radiosonde as case studies and statistical forms. MLH L shows good performance in calculating the convective layer height in the daytime and the residual layer height at night. While MLH L has the potential to describe the stable layer height at night, the performance is limited due to the high range gate of lidar. A nearly 6year climatology for the diurnal cycle of the MLH is calculated for convective and stable conditions using the dataset of MLH L from lidar. The daily maximum MLH L characteristics of seasonal change in Beijing indicate that it is low in winter (1.404 ± 0.751 km) and autumn (1.445 ± 0.837 km) and high in spring (1.647±0.754 km) and summer (1.526±0.581 km). A significant phenomenon is found from 2014 to 2018: the magnitude of the diurnal cycle of MLH L increases year by year, with peak values of 1.291±0.646 km, 1.435±0.755 km, 1.577 ± 0.739 km, 1.597 ± 0.701 km and 1.629 ± 0.751 km, respectively. It may partly benefit from the improvement of air quality. As to converting the column optical depth to surface pollution, the calculated PM 2.5 using MLH L data from lidar shows better accuracy than that from radiosonde compared with observational PM 2.5 . Additionally, the accuracy of calculated PM 2.5 using MLH L shows a diurnal cycle in the daytime, with the peak at 14:00 LST. The study provides a significant dataset of MLH L based on measurements and could be an effective reference for atmospheric models of surface air pollution calculation and analysis.

Introduction
The height of the mixing layer (MLH) is a crucial parameter for near-surface air quality forecasts, pollutant dispersion and the quantification of pollutant emissions (Haeffelin et al., 2012;Seibert et al., 2000;Baars et al., 2008;Liu and Liang, 2010;de Bruine et al., 2017). The pollutants discharged in the boundary layer diffuse vertically under the drive of turbulence (Gan et al., 2011;Monks et al., 2009;Guo et al., 2016) and finally become completely mixed over this layer if sufficient time is given (Emeis et al., 2008). The MLH determines the space in which pollutants diffuse and is thus con-ducive to the estimation of the pollutant concentration near the surface, which might be detrimental to the health of humans and ecosystems Collaud Coen et al., 2014;Singh et al., 2016;Mues et al., 2017).
Within the planetary boundary layer (PBL) the height of the mixing layer (ML) is defined as the height up to which vertical dispersion by turbulent mixing of air pollutants takes place due to the thermal structure of the PBL (Seibert et al., 2000;Schafer et al., 2006;Emeis et al., 2007). The MLH depends largely on the synoptic weather situation (Emeis et al., 2008). The MLH can be estimated by the detection of variance in the mechanical turbulence, the temperature enabling convection or the substance content in the low troposphere. Singh et al. (2016) investigated the evolution of the local boundary layer in the central Himalayan region using a radar wind profiler detecting wind components based on the signal-to-noise ratio profile. Collaud Coen et al. (2014) compared MLH measurements of a microwave radiometer from an atmospheric temperature profile with other measurements in the Swiss Plateau. Mues et al. (2017) used a ceilometer to retrieve the MLH based on the aerosol backscatter signal in the Kathmandu Valley. These measurements are based on different atmospheric parameters, different measuring instruments and various analysis algorithms, leading to MLH results obtained by different methods being inconsistent (Collaud Coen et al., 2014).
In order to move from the general definition to practical measurements, it is necessary to separately consider the structure of the convective boundary layer (CBL) and the stable boundary layer (SBL). In the case of fair weather days, the PBL height has a well-defined structure and diurnal cycle (Collaud Coen et al., 2014;see Fig. S1 in the Supplement). The PBL development under strong convection driven mainly by solar heating is called the CBL (Collaud Coen et al., 2014). The nocturnal SBL shows a more complex internal structure, including a stable layer caused by radiative cooling from the ground, and gradually merges into a neutral layer called the residual layer (RL) (Stull, 1988;Mahrt et al., 1998;Salmond and McKendry, 2005;Collaud Coen et al., 2014). The RL height is the top of the neutral layer and the beginning of the stable free troposphere. The pollutants discharged from the surface at night are restricted to the SBL, while the pollutant emissions from past days tend to stay in the RL. In addition to the dominance of the CBL in the afternoon, the SBL and neutral boundary layer may be formed under certain weather conditions (Stull, 1988;Poulos et al., 2002;Medeiros et al., 2005;Zhang et al., 2018) Since most atmospheric column aerosol particles are usually present in atmosphere below the MLH, the MLH can be used to convert the aerosol optical thickness of the column observed by sun photometers and satellites to the concentration of near-surface pollutants (Sifakis et al., 1998;Dandou et al., 2002;Emeis et al., 2007). Particulate can be used as an important indicator of atmospheric layering because its vertical distribution is strongly affected by the thermal struc-ture of the atmosphere (Neff and Coulter, 1986). Provided the vertical aerosol distribution adapts rapidly to the variational thermal dynamics of the boundary layer, the MLH can be retrieved from the analysis of this aerosol distribution (Emeis et al., 2008).
By the measurement of the profile of aerosol, lidar offers a direct and continuous way to monitor the diurnal cycle of the different layers constituting the PBL (Seibert et al., 2000;De Haij et al., 2006;Emeis et al., 2008;Liu and Liang, 2010;Tang et al., 2016;Su et al., 2017Su et al., , 2019. With recent upgrades of hardware, a ceilometer, also known as an automated lowpower lidar or an automated lidar-ceilometer (ALC), has been demonstrated to capably determine the MLH (Wiegner et al., 2006;Martucci et al., 2007Martucci et al., , 2010Geiß et al., 2017;Kotthaus and Grimmond, 2018;Mues et al., 2017). Recent studies have compared remote sensing measurements (lidar, radar wind profiler, microwave radiometer) with radiosonde (RS) (Wiegner et al., 2006;Milroy et al., 2012;Sawyer and Li, 2013;Cimini et al., 2013;Tang et al., 2016;Singh et al., 2016;Mues et al., 2017;Su et al., 2019), with which convection weather cases have good correlation with differences of 100-300 m, while non-convective weather conditions lead to much larger differences in the MLH estimations if the approaches are intended to measured different structures of the ML such as the CBL, SBL or RL (Collaud Coen et al., 2014). Meteorological radiosondes usually acquire the MLH in the morning (08:00 LT) and at night (20:00 LT), when the diurnal cycle of the ML combined with the stable and convective PBL cannot be well characterized.
In existing studies, numerical simulations, ground-based remote sensing or meteorological radiosondes are used to obtain the characterization of the MLH during short time periods in Beijing, mainly focusing on heavy pollution events (Yang et al., 2005;Quan et al., 2013;Hu et al., 2014;, which underscores the scarceness of continuous high-resolution measurements for a long time period. Depending on the measured atmospheric parameters and observational uncertainties, different measurement approaches may reveal different aspects of PBL structure (Seibert et al., 2000;Seidel et al., 2010;Beyrich and Leps, 2012). Thus, it is of great significance to apply consistent algorithms to consistent types of atmospheric structure parameters when comparing the MLH from different times.
The main aim of this study, therefore, is to present a longterm observed climatology of the MLH diurnal cycle based on lidar observations. For that, the capability of lidar to describe the diurnal evolution of the mixing layer height is evaluated first. The data and methods used are described in Sect. 2. Section 3 contains the results and discussion, which consist of the comparison of lidar-derived MLH with radiosonde measurements, the climatology of the MLH in Beijing, and implications for surface pollution retrieval. Then, it is concluded in Sect. 4. operates at a wavelength of 532 nm with the vertical resolution of 15 m and can detect a long-range profile up to 30 km every second. For the enhancement of the signal-to-noise ratio, 60 profiles are averaged to restore as one, with the time resolution of 1 min. The signal received from the lidar is processed by the subtraction of atmospheric background using the averaged value of the signal received from the height between 22 and 30 km. There are remnants of the previous signal in the system in the absence of optical signal reception, and these signals are called after-pulses. They are monitored every week and removed from the receive signal. Due to the design of the lidar, the received view close to the ground does not completely coincide with the transmitted view. There is a detection blind area in the lidar, and a geometric overlap factor is used to correct the mismatch of the field of view. Then, the correction of the range (range-corrected signal, RCS) is used to retrieve the MLH. The profile of the RCS is expressed as f (z), with z the measurement height. Logarithm calculation of the RCS (expressed in In(β 532 )) is presented in the lidar image (Campbell et al., 2003;Yan et al., 2014;Su et al., 2019). MLH estimation from lidar systems is based on the measurement of the sudden drop in aerosol backscatter at the top of the mixing layer (Seibert et al., 2000). The period of lidar measurements is from 2013 to 2018, nearly 6 years. Except for the lidar data for 2013 mainly in winter and spring (months 1-4 and 11-12), the measurements for 2014-2018 are all annual continued observations.

MLH derived from lidar
Wavelet transforms are commonly used in many studies for MLH determination from lidar observations (Cohn and Angevine, 2000;Davis et al., 2000;Brooks, 2003;De Haij et al., 2006;Baars et al., 2008;Su et al., 2019). When the maximum value of the attenuated backscattering profile convolved with the Haar function is reached, the corresponding height is the MLH. The equation of the wavelet is defined as follows: where b is the transformation of the equation, the equation is cantered and a is the expansion of the equation. The equation of the wavelet covariance transformation W f (a, b), namely the convolution of f (z) with the wavelet function, is defined as follows: where f (z) represents the RCS at different heights, z b denotes the lower limit of the height of the profile and z t represents the upper limit of the height. A valid MLH L is detected, corresponding to the value b when W f (a, b) reaches the biggest local maximum with a coherent scale of a (Brooks, 2003;De Haij et al., 2006;Emeis et al., 2008). In this study, the expansion of a is selected as 420, 435 and 450 m, and the final W f (a, b) is calculated from the averaged corresponding values. Another layer, called MLH L , is detected simultaneously by the first local maximum W f (a, b) from z b , which is assumed to be smaller than or equal to MLH L (De Haij et al., 2006;Mues et al., 2017). Actually, every local maximum corresponds to an aerosol layer and several internal layers appear, making the allocation of a local maximum to an atmospheric feature very difficult (Morille et al., 2007;Geiß et al., 2017;Poltera et al., 2017;Kotthaus and Grimmond, 2018).
Since the absolute maximum in the vertical gradient of the lidar profiles is characterized by a rapid decrease in pollutant concentrations, the MLH L can be associated with the CBL height during daytime (Haeffelin et al., 2012;Poltera et al., 2017) and the RL height during nighttime (Collaud Coen et al., 2014). However, the interpretation of the first local maximum (MLH L ) is critical.
To form a diurnal cycle of the MLH from these several layers, a geodesic approach was applied to pathfinderTURB (Poltera et al., 2017), while COBOLT (Geiß et al., 2017) uses a time-height tracking approach with moving windows. Nevertheless, these methods are all based on the selection of the lowest detected aerosol layer. The height of the lowest detected aerosol layer was regarded as the daytime MLH and the nocturnal stable boundary layer, as reported by Mues et al. (2017) and Kotthaus and Grimmond (2018), respectively. Su et al. (2019) developed a depend time and depend stability (DTDS) algorithm, started with the lowest point and tracked depending on the time and stability, but the nocturnal MLH with SBL height is not evaluated. The detection of nocturnal boundary layer heights, in contrast to the residual layer, is a major challenge (Haeffelin et al., 2012;Lotteraner and Piringer, 2016;de Bruine et al., 2017). Thus, one of the objectives of this study is to investigate the usefulness of MLH L from CE-370 to capture the SBL height over Beijing.
MLH retrievals are eliminated if a cloud flag is marked when the cloud base is found within 6 km of the surface, and a threshold is selected to distinguish between clouds and aerosol layers. To improve the retrieval, a Gaussian filter is applied to retrievals to smooth the temporal variability, and unrealistic outliers are deleted. Due to the limitation of the algorithm and insufficient lidar overlap, the minimum range of the MLH calculation from CE-370 is on the order of 250 m, which is higher than the order of 50 m with a ceilometer. This is due to the optical design of the ceilometer using the same lens for the emitter and the receiver optical paths, which suffers from a low signal-to-noise ratio and provides a lower overlap with the limited power transmitted from the optical design. Detecting significant vertical gradients of attenuated backscatter can be challenging on a clear day when aerosol content is low (Eresmaa et al., 2012;Haeffelin et al., 2012). Compared to CE-370, a ceilometer usually needs a large-scale temporal and vertical average at the cost of a reduction of retrieval in relatively clean atmospheric conditions (de Bruine et al., 2017;Kotthaus and Grimmond, 2018).

MLH from radiosonde
Radiosonde (RS) measurements are one of most widely used methods, especially in China, to derive the SBL height and CBL height due to the ability to characterize the thermodynamic and dynamic states of the boundary layer (Piringer et al., 2007;Seidel et al., 2010;Guo et al., 2019). Meteorological radiosondes perform measurements at the international standard weather station (39.484 • N, 116.282 • E), which was located nearly 11 km from the lidar station. It includes two categories: conventional observations around the year, which are performed at 00:00 UTC (08:00 LST) in the morning and at 12:00 UTC (20:00 LST) in the evening each day, and intensified observations only in summer, which are taken at 06:00 UTC (14:00 LT) in the afternoon. The observed meteorological parameters include atmospheric pressure (P ), temperature (T ), relative humidity (RH), wind speed (WS) and wind direction (WD).
The bulk Richardson number (Ri b ) is a dimensionless parameter combining the thermal energy and the vertical wind shear; it is widely used in MLH climatology (Seidel et al., 2012;Collaud Coen et al., 2014). Ri b is defined as the ratio of turbulence associated with buoyancy to that induced by mechanical shear, which is expressed as where z is the height (z > z s , the subscript "s" denotes the surface), θ characterizes the virtual potential temperature, U and V indicate the two horizontal wind velocity components, and g represents the Earth's gravitational constant. The MLH RS corresponds to the first elevation z, with Ri b greater than a critical threshold taken as 0.25 (Stull, 1988;Seidel et al., 2012;Guo et al., 2016Guo et al., , 2019. In most cases, the exact threshold value has only a small impact on the PBL height due to the large slope of Ri b in this interval (Collaud Coen et al., 2014).

Air pollution model
Data on the MLH are usually combined within an atmospheric model to obtain the surface air pollutant concentration. For example, PM 2.5 remote sensing (PMRS) models, such as that derived by , have the ability to calculate the mass concentration of PM 2.5 above the ground. The PMRS method is designed to employ currently available remote sensing parameters, including aerosol optical depth (AOD), fine-mode fraction (FMF), planetary boundary layer height (PBL height) and atmospheric relative humidity (RH), to derive PM 2.5 from instantaneous remote sensing measurements under different pollution levels Li et al., 2016;Yan et al., 2017). PM 2.5 is calculated following the PMRS model as where AOD indicates aerosol optical depth and FMF represents the fine-mode fraction; VE f is the ratio of the volume and extinction of fine-mode aerosol, which can be calculated from FMF as VE f (FMF) =0.2887FMF 2 -0.4663FMF+0.356. The parameter ρ 2.5,dry indicates the density of dry PM 2.5 , while f 0 (RH) represents the particle hydroscopic growth function, which is f 0 (RH) = (1-RH/100) −1 . The PBL height can be derived from remote sensing and radiosonde measurements. All the parameters are observed by the instruments employed in the same observatory as the lidar. The optical parameters of the column aerosols (AOD and FMF) are obtained by a sky-sun photometer (CE318-DP, CIMEL, France), which is affiliated with the Aerosol RObotic NETwork (AERONET) (Holben et al., 1998;Dubovik et al., 2000). Measurements are automatically scheduled with direct sun irradiance measurements of about 15 min each and angular sky radiance scanning of about 1 h each Che et al., 2014;Wang et al., 2019). Atmospheric meteorological data (relative humidity, wind speed, wind direction, etc.) are obtained by an automatic meteorological monitoring station (BLJW-4). The PM 2.5 mass concentration is obtained by a PM 2.5 monitor (BAM-1020; MetOne, USA), which shows good agreement with the measurements of the national monitoring network near the observatory. All the data are quality-controlled and calculated as 1 h averages, and the measurement period is from 2014 to 2018. The MLH obtained from both lidar and radiosonde within the period is used in the model to calculate the surface PM 2.5 . For convenient comparison with air quality and meteorological parameters, all MLH results are 1 h averaged. 3 Results and discussion

MLH operational measurement
A selection of typical atmospheric conditions included in the dataset of lidar measurements is plotted in Figs. 1-4. The heights of the mixing layer, MLH L and MLH L , are obtained from different criteria using the wavelet covariance transform method. As shown in Fig. 1, the development of a convective mixing layer could clearly be observed, with a sharp decrease in aerosol backscatter between the mixing layer and the free atmosphere. MLH RS is also presented, accompanied by the evolution of MLH L and MLH L .
As shown in Fig. 1a , MLH L and MLH L increase at sunrise. On 2 March 2017, MLH L and MLH L show an obvious diurnal cycle, with a maximum up to 1.0 km. During the evening of 3 March 2017 and the early morning of 4 March 2017, the aerosol layers present two visually obvious layers; MLH L characterizes the first layer height and MLH L represents the upper-layer top. On the next day, due to the existence of cloud, the MLH results are discrete. In the evening and early morning, MLH L deviated from MLH L and approached MLH RS . As shown as the vertical profile from the lidar (RCS, wavelet) and radiosonde (RH, T ) at 20:00 (LST) on 3 March 2017 (Fig. 1b), both MLH L and MLH RS demonstrate 0.53 km, while MLH L shows 1.22 km. Figure 1c indicates that MLH L (0.86 km) approaches MLH RS (0.61 km), although it is a little higher (0.25 km) and much lower than MLH L (1.62 km). In these cases, the result of MLH L , indicated by the first local maximal aerosol gradient, agrees well with MLH RS . However, this is not always true, as shown in Fig. 1d. MLH L (0.752 km) is much higher than MLH RS (0.243 km) but equal to MLH L . It is related to the fact that the stable layer height obtained from radiosonde in this case is out of the range of lidar detection (0.255 km), so MLH L from lidar is disabled to determine the stable layer height. In the summertime (JJA) when the radiosonde is additionally launched at 14:00 LST to detect the convective boundary layer, it can provide the comparison between lidar and RS measurements in the afternoon. As shown in Fig. 2a, MLH L undergoes a rapid increase in the morning and reaches the peak in the afternoon, while MLH L grows with a smaller magnitude. On 25 and 26 August 2014, the aerosol load is relatively low, and MLH L reaches the peak at around 3.0 km in the afternoon, while a lower MLH L peak on 27 August 2014 is observed with a high aerosol content. In the afternoon on these three days, MLH L is consistent with MLH RS , while MLH L is frequently under MLH RS . The measurement of RS in the evening and early morning presents very low values at 0.2-0.3 km. The detailed information represented in Fig. 2b shows that MLH L is equal to MLH RS , which reaches up to 2.95 km, while MLH L is only at 1.24 km. Under clear convective conditions on 25 and 26 August 2014, when there are weak vertical gradients in the aerosol load indicated by a weak RCS, lidar can still obtain good MLH results compared to radiosonde, as shown in Fig. 2c and d with MLH L (2.25 and 2.18 km) approaching MLH RS (2.96 and 2.50 km). The slightly lower value of MLH RS is associ-ated with the fact that aerosol within the mixing layer needs some time to adjust to the thermal structure, and there is a delay in reaching the thermodynamic PBL height (Stull, 1988;Collaud Coen et al., 2014).
As shown in Fig. 3a, the peaks of the three days gradually increase, with values of 1.0, 2.0 and 2.5 km. Due to the high temporal variability of the distribution of aerosol, MLH L presents discontinuity on 17 June 2017. MLH L presents good evolution of the mixing layer height on 15 and 17 June 2017 compared to MLH RS . MLH L corresponds to MLH RS for most of the time in the morning and evening from 15 to 17 June 2017. But when the stable layer height from radiosonde is around or below 0.25 km, for example at 08:00 LST on 16 June 2017 with the RS measurement of 0.27 km (Fig. 3b), MLH L misses the height of the SBL but points to the height of the residual layer (0.61 km). When the stable layer height is higher than 0.25 km, MLH L (0.62 km) tends to approach to MLH RS (0.62 km) (Fig. 3d). However, in the afternoon MLH L was closer to MLH RS compared to MLH L (Fig. 3c). As shown in Fig. 4a , MLH L (MLH L ) presents the diurnal cycle with the maximum of 1.2 km, while the next two days stay stable the whole day, and the height of the SBL is missed by MLH L (Fig. 4b-d).

Intercomparison of different MLH approaches
A comparison of MLH estimated by lidar and radiosonde at 08:00 and 20:00 LST, is shown in Fig. 5. The same observation period of nearly 6 years (2013-2018) is considered, during which the data are continuous except for 2013. As shown in the histograms in Fig. 5, the total column is the annual relative frequency, and the different colors indicate the contribution of each season to the total. There is a wide discrepancy between MLH L and MLH RS at 08:00 and 20:00 LST. The frequency of MLH from a radiosonde lower than 0.25 km is nearly 35 %, whereas there are no data for MLH L from lidar due to the limited detection range. This lower value mainly occurs in winter and autumn, when MLH tends to be lower (Tang et al., 2016). Specifically, the rate of MLH L from lidar lower than 0.5 km is nearly 18 % and 12 % at 08:00 and 20:00 LST, respectively, while the corresponding frequency for radiosonde is beyond 75 % and 66 %. The frequency of the larger MLH L value at 20:00 LST is bigger than that of 08:00 LST from both lidar and radiosonde. It is reasonable that the residual layer has not yet collapsed entirely at 20:00 LST, while the CBL has not developed well in the early morning. As for the MLH L , its distribution trend is more similar to MLH RS than MLH L (see Fig. S2), and the correlation between MLH L and MLH RS is a little higher than that between MLH L and MLH RS , in spite of the fact that it is still not good (see Figs. S3 and S4). This indicates that MLH L has the potential to determine the SBL height as radiosonde does.
As for the seasonal variation of both lidar and RS measurements at 08:00 LST, the frequency of the larger MLH L value in summer is minimal, indicating that the summer MLH is lower than in other seasons. As for radiosonde, MLH L lower than 0.25 km is mostly distributed in winter, with a rate of around 15 % for both 08:00 and 20:00 LST, and the frequency decreases rapidly when MLH L becomes higher than 0.25 km.
The poor agreement between MLH from lidar and MLH RS is also reported in the study of Su et al. (2019), in which it is shown that the correlation of the PBL height measurement between the lidar and radiosonde is 0.14 at 06:30 LST. The significant scatter in the morning and evening is associated with the complicated structure of the boundary layer, as indicated by the existence of a stable boundary layer and residual layer (Su et al., 2019;Tang et al., 2016). In this study, regardless of MLH L and MLH L , more than 35 % of the measurements of the SBL height are not within the scope of the lidar detection. Additionally, in the evening and early morning in some cases, a sufficiently clear variety cannot be found in the backscatter profile at the top of the SBL within the previously well-mixed layer (Russell et al., 1974;Seibert et al., 2000).
The comparison of MLH L and MLH RS at 14:00 LST in summer is presented in Fig. 6, with both mainly indicating the CBL height. MLH L shows very good agreement with MLH RS , with a correlation coefficient of 0.692 and a root mean square error (RMSE) of 0.573 km. It is noted that the slope of the linear fitting line is smaller than the 1 : 1 line, indicating that MLH RS tends to be larger than MLH L in the afternoon, which is consistent with the case study. Although the comparison is only for summer, it can be generally concluded that MLH L from lidar in the afternoon characterizes the CBL height with good accuracy. As shown in Fig. S4, the correlation of MLH L and MLH RS at 14:00 LST is 0.330, and the value of MLH L is generally lower than MLH RS , in-dicating that it is overall unsuitable for MLH L to describe the CBL height in the afternoon.
In fact, complete agreement is not expected between MLH L (MLH L ) derived from lidar and MLH RS from radiosonde for several reasons. First, the two systems measure different atmospheric parameters (aerosol for lidar and temperature, humidity, wind for radiosonde) with varying height resolution and accuracy, and these parameters are influenced in different ways by the processes occurring within the PBL (Seibert et al., 2000). Additionally, it is difficult to identify a clear upper boundary of the mixing layer because the measured parameter is actually not a fixed point but rather a transition layer between two atmospheric states (Stull, 1988;Garratt, 1992;Collaud Coen et al., 2014).
A dataset containing nearly 6 years of measurements in Beijing is used for assessment of the overall performance of the wavelet MLH algorithm (MLH L and MLH L ) with respect to the diurnal availability, as shown in Fig. 7. It can be seen that the expected shape representing the growth of a convective mixing layer is observed. Owing to solar heating of the surface, when the convective layer begins to rise due to upward convection in the early morning and the nocturnal residual layer tends to collapse, MLH L from lidar presents the minimal value. After that, MLH L grows contin-  uously and reaches its maximum height around 15:00 LST, with the value of 1.449, similar to results found for Vienna (Lotteraner and Piringer, 2016) and Berlin (Geiß et al., 2017). The shaded areas indicate the temporal variability as calculated from the standard deviation of MLH L (MLH L ). It is on the order of 600 m for MLH L (300 m for MLH L ) during 09:00-15:00 LST, while the other period is on the order of 700 km (400 m for MLH L ). The larger standard deviation is attributed to the variability of the residual layer.
The diurnal cycles derived from MLH L match the RS results well in the afternoon, but they are larger than MLH RS in the early morning and evening. Contrarily, MLH L tends to approach MLH RS in the early morning and evening but stay far from MLH RS in the afternoon. The difference between MLH L (1.405 ± 0.675 km, mean ± standard deviation of the mean) and MLH RS (1.524 ± 0.582 km) at 14:00 LST is around 0.120 km. This is reasonable considering that RS data are acquired at 14:00 LST only in summertime when MLH L is usually larger throughout the year. However, the discrepancy between MLH L (0.912±0.315 km) and MLH RS (1.524 ± 0.582 km) at 1400 is around 0.612 km. Actually, MLH L is nearly 0.46 km larger than MLH L throughout the day, with a bigger standard deviation. As for the measurement at 08:00 LST, the difference between MLH L (1.196±0.710 km) and MLH RS (0.434±0.364 km) is around 0.762 km, which is larger than the difference between MLH L (0.755±0.334 km) and MLH RS (0.434±0.364 km). The discrepancy between lidar and RS measurements at 20:00 LST is similar.
Overall, MLH L can capture the high SBL height in the nocturnal time, when it is larger than 300 m. The stable layer height detected by MLH L in the nighttime is the layer in which ground-emitted atmospheric pollutants are trapped; it contributes to the assessment of the surface pollutant concentration when there are emissions in the nocturnal time using numerical models (Collaud Coen et al., 2014). Due to incomplete optical overlap, in some cases the point derived from MLH L is the residual layer height rather than the low nocturnal SBL height. And in the daytime, MLH L tends to be lower than the CBL height. In the studies of Mues et al. (2017) and Kotthaus and Grimmond (2018), the MLH in the daytime is usually assigned as the lowest layer detected by a ceilometer. Using higher-power lasers (CE-370) with an increasing signal-to-noise ratio (SNR) and a small gradient detected, the attribution of the lowest layer in the daytime may remain open, since the first local maximum gradient (MLH L ) does not always correspond to the biggest local maximum (MLH L ). Our study indicates that MLH L retrieves the consistent RL height during the night following the CBL diurnal maximal. The RL height corresponds to trapped atmospheric constituents discharged some hours before, which can be employed to convert column-mean optical depths into nearsurface air quality information from remote sensing. And the SBL height provided by radiosonde at 08:00 and 20:00 LST can be considered complementary to the lidar approaches.

Seasonal variation
The seasonal mean diurnal cycle of the MLH L from lidar is shown in Fig. 8. An evident seasonal variation in the magnitude of the diurnal cycle is observed. As shown in Table 1 and Fig. 8, the smallest MLH L magnitude is found in winter, with the peak value of 1.404 ± 0.751 km at 15:00 LST, whereas spring demonstrates the maximum magnitude with 1.647±0.754 km at 15:00 LST. The maximum in summer is 1.526 ± 0.581 km, and the maximum in autumn is 1.445 ± 0.837 km. From the all-day average of the four seasons, the averages in spring, summer, autumn and winter are 1.409, 1.261, 1.297 and 1.228 km, respectively. In summer, MLH L acquired by lidar at 14:00 and 15:00 LST is 1.430 and 1.507 km, respectively, while MLH RS at 14:00 LST is 1.524 km. The measurement of MLH L at 15:00 LST is closer to MLH RS at 14:00 LST. This is consistent with the case study in that it takes some time for aerosol to diffuse upward with the drive of thermal turbulence. As for the statistical variation, the values of autumn MLH L vary most at each hour with the bigger standard deviation, indicating great fluctuations in the long measurement period, while the variation of summer MLH L values for most hours is relatively stable. It should be noted that summer exhibits the biggest amplitude of the diurnal variation of MLH L , with the deepest drop (0.93 km) increasing to the peak value of 1.51 km. Tang et al. (2016) indicate that the lower MLH L value for summer nights and early mornings is attributed to effect of the mountain plain wind. When the local mountain breeze from the northeast in the summer night superimposes the surface cooling, leading to an increase in the thickness of the inversion layer, the height of the mixed layer gradually decreases. After sunrise, with the drive of thermal turbulence, the residual layer height observed by lidar is gradually replaced by a convective boundary layer height, with MLH L increasing rapidly; after 12:00 LT, the plain wind from the southwesterly direction gradually dominates. Previous studies have suggested that the seasonal variation in the MLH L may be associated with radiation flux (Stull, 1988;Kamp and McKendry, 2010;Munoz and Undurraga, 2010), which is consistent with our results. The observational data from Tang et al. (2016) indicate that radiation flux in spring is more than that in summer. The relatively low values in autumn and winter are likely related to the low radiation flux.

Interannual variation
Interannual variations of the MLH L diurnal cycle are investigated in Beijing from 2013 to 2018, as shown in Fig. 9. Diurnal variations of MLH L in different years all have the same patterns but with a different magnitude. Clearly, from 2013 to 2018, the values of the diurnal circle of MLH L increase year by year, including both the RL height at night and the CBL height in the daytime. Since the data for 2013 are mainly from winter and spring, the MLH L seems stable, unlike the amplitude in other years. As shown in Fig. 9 and Table 2, from 2014 to 2018, the MLH L all-day maximum values around 15:00 LST grow year by year, with values of 1.291±0.646 km, 1.435±0.755 km, 1.577±0.739 km, 1.597 ± 0.701 km and 1.629 ± 0.751 km, respectively. The all-day averages of MLH L are 1.110, 1.216, 1.352, 1.391 and 1.502 km, respectively, also showing an increasing trend. This indicates that the volume available for the dispersion of pollutants is extending, which is beneficial to the mitigation of surface pollution. As shown in Fig. S5, from 2014 to 2018, the cumulative increase in the mean MLH L for the whole  day was 0.392 km, and the total increase in the maximum is 0.338 km. As for the annual increase in MLH L , the average of the all-day increments in 2016 is the largest (0.136 km), while the average of the all-day increments in 2017 is the smallest (0.039 km). In particular, the interannual variation of MLH L in the period from 10:00 to 15:00 LST is calculated, when the PBL is characterized by an obvious convective boundary layer. From 2014 to 2018, the average CBL height shows a significant increasing trend of 1.075, 1.212, 1.324, 1.351 and 1.533 km, respectively. The total increase in the average CBL height is 0.458 km. As for the annual increase in CBL height, in 2018 it is the largest (0.182 km), while the average increment for the whole day in 2017 is the smallest (0.027 km).
It is found that, based on measurements from 2014 to 2017, MLH L has a strong negative correlation with AOD (R = −0.41) and with relative humidity (R = −0.21), while MLH L presents a positive correlation with wind speed (R = 0.43) and shows no correlation with temperature (Fig. 10). As shown by the study of Wang et al. (2019), from 2014 to 2017 in Beijing, AOD and surface PM 2.5 have a tendency to decrease year by year, while MLH L increases gradually. The reduction of AOD and surface PM 2.5 is revealed to relate to pollution emission control in recent years (Zhang et al., 2019). Compared with 2016, relative humidity increased in 2017 and wind speed weakened, which is not beneficial for the development of the MLH but is consistent with the small Figure 10. The correlation between MLH L and AOD, relative humidity, wind speed and temperature for measurements from 2014 to 2017. increase in MLH L (0.027 km) in 2017. In addition to the effects of meteorological conditions, the increase in MLH L benefits from the improvement of air quality in Beijing in recent years . Due to the scattering and absorbing of aerosol, solar radiation received from the ground decreases. It is thermal buoyancy generated from surface radiation that drives the PBL to develop. Thus, the development of the MLH is suppressed under a high aerosol load. Hence, with the relief of the radiation effect by aerosol during these years, the turbulence increases, thus leading to a larger PBL height (Ding et al., 2016;Li et al., 2017;. The MLH affects the concentration of pollutants near the surface, while the total radiation of aerosol within the column atmosphere in turn influences the MLH.

Implications for surface pollution retrieval
The vertical structure of the ML is important for pollution concentrations at the surface due to its impact on the volume into which pollutants are mixed . Mues et al. (2017) reported that black carbon concentrations show a clear anticorrelation with MLH measurements. Hu et al. (2014) found a negative correlation between near-surface O 3 and MLH for seven cities in the North China Plain. In the study, as shown in Fig. 11, the correlation between MLH L and observed PM 2.5 data from the same observatory shows a high negative correlation (R = −0.569) with 4 years of measurement (2014)(2015)(2016)(2017). Actually, the pollutant concentration near the surface is affected by the overall effect of local emissions and meteorological conditions, with variation in different spatiotemporal distributions. The MLH is just one of these influencing factors. Geiß et al. (2017) indicated that when the MLH and near-surface concentrations are linked, it is necessary to take the locations and the details, i.e., meteorological conditions and local sources, of the MLH retrieval into account. In fact, all the data used in our study are from the same observatory. And the PMRS model used to calculate the surface PM 2.5 concentration includes the parameters for emissions (AOD) and meteorological conditions (RH).
Due to the difference in the sources of the MLH, lidar and radiosonde, the comparison of derived PM 2.5 _lidar and PM 2.5 _RS with in situ observational PM 2.5 data at 08:00 LST is presented in Fig. 12a and b. MLH from lidar shows reasonably good performance for the retrieval of PM 2.5 in the morning, with a correlation coefficient of 0.741 and RMSE 46.69 µg m −3 . However, the calculated PM 2.5 from MLH RS obviously overestimates the surface pollution, with lower correlation coefficients and a larger standard deviation. The large overestimation should be attributed to the underestimation of the aerosol layer height. In the morning when the PBL is not well developed, above MLH RS there is still a large amount of aerosol, as seen in the lidar images in Figs. 1-4. The discrepancy makes sense given the method using the observed total amount of pollutant in the column atmosphere, including emissions from the surface and the residual aerosol from the day before. Therefore, MLH L from lidar, as a good indicator of the aerosol layer height, is more suitable for estimating surface air pollution from column-mean optical depths.
As presented in Fig. 12c, the calculated PM 2.5 _lidar data for the daytime period (08:00-17:00 LST) show higher correlation (0.846) than those for only the early morning, and the slightly larger RMSE (55.58 µg m −3 ) is associated with the larger number of samples in the statistics. Considering the uncertainty of the series of parameters used in the model, the agreement between calculated PM 2.5 _lidar and in situ measurements is reasonably good. Actually, the accuracy also shows a diurnal cycle, with the peak of correlation coeffi- cients (0.927) at 14:00 LST (Fig. 12d). The correlations at 12:00, 13:00, 14:00 and 15:00 LST were 0.894, 0.922, 0.927 and 0.900, respectively. The higher accuracy may be due to the completed mixing of the aerosol at noon and the vertical distribution of the aerosol tending to be uniform. The correlation between 08:00, 09:00 and 17:00 LST is less than 0.8, which is related to the complex boundary layer structure in the morning and at nightfall. It is difficult to achieve full mixing of the aerosol in the stable boundary layer or the residual layer. The smaller RMSE is related to the limited number of samples. Therefore, the daily variation of the accuracy in the calculated surface pollutant using MLH retrieval by lidar varies with the daily variation of aerosol mixing uniformity at different times during the daytime. Based on the observational data, PM 2.5 tends to peak in the morning and evening. In contrast, afternoon usually witnessed lower mass concentrations due to the rapid vertical diffusion of aerosols . Thus, MLH L from lidar can offer a significant contribution to retrieving the diurnal circle of surface air pollution.

Summary and conclusions
To acquire high-resolution observations of MLH diurnal variation, a study using lidar was performed from January 2013 to December 2018 in the Beijing urban area. The detection of the MLH based on two wavelet methods (MLH L and MLH L ) applied to lidar observations was conducted. The two data results are compared with radiosonde as case studies and statistical forms. The temporal resolution (two or three measurements per day) of PBL detection by RS is not able to provide the mixing layer height diurnal cycle, regardless of its good precision. The MLH shows good performance for the convective layer height in the daytime and the residual layer height at night. MLH L has the potential to describe the stable layer height at night sometimes, even though the capability is limited due to the highly incomplete overlap with the lidar used in the study. The stable layer height detected by MLH L in the nighttime is the layer in which ground-emitted atmospheric pollutants are trapped; it contributes to the assessment of the surface pollutant concentration when there are emissions in the nocturnal time using numerical models. While the residual layer height corresponds to trapped atmospheric constituents discharged some hours before, it can be employed to convert column-mean optical depths into nearsurface air quality information from remote sensing. And MLH L does not always capture the convective layer height as MLH L in the afternoon. Nevertheless, MLH L could be useful as complementary information for the stable layer height in datasets of MLH L .
A nearly 6-year climatology for the MLH L diurnal cycle is calculated for convective and stable conditions. It is true that the height of the mixing layer obtained by different approaches may be different. We focus on the temporal change in the aerosol layer height with a consistent method using a dataset of MLH L . The maximum MLH L characteristics of seasonal change in Beijing indicate that it is low in winter (1.404 ± 0.751 km) and autumn (1.445 ± 0.837 km) and high in spring (1.647±0.754 km) and summer (1.526±0.581 km). A significant phenomenon is found from 2014 to 2018: the magnitude of the diurnal cycle of MLH L increases year by year. The cumulative increase in the mean MLH L for the whole day is 0.392 km, and the total increase in the maximum is 0.338 km. It may partly benefit from the improvement of air quality. As for converting the column optical depth to surface pollution, the calculated PM 2.5 using MLH L data from lidar shows better accuracy than that from radiosonde compared with observational PM 2.5 . Additionally, the accuracy of calculated PM 2.5 using MLH L shows a diurnal cycle in the daytime, with the peak at 14:00 LST. For the operational measurement of PBL height, the MLH from lidar has the capability to mark the diurnal circle of the mixing layer height and can be used as an effective parameter for the vertical distribution of aerosols, providing an important reference to obtain near-ground pollutant concentrations with remote sensing.
Actually, interpreting data from aerosol lidar is often not straightforward because the detected aerosol layers are not always the result of ongoing vertical mixing; they may originate from advective transport or past accumulation processes (Russell et al., 1974;Coulter, 1979;Baxter, 1991;Batchvarova et al., 1999;Tang et al., 2016). Each detection method has good performances only for defined ML structures and under specific meteorological conditions. Therefore, the combination of several methods and instruments may contribute to characterizing the complete diurnal cycle of the complex ML structure (Wiegner et al., 2006;de Bruine et al., 2017;Morille et al., 2007;Kotthaus and Grimmond, 2018).
Data availability. The lidar data used in this study can be requested from the corresponding author (lizq@radi.ac.cn).
Author contributions. ZL conceived and designed the study. JG collected and processed the radiosonde data. YL, HX and PG collected the lidar data. HW improved the MLH retrieval algorithm from lidar. HW processed and analyzed the data and prepared the paper.
Competing interests. The authors declare that they have no conflict of interest.
Special issue statement. This article is part of the special issue "Satellite and ground-based remote sensing of aerosol optical, physical, and chemical properties over China". It is not associated with a conference.