A multi-source marine data fusion method based on multi-scale optimal interpolation
By constructing the background error covariance matrix using a multi-scale optimal interpolation algorithm and Kronecker product, the problems of high computational resource consumption and low efficiency in existing technologies are solved. This enables efficient fusion of multi-source and multi-scale ocean observation data, generating high-resolution data products that are continuous in time and regular in space.
Patent Information
- Application Number
- CN202510567994.3
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-04-30
- Publication Date
- 2026-03-03
- Estimated Expiration
- 2045-04-30
AI Technical Summary
Existing optimal interpolation algorithms cannot effectively handle observation operators when fusing multi-source observation data, especially satellite observation data with different spatial resolutions. They also consume a lot of computational resources, have low computational efficiency, and are difficult to achieve multi-scale data fusion.
A method based on multi-scale optimal interpolation is adopted, and the background error covariance matrix is constructed by Kronecker product. Combined with the mathematical method of Kronecker product, the background error covariance matrix is constructed, which reduces the amount of computation, supports the calculation of non-uniform anisotropic background error matrix, and realizes efficient fusion of multi-source and multi-scale observation data.
It reduces computational load and computer memory requirements, supports the calculation of non-uniform anisotropic background error matrix, and can directly fuse ocean observation data of different resolutions, successively reducing observation errors and achieving effective fusion of multi-source ocean observation data.
Smart Images

Figure CN121009483B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of multi-source marine data fusion, and in particular to a multi-source marine data fusion method based on multi-scale optimal interpolation. Background Technology
[0002] Activities related to national welfare and people's livelihood, such as marine fisheries production, energy development, and shipping, as well as maritime operational safety assurance, search and rescue for catastrophic accidents such as oil spills, and pollution control, all urgently require high-precision and high-resolution sea state data. In the past 20 years, marine observation technology has made breakthroughs, and the amount of observation data has increased explosively. The main observations can be classified into two categories: the first is in-situ observation relying on surface and underwater observation platforms; the second is remote sensing observation relying on satellites, aircraft, and shore-based radar. The sensors of this type of remote sensing observation are far from the seawater, and invert sea surface height, seawater temperature, and salinity by receiving the radiation, reflection, or backscattering signals of electromagnetic waves from the water surface. All of these observations have extremely irregular spatiotemporal distributions, different spatiotemporal resolutions, and different observation errors, making them difficult to directly apply to solving practical problems, let alone to the study of ocean dynamics.
[0003] Satellite altimeter observations are considered revolutionary for the development of ocean science. However, the spatiotemporal distribution of satellite altimeter observations is extremely irregular. They can only acquire sea surface height data at points below the satellite's orbit. Along the satellite's orbit, the altimeter generates an observation every 7 km, but the orbital interval of a single satellite is greater than 100 km, and the orbital period varies from 10 to 30 days. This results in discontinuity and non-uniformity in time and space, greatly limiting the direct application of satellite altimeter data in ocean research. The method of fusing data from multiple satellite altimeters can eliminate this limitation. In particular, there are currently more than ten altimeter satellites in the world. By fusing ocean altimeter data from multiple satellites and eliminating missing values in the satellite altimeter data, regular gridded data products that can be analyzed and used are generated, which is the key to effectively utilizing satellite ocean altimeter data.
[0004] Satellite sea surface temperature (SST) is another important type of ocean observation. Satellite remote sensing can acquire large-scale, near-real-time sea surface temperature data, but it can only acquire SST within the width of the scan swath. Gaps exist between different swaths, and the spatiotemporal resolution of SST acquired by different sensors varies greatly. Currently, satellite sea surface temperature remote sensing mainly uses two bands: microwave and infrared. The spatial resolution of the microwave and infrared bands is extremely different. The microwave band has a longer wavelength, at the centimeter level, with a spatial resolution of 25 kilometers for observed SST. The infrared band has a shorter wavelength, at the micrometer to millimeter level, with a spatial resolution of 1 kilometer or even higher for observed SST. In addition, although infrared SST remote sensing has high resolution, it cannot penetrate clouds, rain, and fog. Areas with clouds, rain, and fog cannot be observed with infrared SST, and infrared SST remote sensing suffers from large errors due to the influence of aerosols and other factors. Therefore, the application potential of raw SST data generated by satellite remote sensing is greatly limited. To fully realize its application potential, it is necessary to generate high-resolution SST data products that are continuous in time and have regular spatial grids.
[0005] Currently, utilizing algorithms to fuse various observational data to generate high-resolution, high-quality analytical fields is one of the fundamental approaches to effectively utilize ocean satellite data. Among these methods, the traditional optimal interpolation (OI) algorithm is the most commonly used. For example, CN116340800A discloses a method and apparatus for assimilating ocean surface data based on ensemble optimal interpolation of seasonal samples. The method includes: acquiring four sets of seasonal ensemble sample data corresponding to ocean surface observation data; calculating the background error of each seasonal ensemble sample data based on the ensemble sample mean of each seasonal sample data; determining the optimal interpolation method for each seasonal ensemble sample data based on the background error; and assimilating the ocean surface observation data based on the optimal interpolation method to obtain assimilated ocean surface data corresponding to each seasonal ensemble sample data. Thus, the optimal interpolation method determined based on the background error of each seasonal ensemble sample achieves accurate assimilation of ocean surface observation data, increasing the quantity and spatiotemporal uniformity of ocean surface observation data. However, this method does not consider data from multiple sources with different resolutions and cannot achieve multi-scale data fusion. Existing optimal interpolation algorithms cannot effectively handle observation operators when fusing multi-source observation data, especially satellite observation data with different spatial resolutions, thus failing to effectively fuse data of different spatial resolutions. Furthermore, existing optimal interpolation algorithms consume significant computational resources when calculating the background field error covariance matrix, resulting in high memory requirements and limited computational efficiency. Summary of the Invention
[0006] The purpose of this invention is to provide a multi-source ocean data fusion method based on multi-scale optimal interpolation, so as to achieve efficient fusion of multi-source and multi-scale observation data.
[0007] The objective of this invention can be achieved through the following technical solutions:
[0008] A multi-source ocean data fusion method based on multi-scale optimal interpolation includes the following steps:
[0009] Step 1) Data Acquisition: Acquire multi-source ocean data, including static data and real-time dynamic data. The static data includes global high-resolution land-sea distribution data and global climate data, while the real-time dynamic data includes multi-source observation data at different resolutions.
[0010] Step 2) Constructing the regional grid: Based on the acquired multi-source ocean data and user requirements, determine the starting and ending latitude and longitude coordinates, set the grid resolution, generate the regional grid, and mark the land and sea distribution of the grid;
[0011] Step 3) Generate the regional climatological background field: Interpolate global climatological data onto the regional grid to generate the regional climatological background field;
[0012] Step 4) Process low-resolution observation data: Calculate the location index of the low-resolution observation data in the regional grid and remove invalid data located in this region;
[0013] Step 5) Calculate the regional low-resolution analysis field: Use a multi-scale optimal interpolation algorithm to fuse the regional climatological background field and low-resolution observation data to generate the regional low-resolution analysis field;
[0014] Step 6) Processing high-resolution observation data: Calculate the location index of the high-resolution observation data in the regional grid and remove invalid data located in this region;
[0015] Step 7) Calculate the high-resolution analysis field of the region: Use a multi-scale optimal interpolation algorithm to fuse the low-resolution analysis field of the region with the high-resolution observation data to generate the high-resolution analysis field of the region;
[0016] The multi-scale optimal interpolation algorithm constructs the background field error covariance matrix based on the Kronecker product, thereby achieving optimal interpolation and simultaneous calculation.
[0017] The real-time dynamic data includes at least low-resolution microwave remote sensing data and high-resolution near-infrared remote sensing data.
[0018] Step 2) specifically includes the following steps:
[0019] Step 2.1) Define the starting and ending latitude and longitude coordinates of the region grid: Set the starting and ending latitude and longitude coordinates of the region grid in degrees;
[0020] Step 2.2) Define the grid resolution: Set the resolution of the regional grid, with the unit being degrees;
[0021] Step 2.3) Generate the regional grid: Using the starting and ending longitude and latitude coordinates of the regional grid defined in Step 2.1) and the grid resolution defined in Step 2.2), generate longitude vectors and latitude vectors. The longitude vector has N columns, and the latitude vector has M rows. Copy the longitude vector M times and copy the latitude vector N columns to generate the regional grid;
[0022] Step 2.4) Mark the land-sea distribution of the regional grid: Use the global high-resolution land-sea distribution data obtained in Step 1) to mark the land-sea distribution in the regional grid, where 0 represents land and 1 represents ocean.
[0023] The specific steps of Step 3) are as follows:
[0024] Step 3.1) Calculate the regional climate state background field: Read the global climate state data and interpolate the ocean elements onto the regional grid generated in Step 2) using the bilinear interpolation method. The specific interpolation method is as follows: Assume that the target point is to the left of the abscissa x, and find the two nearest abscissas x1 and x2, where x1 < x2. According to the relative positions of the abscissa of the target point and the two nearest neighbor points, use the bilinear interpolation formula to calculate the interpolation result of the target point:
[0025] S = (1 - α) * S(x1, y) + α * S(x2, y)
[0026] where S(x1, y) and S(x2, y) are the values of the two nearest neighbor points respectively, and α is the relative position of the target point, and the calculation method is
[0027] Step 3.2) Mark the land-sea distribution of the regional climate state background field: Use the land-sea distribution of the grid marked in Step 2) as the mask variable, and multiply the mask variable by the regional climate state background field to make the ocean element values at the land grid points be 0.
[0028] The specific steps of Step 4) are as follows:
[0029] Step 4.1) Eliminate the small-scale signals of the low-resolution observation data: Calculate the average value of the obtained low-resolution observation data within a preset time interval to eliminate the small-scale signals in the low-resolution observation data;
[0030] Step 4.2) Obtain the regional low-resolution observation data: Intercept the low-resolution observation data located in the regional grid. Use the first variable to represent the longitude coordinates of the regional low-resolution observation data, the second variable to represent the latitude coordinates of the regional low-resolution observation data, and the third variable to represent the ocean element values in the regional low-resolution observation data;
[0031] Step 4.3) Remove invalid values from low-resolution observation data in the region: Clean the third variable to remove outliers and invalid values;
[0032] Step 4.4) Find the nearest grid point of low-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data falling within the regional grid, and use it as the location index.
[0033] Step 5) specifically includes the following steps:
[0034] Step 5.1) Set background error: Specify the error value for the climatological background field data at each grid point in the regional grid;
[0035] Step 5.2) Calculate the observation error: The observation error includes the measurement error and the sampling error. The measurement error of the low-resolution observation data is provided by the dataset, and the sampling error is a preset percentage of the measurement error. The observation error of each observation point is stored in the observation error matrix.
[0036] Step 5.3) Calculate the observation increment: Calculate the difference between the regional low-resolution observation data mapped to the regional grid points and the regional climatological background field to obtain the first observation increment. The calculation formula is as follows:
[0037] d low =y o_low_res -H low x b
[0038] Among them, y o_low_res For the region's low-resolution observation data vector, x b H serves as the regional climatological background field. low For the first observation operator, d low This is the first observation increment;
[0039] Step 5.4) Calculate the regional low-resolution analysis field: Using a multi-scale optimal interpolation algorithm, fuse the regional climatological background field with the regional low-resolution observation data to generate the regional low-resolution analysis field x. a_low_res The calculation formula is as follows:
[0040]
[0041] Among them, R low x is the first observation error covariance; a_low_res For the region's low-resolution analysis field; B low The background field error covariance, representing the regional climatological background field, determines the weighting between observational and background field information.
[0042] Step 6) specifically includes the following steps:
[0043] Step 6.1) Obtain regional high-resolution observation data: Extract high-resolution observation data located in the regional grid, use the fourth variable to represent the longitude coordinates of the regional high-resolution observation data, use the fifth variable to represent the latitude coordinates of the regional high-resolution observation data, and use the sixth variable to represent the ocean element values in the regional high-resolution observation data.
[0044] Step 6.2) Remove invalid values from high-resolution regional observation data: Clean the sixth variable to remove outliers and invalid values;
[0045] Step 6.3) Find the nearest grid point of the high-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data falling within the regional grid, and use it as the location index.
[0046] Step 7) specifically includes the following steps:
[0047] Step 7.1) Set background error: Specify the error value of the low-resolution analysis field data of the region at each grid point in the region grid;
[0048] Step 7.2) Calculate the observation error: The observation error includes measurement error and sampling error. The measurement error of high-resolution observation data is provided by the dataset. The sampling error is a preset percentage of the measurement error. The observation error of each observation point is stored in the observation error matrix in the form of elements.
[0049] Step 7.3) Calculate the observation increment: Calculate the difference between the high-resolution observation data of the region mapped to the regional grid points and the low-resolution analysis field of the region to obtain the second observation increment. The calculation formula is as follows:
[0050] d high =y o_high_res -H high x a_low_res
[0051] Among them, y o_high_res x is the region's high-resolution observation vector. a_low_res For the region's low-resolution analysis field, H high For the second observation operator, d high This is the second observation increment;
[0052] Step 7.4) Calculate the multi-scale regional high-resolution analysis field: Using the regional low-resolution analysis field as the background field, and employing a multi-scale optimal interpolation algorithm, fuse the regional high-resolution observation data to generate the multi-scale regional high-resolution analysis field x. a _high_res The calculation formula is as follows:
[0053]
[0054] Among them, R high x is the covariance of the second observation error; a_high_res For high-resolution analysis of multi-scale regions; B high The background field error covariance, which serves as the background field for the low-resolution analysis field in the region, determines the weighting between the observation and background field information.
[0055] The distance calculation method of the minimum distance method is as follows:
[0056]
[0057] Among them, Dis AB Let A be the distance between points A and B, where the coordinates of point A are (x1, y1) and the coordinates of point B are (x2, y2). The background field error covariance is decomposed as follows during calculation:
[0058] B=ΣCΣ
[0059] Where B is the background field error covariance of the regional climatological background field. low Or, the background field error covariance B is used as the background field in a low-resolution analysis field of the region. high Σ is a diagonal matrix, with diagonal elements representing the root mean square error of the background error at each grid point. C is the correlation coefficient matrix, which propagates the observed information to the surrounding blank area and consists of the error correlation coefficients between grid points in the background field. The formula is as follows:
[0060]
[0061] in, C represents the Kronecker product. x Let C be the correlation coefficient matrix in the x-direction. y Let be the correlation coefficient matrix in the y-direction;
[0062] Assuming the spatial grid size is m×n, then C x With C y The calculation matrix forms are as follows:
[0063]
[0064]
[0065] Among them, C x With C y Represented using a Gaussian function:
[0066]
[0067] Where deltx is the distance between grid points in the x-direction, and delty is the distance between grid points in the y-direction; L x L is the length scale of the correlation coefficient in the x-direction. y The length scale of the correlation coefficient in the y-direction is expressed in degrees. The larger the length scale of the correlation coefficient, the greater the distance that the observation information travels.
[0068] Compared with the prior art, the present invention has the following beneficial effects:
[0069] (1) In the OI algorithm, the background error covariance matrix B accounts for the largest part of the computational load. This invention is based on the multi-scale optimal interpolation algorithm. By introducing the mathematical method of Kronecker product, the background error covariance matrix B is constructed, which makes B a symmetric positive definite matrix. Through matrix transformation, the (MN) in the traditional OI algorithm is reduced. 2 The computational cost of the background error covariance is reduced to M. 2 +N 2 The computational complexity of the size (M is the number of grid cells in the x-direction, and N is the number of grid cells in the y-direction) is reduced, thus decreasing the computational complexity of the OI algorithm and increasing the computational speed.
[0070] (2) Traditional OI algorithms require temporary storage of the background field error covariance matrix B in the computer during calculation. The size of B is (MN). 2 When M and N are too large, computer memory can hardly support the calculation. This invention introduces the mathematical method of the Kronecker product, which only requires temporary storage of a one-dimensional matrix B in the computer. x With B y Their sizes are M 2 With N 2 This effectively solves the huge demand for computer memory in data fusion.
[0071] (3) This method breaks through the existing technology’s description of the ocean as uniform and isotropic, and supports the calculation using a background error matrix with varying non-uniform anisotropy.
[0072] (4) Based on the multi-scale optimal interpolation algorithm, this invention avoids the point-by-point calculation of the existing optimal interpolation algorithm and can complete the interpolation calculation for the set area in one go.
[0073] (5) The present invention can directly fuse ocean observation data of different resolutions, and successively fuse observation data of different resolutions to gradually reduce observation errors and achieve effective fusion of multi-source ocean observation data. Attached Figure Description
[0074] Figure 1 This is a flowchart of the method of the present invention;
[0075] Figure 2 This is a region grid in one embodiment;
[0076] Figure 3 A visualization of the sea surface temperature background field in one embodiment;
[0077] Figure 4 Sea surface temperature as observed by low-resolution remote sensing in one embodiment;
[0078] Figure 5 This is a schematic diagram of the distribution of observation stations in one embodiment;
[0079] Figure 6 This is a low-resolution remote sensing sea surface temperature fusion result in one embodiment;
[0080] Figure 7 Sea surface temperature as observed by high-resolution remote sensing in one embodiment;
[0081] Figure 8 This is a high-resolution remote sensing sea surface temperature fusion result from one embodiment. Detailed Implementation
[0082] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. These embodiments are based on the technical solution of the present invention and provide detailed implementation methods and specific operating procedures. However, the scope of protection of the present invention is not limited to the following embodiments.
[0083] This invention discloses a multi-source ocean data fusion method based on multi-scale optimal interpolation. Its purpose is to fuse observational data from different observation techniques, with varying error characteristics, temporal discontinuities, and spatial irregularities into a high-resolution data product that is temporally continuous and spatially regularly distributed. This method, based on the fully optimal interpolation mathematical formula, introduces the Kronecker product to construct the background error covariance matrix, enabling simultaneous calculation of optimal interpolation over large areas with massive amounts of observational data. It overcomes the limitations of traditional optimal interpolation algorithms, which require enormous computer memory to calculate the background error covariance matrix, leading to point-by-point calculations, huge computational loads, and the need for oversimplification of the fully optimal interpolation formula. Furthermore, this method incorporates multi-scale methods that can consider different spatial resolutions of ocean satellite observation data.
[0084] First, this embodiment explains the calculation and determination of the traditional optimal interpolation algorithm (OI):
[0085] Traditional OI methods can be summarized as taking the feature value T at any spatial point x as an example. xEstimate based on the observations of a limited number of observation points r (r = 1, 2, 3, ..., M) The actual value; recorded here. The observation error is random noise, and its squared error is E; This is the true value plus random error, so we have Then the background field is introduced. have in This represents the background error; traditional OI provides T. x The optimal estimate is:
[0086]
[0087] in, For the analysis field, It is A rq The inverse matrix, A rq The sum of the background error covariance and the observation error covariance between observation points r and q is calculated using the following formula:
[0088]
[0089] Here, angle brackets represent the statistical average, denoted as... It is the background error at point r. It is the background error at point q, Eδ rq For observation error covariance,
[0090]
[0091] C xr Let x be the background error covariance between grid point x and observation point r.
[0092] The aforementioned traditional OI algorithm is derived based on the theory of optimal interpolation. Its solution has a certain degree of optimality and basically achieves the minimum squared error of the final solution. Currently, the most widely used European satellite altimeter data AVISO and the fused sea surface temperature data both adopt this method.
[0093] However, the aforementioned traditional OI algorithm faces significant difficulties in practical applications:
[0094] ① Formula (1) adopts a point-by-point calculation method, which makes the actual calculation process slow. When the amount of data is huge, this point-by-point calculation method is very time-consuming.
[0095] ② When performing calculations, formula (1) requires searching for all observation data within the search radius point by point. Not only is the entire process computationally intensive, but when the search radius is too small, the continuity between observation points is poor, which can easily lead to spatial discontinuity.
[0096] Theoretically, equation (1) is not completely optimal, but rather a simplification of the completely optimal interpolation formula; the completely optimal interpolation formula can be written as:
[0097] x a =x b +BH T (HBH T +R) -1 (y o -Hx b (4)
[0098] Here x b y is the background field, and B is the background error covariance matrix. o The observation vector consists of all observation points, R is the observation error covariance matrix, H is the observation operator that maps the background field values to the observation points, and x a It is the fusion variable analysis value of all regular grid points.
[0099] In equation (4), the background error covariance matrix B is difficult to calculate, where x a Let B be a vector with n elements, representing the total number of grid points in the fused data. Since B is a matrix with n×n elements, for example, the two-dimensional grid points are only 100×100=10. 4 Therefore, the number of elements in B is as high as 10. 8 With so many elements in B, neither the computational cost nor the memory requirements can be directly used for calculation in equation (4). Therefore, equation (1) simplifies equation (4).
[0100] This invention proposes a novel method for calculating the background error covariance matrix B. In the fully optimal interpolation algorithm, B can be expressed as:
[0101] B=ΣCΣ (5)
[0102] Σ: a diagonal matrix, where the diagonal elements are the root mean square of the background error for each grid point; C: the correlation coefficient matrix.
[0103] In equation (3), the background error covariance matrix B is expressed as uniform and isotropic. However, in reality, the ocean has non-uniform characteristics, and it is difficult to accurately calculate the non-uniform and anisotropic background error covariance matrix using equation (3).
[0104] This method overcomes the difficulty of calculating the background error covariance matrix B in equation (4) by introducing the mathematical method of Kronecker product, thus obtaining a new optimal interpolation algorithm.
[0105] Specifically, such as Figure 1 As shown, the method includes the following steps:
[0106] Step 1) Data Acquisition: Acquire multi-source ocean data, including static data and real-time dynamic data. The static data includes global high-resolution land-sea distribution data and global climate data. The real-time dynamic data includes multi-source observation data at different resolutions.
[0107] Step 1) specifically includes the following steps:
[0108] Step 1.1) Obtain global high-resolution land and sea distribution data: Collect GHRSST 1km data, obtain the mask variable, modify its classification criteria, and reclassify the 5 categories (open ocean, land, lake, sea ice, lake ice) into two categories: ocean and land, to generate global high-resolution land and sea distribution data, where 0 represents land and 1 represents ocean, and the modified data is still stored as a mask variable.
[0109] Step 1.2) Download global climate data: Commonly used global climate data include the WOA series datasets. In this embodiment, the WOA2023 climate (1991-2020) dataset is downloaded to obtain the sea surface temperature data.
[0110] Step 1.3) Download real-time dynamic low-resolution observation data: In this embodiment, download AMSR2 data and obtain the sea surface temperature data therein;
[0111] Step 1.4) Download real-time dynamic high-resolution observation data: In this embodiment, VIIRS near-infrared remote sensing data is downloaded to obtain sea surface temperature data.
[0112] Step 2) Constructing the regional grid: Based on the acquired multi-source ocean data and user requirements, determine the starting and ending latitude and longitude coordinates, set the grid resolution, generate the regional grid, and mark the land and sea distribution of the grid.
[0113] Step 2) specifically includes the following steps:
[0114] Step 2.1) Define the start and end latitude and longitude of the region grid: Set the start latitude and longitude coordinates (lon_start, lat_start) and end latitude and longitude coordinates (lon_end, lat_end) of the region grid, in degrees;
[0115] Step 2.2) Define the grid resolution: Set the resolution (θ) of the region grid in degrees;
[0116] Step 2.3) Generate the regional grid: Using the starting and ending longitude and latitude coordinates of the regional grid defined in Step 2.1) and the grid resolution defined in Step 2.2), generate a longitude vector with N columns and a latitude vector with M rows. Duplicate the longitude vector M times and duplicate the latitude vector N columns to generate the regional grid, as Figure 2 shown;
[0117] Step 2.4) Mark the land-sea distribution (mask) of the regional grid: Use the global high-resolution land-sea distribution data obtained in Step 1.1) to mark the land-sea distribution in the regional grid, where 0 represents land and 1 represents ocean.
[0118] Step 3) Generate the regional climatological background field: Interpolate the global climatological data onto the regional grid to generate the regional climatological background field.
[0119] Step 3) Specifically includes the following steps:
[0120] Step 3.1) Calculate the regional climatological background field: Read the global climatological data and interpolate the ocean elements onto the regional grid generated in Step 2.3) using the bilinear interpolation method. The specific interpolation method is as follows: Assume the target point is to the left of the abscissa x, find two nearest abscissas x1 and x2 respectively, where x1 < x2. According to the relative position of the target point's abscissa and the two nearest neighbor points, use the bilinear interpolation formula to calculate the interpolation result of the target point:
[0121] S = (1 - α) * S(x1, y) + α * S(x2, y) (6)
[0122] where S(x1, y) and S(x2, y) are the values of the two nearest neighbor points respectively, and α is the relative position of the target point. The calculation method is
[0123] Step 3.2) Mark the land-sea distribution of the regional climatological background field: Use the land-sea distribution of the grid marked in Step 2.4) as the mask variable, and multiply the mask variable by the regional climatological background field to make the values of ocean elements at land grid points equal to 0, as shown in the following formula:
[0124]
[0125] is the mask, is the climatological sea background field. Multiplying the two can assign the value of ocean elements at land to 0. The result is as Figure 3 shown.
[0126] Step 4) Process the low-resolution observation data: Calculate the position index of the low-resolution observation data in the regional grid and剔除 the invalid data in this region.
[0127] Step 4) specifically includes the following steps:
[0128] Step 4.1) Eliminate small-scale signals in low-resolution observation data: Calculate the 48-hour average (today and the previous day) of the low-resolution observation data (i.e., AMSR2 sea surface temperature observation data) obtained in step 1.3) to eliminate small-scale signals in the low-resolution observation data.
[0129] Step 4.2) Obtain regional low-resolution observation data: Extract low-resolution observation data located in the regional grid, use xob to represent the longitude coordinates of the regional low-resolution observation data, use yob variable to represent the latitude coordinates of the regional low-resolution observation data, and use ho variable to represent the ocean element values in the regional low-resolution observation data.
[0130] Step 4.3) Remove invalid values from low-resolution observation data in the region: Clean the ho variable to remove outliers and invalid values;
[0131] Step 4.4) Find the nearest grid point of low-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data falling within the regional grid, and use it as the location index.
[0132] In one embodiment, the distance calculation method of the minimum distance method is as follows:
[0133]
[0134] Among them, Dis AB Let be the distance between points A and B, where the coordinates of point A are (x1, y1) and the coordinates of point B are (x2, y2).
[0135] In this embodiment, the processing results of low-resolution remote sensing sea surface temperature data are as follows: Figure 4 As shown.
[0136] Step 5) Calculate the regional low-resolution analysis field: Use a multi-scale optimal interpolation algorithm to fuse the regional climatological background field with low-resolution observation data to generate the regional low-resolution analysis field.
[0137] Step 5) specifically includes the following steps:
[0138] Step 5.1) Set background error: Specify the error value for the climatological background field data at each grid point in the regional grid;
[0139] Step 5.2) Calculate the observation error: The observation error includes measurement error and sampling error. The measurement error of low-resolution observation data is provided by the dataset. The sampling error is 20% of the measurement error. The observation error of each observation point is stored in the observation error matrix. The observation error of the commonly used AMSR2 L2 data is 0.6℃.
[0140] Step 5.3) Calculate the observation increment: Calculate the difference between the regional low-resolution observation data mapped to the regional grid points and the regional climatological background field to obtain the first observation increment. The calculation formula is as follows:
[0141] d low =y o_low_res -H low x b (9)
[0142] Among them, y o_low_res For the region's low-resolution observation data vector, x b H serves as the regional climatological background field. low For the first observation operator, d low This is the first observation increment;
[0143] Step 5.4) Calculate the regional low-resolution analysis field: Using a multi-scale optimal interpolation algorithm, fuse the regional climatological background field with the regional low-resolution observation data to generate the regional low-resolution analysis field x. a_low_res The calculation formula is as follows:
[0144]
[0145] Among them, R low x is the first observation error covariance; a_low_res For the region's low-resolution analysis field; B low The background field error covariance, representing the regional climatological background field, determines the weighting between observational and background field information.
[0146] The background field error covariance is decomposed as follows during calculation:
[0147] B=ΣCΣ (11)
[0148] Where B is the background field error covariance of the regional climatological background field. low Or, the background field error covariance B is used as the background field in a low-resolution analysis field of the region. high Σ is a diagonal matrix, with diagonal elements representing the root mean square error of the background error at each grid point. C is the correlation coefficient matrix, which propagates the observed information to the surrounding blank area and consists of the error correlation coefficients between grid points in the background field. The formula is as follows:
[0149]
[0150] in, C represents the Kronecker product. x Let C be the correlation coefficient matrix in the x-direction. y Let be the correlation coefficient matrix in the y-direction;
[0151] Assuming the spatial grid size is m×n, then C x With C y The calculation matrix forms are as follows:
[0152]
[0153] In this embodiment, the observation station is as follows: Figure 5 As shown, red dots represent observation point locations, and blue dashed lines indicate the grid points closest to the observation data. To find the grid point index (I, J) closest to the q high-resolution observation points, where I = 1:q and J = 1:q, then C... x With C y Further expressed as:
[0154]
[0155] Then calculate for:
[0156]
[0157] Among them, C x With C y Represented using a Gaussian function:
[0158]
[0159] Where deltx is the distance between grid points in the x-direction, and delty is the distance between grid points in the y-direction; L x L is the length scale of the correlation coefficient in the x-direction. y The length scale of the correlation coefficient in the y-direction is expressed in degrees. The larger the length scale of the correlation coefficient, the greater the distance that the observation information travels.
[0160] In this embodiment, the low-resolution sea surface temperature results fused in this step are as follows: Figure 6 As shown.
[0161] Step 6) Process high-resolution observation data: Calculate the location index of the high-resolution observation data in the regional grid and remove invalid data located in this region.
[0162] Step 6) specifically includes the following steps:
[0163] Step 6.1) Obtain regional high-resolution observation data: Extract high-resolution observation data located in the regional grid, use xob′ variable to represent the longitude coordinates of the regional high-resolution observation data, use yob′ variable to represent the latitude coordinates of the regional high-resolution observation data, and use ho′ variable to represent the ocean element values in the regional high-resolution observation data.
[0164] Step 6.2) Remove invalid values from high-resolution regional observation data: Clean the ho′ variable to remove outliers and invalid values;
[0165] Step 6.3) Find the nearest grid point of the high-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data falling within the regional grid, and use it as the location index. In this step, the distance calculation of the minimum distance method refers to Equation (8).
[0166] In this embodiment, the result of processing high-resolution remote sensing sea surface temperature data is as follows: Figure 7 As shown.
[0167] Step 7) Calculate the high-resolution analysis field of the region: Use a multi-scale optimal interpolation algorithm to fuse the low-resolution analysis field of the region with the high-resolution observation data to generate the high-resolution analysis field of the region.
[0168] Step 7) specifically includes the following steps:
[0169] Step 7.1) Set background error: Specify the error value of the low-resolution analysis field data of the region at each grid point in the region grid;
[0170] Step 7.2) Calculate the observation error: The observation error includes the measurement error and the sampling error. The measurement error of high-resolution observation data is provided by the dataset. The sampling error is 20% of the measurement error. The observation error of each observation point is stored in the observation error matrix in the form of elements. The observation error of the commonly used VIIRS dataset is 0.7℃.
[0171] Step 7.3) Calculate the observation increment: Calculate the difference between the high-resolution observation data of the region mapped to the regional grid points and the low-resolution analysis field of the region to obtain the second observation increment. The calculation formula is as follows:
[0172] d high =y o_high_res -H high x a_low_res (twenty one)
[0173] Among them, y o_high_res x is the region's high-resolution observation vector. a_low_res For the region's low-resolution analysis field, H high For the second observation operator, d highThis is the second observation increment;
[0174] Step 7.4) Calculate the multi-scale regional high-resolution analysis field: Using the regional low-resolution analysis field as the background field, and employing a multi-scale optimal interpolation algorithm, fuse the regional high-resolution observation data to generate the multi-scale regional high-resolution analysis field x. a _high_res The calculation formula is as follows:
[0175]
[0176] Among them, R high x is the covariance of the second observation error; a_high_res For high-resolution analysis of multi-scale regions; B high The background field error covariance, which serves as the background field for the low-resolution analysis field in the region, determines the weighting between the observation and background field information. In this embodiment, the background field error covariance B... high The calculation can be referred to formulas (11)-(20), which will not be repeated here.
[0177] In the implementation of the above method, after the multi-scale regional high-resolution analysis field of the day is generated by the cold start method for the first time, the multi-scale regional high-resolution analysis field of the next day is generated by the hot start method using the multi-scale regional high-resolution analysis field of the previous day as the regional background field.
[0178] The preferred embodiments of the present invention have been described in detail above. It should be understood that those skilled in the art can make numerous modifications and variations based on the concept of the present invention without creative effort. Therefore, all technical solutions that can be obtained by those skilled in the art based on the concept of the present invention through logical analysis, reasoning, or limited experimentation on the basis of existing technology should be within the scope of protection defined by the claims.
Claims
1. A multi-source ocean data fusion method based on multi-scale optimal interpolation, characterized in that, Includes the following steps: Step 1) Data Acquisition: Acquire multi-source ocean data, including static data and real-time dynamic data. The static data includes global high-resolution land-sea distribution data and global climate data, while the real-time dynamic data includes multi-source observation data at different resolutions. Step 2) Constructing the regional grid: Based on the acquired multi-source ocean data and user requirements, determine the starting and ending latitude and longitude coordinates, set the grid resolution, generate the regional grid, and mark the land and sea distribution of the grid; Step 3) Generate the regional climatological background field: Interpolate global climatological data onto the regional grid to generate the regional climatological background field; Step 4) Process low-resolution observation data: Calculate the location index of the low-resolution observation data in the regional grid and remove invalid data located in this region; Step 5) Calculate the regional low-resolution analysis field: Use a multi-scale optimal interpolation algorithm to fuse the regional climatological background field and low-resolution observation data to generate the regional low-resolution analysis field; Step 6) Processing high-resolution observation data: Calculate the location index of the high-resolution observation data in the regional grid and remove invalid data located in this region; Step 7) Calculate the high-resolution analysis field of the region: Use a multi-scale optimal interpolation algorithm to fuse the low-resolution analysis field of the region with the high-resolution observation data to generate the high-resolution analysis field of the region; The multi-scale optimal interpolation algorithm is based on the Kronecker product to construct the background field error covariance matrix, thereby achieving optimal interpolation and simultaneous calculation. Step 5) specifically includes the following steps: Step 5.1) Set background error: Specify the error value for the climatological background field data at each grid point in the regional grid; Step 5.2) Calculate the observation error: The observation error includes the measurement error and the sampling error. The measurement error of the low-resolution observation data is provided by the dataset, and the sampling error is a preset percentage of the measurement error. The observation error of each observation point is stored in the observation error matrix. Step 5.3) Calculate the observation increment: Calculate the difference between the regional low-resolution observation data mapped to the regional grid points and the regional climatological background field to obtain the first observation increment. The calculation formula is as follows: in, This is a vector of low-resolution observation data for the region. As the regional climatological background field, As the first observation operator, This is the first observation increment; Step 5.4) Calculate the regional low-resolution analysis field: Using a multi-scale optimal interpolation algorithm, fuse the regional climatological background field with the regional low-resolution observation data to generate the regional low-resolution analysis field. The calculation formula is as follows: in, The first observation error covariance; For low-resolution analysis of the region; The background field error covariance, representing the regional climatological background field, determines the weighting between observational and background field information.
2. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, The real-time dynamic data includes at least low-resolution microwave remote sensing data and high-resolution near-infrared remote sensing data.
3. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, Step 2) specifically includes the following steps: Step 2.1) Define the starting and ending latitude and longitude coordinates of the region grid: Set the starting and ending latitude and longitude coordinates of the region grid in degrees; Step 2.2) Define grid resolution: Set the resolution of the region grid in degrees; Step 2.3) Generate the region grid: Using the start and end latitude and longitude coordinates of the region grid defined in Step 2.1) and the grid resolution defined in Step 2.2), generate longitude and latitude vectors. The longitude vector has... Columns, dimension vectors have Copy the longitude vector. Next, copy the dimension vector. Columns generate a region grid; Step 2.4) Mark the land-sea distribution of the region grid: Use the global high-resolution land-sea distribution data obtained in Step 1) to mark the land-sea distribution in the region grid, where 0 represents land and 1 represents ocean.
4. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, Step 3) specifically includes the following steps: Step 3.1) Calculate the regional climatological background field: Read global climatological data and use bilinear interpolation to interpolate ocean elements onto the regional grid generated in Step 2). The specific interpolation method is as follows: Set the target point to be located on the x-axis. To the left of [the target], find the two nearest neighbors with x-coordinates [their x-coordinates]. and ,in, < Based on the x-coordinate of the target point and the relative positions of its two nearest neighbors, the bilinear interpolation formula is used to calculate the interpolation result of the target point: in, and These are the values of the two nearest neighbors. It is the relative position of the target point, calculated as follows: ; Step 3.2) Mark the land-sea distribution of the regional climatological background field: Use the land-sea distribution of the grid marked in Step 2) as a mask variable, and multiply the mask variable with the regional climatological background field so that the value of the marine element at the land grid point is 0.
5. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, Step 4) specifically includes the following steps: Step 4.1) Eliminate small-scale signals in low-resolution observation data: Calculate the average value of the acquired low-resolution observation data within a preset time interval to eliminate small-scale signals in the low-resolution observation data; Step 4.2) Obtain regional low-resolution observation data: Extract low-resolution observation data located in the regional grid, use the first variable to represent the longitude coordinates of the regional low-resolution observation data, use the second variable to represent the latitude coordinates of the regional low-resolution observation data, and use the third variable to represent the ocean element values in the regional low-resolution observation data; Step 4.3) Remove invalid values from low-resolution observation data in the region: Clean the third variable to remove outliers and invalid values; Step 4.4) Find the nearest grid point of low-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data falling within the regional grid, and use it as the location index.
6. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, Step 6) specifically includes the following steps: Step 6.1) Obtain regional high-resolution observation data: Extract high-resolution observation data located in the regional grid, use the fourth variable to represent the longitude coordinates of the regional high-resolution observation data, use the fifth variable to represent the latitude coordinates of the regional high-resolution observation data, and use the sixth variable to represent the ocean element values in the regional high-resolution observation data. Step 6.2) Remove invalid values from high-resolution regional observation data: Clean the sixth variable to remove outliers and invalid values; Step 6.3) Find the nearest grid point of the high-resolution observation data in the region: Use the minimum distance method to find the index value of the nearest grid point of each observation data that falls within the regional grid, and use it as the location index.
7. The multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1, characterized in that, Step 7) specifically includes the following steps: Step 7.1) Set background error: Specify the error value of the low-resolution analysis field data of the region at each grid point in the region grid; Step 7.2) Calculate the observation error: The observation error includes measurement error and sampling error. The measurement error of high-resolution observation data is provided by the dataset. The sampling error is a preset percentage of the measurement error. The observation error of each observation point is stored in the observation error matrix in the form of elements. Step 7.3) Calculate the observation increment: Calculate the difference between the high-resolution regional observation data mapped to the regional grid points and the low-resolution regional analysis field to obtain the second observation increment. The calculation formula is as follows: in, For regional high-resolution observation vectors, For the low-resolution analysis field of the region. For the second observation operator, This is the second observation increment; Step 7.4) Calculate the multi-scale regional high-resolution analysis field: Using the regional low-resolution analysis field as the background field, and employing a multi-scale optimal interpolation algorithm, fuse the regional high-resolution observation data to generate the multi-scale regional high-resolution analysis field. The calculation formula is as follows: in, The second observation error covariance; For high-resolution analysis of multi-scale regions; The background field error covariance, which serves as the background field for the low-resolution analysis field in the region, determines the weighting between the observation and background field information.
8. A multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 5 or 6, characterized in that, The distance calculation method of the minimum distance method is as follows: in, for The distance between two points The coordinates of the point are ( , ), B The coordinates of the point are ( , ).
9. A multi-source ocean data fusion method based on multi-scale optimal interpolation according to claim 1 or 7, characterized in that, The background field error covariance is decomposed as follows during calculation: in, Background field error covariance of the regional climatological background field Or, the background field error covariance of the low-resolution analysis field in the region as the background field. , It is a diagonal matrix, and the diagonal elements are the root mean square error of the background at each grid point. The correlation coefficient matrix, which propagates observation information to the surrounding blank area, consists of the error correlation coefficients between grid points in the background field. The formula is as follows: in, Indicates the Kronecker product. for The correlation coefficient matrix of directions for The correlation coefficient matrix of directions; Assume the spatial grid size is ,but and The calculation matrix forms are as follows: in, and Represented using a Gaussian function: in, for Distance between directional grid points for Distance between directional grid points; for directional correlation coefficient length scale for The length scale of the directional correlation coefficient is measured in degrees. The larger the length scale of the correlation coefficient, the greater the distance that the observation information travels.