A three-dimensional variational assimilation method of inverting vertical velocity from reflectivity factor

By inverting the reflectivity factor of the Doppler weather radar and combining it with the model background field, and using the adiabatic Richardson equation to assimilate the vertical velocity, the problem of inconsistency between dynamics and microphysics in existing technologies is solved, thus improving the accuracy and stability of convective-scale forecasts.

CN116559978BActive Publication Date: 2026-06-26LANZHOU UNIV

Patent Information

Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
LANZHOU UNIV
Filing Date
2023-05-15
Publication Date
2026-06-26

AI Technical Summary

Technical Problem

Existing technologies fail to effectively incorporate vertical velocity information when assimilating the reflectivity factor of Doppler weather radar, resulting in a lack of dynamic-microphysical consistency in the initial field of the model, which affects the effectiveness of convective-scale forecasts.

Method used

By combining Doppler weather radar reflectivity factor observations with the model background field, the proxy vertical velocity is inverted, and the adiabatic Richardson equation is used as the physical transformation observation operator for three-dimensional variational assimilation. The mass conservation, adiabatic and hydrostatic equilibrium relationships are integrated to optimize the initial field of the dynamic field.

Benefits of technology

It improves the accuracy and stability of convective-scale model forecasts, reduces the dynamic spin-up time of model integration, and improves the forecast performance of convective systems.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN116559978B_ABST
    Figure CN116559978B_ABST
Patent Text Reader

Abstract

The present application relates to a kind of three-dimensional variation assimilation method of reflectivity factor inversion vertical velocity, the method includes the following steps: (1) by Doppler weather radar reflectivity factor observation combines mode background field vertical velocity statistical characteristics inversion proxy vertical velocity observation w ; (2) using adiabatic richardson equation as physical conversion observation operator assimilation proxy vertical velocity observation, which is composed of cost function and its gradient value;(3) through three-dimensional variation minimization process constantly repeat step (2), obtain the control variable that makes gradient obtain minimum cv , then according to analysis increment, by the horizontal, vertical and physical conversion of control variable, obtain the analysis increment of this assimilation time, and further obtain the final analysis field.The present application realizes the three-dimensional variation assimilation of vertical velocity, and increases the dynamic information of mode initial field by assimilating the vertical velocity proxy observation of reflectivity factor inversion, so as to improve the convection scale prediction effect of mode.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the field of weather forecasting, and more particularly to a three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion. Background Technology

[0002] In the field of mesoscale numerical weather prediction, higher model resolution (cloud resolution scale) and an initial field containing sufficient small and medium-scale information are important ways to improve forecast performance. As an important source of small and medium-scale information in the initial field of a model, the Doppler weather radar can reduce the model spin-up time to a certain extent through its assimilation, playing an important and irreplaceable role in weather forecasting, especially short-term and nowcasting (Sun Juanzhen, Chen Mingxuan, Fan Shuiyong (2016), Radar Data Assimilation Methods: Retrospect and Prospect, Progress in Meteorological Science and Technology, (3):11).

[0003] Doppler weather radar can acquire radial velocity and reflectivity factor observations. The reflectivity factor contains information on condensate and water vapor. In the three-dimensional variational assimilation system, it can improve model variables such as condensate and cloud water vapor in the initial field through direct assimilation (Xiao, Q., and J. Sun (2007), Multiple-radar data assimilation and short-range quantitative precipitation forecasting of a squall line observed during IHOP_2002. Mon. Wea. Rev., 135, 3381–3404) or indirect assimilation (Wang, H., J. Sun, S. Fan, and X.-Y. Huang (2013), Indirect assimilation of radar reflectivity with WRF 3D-Var and its impact on prediction of four summertime convective events. J. Appl. Meteorol. Climatol., 52, 889–902). It is widely used in operational numerical weather prediction systems to initialize high-resolution convective forecasts. However, if only water condensate and water vapor are adjusted when assimilating weather radar data without adjusting the dynamic field, the initial field of the model will lack microphysical and dynamic consistency among different state variables. This will cause the rainwater formed in the initial field to fall quickly in the early stage of model integration, resulting in no rainwater to fall in subsequent models. Since radial wind contains information about the atmospheric dynamic field and is a wind component in a certain direction, it can be combined with reflectivity factor assimilation to generate a more coordinated dynamic and microphysical analysis field, thereby obtaining better forecast results (Xiao, Q., and J. Sun (2007), Multiple-radar data assimilation and short-range quantitative precipitation forecasting of asquall line observed during IHOP_2002. Mon. Wea. Rev., 135, 3381–3404).Simultaneously, a more harmonious initial field for microphysical and thermodynamic fields can be obtained by assimilating the latent heat obtained from the reflectivity factor (Liu, Y., T.T. Warner, J.F. B. Worers, et al. (2008), The operational mesogamma-scale analysis and forecast system of the USArmy Test and evaluation command. Part I: overview of the modeling system, the forecast products, and how the products are used. J Appl Meteor Climatology, 47: 1077–1092). Furthermore, strong vertical velocity is one of the main characteristics of mesoscale weather systems and is an important factor in the development and maintenance of these systems. In convective-scale numerical forecasting, after assimilating the condensate and water vapor information retrieved from the Doppler weather radar, vertical motion is crucial for the maintenance of clouds formed at the initial moment and for forecasted precipitation, reducing the dynamic spin-up time of simulated convection. For severe weather forecasts caused by deep convection, especially for severe weather forecasts of deep convection with weak exotropic forcing, the spin-up time of the model's initial field is even more important. Therefore, assimilating the vertical velocity in the initial field may help initialize the dynamic field of the mode.

[0004] Since the reflectivity factor is an indicator of the strength of convective weather system development, it contains not only microphysical processes (hydrocondensate information) but also dynamic processes within the cloud (vertical velocity information within the cloud). Based on the method for constructing vertical velocity within the cloud (Yuter, SE, and RAHouze (1995), Three-dimensional kinematic and microphysical evolution of Florida cumulonimbus. Part II: Frequency distributions of vertical velocity, reflectivity, and differential reflectivity. Monthly Weather Review, 123, 1921–1940), vertical velocity observations within the cloud can be further constructed using radar reflectivity factors (Liu Hongya et al. (2010), Three-dimensional variational assimilation of radar data for a case study of heavy rain, Acta Meteorologica Sinica, 68(6):779–789), and then more dynamic field information can be provided to the initial field by assimilating the vertical velocity. However, the vertical velocity inversion formula used by Liu Hongya et al. (2010) assumes that the maximum vertical velocity in the convective region occurs at the same height, i.e., 6 km. However, in reality, the altitude at which the maximum vertical velocity occurs is always different at different stages of the convective system's life cycle (Schumacher, C., S.N. Stevenson, and C.R. Williams (2015), Vertical motions of the tropical convective cloud spectrum over Darwin, Australia. Quarterly Journal of the Royal Meteorological Society, 141, 2277–2288). Furthermore, previous sensitivity experiments on the altitude of the maximum vertical velocity have shown that the impact of vertical velocity assimilation on convective-scale precipitation forecasts is significantly affected by the altitude of the maximum vertical velocity. Considering that the purpose of three-dimensional variational assimilation is to add an increment that satisfies physical constraints to the background field, assimilating observations with statistical characteristics similar to the background field is more likely to produce an equilibrium analytical field. Therefore, determining the parameters in the vertical velocity inversion formula based on the statistical characteristics of the background field is reasonable.

[0005] Three-dimensional variational assimilation methods are widely used in operational weather forecasting due to their low computational cost and ease of adding nonlinear physical constraints. Based on the principles of three-dimensional variation, one approach to directly assimilate vertical velocity is to treat it as a control variable characterizing atmospheric state. This method only requires prior calculation of the background error covariance including the vertical velocity control variable; however, this assimilation method lacks physical constraints on the vertical velocity and other model variables at the analysis time, easily generating high-frequency noise and leading to model instability. Another approach is to use physically constrained observation operators to establish a physical transformation between model variables and vertical velocity observations (e.g., the mass continuity equation and Richardson's equation). Among these, using Richardson's equation, which integrates the mass conservation equation, adiabatic thermodynamic equations, and hydrostatic equilibrium relationships, as the observation operator to assimilate vertical velocity yields a more consistent initial field for the three-dimensional dynamic-thermal field. Summary of the Invention

[0006] The technical problem to be solved by the present invention is to provide a three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor to improve the forecasting effect of convection scale models.

[0007] To address the aforementioned problems, the present invention provides a three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion, comprising the following steps:

[0008] (1) The proxy vertical velocity observation w was retrieved by combining the Doppler weather radar reflectivity factor observation with the statistical characteristics of the vertical velocity in the background field of the model.

[0009] (2) The adiabatic Richardson equation is used as the physical transformation observation operator to assimilate the surrogate vertical velocity observation. This assimilation surrogate vertical velocity observation is composed of the cost function J and its gradient. Value composition;

[0010] (3) By repeatedly performing step (2) through the three-dimensional variational minimization process, the control variable cv that minimizes the gradient is obtained. Then, based on the analysis increment δx, the analysis increment δx at this assimilation moment is obtained through the horizontal, vertical and physical transformations of the control variable, and finally the analysis field is obtained.

[0011] In step (1), the inverted proxy vertical velocity observation w (unit: m / s) is calculated using the following formula:

[0012]

[0013] In the formula: Z is the reflectivity factor, in dBZ; e≈2.71828; λ determines the altitude range of the main distribution of the vertical velocity inversion value, in km; h is the altitude at which the reflectivity factor is observed, in km; h0 is the altitude at which the maximum vertical velocity is located, in km;

[0014] The method for determining h0 is to statistically analyze the average vertical velocity profile of the background field mode grid point combination with a reflectivity factor greater than 35 dBZ at a certain analysis time, and use the height of the maximum vertical velocity of the profile as the parameter h0 for inverting the vertical velocity at that analysis time.

[0015] The method for determining λ is to find the heights at which the average vertical velocity first reaches zero, both upwards and downwards from h0, and then divide the difference between the two heights by 20 to obtain the parameter λ.

[0016] The expression for the cost function J in step (2) is:

[0017]

[0018] In the formula: cv is the control variable; R is the observation error covariance; d is the observation increment; x b For background field; y o For vertical velocity observation; H(x) b The observation operator () is applied to the background field to calculate the vertical velocity of the background field at the observation location; H is the observation operator, which is the algorithm for converting model variables into observations; H p For physical transformation operators; H s This is the interpolation operator; T represents the transpose of the matrix; –1 represents the inverse of the matrix;

[0019] The HUcv term refers to the term obtained by applying the tangent linear operator H of the observation operator to the analytical increment δx; the analytical increment... x represents the analysis field; U h U v U p These represent horizontal transformation, vertical transformation, and physical transformation, respectively; the relationship between U and the background error covariance B is: B = UU T .

[0020] The gradient of the cost function in step (2) The expression is:

[0021]

[0022] Where: H T U is the adjoint of the tangent linear operator H; T is the transpose of matrix U; –1 represents the inverse of the matrix.

[0023] Compared with the prior art, the present invention has the following advantages:

[0024] 1. This invention establishes a proxy vertical velocity inversion (obtained by reflectivity factor inversion) based on the statistical characteristics of the model background field and its observation operator, realizing three-dimensional variational assimilation of vertical velocity, which can effectively increase the dynamic information of the model's initial field, thereby improving the model's convective-scale forecasting effect.

[0025] 2. The vertical velocity observation operator used in this invention is the Richardson equation, which integrates the mass conservation equation, the adiabatic thermodynamic equation, and the static equilibrium relationship, and has the characteristics of high accuracy and good stability.

[0026] 3. The assimilation method of this invention adopts a three-dimensional variational method, which has high applicability in operational weather forecasting. Attached Figure Description

[0027] The specific embodiments of the present invention will be described in further detail below with reference to the accompanying drawings.

[0028] Figure 1 (a) is a distribution map of combined reflectivity factors greater than 35 dBZ at a certain analysis time of the present invention (15 UTC on July 4, 2020); (b) is a distribution map of the maximum vertical velocity (unit m / s) obtained by reflectivity factor inversion; (c) is a schematic diagram of the distribution (scatter points) of vertical velocity inverted from fixed parameters h0 and λ at different heights and the background field average vertical velocity profile; (d) is a schematic diagram of the distribution (scatter points) of vertical velocity inverted from parameters h0 and λ obtained from background field statistical characteristics at different heights and the background field average vertical velocity profile.

[0029] Figure 2 (a) The horizontal wind analysis increment (m / s; vector) at layer 14 (~850 haPa) of the model at a certain analysis time (July 4, 2020, 15 UTC) after assimilation of vertical velocity according to this invention, and the background field 850 hPa rainwater mixing ratio (g / kg; contour lines: ~0.04-0.4 g / kg) at the same time; (b) and (c) are the horizontal wind (m / s; vector) and divergence (×10) at layers 14 and 24 (~500 hPa) of the model at the same time, respectively. -4 / s; contour lines: dashed lines are negative, solid lines are positive) Schematic diagram of the analysis increment.

[0030] Figure 3 This is a schematic diagram comparing the combined reflectance factor 6-hour forecast field of the Assimilation Vertical Velocity Experiment (DA) of this invention (16-21 UTC, July 4, 2020) with the combined reflectance factor observation (OBS) and the combined reflectance factor 6-hour forecast field of the Unassimilated Experiment (CTL). Detailed Implementation

[0031] A three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion includes the following steps:

[0032] (1) The vertical velocity observation w was obtained by inverting the statistical characteristics of the vertical velocity in the background field of the Doppler weather radar by combining the reflectivity factor observation with the model background field.

[0033] The inverted proxy vertical velocity observation w (unit: m / s) is calculated using the following formula:

[0034]

[0035] In the formula: Z is the reflectivity factor, in dBZ; e≈2.71828; λ determines the altitude range of the main distribution of the vertical velocity inversion value, in km; h is the altitude at which the reflectivity factor is observed, in km; h0 is the altitude at which the maximum vertical velocity is located, in km.

[0036] Liu Hongya et al. (Liu Hongya et al. (2010), Case Study of Heavy Rainfall Using Three-Dimensional Variational Assimilation Radar Data, Acta Meteorologica Sinica, 68(6):779–789) used fixed h0 (6km) and λ (0.4). However, for convection at different development stages, the height at which the maximum vertical velocity is located is different (Schumacher, C., S.N. Stevenson, and C.R. Williams (2015), Vertical motions of the tropical convective cloud spectrum over Darwin, Australia. Quarterly Journal of the Royal Meteorological Society, 141, 2277–2288), which leads to differences in λ. Considering that the essence of three-dimensional variational assimilation is to add an analytical increment based on the difference between observation and background to the background field, it is easier to generate an equilibrium analytical field if the observation and background fields have similar statistical characteristics. Therefore, this invention first combines the statistical information of the vertical velocity of the background field at the analysis time to determine the parameters h0 and λ.

[0037] The method for determining h0 is to statistically analyze the average vertical velocity profile of the background field mode grid combination with a reflectivity factor greater than 35 dBZ at a certain analysis time, and use the height of the maximum vertical velocity of the profile as the parameter h0 for inverting the vertical velocity at that analysis time.

[0038] The method for determining λ is to find the heights at which the average vertical velocity first reaches zero, both upwards and downwards from h0, and then divide the difference between the two heights by 20 to obtain the parameter λ.

[0039] (2) This invention uses the adiabatic Richardson equation as the physical transformation observation operator to assimilate the vertical velocity. The variational method for solving the optimal analysis field is to find the minimum value of the cost function (Lorenc, AC (1986), Analysis methods for numerical weather prediction. Quarterly Journal of the Royal Meteorological Society, 112, 1177–1194).

[0040] Specifically, the adiabatic Richardson equation is used as the physical transformation observation operator to assimilate surrogate vertical velocity observations. These surrogate vertical velocity observations are determined by the cost function J and its gradient. Value composition.

[0041] a. The expression for the cost function J is:

[0042]

[0043] In the formula: cv is the control variable; R is the observation error covariance; d is the observation increment; x b For background field; y o For vertical velocity observation; H(x) b ) is the observation operator applied to the background field to calculate the vertical velocity of the background field at the observation position; T represents the transpose of the matrix; –1 represents the inverse of the matrix.

[0044] H is the observation operator, an algorithm that converts model variables into observations. Since model variables reside at model grid points, the observation operator H consists of two parts:

[0045] H = H p H s

[0046] Where: H p For the physical transformation operator (the physical equation for calculating observations from model variables); H s This is an interpolation operator (interpolation from model space to observation point).

[0047] For vertical velocity observation, H in this invention p The Richardson equation (Richardson, LF, 1922: Weather Prediction by Numerical Process. Cambridge University Press, 236pp) combines the mass conservation equation, the adiabatic thermodynamic equation, and the static equilibrium equation.

[0048]

[0049] In the formula: γ is the specific heat of dry air at constant pressure, c p (1004J / (kg K)) and specific heat at constant volume c v (717J / (kg K)) ratio (γ=c) p / c v =1.4); p is air pressure, in hPa; w is vertical velocity, in m / s; z is altitude, in m; V h Horizontal wind speed (including u and v components, unit: m / s); The divergence of the horizontal wind; The pressure gradient is denoted by γ; g is the acceleration due to gravity (9.80665 m / s²). 2 ); ρ is the air density, in kg / m³ 3 .

[0050] The HUcv term refers to the term obtained by applying the tangent linear operator H of the observation operator to the analytical increment δx; the analytical increment... x represents the analysis field; U h U v U p These represent horizontal transformation, vertical transformation, and physical transformation, respectively; the relationship between U and the background error covariance B is: B = UU T .

[0051] H is the tangent linear operator of the observation operator, and the physical transformation operator H for the vertical velocity observation term. p Its tangent linear operator H p The expression is as follows:

[0052]

[0053] Where: p, ρ, w, V h It can be written as In form, A represents the basic state, and A′ represents the increment. The units of the basic state and the increment are consistent with the units of the variables.

[0054] b. The gradient of the cost function is used to determine the descent direction in the minimization process. (Cost function gradient) The expression is:

[0055]

[0056] In the formula: H T U is the adjoint of the tangent linear operator H; T Let U be the transpose of matrix U.

[0057] (3) By repeatedly performing step (2) through the three-dimensional variational minimization process, the control variable cv that minimizes the gradient is obtained. Then, based on the analysis increment δx: δx = U p Uv U h cv, by controlling the horizontal, vertical, and physical transformations of the variables, yields the analytical increment δx at this assimilation moment, and then, based on x = x b +δx yields the final analysis field.

[0058] This analysis field assimilates the vertical velocity and contains more dynamic information, making it suitable as the initial field for subsequent numerical prediction integration.

[0059] Example 1

[0060] When using a three-dimensional variational assimilation system to assimilate vertical velocities to simulate a heavy precipitation event in a certain area, the following steps are performed for calculation. The numerical model used in this example is the Mesoscale Weather Numerical Forecast System of China Meteorological Administration (CMA-MESO).

[0061] Step a: Calculate the agent's vertical velocity in the CMA-MESO coordinate system. Observation: First, obtain the vertical velocity w in the z-coordinate system, and then obtain the vertical velocity in the CMA-MESO vertical coordinate system through coordinate transformation. The specific steps are as follows:

[0062] Step a1: Calculate the vertical velocity w in the z-coordinate system. First, interpolate the radar base data in the polar coordinate system to the z-coordinate system to obtain the three-dimensional (longitude, latitude, altitude) radar reflectivity factor observation in the z-coordinate system.

[0063] Secondly, using formula (1) (Liu Hongya et al. (2010), a case study of heavy rain using three-dimensional variational assimilation radar data, Acta Meteorologica Sinica, 68(6):779–789), the three-dimensional vertical velocity w within the cloud was calculated from the observation of the three-dimensional radar reflectivity factor:

[0064]

[0065] Vertical velocity is calculated only for reflectivity factors with an intensity greater than 35 dBZ (convective cloud region). In the formula: Z is the reflectivity factor, in dBZ; e≈2.71828; h is the altitude at which the reflectivity factor was observed, in km. The parameter h0 (altitude of maximum vertical velocity, in km) and the coefficient λ are determined by combining the statistical characteristics of the vertical velocity in the model background field. The steps are as follows:

[0066] (1) For a certain analysis time, the average vertical velocity profile of the background field mode grid combination with a reflectivity factor greater than 35dBZ is statistically analyzed, and the height of the maximum vertical velocity of the profile is used as the parameter h0 for inverting the vertical velocity at that analysis time.

[0067] (2) Find the heights at which the vertical velocity first reaches zero, both upward and downward, using h0 as the reference. Divide the difference between the two heights (in km) by 20 to obtain the parameter λ.

[0068] like Figure 1 As shown, (a) is the distribution map of the combined reflectivity factor greater than 35 dBZ at 15 UTC on July 4, 2020; (b) is the distribution map of the maximum vertical velocity calculated according to formula (1) (h0 and γ combined with background field statistical characteristics). It can be seen that the maximum vertical velocity has good consistency with the region with a reflectivity factor greater than 35 dBZ. (c) is the vertical velocity (scattered points) at different heights and the background field average vertical velocity profile (black solid line) obtained according to formula (1) (parameter h0 is set to 6 km, λ is set to a fixed value of 0.4). (d) is the vertical velocity (scattered points) at different heights and the background field average vertical velocity profile (black solid line) obtained according to formula (1) (h0 and λ combined with background field statistical characteristics). It can be seen that the vertical velocity distribution obtained by combining background field vertical velocity statistical information has better consistency with the background field average vertical velocity profile (the height of the maximum vertical velocity is the same as the height range of the main vertical velocity distribution). The maximum vertical velocity obtained by inversion can reach ~2.5 m / s, and is around 6 km. Most vertical velocity inversion values ​​are below 5 km, and most values ​​are between 0 and 0.7 m / s.

[0069] Step a2: Calculate the vertical velocity in the CMA-MESO height terrain following coordinate system.

[0070] Because CMA-MESO uses height-terrain-following coordinates in the vertical direction:

[0071]

[0072] Among them, z0 and z t The vertical velocity in this coordinate system represents the terrain height and the model layer top height (in meters). (Unit: m / s) is:

[0073]

[0074] in:

[0075] Δz0=z t -z0(x, y)

[0076] Δz z =zt -z(x, y)

[0077] w0 is the vertical velocity at the Earth's surface.

[0078] Step b uses the proxy vertical velocity obtained from the inversion in step a as the observation y. o Substitute into the three-dimensional variational assimilation system to calculate the cost function J and its gradient. The goal is to find the minimum value of the cost function.

[0079] Step b1: Calculate the cost function.

[0080]

[0081] Where cv is the control variable, R is the observation error covariance, and d is the observation increment: d = H(x b )-y o , where x b For background field; y o For observation purposes, in this invention, y o Represents vertical velocity observation. H(x) b H(x) is obtained by applying the observation operator to the background field and calculating the vertical velocity of the background field at the observation position; then, the difference between this vertical velocity and the observed vertical velocity is taken to obtain H(x). b )-y o The observation operator H is an algorithm that converts model variables into observations. Since model variables reside at model grid points, the observation operator H consists of two parts:

[0082] H = H p H s (5)

[0083] Among them, H s For the interpolation operator (interpolation from model space to observation point), H p For physical transformation operators (physical equations for calculating observations from model variables), H in this invention is used for vertical velocity observations. p The Richardson equation (Equation (6)) is adopted, which combines the mass conservation equation, the adiabatic thermodynamic equation and the static equilibrium equation.

[0084]

[0085] In the above formula, γ is the specific heat of dry air at constant pressure, c. p (1004J / (kg K)) and specific heat at constant volume c v (717J / (kg K)) ratio (γ=c) p / c v =1.4); κ=c p / R≈3.5, (R is the dry air gas constant: 287.05 J / (kg K)). Π is the dimensionless air pressure. p is the air pressure (in hPa, P0 = 1000 hPa); u and v are the horizontal wind components (in m / s).

[0086] (4) The HUcv term in the formula:

[0087] Analyze incremental U h U v U p These represent horizontal transformation, vertical transformation, and physical transformation, respectively; the relationship between U and the background error covariance B is: B = UU T ), where x is the analysis field.

[0088] H is the tangent linear operator of the observation operator, and the physical transformation operator H for the vertical velocity observation term. p Its tangent linear operator H p The expression is given by equation (7), and the tangent linear operator is applied to the analysis increment to obtain the HUcv term.

[0089]

[0090] In the above formula, Π, u and v can be written as In form, A represents the basic state, and A′ represents the increment.

[0091] Substituting the above calculation results into equation (4), we can calculate J.

[0092] Step b2: Calculate the gradient of the cost function.

[0093]

[0094] Among them, H T This is the adjoint of the tangent linear operator H. For the physical transformation of the vertical velocity observation term, the tangent linear operator H... p Its adjoint operator is the adjoint of equation (7), which has been calculated in step b1 of the HUcv term.

[0095] Step c involves repeatedly performing step b through the CMA-MESO three-dimensional variational minimization process to obtain the control variable cv that minimizes the gradient. Then, based on δx = U... p U v U h cv, by controlling the horizontal, vertical, and physical transformations of the variables, yields the analytical increment δx at this assimilation moment. Therefore, based on x = x... b +δx yields the final analysis field.

[0096] Figure 2 (b) and (c) represent the analysis increments of horizontal wind (vector) and divergence (fill plot) at the 14th (~850hPa) and 24th (~500hPa) layers of the model, respectively, at the analysis time of 15 UTC on July 4, 2020. Figure 2 (a) is a distribution map of the horizontal wind increment at layer 14 of the model and the 850 hPa rainwater mixing ratio (contours) of the model background field at that moment. It can be seen that in the middle layer of the model ( Figure 2 (b) There is a certain convergence increase in Cangzhou, Hebei Province and northern Henan Province, while the model background field in this region is associated with a weak rainwater mixing ratio (~0.04-0.4 g / kg), and at the same time in the middle layer of the model ( Figure 2 (c) In the lower-level convergence area of ​​Cangzhou, there is obvious horizontal wind divergence or weak horizontal wind convergence. This configuration of lower-level convergence and middle-level divergence is conducive to the vertical movement of the atmosphere, which leads to the saturation of water vapor as it rises and condenses, increasing the probability of precipitation.

[0097] Figure 3 This is a comparison map of the 6-hour (16 UTC-21 UTC, July 4, 2020) forecast fields for combined reflectance factors. OBS represents the observed combined reflectance factors; both CTL and DA are forecast experiments using CMA-MESO as the forecast model. The simulation area includes Hebei Province, with a 3 km grid spacing. CTL is the control experiment: using NCEP (National Center for Environmental Prediction) GFS (Global Forecast System) data as initial and boundary conditions, forecasting began at 06 UTC on July 4, 2020, without assimilation, directly forecasting to 00 UTC on July 5, 2020. DA is the assimilation experiment: forecasting began on July 4, 2020, forecasting for 6 hours, then assimilating proxy vertical velocity observations at 12 UTC, then forecasting for 1 hour again and assimilating proxy vertical velocity observations, repeating this until 15 UTC, assimilating proxy vertical velocity observations and then forecasting for 6 hours to 00 UTC on July 5. Figure 3As can be seen, the observed strong echo areas are divided into the Cangzhou area of ​​Hebei (dashed box) and the northern Henan strong echo area. For the northern Henan strong echo area, it moved slowly eastward over 6 hours, with its intensity weakening after 20 UTC. Both the CTL and DA experiments showed poor predictability for this strong echo area. Compared to observations, the CTL experiment's predicted intensity was generally weaker and had a shorter lifespan. The DA experiment's predicted intensity was weaker and its location was further north. For the Cangzhou, Hebei strong echo area, it moved slowly eastward over 6 hours, with its intensity remaining relatively stable between 16-19 UTC, gradually weakening after 19 UTC. The CTL experiment failed to predict the observed strong echo area and exhibited false echo predictions to the southwest of the observed strong echo area. Furthermore, in the subsequent 3-hour forecast field, false echoes also developed to the east of the observed strong echo area. Compared to the CTL experiment, the DA experiment predicted a stronger strong echo area further north and closer in location to the observed area. A certain reflectivity factor was predicted within the observed strong echo region, and the false echo predictions present in the CTL experiment were somewhat mitigated. This demonstrates that assimilating the vertical velocity derived from the reflectivity factor can effectively improve convective-scale forecasts.

[0098] Of course, this method of assimilating vertical velocities derived from the reflectivity factor of Doppler weather radar is also applicable to assimilating vertical velocities obtained from other types of observations.

Claims

1. A three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion, comprising the following steps: (1) The proxy vertical velocity observation w is retrieved by combining Doppler weather radar reflectivity factor observations with the statistical characteristics of vertical velocity in the model background field; the retrieved proxy vertical velocity observation w is calculated using the following formula: In the formula: Z is the reflectivity factor, in dBZ; e≈2.71828; λ determines the altitude range of the main distribution of the vertical velocity inversion value, in km; h is the altitude at which the reflectivity factor is observed, in km; h0 is the altitude at which the maximum vertical velocity is located, in km; The method for determining h0 is to statistically analyze the average vertical velocity profile of the background field mode grid point combination with a reflectivity factor greater than 35 dBZ at a certain analysis time, and use the height of the maximum vertical velocity of the profile as the parameter h0 for inverting the vertical velocity at that analysis time. The method for determining λ is to find the heights at which the average vertical velocity first reaches zero, both upwards and downwards from h0, and then divide the difference between the two heights by 20 to obtain the parameter λ. (2) The adiabatic Richardson equation is used as the physical transformation observation operator to assimilate the surrogate vertical velocity observation. This assimilation surrogate vertical velocity observation is composed of the cost function J and its gradient. Value composition; (3) By repeatedly performing step (2) through the three-dimensional variational minimization process, the control variable cv that minimizes the gradient is obtained. Then, based on the analysis increment δx, the analysis increment δx at this assimilation moment is obtained through the horizontal, vertical and physical transformations of the control variable, and finally the analysis field is obtained.

2. The three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion as described in claim 1, characterized in that: The expression for the cost function J in step (2) is: d=H(x b )-y o ;H=H p H s ; In the formula: cv is the control variable; R is the observation error covariance; d is the observation increment; x b For background field; y o For vertical velocity observation; H(x) b The observation operator () is applied to the background field to calculate the vertical velocity of the background field at the observation location; H is the observation operator, which is the algorithm for converting model variables into observations; H p For physical transformation operators; H s This is the interpolation operator; T represents the transpose of the matrix; –1 represents the inverse of the matrix; The HUcv term refers to the term obtained by applying the tangent linear operator H of the observation operator to the analytical increment δx; the analytical increment δx = xx b =U p U v U h cv = Ucv; x is the analysis field; U h U v U p These represent horizontal transformation, vertical transformation, and physical transformation, respectively; the relationship between U and the background error covariance B is: B = UU T .

3. The three-dimensional variational assimilation method for vertical velocity derived from reflectivity factor inversion as described in claim 1, characterized in that: The gradient of the cost function in step (2) The expression is: In the formula: H T U is the adjoint of the tangent linear operator H; T Let U be the transpose of matrix U.