Open access peer-reviewed chapter

Adaptive Speckle Filtering in Radar Imagery

By Edmond Nezry

Submitted: May 6th 2013Reviewed: April 11th 2014Published: June 11th 2014

DOI: 10.5772/58593

Downloaded: 4264

1. Introduction

Historically, foundations in the domain of radar speckle properties have been laid down from the 1940’s to the 1980’s. Decisive theoretical advances were made by teams led by Professor Fawwaz Ulaby at the University of the Michigan (USA), by Professor Christopher Oliver at the Defence Research Agency in Great Malvern (UK), and by Professor Keith Raney at the Canadian Centre of Remote Sensing (Canada) since the 1970’s. Then, the domain of speckle filtering in SAR images matured in the period 1976-2000, mostly under the impulsion of Dr. Jong Sen Lee of the Naval Research Center, Washington D.C. (USA). Since 1986, the team led by Dr. Armand Lopès at the Centre d'Etude Spatiale des Rayonnements in Toulouse (France), has carried out and then inspired the development of the most efficient speckle filters existing today. Since 2000, with speckle filters having reached a satisfactory level of performance, no significant advances have been made. Nevertheless, in this period, the use of speckle filters in a wide range of applications using SAR imagery has become generalized.

A radar wave can be considered, with a good approximation, as plane, coherent and monochromatic. It is emitted by an antenna towards a target. The target backscatters partially the radar wave in the direction of a receiving antenna. In the vast majority of spaceborne Synthetic Aperture Radars (SAR), a single antenna assumes the two functions of emission and reception (monostatic radar).

The complete radar measurement is the combination of the horizontally (H) and vertically (V) linearly polarised radar waves, at emission and at reception after backscattering by the observed target. After signal calibration, this measurement, affected by noise, enables to restitute for each resolution cell a polarimetric backscattering matrix S:


whose coefficients Spq=|Spq|.exp(jpq)are the complex backscattered field amplitudes for p-transmit and q-received wave polarisations. They relate the backscattered wave vector Edto the incident wave vector Ei: Ed = exp(j.k.r).S. Ei. This complete radar measurement is said "fully polarimetric", often abbreviated in "polarimetric".

Let consider the interaction of the radar wave with an extended random surface, i.e.a surface containing a sufficiently great number of scatterers per resolution cell, with no preponderant scatterer with respect to the others. For whatever configuration pqof polarisation, the signal Spqreceived from such a surface by the antenna becomes, after quadratic detection in intensity Ipq:

Ipq = |Spq|2 = Spq.Spq*(where*denotes complex conjugation)E2

This detected signal intensity Iis proportional in average to the radar "backscattering coefficient" σ°. The backscattering coefficient σ°=4.π.|Spq|2is the average radar cross-section per surface unit [1]. σ°, expressed in m2/m2, is a dimensionless physical quantity. It is a physical property of the sensed surface, which depends principally on its roughness, its dielectric properties, its geometry, and the arrangement of its individual scatterers.

Carrying the radiometric information with regard to the sensed target, σ° is a function of the frequency of the radar wave, of its angle of incidence upon the target, and of the configuration of polarisation. In terms of physical meaning, the radar backscattering coefficient is analogous to the bidirectional reflectance in the domain of optical wavelengths: σ° # 4cos2θ,where is the incidence angle of illumination on the sensed target.

Nevertheless, detected radar images look visually very noisy, exhibiting a very characteristic salt-and-pepper appearance with strong tonal variations from a pixel to the next. Indeed, since radar imaging systems are time-coherent, radar measurements over random rough surfaces are corrupted by "speckle" noise due to the random modulation of waves reflected by the elementary scatterers randomly located in the resolution cell. Then, coherent summation of the phases of elementary scatterers within the resolution cell results in a random phase of the complex pixel value.

This speckle "noise" makes both photo-interpretation and the estimation of σ° extremely difficult. Actually, speckle is a physical phenomenon, which is inherent to all coherent imaging systems (radar, lidar, sonar, echography). In most remote sensing applications using radar/SAR imagery, speckle is generally considered a very strong noise that must be energically filtered to obtain an image on which classic and proven information extraction techniques could be further applied, in particular the techniques used for optical imagery acquired in the visible and near-infrared part of the electromagnetic spectrum.

Therefore, speckle filtering and radar reflectivity restoration are among the main fields of interest in radar images processing for remote sensing. Speckle filtering is a pre-processing aiming at the restoration of σ° value in the first place. This pre-processing must account for both the particular properties of the speckle, and those of extended imaged targets (often called "clutter"). It must also account for the radar imaging system that has sensed the target and for the processor that has generated the images. For stationary targets of infinite size, speckle filtering is equivalent to a simple smoothing using a moving processing (averaging) window. An ideal filter must nevertheless avoid image degradation through excessive smoothing of the signal. To this end, it must respect structural image information (road and water networks, etc.), and the contours of radiometric entities. In addition, it must also respect the particular texture of the clutter, in forested or urban environments for example. Last, it must also identify the strong local radiometric variations due to the presence of strong scatterers (often artificial in nature) from those due to spatially rapid speckle fluctuations.

Therefore, an ideal speckle filter must satisfy to the following specifications:

  1. Preserve accurately the local mean value of the radar reflectivity (i.e.the quantity actually measured by the radar, which is proportional to σ°) to enable, for example, the comparison of radar reflectivities in the framework of a multitemporal analysis of radar acquisition series.

  2. Smooth as much as possible homogeneous image areas and therefore reduce the speckle to increase the Equivalent Number of Looks (ENL) of the radar image (cf.§ 2.1.4). The minimum ENL depends on the desired radiometric accuracy. For example, a 1 dB accuracy with 90% confidence level (i.e.less than 25% variation of the radar intensity) requires an ENL value around 230, and a 2 dB accuracy with 90% confidence level (i.e.less than 60% variation of the radar intensity) needs an ENL value around 40.

  3. Preserve texture as much as possible where it exists in the image (forests, non-homogeneous fields, etc.) to avoid confusions among radiometrically similar areas exhibiting different texture. Therefore, a speckle filter must be able to discriminate heterogeneity effects due to texture from those due to speckle.

  4. Both preserve and denoise image structures (contours, lines) as well as the quasi-deterministic responses due to corner reflector effects within strongly textured areas such as urban environments. Indeed, the energy of artificial radar reflectors responses must be preserved to enable radiometric calibration, in particular when calibration targets are dispersed in the radar image.

  5. Minimise, and whenever possible prevent loss in useful spatial resolution during the speckle filtering process.


2. Statistical properties of speckle and texture in radar images

In this section, the statistical properties of speckle in images produced by coherent imaging systems such as imaging radars, lidars or sonars, are exposed. Since a good speckle filter must restore the texture of the scene imaged by the radar, the statistical properties of texture in radar images are examined as well. This analysis intentionally restrains to the first order statistical properties, since only these are generally used by the estimation techniques involved in speckle reduction methods. Explicit use of second order statistical properties of both the speckle and the imaged scene in the filtering process is adressed in Section 4.

2.1. Speckle in SAR images

First, let consider a natural target with radar backscattering coefficient σ° which remains constant from a resolution cell to the next. Such a target, said stationary, "homogeneous", or "textureless", may be found in a medium-resolution radar image within cultivated fields or forest cut areas for example. Spatial variations of the radar signal backscattered by such a target are therefore only due to speckle and the statistical properties of the radar signal are those of the speckle, which depends only on the radar system and the image production system.

2.1.1. First order statistics for "1-look" complex radar images

The radar imaging system is linear, spatially invariant, and can be characterised at each image pixel (x,r), where xis the azimuth and ris the radial distance with respect to the flight path of the sensor’s carrier, by a complex impulse response h(x,r), which defines the resolution cells. The 3dB widths of the impulse response, often used to quantify the spatial resolutions δ(x) in arimuth and δ(r) in radial distance.

If the number of individual scatterers within a resolution cell is large enough and none of them has absolute predominance in scattering the radar wave [2] [3], speckle can be modelled in output of the radar processor by a circular complex gaussian random process uc(x,r) at the corresponding image pixel of coordinates (x,r). The complex amplitude (complex radar signal, in a complex radar image) in output of the radar processor can be expressed as Ac(x,r)=a(x,r)+j.b(x,r), where a(x,r) is proportional to the field backscattered by the surface element located in (x,r). In the literature, the terms a(x,r) and b(x,r) are generally called i(x,r) ("in-phase" term) and q(x,r) ("in-quadrature" or "out-of-phase" term), respectively.

These complex data in output of the radar processor are called a "1-look complex" image or equivalently, a "Single-Look-Complex" (SLC) image.

For a homogeneous area where A(x,r)=A0(constant), the presence of speckle only, the complex amplitude has the form [4]:

Ac(x,r) = [A(x,r).uc(x,r)]*h(x,r) = A0. uc(x,r)*h(x,r) = A0. V(x,r)E3

where A=A0is the amplitude of the backscattered wave and * is the convolution operator. The random process uc(x,r) represents the circular white gaussian (by application of the central limit theorem over a large number of individual scatterers within the resolution cell) complex process responsible for the speckle. Then, V(x,r) is a correlated complex gaussian process characterising the "gaussian fully developed speckle", which behaves as multiplicative "noise" with a very high spatial frequency dependent on h(x,r).

Owing to the particular character of the speckle phenomenon in coherent imagery, which is the case of radar imagery, information extraction from a radar image results in a statistical estimation problem. It is therefore mandatory to have an as complete as possible statistical speckle model. The statistical model of the fully developed speckle has proven perfectly adapted to SAR images, at the spatial resolutions / wavelength frequencies combinations actually used in radar remote sensing. Note that this speckle model emphasises the importance of wave phase information under the condition that the phase has not been lost during signal detection.

Goodman hypotheses [2] [5] enable to calculate speckle statistics before ( a complex radar image) and after radar signal quadratic detection in amplitude Aor in intensity I. Taking this hypotheses into account, one obtains:

E[i]=E[q]=0           and           Var[i]=Var[q]=σi2 =σq2 =σ2E4

where E[.] denotes the mean expectation and Var[.]denotes the variance. Therefore, the mean backscattered intensity Ris expressed as:


R=E[I] is proportional to the radar backscattering coefficient σ°. Rcan be estimated by R≈<I>, where <.> denotes the averaging operator applied in a neighbourhood of the pixel under consideration. Then, from R, it is possible to retrieve the σ° value through the calibration parameters of the radar image. Problems related to radar image calibration are discussed in detail in [6] [7] [8]. R=E[I] can be estimated locally in a radar image through the spatial averaging <I> of the intensities of a number N of image pixels in a spatial neighbourhood of the image pixel of interest [9]. This estimate is the unbiased Maximum Likelihood estimate with minimal variance [6].

Besides, the in-phase and the out-of-phase components iand qof the complex radar signal are decorrelated, E[i.q]=0, and therefore independent of each other. It results of the above considerations that the probability density functions (pdf) of iand qare Gaussian distributions:

P(i) =1/2π.σ2. exp(i2/2σ2)     and    P(q) =1/2π.σ2. exp(q2/2σ2)E6

with a phase of the complex radar signal ϕ=Arctg(q/i), which is uniformly distributed in the interval [0, 2π].

2.1.2. First order statistics for "1-look" detected radar images

Since iand qare independent variables, their joint pdf is P(i,q)=P(i).P(q). The pdf Pu(I)of the speckle in intensity is obtained by a simple change of variable I=a2+b2in Equation (6). Then, Pu(I)results being an exponential pdf:

Pu(I) = (1/E[I]) . exp(I/E[I])        for      I>0E7

It is important to note that, since E[I] is proportional to the radar reflectivity R, Pu(I)is the pdf of Iconditional to the value of R, therefore Pu(I)=Pu(I/R).

The two first first-order moments of this pdf can be expressed as a function of its standard-deviation as follows:

mean   E[I] = 2σ2    and     standard-deviation   σI = 2σ2E8

To characterise the strength of speckle in radar images, it is convenient to consider the normalised coefficient of variation of the intensity, CI, i.e.the standard-deviation of the radar signal intensity Iover its mean value:

E[CI]  =  σI/E[I] =  1E9

The value of the coefficient of variation of speckle only, sometimes called "contrast index", is a constant for every type of radar image product. Equation (9) means that, in a radar image, the dispersion (variance) of radiometry increases as the mean signal backscattered by the target increases. This justifies in part the qualification of "multiplicative noise" given to the speckle.

Clearly, with a signal-to-noise ratio of 1, the radiometric degradation due to speckle makes very difficult the discrimination of two homogeneous targets with very different radar backscattering coefficients. As an example, a theoretical computation [10] demonstrates that two textureless target classes of homogeneous radar reflectivities (i.e.exhibiting only speckle) radiometrically separated by 2.5 dB present a probability of confusion of 40% in a 1-look SAR image.

Therefore, a first radiometric enhancement is needed to achieve a reduction of the coefficient of variation of the speckle over homogeneous areas. It corresponds to an enhancement of the signal-to-noise ratio and to a preliminary speckle "noise" reduction.

2.1.3. First order statistics in "multilook" images

A first method of speckle reduction consists in averaging incoherently Mindependent samples of the intensity I(or "looks") obtained from the same target:

I =(1/M).k=1MIk with each of the Ik distributed according to Equation (7)E10

The goal of this method is to reduce speckle enough to make radar image photo-interpretation possible. Indeed, experience has shown that Mvalues of the order of 3 or 4 enable a photo-interpreter to use a SAR image [1]. Such values have been adopted for most spaceborne SARs (ERS, Almaz, JERS-1, Radarsat, Envisat, ALOS: 3-looks, Seasat, SIR-B: 4-looks), including recent ones (TerraSAR-X, Cosmo-Skymed, Sentinel-1, etc.).

This operation is realised by splitting the Doppler bandwidth of the backscattered signal into Msections. This is equivalent to divide the length of the synthetic antenna into Msub-antennas. Independent processing of each section results in an independent image called a "look". Nevertheless, this operation results in a degradation by a factor Mof the spatial resolution in azimuth of each look [1].

Over the same target, the mean intensity resulting from this operation remains the same mean intensity of each of the individual looks. If the Mindividual looks were independent, the standard-deviation is divided by M. Thus, the coefficient of variation of the speckle measured over a homogeneous area of the intensity multilook image becomes:

CI (homogeneous area)  =  Cu =  1/ME11

It is important to note that multilook radar image formation is at the expense of the spatial resolution in azimuth. In practice, the value of Mdoes never exceed a few units (less than 16 for airborne SARs, and in general only 3 or 4 for spaceborne SARs). This remains insufficient to improve satisfactorily the signal-to-noise ratio of a radar image. Indeed, in our example of two target classes with homogeneous radar reflectivities radiometrically separated by 2.5 dB, the probability of confusion of 40% in a 1-look SAR image reduces only to 33% in a 4-look image.

If the individual looks are uncorrelated with each other, the pdf of the speckle, which is the sum of Mindependent exponential distributions, becomes a χ2 distribution with 2Mdegrees of freedom:

Pu(I) = Pu(I/R) = (M/E[I])M / (M1)! . exp(M.I/E[I]) . I(M1)E12

2.1.4. "Equivalent Number of Looks" of a SAR image

If the Msections of Doppler bandwidth used to produce the individual looks overlap, the averaged samples, and therefore the individual looks, are correlated. The coefficient of variation of the speckle Cu, measured over a perfectly homogeneous image area will therefore be always superior to what could have been expected from Equation (11). In this situation, speckle strength results as if the number of independent looks were equal to:

L = 1/Cu2   <  ME13

The Lvalue, which is generally a non-integer value, is called the "Equivalent Number of Looks", or "ENL", of the multilook radar image.

Hence, the pdf of the speckle in intensity can be approximated, by extension of Equation (12), for whatever value of L, by a Gamma distribution with parameters E[I] and L:

Pu(I) = Pu(I/R) =(L/E[I]) LΓ(L) . exp(L.I/E[I]) . I(L1)E14

which is rigorously equivalent to Equation (12) when L=M(Linteger) for Mindependent looks.

All speckle models for multilook images (Mcorrelated or uncorrelated averaged looks, i.e.Lequivalent looks, with ML) use this approximated distribution. Note that the pdf of the speckle as formulated in Equation (14) is slightly inexact if the Mlooks are correlated, as it is the case in general. Therefore, the pdf of the speckle for the average of Mcorrelated looks is close to, but is not exactly, a Gamma pdf [11].

Let us consider that the Mlooks (I1,I2,...,IM) are correlated to each other with correlation coefficients ρij, and that the Iiare distributed according to the same exponential marginal pdf (1-look, cf.Equation (7)) with same parameter R=E[I]=E[Ii], for whatever i. This last condition can be actually fulfilled by re-weighting the individual looks to make them contain an identical mean energy E[I].

For the intensity image resulting of the averaging of the Mlooks, I=(∑Ii)/M, the set of looks correlations ρijis taken into account to compute the ENL value L(L<M), which is related to the coefficient of variation of the speckle in intensity, Cu, by:

Cu2 =  1/L = [ 1 + (1/M) .  i=1Mj=1i1σij2] / ME15

The exact pdf of the average intensity Iof Mcorrelated looks, is extremely difficult to establish mathematically in a closed-form (for 2-looks, cf.[11]), and, when established, to manipulate. Nevertheless, Equation (14) is a very close and therefore satisfactory approximation when using the appropriate ENL value Lcalculated using Equation (15).

2.1.5. The logarithmic (homomorphic) transformation

Since the radar backscattering coefficient σ° is generally expressed in dB in physics, and with the goal to make speckle an additive noise that would be independent on the level of the radar signal, some authors [12] [13] [14] chose to use a neperian logarithmic (homomorphic) transformation:

D =  Log (I)E16

Arsenault and April [15] [17] [19] demonstrated that, after this transformation, the pdf of the speckle for a multilook radar image, Equation (14), becomes a Fisher-Tippett distribution:

Pu(D) =  LL/ Γ(L) .  exp[L.(DD0)] . exp{L . exp[(DD0)]}E17

with D0=Log(E[I])

The first order statistical moments of this distribution are as follows [15]:

1) Mean <D>:

<D> = D0  Γ'(L)/Γ(L) + Log(L)E18

Equation (18) shows that the logarithmic transformation causes a signal distortion, increasing with decreasing number Lof independent speckle samples taken into account, and with decreasing R(therefore, σ°) values. There is therefore a tendency to systematically underestimate the value of σ° [15] [16] [18].

2) Variance σD2:

σD2 = ( π2/6) . [ γ+ Γ'(L)/Γ(L)]2  +  2  k=1L2 [ (1/(Lk). j=1Lk1(1/j) ]  E19

Lopès [20] has shown that for N independent L-look radar data samples (NLindependent speckle samples), the standard-deviation of the samples averaged after logarithmic transformation, and then retransformed into intensity, is always significantly larger than the standard-deviation σIof the same samples averaged in intensity. Thus, the logarithmic transformation degrades both the measurement accuracy of the backscattering coefficient σ° and its local range of fluctuations due to local scene texture.

2.2. Texture in SAR images

Texture concerns the spatial distribution of grey levels. It is important in the analysis of SAR images for a wide range of applications (e.g.[21] [22] [23] [24] etc.). In most of these applications, not only the radiometric, but also the textural information must be retrievable after adaptive speckle filtering of the SAR images.

2.2.1. Texture of the imaged scene

As seen above, within a homogeneous area image of a detected radar image, one can consider the observed speckle Ias originating from a real non-gaussian random process u(x,r), with unit-mean <u>=1, and proportional to the radar signal I(multiplicative noise). A simplification of Equation (3) enables to write, for every pixel (x,r) located within a homogeneous image area (i.e.R(x,y)=E[R]=R, (x,y)):

I(x,r) = u(x,r) . E[I(x,r)]= u(x,r) . RE20

In most remote sensing applications, it is reasonable to consider the imaged scene as an arrangement of discrete objects superimposed to a background with mean reflectivity E[R]. The imaged scene is composed of classes of non-homogeneous (R(x,y)E[R]) objects characterised by the statistical properties and parameters of the variable R[25]:

If Tjis the random variable which represents the spatial fluctuations of the reflectance (let consider in the following its analog in radar images, the denoised/speckle-filtered radar reflectivity) R(x,r) around its mean expectation E[R(x,r)], in the terrain class jto which the pixel at coordinates (x,r) belongs, the following results are obtained:

  1. R(x,r) is proportional to the σ° value of the resolution cell containing pixel (x,r);

  2. E[R(x,r)] is proportional to σ° in a neighbourhood of pixel (x,r) belonging to class j;

  3. E[Tj]=1 for every class j;

R(x,r) = E[R(x,r)] . TjE21

The spatial structure of actual imaged scenes, called "texture", induces a mesurable spatial structure (after the properties of  Tj) in the images of these scenes. This is the case, for both optical and radar images. Indeed, a simple visual interpretation of a radar image reveals radiometric spatial variations at a longer scale than those due to speckle only. These variations originate from the spatial fluctuations of the radar backscattering coefficient σ° within a given terrain class, and are affected by the presence of speckle.

To characterise a given class j, one must be able to restore  Tj, that means to separate the respective effects of speckle and texture in the spatial variations of the radar intensity signal. To this end, the multiplicative fully-developed speckle model exposed above will be used. It has been shown [26] [3] that, even when the scene is not stationary, speckle still remains multiplicative, eventually correlated, and independent on the imaged scene as long as Goodman’s conditions [2] [5] are realised.

2.2.2. First order statistics of texture in SAR images Spatially uncorrelated speckle

Considering the scene texture model of [25] in Equation (21), it is possible [27] [28] to model the radar image intensity at pixel of coordinates (x,r) by generalising the multiplicative speckle model for a wide-sense stationary scene (Equation (20)) as follows:

I(x,r) = R(x,r) . u(x,r) = (E[R(x,r)].Tj) . u(x,r)E22

This relationship is valid as long as the spatial variations of R(x,r) happen at length scales longer than the size of the resolution cell of the radar imaging system.

The variance of the intensity Iwithin a given target class is computed using:

σI2 = E[I2]E2[I] = E2[R].(E[T2.u2]E2[T.u])E23

Since the fluctuations of σ°, and therefore those of R, are independent on the speckle, one gets [27] [28]:

(σI/E[I])2 = (σR/E[R])2 .(σu/E[u])2 + (σR/E[R])2 + (σu/E[u])E24


E[I] = E[R.u] = E[R].E[u] = E[R]E25

Therefore, E[R] is locally estimated in the radar image by averaging pixel intensities in a neighbourhood of pixel (x,r):


Introducing the coefficients of variation (CI=σI/<I>for the radar intensity, CR=σR/<R>for the imaged scene, Cu=σu/<u>for the speckle), Equations (24) and (25) lead to the important result [27] [28]:

CR2  =  (CI2  Cu2) / (1 + Cu2)E27

which characterises scene texture in terms of heterogeneity, with Cu2 = CuI2 = 1/Lfor a homogeneous/textureless area of a L-look intensity radar image. Indeed, Equation (27) shows that the more heterogeneous is the scene, the easier its texture can be restored.

For a better description of the scene, one must use the pdfs of its diverse classes (distribution of the random variable Ror σ°). Their knowledge may result, either from a priori knowledge, or from direct estimation in the image. If the general form of the pdf is known, its parameters can be estimated from the image data. In all cases, if one can establish the pdf of R, PR(R), for a given terrain class in the scene, the unconditional pdf of the corresponding radar image intensity Iis:

PI(I) =   0Pu(I) . PR(R) . dR   E28

where Pu(I)is the pdf of the speckle for multilook images (cf.Equation (14)). Equation (28) shows, as does also Equation (27), that the respective contributions of scene texture and speckle micro-texture can be separated in the radar image. Spatially correlated speckle

To produce a SAR image, speckle samples at the output of the SAR processor are correlated through the SAR impulse response to obtain a sampling rate (pixel size) of about half the spatial resolution of the radar sensor, thus avoiding aliasing effects. Using the multiplicative noise model and Equations (24) and (25) [27] [28], the correlation coefficient ρI(Δx,Δy)between two pixels of intensities I1 and I2 separated by (Δxy) in the SAR image:

ρI(Δx,Δy) = (<I1(x,y)<I1>> . <I2(x,y)<I2>>) / (σI1 .σI2)E29

can be related to the scene correlation coefficient ρR(Δx,Δy)and to the speckle correlation coefficient ρu(Δx,Δy)by [10]:

<CI>2. ρI(Δx,Δy) = <CR>2.<Cu>2.ρR(Δx,Δy).ρu(Δx,Δy) + <Cu>2.ρu(Δx,Δy) + <CR>2.ρR(Δx,Δy)E30

Therefore, if ρu(Δx,Δy)>0 for Δx>0 and/or Δy>0, which is the case in correctly sampled SAR images, the local statistics estimated in the neighbourhood of a given pixel will have to be corrected of the effect of speckle spatial correlation properties (cf.Section 4).

2.2.3. Local statistics computation and consequences of the speckle model

Using the simplified multiplicative noise model for the speckle u(cf.Equation (22)), the first-order non-stationary statistics of the scene, <R> and σR2, can be deduced locally (cf.Equations (26) and (27)) from those of the radar image intensity, <I> et σI2.

In practice, CIis estimated in the radar image, and its locally estimated value can be inferior to Cu. Since σR2and CR2are positive quantities, finding cases where CI2<Cu2must be attributed to the limited size of the neighbourhood of the pixel under processing over which <I> and σI2are estimated [29]. Independently of the scene model PR(R), the multiplicative speckle model fixes therefore an inferior threshold CImin=Cuto the possible values of the coefficient of variation CI.

Below the CIminvalue, speckle filters based on the use of local statistics, it means all adaptive speckle filters, are no longer valid. In this case, the image area under processing must be considered homogeneous and textureless.

2.2.4. The coefficient of variation as a heterogeneity indicator

The coefficient of variation CIis a heterogeneity measure particularly well appropriate to the case of radar images and for isotropic textures. Nevertheless, if the source of heterogeneity presents a particuliar orientation (contour, line...), additional detectors able to identify or retrieve this orientation are needed.

CIincreases with scene heterogeneity. Depending on the heterogeneity of the area under processing, it enables to discriminate among the following situations:

  1. homogeneous areas with neither texture, nor structures.

  2. areas with low or moderate heterogeneity where the presence of structures is improbable.

  3. strongly heterogeneous areas where one must search for the possible presence of structures (contour, line, strong scatterers).


3. Adaptive speckle filters for single-channel detected radar images

This section is dedicated to a theoretical analysis of the most common and more efficient speckle filtering techniques developed for Synthetic Aperture Radar (SAR) images, and of their corresponding estimation techniques. These filters use variants of the statistical speckle model exposed in the preceding Section 1. They also use diverse statistical estimators to restore the radar reflectivity (in the sense of the CEOS’s (Committee for Earth Observation Satellites) "radar brightness"). These statistical estimators are, either Minimum Mean Square Error estimators (e.g.Lee et al.filter, Kuan et al.filter, etc...), or autoregressive estimators (e.g.Frost et al.filter), or Bayesian estimators (e.g.Gamma-Gamma and DE-Gamma MAP filters). These estimators are discussed, and their behaviours are analysed.

3.1. Wiener Method: the Frost et al.filter (1982) [30]

3.1.1. Theoretical development

The SAR image model considered by [30] is as follows:

I(t) = [R(t) . u(t)] * h(t)E31

where R(x,y) is a wide-sense stationary random process decribing the radar reflectivity of the observed scene at pixel of coordinates t=(x,y) located in a homogeneous area of the radar image. u(t) is the multiplicative noise due to speckle, modelled by a real white stationary random process with a Gamma pdf (cf.Equation (14)). h(t) is the impulse response function of the system. This model is a simplification of the fully developed speckle model exposed above (cf.Equation (3); [31] [32].

To estimate the radar reflectivity R(t) of the noise-free image, Frost et al.apply a linear filter whose impulse response m(t) minimises the mean quadratic error (MMSE estimator):

 ε= E[{ R(t) I(t)*h(t) }2]E32

The MMSE least-squares solution is valid for homogeneous areas for which R(t) can be represented by a stationary random process. To deduce m(t), [30] assume that the transfer function H(f) is constant over a certain bandwidth. This leads to an uncorrelated multiplicative speckle model:

I(t) = R(t) . u(t)E33

where u(t) is a pseudo-white noise, independent of the signal, also used by other authors [33] [34]. The multiplicative speckle model is deduced from Equation (3) assuming that the bandwidth of the signal R(t) is much narrower that that of the linear filter h(t).

The impulse response of the filter is calculated adopting an autoregressive process model for R(t) with an exponential autocorrelation function, a classic hypothesis for this family of models [35] [36]. Later studies showed that the choice of an exponential autocorrelation function is appropriate to SAR scenes [37].

These hypotheses enable to define an optimal MMSE filter (Wiener filter) with impulse response m(t) :

R^=  t  neighborhood of N pixelst=(x,y)m(t) . I(t)E34


m(t)  =  K1 . exp[ K . CI2 . d(t) ]E35

where Kis the filter parameter, and d(t) is the distance between the pixel located in tand the pixel under processing. K1is a normalisation constant granting the preservation of the mean value of the radar reflectivity:

K1 =  (1/N).  t neighborhood of N pixelst=(x,y)m(t)E36

The final Frost et al.filter (Equation (36)) is not exactly a Wiener filter. It is designed using Yaglom’s method [38]. In this method, contrarily to that of Wiener, one searches for the frequential characteristic, not for the impulse transfer function susceptible not to exist. This function is chosen empirically taking into account the hypotheses it must fit.

3.1.2. Behaviour of the Frost et al. filter

The behaviour of the filter depends on the values of the locally observed coefficient of variation, i.e.on the local heterogeneity:

  • Extremely homogeneous areas, for which CI = 0,     then:     R^ = <I>

  • Very strong scatterers (extreme heterogeneity), where CI , then:      R^ = I

  • Between these two situations, in textured natural areas, when CIvalue increases, neighbour pixels located far away from the pixel under processing are given less weight and the neighbourhood on which the weighted mean value is actually estimated narrows.

Considering the requirements for an ideal speckle filter, the filter should:

  • Restore the mean intensity value over homogeneous areas, where CICu. A consequence is that Kmight be close to 0, unless the value of Cutends towards 0, which is obviously not possible for usual radar images.

  • Preserve the observed intensity when the value of CIis very high but obviously always finite. For high CIvalues, the only pixel to weigh might be the pixel under processing, and Kmight therefore tend towards infinity since CIhas always a finite value.

Since it is impossible to fulfil simultaneously these two conditions (cf.[39] [40], one must make a difficult choice between, either an efficient speckle reduction, or the preservation of image structures and texture.

3.2. Locally adaptive linear MMSE estimators: Lee and Kuan filters

3.2.1. Kuan et al. filter (1985) [41] [32]

The radar image I(t) is modelled as a function of the radar reflectivity image R(t) to be restored and of a white uncorrelated noise n(t) with 0 mean and dependent on R(t). This is an additive noise model different of the multiplicative speckle model used in Equation (33):

I(t) = R(t) + n(t)E37

Kuan et al.[41] introduce a scene model R(t) with non-stationary mean and non-stationary variance (NMNV model), where the non-stationary mean describes the general structure of the image, whereas the non-stationary variance characterises texture and the presence of local structures. With this model, the linear filter minimising the quadratic mean error (LMMSE) has the general form [32]:

R^LMMSE =  E[R]  +  CRI. CI1. ( I <I>)E38

where CRIand CIare the non-stationary spatial covariance matrices of the image. These matrices are diagonal if, as it is assumed in the model, n(t) is a white noise and the image model is NMNV. Then, Equation (38) takes a scalar form [31]:

R^ =    E[R] + ( I <I>) . σR2 / (σR2 +σn2 )E39

The local statistics of the scene Rare deduced from those of the intensity I. The Kuan et al.filter is obtained by replacing locally:

E[R]   by  <I>      and      σR2  by  [σR2 <I2>.σn2 ] / [1+σn2]  E40

Since the variance of the additive noise σn2is numerically equal to the coefficient of variation of the speckle Cuin the multiplicative model, Equation (40) is, in practice, exactly equivalent to Equations (26) and (27) established for the multiplicative speckle model.

It is noteworthy that, since the image model does not assume independence of the noise nand the signal R, the filter is theoretically able to take into account the effects of an eventual dependence between speckle and the radar reflectivity, when speckle is no longer fully developed as in the case of strong scatterers.

3.2.2. Lee filter (1980) [33] [42] [43] [44] [45]

Lee uses the unit-mean uncorrelated multiplicative speckle modele (cf.Equation (33)). A linear approximation is done by developing Equation (33) (Ias a function of ufor pixel located in t) in a first-order Taylor series with respect to u(t):

I(t) = R(t). E[u(t)] +  E[R(t)] .[u(t)E[u(t)] ]        with    E[u(t)] = 1E41

This approximation enables to transform Equation (33) into a weighted sum of the signal and of a noise independant on the signal. The linear MMSE estimator, and Equations (26) and (27) that are consequences of the multiplicative speckle model, enable to establish the Lee filter [33], historically the first speckle filter designed to be adaptive to local radar image statistics:

R^ = <I> + ( I <I>). (σI2 σu2)/σI2E42

It is noteworthy that, assuming the independence of noise and signal in the model used by the Kuan et al.filter, this filter become identical to the Lee filter [41].

3.2.3. Behaviour of the Kuan et al. and Lee filters

The linear MMSE speckle filters of Lee and of Kuan et al., based on the use of the local statistics of the observed intensity can be written under the same general form:

R^(t) =  I(t). W(t) + <I(t)>. [1W(t)]E43


W(t) = (1[Cu2/CI2(t)])/(1+Cu2)         for the Kuan filterE44


W(t) = 1  [Cu2 / CI2(t)]    for the Lee filterE45

Both methods perform a linearly weighted average of the local mean of the intensity and of the observed pixel intensity. In both cases, weights depend on the ratio of the coefficients of variation of the observed intensity and of the noise [7] [40]. As for the Frost et al. filter, the locally observed heterogeneity governs the behaviour of these filters:

  • Homogeneous areas where CI=Cu, therefore R^=<I>

  • Very stong scatterers (extreme heterogeneity) where CI→∞, therefore R^=I

  • Between these two extreme situations, the weight of the observed value of the pixel intensity increases with the heterogeneity of its neighbourhood.

Nevertheless, note that the weight on the estimated mean intensity remains significant for high values of the coefficient of variation. Thus, responses of small impulse targets are cut-off by these filters and image structures suffer some amount of smoothing. To correct this drawback, Lopès et al.[40] have proposed a very efficient enhancement to the Lee and Kuan et al.filters.

Inversely, when CI<Cu, some weighting on the observed intensity reappears, resulting in an amplification of the noise [13]. In this case, the local mean of its neighbour pixels must be assigned to the pixel under processing as filtered value.

The linear MMSE filters differ from the Wiener filters by the fact that the A Priori mean and variance, <R> and σR2, are estimated using the local statistics of the original image I, and not implicitly estimated from an assumed autocorrelation model (cf.[30], Equation (35)). The goal of this approach is to obtain a better adaptivity to the local properties of the scene R, which verifies in pratice.

It is remarkable that the often used MMSE estimation ([33] [43] [44] [45] and [41] [32], among others) is nothing else than the mean of the conditional A Posteriori pdf P(R/I) [46]. Note that the complete evaluation of the terms enabling to obtain the A Posteriori pdf P(R/I) requires a description of PI(I). This last pdf, which depends on the non-linear transformation leading from Rto Iis, in general a K-distribution. The solution adopted by Lee and Kuan et al.consists in a forced linearization of the problem assuming a Gaussian PI(I) [33] [42] [44], which is unjustified for low (< 3) ENL values.

As an effect, this linearization restrains the validity of the linear MMSE estimator to situations where the noise level is not too high (multilook images with high values of L), or where the local heterogeneity CIis not too high. These considerations justify one of the modifications proposed to the Lee and Kuan et al.filters par Lopès et al.[40]: an upper heterogeneity threshold CImax, above which the observed intensity of the pixel under processing is preserved.

3.3. Bayesian Maximum A Posteriori (MAP) speckle filters

3.3.1. Originality of the MAP approach in the case of SAR images

The Bayesian MAP approach [47] [48] consists in characterising the imaged scene and the speckle "noise" by their statistical description, using their associated pdfs.

In the Bayesian perspective, the theory of probabilities is extended to the logics of probabilistic inference. Probabilities are seen as a relationship between a formal hypothesis and a possible conclusion. This relationship corresponds to a certain degree of rational credibility and is limited only by the extreme relationships of certitude and impossibility [49]. The classic deductive logics considering only these extreme relationships (à laSherlock Holmes: "When all impossible hypotheses have been excluded, what remains, even if unlikely, corresponds to reality") is nothing but a particular case of this general development [50].

As a general consequence, the theory of probabilities cannot base on the sole concepts of classic logics (frequentist as in the MMSE or in Wiener’s approaches). In particular, the relationship of probability cannot be defined in terms of certitude, since certitude is viewed as a particular case of probability. The frequentist definition of probabilities based on relations of certitude related to the knowledge of a number of parameters (<R> for example, cf.§ 3.3.5 below) is therefore no longer sufficient [50].

In this context, probabilities are used to describe stochastically an incomplete information on a global phenomenon (here, the radar reflectivity and the superimposed speckle), rather than to describe only the noise randomness that corrupts its comprehension. Probability relationships are viewed as being conditional to the context. This way, the pdf of the speckle Pu(I)(cf.§ 2.1.4, Equation (14)) is therefore formulated under the form P(I/R). Taking into account a reasonable A Priori statistical model for the radar reflectivity P(R), one must also take into account – which is new, with respect to the methods precedently exposed – the probability of R, given the information obtained through observation (the intensites Imeasured by the radar), formulated under the form of an A Posteriori pdf P(R/I). This approach has the great advantage to enable the characterisation of speckled radar image formation while easing the description of non-linear effects.

The least error-cost inference mechanism leading from the observed intensity Ito the A Posteriori deduced radar reflectivity Rthrough Bayes’s theorem [47] [48], allows a rigorous combination of the A Priori knowledge P(R) and of the new knowledge provided by the observation I:

P(R/I)  =   P(I/R) . P(R) / P(I)E46

P(R/I) depends on the pdf of Rintroduced as A Priori information about the scene to restore. Thus, the estimate is influenced by the A Priori statistical knowledge about Ror, by default, the hypotheses made about R. The noisiest are the radar data, the less the A Priori information would contribute to estimate R[46].

In theory, the MAP method enables to avoid direct estimation of the mean of the conditional A Posteriori pdf, which is necessary to the MMSE estimation. This feature is of great interest in the resolution of non-linear problems where the evaluation of the conditional mean is, either difficult, or impossible [51].

The Maximum A Posteriori (MAP) estimation of Rcorresponds to the maximum value of P(R/I), the mode of the A Posteriori pdf P(R/I). The mode and the mean of P(R/I) coincide if P(R/I) is symmetrical, in particular in the Gaussian case. Since the available knowledge is:

  • P(I/R), the conditional pdf of speckle, which is a Gamma distribution (cf.Equation (14)),

  • P(R), the A Priori unconditional pdf of the imaged scene, which is an assumed statistical distribution (not necessarily Gaussian),

there is no particular justification to prefer the MMSE estimation by minimisation of the mean quadratic error implicitly assuming a symmetrical P(R/I), rather than to use this more complete available knowledge through the Bayesian estimation technique, which is less subject to arbitrary assumptions [46].

It is of interest to notice that the "Maximum Likelihood" estimation is nothing more than a particular case of Maximum A Posteriori, in the situation where the conditional A Posteriori and A Priori pdf’s are equal: P(R/I)=P(I/R), with P(R)=P(I). Note that the Maximum Likelihood estimation of the scene Ris the radar image Iitself, which is of interest in the case of strong scatterers whose radar response is deterministic.

Once the first-order statistical speckle model P(I/R) is either known, or reasonably approximated, (cf.Equation (14)), the form of the MAP estimate depends principally on the form of P(R).

The NMNV scene model [32] enables to solve locally, in an analytic manner, the estimation of R. The parameters of the pdf P(R) are locally estimated in a neighbourhood of the pixel under processing. This way, the complicated recursive forms of mathematical resolution as exposed by Hunt [52], Kuan et al.[53] and Geman & Geman [54] are avoided. This approach is particularly adapted to the high spatial resolution of airborne and modern spaceborne SARs. Indeed, it allows a good preservation of scene texture in an adaptive way by taking into account the local fluctuations of σ° in the A Priori model, while being at the same time computationally efficient.

3.3.2. Computation of the MAP estimate

Since the logarithm function is a monotonically increasing function, Bayes’s formula ([47] [48]; Equation (46)) can be rewritten as:

Log[P(R/I)] = Log[P(I/R)]  Log[P(I)] + Log[P(R)]E47

which gives the local MAP estimation of Rwhen Log[P(R/I)] is maximum, i.e.when its first-order derivative with respect to Ris locally equal to 0:

/R Log[P(I/R)]+/R Log[P(R)]=/R Log[P(R/I)]+/R Log[P(I)]  E48


/ RLog[P(I)]=0 because P(I) does not depend on the local fluctuations of R.

P(R/I) reaches its maximum when:

/R Log[P(R/I)]= 0   for   R=R^MAP,   where R^MAPis the MAP estimation of RE49

Then, the general equation of the MAP speckle filter becomes, locally:

/R Log[P(I/R)]+/R Log[P(R)]= 0    for   R =  R^MAPE50

The first term of Equation (50), the Maximum Likelihood term, accounts for the effects of the whole imaging system on the radar image and describes the detected radar intensity once the speckle statistical model is known. The second term, the Maximum A Priori term, represents the prior statistical knowledge with regard to the imaged scene.

In the Bayesian approach, probabilities are used to describe incomplete information rather than randomness. As Equation (50) shows, in the Bayesian inference process, induction is influenced by the prior expectations allowed by the prior knowledge of P(R) [46]. In addition, the non-linear system and scene effects are implicitly taken into account by the restoration process. Therefore, MAP speckle filtering can be considered as a controlled restoration of R, where A Priori knowledge controls the inference restoration process and allows an accurate estimation of the radar backscattering coefficient σ°.

The pdf of the speckle in intensity for a L-looks radar image, P(I/R), is a Gamma distribution with parameters Rand L(cf.Equation (14)). The Maximum Likelihood term for a L-looks radar image is then equal to:

/R Log[P(I/R)]  =  L . ( I/R2  1/R )E51

The Maximum A Priori term must be calculated according to the scene model chosen as A Priori knowledge, hypothesis, or belief.

3.3.3. Gaussian distributed imaged scene: The Gaussian-Gamma MAP filter (1987)

The hypothesis of Gaussian-distributed scene has been adopted as a natural hypothesis by a large number of authors who had worked, either on images from optical sensors, or on images from passive/active microwave sensors. These authors are comforted in this hypothesis by both the force of habit and by the mathematical ease in manipulating a Gaussian distribution.

Kuan et al.[32] have developed a MAP filter (Gaussian-Gamma MAP filter, with a Gaussian-distributed scene and Gamma-distributed speckle) for radar images under this hypothesis. Though, the hypothesis of a Gaussian-distributed scene, although widely spread in the literature, is inappropriate. Indeed, this hypothesis assumes implicitly the theoretical possibility of negative σ° values-which has no physical sense-in the extreme case of a large variance with a low mean value, therefore needing a regularisation of the filter behaviour in such a situation.

Therefore, one must preferably take into account a positive pdf as a realistic scene model. For reasons that are both experimental and theoretical, in the case of natural extended targets as it is most often the case dealt with in remote sensing, a Gamma distributed scene model is better appropriate.

3.3.4. Gamma distributed imaged scene: The Gamma-Gamma MAP filter (1990)

Natural textures observed, either by coherent radar imagery, or by incoherent optical imagery are due to a common contribution corresponding to the variability in spatial distribution of the objects within the scene. Even if the interaction mechanisms between the electromagnetic wave and the observed medium are very different in either case, the natural arrangement of the scene makes the second-order statistics very similar in either kind of imagery [55] [56]. At the scale of a large number of resolution cells, the pdfs of the cross-section variables corresponding to either mechanism belong to the same family of distributions, at least for the high radar frequencies in bands Ku, X, and C for which wave penetration into natural media is limited. This point is more arguable for radar bands L and P.

In a wide range of radar backscattering situations, the Gamma distribution is experimentally the one that best fits, not only the distribution of the radar backscattering coefficient [21] [57], but also the distribution of radiometries observed in incoherent optical images [55]. This scene model has been successfully used also for radar images of the sea [58] [59] [60].

The local pdf of a scene statistically described by a Gamma distribution, has the form:

P(R) =  (α/E[R])α/Γ(α) . exp(α.R/E[R]) . Rα1E52

with α= 1/CR2. The parameter is called the "heterogeneity coefficient".

Note that assuming a Gamma-distributed R, and by performing the integration in Equation (28), the pdf P(I)=P(I/R) of the intensity is a K-distribution [57]. Introduced in 1976 [58] by British researchers of the RSRE (later DRA) to describe the non-Gaussian properties of waves backscattered by objects within a radar resolution cell, the K-distribution has been theoretically recognised as the pdf of the intensity Ibackscattered by a rough non-stationary surface [61] such as most natural scenes observed by a radar.

Nevertheless, in an illustration of the so-called "Cromwell’s rule" [62], even hard A Priori conviction that the scene presents a Gamma-distributed texture must not be insensitive to counter-evidence. Therefore, the complete MAP filter is a set of three filters adapted to diverse situations locally encountered in a radar image: the application of this or that filter is decided depending on the degree of heterogeneity of the image part under processing, that is on the locally estimated value of the coefficient of variation CI. This may eventually imply the determination of thresholds on CIcalculated as a function of a user-defined probability of false alarm with respect to the local presence of texture or of strong scatterers (cf.[63] [64] [10]). Note that these considerations apply to all MAP speckle filters, and can be extended to all other adaptive speckle filters. Textured areas

Assuming that the pdf of the scene is a Gamma distribution, the Maximum A Priori term is locally equal to: / RLog[P(R)]=(α-1)/R-α/<R>. Once E[R] is estimated locally by <I> in an ergodic stationary neighbourhood of the pixel under processing in the radar image, the equation of the Gamma-Gamma MAP filter ([65] [10] [63] [64]) is:

α.R2 + (1+Lα).<I>.R  L.I.<I>  =  0E53

This second-degree equation admits only one real positive solution Rin the interval ranging between the mean intensity <I>=<R> and its observed value I. Therefore, the Gamma-Gamma MAP estimate of the radar reflectivity of the pixel under processing, is:

R^=<I>.(αL1) +  <I>2.(αL1)2 + 4α.L.I.<I> 2αE54

The integration of a heterogeneity/texture detector based on the coefficient de variation and of specific detectors (ratio-of-amplitudes – RoA; [66] [7]) for contours, linear structures and strong scatterers in the filtering process is described by the general algorithm presented in [63] [64] [10]. The integration of texture and structure detectors using second-order statistics (autocorrelation functions) of both the speckle and the radar reflectivity of the scene [67] is presented in Section 4. In all cases, mage areas identified as textured are filtered using Equation (54). Homogeneous areas

In the particular case of a perfectly homogeneous (textureless) scene, with CICu, the radar reflectivity Ris a constant (R=E[R]), and can be statistically represented by a Dirac distribution:

P(R)  =  δ(R)E55

This distribution is the limit of Gamma distributions when tends towards+∞. In such a case, the MAP estimate is equal to the local mean intensity:

R^ = <I>E56

This case is taken into account when the local statistics calculated in the neighbourhood of the processed pixel show a nearly perfect homogeneity of the scene. Strong scatterers and impulsional targets

The other extreme case regards strong scatterers, when speckle is no longer fully developed. The considerations that led to consider the Gamma pdf as an A Priori model for a textured natural scene are no longer valid: we no longer have any A Priori information about the scene. In this situation, the information content of every grey level of the image being A Priori the same, the pdf of the imaged scene can be represented by an uniform distribution:

P(R) =  1/(Rmax  Rmin)E57

with undetermined extreme values Rminand Rmax, which we may consider equal to respectively 0 and+∞. Thus, the Maximum A Priori term of the MAP equation becomes / R[Log(P(R))]=0, and the MAP estimate is:


If the resolution cell contains only one isolated strong scatterer, the response (value Iof the pixel) of this scatterer is deterministic and must therefore be preserved. This situation leads to the same conclusion and is therefore treated similarly.

3.3.5. Scene pdf estimated from the data: The Gamma/Distribution-Entropy MAP filter (1998)

SAR images of dense tropical forest, urban areas, or very strong and rapidly varying topography often show very strong or mixed textures. This is also the case of high-and very-high spatial resolution SAR images. In these situations, it may be hazardous to make an assumption about the probability density function of the radar reflectivity.

Indeed, the MAP technique does not account for any uncertainty in the value of the parameters of the A Priori pdf chosen as a Gamma distribution once it has been locally estimated on a given image area. Hence, in the presence of mixed (forests with underlying structures, for example) or rapidly varying (strongly textured area located on strong slopes, or very-high spatial resolution radar images, for example) texture, the MAP estimator will underestimate the variance of the predictive distribution. Indeed, this predictive distribution can hardly take into account the fact that it results of a compound of a mix of different distributions.

In this context, the A Priori knowledge with regard to the observed scene can hardly be an analytical first order statistical model, chosen on the base of prior scene knowledge. However, to retrieve local statistical scene knowledge directly from SAR image data, Datcu & Walessa [68] [69] proposed to introduce the local entropy of the radar reflectivity, S(R), as a measure of local textural disorder. This concept originates from the theory of information. S(R) is estimated on a neighbourhood (Npixpixels) of the pixel under processing:

S(R) = k=1Npix[Rk. log(Rk)]   E59

Because the radar reflectivities Rkare non-negative and exp[S(R)]/Zis normalised, and since S(R) is a measure of the spread/dispersion of the radar reflectivities of the scene, characterising its heterogeneity, Equation (59) can be treated as a pdf whose entropy is S(R) [70]:

P(R) = (1/Z. exp[S(R)] =  (1/Z). exp(k=1Npix[Rk. log(Rk)] )E60

For a single detected SAR image, the conditional pdf of the speckle can be, as long as speckle is fully developed, modelled as a Gamma distribution:

P(I/R) = (L/R)L/Γ(L) . exp(L.I/R) . IL1E61

Incorporating these scene and speckle models, the Gamma/Distribution-Entropy MAP (Gamma-DE MAP) speckle filter for single-channel detected SAR images is the solution of the following equation [71]:

L.I  L.R  R2.k=1Npix[log(Rk)1/Ln(10)]  = 0 E62

The radar reflectivites Rkin the neighbourhood of the pixel under processing are pre-estimated by a first speckle filtering pass.

Note that the local DE MAP estimation of Ris contrained by the value of the entropy S(R) retrieved from image data. Since this fixes upper/lower bounds to the local entropy S(R), the DE MAP speckle filter combines both the Bayesian MAP and the Maximum/Minimum Entropy estimation techniques.

The DE-MAP filters adapt to a much larger range of textures than the other MAP filters ([10] [11] [72] [73]) developed under the assumption of K-distributed SAR intensity (i.e.Gamma-distributed scene). In particular, these filters are of particular interest in the case of very high-resolution SAR images and strongly textured scenes.

Compared to those of the other MAP filters, performances in terms of speckle reduction are identical. However, texture restoration and structures or point targets preservation, identical for moderate textures, are superior in strongly textured areas. These filters have proven a remarkable efficiency in operational remote sensing (cf. [74]).

3.3.6. Behaviour of the MAP filters

The local MAP estimate is the mode of the local A Posteriori pdf of the radar reflectivity R, under the hypothesis made for P(R). The MAP estimation is therefore probabilistic: it corresponds to the most probable value of radar reflectivity. This justifies to chose the local mean observed intensity as filtered value within areas identified as perfectly homogeneous (cf.§ 2.4.4; [64] [10). Excellent speckle reduction and preservation of the mean value of the radar reflectivity are therefore granted over such areas.

In the presence of texture, the estimator takes into account the non-linear effects in getting from Rto Ivia the imaging system. As for the Lee and Kuan et al.filters, the local A Priori mean and variance, <R> et σR2, are estimated using the first-order non-stationary local statistics of the original image I(cf.Equations (26) and (27)), thus allowing adaptivity to the local isotropic properties of the scene R. At the scale of a terrain parcel, the MAP filter must theoretically restore the whole pdf of the radar reflectivity. This includes its mean E[R] and its variance (equivalent to parameter α), as well as the spatial relationships between pixels if one introduces speckle and scene second-order statistics in the filtering process (cf.[10] [75] [67]).

Nevertheless, in the presence of structures, the NMNV model may be used abusively on an inappropriate neighbourhood of the pixel under processing and therefore the local statistics <R> and σRare no longer those of an ergodic process. Thus, in order to use the filter in the conditions required by the NMNV model, the statistical situation it has been designed for, one must introduce, as one might also do for the Lee and Kuan et al.filters, structure detectors. These detectors contribute to select a neighbourhood of the pixel under processing composed of pixels pertaining to the same thematic class. They may be, either geometrical improvements based on the use of ratio (RoA) detectors ([7] [66] [63] [64] [10]), or detectors based on the use of the second-order statistics of both the speckle and the scene ([67]; cf.Section 4). The desired result is a better preservation of image texture, structures, and responses of strong scatterers in the speckle filtered radar image.


4. Using speckle and scene spatial correlation to preserve structural information in SAR images

The assumptions usually made by adaptive speckle filters ([76] [30] [32] [10]) with regard to the first order statistics of the speckle are:

  • the multiplicative speckle model;

  • the pdf of the speckle which reflects only incompletely the spatial correlation of the speckle through the ENL.

  • The assumptions regarding the imaged scene are:

  • Local wide-sense stationarity and ergodicity;

  • A formal Non-stationary Mean Non-stationary Variance (NMNV) model accounting for the spatial variation of the first order statistical properties [32];

  • For the Maximum A Posteriori (MAP) filter, the general form of the pdf of the scene radar cross-section [10].

The main consequences of these assumptions are that:

  • The formal NMNV model justifies local treatment (processing window) using the local statistics;

  • Assuming spatially uncorrelated speckle leads to the design of scalar (single-point) filters;

  • Local adaptivity to tonal properties (mean intensity <I> and scene radar reflectivity R) is achieved by controlling the filtering process through the local coefficient of variation CRof the radar reflectivity Rof the scene, through the local coefficient of variation CIof the original SAR image intensity I[27] [28].

The major drawbacks of these assumptions are that:

  • The coefficients of variation are statistically [10] sensitive to texture and speckle strength [40], but do not provide direct information on spatial correlation properties and texture directionality. Locally estimated, they can also be biased if speckle samples are not independent (cf.§ 2.2.2), as it is often the case in SAR images. Indeed, speckle can have strong spatial correlation properties, a situation the filters are not designed to deal with;

  • Compliance with the ergodicity and stationarity hypothesis requires a preliminary identification of the structural elements of the scene to fully exploit the NMNV model and to correctly estimate the local mean reflectivity around the pixel under processing;

In this Section, local second order properties, describing spatial relationships between pixels are introduced into single-point speckle adaptive filtering processes, in order to account for the effects of speckle and scene spatial correlations. To this end, texture measures originating from the local autocorrelation functions (ACF) are used to refine the evaluation of the non-stationary first order local statistics as well as to detect the structural elements of the scene [75].

4.1. Effects of the spatial correlation on the speckle

In practice, the usual single-point filters do preserve texture, only due to the spatial variation of the local first order statistics, using the NMNV formal model [32]. Fairly good preservation of textural properties and structural elements (edges, lines, strong scatterers) can be achieved by associating constant false alarm rate structure detectors such as the directional Ratio-Of-Amplitudes (RoA) detectors; ([7] [66] [10] [63] [64]) to the speckle filtering process. The combination of detection and speckle filtering allows to preserve scene structures, and to enhance scene texture estimation on a shape adaptive neighborhood of the pixel under consideration. When the conditions for which they have been developed (especially spatially uncorrelated or low-correlated speckle) are fulfilled, the best single-point adaptive speckle filters and their structure retaining associated processes based on RoA detectors retain enough scene texture to allow its use as an additional discriminator ([23] [24] [63] [64]).

However, the performances of the usual single-point filters, even when refined with associated RoA structure detectors, degrade when the actual spatial correlation of speckle samples becomes significant [77], i.e.generally when SAR images are far too much oversampled [31]. In fact, the uncorrelated speckle model approximation is seldom justified, according to sampling requirements of the radar signal and to the SAR system and processor which have produced the image. In multilook SAR images, spatial correlations are introduced by a series of weighting functions: complex coherent weighting related to data acquisition (antenna pattern, pre-sum filter and Doppler modulation in azimuth), coding of transmitted pulses in range, and looks formation (useful Doppler bandwidth selection, extraction window in azimuth). In practice, they are related to sampling and scale requirements in the design of SAR image products, up to the extreme case of signal oversampling [31].

As an example of the effects of speckle spatial correlation, the case of structure detection using the RoA detectors [10] implemented in the Gamma-Gamma MAP filter is illustrated in [67]. Detection is successfully performed when the correlation of speckle is low between pixels. However, the performance of RoA detectors substantially degrades when the spatial correlation of speckle samples becomes significant. It has been shown theoretically in [67] that spatial speckle correlation results in an increasing detection of non-existing structures, the increase of the probability of false alarm (Pfa). The expected 1% Pfa for uncorrelated speckle raises to about 25% if the correlation of adjacent pixels is 0.7 in both range and azimuth (cf.ERS: about 0.55 in range and 0.65 in azimuth). In addition, the spatial speckle correlation also decreases the probability of detection. The combination of these two effects causes artefacts in the filtered image (visual “crumbled paper” effect), making photo-interpretation uneasy and reducing substantially the performance of further automatic image processing like segmentation [78] and classification [23] [24].

Therefore, second order statistical properties must also be considered, including both scene texture and resolution/sampling related properties, to achieve a more complete restoration of the radar reflectivity.

4.2. A possible solution: Spatial vector filtering

A possible solution is the implementation of vector (multiple-points) filters ([32] [79]) where the spatial covariance matrices of the speckle and of the scene are taken into account. The development of a filtering method using second order statistics, estimation of the complete ACF of the scene through the speckled image ACF, results in a vectorial equation giving rise to a multiple-points filter, such as the LMMSE vector filter developed theoretically by Kuan et al.[32]. A practical implementation of multiple-point filters has been first suggested by Lopès et al.[10], and then detailed by Lopès & Séry [79]. Such filters have good spatial memory, but they result in a complicated implementation and very heavy computations.

4.3. Single-point speckle filtering using spatial second order statistics

To avoid the mathematical complexity and the heavy computational burden of multiple-point filters, an alternative solution is to introduce an appropriate description of both speckle correlation properties and spatial relations between resolution cells into a single-point filter.

Second order statistics have explicitly been used in the past, following the scheme proposed by Woods & Biemond [80] implemented later by Quelle & Boucher [81] in the adaptive Frost speckle (single-point) filter [30]. The filter, which belongs to the family of Wiener filters, is established using Yaglom's method [38], where a frequency characteristic is determined, instead of an impulse transfer function, which actually exists only within wide-sense stationary areas. This frequency characteristic is determined in the speckled image at a distance of one pixel only in range and azimuth, thus mainly limiting the description of correlation properties to those of the speckle.

Considering the NMNV model [32], a much better solution is to introduce locally estimated second order statistics to refine the computation of the local NMNV first order statistics. In this way, speckle correlation properties can be taken into account for filtering. In addition, the local mean radar reflectivity E[R] would be estimated, taking better account of the textural properties of the imaged scene (natural variability of the radar reflectivity, including spatial relations between resolution cells). On a given textured class of the scene, when scene correlation length is smaller than the processing window size, the mean radar reflectivity E[R] will then be modulated by the scene ACF, a function of the correlation coefficients of the radar reflectivity in all possible directions. When the correlation length of the imaged scene is longer than the processing window size, the estimate for the non-stationary mean radar reflectivity tends towards to the classical Maximum Likelihood estimate E[R]=<I>.

4.4. Local ACF’s and texture fields

Spatial relations between pixels are well described by the intensity ACF, defined by a set of correlation coefficients z), where Δz=(Δa, Δr), the normalised ACF of an intensity SAR image, {ρI(Δz)}, is a composition of the underlying scene ACF, {ρR(Δz)}, convoluted with an overlap function depending on the point spread function (PSF).

The local estimates for the normalised intensity ACF, {ρI(Δz)}, are, either deduced from the spatial autocovariance, which is computed on the domain of interest D(Npairs of pixels separated from each other by Δz) as follows:

Co^vI(Δz) = (1/N.D[I(z+Δz)<I>] .[I(z)<I>]E63

or, equivalently, directly computed within the domain of interest Das follows:

ρ^I(Δz) =D([I(z+Δz)<I>].[I(z)<I>])D[I(z)<I> ]2E64

These estimates have the attractive property that their mean square error is generally smaller than that of other estimators; considerations on the estimation accuracy can be found in Rignot & Kwok [82].

The contribution of scene texture must be separated from that of the speckle for all displacements Δz: the scene ACF is deduced from the intensity image ACF, using the following transformation [28]:

ρ^R(Δz) = [ 1 + ρ^I(Δz).CI2] [1 + |Gc(Δz)|2/L] E65

where Lis the equivalent number of independent looks, CIis the local coefficient of variation of the image intensity, and Gc(Δz)is the normalised ACF of the individual 1-look complex amplitudes (i.e.the ACF of the SAR imaging system).

Gc(Δz)depends only on the SAR complex PSF [83]. It is realistic to adopt an exponential form for the ACF of the speckle for the separated looks detected in intensity:

|Gc(Δz)| = exp [(Δa.dx/rx+Δr.dy/ry)]E66

where dx, dyare the pixel dimensions, and rx, ryare the spatial resolutions in azimuth and range directions. However, it is preferable to use directly the true ACF of the speckle in the actual multilook SAR image, if it is available, or if its estimation is possible.

The computation of a local non-stationary estimate of E(R) is performed by a simple convolution in the domain of interest around the pixel under consideration of the normalised ACF of the textured scene {ρR(Δz)}with the intensity:

E^(R) = D[ρ^R(Δz).I(z+Δz)]D[ρ^R(Δz)]  E67

for all ρR(Δz)within D, and as long as the following condition applies:

[ρ^R(z) /ρ^R(z=(0,0))]   >  1/eE68

This means that, along a direction z, the pixels whose correlation ρ^R(z+Δz)to the pixel under process is less than 1/e are ignored for the computation of the first order NMNV local statistics through Equation (67). This way, the algorithm achieves powerful speckle reduction in homogeneous areas, preservation of the textural properties (heterogeneity and directionality) wherever they exist, and correct structure detection.

The implementation of this operator, which estimates the local non-stationary E[R] is inspired from the notion of "texture fields" introduced by Faugeras & Pratt [84]. It is important to notice the following properties of the operator:

  1. in homogeneous areas, Equation (67) acts as a smoothing operator: E[R]<I>, which is the Maximum Likelihood estimate of E[R] for a perfectly homogeneous neighbourhood. (The local energy scaling factor in the denominator of Equation (67) ensures the preservation of the absolute level of the radar reflectivity).

  2. for a strong scatterer, it acts as a Laplacian: E[R]I(z).

  3. if a perfect step edge between two homogeneous areas is present in the neighbourhood (correlations=1,0,-1), it acts as a gradient and automatically stops weighting further pixels. This indicates structure detection capability;

  4. in all other cases, i.e.for a textured scene, the E[R] estimate is a weighted average of the neighbourhood, modulated by the scene ACF in all possible directions.

This concept of "texture fields", based on the analysis of second order statistical properties, generalises image processing concepts such as filtering and detection, which are usually considered as distinct. It provides the possibility to perform simultaneously structure detection and image pre-filtering. Structural elements such as edges are detected (value of the estimate of the correlation coefficient becoming less than 1/e) when the distance from the pixel under processing (central pixel of the processing window) to the edge is reached, going in the direction of the edge. Therefore, the position of the edge is known as soon as the edge enters the “field of view” of the processing window.

The local mean radar reflectivity estimation takes now into account the textural properties of the imaged scene (spatial variability of Rbetween resolution cells, directionality) as well as the spatial correlation properties of the speckle. Introduction of the non-stationary estimates of E[R] and CRinto the scalar equation of a single-point speckle filter improves the restoration of Rfluctuations in the filtered radar image [75]. As speckle and scene correlations are explicitly accounted for, this method can be considered as being complete.

The whole process acts as an adaptive focusing process in addition to the filtering process, since emphasis is put on the restoration and the enhancement of the small-spatial-scale fluctuations of the radar reflectivity. Full profit is taken of the useful resolution offered by the compound sensor and processing systems. Thus, speckle related features and thin scene elements (scene short-scale texture and thin structures) are automatically differentiated. The latter are restored with enhanced spatial accuracy, according to the local statistical context.

There is no geometrisation of the scene structural elements, as when using templates-based detectors used in structure retaining speckle filters ([76] [63] [64] [10]). The neighbourhood on which the NMNV local statistics are estimated is delimited by the estimates of ρR(Δz)for all possible Δzorientations. The only limitation to the potential extension of the domain Dis the spatial extent of the available ACF of the speckle.


5. MAP speckle filters for series of detected SAR images

Since the launch of ERS-1 satellite in 1991 temporal series of calibrated SAR images have been made available. This has stimulated the development of multichannel ("vector") filters especially dedicated to filter the speckle in a series of images, taking into account the correlation of the SAR signal between image acquisitions, thus opening the way to furher developments of change detection techniques specific to SAR series of images.

Interest in multichannel adaptive speckle filtering arises from their ability to combine both multi-image diversity through the exploitation of correlation between speckle and scene between images, and spatial diversity through the locally adaptive averaging of pixel values in the spatial domain. The objectives of this combination are the simultaneously achievement of both, a better restoration of existing texture, and a stronger speckle smoothing in textureless areas.

Although numerous work with regard to multi-channel speckle filtering in SAR images had already being done in a remote past, the introduction of A Priori knowledge or A Priori guess which implies the use of Bayesian methods in the processing of multi-SAR’s images, multi-date SAR images or SAR and optical images began only in 1992. Bayesian MAP vector speckle filters developed for multi-channel SAR images incorporate statistical descriptions of the scene and of the speckle in multi-channel SAR images. These models are able to account for the scene and system effects, which result in a certain amount of correlation between the different channels.

To account for the effects due to the spatial correlation of both the speckle and the scene in SAR images, estimators originating from the local ACF’s (cf.Section 4) are incorporated to these filters to refine the evaluation of the non-stationary first order local statistics on an ergodic wide-sense neighbourhood of the pixel under processing. The goal is to improve the restoration of scene textural and structural properties and the preservation of the useful spatial resolution in the filtered SAR image.

5.1. Multichannel vector MAP speckle filtering

In the case of multi-channel detected SAR images, let define the vector quantities of interest: Iis the speckled intensity vector available in the actual SAR data; Ris the radar reflectivity vector, which is the quantity to restore. The MAP filtering method bases on Bayes’s theorem [47] [48] in its matricial form:

P(R/I)= P(I/R) . P(R) / P(I)E69

where P(I/R) is the joint conditional pdf of the speckle.

For each individual detected SAR image i, P(Ii/Ri) is known to be well approximated by a Gamma distribution for multilook intensity images (cf.§ 2.1.4; Equation (14)).

P(R) is the joint pdf of the radar reflectivity, introduced as statistical A Priori information in the restoration process. To describe the first order statistical properties of natural scenes, it has already been shown that, in most situations, for each individual image i, a Gamma distribution (cf.§ 3.3.4) is an appropriate choice for P(Ri). Nevertheless, when in doubt about the reliability of such a choice, P(Ri) can be estimated directly from the data in the corresponding radar image i.

For multi-channel detected SAR images, MAP filtering is a vector filtering method. For every channel i, the posterior probability is maximum if the following condition is verified:

[Ln(P(I/R))]/Ri+[Ln(P(R))]/Ri  = 0      for      Ri= R^i MAPE70

The MAP speckle filtering process acts as a data fusion process, since the information content of the whole image series is exploited to restore the radar reflectivity in each individual image. Among other advantages, this allows a better detection and restoration of thin scene details.

5.2. Multichannel scene models

To describe the first order statistical properties of a natural scene viewed at the spatial resolution and at the wavelength frequency of a radar remote sensing sensor, it has been shown that a Gamma pdf would be a suitable representation (cf.§ 3.3.4).

However, to describe the first order statistical properties of a natural scene as viewed by diverse radar sensors (different physical scene features), or at different dates (scene evolution over time), there is no analytic multivariate Gamma pdf available under closed-form. Therefore, for the sake of mathematical tractability, a multivariate Gaussian pdf is used as analytic "ersatz" of a multi-channel scene statistical model:

P(R) =[(2π)N|CovR|]1/2. exp[(R<R>)t.CovR1.(R<R>)]  E71

CovRis the local covariance matrix of the imaged scene in all image channels. For each pixel location, CovRis estimated from the local covariance matrix of the intensities observed in all image channels CovI, using the multiplicative speckle model (cf.§ 2.2.2; [27] [28]).

This statistical scene model is convenient to preserve the mathematical tractability of the problem. In addition, the Gaussian model, commonly used to describe the statistical properties of natural scenes in processing of optical imagery, is appropriate as a first-order statistical description of moderate, not too heterogeneous, textures.

However, when scene texture becomes very heterogeneous, this multivariate Gaussian scene model is not any longer reliable, and in the absence of specific prior knowledge about the distribution of scene texture, its model has to be directly estimated from the SAR image.

5.3. Speckle uncorrelated between SAR images: The Gamma-Gaussian MAP filter

Let consider the case of a series of Nimages originating from different SAR systems (different frequencies, incidence angles, or spatial resolution for example, but same spatial sampling). In such a case, it is justified to consider that speckle is independent between the Nimage channels.

Under this assumption, since each one of the SARs senses slightly different physical properties of the same scene, the joint conditional pdf of the fully developed speckle, P(I/R), can reasonably be modelled as a set of Nindependent Gamma distributions P(Ii/Ri):

P(Ii/Ri) =(Li /Ri)Li /Γ(Li). exp(LiIi /Ri) .Ii(Li1)E72

where the Liparameters are the ENLs of the individual SAR images.

The first MAP filtering algorithm for multi-channel detected SAR images results in a set of Nscalar equations, with Nindependent speckle models and a coupled scene model (Equation (71)). The coupled scene model is justified by the fact that the physical properties of the scene do contribute, but not in an identical way, to the formation of the images provided by the different sensors. The set of equations describing the Gamma-Gaussian MAP filter for multi-channel detected SAR images is as follows ([73] [85]):

Li(Ii/Ri21/Ri) (1i)t.CovR1.(R<R>)(R<R>)t.CovR1.(1i)  1/2 Tr[CovR1.CovR/Ri]+(R<R>)t.CovR1.CovR/Ri.CovR1.(R<R>) = 0E73

where (1i) is a vector where all components, but the ith, are equal to zero. The set of Equations (73) is solved numerically.

The introduction of coupling between the scene statistical representations leads to a data fusion process taking profit of the correlation between texture as it is observed in all the images in the series. Indeed, replacing speckle noise model by optical noise (or film grain noise, or any other appropriate noise model) in some of the Nimage channels, this filter adapts also easily to the case of multi-channel optical and SAR images.

5.4. Speckle correlated between SAR images: The Gaussian-Gaussian MAP filter

To filter series of images originating, either from the same SAR system operating in repeat-pass mode, or from different SAR systems with relatively close properties, a second speckle filtering algorithm has been developed. The different SARs may be close in terms of frequency, angle of incidence, spatial resolution and sampling, geometry of images, with difference in polarisation configuration only, or small differences in incidence angle, for example. In such cases, speckle correlation between individual SAR images must be taken into account to deal optimally with system effects on the SAR images series.

Taking into consideration speckle correlation between image channels, the joint conditional pdf of the fully developed speckle P(I/R) should theoretically be a multivariate Gamma pdf. Nevertheless, since there is no analytic multivariate Gamma pdf available in closed-form, another reasonable choice for P(I/R) should be done for the sake of mathematical tractability.

To solve this problem, Lee’s assumption [33] is adopted: for multilook SAR images (more than 3 looks), the joint conditional pdf of the speckle can be reasonably approximated by a Gaussian distribution. Therefore, for convenience, P(I/R) is approximated by a Gaussian multivariate distribution in the case of multilook SAR images:

P(I/R) = [(2π)N|CovS|]1/2. exp[(IR)t.CovS1.(IR)]E74

where CovSis the covariance matrix of the speckle between images of the series. In practice, CovSis estimated within the most homogeneous area one can identify in the series of detected SAR images.

It is noteworthy that CovStakes into account possible different Livalues of the respective ENLs of the SAR images in the series as well as energy unbalance between images. This enables to compose the series combining images from different SARs operating at the same wavelength, or acquired by the same SAR at different times or in different configurations of polarization.

At this point, the second MAP filtering algorithm for multi-channel detected multilook SAR images results in a set of Nscalar equations. The set of equations describing the Gaussian-Gaussian MAP filter for multi-channel detected multilook SAR images is ([73] [85]):

t(1i).CovS1.(IR) +(IR)t.CovS1.(1i)  1/2 Tr[CovR1.CovR/Ri]+(R<R>)t.CovR1.CovR/Ri.CovR1.(R<R>)(1i)t.CovR1.(R<R>) t(R<R>).CovR1.(1i) = 0E75

where (1i)is a vector where all components, but the ith, are equal to zero. Again, this set of equations is solved numerically.

Including a coupled speckle model (Equation (71)) and a coupled scene model (Equation (74)), this speckle filter takes profit of both speckle and scene texture correlations, thus restoring the radar reflectivity through a complete data fusion process.

5.5. Prior knowledge gathered from SAR image series: DE-MAP filters

In the presence of very strong or mixed textures, eventually combined with the presence of strong topographic relief, it may be hazardous to make an A Priori assumption about the pdf of the radar reflectivity. This situation is often the case in SAR images of dense tropical forest or of urban environments located in slanted terrain with rapidly varying slopes and counter-slopes.

This is also often the case in high-and very-high spatial resolution SAR images, where strong texture is omnipresent. In addition, at such very fine scale, the textural properties of the scene vary strongly within the timeframe elapsed between successive image acquisitions by the same SAR, and exhibit strong differences when imaged by different SARs or using different SAR configurations in terms of wavelength, polarisation, or angle of incidence, for example.

5.5.1. Speckle correlated between SAR images: The Gaussian-DE MAP filter

As mentioned precedently (cf.§ 5.4), for a series of Nmultilook (Li>3 for each image i) detected SAR images, the joint conditional pdf of the speckle, P(I/R), can be modelled as a multivariate Gaussian distribution [33] when speckle is correlated between image channels with covariance matrix CovS:

P(I/R)=[(2π)N.|CovS|]1/2. exp[t(IR).CovS1.(IR)]E76

The entropy of the scene texture [68] [69] [70] is introduced (on a sample of Npixpixels of an ergodic neighbourhood of the pixel under processing) as follows:

S(Ri) = kNpix[Rik. log(Rik)]     for the ithchannelE77

Because the radar reflectivities Rkare non-negative and exp(S(Ri))/Zis normalized, and since S(Ri) is a measure of the spread/dispersion of radar reflectivities in image channel i, Equation (77) can be treated as a pdf whose entropy is S(Ri) [70]:

P(Ri) = (1/Z. exp(S(Ri)) = (1/Z) . exp(kNpix[Rik.log(Rik)] )E78

Under this assumption, the Gaussian/Distribution-Entropy MAP (Gaussian-DE MAP) filter for multi-channel multilook detected SAR images (Nimage channels) results in a set of Ncoupled scalar equations of the form [71]:

(1i)t.CovS1.(IR)+(IR)t.CovS1.(1i) Ri2.kNpix[log(Rik)1/Ln(10)]= 0E79

where (li) is a vector where all components, but the ith, are equal to zero.

To estimate P(Ri), the radar reflectivities Rikin a neighbourhood of the pixel under processing in the ith SAR image are pre-estimated in a first speckle filtering pass; indeed, the iterative estimation process does already converge after a first application of Equation (79).

Note that this filter does not need to introduce implicitly scene texture correlation. Indeed, the restoration of the radar reflectivity in each image channel bases only on the A Priori knowledge retrieved from the speckled image channel. Nevertheless, the filter takes profit of speckle correlation between image channels. An application of this filter is illustrated in Figure 1.

Figure 1.

TerraSAR-X 4-looks SAR images (HH and VV polarisations) acquired on August 10, 2007, over Oberpfaffenhofen, Germany (Credits: Astrium and InfoTerra; filtered images: Privateers NV, using ®SARScape).

5.5.2. Speckle correlated between SAR images: The Gamma-DE MAP filter

Note that if speckle is not correlated between the Nimage channels, the Gamma/Distribution-Entropy MAP (Gamma-DE MAP) filter for multi-channel detected SAR images(Nchannels) results in the resolution of a set of Nindependent (uncoupled) scalar equations similar to Equation (62) (cf.§ 3.3.5).

5.6. Vector MAP filtering and control systems

A number of advantages and improvements these filters offer is listed in this section. Most of these advantages arise from the use of the covariance matrices of the speckle and of the imaged scene (i.e.from the use of coupled models) in a Bayesian reconstruction of the radar reflectivity.

5.6.1. Imaging system’s effects and preservation of high-resolution

In a series of SAR (but not only...) images, the resolution cells never overlap perfectly between the different individual images. The covariance matrix of the speckle CovScontains the information about these overlaps, and the filtering process accounts for these overlaps in order to preserve, and even improve the spatial resolution in the filtered SAR images, thus allowing us to restore the thinnest scene structures.

In addition, the simultaneous attempt to detect images structures and targets in all radar image channels improves the probability to detect such structures. As shown in Figure 1, such an improvement is already very substantial using only two SAR images.

5.6.2. Potential for change detection

In series of radar images acquired over time, the covariance matrix of the scene CovRcontains the information about the temporal evolution of the imaged area and enables the detection of changes in scene physical properties (scene temporal evolution).

In a series of SAR images acquired by (not too...) different SAR’s, the covariance matrix of the scene enables the detection of more aspects of the scene, as it is viewed by more SAR sensors.

5.6.3. MAP speckle filters and control systems

Let consider the mathematical form of the set of Equations (73), (75) and (79). Indeed, their formulation is that of a control system, since both equations can be rewritten under the form of Riccati’s [86] continuous time algebraic matricial equations:

 A . X  X .tA  Q + X .tC . P1.C . X = 0E80

Equation (10) represents the optimal state controlled reconstruction at constant gain of linear invariant processes (Rand the textures of the channels) perturbed by white noises (speckle, pixel spatial mismatch between channels), from observed evolving state variables (Iand CovI) [87]. The scene A Priori model acts as a command, and the covariance matrices act as controls. It is also noticeable that the MAP elaboration of information generates automatically the feedback process, which enables to control the filtering process.

In addition, it is highly remarkable that Riccati’s theorem [86] stipulates the existence of a unique positive definite solution for Equation (10). Therefore, this property holds also for the MAP Equations (73), (75) and (79).

The extension of this noise filtering technique to the case of SAR and optical images set is straightforward, as mentioned in § 5.3. It is also noticeable that, extending this method to the case of complex SAR images, this technique presents a very interesting potential for super-resolution. Such methodologies can be found in the literature (cf.[46], for example).

A major interest of control systems is that they offer wide possibilities for the choice and design of additional commands (statistical and physical models) for further data exploitation. In this view, speckle filtering should be regarded as a first step of integrated application oriented control systems.


6. MAP speckle filters for complex SAR data

SLC images are complex valued images. Hence, every image pixel is assigned a complex value, its complex amplitude Ac=i+j.q=I.exp(j.ϕ), which corresponds to an intensity I=i2+q2 and a phase ϕ=Arctg(q/i).

Whereas the intensity Icarries information related to the radar reflectivity Rand the backscattering coefficient σ° of the imaged scene, the phase of a complex pixel of a single SLC image is totally random and does not carry any information about the scene when speckle is fully developed. (It is noteworthy that only, phase differences between SLCs acquired in interferometric conditions, or phase differences between configurations of polarisation in SAR polarimetric data, carry information about the imaged scene).

Therefore, the quantity of interest to restore in SLC images through speckle filtering is the radar reflectivity R. Thus, whereas the input of the speckle filters consists in complex radar data, their output consists in radar reflectivity images, i.e.speckle filtered intensity images.

6.1. Single-look complex (SLC) SAR image: The CGs-DE MAP filter

In SLC radar data the spatial correlation of the speckle can be dealt with by taking into account the spatial covariance matrix of the speckle in the filtering process. Indeed, the joint pdf of a local sample X(Npix-dimention vector) of 1-look spatially correlated speckle is given by Equation (7) (cf.§ 2.1.2), can be written as [88] [89] [90] [91]:

Pu(X) = P(X/R) = [1/(πNpix.|CS|)] .exp [Xt*.CS1.X]E81

where CSis the spatial covariance matrix of the complex speckle, estimated within the most homogeneous part of the SLC radar image.

In high and very-high spatial resolution SLC SAR images, which often exhibit strong textural properties, it becomes difficult to consider any theoretical statistical model as a reliable A Priori knowledge about scene texture. In such a situation, it seems reasonable to retrieve statistical scene knowledge directly from the SAR image data [68] [69]Introducing the local entropy of the radar reflectivity as A Priori scene knowledge [70] (cf.§ 3.3.5), the Complex-Gaussian/Distribution-Entropy MAP (CGs-DE MAP) filter for an SLC radar image is the solution of the following equation [71]:

(1/Npix) . tX*.CS1.X R  R2.k=1Npix[log(Rk)1/Ln(10)] = 0 E82

where Npixis the effective number of neighbour pixels taken into account in the computation of local statistics. The radar reflectivites Rkare pre-estimated in this neighbourhood of the pixel under processing by a first pass of speckle filtering in intensity.

An example of application of this filter is shown in Figure 2.

Figure 2.

RADARSAT-2 1-look SAR image (HH polarisation) acquired on April 30, 2009, over Capetown, South Africa (Credits: Canadian Space Agency and MDA; filtered images: Privateers NV, using ®SARScape).

6.2. Separate complex looks, or series of SLC SAR images

Interest in adaptive speckle filtering of SLC image series arises from their ability to combine both temporal or spectral diversity, and spatial diversity. The objective is to obtain a multilook radar reflectivity image, where additional speckle reduction is obtained in the spatial domain by averaging pixel values through locally adaptive filtering.

The final objective is to obtain very strong speckle smoothing in textureless areas and to reach an ENL high enough (i.e.ENL>>30, and preferably ENL>150) to enable information extraction using classic techniques, or combined use with optical imagery, while better restoring existing or characteristic scene texture.

6.2.1. Compount statistical speckle model for series of SLC SAR images

A SAR sensor can produce series of SLC SAR images in three different situations:

  1. Lseparate complex looks are usually extracted from the useful Doppler band of the radar signal, usually to produce a L-looks multilook image, but the individual looks may be conserved at this stage of SAR image production. Let remind that this operation is done at the cost of a loss in azimuth spatial resolution (cf.§ 2.1.3). These looks generally overlap, in the sense that they have in common a part of the Doppler spectrum of the original full-resolution 1-look signal, in order to preserve as completely as possible the total amount of backscattered radar signal. This overlap introduces a certain amount of speckle correlation between looks extracted from adjacent portions of the Doppler spectrum.

  2. In the case of a non-interferometric 1-look complex SAR data set, the individual SAR SLC images acquired independently of each other show low or inexistent speckle correlation between images.

  3. However, in the case of interferometric 1-look complex SAR data, speckle correlation is close to 1 (corresponding to perfect correlation) or very high, since speckle correlation between SAR data acquisitions is exactly the fundamental requirement for interferometry.

As a consequence of these three possible configurations of a series of SLC images, a speckle filter dealing with the whole series must be able to take into account the complex correlation of speckle between the different images. The local complex covariance matrix CovSdescribes completely both, the radar intensities in each SLC image, and the correlation between images.

Considering the whole set of SLC images, the measurement vector for each pixel is X={yn}, where yn=in+j.qn. When speckle is fully developed, the (in, qn) are statistically independent random processes. However, the ynare correlated complex Gaussian random processes with a joint conditional pdf given by [88]:

P(X/CovS) = exp(tX*.CovS1.X) / (πL.|CovS|)E83

where CovS=E[tX. X*] is the complex covariance matrix of the speckle between SLC images, |CovS| is the determinant of CovS, and (CovS-1/|CovS|) is the inverse complex correlation matrix of the speckle. Theoretically, it is possible to compute the correlation matrix of the speckle from the sensor’s and SAR processor’s parameters. But, since these parameters are in general not easily available to SAR images users, CovSis estimated in practice over the most homogeneous area identified in the images.

6.2.2. Gamma distributed scene model: The Complex-Gaussian–Gamma MAP filter

In the case of separate complex looks corresponding to the same SAR data acquisition, or of series of SLC images acquired over time by the same SAR over a time-invariant scene, it is quite reasonable to assume that the textural properties of the scene are similar in all the SLC’s.

Therefore, scene texture may be locally modelled for all SLCs in the series, as well as in the detected multilook image formed by incoherent pixel-per-pixel averaging, by the same Gamma pdf with parameters E[R] and α(cf.§ 3.3.4, Equation (52)). Locally, E[R] is estimated, either by <I> or using Equation (67), and αis estimated by 1/CR2 (cf.§ 2.2.2, Equations (26) and (27)), in an ergodic wide-sense stationary neighbourhood of the pixel under processing in the (less noisy) multilook radar image.

In this situation, the Complex-Gaussian/Gamma MAP (CGs-Gamma MAP) filter for separate complex looks (L1-look SAR images) has the following expression [11]

R^=<I>.(αL1) + <I>2.(αL1)2+ 4α.<I>.(tX*.(CovS1/|CovS|).X) 2αE84

Note that when the SLC images or the separate looks are independent, i.e.when the non-diagonal terms of CovSare equal to zero, then tX*.(CovS-1/|CovS|).X=L.<I>, and the filter is equivalent to the Gamma-Gamma MAP filter for single detected SAR image (cf.§ 3.3.4). Therefore, independent SLC images can be incoherently summed to produce an L-looks intensity image that is filtered using the Gamma-Gamma MAP filter (Equation (54)) with the same α and Rlocal parameters, for the same final result.

Finally, it is noteworthy that the ENL L’of the multilook radar image can be computed from the complex correlation coefficients ρmnbetween images mand n, computed from the covariance matrix of the speckle, i.e.CovSestimated over the most homogeneous area identified in the images [11]:

L = L . (1 +(2/L) .n=1L1m=n+1L|ρmn|2)1E85

This relationship may prove useful if SLC images are incoherently summed to produce a L’-looks intensity image.

6.2.3. Prior knowledge gathered from the data: The Complex-Gaussian-DE MAP filter

High and very-high spatial resolution SLC SAR images often exhibit strong textural properties. As exposed above, in such a situation, statistical scene knowledge is estimated from the radar data [68] [69]Introducing the local entropy of the radar reflectivity as A Priori scene knowledge [70] (cf.§ 3.3.5), the Complex-Gaussian/Distribution-Entropy MAP (CGs-DE MAP) filter for an SLC image series (Limages) is the solution of the following equation [71]:

tX*.CovS1.X L.Ri Ri2.k=1Npix[log(Rik)1/Ln(10)] = 0E86

where CovSis the covariance matrix of the speckle between the different SLC images [11] and Npixis the effective number of neighbour pixels taken into account in the computation of local statistics. The radar reflectivites Rkare pre-estimated in this neighbourhood of the pixel under processing by a first pass of multichannel speckle filtering of the Limages in intensity (using Equation (84) for example).

An application of this filter to a series of 18 (3x6) SLC SAR images is shown in Figure 3.

Figure 3.

Series of 18 (3x6) COSMO-SkyMed 1-look spotlight SAR images (HH and VV polarisations) acquired between February 2 and September 30, 2009, over Perito Moreno in Patagonia, Argentina (Credits: original images: Italian Space Agency; filtered images: Privateers NV, using ®SARScape).


7. MAP speckle filters for polarimetric radar data

7.1. Polarimetric radar data representation, and polarimetric speckle model

A polarimetric radar system produces, for each pixel location, a scattering measurement matrix, which is the scattering matrix Sof the corresponding surface (cf.Section 1, Equation (1)) corrupted by speckle noise.

The first polarimetric speckle filter ever developed [92] resulted in an optimal summation of the intensities IHH, IHV, IVVof the polarisation channels. The resulting image is called the "improved span" image:

ISpan= IHH+(1+|ρHHVV|2) .IHV. E(IHH)/E(IHV)+ IVV.E(IHH)/E(IVV)E87

where ρHHVVis the complex correlation coefficient between SHHand SVV.

However, Equation (87) clearly shows that the ENL achieved in the span image will barely reach 3, which remains largely insufficient to meet the noise reduction requirements of remote sensing applications. Besides, the physical interpretation of ISpan, which is a mixture of HH, VV, HV radar reflectivities is questionable. In addition, scene information contained in the polarimetric diversity of the radar measurement as well as in the complex correlations between polarimetric channels is lost in the process.

Lee et al.[93] developed a polarimetric optimally combining the polarisation channels to restore the radiometric information in the HH, VV, and HV channels. This polarimetric speckle filter has been shown to preserve polarimetric signatures [94].

Nevertheless, the ultimate aim of polarimetric speckle filtering is to restore the full polarimetric scattering matrix S, or at least the radar reflectivities in the HH, HV, VV configurations of polarisation, while achieving a very high ENL value (>100) in the filtered data. To achieve these objectives, fully polarimetric adaptive speckle filters take profit of both the entire polarimetric diversity (i.e.the full scattering matrix) and the spatial diversity (i.e.local statistics around the location under processing) [95].

7.1.1. Polarimetric covariance matrix

Under the assumption of reciprocity, SHV=SVH, the speckle-free (denoised) scattering matrix Scan be transformed into its covariance matrix CS, as follows:


This representation puts emphasis on the radar reflectivities RHH=|SHH|2, RHV=|SHV|2 and RVV=|SVV|2, as well as the covariances between HH, HV and VV configurations of polarisation.

If the polarimetric measurement is made by a monostatic radar system, the covariance matrix ΣS of the actual polarimetric radar measurement is the speckle-corrupted version of CS.

It is noteworthy that it is possible to produce multilook polarimetric data (L-looks) from the original 1-look data. For L-looks polarimetric data obtained by either spatial averaging (generally in azimuth) or spectral multilooking (cf.§ 1.1.4), ΣS multilook is defined by:

ΣS multilook = [ (1/L).1L(Spq.Spqt*) ]      for all combinations p,q of H,VE89

7.1.2. Polarimetric vector, degrees of coherence, phase differences

For convenience, the polarimetric covariance matrix is often expressed, without changing its overall information content, under the form of a real valued vector, called the "polarimetric feature vector", or "polarimetric vector" X:

tX = ( |SHH|2, |SHV|2, |SVV|2, Re(SHH.SVV*), Im(SHH.SVV*),Re(SHH.SHV*), Im(SHH.SHV*), Re(SVV.SVV*), Im(SVV.SVV*) )   E90

Two quantities of great interest in polarimetric SAR applications can be obtained from the polarimetric vector:

1) The Phase Differences ΔφKL-MN between radar returns in configurations of polarisation KL and MN (where K,L,M,N ∈ {H,V} are preferably considered in most applications), which are computed as follows:

ΔφKL-MN= Arg(SKL.SMN*)= Arctan[Im(SKL)/Re(SKL)]Arctan[Im(SMN)/Re(SMN)] E91

2) The (complex) Degrees of Coherences γKL-MN, which measure wave coherence between radar returns in two configurations of polarisation KL and MN, and express textural properties of an extended wide-sense stationary scene, are defined as follows:

γKL-MN= E[Re(SKL.SMN*)]+j. E[Im(SKL.SMN*)](E[|SKL|2].E[|SMN|2] E92

However, in the actual polarimetric radar measurement, the observations ΔφKL-MN of the phases differences, and the observations ρKL-MN of the complex degrees of coherence are the speckle-corrupted versions of Arg(SKL.SMN*) and γKL-MN, respectivly. Therefore, these quantities, which represent important contributions of radar polarimetry must also be speckle filtered, justifying the development of fully-polarimetric speckle filters that do not limit themselves to the restoration of radar reflectivities.

7.1.3. Polarimetric speckle model

In the case of polarimetric radar data, ΣS is the actually observed polarimetric covariance matrix, and CS is the unspeckled polarimetric covariance matrix, i.e.the quantity to restore through speckle filtering. In the reciprocal case, and for low look correlation, the conditional pdf of ΣS is a complex Wishart distribution of the form [72]:

P(ΣS/CS)=(det ΣS)L3.LL3π3.Γ(L).Γ(L1).Γ(L2).(det CS)L . exp[Tr(L.CS1.ΣS)]  E93

where Lis the ENL of the "improved span" image defined by Equation (87) [92].

7.2. Restoration of the whole polarimetric vector / covariance matrix

7.2.1. Gamma distributed scene model: The Wishart-Gamma MAP filter

Using physical backscattering models, assuming (as a rough approximation) that texture is identical in all polarizations, we get the following approximation [72]:

CS=μ . E[CS]E94

where μis the scalar texture parameter equal to the normalized number of scatterers within the resolution cell, and E(CS) is the mean covariance matrix [92].

P(μ) = αα/Γ(α) . exp(α.μ) . μα1     and     E[μ]= 1E95

The characterisation of scene heterogeneity by only one textural random variable μbases on the assumption that textural properties are identical in all configurations of polarisation. However, this has been shown to be inexact, both experimentally [96] and theoretically [97]. Therefore, μshould be regarded as an average textural characterisation of the scene observed in all possible configurations of polarisation. Besides the advantage of estimating α=1/ CR2 (using Equation (27)) in a less noisy multilook image is an additional justification to do it in the "improved span" image defined by Equation (87) [92].

Introducing the first-order statistical models for fully-polarimetric speckle (Equation (93)) and the polarimetric texture parameter μin the MAP equation (Equation (50), cf.§ 3.3.2), the fully-polarimetric Wishart-Gamma MAP filter restores the value of μas:

μ^=(αL1) + (αL1)2+ 4.α.L.Tr(E[CS]1.ΣS) 2αE96

where Tr(.) denotes the trace of a matrix.

Finally, the restored (speckle filtered) version of the covariance matrix CSis obtained using the maximum likelihood estimator described in [72]:

C^S= μ^ .E[CS]E97

7.2.2. Prior knowledge gathered from polarimetric SAR data: The Wishart-DE MAP filter

In high and very-high spatial resolution polarimetric SAR data, strong and/or mixed textural properties justify to estimate statistical scene knowledge from the data themselves, rather than assuming a A Priori theoretical model. Assuming that textural properties are identical in all configurations of polarisation, the entropy constraint on scene texture ([68] [69] [70]) becomes:

P(μ) = (1/μ). exp(k[μk.log(μk)]           and     E[μ] = 1E98

In this case, the complex Wishart/Distribution-Entropy MAP (CW-DE MAP) filter for polarimetric multilook SAR data is expressed as [71]:

L.Tr(E[CS]1.ΣS)  L.μμ2.k[log(μk)1/Ln(10)] = 0E99

E[CS] is obtained using the maximum likelihood estimator (Equation (97)) described in Lopès et al.[72].

Figure 4.

ALOS-PALSAR 6-looks polarimetric SAR imagery acquired on June 30, 2006 over Bavaria, Germany (Credits: JAXA and MITI; filtered images: Privateers NV, using ®SARScape).

The texture parameters μkin the neighbourhood of the pixel under processing are pre-estimated over an ergodic neighbourhood of the pixel under processing by a first speckle filtering pass using the Wishart-Gamma MAP filter (cf.§ 7.2.1).

Although these fully polarimetric MAP filters assume identical texture properties in the HH, HV, and VV channels, which has been shown both experimentally [96] and theoretically [97] inexact, they have nevertheless been shown to preserve polarimetric signatures [98]. Fully polarimetric speckle filtering is illustrated in Figure 4.



ACF: spatial AutoCorrelation Function (of a signal)

DE: Distribution-Entropy (statistical distribution)

ENL: Equivalent Number of Looks (of a radar image)

LMMSE: Linear Minimum Mean Square Error (estimation)

MAP: Maximum A Posteriori (statistical estimation)

MMSE: Minimum Mean Square Error (estimation)

NMNV: Non-stationary Mean Non-stationary Variance (model)

pdf: Probability Density Function (random variable distribution)

Pfa: Probability of False Alarm (detection)

PSF: Point Spread Function (of an imaging system)

RoA: Ratio Of Amplitudes (structure detector)

SAR: Synthetic Aperture Radar (imaging sensor)

SLC: Single-Look-Complex (radar image)

© 2014 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution 3.0 License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

How to cite and reference

Link to this chapter Copy to clipboard

Cite this chapter Copy to clipboard

Edmond Nezry (June 11th 2014). Adaptive Speckle Filtering in Radar Imagery, Land Applications of Radar Remote Sensing, Francesco Holecz, Paolo Pasquali, Nada Milisavljevic and Damien Closson, IntechOpen, DOI: 10.5772/58593. Available from:

chapter statistics

4264total chapter downloads

3Crossref citations

More statistics for editors and authors

Login to your personal dashboard for more detailed statistics on your publications.

Access personal reporting

Related Content

This Book

Next chapter

Large Scale Mapping of Forests and Land Cover with Synthetic Aperture Radar Data

By Josef Kellndorfer, Oliver Cartus, Jesse Bishop, Wayne Walker and Francesco Holecz

Related Book

First chapter

Plates Amalgamation and Plate Destruction, the Western Gondwana History

By Ngako Vincent and Njonfang Emmanuel

We are IntechOpen, the world's leading publisher of Open Access books. Built by scientists, for scientists. Our readership spans scientists, professors, researchers, librarians, and students, as well as business professionals. We share our knowledge and peer-reveiwed research papers with libraries, scientific and engineering societies, and also work with corporate R&D departments and government entities.

More About Us