Ocean gravity inversion method based on low-rank error suppression and frequency-domain regularized spectral gradient

Abstract

The present application relates to the technical field of ocean gravity, and particularly relates to a marine gravity inversion method based on low-rank error suppression and frequency domain regularization spectral gradient, which can effectively improve marine gravity inversion precision, and comprises the following steps: SWOT SSHA data preprocessing; multi-period static residual geoid height field construction; coupling of low-rank error suppression and frequency domain regularization gradient estimation in series; based on the regularization residual gradient field, integral inversion is carried out by using the inverse Vening-Meinesz formula to obtain marine gravity anomaly results; the present application can constrain the propagation of extremely short wave noise in the gradient estimation and gravity anomaly inversion process, reduce the influence of local covariance modeling instability and high frequency noise amplification on the results, thereby improving the stability of the medium and short wave scale gravity anomaly inversion.

Description

Technical Field

[0001] This invention relates to the field of marine gravimetry technology, specifically a marine gravity inversion method based on low-rank error suppression and frequency domain regular spectral gradient that can effectively improve the accuracy of marine gravity inversion. Background Technology

[0002] Ocean gravity anomalies reflect seafloor topography, lithospheric density differences, and seafloor tectonic activity, serving as crucial foundational data for seafloor topography inversion and marine geophysical research. Compared to terrestrial areas, shipborne gravity observations in the ocean are sparsely distributed, unevenly covered, and costly to acquire, making it difficult to meet the needs of large-scale, high-resolution ocean gravity field modeling. Satellite altimetry, by observing sea surface height and its spatial variations, has become an important observational tool for constructing global and regional ocean gravity fields.

[0003] Early studies established a global ocean gravity anomaly inversion framework based on altimetry data from satellites such as Seasonat, Geosat, ERS-1, and TOPEX / POSEIDON, revealing the application value of satellite altimetry gravity anomalies in identifying seafloor topography, seamounts, and fault structures. With the continuous accumulation of altimetry data from multi-source satellites such as CryoSat-2 and Jason-1, the spatial resolution and small-scale structural characterization capabilities of global ocean gravity models have been further improved. Meanwhile, in terms of methodology, existing studies typically use satellite altimetry sea surface height, geoid undulation, or vertical deviation as a basis, combining methods such as least squares configuration, local surface fitting, gradient component estimation, and the inverse Vening Meinesz formula to recover ocean gravity anomalies, and conducting independent verification using shipborne gravity data. However, traditional nadir altimetry mainly provides one-dimensional observations along the orbit, with limited cross-orbit sampling capabilities, and still has certain shortcomings in recovering mid- and short-wave gravity signals and characterizing local details in complex tectonic regions.

[0004] The SWOT (Surface Water and Ocean Topography) mission has provided new observational conditions for high-precision and high-resolution ocean gravity anomaly inversion. However, using real SWOT KaRIn data for ocean gravity anomaly inversion still faces challenges in error control and method adaptation. On the one hand, SWOT wide-swath interferometric altimetry is affected by factors such as attitude, phase, baseline, KaRIn random noise, and calibration residuals. Real SSHA data may contain residual cross-track correlation errors, banded structures, and periodically correlated systematic biases. These errors are amplified during gradient calculation and further propagate to vertical deviation and gravity anomaly estimation results. Therefore, targeted quality control and systematic error correction are needed before geophysical inversion. On the other hand, SWOT data has characteristics such as a wide-swath two-dimensional grid structure, local missing measurements, and spatial correlation errors. Traditional spatial domain gradient estimation and local fitting methods are easily affected by window scale, observation geometry, missing measurement distribution, and noise amplification in practical applications, thus limiting the stability of mid- and short-wave gravity anomaly recovery. Summary of the Invention

[0005] This invention addresses the shortcomings and deficiencies of existing technologies by proposing a marine gravity inversion method based on low-rank error suppression and frequency-domain regularized spectral gradient, which can effectively improve the accuracy of marine gravity inversion.

[0006] This invention achieves its purpose through the following measures: A method for ocean gravity inversion based on low-rank error suppression and frequency-domain regularized spectral gradient, characterized by the following steps: Step 1: SWOT SSHA data preprocessing: Based on gross error removal and low-frequency tilt correction, empirical orthogonal function analysis based on singular value decomposition is used to analyze SVD–EOF, extract and subtract low-rank structured error modes with spatial coherence in multi-period SWOT SSHA, so as to weaken the influence of structured systematic errors. Step 2: Construction of multi-period static residual geoid height field: To address the problem of high-frequency noise being easily amplified in the gradient calculation of SWOT two-dimensional wide-span data, a frequency domain regularized spectral gradient method SRG is constructed. By constraining the high wavenumber gain of the spectral domain differential operator, the amplification of residual high-frequency noise is suppressed and the stability of gradient estimation is improved. Step 3: Couple low-rank error suppression with frequency domain regularized gradient estimation in series: the former weakens the structured error in multi-period SSHA before gradient calculation, while the latter limits the amplification of residual high-frequency error during gradient calculation. The two are connected by the multi-period static residual geoid height field to jointly control the propagation of error from SSHA to residual geoid gradient and ocean gravity anomaly. Step 4: Based on the regularized residual gradient field, integral inversion is performed using the inverse Vening-Meinesz formula to obtain the ocean gravity anomaly results.

[0007] In step one of this invention, let the first... SSHA observations after cross-calibration height correction on a certain band in a given period are: ,in, For observation period index, and Representing the indices along the track and across the track, respectively, a robust baseline along the track is first constructed on each pixel column in the cross-track direction. And define the high-frequency residuals used for the Hampel test: (1), in, For the robust baseline along the track obtained by along-track sliding median filtering, the residuals After removing the low-frequency background, the high-frequency components are then subjected to a local two-dimensional test to determine the consistency between the pixel and its surrounding two-dimensional neighborhood. The local two-dimensional test employs... Pixel neighborhood, threshold set to ,in For robust scaling parameters estimated based on MAD within a local two-dimensional neighborhood; low-pass processing is performed on the residual field only during the tilt parameter estimation stage. Let the SSHA of the i-th period after gross removal be , and firstly, a pixel-wise robust reference field is constructed using multi-period observations: (2), No. The residuals relative to the reference field for each period are: (3), Then, a one-dimensional Gaussian low-pass filter is applied to the residual field along the track direction to obtain the residual field used to estimate the low-frequency tilt error. The low-frequency tilt term is described by a two-dimensional plane within the segmented window along the track. (4), in, and These are the normalized cross-track and along-track coordinates, respectively. , and For the parameters to be estimated, the tilt term is estimated for the left and right halves respectively. The low-frequency tilt term is subtracted from the SSHA after gross error removal to obtain the SSHA after basic correction: (5), in This represents the low-frequency tilt term corresponding to the low-frequency tilt plane obtained by fitting equation (4).

[0008] In step one of this invention, the SVD-EOF method is also used to extract and subtract the low-rank structured error from the multi-period SSHA residuals. Here, EOF represents the empirical orthogonal function, and SVD is a commonly used numerical method for implementing EOF decomposition. Specifically: Let the low-frequency tilt correction be the first The residual field for each period is: (6), in, This is the SSHA after low-frequency tilt correction. This represents the median of the low-frequency tilt correction SSHA over all observation periods. For the structured residual field used in EOF analysis, let the effective spatial pixels participating in the EOF analysis in each period be . The number of observation periods is , will the EOF residual field for each period Flatten into a column vector within the effective cell set: (7), in, Indicates the first pixel in the set of valid pixels. Each two-dimensional spatial cell position, therefore This represents all valid spatial pixels arranged in a fixed order; express 3D real vector space; This represents the column vector obtained by flattening the two-dimensional residual field according to the above effective pixel order, and then constructing the residual matrix: (8), in, They represent the 1st, 2nd, ..., ... The transpose of the EOF residual vector for each observation period, each ; Total number of observation periods; The matrix is ​​a multi-period residual matrix formed by arranging the residual vectors of each period in rows. Each row corresponds to one observation period, and each column corresponds to one effective spatial cell. Perform singular value decomposition: (9), in: (10), in, It is a left singular vector matrix. It is a singular value matrix. Given a right singular vector matrix, in this decomposition, the matrix... column vectors Indicates the first One EOF spatial mode, The The column represents the amplitude variation of the spatial mode over different observation periods, and the singular values. Then characterize the first The contribution of the mode to the total variance, then the The explanatory variance ratio of each modality is defined as: (11), in, This represents the total number of singular values. ; This represents the singular value traversal index during the summation process; Indicates the first A singular value, This represents the variance contribution of the corresponding mode, used to determine the appropriate number of SVD-EOF error modes and analyze the concentration of variance information in the multi-period residual field; for each period, ridge regression is further used to estimate the mode coefficients, assuming the previous... The mode matrix is ​​composed of EOF spatial modes: (12), No. The residual vector for each period is: (13), Then its modal coefficients are: (14), in: For the first EOF mode coefficients for each period, Let be the ridge regression regularization parameter. For the identity matrix, Set to 0.05, and introduce... The aim is to reduce the impact of intermodal correlation and residual noise on coefficient estimation. Based on the estimated modal coefficients, the first... The low-rank structured error vector obtained from the reconstruction over several cycles is: (15), Reconstruct it into a two-dimensional field: (16), in, Represents the low-rank structured error vector The One element, The first one defined by equation (7) The positions of each effective two-dimensional pixel are rearranged according to their positions to obtain the two-dimensional low-rank structured error two-dimensional field. Finally, subtracting this from the SSHA after low-frequency tilt correction yields the SSHA after further correction for system errors: (17), in, Indicates the first System error correction (SSHA) after low-frequency tilt correction and SVD-EOF low-rank structured error subtraction for each observation period; Indicates the first SSHA after low-frequency tilt correction for each observation period; Indicates based on the selected previous The EOF spatial mode reconstruction yields the first A periodic low-rank structured error field.

[0009] Step three of this invention employs the frequency-domain regularized spectral gradient (SRG) method to estimate the residual geoid gradient in the two-dimensional wavenumber domain, thereby reducing the amplification effect of the unregularized spectral domain differential operator on noise at high wavenumbers. Ocean gravity anomaly inversion relies on the gradient information of the residual height field. For the static residual geoid obtained through multi-period fusion... Its eastward and northward gradients can be expressed as: (18), in, and These represent the eastward and northward spatial coordinates, respectively. and Let the eastward and northward gradients of the static residual geoid be represented, respectively; in the frequency domain, let... The two-dimensional Fourier transform is: (19), where, Represents a two-dimensional Fourier transform. and Given the eastward and northward wavenumbers, respectively, the unregularized spectral gradient is: (20), in, The imaginary unit; Before calculating the gradient, a wavenumber-based regularized response function is applied to the residual height field: (twenty one), in: (twenty two), in, It is the two-dimensional radial wavenumber, i.e., the eastward wavenumber. Northbound wavenumber The synthetic wavenumber amplitude, Here, the regularization intensity parameter controls the attenuation intensity at the high wavenumber end. After regularization, the eastward and northward gradient spectra are written as follows: (twenty three), (twenty four), Subsequently, the regularized gradient in the spatial domain is obtained through a two-dimensional inverse Fourier transform: (25), (26), in, This represents the two-dimensional inverse Fourier transform. This represents the residual eastward gradient component in the spatial domain. This indicates the residual northward gradient component in the subsequent spatial domain. and These represent the spectral representations of the eastward and northward gradients in the two-dimensional wavenumber domain after SRG regularization.

[0010] From the perspective of wavenumber response, the gain of the unregularized gradient operator continuously increases with the increase of the two-dimensional radial wavenumber k. However, after introducing a regularized wavenumber response function in the SRG method, the combined gain function of the unregularized spectral domain differential operator and the regularized wavenumber response function is expressed as: (27), When the wave lasts for hours, Approximately Increasing the wavenumber can preserve gradient information of large-scale and mesoscale signals; when the wavenumber is large, Controlled by the regularization term, it no longer increases continuously, thus avoiding excessive amplification of high-frequency noise.

[0011] In step four of this invention, after obtaining a stable residual geoid gradient, it is used as a residual vertical deviation component, and the residual gravity anomaly is inverted using the inverse Vening-Meinesz formula. The calculation points are designed as follows: Contribution points are The residual eastward and northward gradients are respectively and The residual gravity anomaly is then expressed as: (29), in, Normal gravity, for arrive spherical angular distance, It is the azimuth angle. For IVM kernel functions, The area fraction is a unit sphere.

[0012] The invention has been experimentally verified, and the overall RMSE of the SWOT-SRG gravity anomaly is 3.12 mGal, which is about 10.1% and 8.2% lower than that of the SIOV33.1 and DTU21GRA ocean gravity anomaly models, respectively. Spatial distribution and subscale analysis further show that this method can better characterize the local short- to mid-wave gravity anomaly changes under complex topographic backgrounds. Moreover, the accuracy advantage becomes more obvious as relevant wavelength components are gradually included. The results show that the joint processing of low-rank structured error suppression and frequency domain regularized gradient estimation can improve the stability and reliability of real SWOT KaRIn wide-swath data for regional ocean gravity anomaly recovery, and provide a feasible technical approach for the recovery of short- to mid-wave ocean gravity fields in complex tectonic areas. Attached Figure Description

[0013] Figure 1 This is a three-dimensional terrain change map of the study area in this embodiment of the invention.

[0014] Figure 2 These are the transit tracks 312, 006, 284, 562, 256, 534, and 228, which cover the research area in this embodiment of the invention. Figure 2 In the middle (a), it is the SSHA (ssha_karin) of the 10th cycle. Figure 2 In (b), it is the SSHA (ssha_karin) of the 14th cycle. Figure 2 In the middle (c), the SSHA (ssha_karin +height_cor_xover) after cross-calibration height correction is in the 10th cycle. Figure 2 In the middle (d), the SSHA (ssha_karin +height_cor_xover) after cross-calibration height correction is in the 14th cycle. Figure 2 In the middle (e), the uncertainty of the 10th period is (ssh_karin_uncert). Figure 2 In the middle (f), the uncertainty of the 14th period is (ssh_karin_uncert).

[0015] Figure 3 This is a flowchart of the present invention.

[0016] Figure 4 These are the explained variance ratios and cumulative explained variances of the first 10 EOF modes in this embodiment of the invention, wherein... Figure 4 (a) corresponds to pass 562; Figure 4 (b) corresponds to pass 534.

[0017] Figure 5 This refers to the spatial morphology of the modes of the first three EOF structured errors in passes 562 and 534 in the embodiments of the present invention, wherein... Figure 5 In (a), the modal space morphology of the first EOF structured error of pass 562 is shown. Figure 5 (b) shows the modal space morphology of the second EOF structured error in pass562. Figure 5 (c) represents the modal space morphology of the third EOF structured error in pass 562. Figure 5 (d) represents the modal space morphology of the first EOF structured error in pass 534. Figure 5 (e) represents the modal space morphology of the second EOF structured error in pass 534. Figure 5 In (f), the modal space morphology of the third EOF structured error of pass 534 is shown.

[0018] Figure 6 This invention relates to the wavenumber response of the SRG method and its combined gain with the gradient operator in this embodiment: where... Figure 6 (a) corresponds to the regularized magnitude response; Figure 6 In (b), the corresponding gradient operator is multiplied by the regularized result.

[0019] Figure 7 This is a comparison of SSHA error preprocessing before and after steps 562 and 534 in this embodiment of the invention. Figure 7 In the middle (a), the SSHA of pass 562 is not processed. Figure 7 (b) is the SSHA after removing gross errors from pass 562. Figure 7 In the middle (c), the SSHA is the SSHA after systematic error correction for pass 562. Figure 7 In the middle (d), the SSHA is the unprocessed version of pass 534. Figure 7 In the middle (e), the SSHA is the result of removing gross errors from pass 534. Figure 7 In the middle (f), SSHA is the SSHA after system error correction for pass 534.

[0020] Figure 8 is a comparison of the shipborne and altimeter gravity anomaly profiles in an embodiment of the present invention: where Figure 8 In the middle (a), SO23 is represented. Figure 8 (b) is voyage number 85001211.

[0021] Figure 9 This is an evaluation of the spatial plausibility of the SWOT-SRG gravity anomaly results in this embodiment of the invention: Figure 9(a) corresponds to the SIO V33.1 gravity anomaly. Figure 9 (b) shows the SWOT-SRG gravity anomaly. Figure 9 (c) shows the distribution of differences between SIO V33.1 and shipboard gravity data. Figure 9 In the middle (d), the difference distribution between SWOT-SRG and shipboard gravity data is shown. Figure 9 (e) shows the differences between SWOT-SRG and SIOV33.1. Figure 9 The middle (f) corresponds to the seabed topography.

[0022] Figure 10 The following is a comparison of power spectral density in the embodiments of the present invention: (a) voyage SO23; (b) voyage 85001211. Detailed Implementation

[0023] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0024] The study area in this case is located in the central South China Sea—the waters adjacent to the Philippine island arc, with a latitude and longitude range of approximately 113°E–120°E and 10.5°N–14°N. Figure 1 The location of the study area, seabed topography, and distribution of ship-based survey tracks are given, among which... Figure 1 The study area showcases local seafloor topographic features. The western side of the study area borders the South China Sea basin, while the eastern side is adjacent to the Philippine Islands and island arc tectonic zone, resulting in a complex marine-continental transition and tectonic setting. The region exhibits diverse geomorphic units, including continental slopes, basins, seamounts, and island arcs, with significant variations in water depth and intense local seafloor topographic undulations. Complex topography typically corresponds to strong density interface fluctuations and mid-to-shortwave gravity anomalies; therefore, this area is suitable for testing the ability of inversion methods to recover local detail signals.

[0025] The shipborne gravity data used in this example primarily comes from the National Center for Environmental Information (NCEI) Marine Geophysical Database. This data is used to assess the accuracy of the ocean gravity anomaly data in this case. To provide a more reliable independent verification basis, the shipborne survey data underwent standardized formatting and necessary quality control, including organizing track positions, gravity anomalies, and cruise information, removing obvious anomalies, and selecting valid observation points based on the study area. After screening, 18 shipborne cruises were ultimately selected as independent verification data for the accuracy validation discussed later.

[0026] The SWOT KaRIn data used in this example comes from the SWOT Level-2 KaRIn Low Rate Sea Surface Height product (L2_LR_SSH_Basic) released by the AVISO+ / CNES data center. This product provides sea surface height, sea surface height anomalies, uncertainty estimates, quality indicators, and related instrument and environmental correction information from the KaRIn low-rate data stream. This example primarily uses the Sea Surface Height Anomaly (SSHA) data from this product as the inversion input, specifically the PGD0 reprocessed data from Version D. PGD0 stands for reprocessed Version D products, which includes SWOT science cycles 1–31. Following the official SWOT product organization method, L2_LR_SSH data is released in transit orbits as the basic file unit, with each pass corresponding to a half-orbit observation segment, covering the wide swath observation areas on the left and right sides of the KaRIn. This example selects cycle 1–31 data from descending transit orbits of passes 312, 006, 284, 562, 256, 534, and 228, which cover the study area. Due to differences in data availability across different passes during actual observations, with some passes lacking one or two observation cycles, this example uses the actual available cycles for each pass for data selection and fusion during processing. Figure 2 The variables shown are all provided by the L2_LR_SSH_Basic product. Among them, ssha_karin is the KaRIn sea surface height anomaly, and height_cor_xover is the cross-calibration height correction. The sum of the two can be used to obtain the height-corrected SSHA. ssh_karin_uncert is the sea surface height uncertainty estimate provided by the product, which is used to characterize the spatial differences in the quality of KaRIn observations.

[0027] from Figure 2 As can be seen, the SWOT data exhibits a typical two-dimensional wide-span swath structure, with different passes covering the study area in the form of slanted strips. However, significant data quality issues still exist in the original ssha_karin. Figure 2 (a) and Figure 2As shown in (b), there are strip-like discontinuities between different passes, and adjacent stripes are not entirely consistent in terms of amplitude and spatial variation trends; different cycles also exhibit certain spatial structure differences. This indicates that the single-cycle SWOT SSHA not only includes static gravity field-related signals, but is also affected by ocean time-varying signals, pass-level biases, systematic errors, and local outliers. If directly used for gradient calculation, these discontinuities and anomalous signals are easily amplified, thus affecting the gravity anomaly inversion results.

[0028] After introducing cross-calibration height correction, the low-frequency inconsistency between some passes was reduced. Figure 2 (c) and Figure 2 Figure (d) shows that crosspoint correction can improve the overall deviation between some track strips, but it still does not completely eliminate local stripes, half-amplitude differences, and period discontinuities. Therefore, crosspoint correction cannot completely replace subsequent systematic error processing. Figure 2 (e) and Figure 2 As shown in (f), ssh_karin_uncert exhibits significant spatial non-uniformity. High uncertainty is observed at some swath edges, in local anomaly regions, and on individual passes, indicating that the SWOT SSHA error is not uniformly or randomly distributed. To improve the stability and reliability of the inversion results, this example performs gross error removal, systematic error correction, and multi-period robust fusion on the SWOT SSHA in subsequent processing before conducting spectral-domain regularized gravity anomaly inversion.

[0029] The technical workflow in this example mainly includes three stages: SWOT SSHA data preprocessing, multi-period static residual geoid height field construction, and ocean gravity anomaly inversion. Figure 3 As shown.

[0030] Depend on Figure 3 As can be seen, this example first performs gross error removal and low-frequency tilt correction on the SSHA data in the SWOT L2_LR_SSH_Basic product. Then, it employs empirical orthogonal function analysis (SVD-EOF) based on singular value decomposition to extract and subtract low-rank error modes, thereby eliminating local outliers and mitigating the impact of structured systematic errors. Subsequently, a multi-period static residual geoid height field is constructed by combining the mean sea level model, mean dynamic topography, and reference geoid model. Finally, based on the SWOT two-dimensional wide-swath grid data, the residual gradient is estimated using the frequency domain regularized spectral gradient method (SRG), and the ocean gravity anomaly is inverted using the inverse Vening-Meinesz formula.

[0031] SWOT SSHA Data Preprocessing: In the SWOT L2_LR_SSH_Basic product, KaRIn SSHA may exhibit decimeter- to meter-level spikes / jumps in local track segments due to solution failures, geometric degradation, or correction anomalies. To remove these outliers, this example employs a combination of the along-track Hampel test and local 2D anomaly detection for local gross error removal. The Hampel test identifies outliers using the median of the sliding window, while the local 2D anomaly detection estimates the robust dispersion within the neighborhood based on MAD (median absolute deviation) to identify isolated anomalous patches or local abrupt changes in 2D space.

[0032] Let the first SSHA observations after cross-calibration height correction on a certain band in a given period are: ,in and These represent the indices along the orbit and across the orbit, respectively. To avoid interference from large-scale ocean variations on anomaly detection, a robust baseline along the orbit is first constructed on each pixel column in the cross-orbit direction. And define the high-frequency residuals used for the Hampel test: (1), in, Obtained through sliding median filtering along the track. Residual Since the high-frequency components approximately correspond to the low-frequency background after removal, single-point spikes, short-segment jumps, and pulse-type anomalies along the track direction can be identified using the Hampel test. Further local two-dimensional testing is then employed to determine the consistency between the pixel and its surrounding two-dimensional neighborhood. The local two-dimensional test uses... Pixel neighborhood, threshold set to ,in It represents the robust scaling parameter estimated based on MAD within a local two-dimensional neighborhood.

[0033] Even after gross error removal, low-frequency tilt errors may still remain in the SWOT SSHA. To avoid misfitting real mid- and short-wave ocean signals as system tilt, this example only performs low-pass processing along the track on the residual field during the tilt parameter estimation stage. Let the 1st... The periodic SSHA is ,in and These represent the along-track and cross-track direction indices, respectively. First, a pixel-by-pixel robust reference field is constructed using multi-period observations: (2), No. The residuals relative to the reference field for each period are: (3), Then, a one-dimensional Gaussian low-pass filter is applied to the residual field along the track direction to obtain the residual field used to estimate the low-frequency tilt error. In this example, a two-dimensional plane is used to describe the low-frequency tilt term within a segmented window along the track: (4), in, and These are the normalized cross-track and along-track coordinates, respectively. , and These are the parameters to be estimated. Considering that the left and right sides may have different tilt characteristics, in practice, the tilt term can be estimated separately for the left and right halves. Finally, the low-frequency tilt term is subtracted from the SSHA after gross error removal to obtain the SSHA after basic correction: (5), Low-rank structured error subtraction based on SVD-EOF: The low-frequency tilt correction described in this example is mainly used to reduce low-frequency tilt errors in SWOT altimeter data. However, even after the tilt term is removed, some systematic structured errors may still remain in the SSHA. These errors typically have a relatively stable spatial morphology in the pixel coordinate system and exhibit a certain degree of repeatability and spatial coherence in different periods. Therefore, they can be approximated by a few dominant spatial modes, exhibiting certain low-rank characteristics. Based on this, this example uses the SVD-EOF method to extract and subtract low-rank structured errors from the multi-period SSHA residuals. Here, EOF represents empirical orthogonal functions, and its spatial modes can be understood as representative orthogonal spatial bases in the multi-period residual field; SVD is a commonly used numerical method for realizing EOF decomposition. This method characterizes cross-period repeatable structural errors through a few dominant EOF spatial modes and further uses period-by-period regression estimation and subtraction of the corresponding error components. Thus, while reducing the impact of cross-period repeatable structural errors, it can retain the non-repetitive spatial variations in the SSHA related to the medium and short-wave gravity signals.

[0034] Let the low-frequency tilt correction be the first The residual field for each period is: (6), in, This is the SSHA after low-frequency tilt correction. This is the structured residual field used for EOF analysis. Let the effective number of pixels per period be... The number of observation periods is The first EOF residual field for each period Flatten into a vector within the effective cell set: (7), Then construct the residual matrix: (8), in, Each row corresponds to one observation period, and each column corresponds to one effective spatial cell. Next, we will examine the matrix... Perform singular value decomposition: (9), where: (10), in this decomposition, matrix column vectors Indicates the first One EOF spatial mode, The The column represents the amplitude variation of the spatial mode over different observation periods, and the singular values. Then characterize the first The contribution of each mode to the total variance. Let... If is the total number of singular values, then the th The explanatory variance ratio of each modality is defined as: (11), To determine the appropriate number of SVD-EOF error modes and analyze the concentration of variance information in the multi-period residual field, this example uses passes 562 and 534 as examples to statistically analyze the explained variance ratio and cumulative explained variance of the first 10 EOF modes. Figure 4 As shown.

[0035] Depend on Figure 4 It can be seen that the variance information of the residual matrix is ​​mainly concentrated in the first few EOF modes. The variance contribution of the first mode in the two typical passes reaches 57.4% and 64.0% respectively, and the cumulative explained variance of the first three modes reaches 86.9% and 89.8% respectively, indicating that a few low-order modes can characterize the main low-rank structured error components in the residual field. In contrast, the single-mode variance contribution of the fourth and higher-order modes decreases rapidly, and the increase in cumulative explained variance is significantly reduced. The information contained therein is more likely to be mixed with local random noise or mid- and short-wave sea surface changes. Therefore, in this example, the first three EOF modes are uniformly selected for all passes to construct the error basis. This setting can effectively extract the dominant structured error and avoid excessively weakening the effective signal related to the actual mid- and short-wave geophysical changes due to the introduction of too many high-order modes, thus achieving a relatively robust balance between error suppression and signal preservation. In order to further analyze the spatial structure characteristics of these modes, the first three EOF modes of passes 562 and 534 are reconstructed into spatial fields, as shown in the figure. Figure 5 As shown.

[0036] Figure 5This paper presents the spatial morphology of the first three EOF structured error modes obtained from multi-period SSHA residual field decomposition in a typical SWOT KaRIn transit orbit. Each mode exhibits a clear directional structure, span variation, and differences between left and right observation bands within a wide swath observation band, indicating that the residuals after low-frequency tilt correction are not entirely controlled by random noise dispersion, but rather manifest as systematic errors with a repetitive spatial structure. These errors can be approximated by a few low-rank modes, providing a basis for subsequent mode regression subtraction. Furthermore, differences in mode morphology exist between different SWOT transit orbits, indicating that this error has a certain degree of orbital correlation. Therefore, this example constructs error modes separately for each transit orbit to improve the specificity of the systematic error characterization.

[0037] For each cycle, this example further utilizes ridge regression to estimate the modal coefficients. Let's assume... The mode matrix is ​​composed of EOF spatial modes: (12), No. The residual vector for each period is: (13), then its modal coefficients are: (14), in: For the first EOF mode coefficients for each period, Let be the ridge regression regularization parameter. This is the identity matrix. In this example, it will be... Set to 0.05. (Introduce) The aim is to reduce the impact of intermodal correlation and residual noise on coefficient estimation. Based on the estimated modal coefficients, the first... The low-rank structured error vector obtained from the reconstruction over several cycles is: (15) Reconstruct it into a two-dimensional field: (16), Finally, subtracting this from the SSHA after low-frequency tilt correction yields the SSHA after further correction for system errors: (17).

[0038] Construction of a Multi-Period Static Residual Geoid Height Field: Ocean gravity anomalies are mainly controlled by relatively stable seafloor mass distribution and lithospheric structure, while the single-period SWOT SSHA still contains time-varying ocean signals and random observation errors. To reduce the impact of these non-static signals and errors on subsequent gradient estimation, this example constructs a multi-period static residual geoid height field as the input for ocean gravity anomaly inversion, based on the aforementioned error correction. First, the corrected SWOT SSHA is added to the mean sea surface model (MSS) to obtain the SWOT sea surface height field, which is then uniformly interpolated to... A regular grid is then used. Subsequently, the mean dynamic topographic depth (MDT) and reference geoid height are subtracted from the sea surface height field to obtain the residual geoid height field. The reference geoid height is provided by the EGM2008 longwave reference model. This processing reduces the influence of the mean ocean dynamic signal and the longwave background field, ensuring that the residual geoid height field primarily retains the height variations related to medium- and short-wave ocean gravity anomalies. Multi-period averaging helps reduce the influence of random errors and non-static ocean signals, thus providing more reliable input data for subsequent frequency domain gradient calculations.

[0039] Unlike traditional along-track altimeter data, SWOT KaRIn provides two-dimensional wide-swath gridded observations with dense data points and continuous spatial coverage, making it suitable for frequency domain processing on a unified grid. Based on this characteristic, this example employs the frequency-domain regularized spectral gradient (SRG) method to estimate the residual geoid gradient in the two-dimensional wavenumber domain, thereby reducing the noise amplification effect of the unregularized spectral domain differential operator at high wavenumbers. Ocean gravity anomaly inversion relies on the gradient information of the residual height field; for the residual geoid... Its eastward and northward gradients can be expressed as: (18), In the frequency domain, let The two-dimensional Fourier transform is: (19), in, Represents a two-dimensional Fourier transform. and Given the eastward and northward wavenumbers, respectively, the unregularized spectral gradient is: (20), in, It is the imaginary unit.

[0040] As shown in the above equation, the gain of the unregularized spectral domain differential operator is proportional to the wavenumber. Some local interpolation errors, boundary value errors, or random high-frequency noise errors are significantly amplified at high wavenumbers. To suppress the amplification of high-frequency noise in gradient calculation, this example applies a wavenumber-based regularized response function to the residual height field before calculating the gradient: (twenty one), in: (twenty two), Here, is the regularization intensity parameter, used to control the attenuation intensity at high wavenumbers. After regularization, the eastward and northward gradient spectra are written as follows: (twenty three), (twenty four), Subsequently, the regularized gradient in the spatial domain is obtained through a two-dimensional inverse Fourier transform: (25), (26), in, This represents the two-dimensional inverse Fourier transform.

[0041] From the perspective of wavenumber response, the gain of the unregularized gradient operator continuously increases with the increase of wavenumber k, while the combined gain of the SRG method after introducing a regularized wavenumber response function can be expressed as: (27), When the wave number is small Approximately Increasing the wavenumber can preserve gradient information of large-scale and mesoscale signals; when the wavenumber is large, Controlled by the regularization term, it no longer increases continuously, thus preventing excessive amplification of high-frequency noise. Therefore, the SRG method can preserve the target band signal while limiting the contamination of the gradient field by extremely shortwave noise. To intuitively illustrate the constraint effect of the regularization parameter on signals of different wavelengths, this example plots the wavenumber response function under different regularization intensities, as well as the combined gain after combining it with the spectral gradient operator, such as... Figure 6 As shown.

[0042] Figure 6 The wavenumber response of the SRG method and its combined gain when combined with the gradient operator are shown. Figure 6 (a) shows that the regularized response maintains a high throughput within the target band, indicating that the method does not significantly weaken the mid- and short-wave signals associated with ocean gravity anomalies. Figure 6As shown in (b), the gain of the unregularized gradient operator continues to increase at high wavenumbers, while the synthesized gain after regularization is effectively limited. This indicates that the SRG method can achieve a trade-off between preserving mid-to-shortwave structure and suppressing ultra-shortwave noise, making it particularly suitable for wide-swath data like SWOT, which possesses the advantages of two-dimensional sampling but is also susceptible to stripe edge errors. This example uses a regularization parameter... Its characteristic wavelength is approximately 20 km. This setting allows for high throughput of medium and shortwave signals with wavelengths greater than 20 km, while significantly suppressing extremely shortwave components with wavelengths less than 20 km. For actual SWOT wide-swath data, extremely shortwave noise can affect adjacent medium and shortwave scales through spectral leakage and spatial propagation. Therefore, the SRG method improves the stability of ocean gravity anomaly inversion by suppressing the anomalous amplification of extremely shortwave noise in gradient estimation and reducing its propagation impact on the two-dimensional gradient field and gravity anomaly inversion results.

[0043] Traditional least-squares collocation (LSC) typically takes discrete gradient observations as input, selects a local search window around each grid point to be estimated, and constructs covariance matrices between observations and between observations and the estimated quantity based on covariance functions. Its basic form can be expressed as: (28), in, For the gradient components to be estimated, The gradient observation vector within the local window. To observe the covariance matrix, The cross-covariance matrix between the estimator and the observed values ​​is... This is the error variance matrix. In other words, the core of LSC is local statistical estimation based on the covariance model.

[0044] This method has clear statistical significance, but it is susceptible to local data distribution, missing data, covariance model specification, and matrix ill-conditioned problems in wide-span SWOT grid data. The computational cost also increases significantly with the window size and the number of grid points. In contrast, the SRG method does not require constructing a local window for each grid point or solving for local covariance. Instead, it directly performs gradient estimation on a regular grid using a two-dimensional FFT and wavenumber response function. Its computational process is simpler, its frequency response is easier to interpret, and it is more suitable for batch processing of large-scale SWOT grid data.

[0045] Ocean gravity anomaly inversion based on IVM: After obtaining a stable residual geoid gradient, this example uses it as the residual vertical deviation component and employs the inverse Vening-Meinesz formula to invert the residual gravity anomaly. The design calculation points are... Contribution points are The residual eastward and northward gradients are respectively and Then the residual gravity anomaly can be expressed as: (29), in, Normal gravity, for arrive spherical angular distance, It is the azimuth angle. For IVM kernel functions, Let be the area integral of a unit sphere. This formula shows that the residual gravitational anomaly at the calculation point can be obtained by integrating the residual gradient field of the surrounding region.

[0046] Since the residual gradient field in this example is already gridded, the IVM integral can be implemented using FFT to improve computational efficiency. For the singularity of the integrator kernel in the near region, this example uses the interior effect approximation. The interior gravity effect is related to the spatial derivatives of the gradients in the east and north directions, and can be written as: (30), in, and The grid spacing is denoted as . The residual gravity anomaly can be obtained by integrating the outer region FFT and correcting for inner region effects. Finally, the complete ocean gravity anomaly is recovered by adding the residual gravity anomaly obtained from SWOT inversion to the long-wave gravity anomaly provided by the reference model EGM2008.

[0047] Based on the above-described SWOT data preprocessing, static residual height field construction, residual gradient estimation, and IVM inversion process, this example yields two types of SWOT gravity anomaly results. For ease of description, the gravity anomaly obtained by estimating the residual geoid gradient using the SRG method and then performing IVM inversion is referred to as the SWOT-SRG gravity anomaly; the comparative result obtained by estimating the residual gradient using the LSC method and performing IVM inversion within the same framework is referred to as the SWOT-LSC gravity anomaly.

[0048] To examine the effectiveness of the aforementioned gross error removal and systematic error correction methods in processing the true SWOT SSHA, this example selects passes 562 and 534, which have relatively representative error characteristics, and compares the spatial distributions of the original SSHA, the SSHA after gross error removal, and the SSHA after systematic error correction. Figure 7 As shown.

[0049] Depend on Figure 7 (a) and Figure 7As shown in (d), the original SWOT SSHA exhibits a certain degree of spatial non-uniform error characteristics across different passes. Pass 562 shows relatively small overall anomaly amplitude, but patchy anomalies are still visible at the swath edges and in local areas. Pass 534's original SSHA shows more pronounced cross-track tilt and strong local anomalies, with some areas exhibiting large-scale red-blue transitions, indicating that this pass contains strong low-frequency system bias and local outliers. To further quantify these error characteristics, this example presents statistical results before and after gross error removal from three aspects: overall dispersion, robust background fluctuations, and anomaly tails. Specifically, STD reflects the overall dispersion of the data but is easily affected by extreme outliers; Robust σ characterizes the robust fluctuation level of the main data; P99 |r| and Max |r| mainly reflect the anomaly tails and extreme value amplitudes.

[0050] The STD and extreme value indicators in the original data are significantly increased, especially P99 |r| and Max |r| for pass 534, which are much higher than the robust volatility levels. This indicates that the original SSHA not only contains random fluctuations but also a small number of outliers with large amplitudes. This suggests that the error in the original SWOT SSHA is not entirely uniform random noise, but also includes local gross errors, spatial patches, banding structures, and low-frequency system biases.

[0051] This example identifies and removes anomalous observations along the track direction and in a local two-dimensional neighborhood. The processing results are as follows: Figure 7 (b) and Figure 7 As shown in (e), compared with the original SSHA, the gross error removal significantly reduced local mutations, isolated anomalous patches, and extreme outliers, while preserving the main spatial distribution and its mesoscale continuous variations. After gross error removal, STD, P99 |r|, and Max |r| all decreased, indicating a reduction in overall data dispersion and compression of the outlier tails; the extreme outliers in Pass 534 were most significantly weakened. Meanwhile, the Robust σ of both passes remained at approximately 0.01 m before and after processing, indicating that this step mainly affected local outliers and outlier tails, without significantly altering the overall stable background of the SSHA.

[0052] Building upon the gross error elimination, this example further employs a systematic error correction method to suppress structural errors present in the SWOT wide-span observation. Figure 7 (c) and Figure 7 Figure (f) shows the SSHA distribution after systematic error correction. Compared to the result after gross error removal, the discontinuities along and across tracks in the corrected image are further reduced, and the differences between the left and right halves are also decreased. Especially... Figure 7In (f), the significant overall tilt and span gradient differences in the original data are significantly weakened, and the SSHA field changes from a strong unidirectional variation to a more balanced spatial distribution. This indicates that systematic error correction effectively reduces the structured biases related to orbital orientation, observation geometry, and instrument systematic errors in SWOT wide-swath observations, making the processed SSH field closer to the residual height field that can be used for static gravity field inversion.

[0053] This example uses shipboard gravity anomaly as independent verification data, comparing SWOT-SRG gravity anomaly with SIO V33.1, DTU21GRA, and SWOT-LSC gravity anomaly respectively. SIO V33.1 is the free air gravity anomaly product grav_33.1.nc released by SIO. Figure 8 Comparative results of gravity anomaly profiles from two typical expeditions, SO23 and 85001211, are presented.

[0054] The gravity anomaly of expedition SO23 showed significant fluctuations and distinct peak-to-valley transitions, serving as a representative of strongly fluctuating tracks. Expedition 85001211 exhibited relatively smaller overall anomaly amplitudes, but more complex local positive-to-negative anomaly transitions, which can be used to test the response capabilities of different methods to rapidly changing local structures. Overall, the altimeter gravity anomaly results all reflect the main trends of the ship-based survey profiles, but there are differences in the preservation of local peak-to-valley amplitudes. In expedition SO23, the SWOT-SRG curve showed good consistency with the ship-based survey curve in terms of peak amplitude and peak-to-valley positions during the strong positive anomaly near 118.8°E–119.1°E and its rapid descent phase; in expedition 85001211, it also clearly reflected the rapid transitions between positive and negative anomalies near 114.9°E–115.1°E. In contrast, SWOT-LSC exhibits a certain smoothing effect at local peaks and valleys, weakening the anomalous amplitude; while in the extremely shortwave fluctuation range near 115.8°E–116.0°E, the results of several types of altimeter gravity anomalies are generally relatively smooth, indicating that the inversion results still have limited expression of extremely shortwave components due to the influence of observation resolution, spatial regularization and error constraints.

[0055] The RMSE of SWOT-SRG in expedition SO23 was lower than that of SIO V33.1, SWOT-LSC, and DTU21GRA, and also decreased by approximately 12.16%, 19.55%, and 4.13% respectively in expedition 85001211. According to Figure 8The cross-sectional comparison reveals that the advantages of SWOT-SRG are mainly reflected in the stable maintenance of the amplitude of medium- and short-wave gravity anomalies in local peaks and valleys and rapidly changing segments at the medium- and short-wave scales. Furthermore, statistical results from all 18 major expeditions further demonstrate that SWOT-SRG achieved the lowest RMSE in 12 expeditions, with an overall RMSE of 3.12 mGal, lower than SWOT-LSC, SIO V33.1, and DTU21GRA. This indicates that the frequency-domain regularized spectral gradient method can, while constraining high-frequency noise amplification, reduce the attenuation of effective medium- and short-wave signals by traditional local smoothing inversion, enabling SWOT KaRIn wide-swath altimeter data to exhibit better medium- and short-wave structure resolution and result stability in ocean gravity anomaly inversion.

[0056] To evaluate the spatial plausibility of the SWOT-SRG gravity anomaly results, this example presents the SIO V33.1 and SWOT-SRG gravity anomaly fields, the residual distributions of both relative to shipboard gravity data, the differences between the two types of results, and the seafloor topography of the study area, such as... Figure 9 As shown.

[0057] Depend on Figure 9 (a) and Figure 9 As can be seen in (b), the SWOT-SRG gravity anomaly and the SIO model have a high degree of consistency in spatial distribution. Figure 9 (c) and Figure 9 Figure (d) shows the difference distributions of SIO V33.1 gravity anomalies and SWOT-SRG gravity anomalies relative to shipboard gravity data. Comparing the two residual plots reveals that the differences between the SIO V33.1 and SWOT-SRG gravity anomalies and the shipboard data are generally within a small amplitude range, indicating that both can effectively characterize the basic changes in regional gravity anomalies within the study area. However, on local tracks, the SWOT-SRG gravity anomaly exhibits a smaller residual amplitude, especially in the area marked by the red box. Figure 9 In the middle (d) track, the overall color is closer to zero, while Figure 9 In (c), relatively obvious positive and negative deviations still exist at the same location. This indicates that the SWOT-SRG gravity anomaly shows a local improvement over the SIO V33.1 gravity anomaly in some ship survey constraint areas.

[0058] Figure 9 Table (e) presents the differences between SWOT-SRG and SIO V33.1 gravity anomalies. The overall differences are small, but local positive and negative differences are not randomly distributed, but mainly occur near seamount margins, continental slope transition zones, and areas with strong topographic relief. Figure 9As shown in (f), these regions typically have large topographic gradients and complex density interface variations, making them prone to strong medium- and short-wave gravity anomaly responses. Therefore, Figure 9 The local differences in (e) can be understood as the differences in the mid- and short-wavelength structure representations of the two types of gravity anomaly results. Further comparison... Figure 9 (c) and Figure 9 As shown in (d), the ship-based residuals of SWOT-SRG within the red-boxed area are closer to zero, indicating that this method can better preserve the details of gravity anomalies related to local structures in some areas with strong seabed topographic undulations. This phenomenon also reflects the analytical capability of SWOT KaRIn wide-swath altimetry data for medium- and short-wave gravity anomalies caused by seamounts and complex seabed topography, demonstrating its potential for high spatial resolution ocean gravity anomaly inversion. Related studies have also shown that SWOT wide-swath altimetry data can provide more refined spatial observation information in seamount areas, helping to enhance the characterization of local fine-scale sea surface height changes and their related geophysical signals.

[0059] Power spectral density and subscale accuracy analysis: To further analyze the signal representation capabilities of different methods at different wavelength scales, this example selects two typical cruises, SO23 and 85001211, to conduct power spectral density analysis, such as... Figure 10 As shown.

[0060] according to Figure 10 It can be seen that, overall, the PSD of the three gravity anomaly results in both expeditions increases with increasing wavelength, indicating that the energy of the gravity anomaly profile is mainly concentrated in the mid-to-long wavelength range. In the band above 50 km, the PSD curves of SIO V33.1, SWOT-LSC, and SWOT-SRG are generally similar, indicating that the three have good consistency in representing mid-to-long wavelength gravity anomaly signals. This is consistent with... Figure 8 The main peak and valley shapes of the two expeditions correspond, indicating that large-scale fluctuations in the profile are mainly controlled by medium- and long-wave signals. In contrast, the spectral energy differences in the short-wave band are more pronounced, especially in the wavelength range below 20 km. The PSD of SIOV33 is generally higher than that of SWOT-LSC and SWOT-SRG, while SWOT-SRG exhibits lower short-wave spectral energy in expedition 85001211, indicating that frequency domain regularization has a strong constraint effect on high-frequency components. Figure 8 It can be seen that SWOT-SRG can maintain the main structure of the strong positive anomaly and rapid descent phase of voyage SO23, and can also reflect the rapid switching of local positive and negative anomalies in voyage 85001211, but it does not completely follow the extremely short wave oscillations in the ship survey profile.

[0061] Error variations in ocean gravity anomaly data exhibit significant scale differences. A positive improvement rate indicates higher accuracy of SWOT-SRG, while a negative rate indicates a larger RMSE. At the <20 km ultra-shortwave scale, the improvement rate of SWOT-SRG relative to SIO V33.1 is negative, indicating that this method does not have a stable advantage in the ultra-shortwave band. This scale may simultaneously include real shortwave gravity anomalies, ship-based observation noise, and altimeter inversion errors, and the strong regularization constraints may also weaken some effective shortwave signals. As the scale range expands to <50 km, <100 km, and <160 km, the proportion of effective mid- and shortwave signals increases, and the improvement rate of SWOT-SRG relative to both types of reference results turns positive overall, with a more significant improvement compared to SWOT-LSC. This indicates that its accuracy improvement is mainly concentrated in the mid- and shortwave scale range, which is also the main wavelength range affected by local undulations and regional anomaly changes in the profile. In the full-wavelength statistics, SWOT-SRG maintained a stable positive improvement rate, indicating that after integrating signals at different scales, SWOT-SRG achieves more stable accuracy overall, especially compared to SWOT-LSC. This suggests that the spectral domain regularization strategy can reduce the impact of local covariance modeling instability and high-frequency noise amplification on the results by constraining the propagation of extremely shortwave noise in gradient estimation and gravity anomaly inversion, thereby improving the stability of gravity anomaly inversion at medium and short wave scales.

Claims

1. A method for ocean gravity inversion based on low-rank error suppression and frequency-domain canonical spectral gradient, characterized in that, Includes the following steps: Step 1: SWOT SSHA data preprocessing: Based on gross error removal and low-frequency tilt correction, empirical orthogonal function analysis based on singular value decomposition is used to analyze SVD-EOF, extract and subtract low-rank structured error modes with spatial coherence in multi-period SWOT SSHA, so as to weaken the influence of structured systematic errors. Step 2: Construction of multi-period static residual geoid height field: To address the problem of high-frequency noise being easily amplified in the gradient calculation of SWOT two-dimensional wide-span data, a frequency domain regularized spectral gradient method SRG is constructed. By constraining the high wavenumber gain of the spectral domain differential operator, the amplification of residual high-frequency noise is suppressed and the stability of gradient estimation is improved. Step 3: Couple low-rank error suppression with frequency domain regularized gradient estimation in series: the former weakens the structured error in multi-period SSHA before gradient calculation, while the latter limits the amplification of residual high-frequency error during gradient calculation. The two are connected by the multi-period static residual geoid height field to jointly control the propagation of error from SSHA to residual geoid gradient and ocean gravity anomaly. Step 4: Based on the regularized residual gradient field, integral inversion is performed using the inverse Vening-Meinesz formula to obtain the ocean gravity anomaly results.

2. The ocean gravity inversion method based on low-rank error suppression and frequency-domain regularized spectral gradient according to claim 1, characterized in that, In step one, let the first... SSHA observations after cross-calibration height correction on a certain band in a given period are: ,in, For observation period index, and Representing the indices along the track and across the track, respectively, a robust baseline along the track is first constructed on each pixel column in the cross-track direction. And define the high-frequency residuals used for the Hampel test: (1) ， Wherein, is the robust baseline along the track obtained through along-track sliding median filtering, and the residual is... After removing the low-frequency background, the high-frequency components are then subjected to a local two-dimensional test to determine the consistency between the pixel and its surrounding two-dimensional neighborhood. The local two-dimensional test employs... Pixel neighborhood, threshold set to , where is the robust scaling parameter estimated based on MAD in the local two-dimensional neighborhood; low-pass processing along the track is only performed on the residual field during the tilt parameter estimation stage, and is assumed to be the i-th after gross error removal. The periodic SSHA is A pixel-by-pixel robust reference field is constructed using multi-period observations: (2)， No. The residuals relative to the reference field for each period are: (3)， A one-dimensional Gaussian low-pass filter is applied to the residual field along the track direction to obtain the residual field used to estimate the low-frequency tilt error. The low-frequency tilt term is described by a two-dimensional plane within the segmented window along the track. (4)， in, and These are the normalized cross-track and along-track coordinates, respectively. , and For the parameters to be estimated, the tilt term is estimated for the left and right halves respectively. The low-frequency tilt term is subtracted from the SSHA after gross error removal to obtain the SSHA after basic correction: (5)， in This represents the low-frequency tilt term corresponding to the low-frequency tilt plane obtained by fitting equation (4).

3. The ocean gravity inversion method based on low-rank error suppression and frequency-domain regularized spectral gradient according to claim 2, characterized in that, Step one also employs the SVD-EOF method to extract and subtract the low-rank structured error from the multi-period SSHA residuals. Here, EOF represents the empirical orthogonal function, and SVD is a commonly used numerical method for EOF decomposition. Specifically, let the i-th period after low-frequency tilt correction be... The residual field for each period is: (6)， in, This is the SSHA after low-frequency tilt correction. This represents the median of the low-frequency tilt correction SSHA over all observation periods. For the structured residual field used in EOF analysis, let the effective spatial pixels participating in the EOF analysis in each period be . The number of observation periods is , will the EOF residual field for each period Flatten into a column vector within the effective cell set: (7)， Where represents the position of the i-th two-dimensional spatial cell in the effective cell set, therefore This represents all valid spatial pixels arranged in a fixed order; express 3D real vector space; This represents the column vectors obtained by flattening the two-dimensional residual field according to the above effective pixel order, and constructing the residual matrix: (8)， in, They represent the 1st, 2nd, ..., ... The transpose of the EOF residual vector for each observation period, each ; Total number of observation periods; The matrix is ​​a multi-period residual matrix formed by arranging the residual vectors of each period in rows. Each row corresponds to one observation period, and each column corresponds to one effective spatial cell. Perform singular value decomposition: (9)， in: (10), where, It is a left singular vector matrix. It is a singular value matrix. Given a right singular vector matrix, in this decomposition, the matrix... column vectors Indicates the first One EOF spatial mode, The The column represents the amplitude variation of the spatial mode over different observation periods, and the singular values. Then characterize the first The contribution of the mode to the total variance, then the The explanatory variance ratio of each modality is defined as: (11)， in, This represents the total number of singular values. ; This represents the singular value traversal index during the summation process; Indicates the first A singular value, This represents the variance contribution of the corresponding mode; for each period, ridge regression is used to estimate the mode coefficients, assuming the previous... The mode matrix is ​​composed of EOF spatial modes: (12)， No. The residual vector for each period is: (13), Then its modal coefficients are: (14), in: For the first EOF mode coefficients for each period, Let be the ridge regression regularization parameter. For the identity matrix, Set to 0.05, based on the estimated modal coefficients, the first... The low-rank structured error vector obtained from the reconstruction over several cycles is: (15)， Reconstruct it into a two-dimensional field: (16)， in, Represents the low-rank structured error vector The One element, The first one defined by equation (7) The positions of each effective two-dimensional pixel are rearranged according to their positions to obtain the two-dimensional low-rank structured error two-dimensional field. Finally, subtracting this from the SSHA after low-frequency tilt correction yields the SSHA after system error correction: (17), where, Indicates the first System error correction (SSHA) after low-frequency tilt correction and SVD-EOF low-rank structured error subtraction for each observation period; Indicates the first SSHA after low-frequency tilt correction for each observation period; Indicates based on the selected previous The EOF spatial mode reconstruction yields the first A periodic low-rank structured error field.

4. The ocean gravity inversion method based on low-rank error suppression and frequency-domain regularized spectral gradient according to claim 3, characterized in that, Step 3 employs the frequency-domain regularized spectral gradient (SRG) method to estimate the residual geoid gradient in the two-dimensional wavenumber domain, thereby reducing the amplification effect of the unregularized spectral domain differential operator on noise at high wavenumbers. Ocean gravity anomaly inversion relies on the gradient information of the residual height field. For the static residual geoid obtained through multi-period fusion... Its eastward and northward gradients are expressed as follows: (18)， in, and These represent the eastward and northward spatial coordinates, respectively. and Let the eastward and northward gradients of the static residual geoid be represented, respectively. In the frequency domain, let... The two-dimensional Fourier transform is: (19)， in, Represents a two-dimensional Fourier transform. and Given the eastward and northward wavenumbers, respectively, the unregularized spectral gradient is: (20)， in, The imaginary unit; Before calculating the gradient, a wavenumber-based regularized response function is applied to the residual height field: (21)， in: (22)， It is the two-dimensional radial wavenumber, i.e., the eastward wavenumber. Northbound wavenumber The synthetic wavenumber amplitude, Here, the regularization intensity parameter controls the attenuation intensity at the high wavenumber end. After regularization, the eastward and northward gradient spectra are written as follows: (23)， (24)， Subsequently, the regularized gradient in the spatial domain is obtained through a two-dimensional inverse Fourier transform: （25）， (26)， in, This represents the two-dimensional inverse Fourier transform. This represents the residual eastward gradient component in the spatial domain. This indicates the residual northward gradient component in the subsequent spatial domain. and These represent the spectral representations of the eastward and northward gradients in the two-dimensional wavenumber domain after SRG regularization. After introducing a regularized wavenumber response function using the SRG method, the combined gain function of the unregularized spectral domain differential operator and the regularized wavenumber response function is expressed as: (27)， When the wave count is several hours, Approximately Increasing the wavenumber can preserve gradient information of large-scale and mesoscale signals; when the wavenumber is large, Controlled by the regularization term, it no longer increases continuously.

5. The ocean gravity inversion method based on low-rank error suppression and frequency-domain canonical spectral gradient according to claim 3, characterized in that, Step four: After obtaining a stable residual geoid gradient, this is used as the residual vertical deviation component. The inverse Vening-Meinesz formula is then used to invert the residual gravity anomaly. The calculation points are designed as follows: Contribution points are The residual eastward and northward gradients are respectively and The residual gravity anomaly is then expressed as: (29), in, Normal gravity, for arrive spherical angular distance, It is the azimuth angle. For IVM kernel functions, The area fraction is a unit sphere.

Images