A physical constraint-based lidar boundary layer height inversion method and system

Abstract

The application discloses a physical constraint-based laser radar boundary layer height inversion method and system, physical constraints are combined with atmospheric parameters, a decision tree method is used to estimate the current atmospheric boundary layer height, and thus the boundary layer height is accurately obtained. Compared with a traditional gradient method, the application has the characteristics of simple execution, high accuracy and wide applicability, can overcome the influence of complex atmospheric conditions and background noise on accurate searching of the boundary layer height, and effectively improves the precision of laser radar data inversion of the atmospheric boundary layer height.

Description

Technical Field

[0001] This invention belongs to the field of atmospheric sounding technology, specifically relating to a method and system for retrieving boundary layer height from a lidar system based on physical constraints. Background Technology

[0002] The atmospheric boundary layer, the part of the atmosphere closest to the Earth's surface, comprises the near-surface layer, the mixing layer, and the entrainment layer, playing a crucial role in atmospheric science and applications. It provides vital information about pollutant transport and diffusion, contributing to improved air quality and the implementation of effective pollution control measures. Furthermore, the boundary layer height influences the evolution of meteorological phenomena, significantly impacting the accuracy of weather forecasts and meteorological models. Therefore, accurately retrieving the atmospheric boundary layer height is essential for studying weather change and climate models.

[0003] Lidar, capable of providing high temporal and vertical resolution atmospheric observation data, has been widely used for boundary layer height detection. Lidar can invert the atmospheric boundary layer height by analyzing the difference in aerosol concentrations inside and outside the boundary layer. Traditional lidar boundary layer height algorithms include gradient methods, curve fitting methods, and wavelet transform methods. These algorithms all rely on aerosol concentration information reflected in the backscattering data of the lidar system to invert the atmospheric boundary layer height. For example, gradient methods and curve fitting methods invert the atmospheric boundary layer height by finding the maximum gradient position of the curve fitted from the range correction signal and the original lidar data, while wavelet transform determines the boundary layer height by finding the local maxima of the wavelet-transformed signal. All of these methods depend on the vertical distribution information of aerosols. When atmospheric conditions are complex or background noise is strong, multiple aerosol vertical gradients can appear, disrupting the accurate determination of boundary layer height by traditional lidar boundary layer height algorithms. Therefore, a reliable method is needed to overcome the influence of multiple local gradients generated by complex atmospheric conditions and background noise on traditional lidar boundary layer height algorithms. Summary of the Invention

[0004] To address the shortcomings of existing technologies, this invention provides a physical constraint-based lidar boundary layer height inversion method and system. It utilizes a decision tree algorithm to fuse physical constraints and meteorological parameters to find the boundary layer height, thereby improving the problem that the existence of multiple local gradients under complex atmospheric conditions leads to inaccurate boundary layer height inversion by traditional lidar algorithms.

[0005] To achieve the above objectives, this invention provides a physically constrained lidar boundary layer height inversion method, comprising the following steps:

[0006] Step 1: Calculate the distance correction signal data of the lidar;

[0007] Step 2: Establish a physical constraint feature sequence using the heights corresponding to the three smallest local peaks on the distance correction signal gradient profile;

[0008] Step 3: Use the atmospheric boundary layer height derived from the radiosonde data as a reference value sequence, and use the physical constraint feature sequence and meteorological parameters as feature value sequences. Match the reference values ​​and feature values ​​according to time to construct a dataset.

[0009] Step 4: Construct a decision tree using reference values ​​and feature values, and estimate the boundary layer height using the decision tree algorithm.

[0010] Furthermore, in step 1, the distance correction signal is obtained by multiplying both sides of the lidar signal by the square of the distance.

[0011] Furthermore, the method for constructing the physical constraint feature sequence in step 2 is as follows:

[0012] The formula for calculating the local peak value of the lidar range correction signal gradient is:

[0013] GM = peaks(diff(RCS)) (1)

[0014] In the formula, GM represents the local peak of the gradient, and RCS represents the range correction signal of the lidar. The expression represents calculating the gradient of the signal, and peaks() represents finding local peaks.

[0015] To facilitate differentiation, the three smallest local peaks are denoted as GM1, GM2, and GM3, respectively. The dynamic threshold is taken as one-third of the sum of the largest and smallest local peaks of the distance correction signal gradient profile. If the distance correction signal gradient profile below the dynamic threshold has only one local peak, then GM1, GM2, and GM3 are all the smallest local peaks. If the distance correction signal gradient profile below the dynamic threshold has two local peaks, then GM1 and GM2 are the smaller of the two local peaks, and GM3 is the other local peak. If the distance correction signal gradient profile below the dynamic threshold has three or more local peaks, then GM1, GM2, and GM3 are the three smallest local peaks.

[0016] Furthermore, in step 3, the atmospheric boundary layer height derived from the radiosonde data is used as a reference value sequence, and the physical constraint feature sequences GM1, GM2, and GM3, along with meteorological parameters such as temperature, humidity, wind speed, pressure, solar radiation, net radiation, and soil temperature, are used as feature value sequences. Then, the reference values ​​and feature values ​​are correlated over time to establish a dataset.

[0017] Furthermore, step 4 uses ten-fold cross-validation to verify the accuracy of the decision tree algorithm. The dataset is divided into 10 parts, with each part serving as a test dataset in turn, and the other 9 parts serving as the corresponding training datasets. The feature values ​​of each group of data in each test dataset are input into the constructed decision tree to obtain a predicted value. This predicted value is then compared with the reference value for each group of data to obtain the overall prediction accuracy. The decision tree is constructed using the reference value and feature values, starting with randomly selecting data from the training set. Group data and randomly sample from all features Build a decision tree based on the features, and then randomly draw again with replacement. Group data and Each feature is used to construct another decision tree, and this process continues until N decision trees are constructed. Each node in each decision tree, when making a split decision, selects from a randomly sampled feature... Group sample data and The optimal splitting condition is found among these features, and this process is repeated until the number of samples contained in a node is less than the set minimum number of leaf nodes. So far. When predicting the boundary layer height, feature values ​​containing physical constraints and meteorological conditions are passed as input samples to each decision tree. Each decision tree generates a regression prediction value, and the total regression prediction result P is the average of the regression prediction values ​​of all regression decision trees, that is:

[0018] (2)

[0019] In the formula, N is the number of decision trees. This represents the regression prediction value generated by each decision tree, i.e., the predicted boundary layer height.

[0020] The present invention also provides a physical constraint-based lidar boundary layer height inversion system to implement the physical constraint-based lidar boundary layer height inversion method described above.

[0021] Furthermore, it includes a processor and a memory, the memory being used to store program instructions, and the processor being used to call the program instructions in the memory to execute a physical constraint-based lidar boundary layer height inversion method as described above.

[0022] Alternatively, it may include a readable storage medium storing a computer program that, when executed, implements a physically constrained lidar boundary layer height inversion method as described above.

[0023] Compared with the prior art, the present invention has the following advantages:

[0024] This invention comprehensively considers physical constraints and atmospheric parameters, and utilizes a decision tree algorithm to estimate the current atmospheric boundary layer height, thereby improving the accuracy and efficiency of boundary layer height inversion. This invention can overcome the influence of complex atmospheric conditions and background noise, achieving accurate inversion under various complex atmospheric conditions, and can be widely used in environmental protection, weather forecasting, and other related industries. Attached Figure Description

[0025] To more clearly illustrate the technical solutions in this invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of this invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0026] Figure 1 This is a flowchart of the boundary layer height inversion method of the present invention.

[0027] Figure 2 (a)- Figure 2 (d) is a schematic diagram of the construction of physical constraint features in an embodiment of the present invention, wherein Figure 2 (a) is the distance correction signal-distance correction signal gradient image at 18:00 on February 18, 2017. Figure 2 (b) is a temperature-elevation image at 18:00 on February 18, 2017. Figure 2 (c) is the distance correction signal-distance correction signal gradient image at 12:00 on March 3, 2018. Figure 2 (d) is the temperature-elevation image at 12:00 on March 3, 2018.

[0028] Figure 3 This is a comparison chart of the boundary layer height obtained by the method proposed in this invention with traditional gradient methods, decision tree methods, and radiosonde data inversion.

[0029] Figure 4 (a)- Figure 4 (h) is a comparison chart of the inversion accuracy of the method proposed in this invention and the gradient method, wherein... Figure 4 (a)- Figure 4 (d) shows the correlation images between the boundary layer height retrieved by the gradient method and the boundary layer height retrieved from sounding data under all conditions, stable boundary layer conditions, neutral boundary layer conditions, and convective boundary layer conditions. Figure 4 (e)- Figure 4 (h) are the correlation images of the boundary layer height inverted by the method proposed in this invention with the boundary layer height inverted from the sounding data under all conditions, stable boundary layer conditions, neutral boundary layer conditions, and convective boundary layer conditions. Detailed Implementation

[0030] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below in conjunction with the embodiments and accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of this invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without creative effort are within the scope of protection of this invention.

[0031] Example 1

[0032] like Figure 1 As shown, this embodiment of the invention provides a method for inverting the boundary layer height of a lidar system based on physical constraints, comprising the following steps:

[0033] Step 1: Calculate the distance correction signal data of the lidar.

[0034] The distance correction signal is obtained by multiplying both sides of the lidar signal by the square of the distance.

[0035] Step 2: Establish a physical constraint feature sequence using the heights corresponding to the three smallest local peaks on the distance correction signal gradient profile.

[0036] The formula for calculating the local peak value of the lidar range correction signal gradient is:

[0037] GM = peaks(diff(RCS)) (1)

[0038] In the formula, GM represents the local peak of the gradient, and RCS represents the range correction signal of the lidar. The expression represents calculating the gradient of the signal, and peaks() represents finding local peaks.

[0039] For ease of differentiation, the three minimum local peaks are denoted as GM1, GM2, and GM3, respectively. To avoid the influence of noise, one-third of the sum of the maximum and minimum local peaks of the distance correction signal gradient profile is set as the dynamic threshold. If the distance correction signal gradient profile below the dynamic threshold has only one local peak, then GM1, GM2, and GM3 are all minimum local peaks. Figure 2 (a) is the distance correction signal-distance correction signal gradient image at 18:00 on February 18, 2017, where the solid line represents the distance correction signal, the dashed curve represents the distance correction signal gradient, the dashed straight line represents the dynamic threshold, and the circle represents the local peak value of the local distance correction signal gradient. Figure 2 (b) is a temperature-height image at 18:00 on February 18, 2017, where the horizontal line represents the boundary layer height retrieved from the radiosonde data. From Figure 2As can be seen from (a) and 2(b), there is only one valid local peak value less than the threshold at this moment. Therefore, GM1 = GM2 = GM3 = the only valid local peak value. The gradient method defines the height corresponding to the location of the signal gradient peak value as the boundary layer height. At this time, the heights corresponding to GM1, GM2, and GM3 are all the current boundary layer heights, which are consistent with the boundary layer heights retrieved from the radiosonde data. If the number of local peak values ​​of the distance-corrected signal gradient profile below the dynamic threshold is two, then GM1 and GM2 are the smaller of the two local peak values, and GM3 is the other local peak value. If the number of local peak values ​​of the distance-corrected signal gradient profile below the dynamic threshold is greater than or equal to three, then GM1, GM2, and GM3 are the three smallest local peak values, respectively.

[0040] Figure 2 (c) is the distance correction signal-distance correction signal gradient image at 12:00 on March 3, 2018, where the solid line represents the distance correction signal, the dashed curve represents the distance correction signal gradient, the dashed straight line represents the dynamic threshold, and the circle represents the local peak value of the local distance correction signal gradient. Figure 2 (d) is the temperature-height image at 12:00 on March 3, 2018, where the horizontal line represents the boundary layer height retrieved from the radiosonde data. From Figure 2 As seen in (c) and 2(d), the local peak at approximately 0.8 km corresponds to the boundary layer height retrieved from the radiosonde data. However, the gradient method assigns the height corresponding to the smallest local peak at approximately 1.5 km as the current boundary layer height, leading to a significant inversion error. This may be due to the presence of a residual layer, causing the gradient method to define the residual layer height as the current boundary layer height, resulting in an overestimation. This indicates that the presence of multiple aerosol layers caused by complex atmospheric conditions and background noise affects the gradient method's search for the boundary layer height. Since the distance correction signal gradient has three local peaks, the heights corresponding to the three smallest local peaks can be defined as GM1, GM2, and GM3, establishing a physical constraint feature sequence.

[0041] Step 3: Build a dataset, including a reference value sequence and a feature value sequence, and match the two according to time.

[0042] The atmospheric boundary layer height derived from the radiosonde data is used as the reference value sequence R. The heights corresponding to GM1, GM2, and GM3, along with meteorological parameters such as temperature, humidity, wind speed, pressure, solar radiation, net radiation, and soil temperature, are used as the feature value sequence F. Then, the reference values ​​and feature values ​​are correlated over time to establish a dataset. The radiosonde data used in this embodiment include vertical temperature, pressure, and wind speed. The atmospheric boundary layer height retrieval method used is the Liu and Liang method (Liu S, Liang X Z. Observed Diurnal Cycle Climatology of Planetary Boundary Layer Height[J]. Journal of Climate, 2010, 23(21):5790-5809.).

[0043] Step 4: Use the decision tree algorithm to estimate the boundary layer height.

[0044] Decision trees are constructed using reference values ​​R and feature values ​​F. Ten-fold cross-validation is used to divide the dataset into 10 parts, with each part serving as a test dataset in turn, and the other nine parts serving as the corresponding training datasets. This example assumes the dataset contains 1000 data sets, i.e., 1000 sets of (R, F). When using 1-100 as the test set and 101-1000 as the training set, 25 data sets are randomly selected from these 101-1000 sets, and 6 features are randomly selected from all features to construct a decision tree. Then, with replacement, 25 more data sets and 6 more features are randomly selected to construct another decision tree, and so on, until 300 decision trees are constructed. When 101-200, 201-300, ..., 901-1000 are the training sets, the same method is used to construct 300 decision trees. Each node in each decision tree, when making a split decision, finds the optimal splitting condition from 25 randomly sampled data sets and 6 features. This process is repeated until the number of samples contained in a node is less than the set minimum number of leaf nodes (5). When predicting the boundary layer height, feature values ​​containing physical constraints and meteorological conditions are passed as input samples to each decision tree. Each decision tree generates a regression prediction, and the total regression prediction result is the average of the regression prediction values ​​of 300 trees, i.e.:

[0045] (2)

[0046] In the formula, P represents the overall regression prediction result, which is the final boundary layer height.

[0047] When using 1-100 as the test set and 101-1000 as the training set, 300 decision trees are constructed using the 101-1000 data sets. Each decision tree generates a regression prediction. By sequentially inputting the F values ​​from the 1-100 data sets, 100 predicted values, i.e., 100 boundary layer heights, can be obtained. Similarly, when using 101-200, 201-300, ..., 901-1000 as the test set, 100 predicted values ​​will also be obtained for the corresponding groups. Ultimately, 1000 predicted values ​​will be obtained, meaning each F value will have a predicted value. These values ​​are then compared with the R values ​​from the 1000 data sets to obtain the overall prediction accuracy.

[0048] To demonstrate the effectiveness of this invention in retrieving atmospheric boundary layer height under weak aerosol mixing conditions, the atmospheric boundary layer height on May 13, 2017, was calculated, and the boundary layer height retrieved by the method proposed in this invention was compared with the boundary layer height retrieved by the traditional gradient method and radiosonde data. Figure 3 In the diagram, the inverted triangle represents the boundary layer height predicted by the method proposed in this invention, the equilateral triangle represents the boundary layer height retrieved by the traditional gradient method, and the pentagram represents the boundary layer height retrieved from sounding data. Figure 3 As can be seen, the boundary layer height predicted by the method proposed in this invention is consistently consistent with the boundary layer height retrieved from the radiosonde data. At 06:00, the boundary layer height retrieved by the gradient method is much higher than that retrieved from the radiosonde data, while the boundary layer height predicted by the method proposed in this invention is more consistent with the standard value of the radiosonde data. This may be due to the presence of a residual layer caused by weak turbulence under stable boundary layer conditions, which affects the gradient method's boundary layer retrieval. Furthermore, from... Figure 3 As can be seen, a significant residual layer exists before 10:00, indicating that the traditional gradient method has certain limitations in estimating the atmospheric boundary layer height under weak convection conditions. In contrast, the method proposed in this invention is not limited by atmospheric conditions because it considers physical constraints and meteorological parameters, effectively compensating for the shortcomings of the traditional gradient method.

[0049] To demonstrate that the proposed method can effectively reduce inversion errors, based on 1480-hour average data from a micropulse lidar in a cloudless region from January 1, 2017 to December 31, 2019, and meteorological data, the correlation between the standard boundary layer height values ​​obtained by the gradient method and the proposed method under different atmospheric conditions and with sounding data was compared. Figure 4As shown in (a)-(h), under the overall conditions, the correlation coefficient between the gradient method and the sounding standard value is 0.47, while the correlation coefficient between the method proposed in this invention and the sounding standard value is 0.8. The correlation coefficients of the gradient method under stable boundary layer, neutral boundary layer, and convective boundary layer conditions are 0.20, 0.51, and 0.47, respectively, while the correlation coefficients of the method proposed in this invention are 0.67, 0.78, and 0.81, respectively. This indicates that the method proposed in this invention performs better than the traditional gradient method. This is because the gradient method is easily affected by residual layers, multilayer aerosols, etc., while the method proposed in this invention takes meteorological conditions and physical constraints into account, further improving the stability of atmospheric boundary layer estimation and reducing estimation errors.

[0050] Example 2

[0051] Based on the same inventive concept, the present invention also provides a physical constraint-based lidar boundary layer height inversion system, including a processor and a memory. The memory is used to store program instructions, and the processor is used to call the program instructions in the memory to execute the physical constraint-based lidar boundary layer height inversion method described above.

[0052] Example 3

[0053] Based on the same inventive concept, the present invention also provides a physically constrained lidar boundary layer height inversion system, including a readable storage medium on which a computer program is stored. When the computer program is executed, it implements the physically constrained lidar boundary layer height inversion method described above.

[0054] In specific implementation, the method proposed in the technical solution of this invention can be automatically executed by those skilled in the art using computer software technology. System devices for implementing the method, such as computer-readable storage media storing the corresponding computer program of the technical solution of this invention and computer equipment including the computer program running the corresponding computer program, should also be within the protection scope of this invention.

[0055] The specific embodiments described herein are merely illustrative of the spirit of the invention. Those skilled in the art to which this invention pertains may make various modifications or additions to the described specific embodiments or use similar methods to replace them, without departing from the spirit of the invention or exceeding the scope defined by the appended claims.

Claims

1. A method for inverting the boundary layer height of a lidar system based on physical constraints, characterized in that, Includes the following steps: Step 1: Calculate the distance correction signal data of the lidar; Step 2: Establish a physical constraint feature sequence using the heights corresponding to the three smallest local peaks on the distance correction signal gradient profile; The formula for calculating the local peak value of the lidar range correction signal gradient is: GM = peaks(diff(RCS)) (1) In the formula, GM represents the local peak of the gradient, and RCS represents the range correction signal of the lidar. This indicates calculating the gradient of the signal, and peaks() indicates finding local peaks; The three smallest local peaks are denoted as GM1, GM2, and GM3, respectively. One-third of the sum of the largest and smallest local peaks of the distance correction signal gradient profile is taken as the dynamic threshold. If the distance correction signal gradient profile below the dynamic threshold has only one local peak, then GM1, GM2, and GM3 are all the smallest local peaks. If the distance correction signal gradient profile below the dynamic threshold has two local peaks, then GM1 and GM2 are the smaller of the two local peaks, and GM3 is the other local peak. If the distance correction signal gradient profile below the dynamic threshold has three or more local peaks, then GM1, GM2, and GM3 are the three smallest local peaks. Step 3: Use the atmospheric boundary layer height derived from the radiosonde data as a reference value sequence, and use the physical constraint feature sequence and meteorological parameters as feature value sequences. Match the reference values ​​and feature values ​​according to time to construct a dataset. Step 4: Construct a decision tree using reference values ​​and feature values, and estimate the boundary layer height using the decision tree algorithm.

2. The method for inverting the boundary layer height of a lidar system based on physical constraints as described in claim 1, characterized in that: In step 3, the atmospheric boundary layer height derived from the radiosonde data is used as a reference value sequence, and the physical constraint feature sequence and meteorological parameters such as temperature, humidity, wind speed, pressure, solar radiation, net radiation, and soil temperature are used as feature value sequences. Then, the reference values ​​and feature values ​​are mapped to time to establish a dataset.

3. The method for inverting the boundary layer height of a lidar based on physical constraints as described in claim 1, characterized in that: In step 4, the accuracy of the decision tree algorithm is tested using the ten-fold cross-validation method. The dataset is divided into 10 parts, with each part serving as a test dataset in turn, and the other 9 parts serving as the corresponding training datasets. The feature values ​​of each group of data in each test dataset are input into the constructed decision tree to obtain a predicted value. The predicted value is then compared with the reference value in each group of data to obtain the overall prediction accuracy.

4. The method for inverting the boundary layer height of a lidar based on physical constraints as described in claim 3, characterized in that: In step 4, a decision tree is constructed using reference values ​​and feature values. First, random samples are drawn from the training set data. Group data and randomly sample from all features Build a decision tree based on the features, and then randomly draw again with replacement. Group data and Each feature is used to build another decision tree, and so on, until N decision trees are built.

5. The method for inverting the boundary layer height of a lidar based on physical constraints as described in claim 4, characterized in that: In step 4, each node of each decision tree, when making a split decision, selects from a randomly sampled... Group sample data and The optimal splitting condition is found among these features, and this process is repeated until the number of samples contained in a node is less than the set minimum number of leaf nodes. To date, when predicting the boundary layer height, feature values ​​containing physical constraints and meteorological conditions are passed as input samples to each decision tree. Each decision tree generates a regression prediction value, and the total regression prediction result P is the average of the regression prediction values ​​of all regression decision trees, i.e.: （2） In the formula, N is the number of decision trees. This represents the regression prediction value generated by each decision tree, i.e., the predicted boundary layer height.

6. A lidar boundary layer height inversion system based on physical constraints, characterized in that, It includes a processor and a memory, the memory being used to store program instructions, and the processor being used to call the program instructions in the memory to execute the physical constraint-based lidar boundary layer height inversion method as described in any one of claims 1-5.

7. A lidar boundary layer height inversion system based on physical constraints, characterized in that, The device includes a readable storage medium on which a computer program is stored, and when the computer program is executed, it implements a physical constraint-based lidar boundary layer height inversion method as described in any one of claims 1-5.

Images