Neural Comput & Applic DOI 10.1007/s00521-013-1394-y

ORIGINAL ARTICLE

SAR image despeckling using possibilistic fuzzy C-means clustering and edge detection in bandelet domain I. Shanthi • M. L. Valarmathi

Received: 18 January 2012 / Accepted: 21 March 2013  Springer-Verlag London 2013

Abstract This paper aims at edge preservation and despeckling of synthetic aperture radar (SAR) images using a novel algorithm comprising edge detection and possibilistic fuzzy C-means clustering (PFCM) in the translationinvariant second-generation bandelet transform (TIBT) domain. The edges from the SAR image are first removed using a canny operator, and TIBT and PFCM clustering are employed to decompose and despeckle the edge-removed image, respectively. The edges are then added to the reconstructed image to obtain an enhanced version of the despeckled image. The quality of the image outperforms other despeckling methods such as K-means and fuzzy C-means that do not use edge preservation techniques. Thus, the proposed algorithm effectively realizes both despeckling and edge preservation techniques. Keywords Canny operator  Edge detection  Possibilistic fuzzy C-means clustering  Image processing  Speckle  Synthetic aperture radar (SAR)  Translation-invariant bandelet transform

1 Introduction Synthetic aperture radar (SAR) is highly advantageous over other imaging radars due to its capability of imaging at alltime and all-weather conditions and its high spatial

I. Shanthi (&) Sree Sakthi Engineering College, Karamadai, Coimbatore, Tamil Nadu, India e-mail: [email protected] M. L. Valarmathi Govt College of Technology, Coimbatore, Tamil Nadu, India

resolution, but the multiplicative speckles produced by coherent imaging process make the interpretation and processing of SAR images difficult. Hence, in order to process the SAR images, the speckles must be restrained [1, 2]. Classical despeckling algorithms such as the median filter, Lee filter [3], Frost filter, adaptive filter, Kuan filter and other despeckling algorithms effectively despeckle the image, but blur the edges of the SAR images. In these schemes, depending on the window size, there exists tradeoff between the extent of speckle noise suppression and the capability of preserving fine details. Hence, an algorithm which preserves the edges of the image while despeckling is of great need. The two ways of traditional denoising are the threshold-based method [4, 5] and the statistic-based method [6–8]. The inconvenience of the threshold-based method is it requires calculation of a proper threshold failing in which a deviation exists between the estimated threshold and the exact threshold and its denoising result is nonideal. The statistic-based method supposes that the data fit an approximate model, and its denoising result is less accurate because the assumption is not precise but approximate. To avoid the limitations of the traditional denoising methods, possibilistic fuzzy C-means (PFCM) clustering [9, 10] is utilized to remove noise from an SAR image by pattern classification. It yields a better result by using the bandelet coefficients of wavelet transform—here despeckling is regarded as a pattern classification of two classes, that is, signal and noise by utilizing PFCM clustering in the bandelet domain. The merits of TIBT [11] are as follows: Bandelets [12, 13], including orthogonal and nonorthogonal one, have greater adaptability as compared to ridgelet, curvelet, continuous wavelet transform, etc. The number of its optional directions outperforms other multiscale geometric analysis tools, which can be used to effectively

123

Neural Comput & Applic

approximate the edges. The translation-invariant version of orthogonal bandelets produces more redundant data than orthogonal bandelets, and these redundant data contain the whole details which help to retain the details and remove noise. The paper is organized as follows: Section 2 describes the review of the traditional speckle filters. Section 3 presents the review of the despeckling algorithm. Section 4 describes the methodology. Section 5 presents experimental results and quality evaluation metrics used for evaluating the quality of speckle-reduction technique. Section 6 presents the simulation results on SAR color images, and the conclusion is given in Sect. 7.

"

K 1 1X ðXn Þ2 a ¼ K n¼0 2

# ð3Þ

The parameter L is the size of the image, Ym is the value of each pixel in the image, and b is the additive noise variance in Eq. (4). " # L1 1X 2 2 b ¼ ðYm Þ ð4Þ L m¼0 The main limitation of Lee filter is that it ignores speckle noise in the areas closest to edges and lines. 2.3 Kuan filter

2 Review of the speckle filters 2.1 Median filter Median filter [14] is a nonlinear filter, and it uses a kind of smoothing technique to remove noise. Median filter works best with impulse noise (salt and pepper noise) while retaining sharp edges in the image. In this filter, the center pixel or the pixel being considered within the defined window is replaced by the middle value in the numerically ordered set after sorting the neighborhood pixels. Extra computation time needed to sort the intensity value of each set is the main limitation of median filter.

The Kuan filter [16] is more superior to the Lee filter because within the filter window, it does not make an approximation on the noise variance. The Kaun filter uses different weighting function in which multiplicative model of speckle is changed into an additive linear form given in Eq. (5).   Su W ¼ 1 ð5Þ ð1 þ Su Þ Si The weighting function is calculated from the coefficients of estimated noise variation of the image, Su given in Eq. (6) pffiffiffiffiffiffiffiffiffiffi Su ¼ 1= ENL ð6Þ

2.2 Statistic Lee filter

and Si is the coefficient of the variation of the image given in Eq. (7)

The Lee filter [15] adapts the technique that if the variance over an area is high, that is, near the edges, smoothing will not be performed, and if the variance is low or constant over an area, then the smoothing will be performed. The Lee filter is based on the approach that the speckle noise is multiplicative, and then, the SAR image can be approximated by a linear model given in the Eq. (1)

Si ¼ SD=Vm

Dðm; nÞ ¼ Vm þ W  ðMp  Vm Þ

ð7Þ

where SD is the standard deviation in filter window and Vm is mean intensity value within the window. The main disadvantage of the Kuan filter is that for computation, the value of equivalent number of looks (ENL) parameter is required.

ð1Þ 2.4 Frost filter

where D(m,n) is the grayscale value of the pixel at indices m and n after filtering. If there is smoothing, the difference between Mp (center/middle pixel) and Vm is calculated and multiplied with a weighting function W is given in Eq. (2) and then added with Vm. If there is no smoothing, the filter will output only the mean intensity value of the filter window Vm. W ¼ a2 =ða2 þ b2 Þ 2

ð2Þ

where a is the variance of the pixel values within the filter window given in Eq. (3), K is the size of the filter window, and Xn is the pixel value within the filter window at indices n.

123

A Frost filter [17] is based on the approach that the noise variance within the filter window is adopted by applying exponentially weighting factors W as given in Eq. (8). When the variance within the filter window reduces, then weighing factor also decreases.  2 !! SD W ¼ exp  DAMP  A ð8Þ Vm DAMP is a factor that gives the extent of the exponential damping for the image. The value of the DAMP is set to typically 1, and when the damping value is larger, the

Neural Comput & Applic

damping effect is heavier. SD is the standard deviation of the filter window, Vm is the mean value within the window, and A is the absolute value of the pixel distance between the center pixel to the surrounding pixels in the filter window. The value of the filtered pixel is given in the Eq. (9) .X X Dðm; nÞ ¼ Pn  Wn Wn ð9Þ where D(m,n) is the filtered or denoised pixel value and is calculated by weighed sum of Pn which is the pixel value and Wn is the weights of each pixel in the filter window over the total weighted value of the image. In the Frost filter, the parameters are adjusted according to the local variance in each area. If the variance is high in certain areas, little smoothing occurs and the edges are retained, but when the variance is low, then the filtering will cause extensive smoothing. 2.5 Enhanced Frost The Enhanced Frost filter [18] is an extension of the Frost filter. This filter divides the SAR image into homogeneous, heterogeneous and isolated point target areas. It uses a different exponentially weighting factors W in Eq. (10) to filter each region optimally.   W ¼ exp DAMP  ðSi  Su Þ=ðSmax  Si Þ  A ð10Þ Si is the local coefficient of variation of the filter window. Si ¼ SD=Vm

Img ði; jÞ ¼

X

P n  Mn

.X

ð14Þ

Mn

In order to preserve the quality of the image, the speckle is reduced but not removed in this class and is called as heterogeneous class. In the last class, if Si is larger than Smax, the value of the filtered pixel is replaced by center pixel within the filter window and is due to the fact that isolated points with high reflectivity should be considered for analysis. The Enhanced Frost filter preserves the edges and texture of an image when compared to Frost filtering. 2.6 Gamma/MAP filter Gamma or maximum a posteriori (MAP) filter [18] minimizes the loss of texture information in gamma-distributed scenes which are the image of forested areas, agricultural lands and oceans. In this filtering technique, the smoothing process is determined by the coefficient of variation and contrast ratios of probability density functions. Gamma/ MAP filter is better than Frost and Lee filter but similar to the enhanced Frost filter except that if Si falls between Su and Smax, the filter pixel value is based on the gamma estimation of the contrast ratios within the filter window and is given in Eq. (15). pffiffiffi Dðm; nÞ ¼ ððW  ENL  1Þ  Vm þ PÞ=ð2  WÞ ð15Þ where W is a weighting function. W ¼ ð1 þ S2u Þ=ðS2i  S2u Þ

ð16Þ

and P is given as ð11Þ

Su is the speckle coefficient of variation of the image using equivalent number of looks. pffiffiffiffiffiffiffiffiffiffi Su ¼ 1= ENL ð12Þ Smax is the maximum speckle coefficient of variation of the image. rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi   ð13Þ Smax ¼ 1 þ 2=ENL There are three classes to be considered. In the first class, if Si is less than Su, then the value of the filtered pixel is replaced with the intensity mean Vm of the filter window so that the speckle noise is removed. This is the homogeneous or uniform class. In the second class, the value of Si is between the minimum and maximum values of speckle coefficient of variation, and then, the value of filtered pixel is replaced by the total weighted value calculated from the weighted sum of each pixel value Pn and the weights of each pixel Mn in the filter window over the total weighted value of the image given in Eq. (14).

P ¼ Vm2 ðW  ENL  1Þ2 þ 4  W  ENL  Vm  Mp

ð17Þ

Si is the local coefficient of variation of the filter window. Si ¼ SD=Vm

ð18Þ

SD is the standard deviation in the filter window and Vm is the mean intensity value with in the window. Su is the speckle coefficient of variation of the image using equivalent number of looks and Mp is the center pixel. pffiffiffiffiffiffiffiffiffiffi Su ¼ 1= ENL ð19Þ Smax is the maximum speckle coefficient of variation of the image. pffiffiffiffiffiffiffiffiffiffiffiffi Smax ¼ 2  Su ð20Þ The value of output denoised pixel will be replaced by the mean of filter window when Si is less than Su, but it will be replaced by center pixel of the filter window if Si is greater than Smax. The removal of speckle noise without losing any data is very difficult because all of these filters depend on the local statistical data related to the output

123

Neural Comput & Applic

denoised filter and on the occurence of the filter window over an area. When the filter window happens to cover an area that is uniform, then both the reduction in speckle noise and preservation of edges are possible. When the filter window happens to cover an edge, then the value of the output denoised pixel will be replaced by the statistical data from both sides of the edge that is from two different intensity distributions. An alternative technique is to use wavelet or bandelet transform. In wavelet transform, filtered pixel will be the window size which is variable according to the contents of the image, whereas bandelets are an orthogonal basis and are interpreted as a warped wavelet basis.

The flowchart representation is as follows:

Start

Input SAR Image

Remove edges using canny operator

Perform log transform

Compute TIBT

3 Overview of despeckling algorithm The despeckling algorithm used in this paper involves computation in the TIBT domain: A Canny operator is used to detect and remove the edges from SAR image in the spatial domain. A log-transform is applied to an edgeremoved image to convert the multiplicative speckle noise to additive noise; consecutively, TIBT and possibilistic fuzzy C-means clustering are utilized to decompose and despeckle the edge-removed image, respectively, and an inverse TIBT and an exponential operation are employed to reconstruct the image; finally, the removed edges are added to the reconstructed image. Thus, the result realizes both despeckling and preserving the edges simultaneously, and its flowchart is shown in Fig. 1.

Apply PFCM clustering to high subbands

Determine min & max clusters centers

Set the noise coefficients to zero

Preserve signal coefficients

Combine average subbands and perform inverse TIBT

4 Methodology 4.1 Canny operator

Perform exponential transform

SAR image despeckling requires the preservation of the edges while despeckling the image, so that detection and preservation of the edges of the SAR images is of great need before despeckling. The steps involved in edge detection are [19]: 1. 2.

A proper threshold TC = 0.4812 is chosen for the original SAR image ORIMG0. Edges are removed from ORIMG0 by the use of Canny operator and threshold TC. When ORIMG0(m,n) is detected as an edge point; let R(m,n) = 1; otherwise, let R(m,n) = 0 (where m is the row and n is the column). Successively, the edges are separated from ORIMG0 and saved as an edge image E and an edge-removed image ER. Edge image E is kept unchanged during denoising, and the edge-removed image is despeckled by the despeckling algorithm.

123

Add the edges to the reconstructed image

Output denoised image

Stop Fig. 1 Flowchart of proposed despeckling algorithm

3.

The edges are added back to the despeckled image. The threshold for the Canny operator is set at TC = 0.4812 which is ideal since the main edges of the original images are preserved. The threshold value should be approximately set at 0.5 because the strong

Neural Comput & Applic

edges as well as weak edges have to be detected. If the value of TC is very less between 0.1 and 0.2, then only strong edges are detected. When the value of TC is high between 0.7 and 0.8, the weak edges and also noises will be detected as edges. In order to detect only the original edges, it has been set to the value nearly 0.5, that is, 0.4812. The features of the Canny Operator are good detection, good localization and clear response. Good detection means the ability to locate and mark all real edges. Good localization is the minimal distance between the detected edge and real edge, and clear response is defined as only one response per edge. Image smoothing, differentiation, nonmaximum suppression, edge thresholding and aggregation or feature synthesis are the stages of edge detection process. The reason for choosing canny operator over other edge detection operators is due to the fact that canny operator detects the weak edges only if they are associated with the strong edges. Hence the detection of false edges is minimized. 4.2 Log-transform The log-transform is applied to an edge-removed image to convert the multiplicative noise to additive noise. It is represented by ERðm; nÞ1 ¼ C  logð1 þ jERðm; nÞjÞ

ð21Þ

where ER(m,n) is the edge-removed image with indices m as row and n as column, ER(m,n)1 is the image after taking log-transform to the edge-removed image, and C is a constant. 4.3 Mechanism of TIBT Bandelets are an orthogonal and warped wavelet basis that is adapted to geometric boundaries. The transform on functions which is defined as smooth functions on smoothly bounded domains is the main motivation of bandelets. TIBT [20–22] is a 2-D translation-invariant wavelet transform (TIWT) followed by a 1-D wavelet transform. The preservation of the edges and the fine details of the SAR image are achieved due to the fact that TIBT has increased number of redundant coefficients which enhance the relationship between the coefficients. Better denoising results are obtained by TIBT because TIBT divides each translation-invariant subband into subsquares of size 8 9 8. Wavelet transform has only 3 directions—horizontal, vertical and diagonal—whereas every 8 9 8 subsquare of TIBT has 73 optional directions which are more accurate than wavelet transform. Speckle artifacts exist in the TIBT denoising results are less due to the properties of translation-invariant and redundant

coefficient. The initial threshold for TIBT is estimated indirectly through TIWT due to the fact that there is no means for direct computation of the same. ER(m,n) has increased number of redundant coefficients, which is the resultant image after log-transform, to which TIWT is performed, and ER(m,n)1 is divided into HH, HL, LH and LL subbands. The threshold is calculated from the formula given in Eqs. (22) and (23). 

 Median Xm;n  q¼ Xm;n 2 subband HH1 ð22Þ 0:6745 3q is considered as the approximate threshold value to decompose the image ER(m,n)1 by TIBT. The best direction of each subsquare is computed by minimizing a Lagrangian function. The Lagrangian of second-generation translation-invariant bandelets is LG ¼ kl  ln k2 þCT2

ð23Þ

where l is the 1-D wavelet transform coefficients and ln is the coefficients after performing hard thresholding to l, and the hard thresholding function is  l; jlj  T ln ¼ ð24Þ 0; jlj\T where C is the number of coefficients in l and T is the threshold value. The absolute value of these coefficients is bigger than the value of T. 4.4 Clustering techniques Clustering is defined as the classification of data points into homogeneous and non-homogeneous classes or clusters. The grouping of datasets is based on the following conditions. The items in the same class which are similar to each other are known as homogeneous classes and the items in the different classes which are dissimilar to each other are said to be non-homogeneous classes. In hard thresholding which is also known as non-fuzzy clustering, the data belong to crisp clusters in which every data point belongs to only one cluster. In fuzzy clustering, every point is belonging to more than one cluster, and each datum has membership grade which indicates the degree to which the data points belong to the different clusters. 4.4.1 K-means clustering K-means is a simple, unsupervised learning algorithm used to classify the data points into a certain number of clusters (assume k clusters). Initially, k centroids are defined, one for each cluster. These centroids should be placed far away from each other as much as possible. Then each point belonging to a given dataset is associated with the nearest

123

Neural Comput & Applic

centroid. Again it is needed to recalculate k new centroids, and a new binding is done between the same data points and the nearest new centroid. This loop is continued as centroids change their location step-by-step until no more changes are done and centroids do not move any more. Finally, this algorithm minimizes an objective function, which is a squared error function given by g X h X Q¼ kK n  Vm k 2

ð25Þ

m¼1 n¼1

where kKn  Vm k2 is a chosen distance measure between a data point Kn and the cluster center Vm is an indicator of the distance of the n data points from their respective cluster centers. The flowchart representation is as follows:

Q¼

g X h X

Fuzzy C-means clustering segregates and groups the data points automatically by the use of self-organization properties [23]. When there are no training samples and number of clusters is known, FCM is suitable. Since SAR image despeckling has no training samples and it is a blind restoration or unsupervised process, this can be considered as a pattern classification technique with two classes, that is, signals and noise. FCM clustering is used to suppress the speckles since it avoids the complicated calculation of the proper threshold used in threshold-based denoising methods and statistic-based denoising methods. Since the approximation information of the image is available in low-frequency subband of TIBT, clustering need not be performed on it. The edges and noise are available in the other fine subbands, and all these subbands should be clustered and denoised. Two classes with the biggest and the smallest clustering center values are considered as the signal and preserved to reconstruct the edge-removed image, and the other classes are considered as noise and set to zero to perform denoising. Plenty of redundant data are produced by TIBT, and the above-mentioned rule is used to differentiate the signal and noise. The validation of the rule is proved by the results. Let K (K = K1, K2, K3,…, Km, …, Kh) be the dataset where h is the number of data, Km = (km1, km2, km3,…kms) and s is the feature number of datum Km. The dataset K is partitioned into g classes by using fuzzy memberships P and designate a class for each datum is the aim of FCM clustering. The iterative optimization that minimizes the objective function in FCM. FCM Clustering is based on minimization of the following objective function:

123

ð26Þ

which is subject to 0  Pm;n  1;

g X

Pm;n ¼ 1;

m¼1 h X

Pm;n [ 0;

ð27Þ

1  m  g; 1  n  h;

n¼1

where Pm;n ¼

1 Pg hkKn Vm ki2=ðb1Þ

1  m  g; 1  n  h;

a¼1 kKn Va k

ð28Þ Ph

4.4.2 Fuzzy C-means clustering

Pbm;n kKn  Vm k2

m¼1 n¼1

b

n¼1 Vm ¼ P h

ðPm;n Þ Kn

n¼1 ðPm;n Þ

b

1mg

ð29Þ

where g is the number of clusters and h is the number of data points. Pm,n is the membership of Kn that belongs to the mth class, Vm (m = 1, 2, 3, …, g) is the clustering center with s-dimensional, b is the quantity controlling clustering fuzziness b e (1,?) is the fuzzy index, and kk denotes the euclidean distance. FCM is implemented by initializing the clustering centers Vm,l (m = 1,2, …, g; l = 1, 2, …, s) randomly. To get the best values for P and V, the objective function is minimized and an iterative process is imposed on Pm,n and Vm, until the termination condition is satisfied. Finally, each datum is designated to a category of class, and clustering result is used to distinguish the signal and noise. 4.4.3 Possibilistic fuzzy C-means clustering Fuzzy C-means (FCM) algorithm is one of the most widely used fuzzy clustering techniques. In FCM algorithm, memberships are assigned to which are inversely related to the relative distance to the cluster centers. If it is equidistant from two prototypes, the membership in each cluster will be the same. The membership in each cluster does not depend on the absolute value of the distance from the two centroids as well as from the other points in the data. Due to this, noise points are created far but equidistant from the central structure of the two clusters. On the other hand, if the membership in both prototypes is equal, then such points be given very low or even no membership in either cluster. To overcome this problem, a new clustering model named possibilistic C-means (PCM) relaxes the column sum constraint so that the sum of each column satisfies the looser constraint. In other words, each element of the kth

Neural Comput & Applic

column can be any number between 0 and 1, as long as at least one of them is positive. In this case the value is known as the typicality of relative to cluster (rather than its membership in the cluster). In PCM algorithm, the optimization of the PCM objective function sometimes helps to identify outliers (noise points). However, PCM algorithm ignores noise points because PCM is very sensitive to initializations, and it sometimes generates coincident clusters. Moreover, typicalities can be very sensitive to the choice of the additional parameters needed by the PCM model. To avoid the coincident cluster problem of PCM, objective function is modified by adding an inverse function of the distances between cluster centers. This extra term acts as a repulsive force and keeps the clusters separate and thus avoids coincident clusters. Objective functions used in clustering algorithm, consider only the typicalities (possibilities), attempt to exploit the benefits of both fuzzy and possibilistic clustering. The need for both possibility (i.e., typicality) and membership values is used as a model and companion algorithm to optimize it, and this algorithm is known as fuzzy possibilistic C-means (FPCM). FPCM normalizes the possibility values, so that the sum of possibilities of all data points in a cluster is 1. Although FPCM is much less prone to the problems of both FCM and PCM, the possibility values are very small when the size of the dataset increases. A model that hybridizes FCM and PCM enjoys the benefits of both models and eliminates the problem of FPCM is called as possibilistic fuzzy C-means (PFCM) [23–25] The function Q is minimized by Qa;b ðP; S; V; K Þ ¼

g X h X ðPamn þ Sbmn ÞL2mnA

ð30Þ

m¼1 n¼1

Subject to the constraints a [ 1 and b [ 1, 0 B Pmn, Smn B 1, LmnA = ||Kn - Vm||A and g X

5.1 Quality evaluation metrics The proposed algorithm has been implemented using MATLAB. The performance of the algorithm is evaluated quantitatively for SAR image with speckle noise using the quality metrics like mean square error (MSE), root mean square error (RMSE), peak signal-to-noise ratio (PSNR), signal-to-noise ratio (SNR), equivalent number of looks (ENL), SSI, SMPI and edge save index (ESI). Let G and R denote the original and the denoised image, respectively. 5.1.1 Mean square error (MSE) Mean square error (MSE) for two P 9 Q monochrome images (G and R) where one of the images is considered a noisy approximation of the other is defined as: P1 X 1 X ½Gðm; nÞ  Rðm; nÞ2 PQ m¼0 n¼0 Q1

MSE ¼

ð33Þ

5.1.2 Root mean square error (RMSE) The root mean square error (RMSE) is the square root of the squared error averaged over P * Q window vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u Q1 P1 X u 1 X RMSE ¼ t ½Gðm; nÞ  Rðm; nÞ2 ð34Þ PQ m¼0 n¼0 5.1.3 PSNR The PSNR is most commonly used as a measure of the quality of despeckled image. The PSNR is defined as:   MAXi PSNR ¼ 20 log10 pffiffiffiffiffiffiffiffiffiffi ð35Þ MSE

Pmn ¼ 1

8 n i:e:; P 2 Mfgh and

ð31Þ

where MAX2i is the maximum intensity in the unfiltered image. A higher PSNR would normally indicate that the reconstruction is of higher quality.

Smn ¼ 1

8 i i:e:; S 2 Mfhg

ð32Þ

5.1.4 SNR

m¼1 h X

5 Experimental results

n¼1

The transpose of admissible S’s are members of the set Mfhg and S is considered as a typicality assignment of the h objects to the g clusters. The possibilistic term b Pg Ph b 2 n¼1 Smn LmnA will distribute the Smn with respect m¼1 to all h data points, but not with respect to all g clusters. Under the usual conditions placed on C-means optimization problems, the first order necessary conditions for extrema of Q and b are obtained which is stated as a theorem.

Signal-to-noise ratio (SNR) compares the level of the desired signal to the level of background noise. Larger SNR values correspond to good quality image. PP PQ 2 2 m¼1 n¼1 ðGm;n þ Rm;n Þ SNR ¼ 10 log10 PP PQ ð36Þ 2 m¼1 n¼1 ðGm;n  Rm;n Þ 5.1.5 Equivalent number of looks (ENL) This index is calculated using the following

123

Neural Comput & Applic

 ENL ¼

Mean SD

2 ð37Þ

where SD is the standard deviation. The higher the ENL value for a filter [22], the higher is the efficiency in smoothing speckle noise over homogeneous areas. 5.1.6 Edge save index (ESI) ESI reflects the edge save ability in horizontal (ESI_H) or vertical (ESI_V) direction of the despeckling algorithm. The bigger the ESI, the stronger is the edge save ability. The computation formulae are as follows:  PP PQ1   m¼1 n¼1 Rm;nþ1  Rm;n  ESI H ¼ PP PQ1  ð38Þ   m¼1 n¼1 Gm;nþ1  Gm;n  PQ PP1  Rmþ1;n  Rm;n  n¼1 m  ESI V ¼ PQ PP1  ð39Þ Gmþ1;n  Gm;n  n¼1

m

where G is the original image, R is the reconstructed image, P is the row number of the image, and Q is the column number of the image. 5.1.7 Speckle suppression index (SSI) This index is based on the equation as follows: pffiffiffiffiffiffiffiffiffiffiffiffiffi varðRÞ MeanðGÞ SSI ¼  pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Mean(RÞ Var(GÞ

ð40Þ

5.1.8 Speckle suppression and mean preservation index (SMPI) ENL and SSI are not reliable when the filter overestimates the mean value. Hence an index called speckle suppression and mean preservation index (SMPI) [26] is used. The equation of this index is as follows: sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi VarðRÞ SMPI ¼ Q  ð41Þ VarðGÞ where

Q ¼ K þ jMean(RÞ  Mean(GÞj

MaxðMeanðRÞÞMinðMeanðRÞÞ : Mean(GÞ

and

K¼

Lower values of SMPI indicate

better performance and good noise reduction capability. 5.2 Results and discussion Figure 2a shows the speckled SAR image of the terrain of Bedfordshire, southeast England, with 1-m resolution and a speckle variance of 0.04. Figure 2b shows a 1-m resolution

123

Fig. 2 Real SAR images with homogeneous regions considered for evaluation

Ku-band SAR image over Horsetrack near Albuquerque, NM. Figure 2c is the NASA JPL AIRSAR P-band polarimetric SAR image of LesLandes Forest, France. Figure 2d shows PALSAR quad-polarization SAR image over Tomakomai, Japan. In these images, two homogeneous regions I and II are considered for the calculation of quality evaluation metrics. These homogeneous regions have uniform intensity in all the places throughout the region. Even though the regions seem to be with uniform intensity level, all the pixels will not have the same intensity. PFCM is very sensitive to even very small changes in the intensity value. Because of this reason, ESI-V is lower, and MSE and RMSE are higher in PFCM. Homogeneous regions are considered for evaluation because if the algorithm is efficient in homogenous region then the overall image can also be despeckled efficiently [27]. The proposed algorithm PFCM clustering technique with edge detection in the TIBT domain is compared with other two clustering techniques—K-means clustering and Fuzzy C-means clustering. Hierarchical clustering methods have the following disadvantages: (1) Inability to scale well—the time complexity of hierarchical algorithm is at least O(m2) where m is the total number of instances which is nonlinear with the number of objects using a hierarchical algorithm is also characterized by huge I/O cost. (2) Hierarchical methods can never undo what was done previously, that is, there is no back-tracking capability. Possibilistic clustering is superior to hierarchical clustering because of these two main disadvantages, that is,

Neural Comput & Applic

ability to scale well and I/O cost is less for large number of objects (size of images) [28]. Average increase in PSNR of the PFCM is approximately 0.6. From the tables, it is observed that the edge preservation and despeckling of possibilistic fuzzy C-means yield better results. The following results are obtained after simulation of the speckled image. In Figs. 3b, 4b, 5b and 6b the edges are removed using canny operator and preserved to be added back to the final reconstructed image and the edgeremoved image is taken as the input for log transformation. In Figs. 3c, 4c, 5c and 6c the multiplicative speckle noise in converted to additive noise by log transformation. Figures 3d–6d and 3e–6e show the image after wavelet transformation and bandelet decomposed image, respectively. Then the output is used for performing PFCM clustering, inverse TIBT and exponential transform. The final reconstructed image using K-means clustering, FCM clustering and PFCM clustering are shown in Figs. 3h–6h, 3i–6i, and 3j–6j, respectively.

6 Simulation results 6.1 Simulation results of Bedfordshire See Fig. 3 and Tables 1, 2. 6.2 Simulation results of Horsetrack See Fig. 4 and Tables 3, 4. 6.3 Simulation results of LesLandes forest See Fig. 5 and Tables 5, 6. 6.4 Simulation results of Tomakomai, Japan See Fig. 6 and Tables 7 and 8.

7 Conclusion Thus, the proposed despeckling algorithm based on edge detection and PFCM in the TIBT domain can preserve edges and fine details of the SAR image and enhances the performance. The effective despeckling of the algorithm lies in the fact that TIBT and possibilistic fuzzy C-means clustering are combined. The signal coefficients of the edge-removed image in the TIBT domain can be easily separated from the noise coefficients. PFCM clustering regards denoising as a pattern classification with two

Fig. 3 Despeckling and enhancement results of Bedfordshire. a unit edge matrix, b edge-removed matrix, c after log transformation, d wavelet-transformed image, e bandelet-transformed image, f after inverse bandelet transformation, g after exponential transformation, h final reconstructed image-K-means, i final reconstructed imageFCM, j final reconstructed image-PFCM

123

Neural Comput & Applic

Fig. 4 Despeckling and enhancement results of Horsetrack. a unit edge matrix, b edge-removed matrix, c after log transformation, d wavelet-transformed image, e bandelet-transformed image, f after inverse bandelet transformation, g after exponential transformation, h final reconstructed image-K-means, i final reconstructed imageFCM, j final reconstructed image-PFCM

123

Fig. 5 Despeckling and enhancement results of LesLandes forest. a unit edge matrix, b edge-removed matrix, c after log transformation, d wavelet-transformed image, e bandelet-transformed image, f after inverse bandelet transformation, g after exponential transformation, h final reconstructed image-K-means, i final reconstructed imageFCM, j final reconstructed image-PFCM

Neural Comput & Applic Table 1 Comparison of performance of K-means and FCM and PFCM clustering for Region I of Bedfordshire Parameters

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

5.3429

5.4234

SSI

0.7948

0.8596

0.9856

19.6633 0.0133

16.7886 0.0085

11.554 0.0024

SMPI ESI-H ESI-V

0.0055

MSE (103) RMSE

0.0016

14.1252

0.0022

92.594

105.52

149.29

304.293

324.842

386.387

SNR

0.0435

0.0408

0.0344

PSNR

27.2293

27.7854

29.2639

Bold values indicate better results than other filtering methods

Table 2 Comparison of performance of K-means and FCM and PFCM clustering for Region II of Bedfordshire Parameters

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

5.3270

5.4237

SSI

0.8673

0.8596

0.9856

14.2285 0.3869

16.7887 0.3597

11.515 0.9744

SMPI ESI-H ESI-V

0.3241

MSE (103) RMSE

0.3667

14.1252

0.7659

70.231

66.117

89.056

265.011

257.132

298.4225

SNR

0.0385

0.0396

0.0344

PSNR

28.0296

28.0367

29.2729

Bold values indicate better results than other filtering methods

Table 3 Comparison of performance of K-means and FCM and PFCM clustering for Region I of Horsetrack Parameters

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

4.4323

4.3118

12.9612

SSI

0.8041

0.8153

0.9838

SMPI ESI-H

8.5208 0.2194

7.6876 0.1793

1.7507 0.0336

ESI-V

1.4987 3

MSE (10 ) RMSE

Fig. 6 Despeckling and enhancement results of Tomakomai, Japan. a unit edge matrix, b edge-removed matrix, c after log transformation, d wavelet-transformed image, e bandelet-transformed image, f after inverse bandelet transformation, g after exponential transformation, h final reconstructed image-K-means, i final reconstructed imageFCM, j final reconstructed image-PFCM

K-means clustering

12.404

1.8116

12.563

13.122

13.324

111.8419

114.5515

115.4186

SNR

0.0341

0.0331

0.0329

PSNR

29.3416

29.5966

29.6587

Bold values indicate better results than other filtering methods

classes, so the combination of TIBT with PFCM can effectively remove the noise existing in the edge-removed image. The comparison between the K-means, fuzzy

123

Neural Comput & Applic Table 4 Comparison of performance of K-means and FCM and PFCM clustering for Region II of Horsetrack

Table 7 Comparison of performance of K-means and FCM and PFCM clustering for Region I of Tomakomai, Japan

Parameters

Parameters

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

4.4328

4.3118

12.9612

ENL

5.3065

4.2257

12.2834

SSI

0.8040

0.8153

0.9838

SSI

0.6433

0.7209

0.9807

8.3461 13.8606

7.6876 13.8575

1.7508 14.1948

SMPI ESI-H

6.3832 2.4601

6.4403 0.7242

3.1762 0.3635

0.1591

ESI-V

SMPI ESI-H ESI-V

0.0965

MSE (103) RMSE

11.843 108.82

0.5238 12.592

12.794

112.214

113.114

6.0119

MSE (103) RMSE

6.083

28.97

30.94

8.2062 32.054

170.2071

175.9013

SNR

0.0341

0.0329

0.0326

SNR

0.0340

0.0322

178.91 0.0318

PSNR

29.3433

29.6637

29.7367

PSNR

29.3738

29.8533

29.9635

Bold values indicate better results than other filtering methods

Bold values indicate better results than other filtering methods

Table 5 Comparison of performance of K-means and FCM and PFCM clustering for Region I of LesLandes forest

Table 8 Comparison of performance of K-means and FCM and PFCM clustering for Region II of Tomakomai, Japan

Parameters

Parameters

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

3.3882

4.3592

12.2353

ENL

5.3696

4.2257

12.2834

SSI

0.7681

0.6772

0.9457

SSI

0.6395

0.7209

0.9807

SMPI ESI-H

5.8907 0.1641

8.044286 0.1358

1.6195 0.3723

SMPI ESI-H

5.7561 6.0917

6.4403 1.4044

3.0977 2.0778

1.5421

1.1085

3.3238

ESI-V

3.3192

3.5807

5.6165

ESI-V 3

MSE (10 ) RMSE

11.81

10.477

12.48

108.667

102.358

111.727

3

MSE (10 ) RMSE

27.715

31.901

32.14

166.477

178.608

179.712

SNR

0.0311

0.0336

0.0298

SNR

0.0339

0.0309

0.0306

PSNR

30.1462

29.4711

30.5265

PSNR

29.3882

30.1997

30.2873

Bold values indicate better results than other filtering methods

Bold values indicate better results than other filtering methods

Table 6 Comparison of performance of K-means and FCM and PFCM clustering for Region II of LesLandes forest Parameters

K-means clustering

Fuzzy C-means clustering

Possibilistic fuzzy C-means clustering

ENL

3.3397

4.3572

SSI

0.7737

0.6773

0.9457

SMPI ESI-H

3.7144 2.0069

8.0457 2.6254

1.6193 3.0107

12.2354

ESI-V

3.5527

3.0358

3.8448

MSE (103)

8.375

8.473

9.067

RMSE

91.5159

92.0499

95.2206

SNR

0.0332

0.0338

0.0312

PSNR

29.5824

29.4252

30.1212

Bold values indicate better results than other filtering methods

123

Fig. 7 Filter performance in terms of ENL, SSI, SMPI and PSNR for Region I of Bedfordshire

Neural Comput & Applic

Fig. 8 Filter performance in terms of ENL, SSI, SMPI and PSNR for Region II of Bedfordshire

C-means and possibilistic fuzzy C-means clustering techniques prove that the possibilistic fuzzy C-means clustering is more efficient in despeckling as indicated by the parameters computed (Figs. 7, 8).

References 1. Fattahi H, Zoej M, Mobasheri M, Dehghani M, Sahebi MR (2009) Windowed Fourier transform for noise reduction of SAR interferograms. IEEE Geosci Remote Sens Lett 6(3):418–422 2. Zha X, Fu R, Dai Z, Liu B (2008) Noise reduction in interferograms using the wavelet packet transform and wiener filtering. IEEE Geosci Remote Sens Lett 5(3):404–408 3. Lee G (1981) Refined filtering of image noise using local statistics. Comput Graph Image Process 15(4):380–389 4. Donoho DL (1995) Denoising by soft-thresholding. IEEE Trans Inf Theory 41(3):613–627 5. Pan Q, Zhang L, Dai G, Zhang H (1999) Two denoising methods by wavelet transform. IEEE Trans Signal Process 47(12):3401– 3406 6. Crouse M, Nowak RD, Baraniuk RG (1998) Wavelet-based signal processing using hidden Markov models. IEEE Trans Signal Process 46(4):886–902 7. Chang SG, Yu B, Vetterli M (2000) Adaptive wavelets thresholding for image denoising and compression. IEEE Trans Image Process 9(9):1532–1546 8. Portilla J, Strela V, Wainwright MJ, Simoncelli EP (2003) Image denoising using scale mixtures of Gaussians in the wavelet domain. IEEE Trans Image Process 12(11):1338–1351

9. Bezdek JC (1981) Pattern recognition with fuzzy objective function algorithms. Plenum, New York 10. Gao X, Xie W (2000) Advances in theory and applications of fuzzy clustering. Chin Sci Bull 45(11):961–970 11. Peyre´ G. Toolbox bandelets. [Online] http://www.mathworks. com/matlabcentral/fileexchange/7914 12. Mallat S, Peyre´ G (2007) A review of bandelets methods for geometrical image representation. Numer Algorithms 44:205– 234 13. Peyre G, Mallat S (2005) Discrete bandelets with geometric orthogonal filters. In: Proceedings of ICIP, Sept 2005, vol 1, pp 65–68 14. Fisher R et al. MedianFilter. http://www.dai.edu.ac.uk/HIPR2/ median.htm. Current 4 Sept 2001 15. Lee JS (1986) Speckle suppression and analysis for synthetic aperture radar. Opt Eng 25(5):636–643 16. Kuan DT et al (1987) Adaptive restoration of images with speckle. IEEE Trans ASSP 35(3):373–383 17. Frost VS et al (1982) A model for radar images and its application to adaptive digital filtering of multiplicative noise. IEEE Trans Pattern Anal Mach Intell 4(2):157–166 18. Lopes A et al (1990) Adaptive speckle filters and scene heterogeneity. IEEE Trans Geosci Remote Sens 28(6):992–1000 19. Canny J (1986) A computational approach to edge detection. IEEE Trans Pattern Anal Mach Intell PAMI-8(6):679–698 20. Coifman RT, Donoho DL (1995) Translation invariant denoising. Wavelet and statistics. Springer, New York, pp 125–150 21. Zhang WG, Zhang Q, Yang CS (2012) Edge detection with multiscale products for SAR image despeckling. IEEE Xplore Digit Libr Electron Lett 48(4):211–212 22. Roomi SMM, Kalaiyarasi D (2012) Despeckling of SAR images by optimizing averaged power spectral value in curvelet domain. Int J Inf Sci Tech 2(2):45–56 23. Zhang W, Liu F, Jiao L, Hou B, Wang S, Shang R (2010) SAR image despeckling using edge detection and feature clustering in bandelet domain. IEEE Trans Geosci Remote Sens Lett 7(1):131– 135 24. Krishnapuram R, Keller JM (1996) The possibilistic C-means algorithm: insights and recommendations. IEEE Trans Fuzzy Syst 4(3):385–393 25. Pal NR, Pal K, Keller JM, Bezdek JC (2005) A possibilistic fuzzy c-means clustering algorithm. IEEE Trans Fuzzy Syst 13(4):517–530 26. Shamsoddini A, Trinder JC (2010) Image texture preservation in speckle noise suppression. In: Warner N, Szekely B (eds) ISPRS VII symposium, 2010, IAPRS, Vienna, Austria, pp 239–243 27. Li Y, Gong H, Feng D, Zhang Y (2011) An adaptive method of speckle reduction and feature enhancement for SAR images based on curvelet transform and particle swarm optimization. IEEE Trans Geosci Remote Sens 49(8):3105–3116 28. Rokach L, Maimon O. Clustering methods. In: Data mining and knowledge discovery hand book, pp 321–352

123