1. Introduction
Synthetic aperture radar (SAR) is a coherent imaging system [
1]. Each pixel in SAR images represents the coherent addition of scatterers from a corresponding resolution cell. These scatterers interfere, either constructively or destructively, depending on the phase of the scatterers. As such, the resulting images exhibit bright and dark pixels and are uneven, even for homogeneous regions. This phenomenon is called speckle noise and it often reduces the quality of images and complicates image interpretation [
1,
2]. This study proposes a novel speckle removal algorithm to not only suppress speckle noise but also preserve edges and textures.
The simplest speckle removal approach is spatial multi-looking [
3], which efficiently suppresses speckle noise at the cost of resolution loss. Three types of non-multi-looking processing methods have been proposed to balance spatial resolution and speckle removal performance.
The first is a local spatial filtering method proposed by Lee [
4,
5,
6,
7,
8]. Representative algorithms include Kuan [
9] and maximum a posteriori (MAP) filtering [
10]. Such methods have been implemented in the spatial domain based on Bayesian criteria and a speckle model. Although resolution is well-preserved and speckle noise is suppressed, the edges and textures are not maintained because the speckle model is unsuitable for filtering areas containing strong scattering points.
The second approach involves transform-domain filtering methods, such as linear minimum mean-square error (LMMSE) estimation in the wavelet domain [
11]. These methods perform multi-scale decomposition on the image, implement filtering to each decomposition image, and reconstruct the despeckling result by fusing sub-images. Since transform domain methods can distinguish edges from homogeneous areas, these techniques can more accurately preserve edges and textures compared to spatial filtering algorithms. However, these techniques are often worse for de-noising homogeneous areas than the following approach.
The third approach is adaptive filtering, which includes methods based on partial differential equations (PDEs) [
12] and non-local approaches [
13]. This PDE-based approach gradually suppresses speckle noise during iterative processing and is sensitive to edge preservation. However, repeated iterations tend to diminish texture, particularly in SAR images. The non-local methods exploit similar pixels or blocks in images to implement filtering. It obtains the most comprehensive performance in speckle suppression and preservation of edges and textures. The probabilistic patch-based (PPB) algorithm is a representative of nonlocal methods. It was proposed by Deledalle et al. in 2009 [
14]. In 2015, they proposed a unified nonlocal framework where bias-reduction was introduced to reduce the spreading of bright structures [
15].
Compared with the conventional PPB, the proposed algorithm achieves a more accurate weighting and homogeneous factor to improve the performance of speckle suppression, with a modified bias-reduction method to further balance speckle suppression with the correction of bright structure spreading.
This paper is structured as follows. The conventional PPB algorithm is introduced and analyzed in
Section 2. The three-step algorithm is then proposed to compensate for the limitations of these existing techniques in
Section 3.
Section 4 presents and analyzes corresponding results by comparing the proposed algorithm with conventional PPB, and
Section 5 concludes the paper.
2. Conventional PPB Algorithm
As illustrated in
Figure 1, P
s represents a pixel to be processed in the SAR image. A search window (centered on P
s) is defined to estimate the intensity of P
s, as represented by the pink rectangle in the figure. The conventional PPB algorithm calculates the weight
between P
s and the pixel (P
i) in the search window and replaces the intensity of P
s with [
14]:
where
Ds represents a set composed of pixels in the search window and
denotes the original intensity of P
i.
This patch region is represented by the cyan rectangle in
Figure 1. The weight
can be calculated as [
14]:
where
As,k and
Ai,k are the amplitudes of the
kth pixels in the two patches centered on P
s and P
i, respectively. The greater the weight, the more similar P
s and P
i. The term
L is the equivalent number of looks and
h is defined as [
14]:
and
q is given by:
where
and
denote the expectation and cumulative distribution functions, respectively. A bias reduction method was developed to reduce the spreading of bright structures and the intensity of P
s was modified as follows [
15]:
where
is the intensity after applying bias reduction and
is the intensity of P
s in the raw SAR image. The homogeneous factor (
) corresponding to P
s is given by
and
where
is defined on the interval [0, 1]. If P
s is in the completely homogeneous area,
equals 0. If P
s is in the bright structures,
tends to 1. A TerraSAR-X image with the resolution of one meter and the processing results acquired by applying the conventional PPB algorithm are shown in
Figure 2.
Figure 2a displays a raw unquantized single-look image, where the maximum and minimum intensities are 3.68 × 10
7 and 0, respectively.
Figure 2b shows the result processed by Equation (1), and
Figure 2c shows the result processed by Equations (1) and (5).
The comparison between
Figure 2a and
Figure 2b demonstrates the extent of speckle noise suppression achievable with Equation (1). However, the high intensity of the strong scattering targets present in the patches negatively affect the estimation using Equation (1).
Figure 2a includes three patches (centered at P
1, P
2, and P
3) with intensities of
,
, and
, respectively. The corresponding weights were
and
by applying Equation (2). It is worth noting that
, which indicates that P
2 is much more similar to P
1, whereas the product terms satisfy
. As a result, the contribution of the dissimilar point (P
3) is higher when estimating the intensity of P
1 in Equation (1). This improves the filtering result, which degrades speckle suppression performance. This effect is evident near bright structures, and widens edges and increases the size of strong scattering targets. This effect is referred to as the spreading of bright structures and can be seen in
Figure 2b. The performance of speckle suppression can be further improved by considering the impact of bright structures, which will be discussed in
Section 3.1.
Equation (5) was used to correct for the spreading of bright structures by moderately restoring the original intensities of pixels according to the factor
, as shown in
Figure 2c. However, speckle noise was also restored, particularly near bright structures. This occurred because of the inverse relationship between speckle suppression and the spreading correction in Equation (5). The value of
obtained from Equation (6) w typically close to 1 for pixels near bright structures. As such,
tends to
in Equation (5), which indicates the processed results are similar to the original image and the speckle remains mostly unaffected. We investigated this limitation using two approaches.
The first approach involved calculation of a homogeneous factor
. There are three search windows centered on P
1, P
4, and P
5 in
Figure 2a. These three points were located in homogeneous areas, and we set the size of the search window to 25 results in
,
, and
. As the size of the search window decreased, a sudden decrease occurred in the homogeneous factor, as shown in
Figure 3. For example, as the size of the search window centered on P
4 decreased from 17 to 15, the homogeneous factor decreased from 0.9032 to 0.1942. This occurred because bright structures were excluded from the search window, as shown in
Figure 4. Therefore, a more accurate value of
could be determined by choosing an appropriately-sized search window to avoid interfering with bright structures. This process is discussed further in
Section 3.2. The second approach involves modifying the form of Equation (5) to balance speckle suppression with the correction of bright structure spreading, which will be discussed in
Section 3.3.
3. Three-Step Algorithm for Speckle Suppression
Figure 5 compares the proposed three-step algorithm with the conventional PPB algorithm. The conventional PPB algorithm applies Equations (1) and (5) to the raw image. In the proposed algorithm, the first step improves the calculation accuracy of the weight by pre-processing speckle noise and reducing the effects of bright structures, and better effect of speckle suppression can be obtained using Equation (1). In the second step, an iterative method is utilized to obtain a more accurate value of
by adaptively changing the size of the search window. The final step corrects for spreading and blurring of bright targets using a modified bias-reduction method.
3.1. Speckle Pre-Processing and Weight Correction
The primary objective of speckle pre-processing is to suppress speckle noise in homogeneous areas without losing edge and texture details, which reduces the influence of speckle noise on weight calculation. This study adopts the linear minimum mean-square error (LMMSE) filter for pre-processing [
11]. Although the denoising results produced by this algorithm are not ideal, it is highly suitable for preserving edges and textures.
Then, a threshold was set, which was 25 dB higher than the average intensity of the search window [
16]. Any pixels with an intensity exceeding this threshold were considered to be strong scattering points. The influence of these points on weight calculation was then considered in four cases, as demonstrated in
Figure 1.
Case 1: Patches centered on Ps and Pi do not contain any strong scattering points, which indicates that an influence of strong scattering points on weight calculation does not exist. In this case, the weight was calculated using Equation (2).
Case 2: Both Ps and Pi are strong scattering points. It was assumed that these two points are likely similar. The weight was then calculated using Equation (2).
Case 3: Either Ps or Pi was a strong scattering point (not both). In this instance, the two patches centered on Ps and Pi were thought to be completely different and the weight was accordingly set to 0.
Case 4: The patches centered on Ps or Pi contained strong scattering points, none of which were Ps or Pi. In order to reduce the impact of the strong scattering points, the weight was then determined from Equation (2), in which all intensities for strong scattering points were replaced by the average intensity of the patch.
3.2. Iteratively Calculating the Homogeneous Factor
As illustrated in
Figure 3, the homogeneous factor
is dependent on the size of the search window. Therefore, the simplest approach to improving the accuracy of
was to reduce the window size. However, this also reduces the number of similar pixels and degrades speckle suppression performance in homogeneous areas. As such, an iterative method was developed to adaptively maximize the search window without affecting the accuracy of
. The details of this process are as follows.
(1) The initial side length of the search window centered on Ps is set to and the corresponding homogeneous factor is calculated using Equation (6). Bright structures have little effect on this calculation. If , the estimation of the homogeneous factor for Ps is less than 0.5, which is an empirical threshold. In this case, is the final estimation. Otherwise, the iteration continues. This step can reduce the computational complexity by identifying pixels that require homogeneous factor correction.
(2) Let (i = 1, 2, …), and the corresponding homogeneous factor is determined using Equation (6). If is less than the minimal size of the search window (i.e., 3 × 3), the iteration terminates and represents the final estimation. Otherwise, the process continues to step (3).
(3) The ratio
is calculated as:
The value of
decreases dramatically if the region does not contain any bright structures, as illustrated in
Figure 3. Therefore, if
is less than 0.5, indicating the homogeneous factor is less than half the previous value, the iteration terminates and
is the final estimation. Otherwise, the process continues to step (2) when
i equals 1, or step (4) when
i is greater than 1.
(4) The ratio
is calculated as:
The iteration terminates if is less than 0.5 and becomes the final estimation, as in the previous step. Otherwise, the process returns to step (2).
3.3. Correcting the Spreading and Blurring of Bright Targets
As aforementioned, the minimum size of the search window was set to 3 × 3. Therefore, homogeneous factors were updated by applying the methods proposed in
Section 3.2., with the exception of 3 × 3 regions surrounding bright structures. A modified bias-reduction method is proposed to reduce the spreading of these bright structures.
A new ratio
can be defined as:
which indicates whether significant spreading occurs or not. Equation (5) can be modified to balance speckle suppression with the correction of bright structure spreading as follows:
where
satisfies the following conditions:
Equations (9) and (10) demonstrate that equals in areas that do not exhibit bright structure spreading. The level of speckle suppression is maintained in such areas.
(2) When
, indicating the presence of spreading, the following condition is satisfied:
where
when
,
is less than
, as illustrated in
Figure 6, where
n is a parameter to balance speckle suppression with the correction of bright structure spreading. In the conventional PPB algorithm,
, which corrects for the spreading of bright structures but degrades speckle noise suppression, as discussed in
Section 2. In contrast, for
,
tends to
, which improves the performance of speckle suppression but induces obvious bright structure spreading. Equation (11) makes
, which results in more balanced performance with some suppression of both spreading and speckle.
Figure 7 demonstrates the impact of
n on the speckle suppression and correction for the spreading of bright structures. From left to right, the values of
n for these images are 1, 5, 10, 20, and 50. When
n = 1, the speckle noise in the corresponding image was the most serious. As
n increased, the speckle noise was more effectively suppressed, while the spreading of bright structures worsened. When the value of
n was between 5 and 10, a more balanced performance was obtained. In this study, the value of
n was set to 5.
In the filtering process described by Equation (9), the bright structures are also suppressed and blurred. A matrix denoted by was developed to recover these structures. This matrix is the same size as the image, and each element in this matrix corresponds to the homogeneous factor of a pixel. Bright structures in SAR images can be positioned from the matrix by the canny operator, after which is directly set to the original intensity of these bright structures.
4. Experimental Results and Analysis
Four TerraSAR-X images were used to validate the proposed algorithm, as illustrated in the first column of
Figure 8. Among these,
Figure 8a1 and
Figure 8b1 exhibit clear edges and uniform backgrounds, whereas
Figure 8c1 and
Figure 8d1 include complex structures. All these images contain strong scattering points. These characteristics help demonstrate the comprehensive performance of the proposed technique. The results achieved using the fast non-local means algorithm [
17], conventional PPB algorithm, and proposed three-step algorithm are shown in the second, third, and fourth columns of
Figure 8, respectively.
Several quantitative metrics were used to evaluate
Figure 8: the equivalent number of looks (ENL) [
18], the edge preservation index (EPI) [
19], the mean
, and the standard deviation
of the ratio image [
20,
21]. The results of this evaluation are presented in
Table 1. The terms ENL
1 and ENL
2 were calculated using the areas enclosed by the red frames, labeled 1 and 2, respectively.
Processing and evaluation results indicated that all three algorithms significantly suppress speckle noise. The fast non-local and conventional PPB algorithms have basically the same ability in speckle suppression, which is indicated by the ENL value. The fast non-local algorithm performed the worst in edge preservation. The proposed algorithm produced the highest ENL and EPI values, indicating that it was most successful in both preserving edges and suppressing speckle.
A point-to-point comparison of the texture preservation results is shown in
Figure 9. These images were produced using the ratio between raw and de-speckled data, with corresponding evaluation results shown in the last two columns of
Table 1. The application of an ideal despeckling algorithm would produce a ratio image containing only speckle points, indicating that the mean and standard deviation of the ratio image would be 1 and
, respectively, for an
L-look raw image [
3]. As all raw SAR images in this study were single-look complex images, the ideal mean and standard deviation were both one. As shown in the second column of
Figure 9, the ratio images obtained by the fast non-local algorithm contained bright structures, so the mean and standard deviation of the ratio images were far from one. Ratio images corresponding to the conventional PPB algorithm are shown in the third column of
Figure 9. They contain obvious geometric structures related to the original images, indicating that not only speckle noise but also textures were removed by the conventional PPB algorithm. In contrast, the ratio images produced using the proposed technique exhibited much weaker geometric structure, as shown in the fourth column. This indicates that the proposed algorithm can preserve texture details, with a mean and standard deviation of ratio images closer to one compared with the conventional PPB algorithm. These results demonstrate the superior performance of the proposed method.