A GNSS-MET-based water vapor tomography modeling method and device

Abstract

The application discloses a kind of water vapor tomography modeling method and device based on GNSS-MET, belong to atmospheric water vapor information inversion technical field.The method obtains satellite signal and meteorological parameter by GNSS-MET receiver, and solves out slant path water vapor content (SWV) as tomography input;Select dense and evenly distributed station to construct tomography area, and divide voxel grid;Utilize Gaussian distance weighting to construct horizontal constraint, combined with water vapor vertical index distribution characteristics to construct vertical constraint, and introduce ERA5 reanalysis data as priori constraint, solve the limitation of traditional sounding data constraint;Tomography equation is solved by singular value decomposition (SVD) method, and three-dimensional water vapor density distribution is obtained.The present application utilizes the high space-time resolution characteristics of ERA5 data, significantly improves the accuracy and reliability of water vapor tomography, and is suitable for extreme weather prediction, meteorological monitoring and other fields.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of atmospheric water vapor information inversion, and particularly relates to a water vapor tomography modeling method and device based on GNSS-MET. BACKGROUND

[0002] More than 90% of atmospheric water vapor is collected in the troposphere, and its changes are closely related to extreme weather events and meteorological disasters. High temporal and spatial resolution water vapor monitoring is of great significance for weather forecasting, climate research and disaster warning. Traditional water vapor detection methods mainly include radio sonde, microwave radiometer, weather radar and optical radiometer, etc. However, these methods have obvious limitations: radio sonde observation is sparse and costly, microwave radiometer is easily disturbed by clouds and rain, weather radar is difficult to detect non-precipitation cloud data, and optical radiometer cannot work all-weather.

[0003] In recent years, GNSS-based water vapor inversion technology has been widely used due to its all-weather, low-cost, high-precision and other advantages. GNSS signals will be delayed when passing through the atmosphere, and by solving this delay, the atmospheric water vapor content can be inverted. Combined with the GNSS-MET receiver of the meteorological observation system (MET), satellite signals can be obtained, and meteorological parameters such as temperature, humidity and pressure can be measured synchronously, providing more comprehensive data support for water vapor inversion. However, traditional GNSS water vapor inversion mainly obtains precipitable water vapor (PWV) in the zenith direction, which can only reflect the vertical integral total amount of water vapor above the station, and cannot reveal the three-dimensional spatial distribution characteristics of water vapor.

[0004] To solve this problem, water vapor tomography technology has emerged. This technology establishes a mathematical relationship between the GNSS signal slant path water vapor content (SWV) and the three-dimensional voxel grid, constructs a tomography equation to invert the three-dimensional water vapor density distribution in the region. However, due to the uneven distribution of GNSS signal rays, the tomography equation often has a rank deficiency problem, leading to unstable solution. Existing methods usually introduce horizontal constraints, vertical constraints and prior constraints to improve the solution effect. Among them, the prior constraint is crucial to the accuracy of the tomography result. Traditional methods often use radio sonde data as prior constraints, but sonde data is only valid at specific sites and cannot reflect the horizontal variation of water vapor in the tomography region. Moreover, the temporal and spatial resolution is low, which makes it difficult to meet the high-precision tomography requirements.

[0005] The ERA5 reanalysis data provided by the European Centre for Medium-Range Weather Forecasts (ECMWF) provides a new way to solve the above problems. ERA5 has high temporal and spatial resolution (hourly temporal resolution, 37 pressure layers in the vertical direction), and can provide gridded data of temperature, humidity, pressure and other meteorological parameters on a global scale, which can more accurately describe the spatial distribution characteristics of water vapor. Compared with sonde data, ERA5 is more suitable as a prior constraint background field for water vapor tomography. SUMMARY

[0006] To solve the above technical problems, the application provides a GNSS-MET-based water vapor tomography modeling method and device, which uses ERA5 data to construct prior constraints, combines horizontal and vertical constraints to optimize the tomography equation, and solves the three-dimensional water vapor density distribution through singular value decomposition (SVD).

[0007] To achieve the above purpose, the technical scheme adopted by the application is as follows:

[0008] A GNSS-MET-based water vapor tomography modeling method, the method comprising:

[0009] Step S1, acquiring satellite signals and meteorological parameters based on GNSS-MET receivers respectively, and solving slant path water vapor content SWV;

[0010] Step S2, obtaining historical observation data and performing three-dimensional voxel grid division, converting the slant path water vapor content SWV into an integral form corresponding to the voxel grid, and constructing an observation equation of the tomography model;

[0011] Step S3, applying horizontal constraints, vertical constraints and prior constraints to the observation equation of the tomography model to obtain a joint observation equation; wherein the horizontal constraints are based on a Gaussian distance weighting method, the vertical constraints are based on the exponential distribution characteristics of water vapor, and the prior constraints are a water vapor density background field calculated by a reanalysis data set;

[0012] Step S4, solving the joint observation equation through singular value decomposition to obtain a three-dimensional water vapor density distribution.

[0013] On the other hand, the application provides a GNSS-MET-based water vapor tomography modeling device, comprising a solving module, a voxel division module, a constraint module and an output module, wherein,

[0014] The solving module is used for acquiring satellite signals and meteorological parameters based on GNSS-MET receivers respectively, and solving slant path water vapor content SWV;

[0015] The voxel division module is used for obtaining historical observation data and performing three-dimensional voxel grid division, converting the slant path water vapor content SWV into an integral form corresponding to the voxel grid, and constructing an observation equation of the tomography model;

[0016] The constraint module is used for imposing horizontal constraint, vertical constraint and prior constraint on the observation equation of the tomographic model to obtain a joint observation equation; wherein the horizontal constraint is based on a Gaussian distance weighting method, the vertical constraint is based on the exponential distribution characteristics of water vapor, and the prior constraint adopts a water vapor density background field calculated by a reanalysis data set;

[0017] The output module is used for solving the joint observation equation by a singular value decomposition method to obtain a three-dimensional water vapor density distribution.

[0018] In a third aspect, the present application provides an electronic device, comprising: one or more processors; a memory for storing one or more programs; wherein when the one or more programs are executed by the one or more processors, the one or more processors implement the aforementioned GNSS-MET-based water vapor tomographic modeling method.

[0019] In a fourth aspect, the present application provides a computer-readable storage medium having stored executable instructions, which, when executed by a processor, enable the processor to implement the aforementioned GNSS-MET-based water vapor tomographic modeling method.

[0020] The beneficial effects of the present application are:

[0021] The present application introduces high spatiotemporal resolution ERA5 reanalysis data as prior constraint, which significantly improves the accuracy and reliability of water vapor tomography compared with traditional sounding data constraint method, and can more accurately reflect the three-dimensional distribution characteristics of water vapor. Secondly, the combination of GNSS-MET receiver network observation and intelligent constraint system design effectively solves the rank deficiency problem of tomographic equation, making the inversion result more stable.

[0022] In terms of practicality, the present application realizes 1-hour time resolution water vapor monitoring, greatly improving the timeliness of observation and providing strong support for short-weather forecasting. At the same time, based on the widely laid GNSS-MET network and free ERA5 data, the observation cost is greatly reduced, and the convenience of data acquisition is improved.

[0023] In addition, the present application is not limited by weather conditions and can realize all-weather continuous monitoring. The obtained three-dimensional water vapor products can be directly applied to numerical weather prediction, disaster warning and other fields, and have important business application value. Through multi-source data cross-validation, the reliability of the inversion result is ensured, providing high-quality data support for meteorological research and business application. BRIEF DESCRIPTION OF DRAWINGS

[0024] Figure 1 The flowchart of the GNSS-MET-based water vapor tomographic modeling method of the present application is shown in the figure;

[0025] Figure 2 The three-dimensional water vapor density distribution diagram is shown in the figure.

[0026] Figure 3 The water vapor profile result of the method is compared with a water vapor profile of a sounding balloon. DETAILED DESCRIPTION

[0027] The application will be further described below in combination with the drawings and examples.

[0028] The application provides a water vapor tomography modeling processing method based on GNSS-MET, extracts high-precision water vapor parameter information, obtains a basic observation equation through network observation, constructs a complete tomography equation based on sounding and reanalysis data and a weight matrix, and obtains water vapor density to generate a water vapor profile, and the like.

[0029] As shown in the figure, the water vapor tomography modeling method based on GNSS-MET specifically comprises the following steps: Figure 1

[0030] In step S1, satellite signals and meteorological parameters are obtained based on GNSS-MET receivers, and slant path water vapor content SWV is solved.

[0031] Water vapor parameter estimation results in the current and a period of time in the past are obtained, specifically, navigation satellite system signals are received and processed, information received by a station is saved as an observation file and a meteorological file, temperature and humidity pressure data of the station are obtained and zenith troposphere wet delay gradients in the east-west and south-north directions are inverted, and zenith wet delay ZWD obtained is converted into slant path wet delay SWD through a wet projection function and gradient values:

[0032] (1)

[0033] In the formula, represents a wet mapping function, ele and azi represent a satellite cutoff height angle and an azimuth angle respectively, and and are wet delay horizontal gradient terms in the south-north and east-west directions caused by anisotropy of the atmosphere, and represents a posteriori residual correction of an observation value.

[0034] Slant path wet delay SWD is further converted into slant path water vapor content SWV as a tomography input:

[0035] (2)

[0036] In the formula, represents a wet mapping function, ele and azi represent a satellite cutoff height angle and an azimuth angle respectively, is an atmospheric water vapor conversion coefficient:

[0037] (3)​

[0038] wherein, is the specific gas constant for water vapor, ; is the density of liquid water; is the weighted mean temperature of the atmosphere; is the adjusted wet term refractivity coefficient ; and for calculating and are the atmospheric refractivity constants; , , ; and are the molar masses of water vapor and dry air respectively , .

[0039] Step S2, selecting a GNSS-MET receiver station network with uniform distribution, obtaining its historical observation data and performing three-dimensional voxel grid division, converting the slant path water vapor content SWV into an integral form corresponding to the voxel grid, and constructing the observation equation of the tomography model.

[0040] First, the single-point observation of the GNSS-MET receiver station is obtained, and the satellite and meteorological observation information is obtained. When performing three-dimensional water vapor tomography, part of the stations are selected for network observation to construct a certain size voxel combination for the tomography area range, in order to construct the tomography model equation subsequently.

[0041] Specifically, a certain number of GNSS-MET receiver stations with relatively dense and uniform distribution are selected to obtain observation information, and at least one sounding station is included to perform vertical constraint and subsequent accuracy verification, and a tomography area including these stations is constructed. When performing voxel-based three-dimensional tomography, the tomography area is first divided into voxel blocks in a certain way, and the water vapor density or wet refractivity obtained by tomography finally corresponds to these voxel blocks. The horizontal division of voxels, i.e. the horizontal resolution of tomography, can be adjusted according to the distribution density of stations in the tomography area, while the vertical resolution adopts a non-uniform layering method, taking into account the need to introduce ERA5 as a prior constraint in the subsequent process, and the layering data of ERA5 is only near 37 height layers. Considering that the water vapor distribution should meet the condition that the water vapor content is large at low layer and the vertical resolution is high, and the water vapor content is small at high layer and the vertical resolution is low, the voxel is divided in the vertical direction based on the above two points. The obtained voxels are n in total, wherein , , and represent the number of rows, columns and layers of voxel blocks in the tomography area respectively.​

[0042] In this example, satellite and meteorological observation information of a certain area is acquired, 12 stations and a 45004 sounding station are selected for three-dimensional tomography, the tomography area is divided into voxel blocks in a certain way, and the water vapor density or wet refractive index obtained by tomography finally corresponds to these voxel blocks. The tomography range is selected as 22.19°~22.54° in the latitude direction, with a step of 0.05°, 113.87°~114.35° in the longitude direction, with a step of 0.06°, and in the vertical direction, considering the characteristics of the water vapor density being lower in the upper part and higher in the lower part, the height intervals are set as 300, 300, 300, 300, 500, 800, 1200, 1700, 2200, and 2800 m, and the research area is finally divided into 7 8 10 tomography grids. Tomography is performed in May 2013, and Figure 2 and Figure 3 are the three-dimensional water vapor density distribution at 0 o'clock on May 1 and the water vapor profile compared with the sounding station, respectively.

[0043] The slant path water vapor content SWV is the integral of the water vapor density on the signal ray path, which can be expressed as the sum of the intersection of each grid through which the satellite signal passes and the atmospheric water vapor density in the grid:

[0044] (4)

[0045] In the formula, represents the intersection of the ray in the grid, represents the estimated water vapor density value in the grid, and respectively represent the sequence numbers in the latitude, longitude, and elevation three-dimensional directions. All signal rays in the research area are expressed in the form of the above formula, and the observation equation of the tomography model can be obtained:

[0046] (5)

[0047] In the formula, y is a column vector composed of the water vapor content SWV on the ray passing out from the top of the research area, A is the coefficient matrix of the observation equation, which is generated by the grid numbers and intersections through which the signal ray passes, and x is a column matrix composed of unknown parameters water vapor density.

[0048] Step S3, impose horizontal constraint, vertical constraint and prior constraint on the observation equation of the tomographic model; construct the horizontal constraint according to the Gaussian distance weighting method according to the voxel division mode, construct the vertical constraint according to the exponential distribution of water vapor in the vertical direction, select the appropriate time ERA5 to calculate the water vapor density of each tomographic voxel block as the background field to impose the prior constraint. And calculate the water vapor density of the layer top from the historical ERA5 data as the boundary constraint, and integrate it into the vertical constraint, and finally construct the tomographic model equation in the networking tomographic area.

[0049] The water vapor density can be calculated from the historical reanalysis data ERA5 or sounding data. The difference lies in that the sounding data can directly provide water vapor pressure and temperature, while the ERA5 needs to obtain the water vapor pressure through temperature and relative humidity , unit or :

[0050] (6)

[0051] In the formula, represents the temperature, unit , and then the saturation water vapor pressure is combined with the relative humidity to calculate the water vapor pressure:

[0052] (7)

[0053] In the formula, represents the relative humidity, represents the water vapor pressure. With the water vapor pressure and temperature, the water vapor density of each layer can be further obtained:

[0054] (8)

[0055] In the formula, is the water vapor density of each layer height point, is the water vapor pressure of each layer height point, unit or .

[0056] According to the characteristics that the water vapor density has correlation in the horizontal direction and the closer the distance, the stronger the correlation, the Gaussian weighting function method is adopted to determine the coefficient value, and the closer the distance, the greater the weight value. The equation of imposing horizontal constraint condition is:

[0057] (9)

[0058] In the formula, and are the coefficient and value matrix of horizontal constraint respectively. According to the characteristics of the negative exponential correlation of water vapor density in the vertical direction, the coefficient value is determined according to the relationship of water vapor density between adjacent two layers of grid, and the vertical constraint condition equation is imposed on it:

[0059] (10)

[0060] wherein, and are the coefficient and value matrix of vertical constraint respectively. According to the sounding data of the previous few days before tomography, the boundary constraint can be imposed on the top layer of the tomography area, and the water vapor density and coefficient obtained by sounding are used to replace the last layer in the original vertical constraint. According to the ERA5, the water vapor density near each layer is calculated, and the PWV result of ERA5 at any time in the past is obtained. Match the GNSS-MET station to select several discontinuous PWV time points with the closest water vapor content, and use the average of the ERA5 stratified water vapor density at these time points to impose the prior constraint condition equation:

[0061] (11)

[0062] wherein, I is the prior constraint coefficient unit matrix, and C is the water vapor density value obtained by ERA5. Then, the combined observation equation after imposing various constraints is:

[0063] (12)

[0064] Step S4, the tomography area three-dimensional water vapor density of the combined observation equation is solved by using singular value decomposition SVD method, and the root mean square error, correlation coefficient, water vapor profile, etc. compared with the reanalysis data and sounding data are analyzed.

[0065] The tomography model equation is solved by SVD method, specifically, the coefficient matrix B on the left side of the combined observation equation is first decomposed by SVD:

[0066] (13)

[0067] wherein, is the matrix containing singular values, is the coefficient matrix of the observation equation is the square root of the eigenvalue, is the rank of the coefficient matrix B of the combined observation equation, and are the left singular and right singular vectors respectively. The generalized inverse of B is defined as:

[0068] (14)

[0069] Then the column matrix x of the unknown parameter water vapor density composition obtained by chromatography is combined with B and the observation value matrix L on the right side of the joint observation equation, and can be expressed as:

[0070] (15)

[0071] The RMS and STD at different height points at a certain latitude and longitude point are verified by comparison with the sounding and ERA5 data, and the calculation methods are as follows:

[0072] (16)

[0073] In the formula, refers to the number of layers of chromatography, represents the water vapor density value obtained by chromatography inversion, represents the water vapor density value calculated by sounding data or ERA5. Within the range of the voxel block where the sounding station is located, the water vapor profiles of the three are compared at the height of the latitude and longitude of the sounding station. Since ERA5 can obtain the water vapor density of each chromatography voxel block, the correlation between chromatography and ERA5 results can be analyzed. Selecting the time as May 1, 2013, the water vapor density values of all chromatography grids on May 1, 2013 can be obtained.

[0074] On the other hand, the application provides a GNSS-MET-based water vapor chromatography modeling device, which comprises various modules capable of realizing the various steps of the aforementioned method. Specifically, it comprises a calculation module for obtaining satellite signals and meteorological parameters based on GNSS-MET receivers and calculating slant path water vapor content (SWV); a voxel division module for obtaining historical observation data and performing three-dimensional voxel grid division, converting the slant path water vapor content into an integral form corresponding to the voxel grid to construct an observation equation of the chromatography model; a constraint module for applying horizontal constraints, vertical constraints and prior constraints to the observation equation to generate a joint observation equation; wherein the horizontal constraints are based on the Gaussian distance weighting method, the vertical constraints are based on the exponential distribution characteristics of water vapor, and the prior constraints are based on the water vapor density background field calculated from ERA5 reanalysis data; and an output module for solving the joint observation equation by singular value decomposition method to obtain a three-dimensional water vapor density distribution result.

[0075] In a third aspect, the application provides an electronic device, comprising: one or more processors; a memory for storing one or more programs; wherein when the one or more programs are executed by the one or more processors, the one or more processors implement the aforementioned GNSS-MET-based water vapor chromatography modeling method.

[0076] In a fourth aspect, the present application provides a computer readable storage medium having stored thereon executable instructions that, when executed by a processor, enable the processor to implement the above-mentioned GNSS-MET-based water vapor tomography modeling method.

[0077] The above-described specific embodiments further explain the purpose, technical solutions and beneficial effects of the present application. It should be understood that the above-described specific embodiments are merely examples of the present application and are not intended to limit the present application. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present application shall be included in the protection scope of the present application.

Claims

1. A GNSS-MET based water vapor tomography modeling method, characterized in that, The method comprises: Step S1, acquiring satellite signals and meteorological parameters based on a GNSS-MET receiver respectively, and solving slant path water vapor content SWV; Step S2, acquiring historical observation data and performing three-dimensional voxel grid division, converting the slant path water vapor content SWV into an integral form corresponding to the voxel grid, and constructing an observation equation of a tomographic model; Step S3, applying horizontal constraint, vertical constraint and prior constraint to the observation equation of the tomographic model to obtain a joint observation equation; wherein the horizontal constraint is based on a Gaussian distance weighting method, the vertical constraint is based on the exponential distribution characteristics of water vapor, and the prior constraint adopts a water vapor density background field calculated by a reanalysis data set, including: calculating the water vapor density near each layer based on ERA5 to obtain the PWV result of ERA5 at any time in the past, matching the PWV result with the GNSS-MET station, selecting a plurality of discontinuous PWV moments closest in water vapor content, and applying the ERA5 layered water vapor density average of the plurality of discontinuous PWV moments as the prior constraint condition; Step S4, solving the joint observation equation by a singular value decomposition method to obtain a three-dimensional water vapor density distribution.

2. The GNSS-MET based water vapor tomography modeling method according to claim 1, characterized in that, The solving of the slant path water vapor content SWV in the step S1 specifically comprises: converting zenith wet delay ZWD into slant path wet delay SWD by a wet projection function, and then converting the SWD into SWV in combination with an atmospheric water vapor conversion coefficient.

3. The GNSS-MET based water vapor tomography modeling method according to claim 1, characterized in that, The three-dimensional voxel grid division in the step S2 specifically comprises: adjusting the resolution according to the station distribution density in the horizontal direction, and adopting a non-uniform layering mode of lower density and higher density in the vertical direction.

4. The GNSS-MET based water vapor tomography modeling method according to claim 1, characterized in that, The specific implementation mode of the horizontal constraint in the step S3 is: based on the correlation of water vapor density in the horizontal direction, a Gaussian weighting function is used to determine the constraint coefficient, and the closer the grid points, the greater the weight.

5. The GNSS-MET based water vapor tomography modeling method according to claim 1, characterized in that, The specific implementation mode of the vertical constraint in the step S3 is: based on the characteristic that water vapor is negatively exponentially distributed in the vertical direction, a relationship equation of water vapor density between adjacent layers is established, and sounding data are introduced as top boundary constraints.

6. The GNSS-MET based water vapor tomography modeling method according to claim 1, characterized in that, The step S4 further comprises comparing and verifying the obtained three-dimensional water vapor density distribution with the reanalysis data set and the sounding data, calculating the root mean square error and the correlation coefficient, and drawing a water vapor profile.

7. A GNSS-MET based water vapor tomography modeling apparatus, characterized by, Comprise: A solving module configured to acquire satellite signals and meteorological parameters based on a GNSS-MET receiver respectively, and solve slant path water vapor content SWV; A voxel division module configured to acquire historical observation data and perform three-dimensional voxel grid division, convert the slant path water vapor content SWV into an integral form corresponding to the voxel grid, and construct an observation equation of a tomographic model; The constraint module is used for applying horizontal constraint, vertical constraint and prior constraint to the observation equation of the tomographic model to obtain a joint observation equation; wherein the horizontal constraint is based on a Gaussian distance weighting method, the vertical constraint is based on the exponential distribution characteristics of water vapor, and the prior constraint adopts a water vapor density background field calculated by a reanalysis data set, including: calculating the water vapor density near each layer by ERA5 to obtain the PWV result of ERA5 at any time in the past, matching the PWV result with GNSS-MET stations, selecting a plurality of discontinuous PWV time points with the closest water vapor content, and applying the ERA5 layered water vapor density average of the plurality of discontinuous PWV time points as the prior constraint condition; The output module is used for solving the joint observation equation by a singular value decomposition method to obtain a three-dimensional water vapor density distribution.

8. An electronic device, comprising: It comprises: one or more processors; a memory for storing one or more programs; wherein when the one or more programs are executed by the one or more processors, the one or more processors implement the GNSS-MET-based water vapor tomographic modeling method of any one of claims 1-6.

9. A computer-readable storage medium, characterized in that, a processor, and a memory having stored executable instructions which, when executed by the processor, cause the processor to implement the GNSS-MET-based water vapor tomographic modeling method of any one of claims 1-6.

Images