Evaluation of ceilometer attenuated backscattering coefficients for aerosol profile measurement

Abstract. Aerosol observations with ceilometers have been made worldwide recently. To use ceilometer data to retrieve aerosol profiles, raw signals should be accurately converted to the attenuated backscattering coefficient. Hence, the calibration coefficient for the system constant has to be determined correctly. We conducted a ceilometer–lidar comparative experiment to evaluate the Lufft CHM15k Nimbus product. The attenuated backscattering coefficient using CHM15k was smaller by a factor of 1.48 compared to that of lidar. The calibration coefficient should be periodically corrected using the ceilometer signal itself since lidar data are generally unavailable in the field observations. We recalibrated the product using both Rayleigh fitting and cloud attenuation methods. The correction factor, determined from the recalibration, was 15% (9%) smaller when using the Rayleigh fitting (cloud attenuation) method than the factor determined from lidar. Uncertainties from backscattering ratios at the reference height and the lidar ratio can cause systematic errors in the correction factor determined from the Rayleigh fitting method. Uncertainties due to the multiple scattering factor contribute to systematic errors for the cloud attenuation method. We propose a calibration method using depolarization ratios for future polarization-sensitive ceilometers, which can estimate the calibration coefficient without multiple scattering factors.


Evaluation of ceilometer attenuated backscattering coefficients for aerosol profile measurement
Retrieval of aerosol backscattering and extinction coefficients from ceilometer data is an important issue for the application to aerosol studies. Aerosol profiling data are useful for validating chemical transport models 8 and improving numerical simulations via data assimilation. 9 Methods to retrieve aerosol properties from ceilometer data have been discussed for several years. Jin et al. 10 demonstrated extinction coefficient retrievals for a dense dust case using the Klett method. 11 Heese et al. 12 retrieved backscattering coefficient from CHM15k signals using the backward inversion method. 13 Since there are two unknown parameters (backscattering and extinction coefficients) in the lidar equation, lidar ratio (extinction-to-backscatter ratio) is assumed in the inversion. The lidar ratio typically varies from 10 to 100 sr depending on aerosol types and therefore the assumption can cause uncertainty in the retrieved properties. In combination with sun photometer data, Wiegner and Geiß 14 estimated lidar ratios for each season in Munich, Germany. Cazorla et al. 15 developed a near real-time processing algorithm to retrieve backscattering coefficients from ceilometer and sun photometer data. Wiegner et al. 16 suggested that the forward inversion method 13 is preferable to retrieve backscattering coefficients for ceilometers because low signal-to-noise ratios (SNR) in the free troposphere do not affect the retrieval results for the lower aerosol layers. In addition, the forward inversion method enables to retrieve backscattering coefficients even under cloudy conditions because atmospheric molecular signals are not necessary for the inversion.
In the forward inversion method, correct initial values are needed to retrieve accurate backscattering coefficient. Because measured signals are used for the boundary condition, the signals have to be calibrated correctly. Namely, calibration coefficient, which is the quantity to calibrate the system constant, should be determined correctly. Jin et al. 10 conducted the signal calibration for a CL51 by comparing it to lidar signals measured simultaneously and revealed that calibration coefficients determined at the ceilometer factory are not necessarily correct. Therefore, ceilometer signals should be recalibrated by users to obtain accurate backscattering coefficient. Comparisons with lidar are a straightforward way to calibrate, but opportunities are limited. Periodical calibration with the ceilometer signal itself is desirable to conduct long-term continuous observations. There are two known calibration methods: the Rayleigh fitting method 14 and the cloud attenuation method. 17 The former method is widely used to calibrate lidars and can be applied to photon counting ceilometers. The latter method utilizes the relationship between the integrated attenuated backscattering coefficient and effective lidar ratio for water clouds. Multiple scattering factor in the effective cloud lidar ratio is an unknown parameter and can cause errors in the calibration process.
The correction of the geometric form factor (or overlap function) is also an important issue. Accurate correction of the overlap function is required to obtain actual aerosol profiles. Incomplete overlap correction can cause large errors in the retrieved backscattering and extinction profiles. The overlap function is determined at the factory and the overlap corrected signals are provided, but the correction is not always accurate. Wiegner et al. 16 conducted horizontal measurements with the CL51 and revealed that the overlap function was overestimated between 60 and 500 m. Hervo et al. 18 determined overlap functions for the CHM15k assuming vertically homogeneous aerosol profiles and pointed out that the overlap function varies diurnally depending on the instrument's internal temperature. Ceilometer overlap functions should be validated using lidars that the overlap function is properly corrected for. Using a small telescope, lidar signals with a low full overlap height are obtained. The signals can be used for the overlap correction for lidar signals with a higher full overlap height.
Ceilometer-lidar comparison is practical to evaluate the ceilometer capability. 10,12,19 We had opportunities to compare CHM15k Nimbus with lidars in August and September 2016 at the National Institute for Environmental Studies (NIES) in Tsukuba, Japan. After the comparative experiment, the ceilometer was transferred to Mandalgovi, Mongolia, in April 2017 to observe dust in the Gobi Desert. In this study, we validate the attenuated backscattering coefficient product of the ceilometer via the comparison. We evaluate the correctness of calibration coefficient and overlap function (Sec. 2). In addition, to evaluate the two known calibration methods, we recalibrate the ceilometer signals using the ceilometer signal itself (Sec. 3). For the cloud attenuation method, we use neutral density (ND) filters to avoid cloud signal saturation caused by the CHM15k photon counting system. Correction factors from the recalibration are validated using the lidar data. Systematic errors in the correction factor are also discussed. In Sec. 4, a new calibration method using the depolarization ratio is proposed for upcoming polarization-sensitive ceilometers. The proposed method can solve the problem of multiple scattering factor in the cloud attenuation method. In addition, a method to increase ceilometer signals near the ground is demonstrated.

Ceilometer-Lidar Comparison
The ability of CHM15k to measure aerosols and clouds is investigated by comparison with NIES lidar of the Asian dust and aerosols lidar observation network (AD-Net). The CHM15k and NIES lidar specifications are summarized in Table 1.
The ceilometer employs a transmitter at 1064-nm wavelength; therefore, it is convenient to compare with aerosol lidars equipped with a Nd:YAG laser. Other ceilometers typically use transmitters at 905 to 910 nm, and the Ångström exponent should be known to compare with the lidars. Moreover, signal attenuation due to the water vapor absorption should be considered for ceilometers with 905-to 910-nm lasers, 20 but the CHM15k has no attenuation from water vapor. The transmitter and receiver are placed biaxially. The full overlap height is around 1.5 km. 12 In this study, attenuated backscattering coefficient provided by CHM15k is express as β 0 obs , where the prime denotes the attenuation.
The AD-Net has about 20 lidar stations in East Asia. AD-Net lidars use a Q-switching Nd: YAG laser with a repetition rate of 10 Hz. The lidars measure backscattering intensity at 1064 and 532 nm and depolarization ratios at 532 nm. Signals at 532 nm are calibrated using the Rayleigh fitting method. The calibration coefficient at 1064 nm is determined using water cloud signals. The backscattering coefficient at 1064 nm is almost similar to that at 532 nm for water clouds since the spectral dependency can be negligible. When the spectral difference for the aerosol optical thickness (AOT) from the ground to the cloud bottom is small, the ratio of the range corrected signals at 1064 nm to the attenuated backscattering coefficient at 532 nm can result in the calibration coefficient at 1064 nm. Depolarization channels are calibrated using signals for AE45deg polarization with respect to the polarization plane of the laser by setting a polarization sheet in front of the telescope. The full overlap height is around 500 to 600 m. The overlap function can be determined using signals measured with a small telescope (50 mm diameter). A detailed description of the AD-Net lidar is described by Shimizu et al. 21 The lidar attenuated backscattering coefficient at 1064 nm is expressed as β 0 lid . Figure 1 shows attenuated backscattering coefficients at 1064 nm observed with CHM15k and NIES lidar. CHM15k can measure the distribution of aerosols and clouds as the lidar. CHM15k signals during the nighttime have low noise because of the photon counting system. Cloud signals at high altitudes are clearly visible during the nighttime compared to the lidar signals. Ghost signals appear above low-level clouds, e.g., during the daytime on 28 August, which may have been caused by the response of the avalanche photodiode to signal saturation from the clouds. The signal saturation is clearly observed when signal profiles are compared. Figure 2 shows the profile comparison on 26 August, when signal saturation appears at 1.5 km. By taking into account a factor of 1.48, the attenuated backscattering coefficient profiles are well matched for aerosols. In Fig. 2, CHM15k signals deviate from the lidar signals below 400 m. The CHM15k signals are beforehand corrected using the overlap function provided by the manufacturer. To identify the realistic aerosol profile, we show an enlarged view of the attenuated backscattering coefficient on 26 August in Fig. 3. The CHM15k product has artificial lines below 400 m, whereas the NIES lidar data do not have such artificial distortion. Consequently, the signal difference between CHM15k and lidar in Fig. 2 is due to the incomplete overlap correction for CHM15k signals. The overlap function provided by the manufacturer should be corrected before the recalibration. In this study, we corrected the function using the lidar data simply by calculating the ratio of the lidar attenuated backscattering coefficient (β 0 lid ) to CHM15k attenuated backscattering coefficient (β 0 obs ). Figure 4 shows the β 0 lid ∕β 0 obs profile averaged over 1 month in August 2016. The ratio is normalized at 1 km. CHM15k signals below 150 m are noisy during  the daytime, as shown in Fig. 3(a), and are not used for the calculation. The ratio is around 1.0 to 1.1 from 300 m to 1 km and significantly increases (>1.2) below 250 m. The β 0 obs value at 180 m results in an ∼40% smaller value than the true value. Hereafter, we use the overlap-corrected β 0 obs for the analyses.

Recalibration Using Ceilometer Signals
Because the magnitude of β 0 obs does not match that of β 0 lid , as shown in the previous section, the CHM15k product should be recalibrated by users. However, lidar data are generally unavailable in the field observation. Therefore, the calibration should be performed with the ceilometer signal itself. We conducted the recalibration for β 0 obs with the Rayleigh fitting and cloud attenuation methods. A correction factor is derived from the recalibration and is defined as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 1 1 6 ; 1 4 7 where β 0 true is the true attenuated backscattering coefficient. We assume that the correction factor of 1.48 AE 0.15 derived from the lidar in August 2016 is the true value, and validate correction factors estimated with the two methods using the true correction factor.

Rayleigh Fitting Method
If backscattering coefficient profiles are retrieved by the backward inversion, β 0 obs can be calibrated referring to the retrieved profiles. Problems in the method are errors of the boundary condition and selection of the lidar ratio. We assume that the aerosol backscatter is negligibly small at the reference height where the inversion calculation is started. Because the Rayleigh scattering intensity is inversely proportional to the wavelength to the fourth power, the intensity at 1064 nm is very small compared to that at the visible wavelength, e.g., the intensity at 1064 nm is 1/16 of that at 532 nm. Consequently, the ceilometer signals at the reference height are usually noisy, and particularly during daytime due to background light. In the case of the analog acquisition method, signals at the reference height are strongly affected by noise with low frequency waveforms. For example, complex waveform noise has been reported for CL31 22 and CL51. 10 Waveform noise makes the calibration difficult due to the large systematic error. Because CHM15k utilizes the photon counting method, no waveform noise is included in the signal; therefore, the Rayleigh fitting method can be applied to the calibration.
We estimated the correction factor for the CHM15k data in August 2016. Clear sky data were selected and averaged for each day. To reduce random noise, signals at the reference height were averaged over AE150 m. The averaged signal was used for the boundary condition of the inversion. When the SNR for the averaged data was <3 at 6 km, the data were excluded from the calculation to reduce random errors. We derived 13 profiles after the extraction. The backscattering ratio R [¼ ðβ 1 þ β 2 Þ∕β 2 , where β 1 and β 2 are the backscattering coefficients for aerosols and atmospheric molecules, respectively] is assumed to be 1 at the reference height. The lidar ratio for the inversion calculation was assumed to 50 sr. The estimated correction factor was 1.26 AE 0.14, which is 15% smaller than the true value.
Systematic errors in the correction factor due to the backscattering ratio at the reference height are investigated. Because Rayleigh scattering at 1064 nm is much smaller than that at 532 nm, the backscattering ratio at 1064 nm can significantly change even if few aerosols exist at the reference height. For example, if the backscattering ratio at 532 nm is 1.05 and the background-type aerosol is assumed, the backscattering ratio at 1064 nm is 1.33. 23 We assumed a backscattering ratio of 1 at the reference height, but the estimated correction factor would contain systematic errors if the true backscattering ratio is different from the set backscattering ratio. To investigate the systematic error, simulated lidar signal profiles were created assuming an aerosol layer (from 0 to 2 km) with two AOTs: 0.1 and 0.5. We calculated the correction factor changing the true backscattering ratio R true from 1 to 1.3 for the simulated lidar signal. If R true − R set is 0.3 and AOT is 0.1, the systematic error in the correction factor is approximately −20% [ Fig. 5(a)]. The systematic errors decrease with increasing AOT.
Systematic errors due to the lidar ratio selection are also investigated. We calculated the error changing the true lidar ratio S 1;true from 20 to 80 sr while the setting lidar ratio S 1;set to the fixed value of 50 sr. The results are shown in Fig. 5(b). If the true lidar ratio is smaller than the setting value, the systematic error will be negative. The systematic errors increase with increasing AOT.
We also investigated systematic errors caused by an incomplete overlap correction. The correction factor can vary depending on the overlap function due to changes in the retrieved extinction coefficient. We calculated the correction factor using the CHM15k product without the additional overlap correction from the lidar data. The calculated correction factor was 2% larger than that after the additional overlap correction.

Cloud Attenuation Method
Platt 24 introduced derivation of effective cloud lidar ratios from the lidar equation. By integrating the true attenuated backscattering coefficient for clouds, the integrated value can be expressed by the effective cloud lidar ratio (ηS c ) as follows: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 1 1 6 ; 1 2 8 where γ 0 true is the integrated attenuated backscattering coefficient calculated by the integral of β 0 true from r b to the cloud top height r t . T 2 denotes two-way atmospheric transmittance from the bottom height r 0 to the cloud base height r b . η is the multiple scattering factor, which can change depending on the height, extinction coefficient, and effective radius of the water clouds. Therefore, η is an unknown parameter in the effective cloud lidar ratio. S c is the cloud lidar ratio, and τ c is the cloud optical thickness. Because the typical τ c value for water clouds is 3 to 5, 25 the exponential term in Eq. (2) can be approximated as 0. Using the ηS c − γ 0 true relationship, O'Connor et al. 17 proposed a method to calibrate ceilometer signals assuming that the S c for water clouds is known. S c for water clouds can be theoretically determined; according to Pinnick et al., 26 S c at 1.06 μm for water clouds is 18.2 sr for a wide range of droplet sizes. Because β 0 true ∕β 0 obs equals to γ 0 true ∕γ 0 obs , where γ 0 obs is calculated from β 0 obs , γ 0 true can be expressed as Cγ 0 obs using Eq. (1). Therefore, we can estimate the correction factor using E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 3 ; 1 1 6 ; 3 3 5 The method is typically applied to the analog system. Because photon counting signals saturate for water clouds, as shown in Fig. 2, applying the method for CHM15k is difficult. To solve the saturation problem, we put an ND filter with 1% transmittance on the receiver lens in September 2016. Signal attenuation by the ND filter was corrected for the observed attenuated backscattering coefficient. On September 20, 2016, optically thick water clouds were observed at 1 and 2 km in height [ Fig. 6(a)]. Figure 6(b) shows the time series of the calculated C. We assumed T 2 ¼ 1 because aerosol signals below clouds are very small. η was also assumed to 1. From 19:00 to 19:30, C values are out of range because the signals are not fully attenuated. C tends to decrease with increasing the cloud height, which may be due to the change in η. The average C from 17:00 to 17:30 is 1.34 AE 0.1, which is 9% smaller than the true value.
Systematic errors in the correction factor are discussed here. The term T 2 can cause systematic errors in the correction factor. 27 For example, when the AOT is 0.1, T 2 is 0.82. In this case, the correction factor will be ∼20% smaller than the true value if the term is ignored. We calculated T 2 for the case in Fig. 6 using the lidar data. A lidar ratio of 50 sr was assumed. The calculated T 2 is 0.99 and therefore the effect of T 2 on the correction factor is negligibly small. The term η, which is another source of the systematic errors, can be roughly estimated from the theoretical calculations 28 but has a large uncertainty in the calibration method. If the true value of the correction factor is 1.48, η results in 0.9 for the case in Fig. 6. The problem due to η can be solved using depolarization ratios as discussed in Sec. 4.2.

Improving Signal Intensity at Low Altitudes
As pointed out in Sec. 2, CHM15k signals below 150 m are noisy since the signals near the ground are very weak. Users are difficult to use the signals below 150 m for aerosol studies. The signals near the ground should be improved to make ceilometers more advanced backscatter lidar. To obtain stronger atmospheric scattering signals at the low altitudes, we conducted an experiment using an additional lens in March 2017. We placed a small planoconvex lens (25 mm diameter and 5 m focal length) on the center of the primary receiver lens. Figure 7 shows signal profiles before and after setting the additional lens. The signals remarkably increase between 50 and 100 m after the additional lens is set.
Because the composite focal length is shorter for two lenses, the receiver lens and additional lens, than the focal length of the primary receiver lens, the composite lens can collect low level atmospheric scattering light more effectively while signals from high altitudes decrease. Above 1 km, the signals collected from the composite lens are theoretically 0.1% of the signals from  the primary receiver lens. Consequently, slope changes for high altitude signals due to the composite lens can be ignored. In addition, the effective area of the composite lens accounts for ∼4% of the area of the primary receiver lens. Therefore, the decrease of signal intensity due to the composite lens has little effect on SNR at high altitudes. In practical data analysis, an overlap function is estimated after the additional lens is set. The proposed method can improve SNR at low altitudes and facilitate using ceilometer signals for monitoring aerosols near the ground. The method is simple for users since they only need to set the additional lens for the current system.

Calibration Method Using Depolarization Ratios
The multiple scattering factor η introduces a large uncertainty in the correction factor (or calibration coefficient) in the cloud attenuation method. We present a new calibration method that excludes the η uncertainty using a depolarization ratio measurement. Multiple scattering by water clouds increases the cloud depolarization ratios. Hu 29 reported that η and depolarization ratios are related via the following equation determined by the Monte Carlo simulation: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 4 ; 1 1 6 ; 5 4 3 where δ is the layer integrated depolarization ratio expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 5 ; 1 1 6 ; 4 8 6 δ ¼ where P ⊥ and P k are the perpendicular and parallel components of lidar signals, respectively. The depolarization channels should be calibrated beforehand by known methods (e.g., the AE45deg calibration method 30 ). Therefore, the multiple scattering factor can be represented using the depolarization ratios. In this study, assuming that spectral dependence in multiple scattering is small for water clouds, we used δ at 532 nm instead of δ at 1064 nm because δ at 1064 nm has not been measured at this time. We calculated ηS c using uncalibrated lidar signals at 1064 nm. To estimate the T 2 term in Eq. (2), we used AOT data measured with a Skyradiometer, assuming that only photometer data are available as supplemental data for ceilometer calibration. Measured AOT in August 2016 was 0.09 AE 0.03, resulting in T 2 ¼ 0.83. The AOT was measured under clear sky conditions, but we assumed that the monthly averaged AOT can be used for this method. Figure 8 shows the frequency of occurrence for δ and ηS c in August 2016. The calibration coefficient was determined by fitting the mode of the frequency distribution of ηS c at δ ¼ 0.1 to 0.15 to the theoretical line plotted in Fig. 8. The calibration coefficient is ∼3.3% smaller compared to the value determined by the procedure described in Sec. 2. The difference may come from the spectral dependency of the depolarization ratio for water clouds, incorrect AOT, and random errors. We investigated the impact of the spectral dependency of the depolarization ratio on the calibration coefficient. Changing the spectral ratio of depolarization ratio (1064/532) from 0.9 to 1.1, the error is −5% (7%) when the ratio is 0.9 (1.1).
In Fig. 8, clouds around region A are well matched with the theoretical curve. Region B deviates from the theoretical curve because the cloud signals are not fully attenuated by optically thin clouds. Region C is smaller than the theoretical curve. Clouds around region C have a positive correlation between δ and ηS c ; therefore, these clouds are considered as horizontally oriented ice, 29 which is difficult to discriminate from water clouds without depolarization ratios. In the analysis to estimate the calibration coefficient, clouds around regions B and C should be excluded. The calibration coefficient is determined by matching clouds around region A to the theoretical curve. The main advantage of this method is that the calibration coefficient can be determined without knowing the multiple scattering factor.

Conclusions
CHM15k is a very useful tool to measure aerosol profiles qualitatively, and the attenuated backscattering coefficient product can be used for aerosol studies if the overlap function and calibration coefficient are corrected by users. This study indicates that the CHM15k attenuated backscattering coefficient was smaller by a factor of 1.48 than the lidar data. The Rayleigh fitting method is desirable for recalibrating CHM15k signals because the photon counting acquisition system facilitates measuring clean air signals at high altitudes during nighttime. However, inappropriate selections of lidar ratio and initial backscattering ratio at the reference height can cause systematic errors in the correction factor. The correction factor determined using the Rayleigh fitting method was 15% smaller than that determined from the lidar data.
Recalibration using the cloud attenuation method is also possible for CHM15k if an ND filter with 1% transmittance is used to avoid signal saturation for optically thick clouds. However, calibration times are limited because the ND filter is required. In addition, the multiple scattering factor can induce a large uncertainty in the correction factor. The correction factor using the cloud attenuation method was 9% smaller than the true value. We proposed a method to estimate the correction factor using the relationship between effective cloud lidar ratios and depolarization ratios. The calibration coefficient estimated by the method was ∼3.3% smaller compared to the true value. The method does not require multiple scattering factors in the estimation. The method cannot be applied to the calibration of the current CHM15k system, but it will be useful for future polarization-sensitive ceilometers.
The overlap function below 300 m is incomplete and should be corrected using lidar data. The signals at 180 m were 40% smaller than the lidar data. Due to the narrow FOV, CHM15k signals below 150 m are too noisy to correct. We demonstrated that the signal intensity can drastically increase using a composite lens method. The configuration of the composite lens will be optimized to improve the SNR below 50 m.