# Tensor decomposition based multisensor hyperspectral fusion for spatial super resolution

Fatima Zahrae M'Rani, Matthieu Puigt, Cédric Jamet, Gilles Delmaire, Claire Guilloteau, Gilles Roussel, Vincent Vantrepotte

<figure><img src="image_1.png" /></figure>

To cite this version:

Fatima Zahrae M'Rani, Matthieu Puigt, Cédric Jamet, Gilles Delmaire, Claire Guilloteau, et al. Tensor decomposition based multisensor hyperspectral fusion for spatial super resolution. 15th Workshop on Hyperspectral Image and Signal Processing: Evolutions in Remote Sensing (IEEE WHISPERS 2025), GRSS (Geoscience and Remote Sensing Society); IEEE (Institute of Electrical and Electronics Engineers), Nov 2025, Bellaterra (Barcelona), Spain. hal-05328866

HAL Id: hal-05328866

https://hal.science/hal-05328866v1

Submitted on 23 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
TENSOR DECOMPOSITION BASED MULTISENSOR HYPERSPECTRAL FUSION FOR
SPATIAL SUPER RESOLUTION

Fatima Z. M'Rani<fnref n="1" />,², Matthieu Puigt<fnref n="1" />, Cédric Jamet², Gilles Delmaire<fnref n="1" />,
Claire Guilloteau<fnref n="1" />, Gilles Roussel<fnref n="1" />, and Vincent Vantrepotte²

¹ Univ. Littoral Côte d'Opale, LISIC – UR 4491, F-62219 Longuenesse, France

² Univ. Littoral Côte d'Opale, CNRS, Univ. Lille, IRD, LOG – UMR 8187, F-62930 Wimereux, France

ABSTRACT

In this paper, we address the multi-sharpening problem by simultaneously fusing multiple multispectral images with a single hyperspectral image. We propose two variants of a novel strategy named G-STEREO-1 and G-STEREO-2. Both methods extend the tensor-based multi-sharpening framework named STEREO, by leveraging the complementarity of several multispectral sources to enhance fusion quality. G-STEREO-1 and G-STEREO-2 are based on a joint tensor decomposition model that incorporates a generalized Sylvester equation within a Canonical Polyadic (CP) tensor decomposition scheme. Our approaches overcome the limitations of existing joint tensor-based fusion techniques, which are restricted to fusing only a single multispectral image with a hyperspectral one. Experimental results show that both G-STEREO-1 and G-STEREO-2 consistently outperform these existing methods.

***Index Terms***— Multi-sharpening, Tensor decomposition, Generalized Sylvester equation.

# 1. INTRODUCTION

Monitoring our planet is ensured by many satellites equipped with either multispectral or hyperspectral imagers. Both generate data cubes capturing information two spatial and one spectral dimensions. However, such systems face a fundamental trade-off between spatial and spectral resolution. Physically, the signal-to-noise ratio (SNR) of a multispectral or hyperspectral imaging sensor is proportional to the ratio between the sensor area and the number of observed spectral bands. Therefore, increasing the number of spectral bands in a hyperspectral image (HSI) inherently results in a reduction in spatial resolution. Currently, Earth observation relies on multispectral systems, which provide high spatial resolution but limited spectral content, and hyperspectral systems, which offer rich spectral information at the expense of lower spatial resolution. This results in a limitation that is critical

in many applications. In this context, multi-sharpening—also known as image fusion or spectral super-resolution—has emerged as an active research field, with numerous methods proposed in the literature [1–3]. This study focuses on approaches based on joint tensor decomposition applied to multi-sharpening problems [4]. Some existing methods [5, 6] rely on solving Sylvester equations to update factors derived from tensor-structured data. However, these techniques are limited to processing pairs of images with fixed spatial resolutions. The present work aims to extend [5] to jointly process a set of images, including multiple multispectral images (MSIs) and one HSI, potentially acquired at heterogeneous spatial resolutions. This problem is motivated by the fusion of data from Sentinel-2 (S2) and Sentinel-3 (S3) data<fnref n="1" />. Specifically, analyzing coastal and estuarine environments requires a fine spatial scale combined with high spectral richness, which neither sensor alone can provide. S2 provides 13 spectral bands at three spatial resolutions—i.e., 10, 20, and 60 m—while S3 offers 21 usable bands at a fixed spatial resolution of 300 m. Both instruments cover the visible and near-infrared domains, making their data complementary for land and ocean color observations. This spatial heterogeneity adds a significant layer of complexity to the image fusion task and poses a limitation for many methods proposed in the literature, which typically focus on fusing only a pair of images—one multispectral and one hyperspectral—both with fixed spatial resolutions. When the spectral bands of the different sensors are sufficiently aligned, the fusion obtained with state-of-the-art methods tends to yield consistent results. However, if the spectral discrepancy is too large, the fusion process may introduce distortions or artifacts, ultimately degrading the quality of the reconstructed image. In this work, we propose two novel tensor decomposition methods named G-STEREO-1 and G-STEREO-2, specifically

<footnote>This work was partially supported by Région Hauts-de-France and by the graduate school IFSEA that benefits from a France 2030 grant (ANR-21-EXES-0011) operated by the French National Research Agency.</footnote>

<footnote n="1">¹S2 and S3 are part of the European Copernicus program. S2 is dedicated to monitoring land surfaces and coastal areas, equipped with a multispectral imager (MSI). S3 is designed for oceanic and terrestrial observations. Equipped with the OLCI and other specialized sensors, S3 delivers data with higher spectral resolution but lower spatial resolution than S2. For simplicity, images from S3 will be abusively considered as “hyperspectral” in the remainder of this paper.</footnote>
designed to tackle the multiple MSI and HSI fusion problem. We compare their performance against recent tensor-based multi-sharpening methods [5–7]. The remainder of the paper is organized as follows. We introduce in Sect. 2 the considered multi-sharpening problem. Section 3 briefly recalls the principles of the STEREO method before describing our proposed method. We then investigate their performance in Sect. 4 while we conclude in Sect. 5.

## 2. PROBLEM STATEMENT

This work focuses on the joint analysis of a set of three MSIs—denoted $M_1 \in \mathbb{R}^{I_1 \times J_1 \times K_1}$, $M_2 \in \mathbb{R}^{I_2 \times J_2 \times K_2}$, and $M_3 \in \mathbb{R}^{I_3 \times J_3 \times K_3}$—as well as one HSI $H \in \mathbb{R}^{I_h \times J_h \times K_h}$, all acquired nearly simultaneously over the same area of interest. The spatial dimensions of the images satisfy $I_1 \ge I_2 \ge I_3 \ge I_h$, and $J_1 \ge J_2 \ge J_3 \ge J_h$. Moreover, the number of spectral bands follows the condition $\max\{K_1, K_2, K_3\} \le K_h$. Such a configuration is representative of data provided by S2 and S3 missions. In practice, since S2 provides three spatial resolutions (i.e., 10, 20, and 60 m), while S3 imagery is limited to a single spatial resolution of 300 m, this results in the following spatial scaling relations:

$$I_h = \frac{I_3}{5} = \frac{I_2}{15} = \frac{I_1}{30} \quad (1)$$

and

$$J_h = \frac{J_3}{5} = \frac{J_2}{15} = \frac{J_1}{30}, \quad (2)$$

while the numbers of spectral bands are equal to $K_1 = 4$, $K_2 = 6$, $K_3 = 3$, and $K_h = 21$, respectively.

<figure><img src="image_3.png" /><figcaption>Fig. 1: Framework for the fusion of S2 and S3 images.</figcaption></figure>

The main objective of this paper is to generate a super-resolved data cube $S$ defined as

$$S \triangleq \llbracket X, Y, Z \rrbracket, \quad (3)$$

where $X \in \mathbb{R}^{I_1 \times F}$, $Y \in \mathbb{R}^{J_1 \times F}$, and $Z \in \mathbb{R}^{K_h \times F}$ are matrices obtained by a rank-$F$ canonical polyadic decomposition (CPD) of $S$. It results that the tensors representing the above-

defined MSIs and HSI read

$$\mathcal{M}_1 = \llbracket P_1^{(1)} X, P_2^{(1)} Y, P_3^{(1)} Z \rrbracket, \quad (4)$$

$$\mathcal{M}_2 = \llbracket P_1^{(2)} X, P_2^{(2)} Y, P_3^{(2)} Z \rrbracket, \quad (5)$$

$$\mathcal{M}_3 = \llbracket P_1^{(3)} X, P_2^{(3)} Y, P_3^{(3)} Z \rrbracket, \quad (6)$$

$$\mathcal{H} = \llbracket P_1^{(h)} X, P_2^{(h)} Y, Z \rrbracket, \quad (7)$$

where $P_1^{(k)}, P_2^{(k)}$ (with $k = h, 1, 2$, or 3) are spatial degradation operators while $P_3^{(k)}$ (with $k = 1...3$) denote spectral degradation ones [5]. In this paper, the fusion problem reads

$$\min_{X,Y,Z} L(X, Y, Z), \quad (8)$$

where the cost function is defined as

$$L(X, Y, Z) = \|\mathcal{H} - \llbracket P_1^{(h)} X, P_2^{(h)} Y, Z \rrbracket\|_F^2 + \sum_{k=1}^3 \lambda_k \|\mathcal{M}_k - \llbracket P_1^{(k)} X, P_2^{(k)} Y, P_3^{(k)} Z \rrbracket\|_F^2, \quad (9)$$

where $\|\cdot\|_F$ denote the Frobenius norm and $\lambda_k$ ($k = 1...3$) are weights to balance the influence of each term of the optimization problem, and $P_1^{(1)}$ and $P_2^{(1)}$ are identity operators.

## 3. PROPOSED METHODS

Both proposed approaches to solve Eq. (8) fall under the generalization framework of the STEREO method [5]. The optimisation procedure a cyclic update of the factors $X$, $Y$, and $Z$ as well, i.e.,

$$X \leftarrow \arg \min_{X} \left\| \mathbf{H}^{(1)} - [(Z \odot P_2^{(h)} Y)] (X P_1^{(h)})^{\top} \right\|_F^2 + \sum_{k=1}^{3} \lambda_k \left\| \mathbf{M}_k^{(1)} - [(P_3^{(k)} Z \odot P_2^{(k)} Y)] (X P_1^{(k)})^{\top} \right\|_F^2. \quad (10)$$

$$Y \leftarrow \arg \min_{Y} \left\| \mathbf{H}^{(2)} - [(Z \odot (P_1^{(h)} X)] (Y P_2^{(h)})^{\top} \right\|_F^2 + \sum_{k=1}^{3} \lambda_k \left\| \mathbf{M}_k^{(2)} - [(P_3^{(k)} Z \odot P_1^{(k)} X)] (Y P_2^{(k)})^{\top} \right\|_F^2, \quad (11)$$

and

$$Z \leftarrow \arg \min_{Z} \left\| \mathbf{H}^{(3)} - [(P_2^{(h)} Y \odot P_1^{(h)} X)] Z^{\top} \right\|_F^2 + \sum_{k=1}^{3} \lambda_k \left\| \mathbf{M}_k^{(3)} - [(P_2^{(k)} Y \odot P_1^{(k)} X)] (Z P_3^{(k)})^{\top} \right\|_F^2, \quad (12)$$

where $\mathbf{M}_k^{(j)}$ and $\mathbf{H}^{(j)}$ ($k = 1, ..., 3$, $j = 1, ..., 3$) denote the mode-$j$ unfolding of the tensors $\mathbf{M}_k$ and $\mathcal{H}$, respectively. Let us consider the update of the factor $X$. The gradient of
the loss function (9) with respect to X is given by

$$ \nabla_X L = \sum_{i=1}^{4} A_i X B_i - C, \quad (13) $$

where

$$ A_1 \triangleq P_1^{(h)\top} P_1^{(h)}, \quad B_1 \triangleq U_h^{\top} U_h, \quad (14) $$

$$ A_2 \triangleq \lambda_1 I, \quad B_2 \triangleq U_1^{\top} U_1, \quad (15) $$

$$ A_3 \triangleq \lambda_2 P_1^{(2)\top} P_1^{(2)}, \quad B_3 \triangleq U_2^{\top} U_2, \quad (16) $$

$$ A_4 \triangleq \lambda_3 P_1^{(3)\top} P_1^{(3)}, \quad B_4 \triangleq U_3^{\top} U_3, \quad (17) $$

$$ C = P_1^{(h)\top} \mathbf{H}^{(1)} U_h + \sum_{i=1}^{3} \lambda_i P_1^{(i)\top} \mathbf{M}_i^{(1)\top} U_i, \quad (18) $$

with

$$ U_h \triangleq \mathbf{Z} \odot P_2^{(h)} \mathbf{Y}, \quad U_1 \triangleq P_3^{(1)} \mathbf{Z} \odot \mathbf{Y}, \quad (19) $$

$$ U_2 \triangleq P_3^{(2)} \mathbf{Z} \odot P_2^{(2)} \mathbf{Y}, \quad U_3 \triangleq P_3^{(3)} \mathbf{Z} \odot P_2^{(3)} \mathbf{Y}. \quad (20) $$

By setting this gradient to zero, we obtain a generalised Sylvester equation [8], i.e.,

$$ \sum_{i=1}^{4} A_i X B_i = C. \quad (21) $$

Solving this type of equation proves particularly challenging, as—unlike the classical Sylvester equation—no closed-form solution within this matrix expression is available. A typical approach involves vectorizing each term [9] and solving the resulting least-squares problem which reads

$$ \mathrm{vec}(X) = \left( \sum_{i} (B_{i}^{\top} \otimes A_{i}) \right)^{-1} \cdot \mathrm{vec}(C), \quad (22) $$

where $\otimes$ denotes the Kronecker product. This leads to our first method, which we refer to as G-STEREO-1. In practice—and as in [8]—we replace the least-square solution (22) by a gradient descent. The general process of our generalized approach G-STEREO-1 is illustrated in Algorithm 1.

### Algorithm 1 G-STEREO-1

**Input:**

Initial values for X, Y, and Z

Weight parameters $\lambda_i$ ($i = 1, \dots, 3$)

**Output:**

$$ S = [X, Y, Z] $$

1. repeat

2. Update X using vectorized form Eq. (10)

3. Update Y using vectorized form Eq. (11)

4. Update Z using vectorized form Eq. (12)

5. until Convergence

However, this solution relies on explicit vectorization of the data, which entails the construction and manipulation of large matrices resulting from Kronecker products. These operations are particularly computationally expensive and demand substantial memory. In practice, the Kronecker product leads to a rapid expansion in the size of the matrices involved, while the corresponding matrix inversions further increase the computational load. To address these limitations, we propose a second method—named G-STEREO-2—whose general process is illustrated in Algorithm 2. Instead of relying on explicit vectorization, G-STEREO-2 employs a momentum-based gradient descent strategy to iteratively update the factor matrices, e.g., from the gradient expression (13) to update X. It reduces computational complexity and accelerates convergence.

### Algorithm 2 G-STEREO-2

**Input:**

Initial values for X, Y, and Z

Weight parameters $\lambda_i$ ($i = 1, \dots, 3$)

**Output:**

$$ S = [X, Y, Z] $$

1. repeat

2. Update X from Eq. (10)

3. Update Y from Eq. (11)

4. Update Z from Eq. (12)

5. until Convergence

## 4. EXPERIMENTAL VALIDATION

We evaluate the performance of the proposed G-STEREO-1 and G-STEREO-2 methods using both simulated and real remote sensing data acquired by S2 and S3 sensors. For the simulated experiments, we use the Pavia and Indian Pines hyperspectral datasets<fnref n="2" />. Spectral bands corresponding to those of S2 and S3 are selected. The hyperspectral images are spectrally and spatially degraded to replicate the resolution characteristics of the S2 and S3 sensors, using the spatial resolution ratios defined in Eqs. (1) and (2). This controlled setup enables a quantitative evaluation, as the original HSI serves as ground truth. For the real data experiments, we use actual S2 and S3 images acquired over the same region. Since high-resolution hyperspectral ground truth is not available in this case, the evaluation is qualitative only. In addition to our proposed methods, we consider several existing tensor decomposition-based fusion techniques for comparison, i.e., STEREO [5] (using CPD), SCOTT & B-SCOTT [6]—which are based on Tucker decomposition—and SCUBA [7] which combines CPD and Singular Value Decomposition. Note that B-SCOTT and SCUBA do not require modeling the degradation operators. Moreover, in these techniques, the images

<footnote n="2">²https://lesun.weebly.com/hyperspectral-data-set.html</footnote>
are partitioned into 4 × 4 blocks, each block having a size of 375 × 375 pixels. All the tested methods are evaluated for different values of the tensor decomposition rank *F*, in order to analyze their stability and performance as a function of the model complexity. The state-of-the-art methods are run for 10 outer loop iterations on *M*₁ and *H* while G-STEREO-1 and G-STEREO-2 are run for 5 ones. The weights *λ*ₖ in the loss (9) are set to 1 (and similarly for the 2-image methods). The evaluation relies on four standard quality metrics for multi-sharpening [2], i.e., PSNR (Peak SNR), which measures reconstruction fidelity; SAM (Spectral Angle Mapper), which quantifies spectral distortion; ERGAS (Relative Global Error in Synthesis), a global normalized error; and Q2n, a multi-band image quality index based on structural similarity.

<table><caption>**Table 1:** Average performance across tested images.</caption><thead><tr><th rowspan="2">Method</th><th colspan="4">Avg. over ranks <i>F</i></th><th colspan="4">Average at best rank <i>F</i><sup>*</sup></th></tr><tr><th>PSNR</th><th>SAM</th><th>Q2n</th><th>ERGAS</th><th>PSNR</th><th>SAM</th><th>Q2n</th><th>ERGAS</th></tr></thead><tbody><tr><td>STEREO</td><td>27.76</td><td>7.72</td><td>0.51</td><td>0.80</td><td>26.62</td><td>8.40</td><td>0.51</td><td>0.95</td></tr><tr><td>SCUBA</td><td>22.59</td><td>20.83</td><td>0.38</td><td>3.91</td><td>23.87</td><td>21.29</td><td>0.46</td><td>3.94</td></tr><tr><td>SCOTT</td><td>22.88</td><td>17.78</td><td>0.33</td><td>1.61</td><td>24.15</td><td>12.46</td><td>0.44</td><td>1.33</td></tr><tr><td>B-SCOTT</td><td>20.67</td><td>23.89</td><td>0.24</td><td>1.76</td><td>21.26</td><td>24.07</td><td>0.21</td><td>1.74</td></tr><tr><td>G-STEREO-1</td><td><b>30.74</b></td><td><b>6.52</b></td><td>0.56</td><td><b>0.60</b></td><td><b>33.09</b></td><td><b>6.36</b></td><td><b>0.70</b></td><td><b>0.45</b></td></tr><tr><td>G-STEREO-2</td><td>30.72</td><td>6.59</td><td><b>0.63</b></td><td>0.61</td><td><b>33.09</b></td><td>6.38</td><td><b>0.70</b></td><td><b>0.45</b></td></tr></tbody></table>

<figure><img src="image_4.png" /><figcaption>**Fig. 2:** Performance vs. rank (Pavia).</figcaption></figure>

Table 1 summarizes the comparison by reporting the average value of each metric computed over the tested images and the different tested ranks (left) or the optimal rank (right). The results clearly indicate both proposed methods, i.e., G-STEREO-1 and G-STEREO-2, consistently outperform the other approaches across all evaluation criteria. Furthermore, Figures 2 and 3 show the performance of each method as a function of *F*, on the Pavia and Indian Pines images, respectively. The proposed G-STEREO-1 and G-STEREO-2 approaches achieve similar performance and consistently outperform all competing methods across all values of *F* and

<figure><img src="image_5.png" /><figcaption>**Fig. 3:** Performance vs. rank (Indian Pines).</figcaption></figure>

evaluation metrics. While SCUBA achieves almost the same performance as the proposed approaches for the Pavia image, it provides the worst one on the Indian Pines dataset. The other methods perform significantly worse than our proposed approaches on both images, highlighting the benefit of jointly processing multiple MSIs with one HSI. Figures 4, 5, and 6 illustrate qualitative results for both the simulated and real data scenarios. The results are shown for the optimal tensor rank of each method: $F^* = 110$ for G-STEREO-1, G-STEREO-2, SCUBA, and B-SCOTT; $F^* = 90$ for STEREO; and $F^* = 40$ for SCOTT. The proposed methods yield more accurate and visually consistent reconstructions, with G-STEREO-2 additionally offering improved computational efficiency (two orders of magnitude faster than G-STEREO-1 on average).

<figure><img src="image_6.png" /><figcaption>**Fig. 4:** Reconstructed images at 779 nm (Pavia).</figcaption></figure>

A similar behavior is observed with real S2 and S3 images (Fig. 6), where G-STEREO-1 and G-STEREO-2 maintain higher image quality. In contrast, most other state-of-
<figure><img src="image_7.png" /><figcaption>Fig. 5: Reconstructed images at 779 nm (Indian Pines).</figcaption></figure>

the-art methods show noticeable visual degradation, except
possibly SCUBA.

<figure><img src="image_8.png" /><figcaption>Fig. 6: Reconstructed images at 761 nm (real S2/S3 data).</figcaption></figure>

5. CONCLUSION

In this work, we proposed a generalized extension of the
STEREO method within a joint tensor decomposition frame-
work to generate super-resolved images combining high
spatial and spectral fidelity. By incorporating generalized
Sylvester equations, our approaches effectively handle ten-
sors of varying dimensions in hyperspectral super-resolution
tasks. Future directions will aim to deeper investigate the
performance of the proposed methods, e.g., by investigating
the effect of the values of the weights λk in Eq. (9) or through
comparisons with in situ data (requiring to handle missing
data due to atmospheric correction [10]). Processing S2/S3
time series will be also explored in the future.

6. REFERENCES

[1] L. Loncan, L. B. De Almeida, J. M. Bioucas-Dias, X. Briottet, J. Chanussot, N. Dobigeon, S. Fabre, W. Liao, G. A. Licciardi, M. Simoes, et al., “Hyperspectral pansharpening: A review,” *IEEE Geosci. Remote Sens. Mag.*, vol. 3, no. 3, pp. 27–46, 2015.

[2] N. Yokoya, C. Grohnfeldt, and J. Chanussot, “Hyperspectral and multispectral data fusion: A comparative review of the recent literature,” *IEEE Geosci. Remote Sens. Mag.*, vol. 5, no. 2, pp. 29–56, 2017.

[3] A. Alboody, M. Puigt, G. Roussel, V. Vantrepotte, C. Jamet, and T. K. Tran, “Experimental comparison of multi-sharpening methods applied to Sentinel-2 MSI and Sentinel-3 OLCI images,” in *Proc. IEEE WHISPERS’21*, 2021.

[4] M. Wang, D. Hong, Z. Han, J. Li, J. Yao, L. Gao, B. Zhang, and J. Chanussot, “Tensor decompositions for hyperspectral data processing in remote sensing: A comprehensive review,” *IEEE Geosci. Remote Sens. Mag.*, vol. 11, no. 1, pp. 26–72, 2023.

[5] C. I. Kanatsoulis, X. Fu, N. D. Sidiropoulos, and W.-K. Ma, “Hyperspectral super-resolution: A coupled tensor factorization approach,” *IEEE Trans. Signal Process.*, vol. 66, no. 24, pp. 6503–6517, 2018.

[6] C. Prévost, K. Usevich, P. Comon, and D. Brie, “Hyperspectral super-resolution with coupled tucker approximation: Recoverability and svd-based algorithms,” *IEEE Trans. Signal Process.*, vol. 68, pp. 931–946, 2020.

[7] C. I. Kanatsoulis, X. Fu, N. D. Sidiropoulos, and W.-K. Ma, “Hyperspectral super-resolution: Combining low rank tensor and matrix structure,” in *Proc. IEEE ICIP’18*. IEEE, 2018, pp. 3318–3322.

[8] F. Ding and T. Chen, “Gradient based iterative algorithms for solving a class of matrix equations,” *IEEE Trans. Autom. Control*, vol. 50, no. 8, pp. 1216–1221, 2005.

[9] N. D. Sidiropoulos, L. De Lathauwer, X. Fu, K. Huang, E. E. Papalexakis, and C. Faloutsos, “Tensor decomposition for signal processing and machine learning,” *IEEE Trans. Signal Process.*, vol. 65, no. 13, pp. 3551–3582, 2017.

[10] F. Steinmetz and D. Ramon, “Sentinel-2 MSI and sentinel-3 OLCI consistent ocean colour products using POLYMER,” in *Proc. SPIE “Remote Sensing of the Open and Coastal Ocean and Inland Waters”*, 2018, vol. 10778.