\section{Introduction}
\label{sec:intro}


Data fusion is a topic of interest for the remote sensing community
arising in such diverse formulations as cross-sensor feature extraction \cite{Scarpa2018, Errico2014},
multitemporal analysis \cite{Gaetano2014},
multiresolution fusion \cite{Masi2016, Gargiulo2018a, Scarpa2018a}.
Sentinel-2 images provide global acquisitions of multispectral images with a high revisit frequency, supplying data for services such as risk management (floods, forest fires, subsidence, landslide), land monitoring, food security/early warning systems, water management, soil protection and so forth \cite{Drusch2012}.
Due to technological constraints only four out of thirteen bands are provided at the highest resolution of 10 meters. 
The remaining bands are given at 20 or 60 meters. 
Motivated by the above consideration in \cite{Gargiulo2018a} we proposed a CNN-based data fusion technique
for the super-resolution of the 20-m short wave infrared (SWIR) band using the other higher resolution bands to recover the missing spatial details,
following the rationale behind the more familiar pansharpening problem \cite{Vivone2015, Scarpa2018a}. 

The super-resolution of an image without auxiliary bands can be addressed 
in such diverse manners ranging from a polynomial interpolation to advanced 
CNN-based modeling \cite{Dong2016,Kim2016}.
The availability of high resolution references, or guide, like in pansharpening,
allows to obtain highly accurate super-resolution, 
thanks to the fusion between the spectral information conveyed by the target image
and the spatial content gathered by the guiding band/s.
Classical multiresolution fusion methods are often categorized as component substitution (CS) or multiresolution analysis (MRA),
a survey of which can be found in \cite{Vivone2015}.
Several pansharpening methods,
based on both CS \cite{Shah2008,Tu2001}
MRA \cite{Chavez1991,Ranchin2000}, were adapted to the Sentinel-2/SWIR case and compared
in \cite{Du2016, Wang2016}.

Other very recent approaches,
both model-based \cite{Brodu2017, Lanaras2017} and data-driven CNN-based \cite{Lanaras2018}, 
aimed to super-resolve the whole Sentinel-2 dataset testify 
the relevance of the topic to the community.

In this work we propose an improved version of our previous method \cite{Gargiulo2018a} 
that counts four main integrations, (i) the use of the residual learning strategy \cite{He2016},
(ii) the batch normalization \cite{Ioffe2015}, both (i) and (ii) aimed to speed-up the learning, 
(iii) a high-pass preprocessing of the input which has already 
proven to be effective for pansharpening in \cite{Yang2017},
and, finally, (iv) the extension to all 20-m bands.

The rest of the paper is organized as follows. 
In the next section we summarize the proposed solution. 
Section 3 presents numerical and visual results. 
Finally, Section 4 draws concluding remarks.



\section{Proposed Sentinel-2 image super-resolution}

 \begin{figure*}[h]
\centering
	\pgfuseimage{net2}
	\caption{\small Top-level workflow of the proposed super-resolution method for 20-m bands of Sentinel-2. 
	Only dashed boxes are used in \cite{Gargiulo2018a}.}
	\label{fig:flowchart}
\end{figure*}
The proposed solution is summarized in the top-level flowchart of Fig.\ref{fig:flowchart}. 
At the core is a three-layer CNN block whose internal architecture is nearly the same 
of \cite{Gargiulo2018a}, the main difference being the output shape 
which counts here six bands rather than just the single SWIR.
For the sake of brevity,
we skip the description of the internal setting of the CNN
as the Reader can refer to \cite{Gargiulo2018a} to this purpose, 
focusing on the remaining architectural aspects.

First, 
it can be noticed the additional batch normalization (BN) layer which processes the input stack $(\mathbf{z}_{\rm H},\widetilde{\mathbf{x}}_{\rm H})$
that feeds the network in order to make the learning process robust with respect to the 
statistical fluctuations of the training dataset \cite{Ioffe2015}.

The BN block takes a concatenation stack as input which comprises both the high pass filtered (HPF) 10-m resolution bands of Sentinel-2, $\mathbf{z}_{\rm H}$,
and the upscaled (via bicubic interpolation) HPF version of the target 20-m resolution bands, $\widetilde{\mathbf{x}}_{\rm H}$.
The use of HPF bands in place of the original ones has been introduced in the contex of pansharpening in \cite{Yang2017}
on the basis on the intuition that low spatial frequencies do not carry relevant information about spatial details,
therefore they can be neglected eventually facilitating the training process.

Another relevant trick generally useful when dealing with deep convolutional networks,
and actually effective for super-resolution and pansharpening as well \cite{Scarpa2018a},
is the use of residual learning \cite{He2016}, which is here implemented 
through the skip connection (on top) that directly links the input (upscaled) to the output ($\hat{\mathbf{x}} = \widetilde{\mathbf{x}} + \hat{\mathbf{y}}$).
Intuitively speaking,
since the low-frequency content of the desired output is already comprised in the low resolution
input bands, it is sufficient for the network to learn only how to predict the complementary
detail component $\hat{\mathbf{y}}$ to be combined with the input $\widetilde{\mathbf{x}}$ in order to get the full-resolution product.
The practical consequence of this is that the network learns much faster, 
which is very useful to finetune the network depending on the dataset.

\subsection{Training setting}

In order to train the network's parameters
a sufficiently large number of input-output examples
and the choice of a suitable cost function to minimize on them
are required to run any learning algorithm, 
like for example the Stochastic Gradient Descent
as done in \cite{Gargiulo2018a, Scarpa2018a} and here, as well,
using the $L_1$-norm as cost.
Due to the lack of ideally super-resolved samples to be used to this purpose,
it is becoming a common practice to resort to a self-supervised learning strategy
that sounds like the so-called Wald's protocol, 
a procedure commonly used to compare objectively different pansharpening algorithms
using referenced data \cite{Vivone2015}.
Given a training Sentinel-2 sample $(\mathbf{x},\mathbf{z})$, 
distinguishing between low and high resolution bands, respectively,
a training example is generated as follows.
Both bands subsets undergo a downsampling whose band-wise anti-aliasing filters
mimic the corresponding sensor characteristics (cutoff frequencies)
and the resolution downgrade factor is equal to the scale factor between $\mathbf{z}$ and $\mathbf{x}$,
which is $2\times2$ in our case.
Said $(\mathbf{z}_\downarrow, \mathbf{x}_\downarrow)$ the reduced-resolution sample so obtained, 
the original band set $\mathbf{x}$ can now play the role of ideal output (reference or label)
corresponding to  $(\mathbf{z}_\downarrow, \mathbf{x}_\downarrow)$.
This process is repeated on the whole training image collection whose size 
must be consistent with the number of parameters to train.
The larger the network the larger must be the training dataset to avoid overfitting problems.
Additional details on training can be found in \cite{Gargiulo2018a} as this part has not been modified.

\section{Experimental Results}


\begin{table}
\centering
\def0.8{0.8}
\setlength{\tabcolsep}{3pt}
\footnotesize
\begin{tabular}{l|ccccc}  
\hline
\rule{0mm}{3mm} \textbf{Method} & \textbf{Q-index} & \textbf{ERGAS} & \textbf{HCC} & \textbf{Epochs}   & \textbf{Sec./Ep.}\\

\rule{0mm}{3mm} & (1)  & (0)   & (1) & & \\

\hline
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a} (early stop)  & $ 0.9853 $ & $ 4.950 $ &  $ 0.5686$ & 200  & 3.577 \\
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a}  & $ 0.9908 $ & $ 3.619 $ &  $ 0.8482$ & 1500 & 3.577 \\
\hline 

\rule{0mm}{3mm} Proposed (without $\mathbf{z}$)    & $ 0.9905 $ & $ 3.985 $ &  $ 0.5792$ & 200  & 1.918\\
\rule{0mm}{3mm} Proposed (no HPF)  & $ 0.9848  $ & $ 3.112 $ &  $ 0.8482$ & 200 & 3.613\\
\rule{0mm}{3mm} Proposed & $\mathbf{0.9978 }$ & $\mathbf{1.937}$ & $\mathbf{0.8819} $ & 200 & 3.597\\

\hline
\end{tabular}
\caption{Ablation study}
\label{tab:ablation}
\end{table}

\begin{table}
\centering
\def0.8{0.8}
\setlength{\tabcolsep}{3pt}
\footnotesize
\begin{tabular}{l|ccc}  
\hline
\rule{0mm}{3mm} \textbf{Method} & \textbf{Q-index} & \textbf{ERGAS} & \textbf{HCC}  \\
\rule{0mm}{3mm} & (1)  & (0)   & (1)  \\
\hline 
\rule{0mm}{3mm} bicubic (without $\mathbf{z}$)  & $ 0.9825  $ & $ 5.340 $ & $0.4444 $  \\
\rule{0mm}{3mm} PCA  & $ 0.9596 $ & $ 8.331 $ & $0.7477 $  \\
\rule{0mm}{3mm} IHS  & $0.9434 $ & $ 9.544 $ & $0.7431 $  \\
\rule{0mm}{3mm} HPF  \cite{Chavez1991} & $ 0.9909 $ & $ 3.891 $ & $0.7091 $  \\
\rule{0mm}{3mm} GS2-GLP  & $ 0.9919  $ & $ 3.672 $ & $0.7414 $  \\
\rule{0mm}{3mm} Indusion  & $ 0.9877  $ & $ 4.494 $ & $0.6844 $  \\
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a}  & $ 0.9908 $ & $ 3.619 $ &  $ 0.8482$\\
\rule{0mm}{3mm} Proposed & $\mathbf{0.9978 }$ & $\mathbf{1.937}$ & $\mathbf{0.8819} $\\
\hline

\end{tabular}
\caption{Comparative analysis.}
\label{tab:comparison}
\end{table}


In order to assess the performance of the proposed method a separate test image,
a 2700$\times$2700 scene (about Rome city center),
not used in the training phase,
has been considered.
The numerical accuracy was assessed by means of commonly used indicators like the
quality (Q-)index and ERGAS (see \cite{Vivone2015} for a detailed description of these measures),
in addition to HCC \cite{Gargiulo2018a}, a correlation coefficient computed over detail components.
Comparative methods are the bicubic interpolator (it does not take any auxiliary high-resolution guide $\mathbf{z}$ as input),
just as baseline, the model M5 proposed in \cite{Gargiulo2018a} properly generalized to output all six bands $\hat{\mathbf{x}}$,
plus other classical pansharpening methods (which can be found in \cite{Vivone2015}) suitably adapted to the case of Sentinel-2.

In Tab.\ref{tab:ablation} it is summarized an ablation study restricted to the proposed method and \cite{Gargiulo2018a}.
A first observation is about the training cost in terms of computational time (last two columns). 
The tradeoff between training time and accuracy is highlighted for M5 that requires 
1500 epochs to reach an accuracy closer to that obtained by the proposed model trained with only 200 epochs.
Focusing on the proposal
it can be noticed, as expected, that by using additional guiding high-resolution bands $\mathbf{z}$ is certainly beneficial.
Moreover, the introduction of the high pass filtering provides further improvements, notably in terms of ERGAS.


Tab.\ref{tab:comparison} gathers the numerical comparison between the proposed solution and the reference methods.
The proposal clearly outperforms all comparative solutions with respect to all quality indicators considered here.


\begin{figure}
\centering
\begin{tabular}{cc}
\pgfuseimage{Bx} & \pgfuseimage{Bz} \\
$\widetilde{\mathbf{x}}_\downarrow$ & $\mathbf{z}_\downarrow$ \\[2mm]
\pgfuseimage{Br} & \pgfuseimage{Bhpf} \\
$\mathbf{x}$ (ground-truth) & HPF \cite{Chavez1991}\\[2mm]
\pgfuseimage{Bfive} & \pgfuseimage{Bprop} \\
Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Super-resolution results on a 130$\times$130 tile of the test image (Rome), in the reduced-resolution domain.}
\label{fig:BIG}
\end{figure}

\begin{figure}
\centering
\begin{tabular}{ccc}
\pgfuseimage{Dx} & \pgfuseimage{Dz} & \pgfuseimage{Dr} \\
$\widetilde{\mathbf{x}}_\downarrow$ & $\mathbf{z}_\downarrow$ & $\mathbf{x}$ (ground-truth)\\[2mm]
\pgfuseimage{Dhpf} & \pgfuseimage{Dfive} & \pgfuseimage{Dprop} \\
 HPF \cite{Chavez1991} & Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Zoomed detail in the reduced-resolution domain.}
\label{fig:DET}
\end{figure}

To complete the analysis of the results let us look at some sample result.
For the sake of brevity we only show a limited portion of the Rome image,
selecting a meaningful subset of comparative solutions.
A 130$\times$130 cropped sample is shown in Fig.\ref{fig:BIG}.
Here the reduced resolution behaviour is analysed so that we can benefit 
from the availability of a ground-truth. 
For the purpose of visualization only three out of 6 (4) bands of $\widetilde{\mathbf{x}}_\downarrow$ ($\mathbf{z}_\downarrow$) are shown 
in a false color RGB representation. 
Likewise, truth $\mathbf{x}$ and compared results are shown according with the same restriction.
On the top row is the input, split in low ($\widetilde{\mathbf{x}}_\downarrow$) and high ($\mathbf{z}_\downarrow$) resolution components.
In the middle row are the ideal super-resolution and HPF \cite{Chavez1991}.
On the bottom are gathered M5 \cite{Gargiulo2018a} and the proposal.
From visual inspection it can be appreciated the fidelity of the proposed solution in comparison to the reference methods.
To further highlight the different behaviours a zoomed detail is shown in Fig.\ref{fig:DET} 
whose discussion is skipped leaving the Reader to make his/her own subjective evaluation.

\begin{figure}
\centering
\begin{tabular}{ccc}
\pgfuseimage{Fx} & \pgfuseimage{Fz} &  \\
$\widetilde{\mathbf{x}}$ & $\mathbf{z}$ & \\[2mm]
\pgfuseimage{Fhpf} & \pgfuseimage{Ffive} & \pgfuseimage{Fprop} \\
 HPF \cite{Chavez1991} & Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Zoomed detail in the full-resolution domain.}
\label{fig:FULL}
\end{figure}

Finally, in Fig.\ref{fig:FULL} we show the full-resolution domain results corresponding to the same area
selected in Fig.\ref{fig:DET}. Although we do not have a ground-truth in this case,
the different behaviours among compared solutions can be easily noticed.

\section{Conclusion}
In this work we have proposed a novel CNN-based fusion method designed to double the resolution of the 20-m bands of Sentinel-2
by taking advantage of the companion higher resolution 10-m subset of bands.
The achieved results look quite promising encouraging us to further investigate on this research line in the near future.
Of a particular interest it would be the collection of a much larger training dataset that would enable to reliably train much deeper networks.

\bibliographystyle{IEEEbib}
\newcommand*{\footnotesize}{\footnotesize}
\footnotesize{

\section{Introduction}
\label{sec:intro}


Data fusion is a topic of interest for the remote sensing community
arising in such diverse formulations as cross-sensor feature extraction \cite{Scarpa2018, Errico2014},
multitemporal analysis \cite{Gaetano2014},
multiresolution fusion \cite{Masi2016, Gargiulo2018a, Scarpa2018a}.
Sentinel-2 images provide global acquisitions of multispectral images with a high revisit frequency, supplying data for services such as risk management (floods, forest fires, subsidence, landslide), land monitoring, food security/early warning systems, water management, soil protection and so forth \cite{Drusch2012}.
Due to technological constraints only four out of thirteen bands are provided at the highest resolution of 10 meters. 
The remaining bands are given at 20 or 60 meters. 
Motivated by the above consideration in \cite{Gargiulo2018a} we proposed a CNN-based data fusion technique
for the super-resolution of the 20-m short wave infrared (SWIR) band using the other higher resolution bands to recover the missing spatial details,
following the rationale behind the more familiar pansharpening problem \cite{Vivone2015, Scarpa2018a}. 

The super-resolution of an image without auxiliary bands can be addressed 
in such diverse manners ranging from a polynomial interpolation to advanced 
CNN-based modeling \cite{Dong2016,Kim2016}.
The availability of high resolution references, or guide, like in pansharpening,
allows to obtain highly accurate super-resolution, 
thanks to the fusion between the spectral information conveyed by the target image
and the spatial content gathered by the guiding band/s.
Classical multiresolution fusion methods are often categorized as component substitution (CS) or multiresolution analysis (MRA),
a survey of which can be found in \cite{Vivone2015}.
Several pansharpening methods,
based on both CS \cite{Shah2008,Tu2001}
MRA \cite{Chavez1991,Ranchin2000}, were adapted to the Sentinel-2/SWIR case and compared
in \cite{Du2016, Wang2016}.

Other very recent approaches,
both model-based \cite{Brodu2017, Lanaras2017} and data-driven CNN-based \cite{Lanaras2018}, 
aimed to super-resolve the whole Sentinel-2 dataset testify 
the relevance of the topic to the community.

In this work we propose an improved version of our previous method \cite{Gargiulo2018a} 
that counts four main integrations, (i) the use of the residual learning strategy \cite{He2016},
(ii) the batch normalization \cite{Ioffe2015}, both (i) and (ii) aimed to speed-up the learning, 
(iii) a high-pass preprocessing of the input which has already 
proven to be effective for pansharpening in \cite{Yang2017},
and, finally, (iv) the extension to all 20-m bands.

The rest of the paper is organized as follows. 
In the next section we summarize the proposed solution. 
Section 3 presents numerical and visual results. 
Finally, Section 4 draws concluding remarks.



\section{Proposed Sentinel-2 image super-resolution}

 \begin{figure*}[h]
\centering
	\pgfuseimage{net2}
	\caption{\small Top-level workflow of the proposed super-resolution method for 20-m bands of Sentinel-2. 
	Only dashed boxes are used in \cite{Gargiulo2018a}.}
	\label{fig:flowchart}
\end{figure*}
The proposed solution is summarized in the top-level flowchart of Fig.\ref{fig:flowchart}. 
At the core is a three-layer CNN block whose internal architecture is nearly the same 
of \cite{Gargiulo2018a}, the main difference being the output shape 
which counts here six bands rather than just the single SWIR.
For the sake of brevity,
we skip the description of the internal setting of the CNN
as the Reader can refer to \cite{Gargiulo2018a} to this purpose, 
focusing on the remaining architectural aspects.

First, 
it can be noticed the additional batch normalization (BN) layer which processes the input stack $(\mathbf{z}_{\rm H},\widetilde{\mathbf{x}}_{\rm H})$
that feeds the network in order to make the learning process robust with respect to the 
statistical fluctuations of the training dataset \cite{Ioffe2015}.

The BN block takes a concatenation stack as input which comprises both the high pass filtered (HPF) 10-m resolution bands of Sentinel-2, $\mathbf{z}_{\rm H}$,
and the upscaled (via bicubic interpolation) HPF version of the target 20-m resolution bands, $\widetilde{\mathbf{x}}_{\rm H}$.
The use of HPF bands in place of the original ones has been introduced in the contex of pansharpening in \cite{Yang2017}
on the basis on the intuition that low spatial frequencies do not carry relevant information about spatial details,
therefore they can be neglected eventually facilitating the training process.

Another relevant trick generally useful when dealing with deep convolutional networks,
and actually effective for super-resolution and pansharpening as well \cite{Scarpa2018a},
is the use of residual learning \cite{He2016}, which is here implemented 
through the skip connection (on top) that directly links the input (upscaled) to the output ($\hat{\mathbf{x}} = \widetilde{\mathbf{x}} + \hat{\mathbf{y}}$).
Intuitively speaking,
since the low-frequency content of the desired output is already comprised in the low resolution
input bands, it is sufficient for the network to learn only how to predict the complementary
detail component $\hat{\mathbf{y}}$ to be combined with the input $\widetilde{\mathbf{x}}$ in order to get the full-resolution product.
The practical consequence of this is that the network learns much faster, 
which is very useful to finetune the network depending on the dataset.

\subsection{Training setting}

In order to train the network's parameters
a sufficiently large number of input-output examples
and the choice of a suitable cost function to minimize on them
are required to run any learning algorithm, 
like for example the Stochastic Gradient Descent
as done in \cite{Gargiulo2018a, Scarpa2018a} and here, as well,
using the $L_1$-norm as cost.
Due to the lack of ideally super-resolved samples to be used to this purpose,
it is becoming a common practice to resort to a self-supervised learning strategy
that sounds like the so-called Wald's protocol, 
a procedure commonly used to compare objectively different pansharpening algorithms
using referenced data \cite{Vivone2015}.
Given a training Sentinel-2 sample $(\mathbf{x},\mathbf{z})$, 
distinguishing between low and high resolution bands, respectively,
a training example is generated as follows.
Both bands subsets undergo a downsampling whose band-wise anti-aliasing filters
mimic the corresponding sensor characteristics (cutoff frequencies)
and the resolution downgrade factor is equal to the scale factor between $\mathbf{z}$ and $\mathbf{x}$,
which is $2\times2$ in our case.
Said $(\mathbf{z}_\downarrow, \mathbf{x}_\downarrow)$ the reduced-resolution sample so obtained, 
the original band set $\mathbf{x}$ can now play the role of ideal output (reference or label)
corresponding to  $(\mathbf{z}_\downarrow, \mathbf{x}_\downarrow)$.
This process is repeated on the whole training image collection whose size 
must be consistent with the number of parameters to train.
The larger the network the larger must be the training dataset to avoid overfitting problems.
Additional details on training can be found in \cite{Gargiulo2018a} as this part has not been modified.

\section{Experimental Results}


\begin{table}
\centering
\def0.8{0.8}
\setlength{\tabcolsep}{3pt}
\footnotesize
\begin{tabular}{l|ccccc}  
\hline
\rule{0mm}{3mm} \textbf{Method} & \textbf{Q-index} & \textbf{ERGAS} & \textbf{HCC} & \textbf{Epochs}   & \textbf{Sec./Ep.}\\

\rule{0mm}{3mm} & (1)  & (0)   & (1) & & \\

\hline
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a} (early stop)  & $ 0.9853 $ & $ 4.950 $ &  $ 0.5686$ & 200  & 3.577 \\
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a}  & $ 0.9908 $ & $ 3.619 $ &  $ 0.8482$ & 1500 & 3.577 \\
\hline 

\rule{0mm}{3mm} Proposed (without $\mathbf{z}$)    & $ 0.9905 $ & $ 3.985 $ &  $ 0.5792$ & 200  & 1.918\\
\rule{0mm}{3mm} Proposed (no HPF)  & $ 0.9848  $ & $ 3.112 $ &  $ 0.8482$ & 200 & 3.613\\
\rule{0mm}{3mm} Proposed & $\mathbf{0.9978 }$ & $\mathbf{1.937}$ & $\mathbf{0.8819} $ & 200 & 3.597\\

\hline
\end{tabular}
\caption{Ablation study}
\label{tab:ablation}
\end{table}

\begin{table}
\centering
\def0.8{0.8}
\setlength{\tabcolsep}{3pt}
\footnotesize
\begin{tabular}{l|ccc}  
\hline
\rule{0mm}{3mm} \textbf{Method} & \textbf{Q-index} & \textbf{ERGAS} & \textbf{HCC}  \\
\rule{0mm}{3mm} & (1)  & (0)   & (1)  \\
\hline 
\rule{0mm}{3mm} bicubic (without $\mathbf{z}$)  & $ 0.9825  $ & $ 5.340 $ & $0.4444 $  \\
\rule{0mm}{3mm} PCA  & $ 0.9596 $ & $ 8.331 $ & $0.7477 $  \\
\rule{0mm}{3mm} IHS  & $0.9434 $ & $ 9.544 $ & $0.7431 $  \\
\rule{0mm}{3mm} HPF  \cite{Chavez1991} & $ 0.9909 $ & $ 3.891 $ & $0.7091 $  \\
\rule{0mm}{3mm} GS2-GLP  & $ 0.9919  $ & $ 3.672 $ & $0.7414 $  \\
\rule{0mm}{3mm} Indusion  & $ 0.9877  $ & $ 4.494 $ & $0.6844 $  \\
\rule{0mm}{3mm} M5 \cite{Gargiulo2018a}  & $ 0.9908 $ & $ 3.619 $ &  $ 0.8482$\\
\rule{0mm}{3mm} Proposed & $\mathbf{0.9978 }$ & $\mathbf{1.937}$ & $\mathbf{0.8819} $\\
\hline

\end{tabular}
\caption{Comparative analysis.}
\label{tab:comparison}
\end{table}


In order to assess the performance of the proposed method a separate test image,
a 2700$\times$2700 scene (about Rome city center),
not used in the training phase,
has been considered.
The numerical accuracy was assessed by means of commonly used indicators like the
quality (Q-)index and ERGAS (see \cite{Vivone2015} for a detailed description of these measures),
in addition to HCC \cite{Gargiulo2018a}, a correlation coefficient computed over detail components.
Comparative methods are the bicubic interpolator (it does not take any auxiliary high-resolution guide $\mathbf{z}$ as input),
just as baseline, the model M5 proposed in \cite{Gargiulo2018a} properly generalized to output all six bands $\hat{\mathbf{x}}$,
plus other classical pansharpening methods (which can be found in \cite{Vivone2015}) suitably adapted to the case of Sentinel-2.

In Tab.\ref{tab:ablation} it is summarized an ablation study restricted to the proposed method and \cite{Gargiulo2018a}.
A first observation is about the training cost in terms of computational time (last two columns). 
The tradeoff between training time and accuracy is highlighted for M5 that requires 
1500 epochs to reach an accuracy closer to that obtained by the proposed model trained with only 200 epochs.
Focusing on the proposal
it can be noticed, as expected, that by using additional guiding high-resolution bands $\mathbf{z}$ is certainly beneficial.
Moreover, the introduction of the high pass filtering provides further improvements, notably in terms of ERGAS.


Tab.\ref{tab:comparison} gathers the numerical comparison between the proposed solution and the reference methods.
The proposal clearly outperforms all comparative solutions with respect to all quality indicators considered here.


\begin{figure}
\centering
\begin{tabular}{cc}
\pgfuseimage{Bx} & \pgfuseimage{Bz} \\
$\widetilde{\mathbf{x}}_\downarrow$ & $\mathbf{z}_\downarrow$ \\[2mm]
\pgfuseimage{Br} & \pgfuseimage{Bhpf} \\
$\mathbf{x}$ (ground-truth) & HPF \cite{Chavez1991}\\[2mm]
\pgfuseimage{Bfive} & \pgfuseimage{Bprop} \\
Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Super-resolution results on a 130$\times$130 tile of the test image (Rome), in the reduced-resolution domain.}
\label{fig:BIG}
\end{figure}

\begin{figure}
\centering
\begin{tabular}{ccc}
\pgfuseimage{Dx} & \pgfuseimage{Dz} & \pgfuseimage{Dr} \\
$\widetilde{\mathbf{x}}_\downarrow$ & $\mathbf{z}_\downarrow$ & $\mathbf{x}$ (ground-truth)\\[2mm]
\pgfuseimage{Dhpf} & \pgfuseimage{Dfive} & \pgfuseimage{Dprop} \\
 HPF \cite{Chavez1991} & Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Zoomed detail in the reduced-resolution domain.}
\label{fig:DET}
\end{figure}

To complete the analysis of the results let us look at some sample result.
For the sake of brevity we only show a limited portion of the Rome image,
selecting a meaningful subset of comparative solutions.
A 130$\times$130 cropped sample is shown in Fig.\ref{fig:BIG}.
Here the reduced resolution behaviour is analysed so that we can benefit 
from the availability of a ground-truth. 
For the purpose of visualization only three out of 6 (4) bands of $\widetilde{\mathbf{x}}_\downarrow$ ($\mathbf{z}_\downarrow$) are shown 
in a false color RGB representation. 
Likewise, truth $\mathbf{x}$ and compared results are shown according with the same restriction.
On the top row is the input, split in low ($\widetilde{\mathbf{x}}_\downarrow$) and high ($\mathbf{z}_\downarrow$) resolution components.
In the middle row are the ideal super-resolution and HPF \cite{Chavez1991}.
On the bottom are gathered M5 \cite{Gargiulo2018a} and the proposal.
From visual inspection it can be appreciated the fidelity of the proposed solution in comparison to the reference methods.
To further highlight the different behaviours a zoomed detail is shown in Fig.\ref{fig:DET} 
whose discussion is skipped leaving the Reader to make his/her own subjective evaluation.

\begin{figure}
\centering
\begin{tabular}{ccc}
\pgfuseimage{Fx} & \pgfuseimage{Fz} &  \\
$\widetilde{\mathbf{x}}$ & $\mathbf{z}$ & \\[2mm]
\pgfuseimage{Fhpf} & \pgfuseimage{Ffive} & \pgfuseimage{Fprop} \\
 HPF \cite{Chavez1991} & Model M5 \cite{Gargiulo2018a} & proposed
\end{tabular}
\caption{Zoomed detail in the full-resolution domain.}
\label{fig:FULL}
\end{figure}

Finally, in Fig.\ref{fig:FULL} we show the full-resolution domain results corresponding to the same area
selected in Fig.\ref{fig:DET}. Although we do not have a ground-truth in this case,
the different behaviours among compared solutions can be easily noticed.

\section{Conclusion}
In this work we have proposed a novel CNN-based fusion method designed to double the resolution of the 20-m bands of Sentinel-2
by taking advantage of the companion higher resolution 10-m subset of bands.
The achieved results look quite promising encouraging us to further investigate on this research line in the near future.
Of a particular interest it would be the collection of a much larger training dataset that would enable to reliably train much deeper networks.

\bibliographystyle{IEEEbib}
\newcommand*{\footnotesize}{\footnotesize}
\footnotesize{
