Optimal interpolation assimilation method for self-adaptive observation station space distribution
By dynamically adjusting the influence radius of the weight function of the Cressman interpolation scheme, the problem of poor assimilation effect caused by uneven distribution of observation stations was solved, and high-precision meteorological reanalysis and numerical weather prediction were achieved.
Patent Information
- Application Number
- CN202511313852.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-09-15
- Publication Date
- 2025-12-30
AI Technical Summary
Existing optimal interpolation assimilation methods are not effective in areas where the spatial distribution of observation data is uneven or the data is sparse. In particular, in areas where meteorological observation stations are unevenly distributed, it is difficult to accurately capture local climate change characteristics, resulting in low quality of the analysis field.
By dynamically adjusting the influence radius of the weight function in the Cressman interpolation scheme, and based on the spatial distribution and density of observation stations and their relative positional relationship with the target grid points, an optimal interpolation assimilation method that adapts to the spatial distribution of observation stations is realized. Observation operators are calculated to improve the preservation of local details and the representativeness of the data.
It improves the accuracy of meteorological reanalysis data and numerical weather prediction, and is applicable to the construction or reconstruction of datasets in the fields of meteorology, energy and navigation.
Smart Images
Figure CN121234576A_ABST
Abstract
Description
TECHNICAL FIELD
[0001] The present application relates to the technical field of numerical weather prediction, and particularly relates to an optimal interpolation assimilation method of self-adaptive observation site spatial distribution. BACKGROUND
[0002] Numerical weather prediction is an initial value problem, and the quality of the initial field directly determines the accuracy of the prediction results (Kalnay 2003). Data assimilation is an important means to provide high-quality initial fields for numerical weather prediction models, and it is not only a key factor affecting the level of weather prediction, but also an important way to build high-precision and high-temporal and spatial resolution atmospheric reanalysis data. The development of data assimilation methods has experienced the evolution from simple to complex, from linear to nonlinear, and from deterministic to probabilistic, forming a series of different types of assimilation methods, which can be generally divided into three categories: optimal interpolation (OI), three-dimensional / four-dimensional variational, and ensemble Kalman filter (Mathieu et al. 2008; Samos et al. 2023). Optimal interpolation is a classic data assimilation method based on the least square root error theory, and its core idea is to fuse the observation field (such as spatially discrete distributed meteorological observation station data) and the background field (such as the short-term prediction results on the regular grid points of the atmospheric numerical model) to obtain the analysis field (i.e., the optimal estimation result of the atmospheric state) (Heo et al. 2018); the variational method can fuse more nonlinear information and achieve higher accuracy of assimilation by introducing temporal and spatial constraints through the optimal theory (Zhang et al. 2020). The ensemble Kalman filter simulates the error statistical characteristics through a set of ensemble members, which is suitable for simulating the uncertainty of nonlinear systems (Chang et al. 2024). Although the variational method and the ensemble Kalman filter method are relatively advanced data assimilation methods, they can better handle the assimilation of nonlinear data, but they also have high computational cost, the objective function is easy to fall into local extremum, and the error caused by the divergence of the ensemble sample. The optimal interpolation not only has high computational efficiency and reliable computational accuracy, but also has an essentially equivalent mathematical expression with the variational assimilation method (Eyre et al. 2022; Chang et al. 2024), so it has a significant advantage in the construction of meteorological reanalysis data sets (Lussana et al. 2009; Gusarov et al. 2017; Fiedler et al. 2019). For example, based on the optimal interpolation, the initialization of sea surface temperature and sea ice can be improved, and higher precision and short-term reliability can be optimized in areas with high variability of sea surface temperature (Fiedler et al. 2019). The assimilation of satellite observation data by optimal interpolation can improve the accuracy of numerical weather prediction (Zhang et al. 2017). Optimal interpolation is also used to fuse meteorological observation data from different sources (such as sounding observations and satellite observations) to improve regional climate simulation results (Lindskog et al. 2019).
[0003] However, existing optimal interpolation assimilation methods also have certain limitations, especially in areas with uneven spatial distribution of observation data or sparse data, where the assimilation effect is generally poor (Miatselskaya et al. 2022). This is mainly related to the expression of the observation operator in the optimal interpolation method. Generally, the observation operator uses bilinear interpolation, kriging interpolation, and Cressman interpolation. When different interpolation methods are used as observation operators, the optimal interpolation assimilation method performs differently. Among them, bilinear interpolation is computationally simple, but it cannot fully capture the spatial correlation between observation data (Zhang et al. 2024); kriging interpolation and Cressman interpolation can overcome the shortcomings of bilinear interpolation better (Zhang et al. 2015). Compared with kriging interpolation, Cressman interpolation, by introducing a weighted distance function, can achieve a smoother and more reasonable spatial field estimation, and can balance computational efficiency and error estimation. Therefore, it is widely used as the observation operator for optimal interpolation assimilation methods (Liu et al. 2016). The core of traditional Cressman interpolation is to calculate a weight function based on the spatial distance between observation points and grid points, thereby integrating the observations into the grid. It's important to note that traditional Cressman interpolation relies on the assumption that the radius of influence in the weight function is constant. This assumption holds true in ideal conditions where observation stations are spatially evenly distributed. However, the spatial distribution of meteorological observation stations is actually highly uneven, especially in areas with complex and variable topography and underlying surfaces. Assuming a constant radius of influence can lead to optimal interpolation assimilation methods failing to accurately capture local climate variations (Wang et al. 2023), resulting in poor assimilation performance. In my country, meteorological observation stations are densely distributed in the eastern region but relatively sparse in the western and mountainous areas. This inevitably leads to lower quality analytical fields or meteorological reanalysis datasets obtained using optimal interpolation assimilation methods. Therefore, dynamically adjusting the spatial radius of influence in the weight function of traditional Cressman interpolation based on the spatial distribution, density, and relative position of meteorological observation stations to the target grid points, thereby improving the observation operator, is a potential way to improve optimal interpolation assimilation methods and enhance the quality of meteorological reanalysis datasets. Summary of the Invention
[0004] The technical problem to be solved by this invention is to provide a simple and effective optimal interpolation assimilation method for adaptive observation station spatial distribution of high-precision atmospheric reanalysis data.
[0005] To address the aforementioned problems, the present invention provides an optimal interpolation assimilation method for adaptive spatial distribution of observation stations, comprising the following steps:
[0006] (1) Obtain the background error covariance matrix B, the observation error covariance matrix R, the optimal weight matrix W, and the k-th observation station y. o (k) and background field And calculate the observation increment.
[0007] (2) Perform the first iteration on the background field, with grid point (i,j) as the center and R c (i, j) is the radius pair By performing a weighted average, we obtain the observation operator for the initial iteration at grid point (i,j) in the background field.
[0008]
[0009] In the formula: i, j represent the grid point numbers in the latitudinal and meridional directions, respectively; K represents the total number of observation stations; R c (i, j) represents the influence radius that dynamically changes with the spatial distribution of the observation stations; w ijk It is the weight function in the Cressman interpolation scheme; d ijk This represents the spatial distance between grid point (i,j) and observation station k; Centered on grid point (i,j), with radius R c Number of observation stations within (i, j); initial radius R start =1km; Upper limit radius R max =200km; step length ΔR = 1km; N min This is a preset value for the number of observation stations within the search radius;
[0010] (3) Based on the observation operator of the initial iteration The background field of the first generation iteration was calculated.
[0011]
[0012] This process is repeated γ times until the observation increment is reached. The absolute value is less than or equal to ε and ε = 1 × 10 -6 When the iteration stops, the observation operator that dynamically changes with space is obtained. Should Used for analyzing field x a The calculation.
[0013] The observation increment in step (1) Calculate using the following formula:
[0014]
[0015] In step (2), N minThe value of is determined by comprehensively evaluating and continuously testing the analytical field through statistical measures including absolute error, relative error, spatial correlation coefficient, and temporal correlation coefficient, based on the principle of minimizing errors in time and space and maximizing correlation coefficients.
[0016] Compared with the prior art, the present invention has the following advantages:
[0017] 1. This invention dynamically adjusts the influence radius of the Cressman interpolation scheme weight function based on the spatial distribution and density of observation stations and their relative positional relationship with target grid points in the background field. This enables the observation operator to dynamically change with the spatial distribution of observation stations, forming an optimal interpolation assimilation method that adapts to the spatial distribution of observation stations. This maximizes the preservation of local meteorological information details while increasing data representativeness, thereby improving the effectiveness of the optimal interpolation assimilation method and enhancing the simulation results of numerical models.
[0018] 2. The method of this invention is simple and can effectively construct high-precision atmospheric reanalysis data, improve the accuracy of the initial field of numerical weather prediction models, and is applicable to the construction or reconstruction of datasets in meteorology, new energy, navigation and other fields. Attached Figure Description
[0019] The specific embodiments of the present invention will be described in further detail below with reference to the accompanying drawings.
[0020] Fig. 1 The inner circle represents the average relative error of the old OI assimilation method and the new OI assimilation method of this invention relative to the observations, while the outer circle represents the average relative error of the old method. Wherein: a and ah represent different regions of China, namely a: overall; b: Northeast; c: North China; d: Northwest; e: Southwest; f: Central China; g: South China; h: East China.
[0021] Fig. 2 The average relative errors for (a) relative humidity, (b) precipitation, (c) temperature and (d) wind speed in different months in China; blue represents the original method; red represents the new method of this invention. Detailed Implementation
[0022] An optimal interpolation assimilation method for adaptive spatial distribution of observation stations includes the following steps:
[0023] (1) Obtain the background error covariance matrix B, the observation error covariance matrix R, the optimal weight matrix W, and the k-th observation station y. o (k) and background field And calculate the observation increment.
[0024]
[0025] (2) Perform the first iteration on the background field, with grid point (i,j) as the center and R c (i, j) is the radius pair By performing a weighted average, we obtain the observation operator for the initial iteration at grid point (i,j) in the background field.
[0026]
[0027]
[0028] In the formula: i, j represent the grid point numbers in the latitudinal and meridional directions, respectively; K represents the total number of observation stations; R c (i, j) represents the influence radius that dynamically changes with the spatial distribution of the observation stations; w ijk It is the weight function in the Cressman interpolation scheme; d ijk This represents the spatial distance between grid point (i,j) and observation station k; Centered on grid point (i,j), with radius R c Number of observation stations within (i, j); initial radius R start =1km; Upper limit radius R max =200km; step length ΔR = 1km; N min N is a preset value for the number of observation stations within the search radius; min The value of is determined by comprehensively evaluating and continuously testing the analytical field through statistical measures including absolute error, relative error, spatial correlation coefficient, and temporal correlation coefficient, based on the principle of minimizing errors in time and space and maximizing correlation coefficients.
[0029] (3) Based on the observation operator of the initial iteration The background field of the first generation iteration was calculated.
[0030]
[0031] This process is repeated γ times until the observation increment is reached. The absolute value is less than or equal to ε and ε = 1 × 10 -6 When the iteration stops, the observation operator that dynamically changes with space is obtained. Should Used for analyzing field x a The calculation.
[0032] The specific analytical process of this invention is as follows:
[0033] 1. The objective function of the optimal interpolation assimilation method is defined as follows:
[0034] x a =xb +W[y0-H(x b (1)
[0035] In the formula: x a For the analysis field, x b The background field consists of two m×n matrices (number of spatial latitudinal and meridional grid points); y o The observation data is represented by a vector of length p (the total number of meteorological stations); W is the optimal weight matrix; and H is the observation operator, which transforms the observation data from irregular meteorological stations to regular grid points in the background field through spatial interpolation. The effect of optimal interpolation assimilation, i.e., the analysis field x... a The quality of the observed data y o and background field x b Given a fixed value, it mainly depends on the mathematical expressions of the optimal weight matrix W and the observation operator H, and the optimal weight matrix W is calculated in step 2.
[0036] 2. The optimal weight matrix is defined as follows:
[0037] W=BH T (HBH T +R) -1 (2)
[0038] In the formula: B is the background error covariance matrix; R is the observation error covariance matrix.
[0039] The optimal weight matrix W is jointly determined by the observation operator H, the background error covariance matrix B, and the observation error covariance matrix R. The formulas for calculating the background error covariance matrix B and the observation error covariance matrix R are as follows:
[0040]
[0041] Where: ε b ε0 represents the vector consisting of the variance of the background field grid points themselves and the covariance between grid points over a period of time (e.g., one month); ε0 represents the vector consisting of the variance of the meteorological station itself and the covariance between meteorological stations over a period of time (e.g., one month); the superscript T indicates the transpose of the vector; E represents a two-dimensional matrix. Therefore, in the observed data y o and background field x b Given a fixed set of conditions, the background error covariance matrix B and the observation error covariance matrix R are also fixed over a period of time (e.g., one month). Therefore, the optimal weight matrix W is still primarily determined by the observation operator H. Combined with step 1, it can be seen that the optimal interpolation assimilation effect is ultimately determined by the observation operator H.
[0042] 3. Calculation and improvement of the observation operator H:
[0043] In this invention, the Cressman interpolation scheme is selected as the observation operator for the optimal interpolation assimilation method. The basic principle of the Cressman interpolation scheme is iterative updating, and its calculation steps are as follows:
[0044] (I) Background field (x) b ) as the initial guessing field (denoted as ), calculate the relationship between the k-th observation station and the initial guess field. difference
[0045]
[0046] (II) The Cressman interpolation scheme shown in Equation (6) is adopted as the observation operator, with the grid point (i,j) as the center and r c For radius pair We perform a weighted average to obtain the observation operator for the initial iteration at grid point (i,j) in the background field.
[0047]
[0048] In the formula: i and j represent the grid point numbers in the latitudinal and meridional directions, respectively; K represents the total number of observation stations; and r represents the influence radius in the traditional Cressman interpolation scheme. c w is a constant ijk This is the weight function in the Cressman interpolation scheme, and its expression is as follows:
[0049]
[0050] In the formula: d ijk This represents the spatial distance between grid point (i,j) and observation station k.
[0051] Equation (7) shows that the weight function w ijk It depends on the spatial distance d between the observation station k and the grid point (i,j). ijk With influence radius r c Relationship: If d ijk ≤r c That is, the observation station within the influence radius r c If the weight function w is within the range, then ijk It is a non-zero value, and d ijk The smaller the value, the lower the weight function w. ijk The larger the value of d, the greater the value of d. ijk >r c That is, the observation station within the influence radius r c Outside the range, the weight function w ijk The value of is zero, indicating that the observation station has no effect on the target grid point.
[0052] (III) Using formula (8), the background field of the first iteration is calculated.
[0053]
[0054] Will Substituting into formulas (5) to (7), we obtain the first iteration. By repeating steps (I) to (III) a total of γ times, we can obtain as well as Until the observed increment in formula (5) The absolute value is less than a small quantity (i.e. ε=1×10 -6 At this point, what will be obtained As the final observation operator, substituting it into formulas (1) to (2) in steps 1-2 completes the calculation of the optimal interpolation assimilation method, yielding the final analysis field x. a (Optimal estimate).
[0055] (IV) The radius of influence r in step (II) c The general assumption is that the value is constant. In areas with dense observation data, this can lead to issues with smoothing details in the analysis field, while in areas with sparse observation data (such as mountainous or plateau regions), it can reduce the representativeness of the data. Both of these negatively impact the assimilation results. To overcome this problem, this invention uses the following method to establish the influence radius R, which dynamically varies with the spatial distribution of observation stations. c (i, j), the specific steps are as follows:
[0056] (a) Using formula (9), for background field grid point (i,j), assuming an initial radius R start =1km, step length ΔR =1km, upper limit radius R max =200km, gradually expand the search range, and iterate through the number of observation stations around the search grid point (i,j). when When the radius R is such that the radius R is the influence radius at grid point (i,j), then the radius R is the influence radius at that grid point (i,j). c (i, j).
[0057]
[0058] In the formula: Centered on grid point (i,j), with radius R c The number of observation stations within (i, j); the upper limit radius is set to 200 km to reduce the computational load; N min The preset value for the number of observation stations within the search radius is determined by the following steps (c).
[0059] (b) The influence radius R obtained in (a) varies dynamically with space. c (i, j) replace the fixed influence radius r in step (II) c The new weight function W can be obtained. ijk :
[0060]
[0061] Replace w in formula (6) with formula (10) ijk The γth time obtained according to steps (I) to (III) This refers to the observation operator that dynamically changes with space. In this invention, "dynamically changing with space" refers to R in formula (10). c (i, j) is dynamically changing, while r in formula (7) c It is a fixed constant.
[0062] (c) By using statistical measures such as absolute error, relative error, spatial correlation coefficient, and temporal correlation coefficient, the analytical field is comprehensively evaluated and continuously tested. The optimal N in formula (9) is determined based on the principle of minimizing temporal and spatial errors and maximizing correlation coefficients. min The value of .
[0063] The influence radius R obtained in step (IV) varies with spatial dynamics. c Substituting (i, j) and the observation operator into steps 1-2 yields the optimal interpolation assimilation method for adaptive spatial distribution of observation stations. At this point, the influence radius of the observation operator in the optimal interpolation assimilation method is determined by the number of observation stations surrounding grid point (i, j). Dynamically adjusted: In areas with dense observation stations, R c A smaller (i, j) value improves the assimilation accuracy of local details; in areas with sparse observation stations, a larger R value is used. c The (i, j) values are used to incorporate observation information from further distances (but still subject to weight decay constraints) to correct the background field, avoiding situations where effective correction is impossible due to a lack of local observation data.
[0064] The observation stations in this embodiment are fixed or mobile monitoring points used to directly or indirectly measure atmospheric conditions. In this invention, the observation stations are surface meteorological stations in China provided by the China Meteorological Administration, containing data such as 10m wind speed, 2m air temperature, 2m humidity, and precipitation, which will serve as the observation field in the optimal interpolation assimilation method.
[0065] In traditional optimal interpolation assimilation methods, the influence radius in the weight function of the observation operator (Cressman interpolation scheme) is static and assumed to be a constant. This invention dynamically adjusts the influence radius in the weight function based on the spatial distribution and density of observation stations around the grid points, combined with the error of the analysis field. This allows the observation operator to dynamically change with the spatial distribution of observation stations, thereby improving the effectiveness of the optimal interpolation assimilation method. An assimilation experiment conducted in China from January to December 2013 is used as an example. Figs. 1-2 As can be seen, compared with the traditional optimal interpolation method (static observation operator), the improved optimal interpolation method of this invention (dynamically adjusted observation operator) significantly reduces the simulation error of different regions, times and meteorological elements in China, indicating that this invention has high universality and significant advantages.
Claims
1. A method of optimal interpolation assimilation of adaptive observation site spatial distribution, comprising the steps of: (1) Obtain background error covariance matrix B, observation error covariance matrix R, optimal weight matrix W, the kth observation site y o (k) and background field and calculate the observation increment (2) The first iteration is performed on the background field, and the observation operator at the grid point (i, j) in the background field is obtained by weighted average of the observation operators at the grid points (i, j) and (i, j) with a radius of R c (i, j) as the center The observation operator at the grid point (i, j) in the background field is obtained by weighted average of the observation operators at the grid points (i, j) and (i, j) with a radius of R where i, j represent the grid point number in the latitudinal and longitudinal direction respectively; K represents the total number of observation sites; R c (i, j) is the influence radius which dynamically changes with the spatial distribution of observation sites; w ijk is the weight function in the Cressman interpolation scheme; d ijk represents the spatial distance between the grid point (i, j) and the observation site k; is the number of observation sites within the radius R c (i, j) centered at the grid point (i, j); the initial radius R start = 1 km; the upper limit radius R max = 200 km; the step size ΔR = 1 km; N min is the preset value of the number of observation sites within the search radius; 3. The observation operator according to the initial iteration The background field of the first iteration is computed And so on for γ iterations until the absolute value of the observation increment is less than or equal to ε and ε = 1 x 10 -6 the iteration is stopped, and the observation operator that varies dynamically with space is obtained The is used to analyze the field x a calculation.
2. The optimal interpolation assimilation method of adaptive observation site spatial distribution according to claim 1, characterized in that: The step (1) of observing the increment is calculated by the following equation:
3. The optimal interpolation assimilation method of adaptive observation site spatial distribution according to claim 1, characterized in that: The value of N in step 2 min is determined by statistical quantities including absolute error, relative error, spatial correlation coefficient, and temporal correlation coefficient, through comprehensive evaluation and continuous testing of the analysis field, with the principles of minimum error and highest correlation coefficient in time and space.