SPICE-Based SAR Tomography over Forest Areas Using a Small Number of P-Band Airborne F-SAR Images Characterized by Non-Uniformly Distributed Baselines

Abstract

Synthetic aperture radar tomography (TomoSAR) has been proven to be a useful way to reconstruct vertical structure over forest areas with P-band images, on account of its three-dimensional imaging ability. In the case of a small number of non-uniformly distributed acquisitions, compressive sensing (CS) is generally adopted in TomoSAR. However, the performance of CS depends on the selected hyperparameter, which is closely related to the noise of a pixel. In this paper, to overcome this limitation, we propose a sparse iterative covariance-based estimation (SPICE) approach based on the wavelet and orthogonal sparse basis (W&O-SPICE) for application over forest areas. SPICE is a sparse spectral estimation method that achieves a high vertical resolution, and takes account of the noise adaptively for each resolution cell. Thus, it does not require the user to select a hyperparameter. Furthermore, the used sparse basis not only ensures the sparsity of the forest canopy scattering contribution, but it can also keep the original sparse information of the ground contribution. The proposed method was tested in simulated experiments and the results demonstrated that W&O-SPICE can successfully reconstruct the vertical structure of a forest. Moreover, three P-band fully polarimetric airborne SAR images with non-uniformly distributed baselines were applied to reconstruct the vertical structure of a tropical forest in Mabounie, Gabon. The underlying topography and forest height were estimated, and the root-mean-square errors (RMSEs) were 6.40 m and 4.50 m with respect to the LiDAR digital terrain model (DTM) and canopy height model (CHM), respectively. In addition, W&O-SPICE showed a better performance than W&O-CS, beamforming, Capon, and the iterative adaptive approach (IAA).

1. Introduction

Forests are the largest terrestrial ecosystem on Earth, and they play an irreplaceable role in the carbon cycle and sustainable development [1]. The vertical structure of a forest not only directly reflects the growth and development of the forest, but is also an important and necessary input for estimating above-ground biomass and storage [2,3,4]. Since interferometric synthetic aperture radar (InSAR) can provide penetrability into forest, especially in the long wavelengths such as the L-band and P-band, it has become an invaluable tool for vertical structure reconstruction over forest areas. However, it cannot discriminate the different scatterers within one resolution cell, and the estimated height is the mean of all the scatterers’ heights.

Synthetic aperture radar tomography (TomoSAR), which is a new advanced SAR technique, can solve this problem because it has the vertical resolution to achieve three-dimensional imaging. TomoSAR has been widely applied over forest areas to acquire the vertical structure [5,6]. The idea behind the concept of TomoSAR is that it combines multiple acquisitions to form an additional synthetic aperture along the vertical direction, in addition to the conventional aperture along the azimuth direction [7,8].

In order to obtain a nice tomogram, many SAR images with uniformly distributed baselines are required for tomographic focusing [9,10,11,12,13,14,15,16,17,18]. However, it is usually impossible to obtain enough images with uniformly distributed baselines in every area. Thus, there are two difficulties for TomoSAR [19,20,21]: (1) For a small number of acquisitions, the common Fourier-based TomoSAR focusing approaches bring about some imaging quality problems [22,23]. If we want to get a good tomographic focusing, we have to acquire more images. However, this would require more flights for an airborne sensor or more repeat orbits for a spaceborne platform, which would greatly increase the cost of the data acquisition. (2) The non-uniformly distributed baselines mean that the observations are non-uniformly sampled in space. To address the non-uniform sampling problem, it is common to carry out interpolation to obtain uniformly sampled observations. However, this processing is sensitive to noise and may cause some other error. In summary, there is an urgent need to achieve high-quality tomographic focusing in the case of a small number of non-uniformly distributed acquisitions.

In order to solve this problem, compressive sensing (CS)-based TomoSAR has been proposed [22,23,24,25,26,27,28,29,30,31]. CS, as a sparse estimation technique, breaks the limitation of Shannon’s sampling theorem when the signal is compressible or sparse in a transform domain [22,23,24,25,26,27,28]. It can not only greatly reduce the cost of the data acquisition, but also effectively overcomes the limitation of non-uniformly distributed baselines. Moreover, CS has a very high resolution along the vertical direction. However, CS based on convex optimization solvers such as the CVX solver requires the user to select a hyperparameter that closely depends on the noise of a pixel. And this hyperparameter directly affects the performance of CS. In general, if the value of the hyperparameter is too small, it can lead to overfitting, whereas a too large value can lead to underfitting.

In this paper, to overcome this limitation, we propose the use of the sparse iterative covariance-based estimation (SPICE) approach in TomoSAR for application in forest areas. SPICE, as a new sparse spectral estimation approach, follows the framework of the sparse estimation theorem. Accordingly, SPICE has the same advantages as CS. Moreover, SPICE exploits the minimization of a covariance matrix fitting criterion to address the sparse parameter estimation, taking account of the noise in every pixel adaptively. This successfully avoids the difficulty of selecting a hyperparameter [32,33,34]. However, the backscattered reflectivity along the vertical direction over forest areas is continuous. Thus, the SPICE algorithm cannot be directly applied in tomographic focusing over forest areas. To address this issue, we use the wavelet and orthogonal (W&O) sparse basis to perform sparse expression of the forest signal. The used sparse basis not only ensures the sparsity of the forest canopy scattering contribution, but it can also keep the original sparse information of the ground scattering contribution. For the sake of simplification, this method is named W&O-SPICE in this paper. After making the forest backscattered reflectivity power sparse, the SPICE algorithm can be used in tomographic focusing.

The objective of this paper is to develop a hyperparameter-free sparse spectral estimation method, which can provide high vertical resolution in the case of a small number of non-uniformly distributed acquisitions for the application of TomoSAR over forest areas. The rest of this paper is organized as follows. Section 2 first gives a brief introduction to the TomoSAR imaging model for sparse spectral estimation. The W&O sparse basis and the W&O-SPICE method are then explained. Section 3 describes and analyzes the results of the simulated experiments. Section 4 describes the study area and datasets, presents the tomograms in three polarizations, provides the estimated underlying topography and forest height obtained with three fully polarimetric P-band F-SAR airborne images, and describes the evaluation of the results with LiDAR measurements. A further discussion about the differences between SPICE and other TomoSAR methods is given in Section 5. Finally, our conclusions are drawn in Section 6.

2. Methodology

2.1. The TomoSAR Imaging Model for Sparse Spectral Estimation

After obtaining a stack of N multiple-baseline SAR images over the same area, it is necessary to do some preprocessing steps, including selecting the master image, co-registration, deramping, and phase calibration. The focused complex value $g_n$ at an arbitrary pixel of the nth image can then be expressed as follows [7,8]:

$$g_n={\displaystyle \int }\gamma \left(z\right)\exp \left(j{k}_z\left(n\right)z\right)dz$$

(1)

where $\gamma \left(z\right)$ denotes the complex scattering coefficients along the vertical direction, $k_z\left(n\right)=\frac{4\pi b_n}{\lambda rsin\theta }$ is the vertical wavenumber of the nth image with respect to the master track, $b_n$ is the perpendicular baseline of the nth image, $\mathsf{\lambda }$ is the wavelength, r is the slant range, and $\theta$ is the incidence angle.

Through discretizing the continuous reflectivity function along the vertical direction z by D intervals, the imaging model (Equation (1)) can be approximately written as [7,8]:

$$\mathit{g}=\mathit{A}\mathit{\gamma }+\mathit{e}$$

(2)

where g is the vector of the N observation measurements, $\mathit{\gamma }$ is the unknown discrete reflectivity vector with D elements, and $\mathit{\gamma }=[{\gamma }_1,\cdots ,{\gamma }_D{]}^T$. e is the noise vector containing N elements. A is named the steering matrix with $\mathit{A}=\left[\mathit{a}\left(z_1\right),\cdots ,\mathit{a}\left(z_D\right)\right]$, and the steering vector $\mathit{a}\left(z_d\right)$ given by:

$$\mathit{a}\left(z_d\right)={\left[\exp \left(j{k}_z\left(1\right)z_d\right),\cdots ,\exp \left(j{k}_z\left(N\right)z_d\right)\right]}^T$$

(3)

where T denotes the transpose operator.

However, over forest areas, Equation (2) is not a linear sparse model and we cannot directly apply sparse spectral estimation methods like CS and SPICE to do tomographic focusing.

If we assume that the temporal decorrelation can be ignored, then the covariance matrix of multiple SAR observations g can be expressed as [22,28,29,30]:

$$\mathbf{C}=\mathrm{E}\left(\mathit{g}{\mathit{g}}^H\right)=\mathit{A}E\left(\mathit{\gamma }{\mathit{\gamma }}^H\right){\mathit{A}}^H+{\delta }^2\mathit{I}=\mathit{A}\mathit{P}{\mathit{A}}^H+{\delta }^2\mathit{I}$$

(4)

where $\mathrm{E}\left(·\right)$ is the expectation operator. ${\left(·\right)}^H$ denotes the conjugate transpose. P is the weight matrix, and $\mathit{P}=diag\left(p_1,\cdots ,p_D\right)=diag\left({\left|{\gamma }_1\right|}^2,\cdots ,{\left|{\gamma }_D\right|}^2\right)$. ${\delta }^2$ is the unknown noise power.

In fact, Equation (4) is a linear expression, as follows [22,28]:

$$\mathit{y}=\mathit{B}\mathit{p}+\mathit{ϵ}$$

(5)

where y and $\mathit{ϵ}$ is a vector with $N^2$ elements; B is a $N^2\times D$ coefficient matrix; and p is a vector with D reflectivity powers. Specifically:

$$\begin{gathered} \mathit{y}=vec\left(\widehat{\mathit{C}}\right);{\mathit{y}}_{m+\left(n-1\right)\ast m}={\widehat{\mathit{C}}}_{m,n} \\ {\mathit{B}}_{m+\left(n-1\right)\ast m}={\mathit{A}}_m⨀conj\left({\mathit{A}}_n\right) \\ \mathit{p}=diag\left(\mathit{P}\right) \\ \mathit{ϵ}=vec\left({\mathit{\delta }}^{2}I\right) \end{gathered}$$

(6)

where $1\le m,n\le N$. $\widehat{\mathit{C}}$ is the sample data covariance matrix with $\widehat{C}=\frac{1}{K}{\displaystyle \sum }_{k=1}^{K}g\left(k\right)g^H\left(k\right)$. K is the number of looks. ${\mathit{A}}_m$ is the mth column vector of steering matrix A, and ${\mathit{A}}_n$ is the nth column vector of steering matrix A. ⨀ indicates element-wise multiplication.

When a sparse basis is used to make reflectivity power vector p sparse, Equation (5) becomes a linear sparse model. Thus, we can use spectral estimation methods like CS to do tomographic focusing [22,28,29,30]. Moreover, for the application of TomoSAR over forest areas, the parameters of interest are the reflectivity powers along the vertical direction ${\left\{p_d\right\}}_{d=1}^D$, rather than the complex reflectivity ${\left\{{\gamma }_d\right\}}_{d=1}^D$.

2.2. W&O Sparse Basis

Over forest areas, the reflectivity powers along the vertical direction are not sparse. Thus, we cannot directly apply SPICE estimator to do tomographic focusing. One effective solution is to make the vertical reflectivity power sparse.

The reflectivity power is usually regarded as being contributed from both the ground and canopy, as shown in Figure 1. For the ground backscattering reflectivity power $p_s$, its corresponding spectrum along the vertical direction has an isolated narrow peak and a very small angular spreading, like a point-like signal [35], as shown in the black line in Figure 1. As for the canopy backscattering reflectivity power $p_v$, its corresponding spectrum along the vertical direction has a wide angular spreading, which is much larger than the vertical resolution (the green line in the middle of Figure 1). Thus, the ground backscattering reflectivity power $p_s$ can be considered as sparse, while the canopy backscattering reflectivity power $p_v$ can be considered as non-sparse.

In this paper, the wavelet basis Ψ and unit orthogonal basis I (W&O) are used to perform the sparse expression [30]. The W&O sparse basis can not only ensure the sparsity of the forest canopy scattering contribution, but it can also keep the original sparse information of the ground contribution. This also makes up for the shortcoming of the conventional wavelet basis, i.e., the sparseness of the ground backscattering contribution is weak, especially for the strong ground contribution.

2.3. W&O-SPICE TomoSAR Method

SPICE, as a sparse spectral estimator, has a simple and sound statistical foundation, and takes account of the noise in a natural manner. Furthermore, it does not require the user to make a difficult selection of a hyperparameter [32,33,34].

From the analysis in Section 2.2, Equation (5) can be written as:

$$\mathit{y}={\mathit{B}}_{\mathit{v}}{\Psi }^H{p}_v^{sparse}+{\mathit{B}}_{\mathit{s}}{\mathit{p}}_{\mathit{s}}+\mathit{ϵ}$$

(7)

where ${\mathit{B}}_{\mathit{v}}$ is the coefficient matrix of the canopy contribution and ${\mathit{B}}_{\mathit{s}}$ is the coefficient matrix of the ground contribution $p_v^{sparse}$ is the canopy backscattering reflectivity power in the sparse domain with L elements.

Equation (7) can be rewritten as:

$$\mathit{y}=\mathbf{\Theta }\mathit{x}+\mathit{ϵ}$$

(8)

where $\Theta =\left[\begin{array}{cc}{\mathit{B}}_{\mathit{v}}{\Psi }^H& 0\\ 0& {\mathit{B}}_{\mathit{s}}\end{array}\right]$ is the new mapping matrix. $x=\left[\begin{array}{c}{\mathbf{p}}_v^{sparse}\\ {\mathit{p}}_s\end{array}\right]$. ${\mathit{p}}_v^{sparse}$ is the canopy reflectivity power in the sparse domain and ${\mathit{p}}_{\mathit{s}}$ is the ground reflectivity power. Moreover, a few elements of x are different from zero and $L+D\gg N^2$.

From the above analysis, Equation (8) is a sparse linear model, which is essentially the same as the application model of SPICE in References [32,33,34]. Moreover, y is the vectorization of data sample covariance matrix, which means that y is obtained with the single look. Thus, the SPICE algorithm in the case of single look [32,33,34] can be used to estimate the unknown parameters x in Equation (8).

Equation (8) can be rewritten as follows [32,33,34]:

$$y=\Phi s$$

(9)

where $\Phi =\left[\Theta \hspace{1em}{I}_{N^2\times N^2}\right]$; $s={\left[{\mathit{x}}^{\mathit{T}},{\mathit{ϵ}}^{\mathit{T}}\right]}^T$.

For the case of single look, the estimation metric of SPICE is the following weighted covariance fitting criterion [32,33,34]:

$$\mathrm{f}=\underset{\left\{{\rho }_d\right\}}{\min }\Vert {\mathit{R}}^{-1/2}{\left(\mathit{R}-\mathit{y}{\mathit{y}}^{\mathit{H}}\right)\Vert }^2$$

(10)

where f is the cost function; $\Vert$ $\Vert$ stands for the Frobenius norm of the matrices. ${\mathit{R}}^{-1/2}$ is the Hermitian positive definite square root of ${\mathit{R}}^{-1}$. R is the modeled covariance matrix of y.

According to References [32,33,34], R can be expressed as:

$$R=\Phi \Sigma {\Phi }^H$$

(11)

where $\Sigma$ is the weight matrix, which can be expressed by:

$$\Sigma =\left[\begin{array}{ccccccc}{\rho }_1& 0& \cdots & \cdots & \cdots & \cdots & 0\\ 0& {\rho }_2& 0& \cdots & \cdots & \cdots & 0\\ ⋮& 0& \ddots & ⋮& ⋮& ⋮& ⋮\\ 0& \cdots & \cdots & {\rho }_{L+D}& \cdots & \cdots & 0\\ 0& \cdots & \cdots & \cdots & {\rho }_{L+D+1}& \cdots & 0\\ ⋮& ⋮& ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& \cdots & \cdots & \cdots & \cdots & \cdots & {\rho }_{L+D+N^2}\end{array}\right]$$

(12)

where ${\left\{{\rho }_{d1}\right\}}_{d1=1}^L={\left\{{\left|x_{d1}\right|}^2\right\}}_{d1=1}^L$ are the canopy contributions in the sparse domain and ${\left\{{\rho }_{d2}\right\}}_{d2=L+1}^{L+D}={\left\{{\left|x_{d2}\right|}^2\right\}}_{d1=L+1}^{L+D}$ are the ground contributions. ${\left\{{\rho }_{d3}\right\}}_{d3=L+D+1}^{L+D+N^2}={\left\{{\left|{ϵ}_{d4}\right|}^2\right\}}_{d4=1}^{N^2}$ are the noise contributions.

The solution to the minimization of the cost function f is equal to the minimization of the following equation [34]:

$$\underset{\left\{{\rho }_d\ge 0\right\}}{\min }{\mathit{y}}^H{\mathit{R}}^{-1}\mathit{y} s.t.{\displaystyle \sum }_{d=1}^{L+D+N^2}{w}_d^2{\rho }_d=1$$

(13)

where $w_d=\frac{\Vert {\varphi }_d\Vert }{\Vert y\Vert }$. ${\varphi }_d$ is the dth ($d=1,\cdots ,L+D+N^2$) column vector of matrix $\Phi$.

Through solving Equation (13), the unknown sparse reflectivity power parameters ${\left\{x_d\right\}}_{d=1}^{L+D}$ can be estimated. After that, we can use the sparse basis to obtain the real forest backscattering reflectivity power p in Equation (5), as follows:

$$\mathit{p}=\left[\begin{array}{cc}{\Psi }^H& 0\\ 0& \mathit{I}\end{array}\right]{\left\{x_d\right\}}_{d=1}^{L+D}$$

(14)

The following

Table 1. Pseudocode of the W&O-SPICE algorithm.

Initialization
…	$s_d^0=\frac{{\varphi }_d^H\mathit{y}}{\Vert {\varphi }_d{}^2\Vert }\left(d=1,\cdots ,L+D+N^2\right)$ $w_d=\frac{\Vert {\varphi }_d\Vert }{\Vert \mathit{y}\Vert }$
Iteration
	repeat
	1. ${\rho }_d=\frac{\left|s_d\right|}{w_d}$
	2. $s_d={\rho }_d{\varphi }_d^H{\mathit{R}}^{-1}\mathit{y}$
	3. $\Sigma =\mathrm{diag}\left({\left\{{\rho }_d\right\}}_{d=1}^{L+D+N^2}\right)$
	4. $\mathit{R}=\Phi \Sigma {\Phi }^H$
	Until (convergence)

presents the detailed processes of the SPICE algorithm [34]:

(1)

Calculate the initial value according to $s_d^0=\frac{{\varphi }_d^{H}y}{\Vert {\varphi }_d\Vert {}^2}$ and $w_d=\frac{\Vert {\varphi }_d\Vert }{\Vert y\Vert }\left(\mathrm{d}=1,\cdots ,L+D+N^2\right)$, where ${\varphi }_d$ is the dth column vector of matrix $\Phi$.

(2)

Calculate every unknown parameter’s power ${\left\{{\rho }_d\right\}}_{d=1}^{L+D+N^2}$ in the ith iteration.

(3)

Calculate every unknown parameter ${\left\{s_d\right\}}_{d=1}^{L+D+N^2}$ in the ith iteration.

(4)

Calculate the weight matrix $\Sigma$ in the ith iteration.

(5)

Calculate the covariance matrix R in the ith iteration.

(6)

Calculate the convergence of the last two iterations $\Vert \Sigma \left(i+1\right)-\Sigma \left(i\right)\Vert /\Vert \Sigma \left(i\right)\Vert$. If this value is less than ${10}^{-4}$, the iteration process is terminated; if not, the iteration process continues.

3. Simulation Experiments

In order to demonstrate the feasibility and effectiveness of the proposed method, a set of simulation experiments was carried out. Due to the influence of the imaging wavelength, polarization mode, forest density, and species, different types of forest have different scattering characteristics [35]. Considering this point, we simulated three kinds of forest scattering scenarios: (1) r > 1, where the ground scattering dominates, and r is the power ratio between the ground scattering and the canopy scattering; (2) r = 1, where the ground scattering power is equal to the canopy scattering power; and (3) r < 1, where the canopy scattering is dominant. The wavelength, slant distance, and other parameters used in the simulation experiments are listed in

Table 2. Parameters of the F-SAR airborne system.

Polarization	Wavelength	Incidence Angle	Slant Range	Azimuth Resolution	Range Resolution
HH + HV + VV	0.6897 m	25–55°	6096 m	3.84 m	2.0 m

. The tomographic focusing profiles were acquired from three aspects:

(1)

The reconstruction performance of SPICE based on the wavelet basis and W&O basis, respectively, was investigated in the case of the same angular spreading with different kind of forest scattering scenarios.

(2)

The reconstruction performance of CS based on the W&O basis was investigated with different selected hyperparameter in the case of the same forest scattering scenario (r > 1).

(3)

We then compared the reconstruction performance between CS and SPICE with different ground angular spreading in the case of the same forest scattering scenario (r > 1).

Considering the implementation facility and compactness to solve the L₁-norm minimization [23,24,25,26,27,28,29,30], this paper applies the CVX solver for the CS TomoSAR method. The hyperparameter is set as 1000, which is an appropriate choice.

Based on the results of the above simulations, a number of observations can be made:

(1)

For the three kinds of forest scattering scenarios, wavelet-SPICE fails to reconstruct the vertical structure, while W&O-SPICE successfully obtains the tomogram, as shown in Figure 2.

(2)

For W&O-CS estimator, it is necessary to select the hyperparameter for all the resolution cells. Differences in the selected hyperparameter lead to differences in performance, as shown in Figure 3. When this hyperparameter is too small, the reconstructed profile will miss the backscattering of targets of interest, even bringing about some misestimations (Figure 3a). When this hyperparameter is too large, the reconstructed profile will not only reflect the backscattering of targets of interest, but also contains much noise, as shown by the red circles in Figure 3c. A nice performance can be obtained only through selecting an appropriate hyperparameter (Figure 3b). Moreover, the hyperparameter is closely related to the noise level of a pixel, and the noise level is different in different pixels. Thus, a selected hyperparameter for a resolution cell may not be suitable for another cell. When the selected hyperparameter can be adaptive to the noise in each pixel, CS can have a nice reconstructed profile.

(3)

With the reduction of the ground angular spreading, W&O-SPICE shows a good performance in retrieving the vertical structure of the forest, in all cases, and the result is almost the same as the ground truth (see Figure 4a–c). However, when the ground angular spreading is wide, W&O-CS can reconstruct the vertical structure of the simulated forest signal (Figure 4d). When the ground angular spreading reduces to moderate, W&O-CS can also detect the ground scattering, but has deviation in the canopy amplitude retrieval. When the ground angular spreading is decreased to narrow, just like the point-like signal, there are deviations in both the canopy and ground amplitude reconstruction.

From the above observations, it can be seen that wavelet-SPICE is not able to obtain correct reconstruction results because the wavelet basis only has good sparseness on the smooth signals, such as the backscattering power of the forest canopy. When the smooth signal is mixed with an approximately point-like signal, such as the ground backscattering power, the sparseness of the wavelet basis is seriously reduced. However, the W&O sparse basis can solve this problem, and it can not only ensure the sparsity of the smooth signal, but it can also keep the sparseness of point-like signals. In addition, the SPICE algorithm adaptively takes into account the noise for every resolution cell, which avoids the need to select a hyperparameter.

4. Real-Data Experiments and Results

In order to further demonstrate the feasibility and effectiveness of W&O-SPICE for a small number of non-uniformly distributed baselines, a real airborne SAR dataset was applied in TomoSAR to obtain the vertical structure of a forest.

4.1. Study Area and Dataset

The study area is located at Mabounie (Figure 5), Gabon, west-central Africa, and is 180 km away from the capital airport. It is known as a mining exploration site and is an important test site for tropical rain forests. There are both mature forests and degraded forests in the area. The weather is hot and humid, and the average annual temperature is about 26 °C. Moreover, the annual rainfall ranges between 1500 mm and 3000 mm. The topography ranges between 0 and 170 m [36,37].

Three P-band fully polarimetric SAR images with non-uniformly distributed baselines were obtained over the study area, which had been already preprocessed into tomographic dataset by the German Aerospace Center (DLR). These images were acquired by the F-SAR airborne system in the framework of the AfriSAR campaign, aiming at addressing the planning requirements for the development of the European Space Agency’s (ESA) BIOMASS mission. During the early stage of the BIOMASS feasibility and consolidation study, the vertical structure of the tropical forest had already been reconstructed as part of the TropiSAR 2009 campaign. The objective of TropiSAR 2009 was to study Amazonian forest types characterized by high biomass density (higher than 300 t/ha). However, the tropical forest around the equatorial belt is very different from the Amazonian forest, as their biomass conditions range from 100 t/ha to 300 t/ha [36]. Thus, it was necessary to develop and verify the BIOMASS mission performance on this kind of tropical forest.

The parameters of the collected data are shown in Table 2 and

Table 3. The baseline information for the InSAR pairs.

Identifier	Acquisition Date	Baseline (m)
07-02	11/02/2016	0
07-05		10
07-15		60

. The range resolution is 3.84 m and the azimuth resolution is 2.0 m. The Rayleigh resolution along the elevation is 40 m. The incidence angle varies from 25° to 55° [36].

In order to validate the TomoSAR results, we used the LiDAR data acquired by the LVIS system operated by NASA over this study area in February 2016, in a joint cooperation between the ESA, the DLR, ONERA, the ANPN, NASA, and the University of Maryland. The LVIS system was installed on a NASA Langley B200 plane with a flight height of 7315 m [36]. The LiDAR footprint was about 20 m wide on the ground. Based on these LiDAR data, the digital terrain model (DTM) and canopy height model (CHM) could be generated.

4.2. Results and Analysis

4.2.1. Tomogram of the Selected Profile

An azimuth profile (at the 200th azimuth resolution cell) was selected to perform tomographic focusing (as shown in Figure 6). The size of the estimation window applied in the experiments was 31 × 31 pixels (slant range/azimuth), corresponding to about 180 m × 60 m (ground range/azimuth).

Figure 7 shows tomograms of the selected azimuth profile in the three polarimetric modes. For the sake of a clear contrastive analysis, the LiDAR DTM and CHM were projected into SAR geometric coordinates. A number of observations can be made: (1) The vertical structure of the forest can be reconstructed from each polarimetric dataset. (2) For the HH polarization, the ground backscattering dominates. Moreover, the ground scattering phase centers detected by the W&O-SPICE algorithm show a good consistency with the LiDAR DTM. (3) For the HV polarization, the dominant scattering is the canopy backscattering, and the corresponding detected phase centers have a similar wave trend to the LiDAR CHM. (4) For the VV polarization, the W&O-SPICE algorithm detects some ground scattering phase centers and some canopy scattering phase centers.

The above observations indicate that the backscattering power for the HH polarimetric channel and HV polarimetric channel mainly comes from the ground and forest canopy, respectively. Thus, the underlying topography can be estimated from the results in the HH polarization, and the forest canopy height can be obtained from the results of the HV polarimetric channel.

4.2.2. Underlying Topography Estimation

For one pixel, the first two strongest peaks generally reflect the ground scattering and forest canopy scattering. The heights of these two peaks are respectively the ground height and the forest canopy scattering phase height. Thus, the smaller value between two heights is the ground height. By detecting all the ground peaks of the tomograms in the HH polarimetric channel, the underlying topography was estimated using the W&O-SPICE algorithm, as shown in Figure 8. In order to compare this with the LiDAR DTM, the estimated underlying topography was projected into the geographic coordinate system. Due to the influence of cloud and fog, there are many cavities in the LiDAR DTM, as shown in Figure 8a. The underlying topography estimated by the W&O-SPICE TomoSAR method ranges between 0 and 170 m, which corresponds with the ground truth. Moreover, most of the underlying topography estimated by the W&O-SPICE TomoSAR method is consistent with the LiDAR DTM, especially in the details.

In addition, to evaluate the performance of the W&O-SPICE TomoSAR method, we calculated the mean and the root-mean-square error (RMSE) of the estimated underlying topography with respect to the LiDAR DTM, and the two values are, respectively, 1.05 m and 6.40 m (see

Table 4. The mean and RMSE of the underlying topography estimated by the W&O-SPICE TomoSAR method with respect to the LiDAR DTM.

TomoSAR w.r.t LiDAR	Mean	RMSE
DTM (m)	1.05	6.40

).

4.2.3. Forest Height Estimation

Based on the analysis given in Section 4.2.1, the canopy height was estimated with the results of the HV polarimetric channel. The forest height could then be obtained as the difference between the canopy height and the ground height, as shown in Figure 9. The detailed process of canopy height estimation can be referred in References [10,11,12,28]. Since only three images were applied in the tomographic focusing, it was impossible to provide a correct estimation in some parts (the holes in Figure 9b). However, the estimated forest height is largely similar to the LiDAR CHM. Moreover, the mean and RMSE of the forest height estimation were also calculated with respect to the LiDAR CHM, and the values are, respectively, 0.02 m and 4.50 m (listed in

Table 5. The mean and RMSE of the forest height estimated by the W&O-SPICE TomoSAR method with respect to the LiDAR CHM.

TomoSAR w.r.t LiDAR	Mean	RMSE
Forest height (m)	0.02	4.50

).

The above analysis suggests that the W&O-SPICE algorithm can be used to reconstruct the vertical structure of a forest using only a small number of non-uniformly distributed acquisitions. Although the estimation accuracy of the underlying topography and forest height is not particularly high, it does demonstrate that the W&O-SPICE estimator is a good candidate to estimate the underlying topography and forest height in the case of a small number of polarimetric SAR images with non-uniformly distributed baselines.

5. Discussion

To further verify the performance of the SPICE algorithm based on the W&O sparse basis, the tomographic estimators of CS and the nonparametric spectral estimation methods (beamforming, Capon, and the iterative adaptive approach (IAA)) were also applied in tomographic focusing of the same azimuth profile.

5.1. Comparison between W&O-CS and W&O-SPICE

Figure 10 shows tomograms in the HH and HV polarizations estimated by CS based on the W&O sparse basis. It is found that the backscattering power along the vertical direction from W&O-CS is also related to the polarimetric mode: the results in the HH polarization are suitable to recognize the ground scattering, and those in the HV polarization indicate the forest canopy scattering. However, W&O-CS shows a degraded performance in reconstructing the canopy information, especially the parts marked by the purple circles in Figure 10. This is because the CS algorithm requires the user to select a hyperparameter. If this parameter is too small, it will lead to overfitting, missing some useful information. If this parameter is too large, it will result in underfitting, estimating some unnecessary information besides the parameters of interest. In addition, the selected hyperparameter makes it difficult for CS to adapt to the noise in each resolution cell. Thus, the tomogram estimated by the W&O-CS estimator is clear in some pixels, while it is seriously affected by noise in other pixels. This reduces the accuracy of the underlying topography and forest height estimation.

However, SPICE does not depend on a hyperparameter, and considers the noise of each resolution cell adaptively. Accordingly, the tomograms estimated by W&O-SPICE are almost unaffected by noise. This analysis confirms that SPICE is more suitable for application over forest areas than CS.

5.2. Comparison between the Nonparametric Spectral Estimation Methods and W&O-SPICE

Since nonparametric spectral estimation does not require any prior information, it has been widely used in TomoSAR application over forest areas. Accordingly, the commonly used nonparametric spectral estimators of beamforming, Capon, and IAA [16,38] were applied for the tomographic focusing of the selected azimuth profile in the HH and HV polarizations, as shown in Figure 11. Among the results of these tomographic estimators, beamforming shows the most serious noise in the tomograms. It also fails to separate the ground scattering contribution and the canopy scattering contribution. As for Capon, the noise is much less than for beamforming, but it only detects one scattering phase center in each resolution cell. This means that Capon also cannot separate the ground scattering contribution and the canopy scattering contribution. The reason for this is that both beamforming and Capon follow Shannon’s sampling theorem, and they have a low resolution and serious sidelobes in the case of a small number of acquisitions with non-uniformly distributed baselines. In addition, it can be clearly seen in Figure 11g,h that IAA can distinguish a few ground scattering phase centers and forest canopy scattering phase centers. However, most of the scattering phase centers are still indistinguishable. Moreover, sidelobes exist in the tomograms estimated by IAA. Therefore, these three estimators cannot be used to obtain the underlying topography and forest height in this case study. However, W&O-SPICE can successfully distinguish the ground and forest canopy scattering contributions, and it successfully estimated the underlying topography and forest height over this study area.

The above observations suggests that SPICE has a higher resolution than beamforming, Capon, and IAA in the case of a small number of non-uniformly distributed acquisitions, which is consistent with the findings in References [32,33,34].

From the above analysis, it is found that it makes it possible for TomoSAR to reconstruct the forest vertical structure with a small number of non-uniformly distributed acquisitions. This provides an opportunity for the extraction of underlying topography and forest height in a large scale. Moreover, our future works will also focus on analyzing the feasibility and effectiveness of the proposed approach for polarimetric TomoSAR over forest areas.

6. Conclusions

In this paper, we have proposed the W&O-SPICE TomoSAR method for use with a small number of acquisitions with non-uniformly distributed baselines. The W&O sparse basis is made up of the wavelet basis and orthogonal basis, which can not only ensure the sparsity of the forest backscattering power, but can also keep the sparseness of the ground backscattering power. SPICE, as a sparse spectral estimation technique, can obtain a high vertical resolution, even with a small number of acquisitions with non-uniformly distributed baselines. Moreover, it does not require the user to select a hyperparameter.

Through several sets of simulation experiments, we found that SPICE based on the W&O sparse basis can reconstruct the vertical structure of forest in three kinds of ground-to-canopy scattering power ratios. However, SPICE based on the wavelet basis cannot achieve this. Moreover, W&O-SPICE can still obtain a good reconstruction performance as the ground angular spread becomes narrower and narrower.

Three P-band fully polarimetric airborne SAR images with non-uniformly distributed baselines were applied in tomographic focusing, to demonstrate the feasibility and effectiveness of the proposed method. The results showed that the W&O-SPICE algorithm could successfully reconstruct the vertical structure of the tropical forest in Mabounie, Gabon, for the three polarimetric channels. Based on the results in the HH polarization and HV polarization, both the underlying topography and the forest height were estimated. With respect to the LiDAR DTM and CHM, the estimated RMSEs were 6.40 m and 4.50 m for the underlying topography and forest height, respectively. Furthermore, the TomoSAR methods of W&O-CS, beamforming, Capon, and IAA were also used to obtain tomograms for the selected azimuth profile. Compared to W&O-SPICE, W&O-CS could also recognize the ground and forest canopy scattering phase centers, but it was greatly affected by noise. This is because CS depends on the selected hyperparameter. As for beamforming and Capon, they were both unable to separate the ground scattering contribution and the forest canopy scattering contribution. Although IAA could distinguish a few ground scattering centers and forest canopy scattering centers, it failed to separate most of them.

In conclusion, the W&O-SPICE estimator showed advantages over the TomoSAR estimators of CS, beamforming, Capon, and IAA in this case study. In the case of a small number of non-uniformly distributed acquisitions, W&O-SPICE is better able to reconstruct the vertical structure of forest and estimate the underlying topography and forest height.