Sand-dust vertical flux high-precision inversion method based on laser radar

Abstract

The invention discloses a sand and dust vertical flux high-precision inversion method based on a laser radar. The method comprises the following steps: acquiring multi-wavelength back scattering and polarization information by using a 355 nm, 532 nm and 1064 nm three-wavelength polarization laser radar system; constructing a five-dimensional optical feature vector space, and combining a support vector machine classifier to realize automatic identification of dust particles; a variational data assimilation technology is adopted to fuse radar observation and numerical forecasting information to invert a three-dimensional wind field; inverting sand and dust mass concentration vertical distribution based on the corrected particle spectrum distribution model; flux calculation and uncertainty quantization are realized through adaptive weighted fusion and a Monte Carlo method; the method has the characteristics of high temporal-spatial resolution, high precision, strong adaptability and the like, and can be widely applied to the fields of weather forecast, environment monitoring, climate research and the like.

Description

Technical Field

[0001] The present invention relates to the field of atmospheric environment monitoring technology, and in particular to a high-precision inversion method for vertical dust flux based on laser radar. Background Art

[0002] Dust aerosols are an important component of the global atmospheric composition. The vertical flux of dust is a key parameter that describes the vertical transport process of dust particles in the atmosphere. It is of great significance for understanding the source, sink, and transport mechanism of dust and its impact on the atmospheric radiation balance.

[0003] At present, the measurement of vertical dust flux mainly relies on traditional ground observation stations. These methods have limitations such as poor spatial representativeness, low temporal resolution, and inability to obtain vertical distribution information. Although satellite remote sensing can provide large-scale dust distribution information, its vertical resolution is limited, and it is difficult to accurately obtain detailed information of the near-surface layer. Although numerical models can simulate the dust transport process, the uncertainty of the model parameterization scheme is large, and high-quality observation data is needed for verification and improvement.

[0004] As an active remote sensing method, lidar technology has the advantages of high temporal and spatial resolution and the ability to obtain vertical profile information, providing a new technical approach for the accurate measurement of vertical dust flux. However, existing lidar inversion methods mostly focus on the measurement of dust concentration, and there are still technical challenges in the quantitative inversion of vertical dust flux, including: limited accuracy in dust particle identification, difficulty in obtaining wind field information, and insufficient uncertainty quantification.

[0005] Therefore, there is an urgent need to develop a high-precision inversion method for dust vertical flux based on lidar to meet the actual needs of atmospheric environment monitoring and climate research. Summary of the Invention

[0006] The purpose of the present invention is to provide a high-precision inversion method for dust vertical flux based on lidar, so as to solve the technical problems of low measurement accuracy, insufficient temporal and spatial resolution, and insufficient uncertainty quantification of dust vertical flux in the existing technology, and to achieve high-precision quantitative monitoring of the vertical transmission process of dust particles in the atmosphere.

[0007] In order to achieve the above object, the present invention adopts the following technical solutions: A high-precision inversion method for vertical dust flux based on lidar obtains atmospheric scattering signals through a multi-wavelength polarization lidar system. After signal preprocessing, particle identification, wind field inversion, mass concentration calculation and flux fusion, high-precision quantitative measurement of the vertical transmission flux of dust particles in the atmosphere is achieved.

[0008] The method specifically comprises the following steps: Step S1: Use a three-wavelength polarization lidar system with working wavelengths of 355nm, 532nm and 1064nm to emit laser pulses into the atmosphere and receive backscattered echo signals from atmospheric aerosol particles; perform distance correction, background noise subtraction and signal smoothing on the original echo signals to obtain backscattering coefficient and depolarization ratio profile data at different wavelengths and polarization states.

[0009] Step S2, based on the multi-wavelength backscattering coefficient ratio, depolarization ratio and signal gradient optical characteristic parameters, construct a multi-dimensional feature vector space to distinguish dust particles from non-dust components in the atmosphere, and realize the positioning and boundary identification of the dust layer.

[0010] In step S3, the three-dimensional wind field distribution of the atmosphere is inverted using variational data assimilation technology, combining the lidar Doppler radial velocity observation data and the background wind field information of the numerical weather forecast. The conjugate gradient method is used to solve the cost function containing background field constraints, observation constraints, and physical constraints to obtain a wind speed vector field with high spatiotemporal resolution and the vertical wind speed component.

[0011] In step S4, assuming that the dust particles obey the modified particle spectrum distribution model, the particle spectrum parameters are inverted using the multi-wavelength extinction coefficient observation data, and the number concentration distribution of the dust particles and the vertical profile of the dust mass concentration are calculated by solving the integral equation of the extinction coefficient and the particle spectrum distribution.

[0012] In step S5, the obtained vertical wind speed component is multiplied by the dust mass concentration to calculate the instantaneous dust vertical flux; error propagation analysis and statistical uncertainty quantification are performed, and 10,000 random samplings are performed using the Monte Carlo method to provide a confidence interval for the flux measurement results.

[0013] Furthermore, the signal smoothing process uses a sliding average algorithm for time smoothing with a window length of 30 seconds, and a spatial smoothing with a vertical resolution of 7.5 meters.

[0014] Furthermore, in step S2, dust particle identification is implemented based on a multi-dimensional optical feature vector and a support vector machine classifier; the multi-dimensional optical feature vector is a five-dimensional feature vector , used to characterize the optical properties of particles, ;in, , represents backscattering at near-infrared wavelengths Backscattering with visible wavelengths The ratio of , representing the backscattering of ultraviolet wavelengths Backscattering ratio with visible wavelength; , represents the linear depolarization ratio at 532 nm wavelength; , represents the vertical gradient of the backscattering coefficient at 532 nm with height z, in km; , represents the normalized signal intensity after 532nm wavelength distance square correction, is the emitted laser power at a wavelength of 532 nm.

[0015] The support vector machine classifier uses a support vector machine with a radial basis function kernel for binary classification, and the decision function is: ;in, is the Gaussian radial basis kernel function.

[0016] is the kernel width parameter, , is the Lagrange multiplier; is the training sample label, +1 represents dust particles, and -1 represents non-dust components; is the bias parameter, is the total number of support vectors.

[0017] Furthermore, the training samples of the support vector machine classifier come from lidar observation data under known sandstorm events and clean atmospheric conditions, and are labeled by combining synchronous PM10 ground observations, satellite remote sensing and backward trajectory analysis. The number of sandstorm samples is no less than 5,000, and the number of non-sandstorm samples is no less than 10,000.

[0018] The dust layer boundary identification adopts a threshold judgment method, and when the dust probability of three consecutive height layers is greater than 0.8, it is determined to be the dust layer boundary.

[0019] Furthermore, in step S3, the variational data assimilation technique is used to invert the three-dimensional atmospheric wind field distribution, and the variational assimilation cost function including the background field term, the observation term and the regularization term is minimized based on the wind field inversion. Expressed as: ; in, To find the three-dimensional wind field, are the east-west, north-south and vertical wind speed components, respectively. Prior wind fields provided for numerical weather forecasts, is the background error covariance matrix; is the radial velocity measured by the lidar, is the observation error covariance matrix; is the regularization parameter, is the Laplace operator; is the lidar observation operator, which describes the geometric projection relationship between the three-dimensional wind field and the radar radial velocity: ; is the azimuth of the laser beam relative to the true north; is the elevation angle of the laser beam relative to the horizontal plane; the cost function is solved by the conjugate gradient method, and the convergence criterion is that the gradient modulus is less than 10 -6 or the number of iterations exceeds 1000.

[0020] Furthermore, in step S4, in the modified particle size distribution model, the dust particle number concentration The distribution with radius follows a modified power-law spectrum: ;in, represents the spectral coefficient; represents the spectral index; represents the exponential cutoff radius; and Indicates the lower and upper limits of the particle radius.

[0021] The relationship between the extinction coefficient of each wavelength and the particle spectrum is: ; The dust mass concentration is obtained by integrating the particle volume: ; in, Indicates wavelength The extinction coefficient at ; Represents the particle radius ,wavelength , complex refractive index Mie scattering extinction efficiency factor under conditions; Indicates the density of dust particle material; represents the mass concentration, and the various inversion spectrum parameters are obtained by minimizing the difference between the observed and theoretical extinction coefficients.

[0022] Furthermore, in the calculation of Mie scattering extinction efficiency, the complex refractive index of dust is Using a wavelength-dependent model: ; Real and imaginary parts use empirical formulas respectively; Real wavelength : ; is the real part of the refractive index in the high frequency limit, 1.53, is the dispersion intensity parameter, is the resonant wavelength; the imaginary wavelength : ; is the absorption coefficient at the reference wavelength; As the reference wavelength, take 550nm; is the absorption index, which is taken as 1.2.

[0023] Furthermore, in step S5, an adaptive weighting method based on measurement uncertainty is used to fuse the inversion results; The inversion methods include: a dual-wavelength inversion method based on 355nm and 532nm, a dual-wavelength inversion method based on 532nm and 1064nm, and a three-wavelength joint inversion method, totaling three independent inversion methods. The weight of each inversion method is determined based on the inverse of the uncertainty of each inversion method: ; Uncertainty of each inversion method Contains three components: ;in, is the instrument system error, including laser power stability and detector noise; Processing errors for the algorithm; is the inversion process error; Optimal fusion flux calculation: ;in, For the The mass concentration inverted by the two methods (μg / m³) is the corresponding vertical wind speed (m / s), is the vertical flux (μg / (m²·s)).

[0024] Furthermore, the output is evaluated by randomly perturbing the input parameters to quantify the measurement uncertainty of the vertical dust flux; the input parameter perturbation is set so that the backscattering coefficient perturbation follows a Gaussian distribution with a standard deviation of 5% of the observed value; the extinction coefficient perturbation follows a Gaussian distribution with a standard deviation of 8% of the observed value; the radial velocity perturbation follows a Gaussian distribution with a standard deviation of 0.3m / s; the real part perturbation range of the complex refractive index is ±0.01, and the imaginary part perturbation range is ±0.001.

[0025] Furthermore, the flux uncertainty calculation method is: ;in, For the The flux value of the Monte Carlo experiment, is the number of Monte Carlo trials, which is 10000; is the mean flux.

[0026] Bootstrap resampling is used to estimate the flux distribution: ;in, and They are distributed by Bootstrap and Quantile, is the significance level, corresponding to the 95% confidence interval; the number of bootstrap resampling is 1000 times.

[0027] Compared with the prior art, the present invention has the following beneficial effects: This invention achieves high-precision real-time monitoring of vertical dust flux through the deep integration of multi-wavelength polarization lidar technology and advanced data processing algorithms. Compared with traditional methods, it has significant improvements in measurement accuracy, spatiotemporal resolution and uncertainty quantification: the relative error of measurement is reduced from more than 50% to less than 30%, the time resolution is increased from hours to minutes, the vertical resolution reaches 7.5 meters, and comprehensive uncertainty quantification of flux measurement is achieved, providing reliable scientific data support for sandstorm warning, air quality assessment, climate change research and numerical model verification, and has important scientific value and broad application prospects. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] Figure 1 This is a flow chart of the high-precision inversion method of dust vertical flux based on lidar of the present invention. DETAILED DESCRIPTION

[0029] To make the objectives, technical solutions, and advantages of the present invention more clear, the technical solutions of the present invention are described clearly and completely below. Obviously, the embodiments described are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making any creative efforts are within the scope of protection of the present invention.

[0030] like Figure 1 The figure shows a flow chart of the high-precision inversion method of dust vertical flux based on lidar of the present invention. The atmospheric scattering signal is obtained by a multi-wavelength polarization lidar system. After signal preprocessing, particle identification, wind field inversion, mass concentration calculation and flux fusion, high-precision quantitative measurement of the vertical transmission flux of dust particles in the atmosphere is achieved.

[0031] The specific steps include: Step S1: Use a three-wavelength polarization lidar system with working wavelengths of 355nm, 532nm and 1064nm to emit laser pulses into the atmosphere and receive backscattered echo signals from atmospheric aerosol particles; perform distance correction, background noise subtraction and signal smoothing on the original echo signals to obtain backscattering coefficient and depolarization ratio profile data at different wavelengths and polarization states.

[0032] The three-wavelength polarization lidar system uses an Nd:YAG laser as its light source, generating 355nm ultraviolet light and 532nm visible light through frequency doubling technology, with 1064nm as the fundamental frequency. The system transmits 50mJ at 355nm, 100mJ at 532nm, and 200mJ at 1064nm, with a repetition rate of 20Hz. The receiving system uses a Cassegrain telescope with a 300mm aperture and a 1200mm focal length, equipped with a polarization beam splitter to enable simultaneous detection of parallel and perpendicular polarization components. The distance correction formula is: ,in, is the atmospheric molecular attenuation coefficient.

[0033] The working principle of the three-wavelength polarization lidar system is based on the differences in scattering characteristics of particles of different sizes due to different wavelengths; the 355nm ultraviolet wavelength is sensitive to fine particles and can effectively detect submicron dust particles; the 532nm visible light wavelength is the standard detection wavelength with good signal-to-noise ratio and stability; the 1064nm near-infrared wavelength is sensitive to coarse particles and helps to identify large-sized dust particles; distance correction compensates for the geometric diffusion loss of the laser during propagation and is corrected according to the inverse square law of distance; background noise subtraction estimates the noise level through the statistical characteristics of long-distance signals and uses an adaptive threshold method to eliminate noise interference; the signal smoothing process uses a sliding average algorithm for time smoothing with a window length of 30s, and spatial smoothing uses a vertical resolution of 7.5m; for example, in dust observations in a certain area, the signal-to-noise ratio of the preprocessed signal can be increased by 3-5 times.

[0034] Step S2, based on the multi-wavelength backscattering coefficient ratio, depolarization ratio and signal gradient optical characteristic parameters, construct a multi-dimensional feature vector space to distinguish dust particles from non-dust components in the atmosphere, and realize the positioning and boundary identification of the dust layer.

[0035] Non-dust components include cloud droplets, water vapor, pollution particles or other atmospheric components. Dust particle identification is based on multi-dimensional optical feature vectors and support vector machine classifiers. The multi-dimensional optical feature vector is a five-dimensional feature vector. , used to characterize the optical properties of particles, ;in, , representing backscattering at near-infrared wavelengths Backscattering with visible wavelengths The ratio reflects the particle size information. The value range of this ratio for dust particles is 0.3-0.8; , representing the backscattering of ultraviolet wavelengths The backscattering ratio of visible light wavelength reflects the characteristics of particle composition. The value range of this ratio for dust particles is 0.1-0.5; , represents the linear depolarization ratio at 532 nm wavelength, characterizing the non-spherical shape of the particle. The value is relatively large (0.2-0.4) due to the irregular shape of dust particles; , represents the vertical gradient of the backscattering coefficient at 532 nm with height z, in km⁻²; , represents the normalized signal intensity after distance square correction at 532nm wavelength, where z is the height from the lidar (km), is the emitted laser power (W) at a wavelength of 532 nm.

[0036] The support vector machine classifier uses a support vector machine with a radial basis function kernel for binary classification, and the decision function is: ;in, is the Gaussian radial basis kernel function, is the kernel width parameter, which controls the complexity of the decision boundary, , is the Lagrange multiplier, which only supports the vector corresponding to ; is the training sample label, +1 represents dust particles, and -1 represents non-dust components; is the bias parameter, is the total number of support vectors.

[0037] The training samples of the support vector machine classifier come from lidar observation data under known sandstorm events and clean atmosphere conditions. Sample annotation is carried out in combination with synchronous PM10 ground observations, satellite remote sensing and backward trajectory analysis. The number of sandstorm samples is no less than 5,000, and the number of non-sandstorm samples is no less than 10,000. The support vector machine training adopts 10-fold cross-validation, and the classification accuracy reaches more than 95%.

[0038] The dust layer boundary identification adopts a threshold judgment method, and when the dust probability of three consecutive height layers is greater than 0.8, it is determined to be the dust layer boundary.

[0039] In step S3, the three-dimensional wind field distribution of the atmosphere is inverted using variational data assimilation technology, combining the lidar Doppler radial velocity observation data and the background wind field information of the numerical weather forecast. The conjugate gradient method is used to solve the cost function containing background field constraints, observation constraints, and physical constraints to obtain a wind speed vector field with high spatiotemporal resolution and the vertical wind speed component.

[0040] The background field constraint term in the cost function ensures that the inversion result is meteorologically reasonable, the observation constraint term ensures consistency with the measured data, and the regularization term ensures the smoothness and stability of the solution; the conjugate gradient method converges 3-5 times faster than the traditional gradient descent method and reduces memory consumption by 50%; the vertical wind speed obtained by inversion has an accuracy of up to 0.1m / s, meeting the high-precision requirements of flux calculation.

[0041] The variational data assimilation technique is used to invert the three-dimensional atmospheric wind field distribution. Based on the wind field inversion, the variational assimilation cost function containing background field terms, observation terms and regularization terms is minimized. Expressed as: ; Item 1 is the background field constraint term, To find the three-dimensional wind field, are the east-west, north-south and vertical wind speed components, respectively. Prior wind fields provided for numerical weather forecasts, is the background error covariance matrix, which reflects the credibility of the background field. The diagonal elements are the squares of the standard deviation of the background wind field. The standard deviation of the horizontal wind speed error is 2 m / s, and the standard deviation of the vertical wind speed error is 0.5 m / s.

[0042] Item 2 is the observation constraint, is the radial velocity measured by the lidar, is the observation error covariance matrix, reflecting the radar measurement accuracy. The diagonal elements are the radial velocity observation error variances, which are determined based on the signal-to-noise ratio. The typical experimental value is 0.1-1.0 (m / s)².

[0043] Item 3 Smoothing regularization term, is the regularization parameter, and its value range is 10 -6 -10 -3 , the optimal value is determined by the L-curve method, is the Laplace operator; is the lidar observation operator, which describes the geometric projection relationship between the three-dimensional wind field and the radar radial velocity: ; is the azimuth of the laser beam relative to true north (0°-360°); is the elevation angle of the laser beam relative to the horizontal plane (0°-90°).

[0044] The conjugate gradient method is used to solve the cost function, and the convergence criterion is that the gradient modulus is less than 10 -6 or the number of iterations exceeds 1000.

[0045] In step S4, assuming that the dust particles obey the modified particle spectrum distribution model, the particle spectrum parameters are inverted using the multi-wavelength extinction coefficient observation data, and the number concentration distribution of the dust particles and the vertical profile of the dust mass concentration are calculated by solving the integral equation of the extinction coefficient and the particle spectrum distribution.

[0046] In the modified particle size distribution model, the dust particle number concentration The distribution with radius follows a modified power-law spectrum: ;in, Represents the spectral coefficient in units of , determines the total particle concentration; It represents the spectral index, dimensionless, with a value of 2-4, describing the steepness of the particle size distribution; Indicates the exponential cutoff radius, in units of , controls the cutoff of large particles, with a value of 10-50μm; and Indicates the lower and upper limits of the particle radius, usually 0.1μm and 100μm.

[0047] The relationship between the extinction coefficient of each wavelength and the particle spectrum is: ; The dust mass concentration is obtained by integrating the particle volume: ; in, Indicates wavelength The extinction coefficient at ; Represents the particle radius ,wavelength , complex refractive index Mie scattering extinction efficiency factor under conditions; Indicates the density of dust particles, which is 2.6 g / cm³; Indicates mass concentration in units of , by minimizing the difference between the observed and theoretical extinction coefficients, the various inversion spectrum parameters are obtained: ,in, is the weight factor for each wavelength, which is determined according to the measurement accuracy.

[0048] In the calculation of Mie scattering extinction efficiency, the complex refractive index of dust Using a wavelength-dependent model: ; Real and imaginary parts use empirical formulas respectively; Real wavelength : ; is the real part of the refractive index in the high frequency limit, 1.53, is the dispersion intensity parameter, in units of , the value is 0.01; is the resonant wavelength, in units of , usually in the ultraviolet band; imaginary wavelength : ; is the absorption coefficient at the reference wavelength, ranging from 0.001 to 0.01; As the reference wavelength, take 550nm; is the absorption index, which is taken as 1.2. The wavelength-dependent model takes into account the intrinsic optical properties of dust mineral components (such as quartz and feldspar) and can accurately describe the refractive index changes in the visible-near-infrared band.

[0049] In step S5, the obtained vertical wind speed component is multiplied by the dust mass concentration to calculate the instantaneous dust vertical flux; error propagation analysis and statistical uncertainty quantification are performed, and 10,000 random samplings are performed using the Monte Carlo method to provide a confidence interval for the flux measurement results.

[0050] Adopting an adaptive weighting method based on measurement uncertainty to fuse the inversion results; The inversion methods include: a dual-wavelength inversion method based on 355nm and 532nm, a dual-wavelength inversion method based on 532nm and 1064nm, and a three-wavelength joint inversion method, totaling three independent inversion methods. The weight of each inversion method is determined based on the inverse of the uncertainty of each inversion method: ; Uncertainty of each inversion method Contains three components: ;in, The instrument system error, including laser power stability and detector noise, is about 5-10%; Processing errors for the algorithm; is the inversion process error; Optimal fusion flux calculation: ;in, For the The mass concentration inverted by the two methods (μg / m³) is the corresponding vertical wind speed (m / s), is the vertical flux (μg / (m²·s)).

[0051] The output is evaluated by randomly perturbing the input parameters to quantify the measurement uncertainty of the vertical dust flux; the input parameter perturbations are set to follow a Gaussian distribution for the backscattering coefficient with a standard deviation of 5% of the observed value; the extinction coefficient perturbations follow a Gaussian distribution with a standard deviation of 8% of the observed value; the radial velocity perturbations follow a Gaussian distribution with a standard deviation of 0.3 m / s; the real part of the complex refractive index perturbations range is ±0.01, and the imaginary part perturbations range is ±0.001.

[0052] The flux uncertainty is calculated as: ;in, For the The flux value of the Monte Carlo experiment, is the number of Monte Carlo trials, which is 10000; is the mean flux Bootstrap resampling is used to estimate the flux distribution: ;in, and They are distributed by Bootstrap and Quantile, is the significance level, corresponding to the 95% confidence interval; the number of bootstrap resampling is 1000 times.

[0053] Real-time data quality control module ensures the reliability of inversion results; signal-to-noise ratio quality control establishes data screening criteria based on signal-to-noise ratio: ;in, is the echo signal power at distance z, is the noise standard deviation, is the minimum signal-to-noise ratio threshold (2-3).

[0054] Physical consistency test verifies whether the inversion results meet the mass conservation constraints: Among them, S(z,t) is the source term (such as sand blowing from the ground) and D(z,t) is the sedimentation term.

[0055] Adopt anomaly detection algorithm based on statistical distribution: ;in, is the median, MAD is the absolute median deviation, when When , it is determined to be an outlier and marked or removed.

[0056] This method is applicable to lidar systems on various platforms, including ground-based, vehicle-mounted, ship-mounted, and airborne. It has continuous observation capabilities, high temporal resolution (1-10 seconds), and vertical resolution of 7.5-30 meters. Before multi-platform data fusion, quality unification processing is performed, including spatiotemporal interpolation, resolution matching, and accuracy weighting, to ensure the comparability and compatibility of data from different platforms.

[0057] The dust vertical flux products generated by this method can be applied to business areas such as weather forecasting, environmental monitoring and climate research; the temporal resolution is adjustable from 1 minute to 1 hour, the spatial resolution is 7.5-30 meters vertically and 100 meters to 1 kilometer horizontally, the relative error of flux measurement is less than 30%, and the NetCDF format complies with the WMO and CF standards, and the real-time product delay is less than 10 minutes.

[0058] Establish a comprehensive verification system including satellite remote sensing, ground observations and numerical models.

[0059] Correlation verification: .

[0060] Bias Assessment: .

[0061] The specific implementation methods described above further illustrate the objectives, technical solutions and beneficial effects of the present invention in detail. It should be understood that the above description is only a specific implementation method of the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A high-precision inversion method for vertical dust flux based on lidar, characterized by: Atmospheric scattering signals are acquired by a multi-wavelength polarization lidar system. After signal preprocessing, particle identification, wind field inversion, mass concentration calculation, and flux fusion, high-precision quantitative measurement of the vertical transmission flux of dust particles in the atmosphere is achieved. The method specifically includes the following steps: Step S1: Using a three-wavelength polarization lidar system with operating wavelengths of 355 nm, 532 nm, and 1064 nm, laser pulses are emitted into the atmosphere to receive backscattered echo signals from atmospheric aerosol particles; the original echo signals are subjected to distance correction, background noise subtraction, and signal smoothing to obtain backscattering coefficient and depolarization ratio profile data at different wavelengths and polarization states; Step S2: constructing a multidimensional feature vector space based on the multi-wavelength backscattering coefficient ratio, depolarization ratio, and signal gradient optical characteristic parameters to distinguish dust particles from non-dust components in the atmosphere, thereby locating and identifying the boundaries of the dust layer; Step S3: Combining the lidar Doppler radial velocity observation data and the background wind field information of the numerical weather forecast, the variational data assimilation technique is used to invert the three-dimensional atmospheric wind field distribution; the conjugate gradient method is used to solve the cost function including the background field constraint, the observation constraint, and the physical constraint to obtain the wind speed vector field with high spatiotemporal resolution and the vertical wind speed component; Step S4: Assuming that the dust particles obey the modified particle spectrum distribution model, the particle spectrum parameters are inverted using the multi-wavelength extinction coefficient observation data. By solving the integral equation of the extinction coefficient and the particle spectrum distribution, the number concentration distribution of the dust particles and the vertical profile of the dust mass concentration are calculated; In step S5, the obtained vertical wind speed component is multiplied by the dust mass concentration to calculate the instantaneous dust vertical flux; error propagation analysis and statistical uncertainty quantification are performed, and 10,000 random samplings are performed using the Monte Carlo method to provide a confidence interval for the flux measurement results.

2. The method according to claim 1, characterized in that The signal smoothing process uses a sliding average algorithm for time smoothing with a window length of 30 seconds and a vertical resolution of 7.5 meters for spatial smoothing.

3. The method according to claim 2, characterized in that In step S2, dust particle identification is achieved based on multidimensional optical feature vectors and support vector machine classifiers; the multidimensional optical feature vector is a five-dimensional feature vector , used to characterize the optical properties of particles, ;in, , represents backscattering at near-infrared wavelengths Backscattering with visible wavelengths The ratio of , representing the backscattering of ultraviolet wavelengths Backscattering ratio with visible wavelength; , represents the linear depolarization ratio at 532 nm wavelength; , represents the vertical gradient of the backscattering coefficient at 532 nm with height z; , represents the normalized signal intensity after 532nm wavelength distance square correction, is the emitted laser power at a wavelength of 532 nm; The support vector machine classifier uses a support vector machine with a radial basis function kernel for binary classification, and the decision function is: ;in, is the Gaussian radial basis kernel function, is the kernel width parameter, , is the Lagrange multiplier; is the training sample label, +1 represents dust particles, and -1 represents non-dust components; is the bias parameter, is the total number of support vectors.

4. The method according to claim 3, characterized in that The training samples for the support vector machine classifier come from lidar observation data under known dust events and clean atmospheric conditions. Sample annotation is performed by combining simultaneous PM10 ground observations, satellite remote sensing, and back-trajectory analysis. The number of dust samples should be no less than 5,000, and the number of non-dust samples should be no less than 10,000. The dust layer boundary identification adopts a threshold judgment method, and when the dust probability of three consecutive height layers is greater than 0.8, it is determined to be the dust layer boundary.

5. The method according to claim 4, characterized in that In step S3, the variational data assimilation technique is used to invert the three-dimensional atmospheric wind field distribution, and the variational assimilation cost function including the background field term, the observation term and the regularization term is minimized based on the wind field inversion. Expressed as: ； in, To find the three-dimensional wind field, are the east-west, north-south and vertical wind speed components, respectively. Prior wind fields provided for numerical weather forecasts, is the background error covariance matrix; is the radial velocity measured by the lidar, is the observation error covariance matrix; is the regularization parameter, is the Laplace operator; is the lidar observation operator, which describes the geometric projection relationship between the three-dimensional wind field and the radar radial velocity: ; is the azimuth of the laser beam relative to the true north; is the elevation angle of the laser beam relative to the horizontal plane; the cost function is solved by the conjugate gradient method, and the convergence criterion is that the gradient modulus is less than 10 -6 or the number of iterations exceeds 1000.

6. The method according to claim 5, characterized in that In the step S4, in the modified particle size distribution model, the dust particle number concentration The distribution with radius follows a modified power-law spectrum: ;in, represents the spectral coefficient; represents the spectral index; represents the exponential cutoff radius; and Indicates the lower and upper limits of the particle radius; The relationship between the extinction coefficient of each wavelength and the particle spectrum is: ; The dust mass concentration is obtained by integrating the particle volume: ; in, Indicates wavelength The extinction coefficient at ; Represents the particle radius ,wavelength , complex refractive index Mie scattering extinction efficiency factor under conditions; Indicates the density of dust particle material; represents the mass concentration, and the various inversion spectrum parameters are obtained by minimizing the difference between the observed and theoretical extinction coefficients.

7. The method according to claim 6, characterized in that In the calculation of Mie scattering extinction efficiency, the complex refractive index of dust Using a wavelength-dependent model: ; The real and imaginary parts use empirical formulas respectively; the real part wavelength : ; is the real part of the refractive index in the high frequency limit, 1.53, is the dispersion intensity parameter, is the resonant wavelength; the imaginary wavelength : ; is the absorption coefficient at the reference wavelength; As the reference wavelength, take 550nm; is the absorption index, which is taken as 1.

2.

8. The method according to claim 7, characterized in that In step S5, an adaptive weighting method based on measurement uncertainty is used to fuse the inversion results; The inversion methods include: a dual-wavelength inversion method based on 355nm and 532nm, a dual-wavelength inversion method based on 532nm and 1064nm, and a three-wavelength joint inversion method, totaling three independent inversion methods. The weight of each inversion method is determined based on the inverse of the uncertainty of each inversion method: ; Uncertainty of each inversion method Contains three components: ;in, is the instrument system error, including laser power stability and detector noise; Processing errors for the algorithm; is the inversion process error; Optimal fusion flux calculation: ;in, For the The mass concentration inverted by the two methods (μg / m³) is the corresponding vertical wind speed (m / s), is the vertical flux (μg / (m²·s)).

9. The method according to claim 8, characterized in that The output is evaluated by randomly perturbing the input parameters to quantify the measurement uncertainty of the vertical dust flux; the input parameter perturbations are set to the backscattering coefficient perturbation following a Gaussian distribution with a standard deviation of 5% of the observed value; the extinction coefficient perturbation follows a Gaussian distribution with a standard deviation of 8% of the observed value; the radial velocity perturbation follows a Gaussian distribution with a standard deviation of 0.3 m / s; the real part perturbation range of the complex refractive index is ±0.01, and the imaginary part perturbation range is ±0.

001.

10. The method according to claim 9, characterized in that The flux uncertainty is calculated as: ;in, For the The flux value of the Monte Carlo experiment, is the number of Monte Carlo trials, which is 10000; is the mean flux; Bootstrap resampling is used to estimate the flux distribution: ;in, and They are distributed by Bootstrap and Quantile, is the significance level, corresponding to the 95% confidence interval; the number of bootstrap resampling is 1000 times.

Images