Thunderstorm gale classification forecasting method and system

Abstract

The invention belongs to the technical field of data processing, and particularly relates to a thunderstorm and gale classification forecasting method and system.The method comprises the steps that radar observation data, numerical mode output and ground observation data are obtained, and a unified grid multi-source data set is obtained through preprocessing and space-time registration; extracting radar morphological characteristics through a convolutional neural network, and identifying a thunderstorm and gale overall risk area in combination with environmental factors output in a numerical mode; then extracting local radar features in the region, estimating the cold pool strength and obtaining a time sequence mutation signal; and finally, inputting the three types of information into a classification model, outputting two types of risk judgment results of common thunderstorm gale and downburst, and presenting the results in probability and uncertainty intervals. According to the method, the downburst can be accurately identified, and the refinement level of medium and short time thunderstorm gale forecast is improved.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of data processing, and particularly relates to a thunderstorm gale classification and prediction method and system. BACKGROUND

[0002] Thunderstorm gale is an important disastrous phenomenon of severe convective weather, in which the downburst has the characteristics of strong burst and great destructive power, and often causes serious threats to airport operation, port operation and urban infrastructure.

[0003] The existing short-time thunderstorm gale prediction method is mostly based on radar extrapolation or numerical model results, and it is difficult to accurately identify special types of gales such as downbursts in a short time scale. The conventional method often makes a whole discrimination on storm gale, lacks special prediction means for key subtypes, and leads to generalization of the prediction result and insufficient accuracy.

[0004] Therefore, there is an urgent need for a prediction method capable of effectively identifying and predicting downbursts in a medium-short time scale, so as to improve the pertinence and fine level of disastrous weather defense. SUMMARY

[0005] To this end, the application provides a thunderstorm gale classification and prediction method and system to solve the above technical problems.

[0006] The application provides a thunderstorm gale classification and prediction method, comprising the following method steps:

[0007] Obtaining radar observation data, numerical model output and ground observation data, pre-processing and time and space registration of various types of data to obtain a unified grid multi-source data set;

[0008] Extracting radar shape features from the radar observation data through a convolutional neural network, the radar shape features combining environmental factors based on the numerical model output to identify the overall risk area of thunderstorm gale;

[0009] In the overall risk area, extracting local radar features, estimating cold pool intensity and time series mutation signal;

[0010] Inputting the local radar features, cold pool intensity and time series mutation signal into a classification model to output a downburst risk discrimination result, the risk discrimination result including two types of ordinary thunderstorm gale and downburst, and the output form being probability and uncertainty interval.

[0011] Another aspect of the application also provides a thunderstorm gale classification and prediction system, comprising:

[0012] A multi-source dataset acquisition module is configured to acquire radar observation data, numerical model output and ground observation data, pre-process the various types of data and perform time and space registration to obtain a multi-source dataset of a unified grid;

[0013] A general risk area identification module is configured to extract radar pattern features from the radar observation data by a convolutional neural network, and identify a general risk area of a thunderstorm gale by combining the radar pattern features with environmental factors based on the numerical model output;

[0014] A feature extraction module is configured to extract local radar features, estimate cold pool intensity and capture time series mutation signals in the general risk area;

[0015] A downburst risk identification module is configured to input the local radar features, cold pool intensity and time series mutation signals into a classification model to output a downburst risk discrimination result, the risk discrimination result including two categories of ordinary thunderstorm gales and downbursts, and the output form being a probability and an uncertainty interval.

[0016] Another aspect of the present application also provides an electronic device, which comprises at least one processor and a memory connected in communication with the at least one processor; wherein the memory stores instructions executable by the at least one processor, and the instructions are executed by the at least one processor to enable the at least one processor to perform the thunderstorm gale classification and prediction method as described above.

[0017] Another aspect of the present application provides a computer-readable storage medium having stored thereon computer program instructions executable by a processor to implement the thunderstorm gale classification and prediction method as described above.

[0018] Another aspect of the present application provides a computer program product comprising a computer program, which, when executed by a processor, implements the thunderstorm gale classification and prediction method as described above.

[0019] The present application first extracts radar pattern features by a 3D-CNN and combines environmental factors, and then predicts a general risk area of a thunderstorm gale by a spatio-temporal Transformer model to realize the first-stage precise positioning from global screening to risk focusing; and then extracts local radar features, estimates cold pool intensity and captures time series mutation signals in the risk area, and outputs the probability and uncertainty interval of ordinary thunderstorm gales and downbursts by an XGBoost classifier containing a full connection layer to complete the second-stage fine discrimination. The present application breaks through the limitations of traditional methods in poor generalization and difficulty in distinguishing downbursts, greatly improves the fine and accuracy of medium and short-term prediction, and provides precise support for thunderstorm gale disaster prevention. BRIEF DESCRIPTION OF DRAWINGS

[0020] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following will briefly introduce the drawings needed to be used in the embodiments. It should be understood that the following drawings only show some of the embodiments of the present application, and therefore should not be regarded as a limitation on the scope. For those skilled in the art, other related drawings can also be obtained without creative labor by reading the following detailed description of the non-limiting embodiments with reference to the following drawings:

[0021] Other features, objects and advantages of the present application will become more apparent from the following detailed description of non-limiting embodiments with reference to the drawings:

[0022] Figure 1 A flow chart of a thunderstorm gale classification and prediction method provided for the embodiments of the present application.

[0023] Figure 2 A schematic diagram of a Transformer model architecture provided for the embodiments of the present application.

[0024] Figure 3 A schematic diagram of an LSTM encoder architecture provided for the embodiments of the present application.

[0025] Figure 4 A schematic diagram of a thunderstorm gale classification and prediction system provided for the embodiments of the present application.

[0026] Figure 5 A schematic diagram of an electronic device provided for the embodiments of the present application. DETAILED DESCRIPTION

[0027] In order to make the purposes, technical solutions and advantages of the embodiments of the present application more clear, the following will clearly and completely describe the technical solutions in the embodiments of the present application with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only some of the embodiments of the present application, but not all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the scope of protection of the present application.

[0028] As shown in Figure 1 The embodiments of the present application disclose a thunderstorm gale classification and prediction method 100, comprising the following method steps:

[0029] S101, obtaining radar observation data, numerical model output and ground observation data, pre-processing each type of data and performing time and space registration to obtain a unified grid multi-source data set;

[0030] S102, performing feature extraction on the radar observation data by a convolutional neural network to obtain radar shape features, the radar shape features being combined with environmental factors based on the numerical model output to identify an overall risk area of thunderstorm gale;

[0031] S103, extracting local radar features, estimating cold pool intensity, and time series mutation signals in the overall risk area;

[0032] S104, inputting the local radar features, cold pool intensity, and time series mutation signals into a classification model to output a downburst risk discrimination result, the risk discrimination result including two categories of ordinary thunderstorm gale and downburst, and the output form being probability and uncertainty interval.

[0033] In some embodiments, for S101, for radar observation data of a target area, a volume scan mode (VCP21 or the like) is adopted, and the radar observation data may, for example, include the following parameters without limitation:

[0034] Reflectivity factor (Z, unit: dBZ), reflecting the concentration and phase distribution of precipitation particles;

[0035] Radial velocity (Vr, unit: m / s), revealing the direction and intensity of airflow movement inside the storm;

[0036] Spectrum width (SW, unit: m / s), representing the intensity of turbulence and the degree of wind shear;

[0037] Echo top height (ET, unit: km), reflecting the vertical development intensity of the convective cloud body.

[0038] The data original resolution is 1 km x 1 km horizontally and 5 minutes in time interval.

[0039] It can be understood that, due to problems such as occlusion, ground clutter, non-meteorological echoes (such as insects and birds), and range folding of radar scanning, the radar observation data is preprocessed, which may, for example, include:

[0040] Clutter filtering, using a dual-polarization parameter (if available) or a spectrum width-reflectivity joint threshold method to remove non-meteorological echoes;

[0041] Radial velocity deblocking, using a maximum likelihood method or a neighboring volume scan interpolation method to process velocity folding;

[0042] Interpolation to a uniform grid, interpolating all radar variables to a uniform three-dimensional grid (resolution: 1 km x 1 km x 0.5 km), and the coverage area is usually a circular area with a radius of 200 km;

[0043] Time smoothing, using a sliding average (window = 5) to eliminate noise fluctuations for 6 consecutive frames (30 minutes) of radar data, retain the trend evolution characteristics, and retain the mutation signal for subsequent analysis.

[0044] In one embodiment, for the numerical model output data, a high-resolution regional numerical prediction model (such as WRFv4.3, CMA-MESO3.0 or ECMWF IFS-HRES2.5km product) is selected, and the output frequency is once per hour. Exemplarily, the numerical model output data can include, but is not limited to, the following parameters:

[0045] Temperature (T, K), specific humidity (q, kg / kg), three-dimensional wind field (u, v, w, m / s), vertical velocity (ω, Pa / s), convective available potential energy (CAPE, J / kg), which measures the atmospheric unstable energy, lifting condensation level (LCL, m), which reflects the low-level water vapor saturation degree, 0-6 km vertical wind shear (VS06, m / s), which is calculated as the modulus of the difference between the wind vector at 6 km height and the wind vector near the ground (0-1 km); boundary layer relative humidity (RH BL, %), exemplarily, the average value from 925 hPa to 850 hPa; temperature advection (T_adv, K / h) and humidity advection (q_adv, g / kg·h), which are calculated based on the wind field and the gradient operator, which will not be described here.

[0046] In this embodiment, in order to match the radar grid, the model data is subjected to bilinear interpolation and vertical interpolation processing. Specifically, the 10-30 km resolution data output by the model is interpolated to a 1 km horizontal grid.

[0047] Vertical level resampling is performed using the terrain-following coordinate system (eta coordinate), so that the number of vertical layers is consistent with the radar inversion layer (a total of 20 layers, from the ground to 15 km).

[0048] All variables are aligned to the integral multiple time points of the radar data according to the time stamp (i.e., interpolated every 5 minutes, and the intermediate period is completed by linear interpolation), so as to ensure that the radar and the model data are synchronized in time.

[0049] In one embodiment, for the ground observation data, high-frequency (1 minute interval) observations are provided. Exemplarily, the observation data can include, but is not limited to, the following parameters:

[0050] Air temperature (T_gnd, ℃), dew point temperature (Td, ℃), wind speed (Ws, m / s), wind direction (Wd, °), air pressure (P, hPa).

[0051] Then, the observation data is preprocessed, specifically including:

[0052] Outlier rejection, 3σ criterion and time consistency test (such as wind speed sudden change > 15 m / s / 30 s is considered as an abnormal value);

[0053] Spatial Interpolation, using Kriging or Inverse Distance Weighting (IDW) to interpolate discrete station data to each pixel of the radar grid, generating continuous surface fields (temperature, humidity, wind speed);

[0054] Cold Pool Monitoring Enhancement, calculating the surface temperature drop rate ΔT / Δt (temperature drop slope in the past 15 minutes) at each grid point as an indicator of cold pool development.

[0055] So far, a unified grid multi-source dataset is obtained, forming a set of spatio-temporal unified four-dimensional data cubes (x, y, z, t), and the dimensions are exemplarily:

[0056] Horizontal: 1 km x 1 km, covering an area of about 200 km x 200 km;

[0057] Vertical: 20 layers (0-15 km);

[0058] Time: 1 frame every 5 minutes, historical window 60 minutes (12 frames), future prediction window 0-120 minutes (24 frames).

[0059] In some embodiments, for S102, based on the unified grid data, this embodiment uses a three-dimensional convolutional neural network (3D-CNN) to perform end-to-end feature extraction on the radar data, capturing the morphological features of the storm structure.

[0060] Specifically, as an example, the 3D-CNN architecture design includes,

[0061] Input layer: 12 frames x 4 variables (Z, Vr, SW, ET) x 20 layers x 200 x 200, where the dimensions are (12, 4, 20, 200, 200);

[0062] First layer: 3D convolution kernel (3x3x3), 64 channels, ReLU activation, Batch Normalization (BatchNorm);

[0063] Second layer: 3D convolution kernel (3x3x3), 128 channels, max pooling (2x2x2);

[0064] Third layer: 3D convolution kernel (3x3x3), 256 channels, Dropout = 0.3;

[0065] Fourth layer: Global Average Pooling (GAP), outputting a 128-dimensional feature vector.

[0066] Based on the above architecture model, the extracted radar morphological features include but are not limited to:

[0067] Bow Echo Intensity Index, calculating the ratio of echo front curvature to rear echo gap area by contour fitting;

[0068] Low-level radial divergence (LRD), the divergence of the radial velocity gradient tensor at 1.5 km height, quantifies the outflow burst strength;

[0069] Echo top drop signal (ETDrop), the maximum drop of ET in the past 15 minutes, if >2 km with Z >45 dBZ, is marked as a strong downdraft trigger signal;

[0070] Spectral width anomaly (SWAnomaly), the region with SW >8 m / s and weak Z correlation at the storm rear (<2 km), indicates strong turbulent mixing.

[0071] It can be understood that these morphological features are encoded as a 128-dimensional dense vector, denoted as F_radar∈R 128 , representing the current structure fingerprint of the storm.

[0072] In practical applications, simply relying on radar morphology is prone to "false positive" interference (such as strong squall lines not necessarily producing downbursts), so this embodiment introduces environmental background factors to reflect atmospheric thermal stability and dynamic forcing conditions.

[0073] Specifically, key environmental factors CAPE, VS06, RH_BL, T_adv&q_adv are extracted from numerical model output, the above factors are extracted at each grid point of the radar grid to form a spatial field, and then compressed into a multi-dimensional vector (e.g. 16-dimensional) F_env∈R 16 It can be understood that the dimensions extracted can be set as needed, and this embodiment is not limited.

[0074] Then, the radar morphology features and environmental features are spliced and fused in a unified vector space:

[0075] F_fused=Concat(F_radar,F_env), forming a fused feature vector of length 144.

[0076] To predict the probability of each grid point occurring thunderstorm gale (≥17 m / s) in the next 0-120 minutes, this embodiment uses a spatiotemporal Transformer model instead of traditional extrapolation or statistical regression models.

[0077] Specifically, as an example, as shown in Figure 2 , the Transformer model architecture includes:

[0078] Input layer: historical 12 frames (60 minutes) of fused feature vector sequence F_fused(t-11),..., F_fused(t);

[0079] Position encoding: sinusoidal function to embed time dimension, giving each frame an explicit temporal position;

[0080] Multi-head self-attention (MHA): 8 heads, each head dimension 18, total embedding dimension 144;

[0081] Feed-forward network: two-layer MLP, hidden layer 256, ReLU activation;

[0082] Output layer: sigmoid activation, output a single scalar, the probability P_wind(t+Δt) of the occurrence of thunderstorm gale at this grid at future time t+Δt(Δt=5,10,15,…,120min).

[0083] Optionally, the model training uses a weighted binary cross-entropy loss function to solve the problem of unbalanced positive and negative samples (the proportion of downburst samples is <10%).

[0084] The model output is a probability field with a spatial resolution of 1 km and a time step of 5 minutes, which can generate dynamic probability evolution maps within 2 hours in the future. By setting a threshold (such as P≥60%), high-risk areas can be screened out, i.e., the overall risk area.

[0085] Specifically, the input layer is based on a sequence of 60-minute fused feature vectors, ensuring that the data covers radar morphology and environmental factor information, and the format is adapted to the time series input requirements of the Transformer model. For example, the specific steps are as follows:

[0086] F_fused is formed by concatenating the radar morphology feature vector (F_radar) and the environmental factor feature vector (F_env). For F_radar, it is extracted by 3D-CNN from historical 12 frames (1 frame every 5 minutes, a total of 60 minutes) of radar data (reflectivity factor Z, radial velocity Vr, spectral width SW, echo top height ET), outputting a 128-dimensional vector, covering key morphological features such as arc echo intensity index, low-level radial divergence intensity, and echo top drop signal.

[0087] For F_env, it is extracted from numerical model output, such as convective available potential energy (CAPE), 0-6 km vertical wind shear (VS06), boundary layer relative humidity (RH_BL), temperature advection, and humidity advection, and compressed into a 16-dimensional vector through spatial global averaging and local maximum pooling.

[0088] Concatenation and fusion, such as concatenating F_radar (128-dimensional) and F_env (16-dimensional) in a unified vector space to form F_fused with a dimension of 144, i.e., F_fused∈R 144 , ensuring that a single frame of features contains both storm structure and atmospheric background information.

[0089] For each 1km x 1km grid point, the F_fused of the past 12 frames (t-11 to t moment, a total of 60 minutes) is extracted to form a time sequence {F_fused(t-11), F_fused(t-10),..., F_fused(t)}, and the sequence length is 12;

[0090] The time sequence is converted into a Transformer input format, for example, the dimension is sequence length x feature dimension, that is, 12 x 144, where the sequence length corresponds to 12 time steps, and the feature dimension corresponds to 144-dimensional fusion features, which meets the requirement of the model for two-dimensional input of time sequence and features.

[0091] At the same time, the F_fused sequence of all grid points is Z-score standardized, and the formula is: x'=(x-μ) / σ, where μ is the mean of all training samples F_fused, and σ is the standard deviation; After standardization, the feature mean is 0 and the standard deviation is 1, which avoids the model from overemphasizing certain features due to the difference in magnitude between radar shape features (such as echo top height, unit km) and environmental factors (such as CAPE, unit J / kg), and improves the convergence speed of the model.

[0092] Specifically, the time dimension information is embedded through position encoding, and each frame of feature is given a clear time position information through position encoding, and the embodiment adopts a sine function to achieve this.

[0093] For example, for the pos-th time step (pos=0 to 11, corresponding to t-11 to t moment) in the time sequence and the i-th feature dimension (i=0 to 143), the calculation of the position encoding value PE(pos, i) is divided into two cases:

[0094] When i is even: PE(pos, i)=sin(pos / 10000^(2i / 144))

[0095] When i is odd: PE(pos, i)=cos(pos / 10000^(2i / 144))

[0096] In the formula, "10000^(2i / 144)" is a scaling factor to ensure that the position encoding periods of different feature dimensions are significantly different, so that the model can distinguish between near and long-term time steps (such as t moment and t-11 moment); The alternative use of "sin" and "cos" functions can capture the periodicity and relative position relationship between time steps, and adapt to the time span of 60-minute historical sequence.

[0097] The calculated position encoding matrix (dimension 12x144) is added to the pre-processed fused feature sequence matrix (12x144) element by element to obtain the input matrix X containing time position information, i.e. X = F_fused_sequence + PE;

[0098] The embedding matrix dimension remains unchanged at 12x144, ensuring that the subsequent multi-head self-attention layer can directly receive the input without additional dimension conversion.

[0099] Multi-head self-attention is used to mine the time dependence of the fused feature sequence within the last 60 minutes. For example, 8 heads are set, with a total embedding dimension of 144, and the implementation is as follows:

[0100] S201, attention head splitting and linear mapping, including,

[0101] The total embedding dimension of 144 is evenly distributed among the 8 attention heads, with each head having a dimension of 18 (144 ÷ 8 = 18);

[0102] The input matrix X containing position encoding (12x144) is mapped to the query matrix Q, the key matrix K, and the value matrix V through three independent linear layers, respectively. After mapping, the dimensions of Q, K, and V are all 12x144;

[0103] Q, K, and V are split by attention head, with each head corresponding to Q_head (12x18), K_head (12x18), and V_head (12x18), forming a total of 8 groups of head-level matrices.

[0104] S202, single-head attention calculation, for each group of head-level matrices, the attention weight and output are calculated according to the following steps:

[0105] Similarity calculation: the similarity between Q_head and the transpose matrix of K_head is calculated through dot product operation, formula: Sim = Q_head x K_head^T / sqrt(d_k), where d_k = 18 (single-head dimension), and dividing by sqrt(d_k) can avoid the gradient disappearance after softmax due to too large similarity value;

[0106] Attention weight normalization: apply the softmax function to the similarity matrix Sim (12x12, row corresponding to query time step, column corresponding to key time step) to obtain the attention weight matrix W (12x12), with weight value range [0, 1], the higher the weight of a certain position, the greater the influence of the key time step feature on the query time step;

[0107] Attention output calculation, multiply the weight matrix W with the value matrix V_head, get the single-head attention output O_head (12x18), the formula is: O_head = W x V_head.

[0108] Multi-head attention fusion, for example, the outputs of the 8 attention heads O_head (each 12x18) are spliced by column to form the spliced matrix O_concat (12x144, 18x8 = 144);

[0109] Dimension conversion and feature fusion of O_concat through linear layer, output the final result MHA_out (12x144) of the multi-head self-attention layer, which has captured the time-dependent relationship of the fused features in the past 60 minutes.

[0110] Specifically, the feedforward network is used to perform nonlinear conversion on the features of the multi-head self-attention output, to improve the model's fitting ability for complex meteorological laws. For example, two layers of MLP are set, with a hidden layer of 256 and ReLU activation. The specific implementation is as follows:

[0111] The first layer of MLP is used for feature dimension lifting and nonlinear activation, for example,

[0112] The input is MHA_out (12x144), and the feature dimension is lifted from 144 to 256 through a linear layer, with an output dimension of 12x256;

[0113] The ReLU activation function (f(x) = max(0, x)) is applied to suppress negative feature values and enhance the model's expression ability for nonlinear relationships, while avoiding the problem of gradient disappearance;

[0114] The Dropout layer (dropout_rate = 0.3) is introduced to randomly mask 30% of the neuron outputs, reduce the risk of overfitting, and ensure the model's generalization ability on different grid data.

[0115] The second layer of MLP is used for feature dimension reduction and output preparation, for example,

[0116] The 256-dimensional features are reduced to 144-dimensional features through a linear layer, keeping the same dimension as the input MHA_out (12x144), to ensure the dimension adaptation of the subsequent residual connection;

[0117] Optionally, no activation function and Dropout are introduced to preserve the integrity of the features and prepare for feature fusion after the residual connection.

[0118] Residual connection and layer normalization, for example, the output of the feedforward network and the input of the multi-head self-attention layer MHA_out are connected in residual connection (i.e. FFN_out = FFN(MHA_out) + MHA_out), which alleviates the problem of gradient decay of deep models;

[0119] Layer normalization is applied to the result after residual connection, the formula is: LN(x) = γ × (x - μ_x) / σ_x + β, where μ_x is the feature mean, σ_x is the feature standard deviation, γ and β are learnable parameters, which ensures the stable distribution of the input features of each layer and accelerates the model training.

[0120] After the model training is completed, the future 0-120 minute grid point thunderstorm wind probability is output through forward propagation, and the high-risk area is screened, and the specific implementation is as follows:

[0121] For the input historical 12 frame fusion feature sequence, the model is propagated forward, and the Sigmoid activation function is applied at the output layer to compress the prediction result to the [0, 1] interval, to obtain the thunderstorm wind probability P_wind(t+Δt) of each time step (Δt = 5, 10,..., 120 min) of the grid point in the future;

[0122] For all 1km×1km grid points, P_wind of 24 time steps in the future is calculated according to the same process to form a probability field (dimension: grid point number × 24) of “space × time”;

[0123] Based on the preset threshold (for example, P_wind≥60%), the grid points with probability ≥0.6 in each time step in the future are screened, and all grid points that meet the condition constitute the overall risk area of thunderstorm wind.

[0124] In some embodiments, for S103, in the overall risk area, further extract local fine features with downburst specificity and high discriminability to constitute the input of the classification model.

[0125] Specifically, the local radar morphology feature extraction is performed on the 5×5 neighborhood (5km×5km) around the high-risk grid point as follows:

[0126] The arch echo integrity index is defined as the ratio of the arch front arc length to the complete circumference length, with a value range of 0-1, and the ratio of the arch front arc length to the complete circumference length is calculated, and the downburst corresponds to 0.6-1.0; the downburst often shows strong arch and rear side notch structure;

[0127] Rear inflow notch (RIP) recognition, in the reflectivity image, the concave area of the echo is detected by using morphological opening operation, and the area and maximum depth thereof are calculated; the notch area, unit km 2, a "U-shaped" notch >3km 2 in depth, unit dBZ, the maximum depth of the notch corresponding to >15dBZ echo depression of the downburst;

[0128] Low-level strong divergence core, at 1.0-1.5km height layer, calculate radial velocity divergence field, select the cluster center of divergence <-0.01s-1 as the outflow burst point;

[0129] Echo top drop rate, calculate the maximum drop rate of ET in 10 minutes (ΔET / Δt), if >0.5km / min and accompanied by Z>50dBZ, it is a strong downdraft trigger sign.

[0130] The essence of downburst is that strong downdraft flow impacts the ground to form a cold pool, and its strength directly determines the peak ground wind speed. The present embodiment proposes a double-source cold pool strength evaluation model.

[0131] Specifically, calculate the cold pool temperature difference ΔT_cold, for example, within a radius of 10km around the center point of the high-risk area, calculate the difference between the observed temperature of the ground station and the background temperature output by the model (without storm influence), and take the minimum value as the cold pool core cooling amplitude, if ΔT_cold≤-4℃, it is considered that the cold pool has formed strong cooling effect.

[0132] Calculate the cold pool horizontal extension radius R_cold, for example, use the low-level (925hPa) temperature gradient field Calculate the temperature contour envelope area, convert to equivalent circle radius, the R_cold corresponding to the downburst is usually >15km, while the ordinary gust front is only 5-10km.

[0133] Before the occurrence of downburst, the boundary layer often appears sharp drop of isentropic potential temperature (θ e ) and specific humidity (q), which marks the arrival of downdraft flow to the ground, the present embodiment extracts the time series of isentropic potential temperature (θ e ) and specific humidity (q) from ground observation and model output, for each high-risk grid point, construct θ e (t) and q(t) sequences in the past 60 minutes, input bidirectional LSTM encoder and attention mechanism model. Among them, θ e (t) is calculated in real time from the temperature (T), pressure (P), and specific humidity (q) output by the numerical model at the grid point;

[0134] Among them, the LSTM encoder is double-layer, as an example, each layer has 128 units;

[0135] Attention module, used to calculate the weight α_t of each time step, highlighting the mutation moment;

[0136] Mutation point detection can be, when α_t>0.8 and Δθ e / Δt<-0.5K / min and Δq / Δt<-0.3g / kg / min are simultaneously true, it is determined that a strong sinking arrival signal is reached; the model can capture the cold pool impact precursor in advance, which is significantly better than the traditional threshold method.

[0137] Specifically, θ e is calculated based on the grid data (temperature T, unit K; air pressure P, unit hPa; specific humidity q, unit kg / kg) output by the numerical model, and the time series data of q is extracted from the grid q output by the numerical model and the high-frequency q (1 minute interval, unit g / kg) of the ground observation station, which is mapped to the high-risk grid of the radar unified grid through Kriging interpolation;

[0138] θ e is calculated by using the standard formula of atmospheric thermodynamics, and the formula is:

[0139] θ e = T x (1000 / P) ^ (0.286) x exp[(3.376 / T-0.00254) x q x 1000 / (1+0.61 x q)]

[0140] Wherein, T is the grid temperature (K), P is the grid air pressure (hPa), and q is the grid specific humidity (kg / kg); during the calculation process, the T, P and q output by the numerical model have completed spatial registration (interpolation to 1km x 1km x 0.5km unified grid) and time alignment (matching 5 minute time step), so as to avoid θ e calculation error caused by data misplacement.

[0141] For the θ e sequence and the q sequence of each high-risk grid, Z-score standardization is used to eliminate the dimensional difference, and the formula is: x'=(x-μ) / σ, wherein, μ is the mean of the past 60 minutes sequence of the grid, and σ is the standard deviation of the sequence; the mean of the standardized data is 0, and the standard deviation is 1, which can avoid the excessive bias of the model to one variable caused by the magnitude difference between θ e (unit K) and q (unit g / kg), and improve the convergence speed of the LSTM model.

[0142] The standardized θ e sequence and the q sequence are spliced into a two-channel time series vector, and the input dimension of each grid is "time step x feature number", that is, 12 (time step) x 2 (θ e , q).

[0143] The bidirectional LSTM encoder is used to capture θ eThe long-term time dependence of the q sequence (including the evolution trend of the past time sequence and the possible change characteristics in the future) is enhanced through the double-layer structure to improve the feature extraction capability.

[0144] Specifically, as shown in the figure, Figure 3 the architecture of the LSTM encoder can be:

[0145] The input layer receives a 12x2 two-channel time sequence vector, and the input of each time step is the normalized θ e and q value at this moment;

[0146] The bidirectional LSTM layer (first layer) is set to have 128 hidden units, and uses forward LSTM (from t0to t 11 iterates the sequence to capture the forward time dependence) and backward LSTM (from t 11 iterates the sequence to capture the reverse time dependence) to calculate in parallel, and the output dimension is 12x256 (time step x 2 x hidden unit number, 256=128x2); the activation function of this layer uses tanh, and the forgetting gate threshold is set to the default value. Through the Dropout layer (dropout_rate=0.3), 30% of the neuron outputs are randomly shielded to avoid overfitting;

[0147] The bidirectional LSTM layer (second layer) has an input dimension of 12x256 and is also set to have 128 hidden units. The forward and backward LSTMs are calculated in parallel, and the output dimension is still 12x256. This layer does not introduce Dropout to ensure that key time sequence dependence features are retained and provide high-quality feature vectors for the subsequent attention module;

[0148] The output layer concatenates the time step dimension of the output of the second bidirectional LSTM, and each time step corresponds to a 256-dimensional feature vector (combining the outputs of the forward and backward LSTMs). Finally, a 12x256 time sequence feature matrix is formed, which is used as input to the attention module.

[0149] The attention module is used to assign weights to the 12 time step feature vectors output by the bidirectional LSTM, highlighting the moments corresponding to the sharp drop of θ e and q, which include:

[0150] The additive attention (AdditiveAttention) mechanism is used to learn the similarity between each time step feature vector and the query vector (Query) to calculate the weight α_t of this time step, which includes feature mapping, similarity calculation, and softmax normalization.

[0151] Specifically, the feature mapping layer maps the 256-dimensional time step feature vector output by the bidirectional LSTM and the learnable query vector (dimension 256) to a 128-dimensional space respectively through a linear layer, and the mapping formula is:

[0152] V = W_v x h_t + b_v, Q = W_q x q + b_q

[0153] where h_t is a 256-dimensional feature vector at the t-th time step, q is a query vector, W_v and W_q are linear layer weight matrices, and b_v and b_q are bias terms.

[0154] The similarity calculation layer calculates the similarity between the mapped feature vector and the query vector through a tanh activation function, and the formula is:

[0155] score(t) = tanh(V + Q) x W_s

[0156] where W_s is a 128-dimensional to 1-dimensional weight vector, and score(t) is the similarity score at the t-th time step.

[0157] The softmax normalization layer performs softmax normalization on the similarity scores of the 12 time steps to obtain the weight α_t of each time step, and the formula is:

[0158] α_t = exp(score(t)) / ∑(t = 0 to 11) exp(score(t))

[0159] After normalization, the value range of α_t is [0, 1], and the sum of α_t of all time steps is 1. The higher the weight, the more likely it is that the time step corresponds to θ e and the mutation moment of q.

[0160] Then, combined with the weight α_t output by the attention module and the change rate of θ e and q, the strong downdraft signal judgment of the downburst precursor is realized, and the specific process is as follows:

[0161] The calculation of θ e and q change rate is, for example,

[0162] For each high-risk grid, the θ e sequence and the q sequence, the change rate of adjacent time steps (interval 5 minutes) is calculated, and the formula is:

[0163] Δθ e / Δt = (θ e (t+1) - θ e (t)) / 5, Δq / Δt = (q(t+1) - q(t)) / 5

[0164] where, Δt = 5 minutes (time step interval), and the units of the change rate are K / min (θ e ) and g / kg min (q), respectively; the original data before normalization is used in the calculation process to ensure the physical meaning of the change rate is accurate (e.g., Δθ e / Δt <-0.5 K / min indicates that θ e drops more than 0.5 K per minute).

[0165] The mutation point determination rule is executed, for example,

[0166] When the following three conditions are met simultaneously for each time step t (t = 0 to 10, because the change rate of t and t+1 needs to be calculated), it is determined that a strong sinking arrival signal, i.e., a downburst precursor, has occurred:

[0167] The attention weight α_t>0.8 indicates that this time step is a moment of high weight mutation of the θ e and q sequence;

[0168] The θ e change rate Δθ e / Δt <-0.5 K / min indicates that θ e is in a sharp downward trend, which is consistent with the boundary layer cooling characteristics of a downburst;

[0169] The q change rate Δq / Δt <-0.3 g / kg min indicates that q is in a sharp downward trend, which is a sign of a sudden drop in humidity caused by the sinking of dry and cold air;

[0170] If there is at least one time step in the past 60 minutes that meets the above conditions, the grid point is marked as a potential downburst affected area, and the signal is used as one of the input features of the subsequent XGBoost classifier to assist in the classification and discrimination of downbursts and ordinary thunderstorm gales.

[0171] In this embodiment, the above three types of features, including local radar shape features, cold pool intensity parameters, and time series mutation signals, are combined into a discriminant feature vector and input into a pre-trained XGBoost classifier. The specific output of the classifier is P_DB∈[0,1], which is the probability of the occurrence of a downburst at the grid point.

[0172] Specifically, the above-mentioned quantified feature parameters are spliced in the order of local radar features (4 items), cold pool intensity (2 items), and time series mutation signals (4 items) to form a 10-dimensional initial feature vector. Z-score standardization (formula: x' = (x-μ) / σ, where μ is the mean of the training set features and σ is the standard deviation) is used to eliminate dimensional differences and ensure that the contribution of each feature to the model is balanced, and finally a 10-dimensional input vector after standardization is obtained.

[0173] In this embodiment, the XGBoost classifier is used as the core discriminant model, which has the advantages of high interpretability, strong anti-overfitting ability, and support for mixed type feature input.

[0174] Specifically, the learning rate (eta) is 0.1, which controls the contribution weight of each tree to avoid overfitting; the maximum tree depth (max_depth) is 6, which limits the complexity of a single tree to adapt to the fitting needs of 10-dimensional features; the sub-sample ratio (subsample) is 0.8, which randomly selects 80% of the samples for training each tree to enhance the generalization ability; the column sample ratio (colsample_bytree) is 0.8, which randomly selects 80% of the features for each tree to reduce the impact of feature redundancy; the number of weak learners (n_estimators) is 200, which improves the discriminant accuracy by integrating multiple trees.

[0175] The objective function is set as a binary classification task, with class 0 corresponding to ordinary thunderstorm gales and class 1 corresponding to downbursts; the loss function uses weighted binary cross-entropy (Weighted Binary Cross-Entropy), formula: Loss = -[w0·y0·log(p0)+w1·y1·log(p1)], where w0 and w1 are class weights, and w0 = 0.1 and w1 = 0.9 are set because the downburst sample proportion is <10%, balancing the training weights of positive and negative samples to avoid model bias towards the majority class (ordinary thunderstorm gales).

[0176] The standardized 10-dimensional feature vector is input into the trained XGBoost classifier, and the model predicts the downburst risk probability through the integration of multiple decision trees, with an example calculation as follows:

[0177] The prediction output of a single decision tree, each decision tree traverses from the root node based on the input feature vector, judges through feature threshold values (such as "whether the bow echo integrity index is >0.8" or "whether the cold pool temperature difference is ≤-4℃"), and finally falls into a leaf node, which corresponds to the probability value p_i (i is the tree number, 1≤i≤200) of class 1 (downburst).

[0178] The probability integration of multiple trees takes the average value of p_i output by the 200 trees using soft voting, and obtains the final downburst risk probability P_DB, formula: P_DB = (1 / 200) × ∑(i=1 to 200) p_i, with a value range of [0, 1].

[0179] When P_DB ≥ 0.6, it is determined as a high-risk downburst; when 0.3 < P_DB < 0.6, it is determined as a moderate-risk downburst; and when P_DB ≤ 0.3, it is determined as an ordinary thunderstorm gale, achieving clear differentiation between the two types of risks.

[0180] Optionally, in order to enhance the feature nonlinear fitting capability while retaining the XGBoost's original advantages of overfitting resistance and high interpretability, a fully connected layer is introduced between the XGBoost classifier and the input layer.

[0181] Specifically, two fully connected layers are set to realize feature dimensionality and nonlinear conversion, and the specific parameters are as follows:

[0182] The first fully connected layer has an input dimension of 10 (the dimension of the standardized feature vector), an output dimension of 64, and a ReLU (formula: f(x) = max(0, x)) activation function selected to solve the gradient disappearance problem. A Dropout layer (dropout_rate = 0.3) is also introduced to randomly mask 30% of the neurons to avoid overfitting.

[0183] The second fully connected layer has an input dimension of 64 (the output of the first layer), an output dimension of 16, and a ReLU activation function. No Dropout is introduced to retain key feature information. The final output is a 16-dimensional high-order feature vector, which is used as the input of the XGBoost classifier.

[0184] In addition, a weight regularization term (L2 regularization, formula: Loss_total = Loss_binary + λ × ∑(w 2 ), where λ = 0.01) is included in the loss calculation to constrain the weight scale of the fully connected layer and avoid parameter overfitting.

[0185] After introducing the fully connected layer, the calculation of the downburst risk probability includes, for example:

[0186] The standardized 10-dimensional original feature vector is input into the fully connected layer, and through two layers of ReLU activation and Dropout (Dropout is turned off during inference, and a weight scaling strategy is used, i.e., the Dropout rate of 0.3 during training is converted to a weight coefficient of 0.7 to ensure that the output magnitudes are consistent), a 16-dimensional high-order feature vector is output. This vector has fused nonlinear correlation information of the original features (such as the interactive influence of the bow echo integrity index and the cold pool temperature difference).

[0187] The prediction output of a single decision tree, as well as the probability integration and analogy determination of multiple trees, are the same as in the previous embodiments and will not be repeated here.

[0188] In this embodiment, optionally, to improve the reliability of the prediction results, the Monte Carlo Dropout (MC-Dropout) method is used to estimate the 95% confidence interval [P_low, P_high] to quantify the uncertainty of the prediction results. For example, the specific steps are as follows:

[0189] As in the previous embodiment, a Dropout layer is reserved in the fully connected layer of the XGBoost model, with dropout_rate = 0.3 (i.e. 30% of neurons are randomly turned off each time of forward propagation). For the same input feature vector, 50 predictions are repeated to obtain 50 samples of the downburst risk probability {P_DB 1 ,P_DB 2 ,...,P_DB 50}.

[0190] The 50 probability samples are sorted in ascending order, and the 2.5th percentile (corresponding to 50 x 2.5% = 1.25, rounded up to the 2nd sample) is taken as P_low, and the 97.5th percentile (corresponding to 50 x 97.5% = 48.75, rounded down to the 48th sample) is taken as P_high, forming a 95% confidence interval.

[0191] When the interval width (P_high-P_low) < 0.2, it indicates low uncertainty, and the forecast result is reliable; when 0.2 ≤ interval width < 0.4, it indicates moderate uncertainty, which needs to be further verified in combination with real-time radar observations; when the interval width ≥ 0.4, it indicates high uncertainty.

[0192] Figure 4 A thunderstorm wind classification forecasting system 400 is shown. The device embodiment corresponds to the method embodiment shown in Figure 1 , and specifically includes:

[0193] A multi-source data set acquisition module 401 is configured to acquire radar observation data, numerical model output, and ground observation data, pre-process each type of data, and perform time and space registration to obtain a unified grid multi-source data set;

[0194] A general risk area identification module 402 is configured to extract radar shape features from the radar observation data by a convolutional neural network, and identify a general risk area of thunderstorm wind by combining the radar shape features with environmental factors based on the numerical model output;

[0195] A feature extraction module 403 is configured to extract local radar features, estimate cold pool intensity, and time series mutation signals in the general risk area;

[0196] A downburst risk identification module 404 is configured to input the local radar features, cold pool intensity, and time series mutation signals into a classification model to output a downburst risk discrimination result, the risk discrimination result including two categories of ordinary thunderstorm wind and downburst, and the output form being probability and uncertainty interval.

[0197] The specific limitations of the thunderstorm wind classification and forecasting system can refer to the limitations of the thunderstorm wind classification and forecasting method described above, which will not be repeated here. Each module in the above thunderstorm wind classification and forecasting system can be realized by software, hardware and their combination. The above modules can be embedded in or independent of the processor in the electronic device in hardware format, or stored in the memory in the electronic device in software format, so that the processor calls the operation of each module.

[0198] It should be noted that, in order to highlight the innovative part of the present application, modules not closely related to solving the technical problems proposed in the present application are not introduced in the present embodiment, but this does not mean that there are no other modules in the present embodiment.

[0199] As shown in Figure 5 The electronic device 1 provided by the present application can include a memory 12, a processor 13 and a bus, and can also include a computer program stored in the memory 12 and executable on the processor 13, such as a thunderstorm wind classification and forecasting program.

[0200] The memory 12 includes at least one type of readable storage medium, including flash memory, mobile hard disk, multimedia card, card-type memory, magnetic memory, disk, optical disk, etc. The memory 12 can be an internal storage unit of the electronic device 1 in some embodiments, such as the mobile hard disk of the electronic device 1. The memory 12 can also be an external storage device of the electronic device 1 in other embodiments, such as the plug-in mobile hard disk, smart storage card, secure digital card, flash memory card, etc. equipped on the electronic device 1. Further, the memory 12 can include both the internal storage unit and the external storage device of the electronic device 1. The memory 12 can be used to store application software and various data installed on the electronic device 1.

[0201] The processor 13 executes the operating system of the electronic device 1 and various installed application programs. The processor 13 executes the application program to realize the steps in the above thunderstorm wind classification and forecasting method.

[0202] For example, the computer program can be divided into one or more modules, which are stored in the memory 12 and executed by the processor 13 to complete the present application. The one or more modules can be a series of computer program instruction segments that can complete a specific function, which are used to describe the execution process of the computer program in the electronic device 1.

[0203] The integrated unit implemented in the form of the software function module can be stored in a computer readable storage medium, which can be nonvolatile or volatile. The software function module is stored in a storage medium and includes a plurality of instructions for enabling a computer device (which can be a personal computer, a computer device, or a network device, etc.) or a processor to execute part of the functions of the embodiments of the present application.

[0204] The embodiments of the present application further provide a computer program product, which comprises a computer program. The computer program is executed by a processor to implement the steps in the above embodiments.

[0205] The above embodiments only illustrate the principles and effects of the present application, but are not used to limit the present application. Any modification or change made by any person skilled in the art without departing from the spirit and scope of the present application shall be covered by the claims of the present application.

Claims

1. A method for classifying and forecasting thunderstorms and strong winds, characterized in that, Includes the following steps: Acquire radar observation data, numerical model output and ground observation data, preprocess the various types of data and perform time and space registration to obtain a unified grid multi-source dataset; Radar morphology features are obtained by extracting features from the radar observation data using a convolutional neural network. These radar morphology features are then combined with environmental factors based on numerical model output to identify the overall risk area of ​​thunderstorms and strong winds. Within the overall risk area, local radar features are extracted, cold pool intensity is estimated, and time-series abrupt change signals are analyzed. The local radar features, cold pool intensity, and time-series abrupt change signals are input into the classification model, and the downburst risk discrimination results are output. The risk discrimination results include two categories: ordinary thunderstorm winds and downbursts. The output format is probability and uncertainty interval.

2. The method for classifying and forecasting thunderstorms and strong winds according to claim 1, characterized in that, The method for identifying the overall risk area of ​​thunderstorms and strong winds includes fusing radar morphological features and environmental factor features in a unified vector space to obtain fused features, and inputting the fused features into a Transformer model to predict the probability of thunderstorms and strong winds occurring at each grid point in the future at a preset time.

3. The method for classifying and forecasting thunderstorms and strong winds according to claim 1, characterized in that, Extracting local radar features includes, Extract the bow-shaped echo and the rear inflow gap, analyze the strong divergence characteristics in the low-level radial velocity field, and extract the signal of the sudden drop in echo top height in the vertical profile. Estimating the intensity of the cold pool includes calculating the temperature difference of the cold pool based on the ground-observed temperature and the temperature field around the overall risk area; and estimating the horizontal expansion of the cold pool using the low-level wind field and temperature gradient output by the numerical model. The time-series abrupt change signals include the temporal variations of boundary layer isentropic potential temperature and humidity extracted from numerical model outputs and ground observation data; and the abrupt change points of the temporal variations are identified using an attention mechanism model.

4. The method for classifying and forecasting thunderstorms and strong winds according to claim 1, characterized in that, The process of preprocessing various types of data and performing temporal and spatial registration yields a unified grid multi-source dataset, including: Bilinear interpolation is used to interpolate both the preprocessed radar data and the numerical model output data to a unified grid; time linear interpolation is performed on the numerical model data to match the time step of the radar data resolution; and a moving average filter is applied to smooth all time series data.

5. The method for classifying and forecasting thunderstorms and strong winds according to claim 1, characterized in that, The classification model further includes a fully connected layer. The classification process involves inputting the local radar features, cold pool intensity, and time-series abrupt change signals into the fully connected layer, passing two layers of ReLU activation and Dropout to output a high-order feature vector, and then sequentially performing predictions through a single decision tree and probability integration and analogy judgment of multiple trees to output the prediction result. The Monte Carlo method is used to estimate the 95% confidence interval to quantify the uncertainty of the prediction result.

6. A thunderstorm and gale classification forecasting system, characterized in that, include: The multi-source dataset acquisition module is used to acquire radar observation data, numerical model output and ground observation data, preprocess the various types of data and perform time and space registration to obtain a unified grid multi-source dataset. The overall risk area identification module is used to extract radar morphology features from the radar observation data through a convolutional neural network. The radar morphology features are combined with environmental factors based on numerical model output to identify the overall risk area of ​​thunderstorms and strong winds. The feature extraction module is used to extract local radar features, estimate cold pool intensity, and time-series abrupt change signals within the overall risk area. The downburst risk identification module is used to input the local radar features, cold pool intensity and time series abrupt change signal into the classification model, and output the downburst risk discrimination result. The risk discrimination result includes two categories: ordinary thunderstorm wind and downburst. The output form is probability and uncertainty interval.

7. The method for classifying and forecasting thunderstorms and strong winds according to claim 6, characterized in that, The classification model further includes a fully connected layer. The classification process involves inputting the local radar features, cold pool intensity, and time-series abrupt change signals into the fully connected layer, passing two layers of ReLU activation and Dropout to output a high-order feature vector, and then sequentially performing predictions through a single decision tree and probability integration and analogy judgment of multiple trees to output the prediction result. The Monte Carlo method is used to estimate the 95% confidence interval to quantify the uncertainty of the prediction result.

8. An electronic device, the electronic device comprising: At least one processor; and a memory communicatively connected to the at least one processor; wherein, The memory stores instructions that can be executed by the at least one processor to enable the at least one processor to perform any one of the methods of claims 1-6.

9. A computer-readable medium having stored thereon computer program instructions that can be executed by a processor to implement the method of any one of claims 1-6.

10. A computer program product comprising a computer program that, when executed by a processor, implements the method of any one of claims 1-6.

Images