Western Mediterranean Altimetry and Gravity Covariance Study for DOT Estimation Through MIMOS

1. Introduction

Remote sensing is nowadays an essential component of the Earth Observation System. New satellite missions provide a value added in environmental monitoring, as the remote part of an integrated with in-situ ground observations system, like the EU Copernicus (https://www.copernicus.eu/en). Radar satellite altimetry has been used for over 50 years for monitoring the sea level and gravity field on a global or regional scale [1]. The last 20 years satellite altimetry has been transformed into its operational phase [2] with near-real-time products that can be used in, e.g., oceanography [3,4,5], gravity field reconstruction [6,7,8,9], hydrology [10,11,12] and height system unification [13,14,15,16].

The use of satellite altimetry in mean Dynamic Ocean Topography (DOT) or Mean Dynamic Topography (MDT) estimation requires accurate information on the gravity field [17]. The advent of SWOT satellite mission data in combination with the recent satellites of Copernicus program (Sentinel satellites) as well as other satellites (e.g., CryoSat-2) can ensure the adequate accuracy and resolution of the final product. Global as well as regional models of DOT can be computed using various methods, e.g., by merging information from altimeter data, GRACE and GOCE gravity field and oceanographic in situ measurements through a filtering procedure [18], by minimizing a global cost function [19], using an optimal filtering technique using least-squares from a Mean Sea Surface (MSS) and the geoid information from GOCE [20], applying spherical harmonic expansion algorithm in a global product estimation [21]. Input-Output System Theory was introduced in gravity field modelling by Sideris [22], while the extension to the Multiple-Input-Multiple-Output System (MIMOS) theory in DOT estimation was presented in [23]. Since then, many issues in the application of the method were resolved, e.g., the sensitivity of the method in the input error, the kernel representation, the data resolution and the Global Geopotential Model (GGM) expansion [24], the resolution of the final product [25], the noise consideration in the repeated altimetric observations [23] and the incorporation of GOCE gravity field information in the estimation procedure [26]. In addition, the method has been applied in various areas, e.g., the Atlantic Ocean [23], the Eastern Mediterranean Sea [26] and the Central Mediterranean Sea [27]. The final estimated mean DOT and the respective geostrophic circulation components can be compared with purely oceanographic in-situ data [26,27] for validation purposes. The major disadvantage of MIMOS theory application to DOT estimation is the noise of the gravity field observations. The gravity anomaly error Power Spectral Density (PSD) can be computed from simulated noise fields or following other simulation approaches [22,23,24,25,26,27].

In this study, two areas in the Western Mediterranean Sea are chosen to investigate the signal and the error covariance of the altimetric data (CryoSat-2, SARAL-AltiKa, Sentinel-3A and 3B) from 2017 to 2024. The signal and error PSD of each mission is computed taking advantage from the successive information of the altimetry data. A similar computation scheme is tested for the gravity information. Various analytical models are fitted into the data-estimated covariances and the annual mean DOT is determined through MIMOS. In addition, the first year of SWOT data contribution and the estimation of first year SWOT signal and error PSD and covariance are considered. Special attention is given to the investigation of the effect of gravity signal and its error from PSD in MIMOS annual mean DOT approximation using either estimated or random information.

2. Data Used and Preprocessing

2.1. Areas Under Study

Two areas were chosen in the Western Mediterranean Sea. AREA A is the sea part before the Strait of Gibraltar. This area is at the border between the Mediterranean Sea and the Atlantic Ocean and is characterized by complex geological and geodynamical features [28,29,30]. AREA B covers the main part of the Western Mediterranean Sea between the east coast of Spain, the islands of Corsica and Sardinia, the northern coast of Africa and the south shoreline of Italy and France. The areas under study are depicted in Figure 1 (AREA A: -6° ≤ λ ≤ 0° and 34° 30΄ ≤ φ ≤ 37° 30΄; AREA B: 0° ≤ λ ≤ 9° and 36° ≤ φ ≤ 44°).

2.2. Data

Seaborne gravity anomalies collected during the GEOMED 2 project [31,32] were used for the retrieval of the gravity field information in both areas. These data sets were part of Morelli, BGI and SHOM databases which were processed and edited under GEOMED 2 [32]. 48195 point gravity anomalies for AREA A and 260205 values for AREA B were finally chosen. The data distribution for both areas is depicted in Figure 2a – AREA A and in Figure 2b – AREA B. As it is noted in Tziavos [32], seaborne gravity data of GEOMED 2 project have an a-posteriori sd ranging from 2 to 5 mGal according to numerous editing tests (debiasing, crossover adjustments, comparison with altimetry derived gravity, etc.).

Altimetry data sets from SARAL/AltiKa, CryoSat-2, Sentinel-3A, Sentinel-3B and the first year of SWOT satellite were collected from AVISO website (https://www.aviso.altimetry.fr/en/home.html). The L2P version 4 non-time-critical processed data (Sea Surface Heights – SSHs) [33,34] were confined to the areas under study. The SSHs were computed by the addition of the Mean Sea Surface (MSS) to the Sea Level Anomaly (SLA) product delivered by AVISO. Especially for SWOT, the L3_LR_SSH Expert product [35] was utilized because it is consistent with the other satellites. All altimetric data were referenced to the WGS84 ellipsoid and the zero-tide system.

2.3. Preprocessing

An important step before the PSD signal and error estimation is the smoothing of the data. In gravity field modelling this smoothing can be achieved using long wavelength information from a global geopotential model. The removal of these gravity features lead to smoother potential field, ready for appropriate gridding before spectral calculations. The seaborne gravity anomalies were reduced to the XGM2019 [36] reference geopotential model expanded to degree and order 2190. The statistics of the reduced gravity anomalies after the removal of the contribution of the geopotential model are presented in Table 1 for AREA A and in Table 2 for AREA B.

Figure 3 and Figure 4 confirm the significant smoothing of the gravity anomalies in the preprocessing procedure in both studied areas.

3. Results – Discussion

3.1. Empirical and Analytical Gravity Covariance Models

The empirical covariance of the reduced gravity data was computed using EMPCOV routine from the GRAVSOFT software suite [37]. The resolution for the empirical covariance spacing was chosen equal to 3 arcmin, in accordance to the density of the original data. Four analytical models were tested to fit the empirical values, namely the simple exponential model [38]

C(ψ) = a·exp(-bψ),

(1)

where C(ψ) is the empirical covariance value at ψ spherical distance, and a, b are the unknown coefficients of the analytical model; the planar Gaussian model [39]

C(s) = a·exp(-s²/b²),

(2)

where C(s) is the empirical covariance value at s horizontal distance, and a, b are the unknown coefficients of the analytical model; the 3rd order Auto-Regressive Moving Average (ARMA) model [40]

C(ψ) = a·(1+b₁ψ+b₂ψ²+b₃ψ³)·exp(-cψ),

(3)

where C(ψ) is the empirical covariance value at ψ spherical distance, and a, b and c are the unknown coefficients of the analytical model; the Tscherning – Rapp covariance model [41]

$$\mathrm{C}\left(\psi \right)=\mathrm{a}{\sum }_{\mathrm{n}=2}^{\mathrm{M}}{\mathsf{\epsilon }}_{\mathrm{n}}{\displaystyle \frac{{\left(\mathrm{n}-1\right)}^2}{{\mathrm{R}}^2}}{\left({\displaystyle \frac{{\mathrm{R}}^2}{{\mathrm{r}}_{\mathrm{P}}{\mathrm{r}}_{\mathrm{Q}}}}\right)}^{\mathrm{n}+2}{\mathrm{P}}_{\mathrm{n}}\left(\cos \psi \right)+{\sum }_{\mathrm{n}=\mathrm{M}+1}^{\text{∞}}{\mathrm{c}}_{\mathrm{n}}\left(\mathsf{\Delta }\mathrm{g}\right){\left({\displaystyle \frac{R_B^2}{{\mathrm{r}}_{\mathrm{P}}{\mathrm{r}}_{\mathrm{Q}}}}\right)}^{\mathrm{n}+2}{\mathrm{P}}_{\mathrm{n}}\left(\cos \psi \right),$$

(4)

where C(ψ) is the empirical covariance value at ψ spherical distance between the computation point and the data point, a is a scale factor that needs to be estimated, M is the maximum expansion degree of the geopotential model (2190 in our case), ε_n are the error degree variances of XGM2019, R is a mean earth radius, $R_B$ is the radius of Bjerhammar sphere [42], ${\mathrm{r}}_{\mathrm{P}}$ and ${\mathrm{r}}_{\mathrm{Q}}$ are the geocentric distances of points P and Q, P_n are the fully normalized associated Legendre polynomials of degree n and ${\mathrm{c}}_{\mathrm{n}}\left(\mathsf{\Delta }\mathrm{g}\right)$ are the anomaly degree variances from the Tscherning–Rapp model [41]. The empirical covariance values as well as the analytical model estimated in each test area are depicted in Figure 5 and Figure 6.

As shown in Figure 5 and Figure 6, the covariances of the two study areas have a correlation length of approximately 0.1 degree. This fact confirms the homogeneity of the GEOMED 2 database which is used in both areas. In addition, AREA A is smaller than AREA B and this is the reason of the different shape in the final empirical covariance figure. Another point of interest is the estimated variance in each area: in AREA A the observed variance was equal to 27.19 mGal², while in AREA B was found to be 19.45 mGal².

In order to validate the analytical model fitting procedure, quality measures were used based on the variance of the adjustment (a-posteriori variance), the correlation coefficient [43]

$$corrcoef={\displaystyle \frac{C_{xy}}{\sqrt{C_{xx}{C}_{yy}}}}={\displaystyle \frac{{\sum }_{i=1}^n\left(C_{model}\left(i\right)-{\overline{C}}_{model}\right)\left(C_{empirical}\left(i\right)-{\overline{C}}_{empirical}\right)}{\sqrt{{\sum }_{i=1}^n{\left(C_{model}\left(i\right)-{\overline{C}}_{model}\right)}^2}\sqrt{{\sum }_{i=1}^n{\left(C_{empirical}\left(i\right)-{\overline{C}}_{empirical}\right)}^2}}},$$

(5)

and the misfit coefficient [44]

$$misfitcoef={\displaystyle \frac{1}{C_0}}\sqrt{{\displaystyle \frac{1}{n}}{\sum }_{i=1}^n{\left(C_{model}\left(i\right)-C_{empirical}\left(i\right)\right)}^2},$$

(6)

where $n$ is the number of the empirical model covariance values, $C_{model}\left(i\right)$ is the model covariance value at the $i$ point, $C_{empirical}\left(i\right)$ is the empirical covariance value at the same point $i$, ${\overline{C}}_{model}$ is the mean covariance value from the analytical model, ${\overline{C}}_{empirical}$ is the mean covariance value from the empirical function and $C_0$ is the variance of the fitted analytical model. The results of the computed statistical measures are presented in Table 3 for AREA A and in Table 4 for AREA B.

Considering Figure 4 and Figure 5 as well as Table 3 and Table 4, the best performance of the 3rd order ARMA model is confirmed in AREA A. In AREA B, all tested models performed equally with a slightly best fitting of ARMA and Gaussian models.

3.2. Periodogram and Inverse Correlogram Approach for the Gravity Data

The major difficulty in the application of the MIMOS theory in gravity field modelling is the necessity of the input error PSD of the observations. Since only the variance of the observations is known and not the error itself, the direct computation of the observation error PSD is not straightforward. Given that the observation variances change from point to point, the noise is not part of a stationary stochastic process [22,23,45]. This fact complicates the solution in the frequency domain, where the simple algebraic equations are transformed into integral relations [46] and the simplified spectral character is cancelled. Various solutions were proposed for the observation noise statistical information as the simulation of random noise field [47], stationary simulation models of the input error PSD [48,49] and the PSD estimation from successive measurement information, as in the case of satellite altimetry [23,45,50].

A similar approach was tested for the gravity measurements in both areas under study. Gravity information was divided into ten (10) randomly distributed fields in order to evaluate the PSD estimation from the data (periodogram approach). The random division of the gravity information has of course an inherent simulation assumption. Nevertheless, since the original data records were not available, as well as the measurement procedure and the date of each campaign, the division into ten gravity fields randomly was the only possibility in the periodogram approach. The randomly distributed gravity field selection scheme is depicted in Figure 7.

The randomly divided gravity fields were gridded using adjustable tension continuous curvature splines [51] by surface and blockmean modules of the Generic Mapping Tools version 6 [52]. The grid spacing was chosen equal to 3 arcmin and the tension factor to 0.25. Using the multiple gravity field information, a similar to the successive altimetric data algorithm [23] was followed. Each gravity observation field is divided into the unknown signal and noise information:

$$\begin{array}{c}{\Delta g}_1^{obs}={\Delta g}^{signal}+m_1\\ \begin{array}{c}{\Delta g}_2^{obs}={\Delta g}^{signal}+m_2\\ ⋮\\ {\Delta g}_{k+1}^{obs}={\Delta g}^{signal}+m_{k+1}\end{array}\end{array},$$

(7)

where ${\Delta g}_i^{obs}$ is the gravity anomaly observation at a specific grid point, k+1 is the number of the successive field information (in our case k+1=10), ${\Delta g}^{signal}$ is the unknown gravity anomaly signal at the specific grid point and $m_i$ is the unknown gravity anomaly observation error (i = 1,…,k+1). Transforming the grids into the spectral domain and subtracting each gravity field from the successive one, the gravity anomaly error PSD estimation is derived:

$${\Delta G}_i^{obs}-{\Delta G}_{i+1}^{obs}=M_i-M_{i+1}withi=1,...,k,$$

(8)

where ${\Delta G}_i^{obs}$ and $M_i$ are the spectra of the gravity observations and noise at the i cluster field. The mean gravity anomaly error PSD can be estimated using the spectral properties of the periodogram approach [23,50]:

$$2k{\overline{P}}_{mm}={\sum }_{i=1}^k\left({\Delta G}_i^{obs}-{\Delta G}_{i+1}^{obs}\right){\left({\Delta G}_i^{obs}-{\Delta G}_{i+1}^{obs}\right)}^{*}\Rightarrow {\overline{P}}_{mm}={\displaystyle \frac{1}{2k}}{\sum }_{i=1}^k{P}_{{\Delta g}_i-{\Delta g}_{i+1}},$$

(9)

where ${\overline{P}}_{mm}$ is the mean gravity anomaly error PSD, and * denotes the complex conjugate of the spectrum. In a similar manner, the estimation of the mean gravity signal PSD is computed:

$${\Delta G}_i^{obs}+{\Delta G}_{i+1}^{obs}={2{\Delta G}^{signal}+M}_i+M_{i+1}withi=1,...,k,$$

(10)

where ${\Delta G}^{signal}$ is the spectrum of the unknown gravity anomaly signal. The mean gravity signal PSD can be derived using the complete successive gravity information:

$${\sum }_{i=1}^k{P}_{{\Delta g}_i+{\Delta g}_{i+1}}=4k{\overline{P}}_{\Delta g\Delta g}+2k{\overline{P}}_{mm}\Rightarrow {\overline{P}}_{\Delta g\Delta g}={\displaystyle \frac{1}{4k}}{\sum }_{i=1}^k\left(P_{{\Delta g}_i+{\Delta g}_{i+1}}-P_{{\Delta g}_i-{\Delta g}_{i+1}}\right),$$

(11)

The estimation of the gravity signal and error covariance functions can be derived using the computed PSD information. The inverse correlogram approach (PSD estimation from the covariance using FFT) can be applied to the 2D signal and error gravity PSDs. It is needed to note here that the inverse correlogram approach transforms originally the PSD into the correlation function. Since the gravity anomaly observations were reduced to a global GGM (see Section 2.3), the inverse correlogram approach leads to the estimation of the covariance function [53,54]. The mean gravity signal covariance can be estimated using the inverse FFT operator F^-1:

$${\overline{C}}_{\Delta g\Delta g}={\mathit{F}}^{-1}\left\{{\overline{\mathit{P}}}_{\mathit{\Delta }\mathit{g}\mathit{\Delta }\mathit{g}}\right\},$$

(12)

and the mean gravity error covariance:

$${\overline{C}}_{mm}={\mathit{F}}^{-1}\left\{{\overline{\mathit{P}}}_{\mathit{m}\mathit{m}}\right\}.$$

(13)

The 2D gravity signal and error covariances derived from the inverse correlogram approach were averaged over azimuth in order to estimate the 1D isotropic covariances and test the analytical model fitting procedure. The simple exponential and 3rd order ARMA analytical models from Section 3.2 were used for fitting the empirical data estimated from the periodogram and inverse correlogram approach. The empirical as well as the analytical model for the signal and error gravity anomaly covariances are depicted in Figure 8 for AREA A and Figure 9 for AREA B.

As seen in Figure 8 and Figure 9, the 3rd order ARMA analytical model provided the best fitting results in the signal covariance. In the error covariance fitting, ARMA as well as the simple exponential model gave comparable results. The quality measures for the fitting procedure are presented in Table 5 for AREA A and in Table 6 for AREA B.

3.3. Comments on the Gravity Anomaly Covariance Estimation Schemes

Comparing the empirical values of the covariances derived from the classical gravity anomaly observations cluster procedure using EMPCOV [37] and the data-driven periodogram PSD estimation and inverse correlogram approach various issues can be pointed out. The signal covariance derived from the inverse correlogram approach seems smoother than the classical empirical one derived from the observations. We must note here that the inverse correlogram approach leads to the signal covariance estimation, while the classical estimation computes an observation (signal and error) covariance. Of course, the PSD-derived covariance can be estimated under certain assumptions, e.g., the random observation separation (ignoring the exact campaign date) and the noise stationarity hypothesis. The variance of the classical empirical covariance in AREA A was estimated equal to 27 mGal² while in the inverse correlogram approach was 36 mGal². On the contrary, in AREA B the variances of the classical and the PSD-driven procedure were computed 19 mGal² and 6.2 mGal², respectively. In addition, the correlation length of the PSD-derived covariances is shorter than the classical estimated covariances in both areas and drops from 0.1 degrees to 0.05 degrees. This fact can be attributed to the signal covariance estimation procedure that uses the abovementioned assumptions. The random observations separation and the stationarity assumption tend to a Dirac-type covariance estimation and reduce the correlation length.

3.4. Satellite Altimetry PSDs and Covariance Estimation

Satellite altimetry data from 2017 to 2024 were utilized similarly to the gravity approach for the signal and noise PSD and covariance estimation. The repeated altimetric tracks provide a multiple sample configuration. Using this information, a PSD estimation is possible. Following [50] and [23], the observed corrected Sea Surface Height (SSH) can be derived as

h_i = h + Δζ_i + n_i ,

(14)

where h_i is the observed corrected SSH, i = 1,…,k+1 is the number of repeated tracks, h is the time-invariant part of the altimetric observation equation (marine geoid plus mean DOT), Δζ_i is the time-variant part of the DOT, and n_i is the altimetric observation noise. The environmental, geophysical and orbit corrections described in [34] and [35] were applied in order to estimate the corrected SSH. Assuming that Δζ_i follows a random distribution and no correlation between Δζ_i and n_i exists, then a new random variable e_i can be assigned by

e_i = Δζ_i + n_{i .}

(15)

This new variable incorporates the statistical information of the observation noise, as well as the assumed random time-variant part of the DOT. Of course, as shown in [50] and [23], this assumption is not exactly valid in the case of slow-moving oceanographic features; in this case aliasing effects can be incorporated into the oceanic geoid signal.

Repeated altimetric data from CryoSat-2, SARAL/AltiKa, Sentinel-3A and 3B as well as the first year of SWOT data were used. The estimation of the mean yearly 2D signal and error covariance was based on the procedure presented in [23], i.e., data periodogram and inverse correlogram approach. The final 2D covariances were averaged over azimuth in order to compute the isotropic ones. The same, as in the gravity data processing, analytical models were tested in the altimetry case. The estimated yearly satellite altimetry signal and noise covariances for the year 2017 are depicted in Figure 10. SARAL/AltiKa, CryoSat-2 and Sentinel-3A were active at that time. The SSH signal variance is comparable for all the investigated satellites at the order of 1.8 m². For the error variances, Sentinel-3A data led to 0.33 cm² while SARAL/AltiKa and CryoSat-2 error variance were computed at the 6.5 and 8.2 cm² level, respectively. This difference can be attributed to the better performance of Sentinel-3A instrumentation, something that can be confirmed during the whole study period.

In Figure 11, the yearly signal and error altimetry covariance for the first semester of 2024 are presented. For the first semester of 2024, all altimetry satellites (CryoSat-2, SARAL/AltiKa, Sentinel-3A, Sentinel-3B, SWOT) led to similar SSH signal covariances with a variance of 0.42 m², approximately. A slightly lower variance was estimated for the SWOT data, namely 0.33 cm². Different behavior can be seen in the SSH noise covariance estimation: SARAL/AltiKa and CryoSat-2 presented a variance of 7.5 cm², while the Sentinels led to an order of magnitude less estimation (0.18 cm² for Sentinel-3A and 0.55 cm² for Sentinel-3B). The magnitude of the noise in the Sentinels is consistent over all the period of study (2017 – 2024) and confirms the quality of the altimetric instrumentation and the accurate preprocessing of the data.

Different behavior for the noise is depicted in Figure 11 regarding SWOT altimetric data. The noise variance of the first semester of 2024 for SWOT was estimated to 22 cm², a magnitude larger than the noise variances of SARAL/AltiKa and CryoSat-2. This fact can be attributed to the initial phase of the satellite. In order to confirm this behavior, the whole 1st year of SWOT data signal and noise covariance was analyzed. Cycles 1 to 17 of SWOT altimetry data (July 2023 – June 2024) were processed in order to estimate 2D signal and noise PSD and covariance. In Figure 12 (a) and (b) the signal and noise covariances for the first year of SWOT data in AREA A and AREA B are presented, respectively. The error variance in both areas was computed larger of one order of magnitude than the other satellites confirming the results shown in Figure 11. The different instrumentation, the initial phases of the satellite and the special measurement principle (SAR altimetry [55]) can probably affect the error covariance estimation. This will be the subject of another study, where more than 3 years of SWOT will be available and thus more reliable conclusions could be drawn.

3.5. MIMOS Yearly Mean DOT Estimation

Yearly steady mean DOT was estimated through MIMOS following the algorithm presented by Andritsanos et al. [23]. The altimetry estimated signal and error PSD were used in the procedure. In the case of gravity anomaly data, simulated (random noise of 5 mGal sd) as well as estimated PSD (see Section 3.3) were utilized to examine the effect of the estimated noise to the final yearly mean DOT. The final annual mean DOT for the year 2023 using real estimated gravity PSD is depicted in Figure 13a for AREA A and 13b for AREA B. It is noticed that the values near the coastlines are not accurate due to the lack of enough gravity data and the inherent problems of altimetric data in shallow waters.

The differences in the annual DOT estimation using real estimated gravity PSD versus random gravity noise for the year 2017 and 2024 are depicted in Figure 14 and Figure 15. The differences between the two solutions are minor and the major values are found near the coastlines due to the lack of gravity data. The simulated versus the real data input gravity PSD performed equally in the final mean DOT estimation. The sensitivity of the input data noise to the final MIMOS estimates is studied in [24]. In both areas, the estimated input gravity PSD/covariance led to similar DOT output variance w.r.t. the simulated one.

The estimated 2D and averaged 1D mean DOT output covariances are presented in Figure 16 for 2017 in AREA A and in Figure 17 for 2023 in AREA B.

As shown in Figure 16 the variance of the output error annual mean DOT for the year 2017 in AREA A was estimated 2.6 cm² and 2.8 cm² for simulated input gravity noise and real data estimated PSD, respectively. In addition, considering Figure 17, the output variance was computed 1.7 cm² using simulated gravity noise, while in the real data estimated PSD the value was 1.1 cm². On the other hand, as shown in Figure 16 and Figure 17, the computed output annual mean DOT error covariances showed larger correlation length than the simulated ones. The use of simulated random gravity noise led to a spike-like DOT output error covariance.

4. Summary and Conclusions

Signal and error gravity and altimetry PSD and covariance studies in two areas in the Western Mediterranean Sea were presented in this work. Gravity data from GEOMED-2 project were combined with altimetry data from CryoSat-2, SARAL/AltiKa, Sentinel-3A, Sentinel-3B and the first year of SWOT data through a MIMOS scheme in order to estimate yearly mean DOT.

The gravity anomaly covariance was modelled using the classical empirical estimation and the analytical model fitting. A 3rd order ARMA analytical model yielded the best empirical fitting according to proper statistical measures. In the second test, a successive algorithm for estimating the gravity anomaly signal and error PSD using periodogram approach was used. Due to the lack of gravity campaign and date information, randomly selected fields were constructed and the successive algorithm led to the signal and noise covariances derived directly from the data using inverse correlogram approach. The 3rd order ARMA analytical model led to better fitting than the other models (simple exponential, gaussian, etc.) in the signal covariance modelling, while a simple exponential model seemed to fit better in the noise covariance adjustment. The signal covariance derived from the inverse correlogram approach seems smoother than the classical empirical one derived from the observations. Of course, the PSD-derived covariance can be estimated under certain assumptions, e.g., the random observation separation due to the lack of the campaign date and the noise stationarity hypothesis. In addition, due to the abovementioned assumptions, the correlation length of the PSD-derived covariances is shorter than that of the classical covariances.

Satellite altimetry data from 2017 to 2024 were tested in a similar manner for signal and noise PSD and covariance estimation. Repeated altimetric data from CryoSat-2, SARAL-AltiKa, Sentinel-3A and 3B as well as the first year of SWOT data were used. The estimation of the mean yearly 2D signal and error covariance was based on the data periodogram and inverse correlogram repeated algorithm approach. The signal SSH variance is comparable for all the investigated satellites, while for the error variances the Sentinels performed better than SARAL/AltiKa and CRYOSAT. This fact can be attributed to the better performance of the Sentinels’ instrumentation, something that can be confirmed during the whole study period. An interesting time period is the first semester of 2024, when all tested altimetry satellites were active. The data from the recently launched mission of SWOT led to a slightly lower variance in the SSH signal covariance than the other operated satellites. An interesting part of the study is the estimation of the SWOT SSH noise covariance: the noise variance of the first semester of 2024 was estimated a magnitude larger than the noise variance of SARAL/AltiKa and CryoSat-2. This fact can be attributed to the initial phase of the satellite. In order to confirm this behavior, the whole 1st year of SWOT data signal and noise covariance was analyzed and the abovementioned conclusion on the SWOT data error budget was justified. Of course, more reliable conclusions can be drawn after the 3rd year of complete mission data collection.

Using the statistical information of the gravity and altimetry observations a MIMOS scheme was utilized for estimating yearly mean DOT in the areas under study. In the case of gravity anomaly data, simulated as well as estimated PSD were utilized in order to examine the effect of the estimated noise in the final yearly mean DOT. The differences between the two solutions were minor with the larger values being found near the coastlines due to the lack of gravity data. The simulated versus the real data input gravity PSD performed equally in the final mean DOT estimation. On the other hand, the computed output annual mean DOT error covariances were more sensitive to the use of random or estimated input gravity noise focusing on the different correlation length. The use of simulated random gravity noise led to a spike-like DOT output error covariance.

A final remark goes to SWOT altimetry mission. A thorough study of the newly available SWOT data is needed with special attention to its advanced measurement approach. The nadir along track measurements of the classical altimetric procedure can be combined with the swath-based data towards the signal and error PSD estimation. The generic swath grid can be used in a direct PSD estimation without the smoothing of the gridding procedure of the along track observations.