Three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint

Abstract

The invention relates to a three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint. The method comprises the following steps: S1, continuously collecting radar reflectivity and radial speed data within 2 hours; s2, converting radar reflectivity data into three-dimensional data under a Cartesian rectangular coordinate system, mapping radial speed data into a radial rectangular coordinate system, and constructing a sample set; s3, acquiring middle-layer radial convergence feature information in the radial velocity diagram; s4, building a radar echo extrapolation model; s5, introducing deformation convolution to complete automatic extraction of the dynamic field information of the convection system; s6, dividing the sample set into a training set, a verification set and a test set; and S7, obtaining a reflectivity image in the next one hour. According to the method, the three-dimensional structural characteristics and the dynamic characteristics of the weather system are fully considered, and meanwhile, the complex evolution process of the weather system is accurately simulated and predicted in combination with the multi-source physical constraint information.

Description

Technical Field

[0001] The invention belongs to the technical field of meteorological radars, and relates to technologies related to radar signal processing, echo extrapolation and radar data analysis, and specifically to a three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraints. Background Art

[0002] The development and evolution of weather systems are accompanied by complex atmospheric activities, among which horizontal atmospheric flow dominates advection motion, while vertical convection plays a key role in affecting its life cycle. Existing radar echo extrapolation technology is insufficient in simulating complex changes such as the formation and extinction of weather systems. Traditional methods rely on radar reflectivity to estimate cloud motion vectors, which are only applicable to simple linear models. Although deep learning can capture nonlinear motion, it still needs to be improved in reflecting echo generation and extinction.

[0003] In view of this, in order to improve the accuracy and reliability of weather radar echo extrapolation, it is necessary to analyze the weather system as a three-dimensional whole so that the model can deeply learn the spatiotemporal correlation between different altitude layers. At present, some scholars have conducted research on this, but most models are only based on historical radar echo data, and do not fully consider physical constraint information such as environmental thermal instability, vertical wind shear, and vertical distribution of water vapor, which are crucial to the development and evolution of severe convective weather. Recent studies have begun to incorporate multi-source data, but there is currently a lack of research on integrating the strong downward airflow within the convective system, namely the mid-level radial convergence (MARC) feature, into the extrapolation model. As an important indicator of the dynamic characteristics of the convective system, MARC is of great significance for predicting the motion trajectory and speed of clouds, and taking it into consideration can improve the accuracy of radar extrapolation technology in predicting high-intensity reflection areas.

[0004] In order to further improve the accuracy of radar extrapolation technology in predicting high-intensity reflection areas, a three-dimensional radar echo extrapolation method that integrates radial velocity information is urgently needed. Summary of the invention

[0005] The purpose of the present invention is to overcome the shortcomings of the prior art and provide a three-dimensional radar echo extrapolation method based on deformed convolution dynamic field constraints, so as to fully consider the three-dimensional structural characteristics and dynamic characteristics of the weather system, and at the same time combine multi-source physical constraint information to achieve accurate simulation and prediction of the complex evolution process of the weather system.

[0006] The present invention solves the technical problem by the following technical solutions:

[0007] A three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraints, the method comprises the following steps:

[0008] S1. The radar operates in VCP21 volume scanning mode, that is, it completes 9 scans at different elevation angles every 6 minutes, and continuously collects radar reflectivity and radial velocity data for 2 hours;

[0009] S2. For the data at each radar volume scan time: convert the radar reflectivity data from the polar coordinate system into three-dimensional data in the Cartesian rectangular coordinate system; perform deblurring processing on the radial velocity data and map the radial velocity data into the radial rectangular coordinate system; use the reflectivity and radial velocity images within the first hour as inputs to the extrapolation model, and use the reflectivity images within the second hour as labels to construct a sample set;

[0010] S3, using the mid-layer radial convergence automatic identification algorithm to obtain the mid-layer radial convergence feature information in the radial velocity map;

[0011] S4, build a radar echo extrapolation model that takes the three-dimensional radar reflectivity and radial velocity images within one hour as input and the reflectivity images of the following hour as output;

[0012] S5. Introduce deformation convolution into the radar echo extrapolation model to complete the automatic extraction of the dynamic field information of the convective system;

[0013] S6. Dividing the sample set into a training set, a validation set, and a test set, training the radar echo extrapolation model on the training set, performing hyperparameter tuning on the validation set, and testing the radar echo extrapolation model on the test set;

[0014] S7. Using the radar reflectivity data within any hour, convert the radar reflectivity data at each radar scanning moment from the polar coordinate system into three-dimensional data under the Cartesian rectangular coordinate system using the steps of S2, deblur the radial velocity data at each elevation angle, map the radial velocity data to the radial rectangular coordinate system, and use the above two types of data within the hour as inputs to the model to obtain the reflectivity image for the next hour.

[0015] Moreover, the S2 converts the radar reflectivity data at each radar scanning moment from the polar coordinate system into three-dimensional data in the Cartesian rectangular coordinate system by first converting the radar base data from the polar coordinate system to the Cartesian coordinate system using a bilinear interpolation method in the horizontal dimension, and then interpolating in the vertical dimension to obtain three-dimensional equidistant grid field data;

[0016] The step S2 of performing deblurring on the radial velocity data at each elevation angle refers to using a deblurring algorithm based on minimizing the velocity difference between regions to solve the velocity ambiguity problem, and mapping the deblurred radial velocity data to a radial rectangular coordinate system.

[0017] Moreover, the S4 adopts the strategy of "independent encoding followed by fusion" to process the three-dimensional radar echo images of ten consecutive moments in the past hour. The reflectivity data and the mid-layer radial convergence structure information obtain feature mappings through their respective encoders. The decoder integrates the two types of feature information and decodes and generates the future extrapolated echo sequence based on it. The radar echo extrapolation model completes the effective extrapolation of the three-dimensional radar echo data by jointly using the ConvLSTM network based on 3D convolution and the convolutional neural network CNN.

[0018] The ConvLSTM network based on 3D convolution replaces the two-dimensional convolution operation of the ConvLSTM unit with a three-dimensional convolution operation to process the three-dimensional radar data corresponding to each radar scanning moment, and uses the unique loop structure of the ConvLSTM unit to extract and synthesize the features of the three-dimensional radar data at different times within an hour; each ConvLSTM unit includes three key components: a forget gate, an input gate, and an output gate. The unit receives the data input X at the current time point. t , the state output H at the previous time point t-1 and long-term memory unit C t-1 As input, it generates an updated long-term memory unit C t and the state output H at the current time point t As output, the operations performed by the ConvLSTM cell and the calculation of its output are shown below:

[0019] i t =σ(W xi *X t +W hi *H t-1 +W ci *C t-1 +b i )

[0020] f t =σ(W xf *X t +(W hf *H t-1 +W cf *C t-1 +b f )

[0021] o t =σ(W xo *X t +W ho *H t-1 +W co *C t-1 +b o )

[0022]

[0023] In the formula, * represents the convolution operation, represents the Hadamard product, σ is the activation function, W and b represent the corresponding convolution kernel weight and bias term, tanh is the hyperbolic tangent function, i t 、f t and t They represent input gate, forget gate and output gate respectively.

[0024] Moreover, the deformation convolution of S5 refers to adjusting the convolution kernel shape according to actual needs, further expanding the extraction and representation of effective speed features and their associated information. The two modes of deformation convolution are:

[0025] (1) Based on the maximum statistical distance between positive and negative velocity points, the height K of the convolution kernel in each row of the deformed convolution is H The width of the convolution kernel is fixed to 15. W It is determined dynamically by the following formula according to its specific radial position in the image;

[0026]

[0027] Where r represents the distance from the center of the convolution kernel to the top of the radial velocity image;

[0028] In order to achieve the maximum speed difference dv max and the average speed difference dv mean The quantization of the feature parameters included, the weight distribution inside the convolution kernel is determined by the following formula;

[0029]

[0030] Where i represents the distance between the pixel point and the top of the template, j represents the distance between the pixel point and the left side of the template, (i, j) is the coordinate position of the pixel point in the template, Ω is the relevant area for feature calculation, and K is l (i,j) is the weight at the coordinate (i,j);

[0031] (2) Using the template matching method, we implement speed matching and calculation of different zones, limit the distance between the speed extreme points, and evaluate the relative position between the center of the maximum inflow speed and the center of the outflow speed, so as to accurately identify the key speed features in the airflow field. The airflow field speed feature identification template is in the shape of a regular hexagon, which is specifically used to identify four types of airflow field structures: convergence, divergence, cyclone, and anticyclone.

[0032] Based on the consideration of statistical distance, the height K of the convolution kernel in each row of the deformed convolution used to determine the rotation type is H The width of the convolution kernel is fixed to 11. W It is determined dynamically by the following formula according to its specific position in the image;

[0033]

[0034] Where h∈{-5,-4,-3,…,3,4,5}.

[0035] Moreover, the S6 trains the radar echo extrapolation model using a weighted mean square error (WMSE) as a model loss function to guide the model to pay more attention to the prediction accuracy of the strong echo area. The WMSE is defined as follows:

[0036]

[0037] In the formula, I h,x,y Indicates the reflectivity at the coordinates (h,x,y) in the label image. Represents the reflectivity at (h, x, y) in the predicted image, Weight h,x,y The weight value representing the position can be set according to the different echo strengths to make the model pay more attention to the strong echo area.

[0038] The positive effects that the present invention can produce are:

[0039] The present invention analyzes the weather system as a three-dimensional whole, incorporates the mid-layer radial convergence characteristics within the convective system to improve the accuracy of the prediction of high-intensity reflection areas, uses the 3D convolution ConvLSTM network and the convolutional neural network CNN to build a radar echo extrapolation model, adopts the "independent encoding and then fusion" strategy to encode the reflectivity data and MARC structure information separately and then integrate them by the decoder, and introduces the deformation convolution method to expand the extraction and characterization of effective velocity characteristics and their associated information. The convergence feature validity test and the deformation convolution test verified that the present invention integrates the dynamic field structure constraints into the extrapolation model and applies the deformation convolution to significantly improve the accuracy of radar echo extrapolation prediction. BRIEF DESCRIPTION OF THE DRAWINGS

[0040] Figure 1 is a flow chart of the present invention;

[0041] Figure 2 It is a schematic diagram of a three-dimensional radar echo extrapolation model based on a deformed convolution dynamic field constraint of the present invention;

[0042] Figure 3 It is a flow chart of de-blurring velocity in a region of a single elevation angle radial velocity diagram of the present invention;

[0043] Figure 4 It is a flow chart of deblurring velocity between regions of a single elevation radial velocity diagram of the present invention;

[0044] Figure 5 It is a schematic diagram of the template-defined area in the polar coordinate system of the present invention;

[0045] Figure 6 It is a schematic diagram of the template-defined area in the radial rectangular coordinate system of the present invention;

[0046] Figure 7 It is a schematic diagram of the three-dimensional echo extrapolation model under the MARC constraint of the present invention;

[0047] Figure 8 Schematic diagram of the airflow field velocity feature recognition template of the present invention. DETAILED DESCRIPTION

[0048] The present invention is further described in detail below through specific examples. The following examples are only illustrative and not restrictive, and the protection scope of the present invention cannot be limited thereto.

[0049] A three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraints, the innovation of which lies in: the steps of the method are:

[0050] (1) S-band Doppler weather radar data were collected hourly from May 2015 to October 2016 from 12 radar sites in 15 different locations in North and East China. In VCP21 volume scanning mode, the radar will complete conical scans at 9 preset elevation angles in sequence within 6 minutes, thereby continuously collecting radar reflectivity and radial velocity data for 2 hours.

[0051] (2) The radar reflectivity data at each radar scanning moment is converted from the polar coordinate system to three-dimensional data in the Cartesian rectangular coordinate system, the radial velocity data at each elevation angle is deblurred, and the radial velocity data is mapped to the radial rectangular coordinate system. The radar image data of 20 consecutive moments is defined as a sample. Among them, the first 10 moments cover a period of 1 hour (6min×10), and the corresponding radar image data includes both reflectivity and radial velocity types, which are used as the input of this radar echo extrapolation model, aiming to predict the next 10 consecutive radar images in the next hour. The last 10 radar images in the sample are composed of only one type of reflectivity, which is used as label data for supervising model training and for evaluating the performance of the model. This embodiment includes 9540 samples, and there is no time overlap between the samples, which together constitute a sample set.

[0052] Defuzzification of radial velocity data at each elevation angle refers to defuzzification of radial velocity data using a defuzzification algorithm MVDR based on minimization of inter-region velocity differences. Velocity fuzzification refers to the situation where the true radial velocity exceeds the Nyquist unambiguous interval [-V N ,V N], the measured velocity value will be "folded" and cause deviation, which is manifested as false jumps of positive and negative maximum velocity values ​​in adjacent areas on the radial velocity diagram. The velocity defuzzification algorithm MVDR includes two parts: velocity defuzzification within isolated areas and velocity defuzzification between areas.

[0053] To deblur the speed in a region, we first need to calculate the speed difference between each point in the connected region and its adjacent points to identify possible speed ambiguity points. Based on these points, the region is further divided into multiple sub-regions. Next, deblurring is performed based on the principle of minimizing the speed difference between sub-regions. Let the number of isolated regions be n r The process of deblurring the velocity in the region of the single elevation radial velocity diagram is as follows: Figure 3 shown.

[0054] Inter-region velocity defuzzification is a method for dealing with velocity ambiguity in isolated velocity regions far from the center of the radar. The specific processing steps are as follows: First, start from the region close to the radar and gradually advance the defuzzification operation to the region far from the radar. Second, give greater weight to the region with closer distance and higher credibility to ensure the accuracy of the reference velocity. Assume that the number of regions to be executed is n a , sort these areas in ascending order of distance from the radar, and the inter-area velocity deblurring process of the single elevation radial velocity map is as follows: Figure 4 shown.

[0055] After velocity deblurring, high-resolution radial velocity information at nine elevation angles in the polar coordinate system is obtained. On this basis, the radial velocity data is mapped to the radial rectangular coordinate system, thus overcoming the problem that polar coordinate representation is not conducive to computer processing and the azimuth information is lost in the Cartesian coordinate system. The specific definition of this coordinate system is: in the plane rectangular coordinate system, the horizontal axis is positive to the right, representing the angle information in the polar coordinate system, where the origin 0° represents the due north direction, and the angle resolution is 1°; the vertical axis is positive downward, representing the radial distance of each point to the radar, and the distance resolution is 1 km. In order to prevent the convergence zone near the due north direction from being cut off, the image is expanded by 20°.

[0056] (3) There are two main steps when the MARC automatic recognition algorithm extracts MARC feature information from the radial velocity map. First, the convergence points are identified in the radial velocity map at each elevation angle, and the characteristic parameters such as its convergence intensity are calculated through the template to determine the strong convergence area at each elevation layer. Then, the strong convergence area is matched longitudinally, and the strong convergence area that exists simultaneously at two or more elevation angles is identified as a strong convergence field, and its characteristic vector is given to obtain its three-dimensional structure and position information. The specific steps are as follows:

[0057] 1. Convergence point identification and convergence intensity calculation

[0058] In the reflectivity map under the radial rectangular coordinate system, the area with reflectivity value higher than 35 dBZ is extracted as a mask, and for the radial velocity data under the same coordinate system, only the foreground area in the mask is retained as the target research area in the subsequent process.

[0059] For the processed single elevation radial velocity map, extract the boundaries of all positive velocity regions (v>0) and negative velocity regions (v≤0), determine the adjacent boundary points of the positive and negative velocity regions, and identify the boundary points where the upper and lower velocity values ​​change from positive to negative as convergence points. For each convergence point p, calculate the eigenvector in its related area Ω. The parameters are the maximum velocity difference dv max , average speed difference dv mean and the coordinates of the positive and negative velocity extreme points and The calculation process of each parameter is as follows:

[0060]

[0061] Where n + , and Represents the number of positive velocity pixels, the total positive velocity, and the maximum positive velocity, n - , and Represent the number of negative speed pixels, the sum of negative speeds, and the maximum negative speed, respectively. Ω represents the positive and negative speed areas centered on point p that meet specific spatial and scale constraints. Here, a template centered on point p is used to select the relevant area. Specifically, the spatial constraint requires that the deviation between the connecting line between the maximum positive speed point and the maximum negative speed point and the radial direction should not be too large; the scale constraint stipulates that the radial distance between the positive and negative speed points should not be too far. Based on this, the template is designed as two isosceles triangles with a common vertex, and the area covered by the template is Ω. The template-defined area in the polar coordinate system is as follows: Figure 5 shown.

[0062] In the polar coordinate system, the size of the two isosceles triangle templates is fixed, but they need to rotate with the azimuth angle to keep consistent with the radial direction, which increases the computational complexity. In the radial rectangular coordinate system, the template can be directly aligned with the radial direction, but coordinate system mapping is required. When the radial distance r between the center point of the template and the radar increases, the central angle corresponding to the bottom edge of the template becomes smaller. The template-limited area in the radial rectangular coordinate system is as follows: Figure 6 As shown. If the coordinates of the point to be measured p are (θ, r), the arc length calculation formula can be used to calculate that the arc length covered by the distance h from point p in the radial direction is approximately equal to (4 / 7)|h|×2, and the corresponding central angle range can also be determined as [θ-k h ,θ+kh ],in:

[0063]

[0064] 2. Single elevation angle strong convergence zone positioning

[0065] When the characteristic parameters of an amplitude point satisfy dv max ≥25m / s and dv mean When the speed is ≥12m / s, it is recorded as a "strong convergence point". The judgment of the average speed difference is to prevent individual noise points from generating misleading high values. The threshold of 12m / s comes from the statistical analysis results of the samples. For each strong convergence area, the size of its minimum enclosing rectangle ω, the coordinates of the center point and the characteristic vector of the convergence point are recorded. Finally, by traversing and comparing the feature vectors of all points in all strong convergence areas, the representative point p of each strong convergence area is determined. 0 .

[0066] 3. Longitudinal matching of strong convergence zone

[0067] The method of searching from low elevation layer to upward layer is adopted to match the strong convergence areas between adjacent elevation layers. When the center distance of the outer rectangles of the convergence areas in the adjacent high and low elevation layers is less than 15 kilometers, the two convergence areas are considered to be related. In this process, there are two special cases that need to be explained:

[0068] First, if there is a "many-to-one" or "one-to-many" convergence zone association, this study selects dv max The largest convergence zone is retained and the rest are discarded.

[0069] Second, for the stronger convergence zones that fail to find matching objects at adjacent elevation angles, especially when these unmatched convergence zones appear at higher elevation layers, they may represent newly emerged stronger convergence zones, so they should be retained and continue to be matched with convergence zones at higher elevation layers.

[0070] Once a strong convergence area can obtain at least two correlation matches at different elevation angles in the vertical direction, it can be judged as a significant MARC.

[0071] (4) Build a radar echo extrapolation model that takes the three-dimensional radar reflectivity and radial velocity images within one hour as input and the reflectivity images of the following hour as output. This study adopts a dual encoder input mechanism and adopts the strategy of "independent encoding followed by fusion" to process the data based on the physical properties and semantic differences of the input data. Specifically, the reflectivity data and MARC structure information obtain feature maps through their respective encoders, and the decoder integrates these two types of feature information and decodes them to generate the future extrapolated echo sequence.

[0072] The core building blocks of the model are ConvLSTM and CNN. In view of the characteristics of 3D radar data, the ConvLSTM unit is adapted and replaced with 3D convolution to form a 3DConvLSTM unit, which significantly improves the model's ability to process 3D data. Finally, through the combined use of 3DConvLSTM and CNN units, effective extrapolation of 3D radar echo data is achieved. The structure of the 3D echo extrapolation model under MARC constraints is shown in the figure. Figure 7 shown.

[0073] The reflectivity encoder receives the three-dimensional echo data at each time point as input. The CNN module maps the original image data to the feature space, and then fuses the echo features from different height levels to express the redistribution of echo intensity caused by spatiotemporal atmospheric motion. The reflectivity encoder contains three three-dimensional convolutional layers, among which the first convolutional layer keeps the spatial dimension of the input image unchanged, promoting the effective conversion from image space to feature space; the subsequent two convolutional layers realize the deep fusion of cross-level features through specific parameter configuration. The reflectivity encoder parameter details are shown in Table 1.

[0074] Table 1 Reflectivity encoder parameter information table

[0075]

[0076] The dynamic field encoder processes the radial velocity data at the same time point, extracts the dynamic properties of the convective system and the dynamic behavior of the cloud, including the three-dimensional structure and position information of the MARC features. The encoder is composed of a stack of 3DConvLSTM modules and CNN modules, which deeply calculates the spatiotemporal sequence data, expresses the potential data patterns and dependencies, and thus provides accurate constraints for radar echoes. The parameters of the dynamic field encoder are shown in Table 2.

[0077] Table 2 Power field encoder parameter information table

[0078]

[0079] The decoder integrates the implicit representations of the reflectivity encoder and the dynamic field encoder outputs, and reconstructs the morphology, motion, and evolution of the echo through cyclic decoding operations. The decoder extracts the echo features at future time points from the hidden state output by 3DConvLSTM, and uses deconvolution technology to reversely map these features from the abstract feature space back to the original image space, completes the data upsampling, and finally generates the predicted radar echo sequence. The decoder parameters are shown in Table 3.

[0080] Table 3 Decoder parameter information table

[0081]

[0082] (5) This embodiment uses a deformed convolution method to improve the radar echo extrapolation model, adjusts the convolution kernel shape according to actual needs, simplifies the calculation process, and further expands the extraction and representation of effective velocity features and their associated information to achieve rapid and automatic MARC recognition. Figure 2 Specifically, there are two main methods for deformable convolution:

[0083] 1. Based on the consideration of the maximum statistical distance between positive and negative velocity points, the height K of the convolution kernel in each row of the deformed convolution is H The width of the convolution kernel is fixed to 15. W It is dynamically calculated and determined by the following formula according to its specific radial position in the image.

[0084]

[0085] Where r represents the distance from the center of the convolution kernel to the top of the radial velocity image.

[0086] In order to achieve the maximum speed difference dv max and the average speed difference dv mean The quantization of the feature parameters included, the weight distribution inside the convolution kernel is determined by the following formula.

[0087]

[0088] Where i represents the distance from the pixel to the top of the template, j represents the distance from the pixel to the left of the template, (i, j) is the coordinate position of the pixel in the template, Ω is the relevant area for feature calculation, and its definition process is described in step (3). K l (i,j) is the weight at the coordinate (i,j).

[0089] Through the above weight settings, the deformed convolution can perform the convolution operation along the specified direction, process the input radial velocity data, and finally generate the corresponding feature vector.

[0090] 2. In the complex airflow environment of the convective system, in addition to the airflow convergence structure, the various structural features such as airflow divergence and rotation presented by the relative position of the maximum inflow and outflow velocity centers are important indicators of the potential occurrence of severe convective weather. The present invention adopts a template matching method to implement partitioned velocity matching and calculation, limit the distance between velocity extreme points and evaluate the relative position between the maximum inflow velocity center and the outflow velocity center, and accurately identify key velocity features in the airflow field. The airflow field velocity feature recognition template is as follows: Figure 8 shown.

[0091] The feature recognition template is in the shape of a regular hexagon and is specifically used to identify four types of airflow field structures: convergence, divergence, cyclone, and anticyclone. In the template, 1) isosceles triangles a and b with common vertices and the same scale are responsible for capturing the two types of features: convergence and divergence. Their height is set according to the maximum statistical distance between positive and negative velocity points in the convergence and divergence fields, and the length of the base is determined according to the maximum allowable value of Δθ; 2) rhombuses c and d with the same scale are used to capture the two types of features: cyclone and anticyclone. The spatial constraints are: the offset angle between the line between the maximum positive velocity point and the maximum negative velocity point and the radial direction must be large enough, and the radial distance between the positive / negative velocity center points should not be too large. The above triangle group and rhombus group form a regular hexagonal area with a complementary relationship in spatial distribution. The center point of the area is denoted as p. The convergence, divergence, cyclone, and anticyclone attributes of this point are obtained from the data of the a, b area and the c, d area in the template.

[0092] Based on the consideration of statistical distance, the height K of the convolution kernel in each row of the deformed convolution used to determine the rotation type is H The width of the convolution kernel is fixed to 11. W It is determined dynamically by the following formula based on its specific position in the image.

[0093]

[0094] Where h∈{-5,-4,-3,…,3,4,5}.

[0095] (6) Dividing the sample set into a training set, a validation set, and a test set, training the radar echo extrapolation model on the training set, performing hyperparameter tuning on the validation set, and performing model testing on the test set.

[0096] Each radar completes a volume scan every 6 minutes on average, using 10 consecutive 3D radar reflectivity and radial velocity images within one hour as radar echo extrapolation model input, and 10 consecutive reflectivity images within the next hour as labels, a total of 9540 samples, which are randomly assigned to the training set, validation set, and test set in a ratio of 4:1:1. Finally, 6360 training set samples, 1590 validation set samples, and 1590 test set samples were obtained.

[0097] The radar echo extrapolation model is trained and the weighted root mean square loss WMSE is used as the model loss function to guide the model to pay more attention to the prediction accuracy of the strong echo area. WMSE is defined as follows:

[0098]

[0099] In the formula, I h,x,y Indicates the reflectivity at the coordinates (h,x,y) in the label image. Represents the reflectivity at (h, x, y) in the predicted image, Weight h,x,y The weight value representing the position can be set according to the different echo strengths to make the model pay more attention to the strong echo area. h,x,y The setup is as follows:

[0100]

[0101] (7) Using the radar reflectivity data and radial velocity data within any hour, the radar reflectivity data at each radar scanning moment is converted from the polar coordinate system to three-dimensional data in the Cartesian rectangular coordinate system, the radial velocity data at each elevation angle is deblurred, and the radial velocity data is mapped to the radial rectangular coordinate system. The above two types of data within the hour are used as the input of the model to obtain the reflectivity image of the next hour.

[0102] This embodiment uses the critical success index CSI, Heidke skill score HSS and fair skill score ETS to measure the consistency between the proportion of estimated reflectivity above the threshold and the proportion of actual reflectivity above the threshold under a certain reflectivity threshold. The larger the values ​​of CSI, ETS and HSS, the better the model performance. CSI, HSS and ETS are defined as follows:

[0103]

[0104] The calculation formula for REF is as follows:

[0105]

[0106] Among them, TP, FN, FP, and TN are defined by the following confusion matrix:

[0107] Table 4 Confusion matrix

[0108]

[0109] The experimental results and analysis of this embodiment include three parts. The first two parts are the verification of the effectiveness of the model architecture, including the determination of different deformation convolution strategies and the effectiveness of incorporating dynamic field convergence features. The third part is a comparative test between the model of this embodiment and other typical models.

[0110] 1. Comparative test of different deformation convolution strategies

[0111] The identification of significant MARC requires the maximum speed difference dv to be met at the same time max ≥25m / s and average speed difference dv mean≥12m / s, here, the former is called condition 1 and the latter is called condition 2. Among them, condition 1 is the primary condition and is also an important criterion for determining whether this MARC reaches the "significant level"; condition 2 is a supplementary condition for eliminating false "significant MARC" caused by noise. A statistical analysis was performed on all the results that actually contain significant MARC in the modeling data, and it was found that the proportion of MARC that can meet condition 2 and can also meet the "significant MARC" standard (condition 1) exceeds 95%. Therefore, this embodiment designs a comparative test scheme of "different deformation convolution strategies", that is, based on Figure 2 The following four models were generated successively and the effect comparison tests of these four models were organized:

[0112] Model 1, let Figure 2 The "deformed convolution" in the algorithm uses the convolution kernel calculation method based on the maximum statistical distance between positive and negative velocity points, and uses the above-mentioned MARC identification condition 1 to obtain the convergence-divergence type dynamic field three-dimensional distribution through the encoder, and cooperates with Figure 2 The other modules are generated by training the training sample set and the validation sample set.

[0113] Model 2 only replaces "Condition 1" of MARC identification in Model 1 with "Condition 2" and is trained and generated using the same sample set.

[0114] Model 3, combining condition 1 of model 1 and condition 2 of model 2, integrates two types of convergent-divergent dynamic field three-dimensional distribution, and synergizes Figure 2 The other modules are trained and generated using the same sample set.

[0115] Model 4, select the "MARC condition" used by the model with the best test quality among models 1 to 3, and obtain the three-dimensional distribution of the convergence-divergence dynamic field through the encoder; Figure 2 The "deformation convolution" in the algorithm uses the convolution kernel calculation method based on the maximum statistical distance between positive and negative velocity points, and obtains the three-dimensional distribution of the rotational dynamic field through the encoder; then the above two types of three-dimensional dynamic distribution fields are fused to coordinate Figure 2 The other modules are generated by training with the same sample set.

[0116] The above four models were tested using the test sample set that did not participate in the model training. The prediction results given by the model extrapolated for 1 hour were used to calculate the means of the three indicators of CSI, HSS and ETS. The evaluation results based on all test samples are shown in Table 5.

[0117] ① Comparing Model 2 and Model 4 under different thresholds, it can be seen that in the establishment of the strong convective reflectivity echo extrapolation model, the more types of dynamic field information are incorporated, the better, and the integration of convergence-divergence dynamic field information has the greatest improvement on the radar echo extrapolation effect;

[0118] ② Comparing models 1 to 3 under different thresholds, it can be seen that model 2 performs best. This result is explainable. In the process of acquiring radial velocity data, various errors are inevitable due to the existence of various interferences and the limitations of the equipment itself. Even though the present embodiment has performed velocity deblurring, it is still impossible to completely avoid misleading high values ​​generated by individual noise points. The maximum velocity difference is easily affected by noise, resulting in the application effect in the extrapolation task being inferior to the extrapolation effect based on the average velocity difference.

[0119] In particular, the winning model 2 under “different deformation convolution strategies” is recorded as 3DConvLSTM-Dconv.

[0120] Table 5 Test scores obtained by the models under four deformation convolution strategies

[0121]

[0122] 2. Effectiveness of integration and convergence features and comparative testing of different integration schemes

[0123] In order to verify the effectiveness of the three-dimensional radar echo extrapolation model under dynamic field constraints proposed in this embodiment, this embodiment uses a test sample set that does not participate in model training to Figure 7 The extrapolation model 3DConvLSTM-MARC under the constraint of the middle-layer radial convergence feature field, the extrapolation model 3DConvLSTM-Dconv based on deformable convolution, and the extrapolation model 3DConvLSTM without incorporating dynamic field information are organized for comparative testing. The test results are shown in Table 6.

[0124] Table 6 Scoring results of three extrapolation models with and without dynamic field information

[0125]

[0126] In terms of the three comprehensive evaluation indicators of CSI, HSS and ETS, the 3DConvLSTM model without dynamic field information scored the lowest. After introducing MARC distribution information, the CSI score increased by 0.99% to 1.53%, and the HSS score increased by 0.77% to 1.88%, of which the improvement rates at the 40dBZ threshold were 5.8% and 4.5%, respectively. The improved extrapolation model 3DConvLSTM-Dconv based on deformable convolution has been further improved on the basis of the 3DConvLSTM-MARC model. Compared with the 3DConvLSTM model, the CSI score increased by 1.42% to 1.98%, and the HSS score increased by 1.2% to 2.43%, of which the improvement rates at the 40dBZ threshold were 7.4% and 5.8%, respectively.

[0127] It can be seen that incorporating the dynamic field structure constraints based on radial velocity data into the radar echo extrapolation model can significantly improve the accuracy of the model's extrapolation for areas with strong reflectivity. Compared with the 3DConvLSTM-MARC model, the 3DConvLSTM-Dconv extrapolation model can more comprehensively capture the shape and change patterns of the convective system, especially when the prediction time is long, the extrapolation results it provides show higher accuracy.

[0128] 3. Comparison test with typical extrapolation model based only on reflectivity data

[0129] In order to verify the effectiveness of the extrapolation model proposed in this embodiment, other typical extrapolation models based only on reflectivity data are selected for comparative testing. The comparative models include constant extrapolation Persistence, Optical flow method and ConvLSTM. The CSI, HSS and ETS scoring results of the above three comparison models and the model of this embodiment for the 2km altitude reflectivity echo at different reflectivity thresholds are shown in Table 7.

[0130] Table 7 Scoring results of different models

[0131]

[0132] Under various reflectivity thresholds, the CSI, HSS and ETS scores of the model in this embodiment are the highest. Taking the CSI score as an example, compared with the traditional optical flow algorithm, the result of the deep learning method ConvLSTM is 2.22% to 4.57% higher than it, and the 3DConvLSTM-Dconv model of this embodiment is further improved by 4.74% to 6.47% compared with the ConvLSTM method, especially in the strong reflectivity (40dBZ) area, relative to ConvLSTM, optical flow method and constant extrapolation, the 3DConvLSTM-Dconv model has an improvement rate of 26%, 56% and 179% on CSI respectively.

[0133] In general, the score of the extrapolation model based on deep learning is better than that of the traditional extrapolation algorithm, which reflects its advantages in nonlinear process modeling. The 3DConvLSTM-Dconv model further improves the prediction accuracy based on the ConvLSTM method, confirming that the modeling scheme of increasing the data dimension from 2D to 3D and then incorporating the dynamic field constraint based on radial velocity data is effective and feasible, and the improvement of model quality is significant.

[0134] Although the embodiments and drawings of the present invention are disclosed for illustrative purposes, those skilled in the art will appreciate that various substitutions, changes and modifications are possible without departing from the spirit and scope of the present invention and the appended claims. Therefore, the scope of the present invention is not limited to the contents disclosed in the embodiments and drawings.

Claims

1. A three-dimensional radar echo extrapolation method based on deformed convolution dynamic field constraints, characterized in that: The steps of the method are: S1. The radar operates in VCP21 volume scanning mode, that is, it completes 9 scans at different elevation angles every 6 minutes, and continuously collects radar reflectivity and radial velocity data for 2 hours; S2. For the data at each radar volume scan time: convert the radar reflectivity data from the polar coordinate system into three-dimensional data in the Cartesian rectangular coordinate system; perform deblurring processing on the radial velocity data and map the radial velocity data into the radial rectangular coordinate system; use the reflectivity and radial velocity images within the first hour as inputs to the extrapolation model, and use the reflectivity images within the second hour as labels to construct a sample set; S3, using the mid-layer radial convergence automatic identification algorithm to obtain the mid-layer radial convergence feature information in the radial velocity map; S4, build a radar echo extrapolation model that takes the three-dimensional radar reflectivity and radial velocity images within one hour as input and the reflectivity images of the following hour as output; S5. Introduce deformation convolution into the radar echo extrapolation model to complete the automatic extraction of the dynamic field information of the convective system; S6. Dividing the sample set into a training set, a validation set, and a test set, training the radar echo extrapolation model on the training set, performing hyperparameter tuning on the validation set, and testing the radar echo extrapolation model on the test set; S7. Using the radar reflectivity data within any hour, convert the radar reflectivity data at each radar scanning moment from the polar coordinate system into three-dimensional data under the Cartesian rectangular coordinate system using the steps of S2, deblur the radial velocity data at each elevation angle, map the radial velocity data to the radial rectangular coordinate system, and use the above two types of data within the hour as inputs to the model to obtain the reflectivity image for the next hour.

2. The three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint according to claim 1 is characterized in that: The S2 converts the radar reflectivity data at each radar scanning moment from the polar coordinate system to three-dimensional data in the Cartesian rectangular coordinate system by first converting the radar base data from the polar coordinate system to the Cartesian coordinate system using a bilinear interpolation method in the horizontal dimension, and then interpolating in the vertical dimension to obtain three-dimensional equidistant grid field data; The step S2 of performing deblurring on the radial velocity data at each elevation angle refers to using a deblurring algorithm based on minimizing the velocity difference between regions to solve the velocity ambiguity problem, and mapping the deblurred radial velocity data to a radial rectangular coordinate system.

3. The three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint according to claim 1 is characterized in that: The S4 adopts the strategy of "independent encoding followed by fusion" to process the three-dimensional radar echo images of ten consecutive moments in the past hour. The reflectivity data and the mid-layer radial convergence structure information are respectively obtained through their own encoders to obtain feature mappings. The decoder integrates the two types of feature information and decodes and generates the future extrapolated echo sequence based on it. The radar echo extrapolation model completes the effective extrapolation of the three-dimensional radar echo data by jointly using the ConvLSTM network based on 3D convolution and the convolutional neural network CNN. The ConvLSTM network based on 3D convolution replaces the two-dimensional convolution operation of the ConvLSTM unit with a three-dimensional convolution operation to process the three-dimensional radar data corresponding to each radar scanning moment, and uses the unique loop structure of the ConvLSTM unit to extract and integrate the features of the three-dimensional radar data at different times within an hour; Each ConvLSTM unit contains three key components: forget gate, input gate and output gate. The unit receives the data input X at the current time point. t , the state output H at the previous time point t-1 and long-term memory unit C t-1 As input, it generates an updated long-term memory unit C t and the state output H at the current time point t As output, the operations performed by the ConvLSTM cell and the calculation of its output are shown below: i t =σ(W xi *X t +W hi *H t-1 +W ci *C t-1 +b i ) f t =σ(W xf *X t +W hf *H t-1 +W cf *C t-1 +b f ) o t =σ(W xo *X t +W ho *H t-1 +W co *C t-1 +b o ) In the formula, * represents the convolution operation, represents the Hadamard product, σ is the activation function, W and b represent the corresponding convolution kernel weight and bias term, tanh is the hyperbolic tangent function, i t 、f t and t They represent input gate, forget gate and output gate respectively.

4. The three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint according to claim 1 is characterized in that: The deformable convolution of S5 refers to adjusting the convolution kernel shape according to actual needs, further expanding the extraction and representation of effective speed features and their associated information. The two modes of the deformable convolution are: (1) Based on the maximum statistical distance between positive and negative velocity points, the height K of the convolution kernel in each row of the deformed convolution is H The width of the convolution kernel is fixed to 15. W It is determined dynamically by the following formula according to its specific radial position in the image; Where f represents the distance from the center of the convolution kernel to the top of the radial velocity image; In order to achieve the maximum speed difference dv max and the average speed difference dv mean The quantization of the feature parameters included, the weight distribution inside the convolution kernel is determined by the following formula; In the formula, i represents the distance between the pixel point and the top of the template, j represents the distance between the pixel point and the left side of the template, (i, j) is the coordinate position of the pixel point in the template, Ω is the relevant area for feature calculation, and K is l (i, j) is the weight at the coordinate (i, j); (2) Using the template matching method, we implement speed matching and calculation of different zones, limit the distance between the speed extreme points, and evaluate the relative position between the center of the maximum inflow speed and the center of the outflow speed, so as to accurately identify the key speed features in the airflow field. The airflow field speed feature identification template is in the shape of a regular hexagon, which is specifically used to identify four types of airflow field structures: convergence, divergence, cyclone, and anticyclone. Based on the consideration of statistical distance, the height K of the convolution kernel in each row of the deformed convolution used to determine the rotation type is H The width of the convolution kernel is fixed to 11. W It is determined dynamically by the following formula according to its specific position in the image; Where h∈{-5,-4,-3,...,3,4,5}.

5. The three-dimensional radar echo extrapolation method based on deformation convolution dynamic field constraint according to claim 1 is characterized in that: The S6 trains the radar echo extrapolation model using a weighted mean square error (WMSE) as a model loss function to guide the model to pay more attention to the prediction accuracy of the strong echo area. The WMSE is defined as follows: In the formula, I h，x，y Indicates the reflectivity at the coordinates (h, x, y) in the label image. Represents the reflectivity at (h, x, y) in the predicted image, Weight h，x，y The weight value representing the position can be set according to the different echo strengths to make the model pay more attention to the strong echo area.

Images