Comparison of Accelerated Versions of the Iterative Gradient Method to Ameliorate the Spatial Resolution of Microwave Radiometer Products

Abstract

In this study, the enhancement of the spatial resolution of microwave radiometer measurements is addressed by contrasting the accuracy of a gradient-like antenna pattern deconvolution method with its accelerated versions. The latter are methods that allow reaching a given accuracy with a reduced number of iterations. The analysis points out that accelerated methods result in improved performance when dealing with spot-like discontinuities; while they perform in a similar way to the canonical gradient method in case of large discontinuities. A key application of such techniques is the research on global warming and climate change, which has recently gained critical importance in many scientific fields, mainly due to the huge societal and economic impact of such topics over the entire planet. In this context, the availability of reliable long time series of remotely sensed Earth data is of paramount importance to identify and study climate trends. Such data can be obtained by large-scale sensors, with the obvious drawback of a poor spatial resolution that strongly limits their applicability in regional studies. Iterative gradient techniques allow obtaining super-resolution gridded passive microwave products that can be used in long time series of consistently calibrated brightness temperature maps in support of climate studies.

1. Introduction

Climate studies call for long time series of reliable measurements with a synoptic point of view. In order to make such a measurement, space observational techniques are of paramount interest, and sensors onboard low Earth observing (LEO) satellites provide an essential source of information to measure key geophysical quantities of interest. Among them, microwave radiometers are sensors that are operated for a long time, ensuring a long time series at a resolution scale that is usually meant to be adequate for climate studies. In particular, there is a growing interest to monitor the evolution of fragile ecosystems such as coastal or polar areas. Coastal areas are usually economically relevant areas both for the civil infrastructure built and the presence of a large population fraction. For example, 52% of the U.S. population lives in such areas with a population density that is from 3 to 4 times the average one [1] and there is a major concern about the climate change impact on such communities. The interdisciplinary study detailed in [1] calls for reliable and independent measurements and this role can be played by satellite microwave remote sensing integrating large-scale sensors such as microwave radiometers and small-scale sensors. In view of this, it is important to provide microwave radiometer products at a spatial grid that is comparable with the one resulting from smaller-scale sensors. Another fundamental scenario is the polar regions. For instance, in Antarctica, there is a major concern about some specific regions due to the increase of 5° of the upper ocean temperatures. In fact, it has been established that the Antarctic circumpolar current is warming more rapidly than the global ocean as a whole and this has a direct impact on the polar ice system [2]. Again, the availability of spatial resolution-enhanced products is a key asset to allow a finer spatial resolution analysis of those areas that may reveal details associated with smaller-scale phenomena. As a matter of fact, the intrinsic benefits of microwave radiometer products, i.e., dense temporal and spatial coverage over the years, are strengthened by the capability to generate new spatially enhanced products.

In addition to the data fusion among different sensors, spatial resolution-enhanced radiometer products are also key assets to fully exploit the powerful capabilities of multi-channel microwave radiometers. This category includes instruments that allow collecting measurements exploiting frequency and polarisation diversity. For example, the seven-channel four-frequency Special Sensor Microwave/Imager (SSM/I), the 24-channel (ranging from 19 to 183 GHz) Special Sensor Microwave Imager Sounder (SSMIS), the 12-channel Advanced Microwave Scanning Radiometer 2 (AMSR2), and the 10-channel Microwave Radiation Imager (MWRI). These multi-channel microwave radiometers enable robust and effective estimation of the geophysical parameter of interest by a proper combination of measurements collected at different frequencies and/or polarization. However, since the footprint size changes according to the adopted frequency channel, multi-frequency estimation of geophysical parameters is not at all straightforward when the spatial grid of the finer spatial resolution channel is not degraded to match the coarser resolution grid. An alternative approach consists, even in this case, of enhancing the spatial resolution of the microwave radiometer product. In the last decades, several approaches have been proposed in the literature to achieve microwave radiometer products at an enhanced spatial resolution. Those products can be generated using operational methods that can be roughly classified into two groups: data-fusion approaches and methods based on the inversion of the antenna pattern. The former ones basically consist of introducing additional information mainly coming from complementary observations [3]. For instance, in [4] the multi-channel capability of MWRs is exploited to enhance the spatial resolution of the low-resolution channels using the finer-resolution channels according to a technique called smoothing filter-based intensity modulation (SFIM). On the other side, methods based on the inversion of the antenna pattern allow reconstruction of the finer-resolution product by exploiting the oversampling that comes with the microwave radiometer measurement mechanism. The improvement in spatial resolution is achieved at the expense of noise amplification, leading to noisier products.

In this study, the focus is on the second group of methods, which, due to the large-scale nature of the microwave radiometer measurements, are typically iterative methods, i.e., techniques that enhance the spatial resolution of the product through a series of iterations. Those methods may result in a low convergence rate that triggered the development of acceleration procedures [5,6,7] that help limit the computational burden while processing long time series of radiometer data. Within this context, we first show a typical brightness field measured by a microwave radiometer with the nominal spatial grid; then, the field reconstructed at finer spatial resolution using the canonical iterative gradient method is shown. Finally, we contrast two accelerated versions of the canonical gradient-like method to discuss their reconstruction’s accuracy with respect to convergence speed and noise amplification. The remainder of the article is organized as follows: in Section 2, the iterative schemes used in this study are introduced. Numerical experiments performed over simulated radiometer profiles are discussed in Section 3 using quantitative metrics. In Section 4, the conclusions are drawn.

2. Methodology

In this section, first, the theoretical background that underpins antenna pattern deconvolution for microwave radiometer spatial resolution enhancement is presented. Then, the two accelerated methods proposed to speed up the convergence rate of the canonical gradient-like iterative procedure are discussed.

Neglecting atmospheric effects, the generic t-th measurement provided by a given channel of an MWR can be represented as a Fredholm operator of the first kind:

$$T_{At}={\overline{G}}_t^{-1}{\int \int }_{\mathsf{\Omega }}{G}_t\left(s\right)T_B\left(\mathbf{s}\right)ds$$

(1)

where $T_A(·)$ denotes the known radiometer measurements, $T_B(·)$ is the unknown brightness field to be reconstructed, and ${\overline{G}}_t\equiv {\int \int }_{\mathsf{\Omega }}{G}_t\left(s\right)ds$, is the normalizing constant value with $G_t(·)$ being the smooth integral kernel that depends on the projection of the integrated antenna pattern onto the surface. Here $s=(\gamma ,\theta )$ is the integration variable, where $\gamma$ and $\theta$ are the azimuth and the elevation angle [8]. The linearity of the integral operator and the band-limited nature of the system function $G_t(·)$ allow recasting (1) in a discrete setting as follows [9]:

$$A\mathbf{x}=\mathbf{b}+\delta ,$$

(2)

where $\mathbf{x}\in {\mathbb{R}}^n$ is the unknown brightness temperature, $\mathbf{b}\in {\mathbb{R}}^m$ is the measured antenna temperature, $A\in {\mathbb{R}}^{m\times n}$ is the matrix that describes the projected antenna pattern with $m<n$, and $\delta \in {\mathbb{R}}^m$ is the additive Gaussian noise that affects the measurements.

Antenna pattern deconvolution methods attempt at inverting the linear system (2) to retrieve $\mathbf{x}$ on a finer spatial grid. In general, the discrete problem to be solved can be fully determined, over-determined, or under-determined. The latter is the case of most remote sensing applications and, in particular, it is the case that applies when dealing with the enhancement of the spatial resolution of microwave radiometer measurements, since one aims at reconstructing the unknown function at a finer spatial resolution [10,11]. Since the problem is underdetermined, i.e., there is no unique solution, additional constraints must be applied to reconstruct the right unknown function. In addition, due to the ill-posedness of the continuous problem (1), the corresponding equation (Equation (2)) is an ill-conditioned discrete problem that requires regularization methods to ensure the stability of the solution [9,12,13].

The general approach to solve (2) consists of minimizing the square of the 2-norm of the residual $A\mathbf{x}-\mathbf{b}$, i.e., minimizing the least square functional ${\mathsf{\Omega }}_2$ defined as:

$${\mathsf{\Omega }}_2\left(\mathbf{x}\right)=\frac{1}{2}{\parallel A\mathbf{x}-\mathbf{b}\parallel }_2^2=\frac{1}{2}\sum _{i=1}^m{(A\mathbf{x}-\mathbf{b})}_i^2$$

(3)

In the past years, several regularization methods have been proposed to deal with spatial resolution enhancement of microwave remotely sensed products. These methods can be broadly split into two main categories: direct and iterative. Direct methods, e.g., truncated singular value decomposition [8,12], the Tikhonov method [14], and Backus–Gilbert [15], are based on standard decomposition positions in numerical linear algebra. The key advantages of direct methods rely on their computer-time effectiveness; however, they are severely limited by the size of the problem to be inverted. In particular, they are generally not suitable for large-size problems arising from microwave remote sensing. On the other side, iterative methods are better suitable for large-scale problems since they converge to the derived solution through an iterative scheme where the number of iterations plays the role of regularization parameter. Among the iterative methods, the basic kernel to minimize (3) is provided by the gradient method, whose $k+1$-th iteration is given by:

$${\mathbf{x}}_{k+1}={\mathbf{x}}_k-{\lambda }_k{\mathbf{p}}_k$$

where ${\mathbf{x}}_0$ is an arbitrary initial brightness temperature field, ${\lambda }_k>0$ is the step size, and ${\mathbf{p}}_k$ is an ascent direction at point ${\mathbf{x}}_k$ for ${\mathsf{\Omega }}_2$. A wide set of popular iterative methods for the solution of $\nabla {\mathsf{\Omega }}_2\left(\mathbf{x}\right)=0$ belongs to this class, e.g., the Landweber (LW) method that calls for ${\mathbf{p}}_{k-1}=\nabla {\mathsf{\Omega }}_2\left({\mathbf{x}}_{k-1}\right)$, and the conjugate gradient (CG) method where ${\mathbf{p}}_{k-1}$ is a special linear combination of $\nabla {\mathsf{\Omega }}_2\left({\mathbf{x}}_{k-1}\right)$ and ${\mathbf{p}}_{k-2}$. For the purpose of this study, we focus on the LW whose $k+1$-th iteration reads as follows:

$${\mathbf{x}}_{k+1}={\mathbf{x}}_k-\lambda A\ast (A{\mathbf{x}}_k-\mathbf{b}),$$

(4)

where the so-called step-size $0<\lambda <2/\parallel A\ast {A\parallel }_2=2/{\parallel A\parallel }_2^2$ ensures that ${\mathsf{\Omega }}_2\left({\mathbf{x}}_{k+1}\right)<{\mathsf{\Omega }}_2\left({\mathbf{x}}_k\right)$ as well as the convergence to the generalized inverse ${\mathbf{x}}^{†}\in N{\left(A\right)}^{\perp }$ [16], where here the symbol † denotes the (Moore–Penrose) generalized inverse [17].

The LW method, which presents very good regularization and robustness features, results in a slow convergence rate. In the following subsections, two state-of-the-art kernels introduced to improve the convergence rate of LW are briefly reviewed.

2.1. Tikhonov-like Acceleration Technique (ILW)

In this subsection, the improved LW (ILW) technique introduced in [5] is briefly reviewed. The method consists of modifying the functional (5) by introducing a Tikhonov-like penalty term that aims at de-regularising the gradient method to speed up the convergence rate. Hence, a new functional is considered as follows:

$${\overline{\mathsf{\Omega }}}_2\left(\mathbf{x}\right)=\frac{1}{2}{\parallel A\mathbf{x}-\mathbf{b}\parallel }_2^2+\frac{\alpha }{2}{\parallel \mathbf{x}\parallel }_S^2$$

(5)

where $\alpha$ is a preassigned regularization parameter, and S is a positive semi-definite self-adjoint operator defined as [5]:

$$S\equiv \left(I-\frac{A\ast A}{{\parallel A\parallel }_2^2}\right)$$

(6)

In general, in the literature $\alpha >0$ is adopted to improve the regularisation capabilities of the iterative scheme by penalizing significantly noisy components. In the case of ILW, an opposite paradigm is considered that consists of using a negative $\alpha$ parameter, i.e., the term ${\parallel x\parallel }^2$ acts as a de-regularization term that results in an opposite behavior with respect to Tikhonov regularization boosting the convergence speed. Starting from a null initial guess, i.e., $x_0=0$, the k-th iteration of the ILW method can be written as:

$$\begin{array}{cc}\hfill {\mathbf{x}}_k& ={\mathbf{x}}_{k-1}-\lambda A\ast \left(A{\mathbf{x}}_{k-1}-\mathbf{b}\right)-{\beta }_{k}S{\mathbf{x}}_{k-1}=\hfill \\ \hfill & =(1-{\beta }_k){\mathbf{x}}_{k-1}-A\ast \left(\left(\lambda -\frac{{\beta }_k}{{\parallel A\parallel }_2^2}\right)A{\mathbf{x}}_{k-1}-\lambda \mathbf{b}\right)\hfill \end{array}$$

(7)

where ${\beta }_k=\lambda {\alpha }_k<0$ is the so-called de-regularization parameter, whose absolute value determines the weight of the penalty term ${\parallel \mathbf{x}\parallel }_S^2$ at the k-th iteration. The iteration-dependency is to ensure that ${\beta }_k$ goes to zero asymptotically, in order to make the final result less affected by noise amplification resulting from an excessive de-regularization.

2.2. Preconditioned Acceleration Technique (LW-P)

In this subsection, the preconditioned LW (LW-P) technique introduced in [6] is briefly reviewed. The key idea of this method is to reduce the number of iterations to achieve a given accuracy by using a filtered preconditioner operator [7,18,19]. The latter is an invertible matrix operator $P_f$ that, exploiting the structured nature of the $\mathbf{A}$ matrix, allows clustering at unity the singular values related to the signal subspace only. This limits the amplification of the components mainly corrupted by noise. The regularization behavior of the LW-P method is governed by the choice of the function that is used to build the preconditioner operator.

Starting from the linear system (2), an equivalent system, i.e., a least square preconditioned system, can be obtained by writing:

$$P_f^{-1}\mathbf{A}\ast \mathbf{A}x=P_f^{-1}\mathbf{A}\ast \mathbf{b}.$$

(8)

thanks to the fact that $P_f$ is invertible. Given (8), the k-th iteration of the LW-P method can be written as follows [6]:

$${\mathbf{x}}_k={\mathbf{x}}_{k-1}-\lambda P_f^{-1}A\ast \left(A{\mathbf{x}}_{k-1}-\mathbf{b}\right).$$

(9)

Equation (9) was constructed so that the solution obtained by the LW-P method is still consistent with the generalized solution of (2), as well the non-preconditioned scheme (4). In order to get a better understanding of the regularizing properties of the preconditioner operator $P_f$, it must be pointed out that the latter should fulfill the following rationale: it should approximate the generalized inverse $A\ast A$ in the recovered signal subspace while penalizing the components of the noise subspace. This means that:

$$P_f\simeq \left\{\begin{array}{c}A\ast A\mathrm{signal}\mathrm{subspace}\hfill \\ I\mathrm{noise}\mathrm{subspace}\hfill \end{array}\right.$$

(10)

3. Numerical Experiments

In this section, first, a thought numerical example is presented to demonstrate the ability of the canonical gradient method to reconstruct a brightness field at enhanced spatial resolution. Two real brightness fields, measured by the Special Sensor Microwave/Imager (SSM/I) and the microwave L-band radiometer flying onboard the Soil Moisture Active Passive (SMAP) mission, are considered. Then, the two acceleration procedures discussed in the previous section are contrasted with the canonical gradient method to analyze their convergence speed. This analysis is carried out on simulated 1D brightness profiles to better appreciate the different convergence speeds.

3.1. Spatial Resolution Enhancement with Gradient Methods

In this subsection, the added-value product resulting from the enhancement of the spatial resolution of microwave radiometer products is showcased using two data sets coming from the SSM/I and SMAP sensors. SSM-I is a seven-channel, four-frequency, linearly polarized radiometer that collects both horizontally and vertically polarized radiation at 19.35, 37.0, and 85.5 GHz and vertical only at 22.2 GHz [20]. In this study, we refer to the lowest spatial resolution channel, namely the 19.35 GHz, which calls for a spatial resolution equal to $69\times 43$ km. The SMAP is a US mission equipped with an L-band microwave radiometer whose spatial resolution is equal to $39\times 47$ km. In both cases, the reconstruction problem is formulated according to Equation (2) and reference is made to the real configurations. In the along-scan direction, the SSM/I radiometer performs 64 measurements over the nominal swath of 1400 km, while for the SMAP case, the sensor is supposed to perform uniformly spaced (every 11 km) measurements over a 1000 km swath in the along-scan direction. The aim of this experiment is to show that a super-resolution brightness temperature product can be obtained for remotely sensed MWR measurements by means of iterative gradient methods.

The actual SSM/I brightness field measured in 1988 over the Northern Europe area is shown in Figure 1a while the actual SMAP brightness field, collected in 2021 over an area that includes part of England, Ireland, and Iceland is shown in Figure 2a. As expected, both the radiometer measurements call for very blurred edges and coastlines that are not well defined. Moreover, the two islands enclosed in the red circles of Figure 1a and Figure 2a, i.e., the Island of Fionia for the SSM/I case and the Isle of Man for the SMAP case, are not well separated from the mainland. The reconstructed fields at enhanced resolution obtained using the LW method are shown in Figure 1b and Figure 2b for the SSM/I and SMAP cases, respectively. It can be noted that the LW reconstruction provides sharper edges and better-defined coastlines with respect to the original radiometer measurements. The two islands enclosed in the red boxes of Figure 1b and Figure 2b are better defined and better separated from the mainland. These kinds of super-resolution products are very important for regional studies such as the monitoring of coastal regions and islands, strongly threatened by the effects of climate change. In addition, the LW method allows reconstructing the enhanced resolution product on a given spatial grid that can be chosen in such a way as to generate spatially consistent time series of measurements collected by microwave radiometers calling for different native spatial resolutions. This rationale inspired, for instance, the Equal-Area Scalable Earth (EASE)-Grid 2 map projection [21], allowing us to include measurements obtained by different sensors or by different channels of the same sensor in long time series of consistently calibrated brightness temperature images [22]. Such time series can be employed in generating high-resolution maps of geophysical parameters that are critical in monitoring the effects of global warming and climate change. Unfortunately, these added-value enhanced spatial resolution products call for a processing time that is not very fast due to the low LW convergence speed. This is a key drawback when large time series of measurements are to be processed. Hence, in the following subsection, two child versions of LW that accelerate its convergence speed while maintaining the same reconstruction accuracy are contrasted.

3.2. Comparison between Accelerated Gradient Methods

Once the key advantages of the enhanced spatial resolution added-value product are discussed, in this sub-section we focus on the numerical issues linked to the algorithmic implementation of the iterative gradient method. In fact, although the LW method is reliable and robust, its convergence rate is low. Hence, in this section, we contrast the convergence rate of LW with its improved versions to discuss the advantages and drawbacks of the speed-up methods. It is worth noting that all the methods are able to provide the same reconstruction accuracy with a different number of iterations. Hence, to better discuss these issues, we design a specific 1D simulation framework that mimics a 1D profile measured by the SSM/I radiometer at the 19.35 GHz channel [11,20].

3.2.1. Metrics

To qualitatively and quantitatively discuss the enhancement capabilities of the three methods, four objective metrics are introduced:

• Improvement Factor (IF), defined as the ratio between the -3 dB width of the enhanced and the non-enhanced spike-like profiles, can assume values in the range $[1,\infty )$. Following its definition, the larger the IF is, the higher the enhancement capability achieved by the method.
• Relative Global Error (ERR), defined as:

  $$Er{r}_k={∥\frac{{\mathbf{x}}_{rec,k}-{\mathbf{x}}_{ref}}{{\mathbf{x}}_{ref}}∥}_2^2,$$

  (11)

  where $x_{ref}$ is the value of the reference profile and $x_{rec,k}$ is the reconstruction at the k-th iteration;
• Brightness Temperature Peak Error ($\Delta T_{B,p}$), defined as:

  $$\Delta T_{B,p}\equiv T_{B,p}^{RP}-T_{B,p}^{REC}$$

  (12)

  with $REC$ = (LW, ILW, LW-P), see Figure 3. It quantifies the accuracy of the methods to reconstruct the true brightness level of a spot-like discontinuity.
• Noise Amplification (NA), defined as the RMSE evaluated over a background area of the signal, see Figure 3. It allows a quantitative analysis of the noise amplification due to the reconstruction method.

3.2.2. Experiments

The first experiment aims at estimating the point spreading function resulting from the inversion scheme. To estimate the point spreading function, a Dirac delta function is needed that, according to the guidelines developed in [7,23] we approximate using a non-realistic reference profile (RP) that consists of a 1-pixel Kronecker function calling for a ${10}^6$ K brightness temperature. This allows discussion of the enhancement of the spatial resolution by contrasting the point spreading function in the non-enhanced and enhanced cases. The noisy radiometer measured profile (MP), generated using the direct problem of (2) and adding a zero-mean Gaussian additive noise (AWGN) $\delta$ with ${\sigma }_{noise}=1.06$ K standard deviation, is depicted together with the LW (yellow line), ILW (green line), and LW-P (red line) reconstructions in Figure 4a, using a dB scale. Note that the RP (shown as a cyan line) is interpolated onto the finer spatial grid using conventional linear interpolation. The relative global error is depicted in Figure 4b.

It can be noted that, although all the methods improve the spatial resolution with respect to the radiometer measurements, LW-P performs best since it results in a brightness level that is the closest to the RP one. The quantitative proof of this visual analysis is given in

Table 1. IF values for the first experiment using the LW, ILW, and LW-P methods.

	LW	ILW	LW-P
IF	1.29	1.29	1.57

, where the IF metric is listed for the three methods. The latter confirms that LW-P (1.57) outperforms both LW and ILW (1.29). Due to the iteration-dependency of the ${\beta }_k$ parameter in (7), the LW and ILW methods achieve the same IF value. In fact, by visually analyzing the zoomed-in blue panel in Figure 4b, it can be noted that the ILW method achieves a better performance in the first 100 iterations, but in the subsequent iterations, its accuracy overlaps with LW method’s one.

The second and third experiments consist of reconstructing two realistic RPs, see the blue dotted lines in Figure 5a,b. The RP of Figure 5a (dotted line) mimics a spot-like discontinuity that can be observed in real environments when targeting small islands. The RP of Figure 5b (dotted line) mimics a large discontinuity that can be observed where land/sea transitions are present. Again, the corresponding radiometer noisy MPs (see Figure 5a,b, cyan lines) are generated using the direct problem of (2) with a zero-mean Gaussian additive noise $\delta$ characterized by a ${\sigma }_{noise}=1.06$ K standard deviation to match the one that identifies SSM/I measurements. Moreover, since the aim of this study is to compare the performance of the three methods at a given iteration number, three ERR thresholds are selected:

$$\begin{array}{cc}\hfill \left[\mathrm{ET}1\right]:ERR& =0.70\hfill \\ \hfill \left[\mathrm{ET}2\right]:ERR& =0.60\hfill \\ \hfill \left[\mathrm{ET}3\right]:ERR& =0.45\hfill \end{array}$$

(13)

The three methods are contrasted in correspondence with the three ERR values by analyzing the number of iterations to achieve such an ERR value, the metric NA, and the metric $\Delta T_{B,p}$. In addition, the metrics NA and $\Delta T_{B,p}$ are also evaluated at those ERR values.

The second experiment refers to a single spike-like signal (with a width of 50 samples), see the dotted blue line in Figure 5a. The LW (yellow line), ILW (green line), and LW-P (red line) methods’ reconstructions are shown in Figure 6. The reconstructions refer to ERR equal to ET1 (a), ET2 (b), and ET3 (c). Panel (d) displays the reconstruction achieved at k = 1000 iterations. The ERR values are depicted in Figure 7, where the three error thresholds defined in Equation (13) are annotated as dotted horizontal black lines and a zoomed-in plot of the first 100 iterations is provided in the blue panel on the right side of the figure. To quantitatively analyze the reconstruction’s performance, the number of iterations required to reach the three error thresholds ET1, ET2, and ET3, and the metrics are listed in

Table 2. Using two different values for the standard deviation of the noise $\delta$ in (2), the table shows the number of iterations required to reach the three error thresholds ET1, ET2, and ET3, and metrics related to the reconstruction of Figure 6. The values of the metrics at the 1000th iteration are also shown.

${\mathit{\sigma }}_{\mathbf{noise}}=1.06$
ET	LW			ILW			LW-P
	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$
ET1	41	0.32	185.34	22	0.32	190.60	31	0.32	190.13
ET2	85	0.47	131.10	52	0.45	140.85	52	0.47	144.71
ET3	585	0.91	43.72	528	0.92	43.80	169	1.02	53.13
	1000	1.14	29	1000	1.17	28	1000	3.66	−2
${\mathbf{\sigma }}_{\mathbf{noise}}=\mathbf{2}$
ET	LW			ILW			LW-P
	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$
ET1	42	0.69	185.06	23	0.89	190.01	32	0.84	189.07
ET2	88	1.17	131.14	55	1.22	139.76	54	1.27	143.46
ET3	641	2.52	45.24	584	2.53	45.29	182	2.73	52.48
	1000	2.98	34	1000	3.05	32	1000	7.29	−2

. As expected, the ILW method requires fewer iterations to reach the accuracy level of the first threshold ET1 (23) with respect to LW(42) and LW-P(32), while after the iteration $k=50$ the LW-P starts outperforming the other two methods. By analyzing Figure 6d, i.e., the reconstructions at the 1000th iteration, it can be clearly noted that LW-P reconstruction performs best, achieving a remarkable accuracy in reconstructing the top of the spike function ($\Delta T_{B,p}=-2$), while the LW and ILW methods result in similar performances ($\Delta T_{B,p}=34$ and $\Delta T_{B,p}=32$, respectively). This result is achieved at the expense of an increase in the fluctuations over the signal’s background. In fact, after 1000 iterations, LW-P calls for a noise amplification value of $\mathrm{NA}=7.29$ which is significantly larger with respect to the LW and ILW cases, i.e., $\mathrm{NA}=2.98$ and $\mathrm{NA}=3.05$, respectively. In order to test the consistency of the obtained results in presence of a higher noise level, the same experiment is repeated using ${\sigma }_{noise}=2$ K, see Table 2. Experimental results confirm that the regularization properties of the three methods can deal with a higher noise level even if the accelerated and preconditioned methods ILW and LW-P exploit de-regularization in order to outperform LW. Hence, for the ${\sigma }_{noise}=2$, one can apply the same comments of the ${\sigma }_{noise}=1.06$ case (SSM/I nominal AWGN), with the obvious difference that an higher RMSE is present.

The third experiment refers to a pulse-like signal calling for a brightness temperature $T_B=300$ K and a width of 600 samples, see the dotted blue line in Figure 8b. The LW (yellow line), ILW (green line), and LW-P (red line) methods’ reconstructions are shown in Figure 8, corresponding to ERR values of (a) ET1, (b) ET2, and (d) ET3. Figure 8d shows the reconstructions after 1000 iterations. The ERR values are depicted in Figure 9, where the three error thresholds defined in Equation (13) are depicted as dotted horizontal black lines and a zoomed-in plot of the first 100 iterations is provided in the blue panel on the right side of the figure. To quantitatively analyze the reconstruction’s performance, the number of iterations required to reach the three error thresholds ET1, ET2, and ET3, and the metrics are evaluated (

Table 3. Number of iterations required to reach the three error thresholds ET1, ET2, and ET3, and metrics related to the reconstruction of Figure 8. In the last row, the values of the metrics at the 1000th iteration are shown.

ET	LW			ILW			LW-P
	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$	Iter #	RMSE	$\Delta T_{B,p}$
ET1	13	0.18	202.92	13	0.57	187.99	14	0.27	0.27
ET2	19	0.26	166.57	18	0.60	160.47	20	0.37	0.37
ET3	29	0.38	119.32	28	0.67	115.74	31	0.54	0.54
	1000	1.65	−23	1000	1.66	−24	1000	3.53	−25

). In this case, the results are different from the single isolated spike case. In fact, all three methods manage to achieve an enhancement with respect to the blurred radiometer measurements with similar performance and accuracy. By analyzing the results in Table 3, it can be noted that the difference in reaching a given error threshold in terms of iteration number is negligible with respect to the second experiment’s case. Moreover, all three methods achieve similar reconstruction performance after 1000 iterations, with LW-P suffering from the stronger noise amplification effect on the background area of the signal. This implies that, in case of large discontinuities, e.g., land/sea transitions, there is no benefit in accelerating gradient-like methods.

4. Conclusions

A study aimed at contrasting three different iterative gradient schemes is presented: the classical LW method and two accelerated versions, i.e., the ILW and the LW-P methods. The latter aim at ameliorating the convergence rate of the LW method, which is intrinsically slow, and hence could be exploited when processing long time series of microwave passive data.

The qualitative and quantitative analysis of the reconstruction accuracy with respect to the number of iterations and the noise amplification shows that all three methods successfully provide enhanced resolution products. However, the LW-P method allows for obtaining the highest improvement factor, i.e., IF = 1.57. When dealing with a single hot spot that resembles an isolated island, the ILW method outperforms the other two methods within the first 100 iterations, while the LW-P method achieves the best results in terms of accuracy of the reconstructions, at the expense of stronger noise amplification and a higher number of iterations with respect to the ILW case. When dealing with a pulse-like signal, the difference between the three methods is negligible, and LW-P should be avoided due to a higher noise amplification with respect to LW and ILW.

From an operational viewpoint, these speed-up methods are very important since they can be used to process the long time series of measurements collected at different spatial resolutions by microwave radiometers to enable added-value products for climate studies.