CN112526548B - Rainfall identification method and device based on wind-measuring laser radar - Google Patents
Rainfall identification method and device based on wind-measuring laser radar Download PDFInfo
- Publication number
- CN112526548B CN112526548B CN202011291561.3A CN202011291561A CN112526548B CN 112526548 B CN112526548 B CN 112526548B CN 202011291561 A CN202011291561 A CN 202011291561A CN 112526548 B CN112526548 B CN 112526548B
- Authority
- CN
- China
- Prior art keywords
- power spectrum
- rainfall
- signal
- peak
- spectrum
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 51
- 238000001228 spectrum Methods 0.000 claims abstract description 210
- 239000000443 aerosol Substances 0.000 claims abstract description 44
- 230000002902 bimodal effect Effects 0.000 claims abstract description 44
- 238000007689 inspection Methods 0.000 claims description 19
- 238000001514 detection method Methods 0.000 claims description 16
- 230000003595 spectral effect Effects 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 claims description 7
- 238000000926 separation method Methods 0.000 abstract description 8
- 238000009792 diffusion process Methods 0.000 abstract description 4
- 230000009286 beneficial effect Effects 0.000 abstract 1
- 230000001427 coherent effect Effects 0.000 description 3
- 238000005259 measurement Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000001556 precipitation Methods 0.000 description 2
- 238000010998 test method Methods 0.000 description 2
- 229920002430 Fibre-reinforced plastic Polymers 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 239000011151 fibre-reinforced plastic Substances 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000009966 trimming Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
- G01S17/88—Lidar systems specially adapted for specific applications
- G01S17/95—Lidar systems specially adapted for specific applications for meteorological use
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01W—METEOROLOGY
- G01W1/00—Meteorology
- G01W1/10—Devices for predicting weather conditions
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Environmental & Geological Engineering (AREA)
- Remote Sensing (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Electromagnetism (AREA)
- Computer Networks & Wireless Communication (AREA)
- Life Sciences & Earth Sciences (AREA)
- Atmospheric Sciences (AREA)
- Biodiversity & Conservation Biology (AREA)
- Ecology (AREA)
- Environmental Sciences (AREA)
- Optical Radar Systems And Details Thereof (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
The invention discloses a rainfall identification method and device based on a wind lidar. The invention can identify rainfall and separate wind speed and rain speed by using the skewness coefficient of the power spectrum of the wind-measuring laser radar. The basic principle is that under the rainfall condition, aerosol and raindrop signals exist simultaneously, the power spectrum has obvious skewness, and the rainfall power spectrum is firstly identified according to the skewness. The identified rainfall power spectrum is then subjected to a bimodal test. And (3) directly performing bimodal fitting on the power spectrum with the obvious bimodal structure, and performing bimodal fitting again on the power spectrum without the obvious bimodal structure through fewer fitting parameters. Finally, the separation of wind speed and rain speed is realized. The method is applied to guaranteeing aviation safety, is beneficial to improving a climate model, a weather forecast model and an atmospheric pollution diffusion model, and has the advantages of high flexibility, high space-time resolution and high accuracy.
Description
Technical Field
The invention relates to the field of wind lidar and atmosphere detection, in particular to a rainfall identification method and a wind speed and rain speed separation method based on the wind lidar.
Background
High spatial-temporal resolution rainfall detection engineering applications are critical, for example, frozen precipitation may reduce the aerodynamic performance of an aircraft. Accurate measurement of rainfall enables evaluation of the effect of the artificial rainfall operation. In addition, rainfall detection plays an important role in understanding the microscopic physical process of an atmospheric boundary layer, improving a meteorological model, forecasting weather and the like.
Compared with the traditional rainfall detection equipment, the wind-measuring laser radar has the advantages of multiple scanning modes, high space-time resolution and high precision. Lidar still has some problems in rainfall identification. Under rainfall conditions, errors are generated in wind speed inversion due to interference of rain spectrums. At present, rainfall can be judged through the increase of the spectrum width, but wind shear, turbulence and multiple scattering of cloud layers can all cause the increase of the spectrum width, and noise interference is identified to the rainfall. When the aerosol signal is covered by the raindrop signal (or when the raindrop signal is covered by the aerosol signal), and the bimodal structure is not obvious, the spectrum width is not obviously increased, and the rainfall is difficult to identify.
Disclosure of Invention
In order to solve the problems, the invention aims to provide a method and a device for identifying rainfall based on laser radar. On the basis of using only one coherent wind lidar, according to the power spectrum of the wind lidar, the rainfall identification and the wind speed and rain speed separation are realized. Has important significance for pollution diffusion modeling, atmospheric weather mode, aviation safety and the like.
The invention is realized in the following way: a rainfall identification method based on wind lidar comprises the following steps:
acquiring original echo data of a laser radar;
performing fast Fourier transform on the original echo data to obtain a power spectrum;
calculating the skewness coefficient of the power spectrum;
when the absolute value of the skewness coefficient is larger than a preset first threshold value, identifying the absolute value as a rainfall power spectrum;
performing double-peak inspection and spectrum width inspection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to an inspection result; classifying the power spectrum with obvious double-peak structure into a type A; classifying power spectrums without a bimodal structure and with the spectrum width obviously larger than a threshold value into a class, and marking the class as a class B; classifying the power spectrum without a bimodal structure and with the spectrum width smaller than a threshold value into a class, and marking the class as a class C;
and respectively carrying out double-peak fitting on various rainfall power spectrums, and separating the wind speed and the rain speed.
Further, the calculating the skewness coefficient of the power spectrum includes:
the skewness coefficient Sk of the power spectrum is the skewness degree of the strongest position of the relative signal, and the calculating method comprises the following steps:
wherein v is i For Doppler velocity of the power spectrum, P i For the signal strength of the power spectrum,the speed corresponding to the strongest signal position; sk (Sk) 0 Is the initial skewness.
Further, the method comprises the steps of: calculating a degree of skew of the signal relative to the mean position, calculating a pearson first skewness coefficient, calculating a pearson second skewness coefficient, calculating a medcoupling skewness, calculating one of Groeneveld and Meeden coefficients.
Further, after performing the fast fourier transform on the raw echo data to obtain a power spectrum, the method further includes:
and carrying out spline interpolation on the obtained power spectrum.
Further, the method for carrying out bimodal inspection on the rainfall power spectrum according to the preset method comprises the following steps:
1) Performing cubic spline interpolation on the power spectrum to obtain a smooth power spectrum curve and a normalized power spectrum f (x);
2) Calculating a first derivative f' (x) of the calculated power spectrum curve, and normalizing;
3) Finding the zero point x of the first derivative f' (x) i If f' (x i-1 )>0, and f' (x) i+1 )<0, and signal strength f (x i ) If the signal is larger than the preset minimum signal threshold value, judging x i There are distinct peak points; when the zero value points are more than 1, the rainfall power spectrum is judged to have an obvious multimodal structure.
Further, the method for carrying out bimodal inspection on the rainfall power spectrum according to the preset method comprises the following steps:
if Sk is min <|Sk|<Sk min2 Judging that the rainfall power spectrum does not have an obvious double-peak structure; if |Sk|>Sk min2 Judging that the rainfall power spectrum has an obvious double-peak structure; where Sk is the skewness coefficient of the power spectrum, sk min For a preset first bias threshold Sk min2 Is a preset second bias threshold.
Further, respectively carrying out double-peak fitting on various rainfall power spectrums, separating wind speed and rain speed, and comprising the following steps:
if the rainfall power spectrum is of type A, adopting a double Gaussian model to perform double peak fitting on the rainfall power spectrum:
wherein I is a For aerosol signal intensity, I r Is the raindrop signal intensity, f is the Doppler frequency, f a Is the Doppler frequency of the aerosol, f r For Doppler frequency, sigma, of raindrops a Is aerosol spectrum wide, sigma r The rain drop spectrum is wide.
Further, respectively carrying out double-peak fitting on various rainfall power spectrums, separating wind speed and rain speed, and comprising the following steps:
if the rainfall power spectrum is of type B:
selecting adjacent rainfall power spectrums with the same height, fixing (f r -f a ) Performing double-peak fitting;
alternatively, the spectral width σ of the aerosol is fixed a Bimodal fitting was performed.
Further, respectively carrying out double-peak fitting on various rainfall power spectrums, separating wind speed and rain speed, and comprising the following steps:
if the rainfall power spectrum is of type C, then:
when Sk is>At 0, the main peak of the signal is an aerosol signal S a (f)，S a (f) From the strongest signal position x m Signal component on left side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the raindrop signal S r (f) S, i.e r (f)＝S(f)-S a (f)；
When Sk is<At 0, the main peak of the signal is raindrop signalS r (f)，S r (f) From the strongest signal position x m Signal component on right side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the aerosol signal S a (f) S, i.e a (f)＝S(f)-S r (f)。
It should be noted that the invention also includes the steps of calculating the wind speed and the rain speed.
Respectively carrying out double-peak fitting on various rainfall power spectrums, and after separating wind speed and rain speed, further comprising: and calculating the wind speed and rain speed information according to the position information of the aerosol signals and the rain drop signals corresponding to the power spectrum and combining the Doppler principle.
Correspondingly, the invention also provides a rainfall identification device based on the wind-measuring laser radar, which comprises:
the original data acquisition module is used for acquiring original echo data of the laser radar;
the power spectrum acquisition module is used for performing fast Fourier transform on the original echo data to obtain a power spectrum;
the skewness coefficient calculation module is used for calculating the skewness coefficient of the power spectrum;
the rainfall spectrum pre-judging module is used for identifying a rainfall power spectrum when the absolute value of the skewness coefficient is larger than a preset first threshold value;
the detection classification module is used for carrying out double-peak detection and spectrum width detection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to detection results; classifying the power spectrum with obvious double-peak structure into a type A; classifying power spectrums without a bimodal structure and with the spectrum width obviously larger than a threshold value into a class, and marking the class as a class B; classifying the power spectrum without a bimodal structure and with the spectrum width smaller than a threshold value into a class, and marking the class as a class C;
and the double-peak fitting module is used for respectively carrying out double-peak fitting on various rainfall power spectrums and separating the wind speed and the rain speed.
In summary, the invention provides a rainfall identification method and device based on a wind-finding laser radar, which are characterized in that original echo data of the laser radar are obtained, and fast Fourier transform is performed on the original echo data to obtain a power spectrum; calculating the skewness coefficient of the power spectrum; when the absolute value of the skewness coefficient is larger than a preset first threshold value, identifying the absolute value as a rainfall power spectrum; and (3) carrying out bimodal inspection and spectrum width inspection on the rainfall power spectrum according to a preset method, classifying the power spectrum according to an inspection result, and respectively carrying out bimodal fitting on various rainfall power spectrums to realize separation of wind speed and rain speed. On the basis of using only one coherent wind-finding laser radar, the method can accurately identify, separate and calculate the wind speed and the rain speed by only one measurement of the laser radar and analysis and calculation of the power spectrum of the laser radar echo signal, and has important significance for pollution diffusion modeling, atmospheric weather mode, aviation safety and the like.
Drawings
In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings that are needed in the description of the embodiments will be briefly described below, it being obvious that the drawings in the following description are only some embodiments of the present invention, and that other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of a rainfall identification method based on a wind lidar according to an embodiment of the present invention;
FIG. 2 is a graph showing typical power spectrum shapes and wind speed and rain speed separation results under rainfall conditions provided by the embodiment of the invention;
fig. 3 is a flowchart of a rainfall identification method based on wind lidar according to an embodiment of the present invention.
Fig. 4 is a block diagram of a rainfall identification device based on a wind lidar according to an embodiment of the present invention.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and fully with reference to the accompanying drawings, in which it is evident that the embodiments described are only some, but not all embodiments of the invention. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to fall within the scope of the invention.
Examples:
fig. 1 is a flowchart of a rainfall identification method based on a wind lidar according to an embodiment of the present invention, and fig. 2 is a typical power spectrum shape and a wind speed and rain speed separation result under rainfall conditions provided by an embodiment of the present invention; fig. 3 is a flowchart of a rainfall identification method based on wind lidar according to an embodiment of the present invention.
In connection with fig. 1 to 3, the following is implemented: the invention discloses a rainfall identification method based on a wind lidar, which comprises the following steps:
s1, acquiring original echo data of a laser radar.
The method only needs to acquire the original echo data of the single laser radar.
S2, performing fast Fourier transform on the original echo data to obtain a power spectrum.
S3, calculating the skewness coefficient of the power spectrum obtained in the step 2.
S4, judging whether the absolute value of the skewness coefficient of the power spectrum is larger than a preset first threshold value, and recognizing the power spectrum as the rainfall power spectrum when the absolute value of the skewness coefficient is larger than the first threshold value.
The basic principle of the invention is that when only an aerosol signal exists, the power spectrum becomes a unimodal structure, as shown in fig. 2 (a), the spectrum width is smaller, and the deviation is close to 0. The cloud layer multiple scattering and turbulence can both cause the spectrum width to be increased, but only aerosol signals are generated, and the deflection change cannot be obviously changed. The single pair of skewness has less effect, as shown in fig. 2 (b) and (c), and when both a rain signal and an aerosol signal are present, the skewness of the strongest position of the relative signal changes significantly, as shown in fig. 2 (d) - (l). Although turbulence and wind shear can also lead to an increase in the spectral width, only the aerosol signal, and the deflection change, will not change significantly.
Preferably, the radar power spectrum skewness Sk calculation method is the skewness degree of the strongest position of the relative signal, namely:
wherein v is i And P i Representing the doppler velocity and signal strength of the power spectrum,representing the velocity corresponding to the strongest position of the signal. Sk (Sk) 0 The initial bias is the result of the radar system itself (e.g., the pulse shape of the laser). Sk (Sk) 0 The power spectrum is estimated under the clear sky condition of weak turbulence and weak wind shear.
When |Sk| is greater than the threshold Sk min If so, it is determined that there is rainfall. Sk in the present embodiment min Taking 0.2.
The degree of deviation is close to 0 under the condition of no rainfall. The speed towards the lidar is defined as positive and vice versa. Under rainfall conditions, the raindrop signal is always on the right side of the aerosol signal, as the raindrop drops faster than the aerosol. When Sk is>Sk min Indicating the presence of rainfall and stronger aerosol signal as shown in fig. 2 (f), (h), (j) and (l). When Sk is<-Sk min When rainfall is present, and the raindrop signal is stronger, as shown in fig. 2 (d), (j), (i) and (k).
The definition of the degree of skewness is in addition to the degree of skewness described above with respect to the strongest signal location. Other classical definition methods of statistics, i.e. the degree of skew of the signal with respect to the mean position, may also be used. Further alternatively, other measures of the degree of deflection of the curve, such as pearson first skewness coefficient (modal skewness), pearson second skewness coefficient (median skewness), medcoupling skewness, groeneveld & Meeden coefficients, and so forth.
When the power spectrum is sparse, spline interpolation can be performed under the condition of keeping the original power spectrum shape in order to improve the accuracy of the skewness, and sampling points participating in the skewness calculation are increased.
S5, carrying out double-peak inspection and spectrum width inspection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to an inspection result. The power spectra with a distinct bimodal structure are classified as type a, as shown in fig. 2 (d), (e), (f), (i), (j); the power spectra without a bimodal structure but with a spectral width significantly greater than the threshold are classified as type B, as shown in fig. 2 (g) and (h); the power spectra without a bimodal structure and with a spectral width less than the threshold are classified as type C, as shown in fig. 2 (k) and (l).
A bimodal test is first performed on the power spectrum identified as rainfall.
An alternative bimodal test method is as follows:
1) And (3) performing cubic spline interpolation on the power spectrum to obtain a smooth power spectrum curve and a normalized power spectrum f (x).
2) Calculating the first derivative f of the calculated power spectrum curve ′ (x) And normalizing.
3) Finding the zero point x of the first derivative f' (x) i If f' (x i-1 )>0, and f' (x) i+1 )<0, and signal strength f (x i ) Greater than the threshold of the minimum signal, x i There is a distinct peak point. When more than 1 such zero point exists, then there is a multimodal structure that proves significant in the power spectrum.
Optionally, the power spectrum double peak test method comprises the following steps: when SK min <|Sk|<Sk min2 When the power spectrum does not have a distinct bimodal structure. When |Sk|>Sk min2 When the power spectrum has a distinct bimodal structure. The principle is as follows: the more pronounced the bimodal structure, the greater the corresponding skewness.
S6, respectively carrying out double-peak fitting on various rainfall power spectrums, and separating the wind speed and the rain speed.
The power spectrum (type a) with a distinct bimodal structure for type a can be directly subjected to a bimodal fit. Bimodal fits may employ a double gaussian model fit:
wherein I and f, σ are signal strength, doppler frequency and spectral width, respectively. Subscripts a and r represent aerosol and raindrops, respectively. The aerosol and raindrop signal classification results are shown in fig. 2 (d), (e), (f), (i), (j).
For type B, the power spectrum does not have a single spectrum width of a bimodal structure and is obviously increased, which indicates that the rainfall power spectrum at the moment has equivalent aerosol signals and raindrop signals of intensity, and the direct bimodal fit can cause larger error due to more parameters of the bimodal fit and no obvious bimodal structure. The following two methods can be adopted to reduce fitting parameters and improve the accuracy of the double-peak fitting. The aerosol and raindrop signal classification results are shown in fig. 2 (i) and (j).
Method 1: selecting the adjacent power spectrums with the same height, fixing (f r -f a ) A bimodal fit is performed simultaneously. The principle is as follows:
for the fixed vertical scanning mode, the radar can only detect Doppler frequency in the vertical direction, the vertical wind speed and the vertical rain speed in adjacent time are almost unchanged at the same height, and the (f) of the adjacent power spectrum at the same height corresponds to the same height r -f a ) Almost unchanged.
For cone scanning mode, under precipitation conditions, radial velocity of aerosol and raindropsAnd->The method comprises the following steps:
wherein V is || ,θ 0 And V ⊥ Respectively horizontal speed, horizontal speed direction and vertical speed. Alpha i Andis the azimuth and elevation of the radar. Aerosol and raindrops in the horizontal direction due to wind speed dragging effectIs nearly the same, the radial velocity difference dV between the two i The method comprises the following steps:
i.e. the difference in velocity of aerosol and raindrops is mainly represented by the difference in vertical velocity when the elevation angle is constant. The vertical wind speed and the vertical rain speed are almost unchanged, and the corresponding (f) of the adjacent power spectrums with the same height r -f a ) Almost unchanged.
Thus, in a bimodal fit, the same highly adjacent power spectrum corresponds to (f r -f a ) And the double-peak fitting is performed simultaneously, so that the degree of freedom of the double fitting is reduced, and the fitting accuracy is improved.
Method 2: fixing the spectral width sigma of the aerosol a And (5) performing bimodal fitting.
The spectrum width of the aerosol is almost unchanged under the conditions of weak wind trimming and weak turbulence, and can be obtained by power spectrum estimation under the condition of sunny days.
The above two methods for reducing the fitting parameters can be used simultaneously according to the actual weather conditions or one of them can be selected for use.
For type C, the power spectrum does not have a distinct bimodal structure and the spectral width increase is not distinct, with one signal component dominant. The other signal has a significantly weaker dominant signal. The wind speed and rain speed separation can be carried out by adopting a following method.
When Sk is>When 0, the power spectrum is right-biased, the aerosol signal is dominant, i.e. the main peak is the aerosol signal, the aerosol signal S a (f) Can be represented by the strongest signal position x m Signal component on left side relative x=x m Axisymmetric. Raindrop signal S r (f) Subtracting the main peak from the power spectrum to obtain a residual peak, namely: s is S r (f)＝S(f)-S a (f) A. The invention relates to a method for producing a fibre-reinforced plastic composite The aerosol and raindrop signal classification results are shown in fig. 2 (l).
When Sk is<When 0, the power spectrum is left-biased, the rain drop signal is dominant, i.e. the main peak is the rain drop signal, the rain drop signal S r (f) Can be represented by the strongest signal position x m Signal component on right side relative x=x m Axisymmetric. Aerosol signal S a (f) Subtracting the main peak from the power spectrum to obtain a residual peak, namely: s is S a (f)＝S(f)-S r (f) The aerosol and raindrop signal classification results are shown in fig. 2 (k).
Fig. 4 is a block diagram of a rainfall identification device based on a wind lidar according to an embodiment of the present invention. As shown in fig. 4, the rainfall identification device based on wind lidar of the present invention further includes:
the original data acquisition module is used for acquiring original echo data of the laser radar;
the power spectrum acquisition module is used for performing fast Fourier transform on the original echo data to obtain a power spectrum;
the skewness coefficient calculation module is used for calculating the skewness coefficient of the power spectrum;
the rainfall spectrum pre-judging module is used for identifying a rainfall power spectrum when the absolute value of the skewness coefficient is larger than a preset first threshold value;
the detection classification module is used for carrying out double-peak detection and spectrum width detection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to detection results; classifying the power spectrum with obvious double-peak structure into a type A; classifying power spectrums without a bimodal structure and with the spectrum width obviously larger than a threshold value into a class, and marking the class as a class B; classifying the power spectrum without a bimodal structure and with the spectrum width smaller than a threshold value into a class, and marking the class as a class C;
and the double-peak fitting module is used for respectively carrying out double-peak fitting on various rainfall power spectrums and separating the wind speed and the rain speed.
In summary, the invention provides a rainfall identification method and device based on a wind-finding laser radar, which are characterized in that original echo data of the laser radar are obtained, and fast Fourier transform is performed on the original echo data to obtain a power spectrum; calculating the skewness coefficient of the power spectrum; when the absolute value of the skewness coefficient is larger than a preset first threshold value, identifying the absolute value as a rainfall power spectrum; and (3) carrying out bimodal inspection and spectrum width inspection on the rainfall power spectrum according to a preset method, classifying the power spectrum according to an inspection result, and respectively carrying out bimodal fitting on various rainfall power spectrums to realize separation of wind speed and rain speed. On the basis of using only one coherent wind-finding laser radar, the method can accurately identify, separate and calculate the wind speed and the rain speed by only one measurement of the laser radar and analysis and calculation of the power spectrum of the laser radar echo signal, and has important significance for pollution diffusion modeling, atmospheric weather mode, aviation safety and the like.
The foregoing is only a preferred embodiment of the present invention, but the scope of the present invention is not limited thereto, and any changes or substitutions easily contemplated by those skilled in the art within the scope of the present invention should be included in the scope of the present invention. Therefore, the protection scope of the present invention should be subject to the protection scope of the claims.
Claims (6)
1. A rainfall identification method based on a wind-measuring laser radar is characterized by comprising the following steps of: the method comprises the following steps:
acquiring original echo data of a laser radar;
performing fast Fourier transform on the original echo data to obtain a power spectrum;
calculating the skewness coefficient of the power spectrum; the calculating the skewness coefficient of the power spectrum comprises the following steps:
the skewness coefficient Sk of the power spectrum is the skewness degree of the strongest position of the relative signal, and the calculating method comprises the following steps:
wherein v is i For Doppler velocity of the power spectrum, P i For the signal strength of the power spectrum,the speed corresponding to the strongest signal position; sk (Sk) 0 Is the initial skewness; when the absolute value of the skewness coefficient is larger than a preset first threshold value, identifying the absolute value as a rainfall power spectrum;
performing double-peak inspection and spectrum width inspection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to an inspection result; classifying the power spectrum with obvious double-peak structure into a type A; classifying power spectrums without a bimodal structure and with the spectrum width obviously larger than a threshold value into a class, and marking the class as a class B; classifying the power spectrum without a bimodal structure and with the spectrum width smaller than a threshold value into a class, and marking the class as a class C;
respectively carrying out double-peak fitting on the rainfall power spectrum type A and the rainfall power spectrum type B, and separating the wind speed and the rain speed;
if the rainfall power spectrum is of type A, adopting a double Gaussian model to perform double peak fitting on the rainfall power spectrum:
wherein I is a For aerosol signal intensity, I r Is the raindrop signal intensity, f is the Doppler frequency, f a Is the Doppler frequency of the aerosol, f r For Doppler frequency, sigma, of raindrops a Is aerosol spectrum wide, sigma r The rain drop spectrum is wide;
if the rainfall power spectrum is of type B:
selecting adjacent rainfall power spectrums with the same height, fixing (f r -f a ) Performing double-peak fitting;
alternatively, the spectral width σ of the aerosol is fixed a Performing double-peak fitting;
if the rainfall power spectrum is of type C, then:
when Sk > 0, the main signal peak is the aerosol signal S a (f)，S a (f) From the strongest signal position x m Signal component on left side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the raindrop signal S r (f) S, i.e r (f)＝S(f)-S a (f)；
When Sk < 0, the main peak of the signal is the raindrop signal S r (f)，S r (f) From the strongest signal position x m Signal component on right side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the aerosol signal S a (f) S, i.e a (f)＝S(f)-S r (f)。
2. The method of claim 1, wherein calculating the skewness factor of the power spectrum comprises: calculating a degree of skew of the signal relative to the mean position, calculating a pearson first skewness coefficient, calculating a pearson second skewness coefficient, calculating a medcoupling skewness, calculating one of Groeneveld and Meeden coefficients.
3. The method of claim 1, wherein after performing a fast fourier transform on the raw echo data to obtain a power spectrum, further comprising:
and carrying out spline interpolation on the obtained power spectrum.
4. The method of claim 1, wherein the bimodal inspection of the rainfall power spectrum according to a preset method comprises:
1) Performing cubic spline interpolation on the power spectrum to obtain a smooth power spectrum curve and a normalized power spectrum f (x);
2) Calculating a first derivative f' (x) of the calculated power spectrum curve, and normalizing;
3) Finding the zero point x of the first derivative f' (x) i If f' (x i-1 ) > 0, and f' (x) i+1 ) < 0, and signal strength f (x i ) If the signal is larger than the preset minimum signal threshold value, judging x i There are distinct peak points; when the zero value points are more than 1, the rainfall power spectrum is judged to have an obvious multimodal structure.
5. The method of claim 1, wherein the bimodal inspection of the rainfall power spectrum according to a preset method comprises:
if Sk is min ＜|Sk|＜Sk min2 Judging that the rainfall power spectrum does not have an obvious double-peak structure; if |Sk| > Sk min2 Judging that the rainfall power spectrum has an obvious double-peak structure; where Sk is the skewness coefficient of the power spectrum, sk min For a preset first bias threshold Sk min2 Is a preset second bias threshold.
6. Rainfall recognition device based on wind lidar, characterized by comprising:
the original data acquisition module is used for acquiring original echo data of the laser radar;
the power spectrum acquisition module is used for performing fast Fourier transform on the original echo data to obtain a power spectrum;
the skewness coefficient calculation module is used for calculating the skewness coefficient of the power spectrum; the calculating the skewness coefficient of the power spectrum comprises the following steps:
the skewness coefficient Sk of the power spectrum is the skewness degree of the strongest position of the relative signal, and the calculating method comprises the following steps:
wherein v is i For Doppler velocity of the power spectrum, P i For the signal strength of the power spectrum,the speed corresponding to the strongest signal position; sk (Sk) 0 Is the initial skewness; the rainfall spectrum pre-judging module is used for identifying a rainfall power spectrum when the absolute value of the skewness coefficient is larger than a preset first threshold value;
the detection classification module is used for carrying out double-peak detection and spectrum width detection on the rainfall power spectrum according to a preset method, and classifying the power spectrum according to detection results; classifying the power spectrum with obvious double-peak structure into a type A; classifying power spectrums without a bimodal structure and with the spectrum width obviously larger than a threshold value into a class, and marking the class as a class B; classifying the power spectrum without a bimodal structure and with the spectrum width smaller than a threshold value into a class, and marking the class as a class C;
the double-peak fitting module is used for respectively carrying out double-peak fitting on the rainfall power spectrum type A and the rainfall power spectrum type B and separating the wind speed and the rain speed;
if the rainfall power spectrum is of type A, adopting a double Gaussian model to perform double peak fitting on the rainfall power spectrum:
wherein I is a For aerosol signal intensity, I r Is the raindrop signal intensity, f is the Doppler frequency, f a Is the Doppler frequency of the aerosol, f r For Doppler frequency, sigma, of raindrops a Is aerosol spectrum wide, sigma r The rain drop spectrum is wide;
if the rainfall power spectrum is of type B:
selecting adjacent rainfall power spectrums with the same height, fixing (f r -f a ) Performing double-peak fitting;
alternatively, the spectral width σ of the aerosol is fixed a Performing double-peak fitting;
if the rainfall power spectrum is of type C, then:
when Sk > 0, the main signal peak is the aerosol signal S a (f)，S a (f) From the strongest signal position x m Signal component on left side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the raindrop signal S r (f) S, i.e r (f)＝S(f)-S a (f)；
When Sk < 0, the main peak of the signal is the raindrop signal S r (f)，S r (f) From the strongest signal position x m Signal component on right side relative x=x m Axisymmetric to obtain; the residual peak obtained by subtracting the main peak from the power spectrum is the aerosol signal S a (f) S, i.e a (f)＝S(f)-S r (f)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011291561.3A CN112526548B (en) | 2020-11-18 | 2020-11-18 | Rainfall identification method and device based on wind-measuring laser radar |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011291561.3A CN112526548B (en) | 2020-11-18 | 2020-11-18 | Rainfall identification method and device based on wind-measuring laser radar |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112526548A CN112526548A (en) | 2021-03-19 |
| CN112526548B true CN112526548B (en) | 2024-02-23 |
Family
ID=74981140
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011291561.3A Active CN112526548B (en) | 2020-11-18 | 2020-11-18 | Rainfall identification method and device based on wind-measuring laser radar |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112526548B (en) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113176586B (en) * | 2021-03-25 | 2023-05-02 | 北京大学 | Raindrop measuring method and solid-state area array laser raindrop spectrometer |
| CN114637012A (en) * | 2021-10-26 | 2022-06-17 | 北京理工大学 | non-Gaussian weather radar target rapid identification method based on generalized PPP |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2012103050A (en) * | 2010-11-08 | 2012-05-31 | Japan Aerospace Exploration Agency | Method for remote detection of air turbulence, and execution device thereof |
| CN104316929A (en) * | 2014-09-12 | 2015-01-28 | 何平 | Power spectrum identification method and apparatus |
| CN109324328A (en) * | 2018-12-17 | 2019-02-12 | 北京无线电测量研究所 | The method and device of wind profile radar vertical beam spectrum of turbulence when extracting precipitation |
| CN109343063A (en) * | 2018-12-17 | 2019-02-15 | 中国气象局气象探测中心 | A kind of millimeter wave cloud metre clear air echo automatic identifying method and system |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| GB0411097D0 (en) * | 2004-05-19 | 2004-06-23 | Qinetiq Ltd | Laser radar device and method |
| JP6413457B2 (en) * | 2014-08-08 | 2018-10-31 | 株式会社デンソー | Precipitation determination device |
-
2020
- 2020-11-18 CN CN202011291561.3A patent/CN112526548B/en active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2012103050A (en) * | 2010-11-08 | 2012-05-31 | Japan Aerospace Exploration Agency | Method for remote detection of air turbulence, and execution device thereof |
| CN104316929A (en) * | 2014-09-12 | 2015-01-28 | 何平 | Power spectrum identification method and apparatus |
| CN109324328A (en) * | 2018-12-17 | 2019-02-12 | 北京无线电测量研究所 | The method and device of wind profile radar vertical beam spectrum of turbulence when extracting precipitation |
| CN109343063A (en) * | 2018-12-17 | 2019-02-15 | 中国气象局气象探测中心 | A kind of millimeter wave cloud metre clear air echo automatic identifying method and system |
Non-Patent Citations (5)
| Title |
|---|
| Simultaneous wind and rainfall detection by power spectrum analysis using a VAD scanning coherent Doppler lidar;TIANWEN WEI et al.;《Optics Express》;全文 * |
| 基于双高斯拟合的风廓线雷达反演雨滴谱;何越;何平;林晓萌;应用气象学报(第005期);全文 * |
| 基于相干多普勒激光雷达的风雨同时探测;魏天问;夏海云;红外与激光工程(第0s2期);全文 * |
| 降水回波谱参数估算雨滴谱参数的算法研究;阮征;刘褚D;马建立;金龙;;高原气象(第04期);全文 * |
| 风廓线雷达探测降水云体中雨滴谱的试验研究;王晓蕾;阮征;葛润生;陈钟荣;;高原气象(第02期);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112526548A (en) | 2021-03-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| WO2022016884A1 (en) | Method for extracting sea surface wind speed on basis of k-means clustering algorithm | |
| CN111583214A (en) | Method for inverting sea surface wind speed by using marine radar image based on RBF neural network | |
| KR101531224B1 (en) | Quantitative precipitation estimation system based dual polarization radars and method thereof | |
| CN112526547B (en) | Atmospheric boundary layer classification method and device based on wind-measuring laser radar | |
| CN112526548B (en) | Rainfall identification method and device based on wind-measuring laser radar | |
| CN110261857B (en) | Spatial interpolation method for weather radar | |
| CN109901130B (en) | Rotor unmanned aerial vehicle detection and identification method based on Radon transformation and improved 2DPCA | |
| CN105093215A (en) | Doppler information based method for tracking low-altitude low-speed small target through radar | |
| CN111366931A (en) | Cloud precipitation refinement inversion method based on cloud radar | |
| CN112881993B (en) | Method for automatically identifying false flight path caused by radar distribution clutter | |
| CN109116359B (en) | Method for estimating low-altitude wind shear wind field echo wind speed of airborne radar | |
| CN108459323B (en) | Strong-radial-closure-field identification and speed correction method based on Doppler weather radar data | |
| CN104199022A (en) | Target modal estimation based near-space hypersonic velocity target tracking method | |
| CN111142085B (en) | External radiation source radar target classification and identification method based on track feature extraction | |
| CN113808174A (en) | Radar small target tracking method based on full convolution network and Kalman filtering | |
| CN116401932A (en) | Sea fog dissipation time estimation method based on laser radar and millimeter wave radar | |
| JP7234803B2 (en) | Target signal separation device, passive radar device and target signal separation method | |
| CN107783090B (en) | Millimeter wave radar-based radar signal processing method for collision avoidance system of fixed-wing unmanned aerial vehicle | |
| CN111948638B (en) | Target length estimation method based on high-resolution range profile and application thereof | |
| CN114739669A (en) | Rolling bearing state monitoring method and device based on terahertz radar | |
| Kida et al. | Improvement of rain/no-rain classification methods for microwave radiometer observations over the ocean using a 37 GHz emission signature | |
| CN108414997B (en) | Boundary layer height inversion method based on particle characteristic difference | |
| CN116879899A (en) | Method based on aerial precipitation particle spectrum inversion | |
| CN106443678B (en) | Utilize the atmospheric duct monitoring method of wind profile radar and RASS | |
| Li et al. | Improved Moments Estimation for Ground-Based K-Band Doppler Radar Using Cost Function Method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| TA01 | Transfer of patent application right | ||
| TA01 | Transfer of patent application right |
Effective date of registration: 20210812 Address after: 210044 915, meteorological building, Nanjing University of information engineering, Nanjing, Jiangsu Province Applicant after: Zhang Zhen Address before: 215000 room 104, building 3, fengzhouyuan community, Tayuan Road, Huqiu District, Suzhou City, Jiangsu Province Applicant before: Dong Jingjing |
|
| TA01 | Transfer of patent application right | ||
| TA01 | Transfer of patent application right |
Effective date of registration: 20210827 Address after: 211500 room 407-1, building A2, No. 9 Kechuang Avenue, intelligent manufacturing Industrial Park (Zhihe Park), Jiangbei new area, Nanjing, Jiangsu Province Applicant after: Nanjing taiaixin Technology Co.,Ltd. Address before: 210044 915, meteorological building, Nanjing University of information engineering, Nanjing, Jiangsu Province Applicant before: Zhang Zhen |
|
| GR01 | Patent grant | ||
| GR01 | Patent grant |