Journal of Systems Engineering and Electronics Vol. 19, No. 6, 2008, pp.1121–1126
Polarimetric whitening filter for POLSAR image based on subspace decomposition∗ Yang Jian, Deng Qiming, Huangfu Yue & Zhang Weijie Dept. of Electronics Engineering, Tsinghua Univ., Beijing 100084, P. R. China (Received September 21, 2007)
Abstract: Speckle filtering is an indispensable pre-processing step for applications of polarimetric synthetic aperture radar (POLSAR), such as terrain classification, target detection, etc. As one of the most typical methods, the polarimetric whitening filter (PWF) can be used to produce a minimum-speckle image by combining the complex elements of the scattering matrix, but polarimetric information is lost after the filtering process. A polarimetric filter based on subspace decomposition which was proposed by Gu et al specializes in retrieving principle scattering characteristics, but the corresponding mean value of an image after filtering is not kept well. A new filter is proposed for improving the disadvantage based on subspace decomposition. Under the constraint that a weighted combination of the polarimetric SAR images equals to the output of the PWF, the Euclidean distance between an unfiltered parameter vector and a signal space vector is minimized so that noises can be reduced. It is also shown that the proposed method is equivalent to the subspace filter in the case of no constraint. Experimental results with the NASA/JPL airborne polarimetric SAR data demonstrate the effectiveness of the proposed method.
Keywords: speckle filtering, synthetic aperture radar, polarimetric, polarimetric whitening filter, subspace decomposition.
1. Introduction It is well known that a coherent synthetic aperture radar (SAR) produces images with considerable speckles[1] . The presence of speckle in SAR images causes degradation and makes terrain classification and target detection difficult[2] . Polarimetric SAR can provide multi-dimensional data by transmitting and receiving electromagnetic waves with different polarizations. The availability of polarimetric SAR data has made it possible to process polarimetric radar data to reduce the image speckle, which can benefit its applications greatly. Several filters have been proposed for speckle reduction in polarimetric images. Novak and Burl[3−4] derived the polarimetric whitening filter (PWF) by optimally combining the elements of the polarimetric covariance matrix, which generates a minimumspeckle intensity image from single-look polarimetric SAR data. Although the PWF filter is theoretically
the optimal method for processing the coherent scattering matrix into pixel intensity to minimize speckle, since only one image is produced after filtering, polarimetric information can not be totally preserved. An extension of PWF for multilook case was formed as MPWF by Liu et al[5] , which processed the multilook covariance matrix to generate a minimum-speckle intensity image. But the same problem on losing polarimetric information was not solved either. Recently, a new approach to smooth speckles was presented by Gu et al[6] , based on subspace decomposition. In this approach, a parameter space consists of two orthogonal subspaces: signal subspace and noise subspace. The full polarimetric information can be derived from the signal subspace. For this method, the knowledge of the statistical property of polarimetric data is unnecessary. And above all, full polarimetric information can be preserved after filtering, which is very helpful for SAR image segmentation. However, the mean value can not be maintained well with this
* This project was supported by the National Natural Science Foundation of China (40571099) and the Research Fund for the Doctoral Program of Higher Education of China.
1122
Yang Jian, Deng Qiming, Huangfu Yue & Zhang Weijie
method, so the image after filtering becomes slightly darker than the original one. In this paper, we provide a new polarimetric filter based on Novak’s PWF and Gu’s subspace filter. The basic idea is to construct a new optimization function and constraint conditions. Under the constraint that a combination of the polarimetric SAR images processed by the proposed method equals to the output of the PWF method, the Euclidean distance between an unfiltered parameter vector and a signal space vector is minimized. This algorithm overcomes the disadvantages of both subspace and PWF methods. The full polarimetric information of an SAR image is preserved while noise is removed, hence various scattering characteristics and detail structures are kept. It is also shown that without constraint, our method is equivalent to the subspace filter. We use L-band, fourlook polarimetric SAR data from the NASA/JPL airborne to demonstrate the effectiveness of the proposed method.
2
Polarimetric whitening filte
Denote the complex scattering matrix S of a target in a linear horizontal (H) and vertical (V ) basis as SHH SHV S= (1) SV H SV V For the reciprocal case, holds. The three elements SHH , SHV and SV V can be used to define a complex ⎛ ⎞ vector SHH ⎜ ⎟ ⎟ (2) X=⎜ ⎝ SHV ⎠ SV V For a homogeneous scene, the vector X is joint complex Gaussian-distributed with a covariance matrix C = E{XX H}. Here E{} denotes the expectation operator and the superscript H stands for complex conjugate transpose. The general form of the covariance matrix C can be⎛written as ⎞ √ 1 0 ρ γ ⎟ ⎜ C = σHH ⎜ (3) ε 0 ⎟ ⎠ ⎝ 0 √ ρ∗ γ 0 γ where σHH = E{|SHH |2 } ε=
E{|SHV |2 } E{|SHH |2 }
E{|SV V |2 } E{|SHH |2 } E{SHH · SV∗ V } ρ=
(4) E{|SHH |2 }E{|SV V |2 } where the superscript ∗ indicates complex conjugate. In order to process the three complex measurements SHH , SHV and SV V into an intensity image having the least speckle, Novak and Burl[3] constructed an image from the quadratic γ=
y = X H AX
(5)
where A is a weighting matrix that is assumed to be Hermitian symmetric and positive definite. Then the ratio of the standard deviation of y to the mean of y is taken as a performance measure and is minimized to remove speckles. Through a series of mathematic deduction, the matrix A must be one that yields three equal eigenvalues for the matrix C ·A. Therefore, A = C −1 is a whitening filter called PWF, and the minimum-speckle image y can be obtained. For the multilook case, let Y be the N -look polarimetric covariance matrix defined as N 1 Y = Xi XiH (6) N i=1 where the subscript i means the ith look’s samples. Liu[5] demonstrated that the minimum-speckle image for multilook polarimetric SAR data can be constructed from the combination of the multilook polarimetric covariance matrix, which is expressed as follows Y22 Y11 + + σHH (1 − |ρ|2 ) σHH ε Y33 2Re (Y13 ρ∗ ) − √ 2 σHH (1 − |ρ| )γ σHH (1 − |ρ|2 ) γ yml =
(7)
where Yij is the ith row and jth column element of Y , and the subscript ml means multilook. From (7), one can find that the filtered image is an optimal combination of the elements of Y . Since all elements of polarimetric covariance matrix are linear combinations of the elements of Mueller matrix, we can also express the filtered image as follows ⎛ ⎞ m11 ⎜ ⎟ ⎜ m22 ⎟ ⎜ ⎟ ⎟ (8) yml = wT ⎜ ⎜ m12 ⎟ ⎜ ⎟ ⎝ m33 ⎠ m34
Polarimetric whitening filter for POLSAR image based on subspace decomposition where mij is the ith row and jth column element of the Mueller matrix and w is the weighted vector shown as below ⎛
γ+1 Re(ρ) + + ⎜ 2γ(1 − |ρ|2 ) (1 − |ρ|2 )√γ ⎜ ⎜ ⎜ 1 1 ⎛ ⎞ ⎜ ⎜ 2 − w1 ⎜ (1 − |ρ| ) γ(1 − |ρ|2 ) ⎜ ⎟ ⎜ ⎜w ⎟ ⎜ ⎜ 2⎟ ⎜ γ+1 Re(ρ) ⎜ ⎟ ⎜ ⎟ ⎜ − − w=⎜ ⎜ w3 ⎟ = ⎜ 2γ(1 − |ρ|2 ) (1 − |ρ|2 )√γ ⎜ ⎟ ⎜ ⎜ w4 ⎟ ⎜ ⎝ ⎠ ⎜ ⎜ 2Re(ρ) ⎜ − w5 √ ⎜ (1 − |ρ|2 ) γ ⎜ ⎜ ⎜ ⎝ 2Im(ρ) 2 √ (1 − |ρ| ) γ
⎞
1 2ε ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ 1 ⎟ ⎟ 2ε ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (9)
where ε, γ, ρ are defined as (4) and
m11 + m22 σHH = E + m12 2
m11 + m22 − m12 σV V = E 2
m11 − m22 σHV = E 2
m33 − m44 ∗ + i · m34 E{SHH · SV V } = E 2
3
(10)
Parameter verctors and parameter space
From the Mueller matrix, one can construct a parameter vector T p = (p1 , p2 , . . . , pn ) (11) There are many choices of the parameters pi such as the elements of a scattering matrix or a covariance matrix, the polarimetric entropy, the similarity parameters and so on. In this article, we define the parameters as shown in Table 1. Here, pi is the normalized element of a Mueller matrix and the normalizing factor F is the root square of the span F =
√
2m11
(12)
The parameter vector lies in an n-dimension linear space which is called the parameter space. Inversely,
Table 1
1123
Definition of parameters
Parameter
Definition
p1
m11 / F
p2
m12 / F
p3
m13 / F
p4
m14 / F
p5
m22 / F
p6
m23 / F
p7
m24 / F
p8
m33 / F
p9
m34 / F
one may obtain a Mueller matrix if the parameter vector is known. Therefore, denoising of a Mueller matrix can be carried through the filtering process of a parameter vector in the n-dimension linear space. Under the ideal condition where noise is not in existence, if there are m (m < n) different scattering components in the moving window, there should be m independent directions in the parameter space, and a parameter vector corresponding to a pixel can be treated as a combination of these typical directions. However, the speckle effect would cause a slight bias on that parameter vector. This biased parameter vector can be decomposed as two orthogonal parts—the signal vectors corresponding to m large eigenvalues, and the noise vectors corresponding to the other (n − m) small ones. Consequently, the n-dimensional linear parameter space can be divided into two orthogonal subspaces—the signal space spanned by m eigenvectors with large norms, and the noise space spanned by the others with small norms. For a parameter vector p, with a mean vector as μp , its covariance matrix is T Cp = E (p − μp ) (p − μp ) (13) which can be decomposed as Cp = QΛQT
(14)
where Λ = diag(λ1 , λ2 , . . . , λn ) and λi (i = 1, 2, . . . , n) are the eigenvalues of Cp . Obviously, Cp is positive definite, so the eigenvalues are real. Arrange λi in decreasing order, as shown by λ1 λ2 . . . λn
(15)
1124
Yang Jian, Deng Qiming, Huangfu Yue & Zhang Weijie
The columns of the n-by-n matrix Q in (14) consist of n orthogonal eigenvectors q1 , q2 , . . . , qn of Cp , that is Q = (q1 , q2 , . . . , qn )
(16)
Assuming there are m large eigenvalues λ1 , λ2 , . . . , λm , (n−m) small ones, where λm λm+1 , then Q can be rewritten as Q = (QS |QN )
(17)
where the n-by-m matrix QS is composed of the eigenvectors q1 , q2 , . . . , qm , corresponding to the m larger eigenvalues, and the n-by-(n − m) matrix QN is composed of the other (n − m) eigenvectors qm+1 , qm+2 , . . . , qn , corresponding to the other (n−m) smaller eigenvalues. Practically, one can accumulate the eigenvalues until their summation reaches a certain threshold, and then make the number of eigenvalues summed as the dimension of the signal subspace. This threshold may be defined as some percentage of the summation of all the eigenvalues, e.g., 80%.
4
Refined subspace filtering algorithm
Denote the unfiltered and filtered parameter vector with p and pˆ respectively. Because polarimetric information is mainly contained in the signal subspace, pˆ should be reconstructed from the sum of the mean vector and the weighted combination of the eigenvectors in the signal subspace. Define the weighted coefficient vector as x = (x1 , x2 , . . . , xm )T , then pˆ can be expressed as pˆ =
m
xi qi + μp = QS x + μp
a minimum-speckle intensity image. Therefore, if we add such a constraint that a weighted combination of pˆ’s elements with the same form as expression (8) equals to the output of the PWF, the disadvantages of the subspace filter and the PWF can both be compensated and promising results can be expected. Mathematically, all this can be expressed as the following constrained optimization problem 2
(19a)
s.t. wT pˆ = yml
(19b)
where 2 is the 2-norm of a vector, w is the weighted vector for pˆ, and yml is the filtering result of PWF divided by the normalizing factor. Noticing from (18) that pˆ can be expressed by means of x, we may rewrite (19a) and (19b) as follows min QS x + μp − p22
(20a)
s.t. k T x + c = yml
(20b)
T where k is QT S w and c is w μp . It is easy to find that if the constraint of (20b) is omitted, the optimization function (20a) would reach its minimum when
x = QT S (p − μp )
(21)
It is the same as that in (18) implying that Gu’s subspace filter[6] is gotten for no-constraint case. When (20b) is added, the whole problem is to seek the minimum value of a quadratic form under linear constraints. Using the Lagrange multipliers, we derive the solution as
(18)
x = QT S (p − μp ) −
i=1
where qi (1 i m) are the eigenvectors belonging to the signal subspace, μp is the mean vector of the parameter vector. Generally speaking, filtered SAR data should be close to the original in the local area, which means the 2-norm of the vector pˆ − p is to be minimized. It should be noted that because the noise part is discarded in the reconstructed vector p, ˆ the overall intensity of the image after filtering is decreased a little. On the other hand, the PWF method mentioned above is the optimal polarimetric filter in terms of providing
min ˆ p − p2
5
T k T QT S (p − μp ) k − (yml − c)k T k k (22)
Experimental results
To evaluate the proposed algorithm, we use a set of four-look polarimetric data of the San Francisco Bay area measured by the NASA/JPL airborne L-band polarimetric SAR. Fig. 1(a) shows the original span image, and Fig. 1(b) shows the resulting image of the subspace filter. Fig. 1(c) and Fig. 1(d) illustrate the minimum-speckle intensity and span images of applying the proposed method, respectively. For comparison, we apply 7 by 7 moving windows to both of
Polarimetric whitening filter for POLSAR image based on subspace decomposition
Fig. 1
1125
Original span and filtering results
the filters. It can be seen that Fig. 1(c) has a higher contrast, where detail structures and edges are well preserved, especially the point target in the sea (as shown in the blocked region). Moreover, Fig. 1(d) has an overall intensity more consistent with Fig. 1(a), while Fig. 1(b) appears a little darker. For quantitatively comparing the speckle reducing performance, we calculated the mean values of three different areas in the original image, the resulting image from the subspace filter and from the proposed method, respectively. The selected three areas contain: urban, grass and ocean. The corresponding results in the testing areas are shown in Table 2. Ob-
viously, the mean values of the proposed method are closer to the original ones for all areas. An important application of speckle filtering is to enhance the separation of all polarized channels for terrain classification. In order to demonstrate the effectiveness of the proposed method, we calculated the polarimetric entropy and alpha angle[7] for the original data, and the filtered data by the propose method, respectively. The corresponding results are plotted in the H/α plane as shown in Fig. 2. Obviously, the entropy extracted from the data processed by the proposed method is much more separable for different classes.
1126
Yang Jian, Deng Qiming, Huangfu Yue & Zhang Weijie
Table 2
Experimental results of subspace filters and proposed method
Mean value Urban Grass Ocean Overall
Original span 0.537 0.169 0.019 0.283
0 8 5 6
Span of subspace 0.401 0.150 0.017 0.209
3 7 5 7
Span of our method 0.535 0.170 0.019 0.286
4 0 6 7
After filtering, the full polarimetric information of an SAR image is preserved, hence scattering characteristics and edges or other detail structures are kept. The experimental results with the NASA/JPL airborne polarimetric SAR data demonstrated that the proposed method is more helpful for terrain classification. However, it should be pointed out that the calculation cost of this method is larger than the others.
References [1] Lee J S. Speckle analysis and smoothing of synthetic aperture radar images. Computer Graphics and Image Processing, 1981, 17(1): 24–32. [2] Lee J S, Grunes M R, De Grandi G. Polarimetric SAR speckle filtering and its impact on classification.
IEEE
Trans. on Geoscience and Remote Sensing, 1999, 37(5): 2363–2373. [3] Novak L M, Burl M C. Optimal speckle reduction in polarimetric SAR imagery. IEEE Trans. on Aerospace Electronic System, 1990, 26(2): 293–305. [4] Novak L M, Burl M C, Irving W W. Optimal polarimetric processing for enhanced target detection. IEEE Trans. on Aerospace Electronic System, 1993, 29(1): 234–243. [5] Liu G, Huang S, Torre A, et al. The multilook polarimetric whitening filter (MPWF) for intensity speckle reduction in polarimetric SAR images. IEEE Trans. on Geoscience and Remote Sensing, 1998, 36(3): 1016–1020. [6] Gu J, Yang J, Zhang H, et al. Speckle filtering in polarimetric SAR data based on the subspace decomposition. IEEE Trans. on Geoscience and Remote Sensing, 2004, 42(8): 1635–1641. Fig. 2
6
Scatter plots for alpha and entropy
Conclusion
A new method for speckle filtering of polarimetric SAR data has been proposed. By minimizing the Euclidean distance between an unfiltered parameter vector and a signal space vector under the constraint that the optimal combination of polarimetric images equals to the PWF output, speckles are reduced.
[7] Cloude S R, Pottier E. An entropy based classification scheme for land applications of polarimetric SAR. IEEE Trans. on Geoscience and Remote Sensing, 1997, 35(1): 68–78.
Yang Jian was born in 1965. He is a Ph. D. and a professor. His research interests include microwave remote sensing, PolSAR, and InSAR. E-mail:
[email protected].