Innowcasting method for convective storm strong wind

Abstract

The invention discloses a nowcasting method for convective storm strong wind. According to the method, the radar data is subjected to electromagnetic interference echo removal and deblurring processing, data points influencing a prediction result are eliminated, the reliability of the radar data is enhanced, and the accuracy of a strong storm area is improved; performing initial region identification by using a region growing algorithm, extracting severe convection structure information in the radar image, and removing a weak region which is too small in area and may be a misidentification result; the vertical direction in the strong storm area is predicted by combining an L-K optical flow method and a Bhatpacarya coeffcient method, the overall dynamic analysis is performed on the radar echo area, and the potential strong wind activity is finely identified and effectively predicted; and weighted voting is carried out on prediction results corresponding to the SVM, the CNN and the U-Net, so that the prediction dimension is widened, the recognition performance is enhanced, the robustness of the integrated model is improved, and the accuracy of storm prediction in the horizontal direction in the strong storm region is improved.

Description

Technical Field

[0001] The invention relates to the technical field of weather forecasting, and in particular to a nowcasting method for convective storm gale. Background Art

[0002] Thunderstorms, hail, heavy rainfall, etc. caused by severe convection are the main causes of local meteorological disasters. Thunderstorms are a typical type of severe convective weather. Although they do not occur frequently, they can cause strong destructiveness and have an immeasurable impact on human life safety, production economy, national defense and military. In recent years, thunderstorms have occurred frequently, causing serious losses of life and property. Therefore, monitoring such weather phenomena has very important research and social value.

[0003] Because the occurrence of thunderstorms and gales is highly random and uncertain, it is difficult to forecast thunderstorms and gales. The extrapolated nowcasting technology based on weather radar has become an important means of predicting short-term weather changes. This technology can be mainly divided into two categories: one is an algorithm based on radar echo monomer centroid tracking, which mainly relies on the continuous tracking of the centroid of a single radar echo to predict its future motion trajectory and change trend. The second category is the radar echo area tracking algorithm, which pays more attention to the overall dynamic analysis of the radar echo area to infer its future evolution. However, the structure of severe convective weather systems is often extremely complex, and its generation and evolution process is difficult to accurately describe with a simple mathematical model. This makes the timeliness of radar-based nowcasting usually limited when simulating the nonlinear motion of convective systems, generally within 1 hour, and no more than 3 hours at most. Although the mesoscale numerical atmospheric model has been rapidly developed and improved in recent years, it has shown higher potential in predicting the development trend of mesoscale weather events. These models can achieve relatively reliable forecasts for the next 3-6 hours by introducing a variety of physical and chemical parameter schemes, flexible spatial scales, and data assimilation technology to improve the accuracy of simulation and forecasting of local meteorological elements. Whether it is the lack of high-resolution observation data during the cold start process in a cloudless state, or the addition of initial fields in a hot start in a cloudy state, the forecast deviations of these numerical models in the initial stage are often large, and their forecast effects are often inferior to the extrapolated forecasts of radar. In addition, the continuous advancement of science and technology has made the public and meteorologically sensitive industries have higher requirements for the timeliness, effectiveness and accuracy of weather forecasts, and it is particularly important to provide more accurate and efficient short-term and nowcasts. Summary of the invention

[0004] In view of the above-mentioned deficiencies in the prior art, the present invention provides a nowcasting method for convective storm gale, which solves the problem that the occurrence of thunderstorm gale has strong randomness and uncertainty, resulting in great difficulty in forecasting thunderstorm gale and failure to achieve the expected effect, and failure to meet the market's forecasting requirements for thunderstorm gale.

[0005] In order to achieve the above-mentioned object of the invention, the technical solution adopted by the present invention is:

[0006] A method for nowcasting convective storm gale is provided, which comprises the following steps:

[0007] S1. Acquire multi-source observation data of convective storm gale and perform preprocessing to obtain preprocessed multi-source observation data;

[0008] S2, using the region growing algorithm to identify the region of the preprocessed multi-source observation data to obtain the strong storm area;

[0009] S3, using the LK optical flow method and the Bhattacharyya coefficient method to perform storm prediction on the preprocessed multi-source observation data in the vertical direction in the strong storm area, and obtain the corresponding vertical prediction results;

[0010] S4, constructing a horizontal storm prediction model; using the horizontal storm prediction model; performing storm prediction on the preprocessed multi-source observation data in the horizontal direction of the strong storm area to obtain corresponding horizontal prediction results;

[0011] S5. Based on the vertical prediction results and the horizontal prediction results, the final prediction results are determined, that is, the nowcast results of the convective storm gale are obtained.

[0012] The beneficial effects of the present invention are as follows: electromagnetic interference echo and deblurring processing are performed on radar data, data points that affect prediction results are eliminated, and the reliability of radar data is enhanced, so as to improve the accuracy of strong storm areas and prediction results; the regional growing algorithm is used for initial region recognition, so as to effectively extract strong convective structure information in radar images, remove weak areas that are too small and may be misidentified results, and accurately display the recognition effect; the LK optical flow method and the Bhattacharyya coefficient method are combined to predict the vertical direction in the strong storm area, focus on the overall dynamic analysis of the radar echo area, and finely identify and effectively predict potential strong wind activities; feature parameters are selected and predictions are performed using SVM, CNN, and U-Net, and weighted voting is performed on the corresponding prediction results to broaden the prediction dimension, enhance the recognition performance, overcome the limitations of a single model, improve the robustness of the model, and improve the accuracy of storm prediction in the horizontal direction in the strong storm area. BRIEF DESCRIPTION OF THE DRAWINGS

[0013] Figure 1 is a method flow chart;

[0014] Figure 2 It is a schematic diagram of bilinear interpolation method;

[0015] FIG3 is a schematic diagram of removing isolated point interference echoes;

[0016] FIG4 is a schematic diagram of removing the banner-shaped interference echo;

[0017] Figure 5 Identify renderings for severe storm areas;

[0018] Figure 6 This is the core reflectivity factor profile of the downburst storm cell in Experiment 1;

[0019] Figure 7 Schematic diagram of the SVM prediction results of Experiment 1;

[0020] Figure 8 Schematic diagram of CNN prediction results of Experiment 1;

[0021] Fig. 9 This is a schematic diagram of the U-Net prediction results of Experiment 1;

[0022] Fig.10 This is a schematic diagram of the horizontal prediction results of Experiment 1;

[0023] Fig.11 This is the core reflectivity factor profile of the downburst storm cell in Experiment 2;

[0024] Fig.12 Schematic diagram of the SVM prediction results of Experiment 2;

[0025] Fig.13 This is a schematic diagram of the CNN prediction results of Experiment 2;

[0026] Fig.14 This is a schematic diagram of the U-Net prediction results of Experiment 2;

[0027] Fig.15 Schematic diagram of the horizontal prediction results of Experiment 2. DETAILED DESCRIPTION

[0028] The specific implementation modes of the present invention are described below so that those skilled in the art can understand the present invention. However, it should be clear that the present invention is not limited to the scope of the specific implementation modes. For those of ordinary skill in the art, as long as various changes are within the spirit and scope of the present invention as defined and determined by the attached claims, these changes are obvious, and all inventions and creations utilizing the concept of the present invention are protected.

[0029] like Figure 1 As shown, a method for nowcasting convective storm gale includes the following steps:

[0030] S1. Acquire multi-source observation data of convective storm gale and perform preprocessing to obtain preprocessed multi-source observation data;

[0031] The multi-source observation data in step S1 includes weather radar data, automatic station data, wind profiler radar data and sounding radar data of multiple time periods;

[0032] Weather radar data includes the differential reflectivity factor Z DR , correlation coefficient ρ HV , differential phase shift φ DP , differential phase shift rate K DP , depolarization ratio LDR, echo top height, echo centroid height, vertically integrated liquid water content, vertically integrated liquid water content density, lower layer radial velocity, differential reflectivity and differential reflectivity column;

[0033] Automatic station data include air temperature, pressure, rainfall, wind speed, wind direction, and relative humidity;

[0034] Wind profile radar data includes overhead wind direction and wind speed;

[0035] Sounding radar data include temperature, pressure, humidity, and wind at all levels within the vertical range of 0 to 30 km in the atmosphere;

[0036] Among them, weather radar data is used to obtain initial horizontal prediction results and initial vertical prediction results; wind profiler radar data and sounding radar data are used to assist in analyzing the environment corresponding to the initial horizontal prediction results and initial vertical prediction results, and to analyze the environmental conditions suitable for the occurrence of thunderstorms and gale winds; automatic station data are used to verify the analysis results and prediction results, and further determine the final area where the storm will occur.

[0037] The present invention uses an S-band dual linear polarization radar, which has a scanning time of completing a volume scan every 6 minutes, an operation mode of 24-hour continuous observation, and operates in VCP21 mode, with a time resolution of 6 minutes and a distance resolution of 250m. The parameters of the S-band dual linear polarization radar are shown in Table 1.

[0038] Table 1 SA dual linear polarization weather radar parameters

[0039]

[0040]

[0041] The automatic station data are obtained through the automatic weather station in the area to be measured, with an observation interval of 10 minutes and an intensive observation of 5 minutes.

[0042] Wind profiler radar data is obtained through wind profiler radar. Wind profiler radar uses the Doppler effect to detect the changes of meteorological elements such as wind direction and wind speed above it with altitude, and can provide three-dimensional wind field information. It has the advantages of high detection time and space resolution and high degree of automation. It can continuously and in real time monitor the atmospheric conditions from 100 meters above the ground surface to 6,000 meters or even higher with high time and space resolution, measure the changes of horizontal wind speed and wind direction at different altitudes within this range, and can also capture the distribution characteristics of key meteorological parameters such as atmospheric refractive index structure constant in detail along the vertical direction. Wind profiler radar observes 12 times per hour, and the operation mode is continuous observation 24 hours a day.

[0043] Sounding radar data is obtained through sounding radar. Sounding radar uses sounding balloons carrying radio transponders to observe and can obtain meteorological element data such as temperature, pressure, humidity, wind, etc. at all levels within the vertical range of 0 to 30 km.

[0044] The preprocessing in step S1 includes the following steps:

[0045] S1-1. Gridding the weather radar data using a bilinear interpolation algorithm to obtain gridded weather radar data;

[0046] Construct a three-dimensional coordinate grid system, with the radar site at the grid center and the latitude and longitude grid spacing adjusted according to the resolution requirements; the grid spacing is generally set to 1 km, 2 km or 5 km, which is converted into longitude and latitude at intervals of 0.01°, 0.02° and 0.05°, and for vertical height, a series of unequally spaced height levels are selected, such as 500 meters, 1 km, 1.5 km, 2 km, 3 km, 4 km, 5 km, 7 km, 10 km, 15 km, etc.

[0047] Apply the bilinear interpolation algorithm to each node in the latitude, longitude and altitude grid to fill in the value. The corresponding process is:

[0048] Based on the calculation rules of the earth coordinate system, the azimuth and horizontal distance of each grid point relative to the radar station are determined. Combined with the actual altitude data, the radial distance between the corresponding grid point and the radar and its corresponding elevation angle are obtained through trigonometric calculations. The corresponding formula is:

[0049] dlon=lon2-lon1

[0050] dlat=lat2-lat1

[0051]

[0052]

[0053] Where (lat1, lon1, H1) and (lat2, lon3, H2) represent the coordinates of two grid points, dlon and dlat represent the latitude and longitude differences of the two grid points, dis represents the horizontal distance, R represents the radius of the earth, arcsin(·), sin(·), cos(·), and arctan(·) represent the inverse sine function, sine function, cosine function, and inverse tangent function, δ represents the angle, θ represents the azimuth, α represents the elevation, Rd represents the radial distance, H represents the altitude, including H1 and H2, and π represents pi;

[0054] Based on the radial distance between each grid point and the radar and its corresponding elevation angle, the bilinear interpolation technique is used to accurately estimate and assign the corresponding radar reflectivity data to each grid point. When the longitude and latitude coordinates of the radar are known, the relative position of any grid point to the radar can be obtained. The coordinates projected to the radar coordinate system are (e, a, r), where e is the elevation angle, a is the azimuth angle, and r is the radial distance (radial bin number). Therefore, two horizontal and vertical lines can be drawn from this point, intersecting with the two adjacent elevation angles to obtain four coordinate points, such as Figure 2 As shown, the grid point value is f a (e,a,r), and the corresponding observation values ​​at each point where the two elevation angles intersect are f a (e1,a,r),f a (e1,a,r1),f a (e2,a,r),f a (e2,a,r2), the corresponding formula is:

[0055]

[0056] in,

[0057] S1-2, using the improved interference echo automatic identification algorithm to remove electromagnetic interference echoes from the gridded weather radar data, to obtain weather radar data after the electromagnetic interference echoes are removed;

[0058] The specific process of the improved interference echo automatic identification algorithm in step S1-2 is as follows:

[0059] S1-2-1. In order to process the non-fixed radial number of each PPI scanning surface and the initial scanning azimuth angle in the gridded weather radar data, it is necessary to determine whether the radial number of each PPI scanning surface of the gridded weather radar data is 360; if so, the PPI scanning surface is retained; otherwise, the PPI scanning surface is removed to obtain the determined weather radar data; wherein, each PPI scanning surface of the gridded weather radar data includes the corresponding radial and elevation angles;

[0060] S1-2-2. Because there are isolated point interference echoes in the gridded weather radar data, and the isolated point interference echoes appear as isolated points or thin lines in the figure, the radar echo image corresponding to the judged weather radar data can be drawn first, and an area range can be set to calculate the number of valid numerical points around the center point within the area range; the area range is generally set to 5*5;

[0061] S1-2-3. Determine whether the percentage of valid numerical points around the center point within the area exceeds the threshold; if so, the corresponding center point is taken as a valid value; otherwise, the corresponding center point is taken as an invalid value and removed; and the filtered weather radar data is obtained; wherein the threshold value is 55%, and the corresponding formula is:

[0062]

[0063] Where ∑(·) represents the summation function, i represents the orientation, j represents the distance library, and R z represents the total number of effective observations of differential reflectivity factors in a 5×5 area centered on (i, j), P i,j Represents the total number of libraries and R z The percentage of the total number of windows, that is, the percentage of valid numerical points around the center point in the area. i,j When it is less than 55%, the pixel corresponding to (i, j) is considered as an interference echo and is removed. In order to balance the purpose of retaining effective information and filtering out clutter, it is necessary to carefully select the regional threshold to avoid excessive removal of isolated interference echoes and cause the loss of effective meteorological information. At the same time, the effect of removing interference echoes should not be reduced in order to retain more effective information. Therefore, multiple tests are conducted every 5% between 25% and 75% to determine the optimal threshold. The results of multiple tests show that setting the threshold to 55% is a more appropriate choice, which can achieve a better balance between retaining effective information and removing clutter.

[0064] As shown in FIG3 , it can be seen that the isolated point echo is removed more cleanly after the electromagnetic interference removal operation.

[0065] Remove banner-shaped interference echoes. Banner-shaped interference echoes often appear as thin or thick strips from the center to the edge in the radar echo image. Such interference may appear in multiple directions, so multiple banner-shaped interference echoes may appear in the same image. Compared with isolated point interference, removing banner-shaped interference is more challenging, so it is necessary to process steps S1-2-4 to S1-2-11 until all banner-shaped interference echoes are removed.

[0066] S1-2-4. Count the total library number and the number of continuous echo libraries in each radial direction of the filtered weather radar data, and record the maximum and minimum values ​​of the total library number;

[0067] S1-2-5, the total library number, the maximum value, the minimum value and the continuous echo library number of each radial direction, and determine the corresponding radial direction with the maximum total library number of valid echoes;

[0068] S1-2-6, according to the formula:

[0069]

[0070] Get the radial echo continuity index R corresponding to each radial direction t ; Among them, R max Indicates the maximum total number of libraries, R min Indicates the minimum value of the total library number;

[0071] S1-2-7, Based on the judgment criteria:

[0072]

[0073] Determine whether each radial direction is an interference radial direction; if so, determine the corresponding radial direction as the initial interference radial direction and enter step S1-2-8; otherwise, determine the corresponding radial direction as the initial normal radial direction and enter step S1-2-10;

[0074] S1-2-8, according to the formula:

[0075] R d =|R N -R i |,i=N±20

[0076] Get the absolute value R of the total distance library of effective echoes on the initial interference radial direction and its nearby radial directions d ; Among them, R N Indicates the total distance bin number of effective echoes on the initial interference radial direction, R i represents the total distance bin number of effective echoes in the radial direction near the initial interference radial direction, and |·| represents the absolute value;

[0077] S1-2-9. Determine the absolute value R of the total distance library of effective echoes d Is it greater than 20? If so, the initial interference radial direction is determined to be the initial normal radial direction; otherwise, the initial interference radial direction is determined to be the final interference radial direction and is removed; proceed to step S1-2-10;

[0078] S1-2-10, determine the starting radial direction N1 and the ending radial direction N2 of each initial normal radial direction N, and perform edge effective echo detection on the area [N1-1, N2+1]. If there is no effective echo, remove the interference echo in this area;

[0079] S1-2-11, repeat step S1-2-10 until a valid echo is obtained, and use the corresponding initial normal radial direction as the normal radial direction;

[0080] The weather radar data after removing the electromagnetic interference echo includes the filtered weather radar data corresponding to the normal radial direction.

[0081] As shown in FIG. 4 , it can be seen that the strip-shaped electromagnetic clutter has also been removed after the electromagnetic interference removal operation. It can be seen that the improved interference echo automatic identification algorithm in step S1-2 performs well in removing isolated point interference echoes, and can also better deal with strip-shaped interference echoes, showing a satisfactory effect.

[0082] S1-3. When using the Doppler radar velocity field, we often encounter the problem of velocity folding (ambiguity). This phenomenon is caused by the need to reduce the maximum detection speed in order to ensure a larger detection distance under the same radar and the same pulse repetition frequency. When the phase shift detected between consecutive radar pulses exceeds the maximum phase shift that can be clearly detected by the radar pulse repetition frequency and wavelength, velocity ambiguity will occur. Therefore, it is necessary to use the velocity deambiguation algorithm to process the weather radar data after removing the electromagnetic interference echo to obtain the weather radar data after velocity deambiguation. The corresponding specific process is:

[0083] First, the initial reference radial is searched. By traversing all radials in a specific elevation layer in the weather radar data after removing the electromagnetic interference echo, a "good" radial with smaller wind shear is selected as the initial reference. Among them, a good radial positioning has no large wind shear between any adjacent distance libraries, that is: |V i -V i-1 |<αV N ; Among them, V N represents the Nyquist speed, α represents an adjustable parameter, |·| represents the absolute value, V i 、V i-1 They represent the i-th radial velocity and the i-1-th radial velocity respectively; where α is taken as 0.8.

[0084] In the screening process, according to the characteristic that there will be no velocity ambiguity near the area where the radar radial velocity is zero, the weak wind shear area is determined by using the characteristic that the radial velocity is zero when the true wind direction is perpendicular to the radial wind direction. In the weak wind shear area, three appropriate adjacent radials are selected as the initial reference radials according to the radial judgment standard. Among them, the radial judgment standard is: in each good radial, according to the formula:

[0085]

[0086] Calculate the average velocity V of the radial distance library along the radial direction M1; Where N1 indicates that the radial velocity value is within the interval ±βV N The total number of reservoirs within, β represents an adjustable parameter, and the corresponding value is 0.3; if the average velocity of the distance reservoir V M1 When the sign changes from two consecutive radial directions to the next two consecutive radial directions, the one with the larger total number of libraries N1 in the middle two radial directions is selected as the initial radial direction.

[0087] If the initial radial cannot be found through the above process, the second iterative search is performed. In the second iterative search, a good radial can be determined if the three radial conditions are met. The three radial conditions are: ① There is no wind shear greater than the wind shear threshold between adjacent measurement point pairs in the radial; ② The number of non-missing velocity measurement points with a total library number N0 of not less than 40 is greater than the set non-missing velocity measurement point number threshold; ③ The average value V of all non-missing velocity measurement points M0 are all less than the average threshold, that is, Among them, the thresholds involved in the three radial conditions are all set according to actual conditions.

[0088] If multiple good radials are found, the average value V is selected M0 The smallest radial direction is used as the initial radial direction. If the initial radial direction cannot be found, relax radial condition ② (total library number N0 is not less than 35) and start the next search. Repeat the above search process until the initial radial direction is found or the total library number N0 changes to less than 5. If the measurement speed N0 is greater than 5 and V M0 <βV N If the image has a good radial direction, no deblurring is performed.

[0089] Each distance library is marked with four different flags to record its processing process. The values ​​of these flags range from 0 to 3, representing: 0 for unprocessed points; 1 for processed points with high confidence; 2 for processed points with low confidence; 3 for points with velocity collapse or missing values. Different flags are handled differently during the processing and error checking stages. Once a distance library is marked as "1", it is considered to be a correctly processed (i.e. good) library and no further processing is required to improve the overall algorithm efficiency and save computing time.

[0090] Secondly, the three initial reference radials are used to deblur the velocity. The velocity averages of these radials are calculated and used as reference velocities for denoising. At the same time, based on the azimuth continuity constraint, the adjacent radials are processed to ensure the continuity and accuracy of the algorithm. That is, after finding a suitable reference radial, the algorithm first uses the average velocity of the distance library V M1 (or use the average value V in specific cases M0As the reference velocity, all points in the initial radial are denoised. Taking the initial radial as the benchmark, the two adjacent radials are further processed. In this process, the algorithm makes full use of the azimuth continuity constraint, finely processes the velocity data of each distance library, and compares it with the library of the adjacent radial in the same range. If there is an available reference velocity in the initial radial, the algorithm will process the corresponding library in the adjacent radial accordingly and assign a mark of "1". If this condition is not met, the library in the adjacent radial will be marked as "0" (i.e., unprocessed state) and left for subsequent processing. The initial radial and its two adjacent radials (a total of three) together constitute the reference data for the first round of radial-by-radial denoising. Since this process uses the continuity constraint in the azimuth direction, only one radial is required as a reference, thereby reducing the potential impact of random errors in the reference velocity.

[0091] Next, the first round of velocity deblurring is carried out radially by radial. Each radial is traversed in clockwise and counterclockwise directions, and the velocity of each distance library is compared with the reference velocity to determine whether deblurring is required. In addition, strict conditions are set to identify the initial distance library to ensure the accuracy and stability of the deblurring process.

[0092] All distance bins within 180° are traversed in clockwise and counterclockwise directions respectively. During the traversal process, the algorithm will compare the speed data of each distance bin with the reference speed V R The reference speed VR is determined based on the average value of the radial speed (clockwise or counterclockwise) in the first three same ranges. R The validity of V must meet the following conditions: the processing flags of the first three distance libraries must all be 1, and any two adjacent speed pairs must satisfy a specific mathematical relationship. R The condition is not met or the difference from the given speed exceeds the threshold αV N , then the distance library will be marked as unprocessed (0).

[0093] The algorithm then turns to the deambiguation process in the radial direction. Prior to this, the initial distance library needs to be identified. The identification of the initial distance library is based on a series of strict conditions: in the radial direction, there are at least 10 kilometers of continuous non-missing adjacent distance libraries on each side; the processing flag of these distance libraries should be 1, and each side contains at least two adjacent distance libraries with flag 1; at the same time, in the first few radials in the same range, there should be at least three adjacent distance libraries with flag 1. Starting from the outermost distance library in the radial direction, the initial distance library that meets the conditions is searched inward step by step. Once the initial distance library is found, the distance library is deambiguated one by one along the radial direction. At this time, the average of the first three good speeds on the same radial direction is used as the reference speed. If data interruption is encountered during the processing, such as a distance library with missing values ​​or velocity ambiguity, the new initial distance library is searched and determined again, and then the search and deambiguation process are repeated until the end of the radial is traversed.

[0094] Finally, the algorithm assigns different flags to each distance library during the processing to record its processing status.

[0095] Deblurring is further performed on all previously unprocessed distance bins. In the currently investigated radial (denoted as radial 0), the first three adjacent radials (labeled as "-1", "-2", and "-3") are first identified, and three adjacent distance bins in these three radials are selected. If the following conditions are met, then according to the formula:

[0096] |V i -V i-1 |<αV N ,i=-1,-2

[0097]

[0098] Defuzzify the distance library corresponding to the observed speed V0: If the above conditions are not met, the algorithm continues to check one radial direction forward and re-check the three speed distance libraries in the same range as the observed speed V0. This process is repeated M times (M = 3 by default) until a suitable VR is found. If no reference speed V that meets the conditions is found after all steps are completed, R , then stop the deblurring process for the distance library. In addition, in order to optimize the processing effect and reduce the influence of random noise, when searching for unprocessed distance libraries along the radial direction, the search is preferentially started from the middle distance starting point. For each unprocessed distance library, search the nearest three consecutive velocity-deblurred libraries along the radial direction until a library that meets the requirements is found or the set range limit is exceeded.

[0099] like Figure 5As shown, by using this algorithm, according to the threshold and echo area limit, the strong storm area is well identified. It can be seen that the strong storm area is mainly located in Figure 5 The southern and central parts of the central region, but the strongest area is still in the south, which can effectively improve the recognition of thunderstorm winds in the vertical and horizontal directions. It can be seen that the de-velocity fuzzification processing in step S1-3 can significantly improve the accuracy and reliability of weather radar data.

[0100] The preprocessed multi-source observation data include automatic station data, wind profiler radar data, sounding radar data and weather radar data after velocity blurring.

[0101] S2, using the region growing algorithm to identify the region of the preprocessed multi-source observation data to obtain the strong storm area;

[0102] Step S2 includes the following steps:

[0103] S2-1, sampling the weather radar data after de-blurring the velocity to obtain spatially continuous time series data;

[0104] S2-2, smoothing the spatially continuous time series data by using a 3×3 Gaussian filter kernel to obtain smoothed time series data;

[0105] S2-3, set threshold;

[0106] S2-4, constructing a data grid corresponding to the smoothed time series data and obtaining corresponding image size information;

[0107] S2-5, traverse the data grid and image size information to construct a feature map array of the same dimension;

[0108] S2-6, initialize the seed point queue structure;

[0109] S2-7, scan the data grid in order from left to right and from top to bottom, and determine whether the intensity in each grid point exceeds the threshold set in step S2-3; if so, determine the corresponding grid point as a strong echo point and add it to the seed point queue structure; otherwise, determine the corresponding grid point as a weak echo point;

[0110] S2-8, using the queue rule, gradually add pixels that are adjacent to the seed point and meet the conditions to the current area until no new pixels can be added, and obtain the corresponding strong echo area; if there are still pixels that have not been assigned to any area, select a new seed point and repeat the above steps until all pixels are assigned to a certain area;

[0111] S2-9, counting the total number of pixels corresponding to each strong echo region to obtain the area of ​​each strong echo region;

[0112] S2-10, setting an area threshold; judging whether the area of ​​each strong echo region is greater than the area threshold; if so, retaining the corresponding strong echo region; otherwise, removing the corresponding strong echo region;

[0113] S2-11. Using the feature map array, store the strong echo area retained in step 2-10, that is, the strong storm area.

[0114] S3. Use the LK optical flow method and the Bhattacharyya coefficient method to perform storm prediction on the preprocessed multi-source observation data in the vertical direction in the strong storm area to obtain corresponding vertical prediction results; wherein the LK optical flow method is the LK pyramid optical flow method.

[0115] LK Pyramid Optical Flow uses optical flow analysis technology to conduct a detailed study of the optical flow field of the vertical profile of the reflectivity factor, clarifying the law of the storm core changing over time. The Lagrangian mechanics model is used to accurately fit the function of the storm core top height drop process, and the histogram and Bhattacharyya coefficient statistical analysis method is used to effectively match and identify the "positive and negative velocity pairs" images in the middle-layer radial velocity field of the storm core. Combined with a series of threshold judgments, the intelligent prediction of downbursts is successfully achieved, that is, step S3 includes the following steps:

[0116] S3-1. Search for the differential reflectivity factor Z in each severe storm area DR Maximum value;

[0117] S3-2. Determine each differential reflectivity factor Z DR Whether the maximum value is greater than the threshold TH1 and the area of ​​each strong storm area is greater than TH2; if so, enter step S3-3; otherwise, end the operation;

[0118] S3-3, select the corresponding de-velocity blurred weather radar data at a time and a time before and after it, and convert the differential reflectivity factor Z DR The maximum value is used as the section, the corresponding section array is generated, and the section image is drawn;

[0119] S3-4, using the LK optical flow method to analyze the profile arrays at three moments, and obtain the corresponding optical flow fields;

[0120] Step S3-4 includes the following steps:

[0121] S3-4-1, converting the dimensions of the profile arrays at three moments to obtain corresponding converted profile arrays;

[0122] S3-4-2, respectively perform LK optical flow calculation on the converted profile arrays at two adjacent moments, and obtain the initial optical flow vector field between the two adjacent moments;

[0123] The optical flow window size of the LK optical flow calculation in step S3-4-2 is 25×25, the maximum number of pyramid layers is set to 4 layers, the maximum number of iterations is 10 times, and the error threshold is 0.03.

[0124] S3-4-3, decomposing the velocity vectors of the two initial optical flow vector fields into two components parallel to the x-axis and the y-axis; determining the median of the components in two directions of each initial optical flow vector field, and merging the two medians to obtain the corresponding initial velocity vector;

[0125] S3-4-4, combine the two initial velocity vectors to obtain the final velocity vector;

[0126] S3-4-5, using the final velocity vector to perform extrapolation movement on the cross-section data corresponding to a moment selected in step S3-3, to obtain a corresponding moved image;

[0127] S3-4-6, calculating the initial optical flow field of the original image composed of the moved image and the cross-sectional data corresponding to a moment selected in step S3-3;

[0128] S3-4-7. Perform matrix addition on the initial optical flow field and the corresponding two initial optical flow vector fields in S3-4-2 to obtain the corresponding intermediate optical flow field; perform weighted addition on the two intermediate optical flow fields to obtain the optical flow field.

[0129] S3-5, calculate the storm core drop rate of the light process; determine whether the storm core drop rate is greater than the threshold TH3; if so, proceed to step S3-6; otherwise, terminate the operation;

[0130] S3-6, read the radar radial velocity low elevation angle volume scan data corresponding to the radar data, and use the value TH4 as the radius range, count the positive and negative velocity area distribution within the radius range, and obtain the corresponding histogram;

[0131] S3-7, converting the histogram into probability density distribution of different speeds in the positive and negative speed areas, and drawing the corresponding probability density distribution histogram; the probability density distribution histogram uses the threshold TH5 as the vertical coordinate accuracy and the threshold TH6 as the horizontal coordinate.

[0132] S3-8, according to the Bhattacharyya coefficient formula:

[0133]

[0134] Get the similarity ρ(ρ1,ρ2); where ρ1 represents the percentage of the speed value of each region of positive speed in the total, ρ2 represents the percentage of the speed value of each region of negative speed in the total, Σ(·) represents the summation function, TH6 represents the total number of regions of positive and negative speeds, ρ1(i) and ρ2(i) represent the percentage of the speed value of each region of the i-th positive speed in the total, and the percentage of the speed value of each region of the i-th negative speed in the total, respectively;

[0135] S3-9, determine whether the similarity ρ(ρ1, ρ2) is greater than the threshold TH7; if so, determine that a downburst event is very likely to occur in the strong storm area and enter step S3-10; otherwise, terminate the operation;

[0136] S3-10. Combined with the automatic station data, wind profiler radar data, and sounding radar data, the optical flow field is calculated using the linear extrapolation method, and the Lagrangian mechanics model is used for extrapolation to obtain the downburst storm monomer core reflectivity factor profile, that is, the vertical prediction result. Among them, the environmental conditions corresponding to the initial profile of the downburst storm monomer core reflectivity factor are analyzed according to the wind profiler radar data and the sounding radar data, and whether the environmental conditions meet the environment for the occurrence of thunderstorms and gale is analyzed. Based on the analysis results, the thunderstorm gale occurrence area in the initial profile of the downburst storm monomer core reflectivity factor is confirmed; then the confirmed thunderstorm gale occurrence area is verified according to the automatic station data to ensure that the thunderstorm gale occurrence area meets the thunderstorm gale occurrence conditions, and finally the downburst storm monomer core reflectivity factor profile is obtained, that is, the vertical prediction result. The analysis of the downburst storm monomer core reflectivity factor profile using the automatic station data, wind profiler radar data, and sounding radar data can further improve the accuracy of storm prediction in the vertical direction.

[0137] The threshold and value ranges designed in step S3 are shown in Table 2.

[0138] Table 2 Threshold reference range

[0139] Threshold Name Threshold Range <![CDATA[TH1]]> 45-55(dBZ) <![CDATA[TH2]]> <![CDATA[50-200(km 2 )]]> <![CDATA[TH3]]> 0.6-2.5 (km / 6min) <![CDATA[TH4]]> 20-50 <![CDATA[TH5]]> 0.1-1 <![CDATA[TH6]]> 20-33 <![CDATA[TH7]]> 0.1-0.9

[0140] S4, constructing a horizontal storm prediction model; using the horizontal storm prediction model; performing storm prediction on the preprocessed multi-source observation data in the horizontal direction of the strong storm area to obtain corresponding horizontal prediction results;

[0141] The horizontal storm prediction model in step S4 includes SVM support vector machine, CNN convolutional neural network, and U-Net sub-model.

[0142] Step S4 includes the following steps:

[0143] S4-1. Determine characteristic data based on the weather radar data after velocity blurring; wherein the characteristic data includes the maximum differential reflectivity factor, echo top height, maximum differential phase shift rate, echo centroid height, vertically integrated liquid water content, vertically integrated liquid water content density, low-layer radial velocity, maximum differential reflectivity and differential reflectivity column;

[0144] S4-2. Input the feature data into the SVM support vector machine to obtain the corresponding SVM prediction result; wherein, the kernel function of the SVM support vector adopts the Gaussian kernel function, which can map the data to an infinite-dimensional feature space, making it easier to model nonlinear problems in high-dimensional space.

[0145] S4-3. Input the feature data into the CNN convolutional neural network to obtain the corresponding CNN prediction results; wherein, the CNN convolutional neural network includes a convolutional layer, a maximum pooling layer, a flattening layer, two fully connected layers and an output layer connected in series; the convolutional layer is composed of multiple convolutional units and adopts the ReLU function, and the parameters of the convolutional units are optimized by the back propagation algorithm.

[0146] S4-4. Input the feature data into the U-Net sub-model to obtain the corresponding U-Net prediction results; the U-Net sub-model includes an encoder (downsampling path) and a decoder (upsampling path); the encoder contains a convolution layer to extract features, and a maximum pooling operation to reduce the data dimension. The decoder uses a transposed convolution layer to upsample the feature map, and connects the feature maps of the corresponding layers of the encoder and decoder through a jump connection to restore the detailed information of the image. Finally, the feature map is mapped to the required number of categories through a 1x1 convolution kernel and output using an appropriate activation function.

[0147] S4-5, using weighted voting to process the SVM prediction results, CNN prediction results, and U-Net prediction results to obtain the final prediction results, that is, to obtain the initial level prediction results;

[0148] The ensemble model of SVM support vector machine, CNN convolutional neural network and U-Net sub-model fusion can broaden the prediction dimension and enhance the recognition performance, overcome the limitations of a single model, improve the robustness of the model, and improve the accuracy of storm prediction in the horizontal direction in strong storm areas.

[0149] S4-6. Analyze the initial horizontal prediction results based on the automatic station data, wind profiler radar data, and sounding radar data to obtain horizontal prediction results. That is, analyze the environmental conditions corresponding to the initial horizontal prediction results based on the wind profiler radar data and sounding radar data, analyze whether the environmental conditions meet the environment for thunderstorms and gale, and confirm the thunderstorm and gale occurrence area in the initial horizontal prediction results based on the analysis results; then verify the confirmed thunderstorm and gale occurrence area based on the automatic station data to ensure that the thunderstorm and gale occurrence area meets the thunderstorm and gale occurrence conditions, and finally obtain the horizontal prediction results. Analyze the initial vertical prediction results using the automatic station data, wind profiler radar data, and sounding radar data to further improve the accuracy of storm prediction in the vertical direction.

[0150] S5. Based on the vertical prediction results and the horizontal prediction results, the final prediction results are determined, that is, the nowcast results of the convective storm gale are obtained.

[0151] The final prediction result in step S5 is the overlapping area of ​​the vertical prediction result and the horizontal prediction result.

[0152] In the embodiment of the present invention, the thunderstorm gale weather process case data is selected as radar data of the occurrence, development and extinction stages of the thunderstorm gale; the hit rate POD, false alarm rate FAR and critical success index SCI are selected as test indicators, wherein:

[0153] The formula corresponding to the hit rate POD is:

[0154]

[0155] Among them, TP represents true positive examples and FN represents false negative examples;

[0156] The formula corresponding to the false alarm rate FA is:

[0157]

[0158] Among them, FP means false positive;

[0159] The formula corresponding to the critical success index SCI is:

[0160]

[0161] The accuracy, precision, recall and F1 Score are selected as the measures of the matching degree between the model prediction results and the real labels.

[0162] Experiment 1: The individual data of thunderstorm and gale weather process in region 1 were selected as weather radar data. The automatic station data were hourly data of the same period as the S-band dual linear polarization radar. The wind profile radar data were hourly data of the same period as the S-band dual linear polarization radar. The sounding radar data were data of a similar period as the S-band dual linear polarization radar.

[0163] According to the probability density data in the probability density distribution histogram of area 1 and the calculation formula of the Bhattacharyya coefficient method, the similarity ρ(ρ1, ρ2) is calculated to be 0.9348, which meets the condition of being greater than the threshold TH7. Therefore, it can be determined that there is a high probability that a downburst event has occurred in this area at this moment. At the same time, combined with the radar data and the automatic station data, the maximum wind speed at this time is 15.5m / s, indicating that the result is correct.

[0164] Then the optical flow field is calculated using the linear extrapolation method, with ref_1446 as the reference and V as the moving vector, and the vertical profile is extrapolated to the next 12 minutes, as shown in Figure 6 As shown, it can be seen that the flow field of the storm is mostly moving downward. At the same time, it can be observed that the shape of the strong echo center and the moving trend in the outward figure are roughly similar to the corresponding areas, covering areas 1, 2, and 3 in area 1.

[0165] like Figure 7 As shown, the dark part is the area where SVM predicts that thunderstorms and gale may occur, which is roughly composed of a bow echo and a block echo, and has a good correspondence with the strong storm area. It can be seen that the main thunderstorm and gale areas are area 1, area 2, and area 3 in area 1.

[0166] like Figure 8 As shown in the figure, the darker areas indicate that the probability of thunderstorms predicted by CNN is higher, and vice versa, the lighter areas indicate that the probability of thunderstorms predicted by CNN is lower. From the identified probability map, it can be seen that the areas most prone to thunderstorms are area 1 and area 3 of region 1, and the darker areas of these two parts are the most obvious.

[0167] like Fig. 9 As shown, darker colors indicate that U-Net predicts a higher probability of thunderstorms, and lighter colors indicate that U-Net predicts a lower probability of thunderstorms. The most likely location for thunderstorms is region 1 in region 1, which has the darkest color and the largest area.

[0168] like Fig.10 As shown, the CNN prediction results, U-Net prediction results, and CNN prediction results are weighted voted, and the horizontal prediction results show that the locations where thunderstorms and gale are most likely to occur are area 1, area 2, and area 3 in area 1.

[0169] Combining the vertical prediction results and horizontal prediction results of Region 1, the locations where thunderstorms and gale are most likely to occur in Region 1 are Area 1, Area 2, and Area 3 of Region 1.

[0170] Experiment 2: The individual data of thunderstorm and gale weather process in Region 2 were selected as weather radar data. The automatic station data were minute-level data in the same period as the S-band dual-linear polarization radar. The wind profile radar data were hour-level data in the same period as the S-band dual-linear polarization radar. The sounding radar data were data in a similar period as the S-band dual-linear polarization radar.

[0171] According to the probability density data in the probability density distribution histogram of area 2 and the calculation formula of the Bhattacharyya coefficient method, the similarity ρ(ρ1,ρ2) is calculated to be 0.9439, which meets the condition of being greater than the threshold TH7. Therefore, it can be determined that there is a high probability that a downburst event has occurred in this area at this moment. At the same time, combined with the maximum wind speed at this time from the radar data and the automatic station data, the result shows that it is correct.

[0172] Then the optical flow field is calculated using the linear extrapolation method, with ref_1854 as the reference and V as the moving vector, and the vertical profile is extrapolated to the next 12 minutes, as shown in Fig.11 As shown, it can be seen from the flow field diagram of the storm that the storm is moving downward. At the same time, it can be observed that the shape of the center of the strong echo and the moving trend in the outward diagram are partially similar to the corresponding areas, covering areas 4 and 5 of area 2.

[0173] like Fig.12 As shown, the dark part is the area where SVM predicts that thunderstorms and gale may occur, which has a good correspondence with the strong storm area. It can be seen that the main thunderstorm and gale areas are area 4 and area 5 in area 2.

[0174] like Fig.13 As shown in the figure, the darker areas indicate that the probability of thunderstorms predicted by CNN is higher, and vice versa, the lighter areas mean that the probability of thunderstorms predicted by CNN is lower. The thunderstorm area is also located in area 4 of area 2, but there are some gaps in the middle of the CNN recognition area. This may be because this area is close to the radar. Although interpolation is performed, some data may still be unreasonable. However, the overall echo is similar to that of the strong storm area.

[0175] like Fig.14 As shown, darker colors indicate that U-Net predicts a higher probability of thunderstorms and gale, and lighter colors indicate that U-Net predicts a lower probability of thunderstorms and gale. The U-Net prediction area is closer to a whole and more complete. The location where thunderstorms and gale are most likely to occur is area 4 in area 2.

[0176] like Fig.15 As shown in the figure, the CNN prediction results, U-Net prediction results, and CNN prediction results are weighted voted. The area where thunderstorms and gale may occur is still in area 4 of area 2, but it can be clearly seen that a strip echo extends to area 6 of area 2. Therefore, thunderstorms and gale may occur in both area 4 and area 6 of area 2.

[0177] Combining the vertical prediction results and horizontal prediction results of Region 2, the locations where thunderstorms and gale are most likely to occur in Region 2 are Area 4 and Area 6 of Region 2.

[0178] As shown in Table 3, the hit rate (POD) of region 1 is as high as 0.862, indicating that the model can effectively capture the actual thunderstorm and wind events; the precision rate (Precision) reaches 0.937, indicating that the model has a high accuracy in identifying thunderstorm and wind events; the recall rate (Recall) is 0.892, indicating that the model can better capture the real thunderstorm and wind events; the accuracy rate (ACC) is 0.926, indicating that the model has a high accuracy rate in the overall sample classification. In terms of the critical success index (CSI), it is 0.799, which also shows that there is a certain degree of accuracy in predicting thunderstorm and wind events. The false alarm rate (FAR) for region 2 is 0.0403, which means that 4.03% of the predictions were incorrectly labeled as high winds; the hit rate (POD) is 0.833, indicating that the model successfully detected 83.3% of the high wind events; the accuracy rate (ACC) is 0.884, which shows that the model has high accuracy overall; the precision rate (Precision) is 0.960, indicating that when the model predicts high winds, it has an accuracy rate of 96.0%; the recall rate (Recall) is also 0.833, which is the same as the detection rate, indicating the ratio of high wind events that are successfully detected; the F1 score is 0.863, which is the weighted average of the precision rate and the recall rate, which comprehensively considers the performance of the model; the critical success index (CSI) is 0.791, indicating that the model has a relatively high ability to successfully predict.

[0179] Table 3 Vertical prediction result index data

[0180] Indicator data Region 1 Region 2 Hit Rate 0.86225 0.83298 False alarm rate 0.02370 0.04031 Critical Success Index 0.79942 0.79073 Accuracy 0.91628 0.88401 Accuracy 0.89657 0.90677 Recall 0.89221 0.83298 F1 value 0.86877 0.86282

[0181] As shown in Table 4, for region 1, SVM, CNN, U-Net, and ensemble models all performed well, and the corresponding index data were all higher than 0.75, and except for the accuracy and F1 value of SVM, which were lower than 0.8, the other index data were all higher than 0.8. For region 2, SVM, CNN, U-Net, and ensemble models all performed well, except for the accuracy of SVM, which was 0.7914, which was lower than 0.8, the other index data were all higher than 0.8. It can be seen that the prediction results of the two embodiments of the present invention are close to the true value, which proves the superiority of the present invention.

[0182] Table 4 Horizontal prediction result indicator data

[0183]

[0184]

[0185] In summary, the present invention performs electromagnetic interference echo removal and deblurring processing on radar data, eliminates data points that affect prediction results, enhances the reliability of radar data, and can improve the accuracy of strong storm areas and prediction results; uses the region growing algorithm to perform initial region recognition, can effectively extract strong convective structure information in radar images, removes weak areas that are too small and may be misidentified results, and accurately displays the recognition effect; combines the LK optical flow method and the Bhattacharyya coefficient method to predict the vertical direction in the strong storm area, focuses on the overall dynamic analysis of the radar echo area, and finely identifies and effectively predicts potential strong wind activities; selects feature parameters and uses SVM, CNN, and U-Net for prediction, and performs weighted voting on the corresponding prediction results, broadens the prediction dimension, enhances the recognition performance, overcomes the limitations of a single model, improves the robustness of the model, and improves the accuracy of storm prediction in the horizontal direction in the strong storm area.

Claims

1. A method for nowcasting convective storm gale, characterized by: The following steps are involved: S1. Acquire multi-source observation data of convective storm gale and perform preprocessing to obtain preprocessed multi-source observation data; S2, using the region growing algorithm to identify the region of the preprocessed multi-source observation data to obtain the strong storm area; S3, using the LK optical flow method and the Bhattacharyya coefficient method to perform storm prediction on the preprocessed multi-source observation data in the vertical direction in the strong storm area, and obtain the corresponding vertical prediction results; S4, constructing and using a horizontal storm prediction model to perform storm prediction on the preprocessed multi-source observation data in the horizontal direction of the strong storm area, and obtaining corresponding horizontal prediction results; S5. Based on the vertical prediction results and the horizontal prediction results, the final prediction results are determined, that is, the nowcast results of the convective storm gale are obtained.

2. The method for nowcasting convective storm gale according to claim 1, characterized in that: The multi-source observation data in step S1 includes weather radar data, automatic station data, wind profiler radar data and sounding radar data of multiple time periods; Weather radar data includes the differential reflectivity factor Z DR , correlation coefficient ρHV, differential phase shift φ DP , differential phase shift rate K DP , depolarization ratio LDR, echo top height, echo centroid height, vertically integrated liquid water content, vertically integrated liquid water content density, lower layer radial velocity, differential reflectivity and differential reflectivity column; Automatic station data include air temperature, pressure, rainfall, wind speed, wind direction, and relative humidity; Wind profile radar data includes overhead wind direction and wind speed; Sounding radar data include temperature, pressure, humidity, and wind at all levels within the vertical range of 0 to 30 km in the atmosphere; The pre-processing in step S1 comprises the following steps: S1-1. Gridding the weather radar data using a bilinear interpolation algorithm to obtain gridded weather radar data; S1-2, using the improved interference echo automatic identification algorithm to remove electromagnetic interference echoes from the gridded weather radar data, to obtain weather radar data after the electromagnetic interference echoes are removed; S1-3, using a velocity defuzzification algorithm to process the weather radar data after removing the electromagnetic interference echo, to obtain the velocity defuzzified weather radar data; The pre-processed multi-source observation data includes automatic station data, wind profiler radar data, sounding radar data and weather radar data after velocity blurring.

3. The method for nowcasting convective storm gale according to claim 2, characterized in that: The specific process of the improved interference echo automatic identification algorithm in step S1-2 is as follows: S1-2-1, determine whether the radial number of each PPI scanning surface of the gridded weather radar data is 360; if so, retain the PPI scanning surface; otherwise, remove the PPI scanning surface to obtain the determined weather radar data; wherein each PPI scanning surface of the gridded weather radar data includes a corresponding radial angle and elevation angle; S1-2-2, draw the radar echo image corresponding to the determined weather radar data, set an area range, and calculate the number of valid numerical points around the center point within the area range; S1-2-3, determine whether the percentage of valid numerical points around the center point within the area exceeds the threshold; if so, the corresponding center point is taken as a valid value; otherwise, the corresponding center point is taken as an invalid value and removed; and the filtered weather radar data is obtained; S1-2-4. Count the total library number and the continuous echo library number in each radial direction of the filtered weather radar data, and record the maximum and minimum values ​​of the total library number; S1-2-5, the total library number, the maximum value, the minimum value and the continuous echo library number of each radial direction, and determine the corresponding radial direction with the maximum total library number of valid echoes; S1-2-6, according to the formula: Get the radial echo continuity index R corresponding to each radial direction t ; Among them, R max Indicates the maximum total number of libraries, R min Indicates the minimum value of the total library number; S1-2-7, determine whether each radial direction meets the judgment standard conditions; if so, determine the corresponding radial direction as the initial interference radial direction and enter step S1-2-8; otherwise, determine the corresponding radial direction as the initial normal radial direction and enter step S1-2-10; wherein, the judgment standard conditions are: when the radial echo continuity index R is met at the same time t Less than 20%, the maximum value of the total number of libraries from the 1st to the 3rd layer elevation angle R max The maximum value of the total number of libraries with elevation angles greater than 300, from the 4th to the 6th layer, R max The maximum value of the total number of libraries with elevation angles greater than 250 from the 7th to the 9th floor is R max Greater than 200 and the minimum value of the total library number R min When it is greater than 20, the radial direction is the initial interference radial direction; S1-2-8, according to the formula: R d =|R N -R i |,i=N±20 Get the absolute value R of the total distance library of effective echoes on the initial interference radial direction and its nearby radial directions d ; Among them, R N Indicates the total distance bin number of effective echoes on the initial interference radial direction, R i represents the total distance bin number of effective echoes in the radial direction near the i-th initial interference radial direction, |·| represents the absolute value, and N represents the initial interference radial direction; S1-2-9. Determine the absolute value R of the total distance library of effective echoes d Is it greater than 20? If so, the initial interference radial direction is determined to be the initial normal radial direction; otherwise, the initial interference radial direction is determined to be the final interference radial direction and is removed; proceed to step S1-2-10; S1-2-10, determine the starting radial direction N1 and the ending radial direction N2 of each initial normal radial direction N, and perform edge effective echo detection on the area [N1-1, N2+1]. If there is no effective echo, remove the interference echo in this area; S1-2-11, repeat step S1-2-10 until a valid echo is obtained, and use the corresponding initial normal radial direction as the normal radial direction; The weather radar data after removing electromagnetic interference echoes includes weather radar data after screening corresponding to a normal radial direction.

4. The method for nowcasting convective storm gale according to claim 2, characterized in that: The step S2 comprises the following steps: S2-1, sampling the weather radar data after de-blurring the velocity to obtain spatially continuous time series data; S2-2, smoothing the spatially continuous time series data by using a 3×3 Gaussian filter kernel to obtain smoothed time series data; S2-3, set threshold; S2-4, constructing a data grid corresponding to the smoothed time series data and obtaining corresponding image size information; S2-5, traverse the data grid and image size information to construct a feature map array of the same dimension; S2-6, initialize the seed point queue structure; S2-7, scan the data grid in order from left to right and from top to bottom, and determine whether the intensity in each grid point exceeds the threshold set in step S2-3; if so, determine the corresponding grid point as a strong echo point and add it to the seed point queue structure; otherwise, determine the corresponding grid point as a weak echo point; S2-8, using the queue rule, gradually add pixels that are adjacent to the seed point and meet the conditions to the current area until no new pixels can be added, and obtain the corresponding strong echo area; if there are still pixels that have not been assigned to any area, select a new seed point and repeat the above steps until all pixels are assigned to a certain area; S2-9, counting the total number of pixels corresponding to each strong echo region to obtain the area of ​​each strong echo region; S2-10, setting an area threshold; judging whether the area of ​​each strong echo region is greater than the area threshold; if so, retaining the corresponding strong echo region; otherwise, removing the corresponding strong echo region; S2-11. Using the feature map array, store the strong echo area retained in step 2-10, that is, the strong storm area.

5. The method for nowcasting convective storm gale according to claim 4, characterized in that: The step S3 comprises the following steps: S3-1. Search for the differential reflectivity factor Z in each severe storm area DR Maximum value; S3-2. Determine each differential reflectivity factor Z DR Whether the maximum value is greater than the threshold TH1 and the area of ​​each strong storm area is greater than TH2; if so, enter step S3-3; otherwise, end the operation; S3-3, select the corresponding de-velocity blurred weather radar data at a time and a time before and after it, and convert the differential reflectivity factor Z DR The maximum value is used as the section, the corresponding section array is generated, and the section image is drawn; S3-4, using the LK optical flow method to analyze the profile arrays at three moments, and obtain the corresponding optical flow fields; S3-5, calculate the storm core drop rate of the light process; determine whether the storm core drop rate is greater than the threshold TH3; if so, proceed to step S3-6; otherwise, terminate the operation; S3-6, read the radar radial velocity low elevation angle volume scan data corresponding to the radar data, and use the value TH4 as the radius range, count the positive and negative velocity area distribution within the radius range, and obtain the corresponding histogram; S3-7, converting the histogram into probability density distribution of different speeds in positive and negative speed regions, and drawing the corresponding probability density distribution histogram; S3-8, according to the Bhattacharyya coefficient formula: Get the similarity ρ(ρ1,ρ2); where ρ1 represents the percentage of the speed value of each region with positive speed in the total, ρ2 represents the percentage of the speed value of each region with negative speed in the total, ∑(·) represents the percentage, TH6 represents the total number of regions with positive and negative speeds, ρ1(i) and ρ2(i) represent the percentage of the speed value of each region with the i-th positive speed in the total, and the percentage of the speed value of each region with the i-th negative speed in the total, respectively; S3-9, determine whether the similarity ρ(ρ1, ρ2) is greater than the threshold TH7; if so, determine that a downburst event is very likely to occur in the strong storm area and enter step S3-10; otherwise, terminate the operation; S3-10. Combining the automatic station data, wind profiler radar data, and sounding radar data, the optical flow field is calculated using the linear extrapolation method, and the Lagrangian mechanics model is used for extrapolation to obtain the core reflectivity factor profile of the downburst storm cell, that is, the vertical prediction result.

6. The method for nowcasting convective storm gale according to claim 5, characterized in that: The step S3-4 comprises the following steps: S3-4-1, converting the dimensions of the profile arrays at three moments to obtain corresponding converted profile arrays; S3-4-2, respectively perform LK optical flow calculation on the converted profile arrays at two adjacent moments, and obtain the initial optical flow vector field between the two adjacent moments; S3-4-3, decomposing the velocity vectors of the two initial optical flow vector fields into two components parallel to the x-axis and the y-axis; determining the median of the components in two directions of each initial optical flow vector field, and merging the two medians to obtain the corresponding initial velocity vector; S3-4-4, combine the two initial velocity vectors to obtain the final velocity vector; S3-4-5, using the final velocity vector to perform extrapolation movement on the cross-section data corresponding to a moment selected in step S3-3, to obtain a corresponding moved image; S3-4-6, calculating the initial optical flow field of the original image composed of the moved image and the cross-sectional data corresponding to a moment selected in step S3-3; S3-4-7. Perform matrix addition on the initial optical flow field and the corresponding two initial optical flow vector fields in S3-4-2 to obtain the corresponding intermediate optical flow field; perform weighted addition on the two intermediate optical flow fields to obtain the optical flow field.

7. The method for nowcasting convective storm gale according to claim 2, characterized in that: The optical flow window size of the LK optical flow calculation in step S3-4-2 is 25×25, the maximum number of pyramid layers is set to 4 layers, the maximum number of iterations is 10 times, and the error threshold is 0.

03.

8. The method for nowcasting convective storm gale according to claim 2, characterized in that: The horizontal storm prediction model in step S4 includes an SVM support vector machine, a CNN convolutional neural network, and a U-Net sub-model.

9. The method for nowcasting convective storm gale according to claim 8, characterized in that: The step S4 comprises the following steps: S4-1. Determine characteristic data based on the weather radar data after velocity blurring; wherein the characteristic data includes the maximum differential reflectivity factor, echo top height, maximum differential phase shift rate, echo centroid height, vertically integrated liquid water content, vertically integrated liquid water content density, low-layer radial velocity, maximum differential reflectivity and differential reflectivity column; S4-2, input the feature data into the SVM support vector machine to obtain the corresponding SVM prediction result; S4-3, input the feature data into the CNN convolutional neural network to obtain the corresponding CNN prediction result; S4-4, input the feature data and into the U-Net sub-model to obtain the corresponding U-Net prediction result; S4-5, using weighted voting to process the SVM prediction results, CNN prediction results, and U-Net prediction results to obtain the final prediction results, that is, to obtain the initial level prediction results; S4-6. Based on the automatic station data, wind profiler radar data, and sounding radar data, the initial horizontal prediction results are adjusted to obtain horizontal prediction results.

10. The method for nowcasting convective storm gale according to claim 2, characterized in that: The final prediction result in step S5 is the overlapped area of ​​the vertical prediction result and the horizontal prediction result.

Images