VIEWS: 8 PAGES: 6 POSTED ON: 12/11/2011 Public Domain
A Wavelet Multiresolution Technique for Texture Analysis and Segmentation of SAR Images G. De Grandi (1), J. Fortuny(1), J.S. Lee(2) , D.L. Schuler(2) (1) European Commission – Joint Research Centre 21027 Ispra (VA), Italy Email:frank.de-grandi@jrc.it (2) Naval Research Laboratory Remote Sensing Division Code 7263 Washington, DC 20375-5351, USA Email:lee@ccf.nrl.navy.mil Abstract A technique for multi-scale texture analysis and segmentation of SAR images is presented. Textural features are extracted using a decomposition based on a wavelet frame. The feature vector is composed of local estimates of normalized variances of the smooth image and of the wavelet coefficients. The technique is also extended to the polarimetric case. To the purpose polarimetric synthesis is performed in the wavelet domain. In the synthesis, selected antenna configurations that are known to be optimal with respect to texture diversity are used to construct the feature vector. Controlled experiments, based on Monte Carlo simulations, are set up to assess the performance of the technique. Finally, examples are given of the application of the technique to real SAR data (ERS-1 and SIR-C) related to tropical swamp forests mapping and oil slick detection. Keywords: texture, segmentation, wavelet frame, multi-scale representation, SAR, polarimetry. 1. INTRODUCTION The notion of texture is well defined at the level of our visual perception. The human eyes and brain are able to delineate areas that exhibit common spatial patterns. However this notion which is fairly simple fro m the intuitive point of view, eluded until now all efforts to cast it into a universal mathematical definition. Despite this, several approaches in image processing have been devised that work reasonably well in certain specific areas of application. Narrowing our scope to hit a concrete target, we concern ourselves here with textures that are somehow associated with the “roughness” of the – eventually multi-dimensional - surface that describes the backscatter as a function of space in a radar image. From the statistical point of view this surface is considered a random process in two dimensions and analyzed using local estimates of a second-order one- point statistics (variance) and of a second-order two-points statistics, the structure function (field of increments of a random process over a distance). The roughness associated with a feature of interest (e.g. a natural target such as a forest) develops at some specific scale multiple of the resolution element, a fact that suggests that measures must be taken at several scales to characterize and discriminate different targets. Questions are then: how do we measure such textures, how do we take scale into account, how do we detect structures or segment homogeneous regions provided by these measures? Additional questions related in particular to Synthetic Aperture Radar (SAR) imagery are: how do we take speckle into account and what is the role of polarimetry in texture segmentation? We will try here to give one possible answer to those questions. 2. THEORY 2.1. Texture Measures Our approach is heavily inspired by by a texture segmentation and feature reduction technique presented in [1][2]. The underpinning mechanism for the extraction of texture measures is a multi-resolution image representation provided by a wavelet frame that acts as a differential operator [3]. Text ural features are captured by spatially local estimates of: i) variance of the smooth signal, ii) variance of the wavelet coefficients. Intuitively, these measures are sensitive to sharp transitions (contours of areas with different reflectivity, point targets), and to edge density, smoothness and swing (at a given scale) within extended targets. In more detail a feature vector holding the texture measures is constructed as follows. We will consider at first the case of scalar radar. The wavelet frame decomposition is imp lemented numerically using an “á trous” algorithm and a filter bank structure. The low-pass and high-pass filter coefficients are given in [3]. Let S j be the image at the output of the low-pass filter (smooth image), Wjx and Wjy the wavelet coefficients in the image row and column directions (output of high pass filters) at scale s = 2 . We assume the convention that the j index j increases with scale. The following steps are performed: i) normalization of the wavelet coefficients by the smooth signal at the same scale to cancel speckle influence; ii) introduction of non-linearity; iii) smoothing by convolution with a separable discrete cubic B-spline kernel enlarged by a factor m . The support of the discrete cubic B-spline is 4. Therefore averaging will take place in a neighborhood of m × 4 pixels. The smoothing window is fixed, and therefore in this respect the approach is not multi-scale. However, within the averaging window measures of textures that occur at different dyadic scales are taken into account by the frame decomposition. In this sense the approach is multi-scale. In summary, for each of the two scales in the decomposition the feature vector will be composed of: (S ) p j−1 2 Wjx 2 W jy 2 v1j = v2 j = v3 j = (1) 2 Sj Sj S p−1 j The same pixel size (equal to the original image pixel size) is kept for all scales and features, resulting in over-sampling. 2.2 Feature Reduction Feature reduction is obtained by an approximate solution of the multiple discriminant transform (MDT) as proposed in [2]. This estimate of the MDT solution does not require a-priori knowledge about the regions (e.g. the within-region and between region scatter matrices) and relies on the availability of local variance estimates at two scales. On the other hand, for this reason, the reduced features are only available at scale 2m. For trials where the characteristic scale of texture is too small to be handled at resolution 2m using the approximate MDT, a Karhunen-Loève (KLT) transform is used instead. 2.3 Extension to the Polarimetric Case We are also interested in dealing with texture that arises from a mixture of polarimetrically different scattering mechanisms [4]. The polarimetric information is included by augmenting the feature vector with local estimates derived from the wavelet decomposition of the backscattered intensity synthesized at a number of antenna configurations. An analytical tool presented in [4] - the polarimetric texture signature – indicates that optimal polarization states for enhancing differences in polarimetric text ure are the linear cross-polarized states with orientation angle ψ = 0 (HV) and ψ = 45 . The latter o o configuration will henceforth be referred to as xpol45. Accordingly the feature vector is constructed by stacking the local estimates at two scales as in (1) derived from the synthesized power at the two optimal polarization states. We have also experimented with a second alternative, where texture measures derived from the Cloude- Pottier scattering mechanisms decomposition - in particular entropy - are used as entries in the feature vector. The frame decomposition of the synthesized power images at different polarization states is performed using the following technique. The polarimetric power synthesis and the wavelet fra me operators commute. I = T(ψ T , χ T , ψ R , χ R ) ⋅ M PW ⋅ I = PW ⋅ T(ψ , χ) ⋅ M = T( ψ, χ) ⋅ PW ⋅ M (2) where I is the synthesized backscattered power, PW is the frame operator, M is the Stokes operator and T is the linear transformation holding information on the transmit and receive antenna configurations (Stokes vectors). Therefore the wavelet transform of the backscattered intensity at any polarization state can be derived from the wavelet transform of the Stokes operator, or, in other words, polarimetric synthesis is performed in the wavelet domain. The advantage of this approach is that the wavelet transform, which is the computationally heavy part, needs to be performed only once. Once the wavelet coefficients for the Stokes operator have been computed then the coefficients of the backscattered intensity for all desired polarization states can be derived by applying the synthesis operator, which is computationally quite lighter. 3. CONTROLLED EXPERIMENTS USING SIMULATED IMAGES A number of controlled experiments based on simulated polarimetric SAR imagery were performed. These experiments allow us to explore the potential of the technique with respect to different detection and segmentation problems and to assess the performance of alternatives in the processing method. In this short communication, we report next about one case concerning the simulation of a fragmented forest. 3.1 Fragmented Forest In this experiment we simulate the polarimetric scattering from a natural extended target – a fragmented forest. In this case we assume that spatial modulation of the radar return is caused by openings in a dense canopy where the electromagnetic wave can penetrate down to a bushy understory. The model consists therefore of two vegetation classes – referenced as forest and understory - with different structural parameters. Scattering at P-band for each class is characterized by a polarimetric covariance matrix that is computed using the University of Texas at Arlington (UTA) wave scattering model. Details on the model parameters can be found in [4]. We expect to generate polarimetric texture by a mixture of the returns from the two classes because different scattering mechanisms come into play. For instance, for the forest class the VV backscatter is the result of both direct scattering from the canopy and indirect scattering branch-soil, while at HH polarization the dominating mechanism is indirect branch-soil scattering. This creates polarimertric diversity between the two channels. For the understory class, direct and indirect scattering from the branches and direct scattering from the surface affect the HH and VV channel equally. a) b) c) f) d) e) Fig. 1. Simulated fragmented forest data set. Mask used to define regions holding a radiometrically “pure” class (closed canopy forest) and a mixture of two classes (fragmented forest) (a). HV amplitude image (b). Feature vector constructed from the synthesized xpol45 image: Variance of smooth image at scale 22 (c). Textural segmentation using K -means clustering with 5 categories. Results when clustering is applied to the first KLT transform in the AMDT procedure: Clusters belonging to high textures related to fragmented forest (grey to black) (d). Clusters belonging to low textures related to closed-canopy forest (grey to white) (e). Segmentation when clustering is applied using the full AMDT procedure (f). A simulated one-look scattering matrix image for this trial is constructed in the following way. Random deviates of scattering matrix samples, corresponding to the covariance matrices computed by the UTA model for the two vegetation classes, are generated by a Monte Carlo technique. A square array of these samples is partitioned into two sub-sets by a curvilinear mask. The first subset is filled with scattering matrix samples belonging to a pure class (either forest or understory). The second set is filled with a mixture of samples belonging to the two classes. The mixture is laid out spatially simply alternating pixels of each class. For the sake of brevity we report here only the case where the pure class is forest. Thus the segmentation problem concerns the ability to discriminate pristine forest from fragmented forest. Since analysis using the polarimetric texture signature indicates that the xpol45 configuration should give maximum of the normalized variance of backscatter for this type of polarimetric mixture [4], while the HV configuration should give a minimum, we build feature vectors corresponding to these polarization states. (a) The mask used for mixing the two classes in the simulated image is shown in Fig. 1 and the HV backscatter image is reported in Fig. 1(b). This image suggests the difficulty of the segmentation problem, in particular if only one polarization were available and only backscatter were used in the classification. In this experiment we use 2 scales and a smoothing kernel with m = 4 . The feature vector component 2 related to the variance of the smooth image at scale 2 is shown in Fig. 1(c). Both the wavelet coefficients and the variance of the smooth image feature a significant between-regions variance. At the coarser scale the smooth image gives mostly information on the edge between regions, a fact which points at the role of scale in texture analysis. A segmentation trial has been set up, using the approximated multiple discriminant transform (AMDT) for feature reduction (see section 2.2) and a K-means clustering procedure with 5 classes. At first we analyze the results when clustering is applied to the output of the first KLT transform (see Fig. 1(de)). The fragmented forest region is well segmented by the categories corresponding to high textures. Pixels corresponding to these clusters are shown in Fig. 1(d ), while pixels belonging to the low-texture categories are shown in Fig. 1(e). There are some commissions that are tantamount to 4% of the class population. Omissions are 2.8 %. Segmentation results for the same 5 categories problem but when the full AMDT procedure is used are shown in Fig. 1(f). Here the segmentation accuracy is much higher (no omission or commis sion errors) at the expense of loss in resolution due to the additional smoothing at a scale factor 2 m which is a requisite for obtaining an approximate MDT solution. 4. APPLICATIONS USIN G REAL SAR DATA The question whether texture can play a significant role in thematic information extraction from SAR data is on the table of the SAR community since a long time. Indeed most of the applications reported in the literature rely on classifications algorithms that assume a piece-wise continuous model of the underlying radar reflectivity. However some natural targets do not exhibit radiometrically “pure” backscatter, and therefore difficulties arise in mapping these targets into a pure category unless texture is taken specifically into account. A case in point is provided by the closed canopy dense tropical forests in Central Africa. A second case, in an entirely different application area, is provided by surfactant slicks in the ocean. In the following we illustrate some results related to those two cases, thus hoping to contribute an answer to the question about the relevance of texture. 4.1 Swamp Forest Mapping in the Congo River Floodplain In tropical Africa and specifically in the Congo River floodplain, the swamp forest upper canopy appears as a smooth surface when imaged by a radar with a resolution of some 30 m. On the other hand the primary lowland rain forest’s canopy is much more inhomogeneous, being composed by many and diverse competing species. Therefore it appears as a more textured surface at the resolution of the radar. This effect is enhanced in imagery acquired by a C-band radar, which, due to lower penetration, senses only a shallow slice of the upper canopy. This difference in backscatter statistical properties can be exploited to discriminate in a very effective way the two forest classes based on texture measures. The physics behind this effect is explained in [5]. In Fig. 2(a) an ESA ERS backscatter image of an area in the Congo floodplain is shown. The image was resized by pixel aggregation to 100 m pixel spacing. The approximate multi discriminant transform reduces the feature vector to one component, which is shown in Fig. 2(b). Notice how in the backscatter image the differences in smoothness between the two types of forest are barely vis ible, whereas the texture image captures neatly these differences (see dark ribbon along the river, and the fern-like patterns spreading from the river bed into the forest domain).The ratio of between- to within-regions variances in the texture image is high enough to allow us to use this image as input product for the segmentation stage. In our trial a K-Means algorithm with four clusters produces the result shown in Fig. 2 (c). Notice how the swamp forest fine structures along the water drainage patterns are well preserved in the map. The swamp forests of the Congo basin were previously mapped in connection with the Global Land Cover Map 2000, using texture features based on block window estimates of the K distribution order parameter. Comparison of the relative accuracy of the two approaches is under way. 4.2 Oil Slicks Detection It is well known that detection of oil slicks based on one-point statistics of SAR backscatter values only is marred by the fact that several other features, such as waves, ship wakes, floor topography, present similar range of values. To achieve reliable detection the classification algorithm must usually be compounded with contextual information (e.g. two-point statistics, morphological operators) and/or take into account ancillary data, such as ship routes and wind data. a) b) c) Fig. 2. Swamp forest mapping in the Congo basin. The original ERS -1 backscatter image (down -sampled to 100 m) (a). The texture image at the output of the approximate multi discriminant transform (b). Segmentation into 4 clusters using K-Means clustering (c). We want here to explore the potential of a simpler approach based on the enrichment of the feature vector by texture measures and polarimetric target decomposition. For the experiment, a fully polarimetric SIR-C image of the English Channel, C-band, is used (courtesy JRC-EMSL). The synthesized HH amplitude image is shown in Fig. 3 (a). We can identify by visual inspection a fairly big oil slick (highlighted area A), features associated with waves (highlighted area B), and several strong point targets related to ships (highlighted area C). In a first trial the feature vector was constructed using the xpol45 configuration. More precisely, the original covariance matrix data (multi-look complex) of the JPL distribution kit were averaged by 4x4 aggregation, and converted to the Stokes operator representation. The wavelet frame operator was then applied to the Stokes data (see section 2.3). Finally by polarimetric synthesis the xpol45 amplitude image was generated. The pixel spacing is 50 m. In the feature vector construction a smoothing kernel dilated by m = 4 is used. The variance of the wavelet coefficients along the image row direction W2x at scale 2 2 is shown in Fig. 3(b) and the component along the column direction W2y in Fig. 3(c). The point targets are very well captured, and the contour of the slick is clearly identified. However if we look at W2y we see that other linear structures, like the wave front corresponding to area B, give similar response. Therefore identification of the slick would require a second level of pattern recognition. As a second trial we added to the feature vector texture meausures derived from the Cloude-Pottier decomposition based on the coherency matrix. In particular we us the entropy image. It turns out that entropy increases significantly in the presence of surfactant slicks due to dampening of resonant Bragg scattering that dominates in clear waters [6]. The wavelet coefficients related to entropy allow now for the slick detection with a good signal to nois e ratio (SNR = 78.4 dB) and can be captured unequivocally by simple image thresholding (Fig. 4(a )). The other side of the coin is that entropy gives a worse signal to noise for the point scatterers. The principal component of the composite feature vector (entropy and the xpol 45 intensity) achieves optimal detection of both features (see Fig. 4(b)) with SNR=31.1 dB for the edges related to the oil slicks. Results of this experiment seem to indicate that the combination of polarimetry, coherent decompositions and texture measures provide a robust vehicle for the oil-slick detection problem. This qualitative statement will have to be of course quantified by a detection performance analysis based on a rigorous validation exercise. 5. CONCLUSIONS We have presented a simple but effective technique for textural segmentation of SAR imagery. The technique exploits the ability of a wavelet frame to capture “correlation“ properties of pixels in neighborhoods of different sizes, the aggregation of features extracted from SAR backscatter at different polarizations, and a mechanism to represent the feature vector in a space that approximately optimizes the between-regions variance to the within-regions variance. The merit of the approach lies on: i) simplicity of the basic mechanism for constructing texture measures, which does not rest on a specific image model for the underlying radar reflectivity; ii) efficient algorithms for the numerical imple mentation; iii) ease of use (only a few free parameters to be set). The combination of these aspects make the approach flexible enough to cope with different thematic applications and computationally suitable also for large scale real- life radar mapping applications. B A C a) b) c) Fig. 3. SIR -C C-band multi-look HH amplitude image acquired over the English Channel (a) (courtesy JRC-EMSL). Variance of the wavelet coefficients along the image row direction, estimated using the xpol45 and a smoothing factor m=4 (b). Variance of the wavelet coefficients in the image column direction (c). a) b) Fig. 4. Segmentation by image thresholding of the feature vector component holding the normalized variance of entropy (a). First component if the KLT transform of the feature vector constructed using both the xpol45 synthesized amplitude and the entropy (b). 6. REFERENCES [1] M. Unser, “Texture classification and segmentation using wavelet frames,” IEEE Trans. Image Processing, vol. 4, no. 11, pp. 1549-1560, Nov. 1995. [2] M. Unser and M. Eden, “Multiresolution feature extraction and selection for texture segmentation”, IEEE Trans. Image Processing, vol. 11, no. 7, pp. 717-728, July 1989. [3] S. Mallat and S. Zhong, “Characterization of signals from multiscale edges,” IEEE Trans. Pattern Anal. Machine Intell., vol. 14, no. 7, pp. 710-732, 1992. [4] G.F. De Grandi, J-S. Lee, D. Schuler, and E. Nezry, “Texture and speckle statistics in polarimetric SAR synthesized images,” IEEE Trans. Geosci. And Remote Sensing, vol. 41, no. 9, pp. 2070-2088, Sept. 2003. [5] M.L. Williams, S. Quegan, D. Blacknell, “Intensity statistics in the distorted Born approximation: application to C-band images of woodland”, Waves in Random Media, vol. 7, no. 4, 1997. [6] D.L. Schuler, J.S. Lee, and G. De Grandi, “Spiral Eddy Detection using Surfactant Slick Patterns and Polarimetric SAR Image Decomposition Techniques”, Proc. Geoscience and Remote Sensing Symposium ,IGARSS'04, IEEE International, vol 1, 20-24 Sept. 2004.