M A S T E R’S T H E S I S
DANIEL SABEL
Scaling Information for Scatterometer Derived Soil Moisture Estimates
MASTER OF SCIENCE PROGRAMME Space Science
Luleå University of Technology
Department of Applied Physics and Mechanical Engineering Division of Physics
Scaling Information for
Scatterometer Derived Soil Moisture Estimates
Daniel Sabel
Master of Space Engineering Programme Lule˚a University of Technology, Sweden
October, 2006
Carried out at the
Institute of Photogrammetry and Remote Sensing, Vienna University of Technology, Austria
under the supervision of Prof. Wolfgang Wagner
and Dr. Klaus Scipal
This study is a part of the development within the SHARE1 project, which is one of the European Space Agency’s DUE Tiger projects. SHARE aims at enabling an operational soil moisture monitoring service for the southern and central regions of the African continent.
All data processing routines were implemented with ITT’s IDL (Interactive Data Language).
Geospatial analysis and production of figures were carried out with IDL, ITT’s ENVI (Envi- ronment for Visualising Images) and ArcMap from ESRI (Environmental Systems Research Institute).
I would like to express my gratitude to Prof. Wolfgang Wagner for giving me the opportunity to carry out this study at the
Institute of Photogrammetry and Remote Sensing (I.P.F.), Vienna University of Technology, Austria.
I am also greatful to Prof. Sverker Fredriksson at the
Department of Applied Physics and Mechanical Engineering, Lule˚a University of Technology, Sweden,
for his work as my examinator and for the excellent feedback on my writing. Special thanks go to my colleagues of the I.P.F. Microwave Remote Sensing Group2, in particular Dr. Klaus Scipal and Dipl.-Geogr. Carsten Pathe, for their support and confidence throughout the work.
Special thanks also go to M.Sc. Zoltan Bartalis, whithout whom I probably would not have ended up in Vienna in the first place.
I wish to thank my mother Anne, my father Sven-Olof and my sisters Emma and Linnea for always believing in me and supporting me, whatever the quest.
1http://www.ukzn.ac.za/sahg/share/
Soil moisture information is an important factor in the fields of hydrology, mete- orology and climatology. For hydrological applications, soil moisture data with a spatial resolution of 1 km are often requested. Current remote sensing methods are limited to a spatial resolution of 25 – 50 km. This study was motivated by the need to improve the resolution of soil moisture information. The aim was to produce and make a first evaluation of scaling information that could be used for validating spatial downscaling of scatterometer derived soil moisture information to achieve a resolution of 1 km.
The backscatter intensity from surface scattering is largely determined by the moisture content in the upper few centimetres of the soil. Based on the assump- tion that soil moisture is mainly driven by large scale atmospheric forcing, the local and regional soil moisture content is likely to be well correlated. Therefore, also the backscatter intensity on the two scales is expected to show similar scal- ing properties. By using time series of ENVISAT ASAR data, a measure of the temporal correlation between the backscatter intensities on the two scales can be calculated. For locations where a high temporal correlation is achieved, the regional soil moisture estimates can be used directly on the 1 km scale. Scal- ing information for a large part of the African continent was produced and the influence of land cover was analysed.
It is concluded that there is a high potential for downscaling scatterometer derived soil moisture information for agricultural regions, including cropland and pasture, and over range land, which includes savannas and land covered with grass or shrub.
Preface i
Abstract iii
Contents v
List of Acronyms viii
1 Introduction 1
1.1 Soil moisture information . . . 1
1.2 Purpose . . . 2
2 Theoretical background 3 2.1 Microwaves . . . 3
2.1.1 Definition . . . 3
2.1.2 Surface interaction . . . 4
2.2 Radar . . . 6
2.2.1 Applications . . . 6
2.2.2 Basic operation . . . 7
2.3 Synthetic Aperture Radar . . . 11
2.3.1 SAR focusing . . . 12
2.3.2 Speckle . . . 17
2.3.3 SAR geocoding . . . 19
2.3.4 The ScanSAR technique . . . 20
2.4 ENVISAT’s Advanced Synthetic Aperture Radar . . . 21
2.4.1 Platform . . . 21
2.4.2 Sensor . . . 23
2.5 Scales and scaling . . . 25
2.5.1 Definition of scale . . . 25
2.5.2 Process, measurement and model scale . . . 25
2.5.3 Definition of scaling . . . 26
2.6 Soil moisture . . . 27
2.6.1 Process . . . 27
2.6.2 In-situ measurement methods . . . 27
2.6.3 Observations from space . . . 28
3 Data sets and study area 31 3.1 Study area . . . 31
3.1.1 Extent . . . 31
3.1.2 Land cover distribution . . . 31
3.2 ASAR Global Monitoring Mode . . . 32
3.2.1 General description . . . 32
3.2.2 Geocoding quality . . . 34
3.3 USGS Global Land Cover Characterization . . . 35
3.4 Global Lakes and Wetlands Database . . . 37
4 Resampling 39 4.1 Purpose . . . 39
4.2 Choice of coordinate system . . . 39
4.3 Data structure and implementation . . . 39
5 Method 41 5.1 Data setup . . . 41
5.2 Correlation . . . 43
5.3 Radiometric normalisation . . . 44
5.4 Line fitting and residuals . . . 49
5.5 Land cover analysis . . . 50
6 Results 51
6.1 Correlation . . . 51 6.2 Influence of land cover . . . 60
7 Discussion 65
8 Conclusions 69
References 71
AMI Active Microwave Instrument
ASAR Advanced Synthetic Aperture Radar ENVISAT Environmental Satellite
ERS European Remote Sensing Satellite ESA European Space Agency
GM Global Monitoring Mode
GLCC Global Land Cover Characterization GLWD Global Lakes and Wetlands Database Landsat Land Remote-sensing Satellite
NASA National Aeronautics and Space Administration SADC Southern African Development Community SAR Synthetic Aperture Radar
SWI Soil Water Index
TM Thematic Mapper
USGS United States Geological Survey
1.1 Soil moisture information
Despite the fact that the water bound in the soil represents only a fraction of the total water reservoir, it is an imperative factor in a range of disciplines. In hydrology, soil moisture affects rainfall runoff, water accumulation and floods (Aubert et al. 2003). In agriculture it is an essential factor in the success of crop growth and prediction of crop yield. The influence on the global energy cycle is prominent through the energy transport between landmass and atmosphere through the processes of evaporation and thermal heat exchange. Thus, soil moisture information is an important ingredient in climate and weather prediction models.
Historically, soil moisture has been measured on a point basis with various in-situ methods. If data over large areas are required, these methods are expen- sive and time consuming. Furthermore, if regional estimates is the objective, the intricate problem of upscaling the point measurements arises. With remote sens- ing techniques, vast areas can be covered and continuously monitored, inherently supplying regional estimates.
At present there are no well established globally applicable methods for retriev- ing absolute volumetric soil moisture information by means of remote sensing.
Wagner (1998) and Scipal (2002) have developed a method to retrieve top soil layer and profile volumetric moisture content from ERS (European Remote Sens- ing Satellite) scatterometer data with an r.m.s. error of 4.9% and a resolution of 50 km. The method is foremost used to provide soil moisture change detec- tion and is not restricted to any particular geographical region. While providing
valuable information on the regional scale, there is a need to improve the spatial resolution for hydrological applications, where the main focus for soil moisture modelling in general is in the range of 1 – 10 km (Bl¨oschl 1996).
1.2 Purpose
The purpose of this study is to produce and evaluate scaling information that can be used for spatial downscaling of scatterometer derived soil moisture estimates to achieve a resolution of 1 km. The scaling information is based on ENVISAT Advanced Synthetic Aperture Radar data.
This chapter is an introduction to the physical mechanisms and technologies governing the acquisition of the radar data used in the analysis of this study.
Furthermore, concepts of importance for the work are defined and described.
The chapter begins with an introduction to microwaves, which is the type of electromagnetic waves used in radar systems, followed by a description of basic radar technology. The ENVISAT satellite and the ASAR sensor is presented.
Finally, the concepts of scale and scaling are defined and a description of soil moisture properties and measuring techniques is given.
2.1 Microwaves
2.1.1 Definition
The microwave part of the electromagnetic spectrum stretches between ∼0.3 GHz and 300 GHz, and is divided into bands. The data analysed in this report are acquired at 5.333 GHz (approximately 5.7 cm wavelength), situated in the C- band (∼4 – 6 GHz) (Ulaby et al. 1981). Microwaves have the advantage over visible light and other parts of the electromagnetic spectrum that the interaction with the atmosphere is weak. In the lower end of the microwave frequency range, up to about 20 GHz, signals can pass without significant distortion or attenuation through layers of clouds and rain, making the C-band frequencies beneficial for remote sensing (Henderson & Lewis 1998).
The concept of coherent electromagnetic waves is essential in the description of synthetic aperture radars. A beam of light is coherent if its waves exhibit approximately the same wavelengths and if the phase difference between locations on the beam is constant in time. A good example of nearly perfectly coherent
light is that produced by a laser. It is nearly monochromatic, i.e., consisting of a narrow wavelength interval, and in phase. The light from the sun is a mixture of a continuum of frequencies and waves with different phases and is therefore highly incoherent (Ulaby et al. 1981).
2.1.2 Surface interaction
When an electromagnetic wave arrives at the Earth’s surface, it undergoes trans- mission, reflection and absorption. The characteristics of light-matter interactions are determined by the geometric and dielectric properties of the matter and of the frequency and polarisation of the wave (Ulaby et al. 1982). The essential material properties are described by the electric permittivity ǫ and the magnetic permeability. The magnetic interaction is weak in comparison with the electric interaction and is usually neglected in remote sensing applications. The permit- tivity of a material can be written as
ǫ = ǫrǫ0, (2.1)
where
ǫr = ǫ′r− iǫ′′r (2.2)
is the relative permittivity and ǫ0 the permittivity of vacuum. The real part of ǫr determines the strength of the electric interaction with the material, while the imaginary part determines the attenuation of the wave as it penetrates the medium. A decrease in ǫ′r leads to a decrease in the reflected wave amplitude, while a decrease in ǫ′′r results in a greater penetration depth into the material. It should be noted that ǫr is frequency dependent (Rees 2001).
Assume that a wave impinges on a surface with an angle θi from the surface normal. This angle is termed the incident angle. Upon an ideally smooth surface,
a b
i r
Figure 2.1: a) Specular reflection of an electromagnetic wave on a smooth surface. b) Diffuse reflection on a rough surface.
the reflected part of the wave undergoes specular reflection and leaves the surface with an angle θr identical with θi but on the opposite side of the normal, as can be seen in figure 2.1a. If the surface height variations are in the order of, or smaller than, the wavelength of the wave, then the surface cannot be considered smooth and the wave will undergo diffuse reflection, as illustrated in figure 2.1b.
In this case, the reflected wave is spread depending on the nature of the surface roughness. A part of the wave will be transmitted into the ground. While being absorbed by the material, the wave will be scattered within the volume if incident on any discontinuities. The scattered waves within the volume will spread in all directions and might be emitted back to the air, leaving the volume. The strength of the volume scattering is determined by the degree of dielectric heterogeneity of the volume. The angular pattern of the volume scattering is determined by the surface roughness, the average relative permittivity of the volume and in the case of non-continuous scatterers, also the geometric size of the scatter elements (Ulaby et al. 1982).
In general, the reflected signal is a combination of surface and volume scat- tering, with contributions not only from the ground surface, but also from other elements, such as vegetation canopies and man-made structures. Furthermore, a wave could experience surface and volume scattering followed by secondary scattering on adjacent targets, e.g., reflected over an open field in the first step, followed by scattering on nearby trees in the second. The modelling of scattering is complex and under continuous research.
It is well documented, as in the work of Hallikainen et al. (1985), that the amount of water in the soil has a strong influence on the dielectric constant.
Therefore, the strength of the backscattered radiation is strongly influenced by the soil’s water content. This fact is exploited for soil moisture retrieval with passive and active microwave sensors.
2.2 Radar
2.2.1 Applications
Radar, or radio detection and ranging, is a technique that was invented to detect distance and speed of objects by sending out an electromagnetic pulse and detect- ing its echo. By measuring the time delay between emission of the pulse and the reception of the echo, the distance can be determined. The relative speed of the object can be established from the frequency shift of the return signal through the use of the Doppler equation (Lillesand et al. 2004).
Today, radars are used in a wide range of applications (Toomay & Hannen 2004), e.g., as weather radars, in missile detection systems, in air traffic con- trol and navigation and in airborne or space-borne systems, such as altimeters, scatterometers and synthetic aperture radars. Altimeters are used for accurate topographic mapping of land or ocean surfaces, with height resolutions of a few
centimetres. The scatterometer measures surface reflectivity and was originally intended for measuring surface winds over the oceans. It has been found ex- tremely useful in other areas, such as ice and land applications. The operation of synthetic aperture radar will be introduced in section 2.3.
2.2.2 Basic operation
As previously mentioned, the basic operation of a radar is to measure the time delay between transmission of an electromagnetic pulse and reception of its echo.
The pulse is generated by electric currents in an antenna. The return pulse is measured by the current induced from the electromagnetic echo and is subse- quently amplified and processed in preparation of information extraction. Most often, a correlation between the transmitted and received signals is needed for the derivation of desired parameters. Many radar systems use the same antenna for transmission and reception. In this case, the antenna needs to switch rapidly between transmission mode and receiver mode, which is done by means of a transmit-receive switch.
Signal modulation
The basic radar signal is a short pulse. The time delay between transmission and reception is given by (Toomay & Hannen 2004)
∆t = 2R
c , (2.3)
where c is the speed of light, and R is the distance between the antenna and the target, called the slant range. The range resolution is given by the smallest dis- tance for which two target echoes can be detected separately and is determined by the time duration of the transmitted pulse. The criterion for two targets to
τ < 2∆R
c , (2.4)
where τ is the duration of the pulse and ∆R is the distance between the targets.
The resolution in the range direction follows directly from equation (2.4) and is given by
rR = cτ 2 = c
2B, (2.5)
where B is the bandwidth, defined as B = 1/τ . Thus, to achieve a high spatial resolution, a short pulse length is desirable. If the pulse length is shortened, the average power of the pulse must be increased so that the received power exceeds that of the background noise by the required amount. For very short pulses, the average power of the pulse is only a small fraction of its peak power, requiring exceedingly powerful, and expensive, antennas. The solution to this problem is to use the so-called chirp technique, in which a longer pulse is “coded” by a sweep of the frequency of the pulse during transmission. In an up-chirped waveform, the part of the pulse first transmitted has a frequency f (t1). During the time of transmission of the pulse (∆t), the frequency is increased linearly so that the last part of the pulse has a frequency f (t2) = f (t1) + ∆f , as illustrated in figure 2.2.
When the echo is received, the transmission time delay of each constituent of the echo relative its initial frequency f (t1) is identified by its frequency. This infor- mation can then be used to compress the pulse with a shift of each constituent by the corresponding time interval, resulting in a short pulse with a high peak value and a high signal-to-noise ratio, without the need of a higher peak output.
The across-track resolution with this technique is given by
Figure 2.2: Chirp frequency modulation. The frequency is increased linearly during pulse transmission so that the received echo can be compressed, resulting in improved spatial range resolution.
rR,chirp = c
2B, (2.6)
where the bandwidth B is equivalent to ∆f in figure 2.2. The chirp technique is used in synthetic aperture radars to achieve high spatial range resolution (Henderson & Lewis 1998).
The radar equation
The energy transfer during signal transmission, scattering and reception of the radar signal is described by the radar equation (Curlander & McDonough 1991),
which can be written as
Pr = PtGtAr
(4π)2R2tR2rσ. (2.7)
Here, Pr is the total power received, Ptthe total power transmitted, Gtthe trans- mitting antenna gain and Ar the effective area of the receiving antenna. The gain is directionally dependent and is the ratio of the antennas radiative power to that of an isotropic antenna, i.e., an antenna that radiates equally in all directions.
The denominator is the product of two attenuation factors that arise from the fact that the signal spreads out spherically during transfer from the transmitting antenna to the target (4πR2t) and again from the target to the receiving antenna (4πR2r). In the monostatic case, when the transmitting antenna coincide with the receiving antenna, we can set
Rr = Rt = R, Gt= G, At= A, (2.8)
so that equation (2.7) simplifies to
Pr = PtGA
16π2R2σ. (2.9)
Here, σ is the radar cross-section. It is a collective term for the interaction of the signal with the target, and is given by
σ = As(1 − fa) Gs, (2.10)
where As is the effective area, Gs the gain and fa the fraction of absorbed radia- tion of the target (Toomay & Hannen 2004). The radar cross-section is a measure
of the ability to reflect the incoming radiation. It varies with the targets geome- try and dielectric properties as well as the direction of observation (through the gain, just as in the case of the antenna). The backscattering coefficient (σ0) is defined as the radar cross-section per unit differential area,
σ0 = dσ
dA, (2.11)
and is a natural parameter, i.e., it is independent of the system that is used to measure it (Curlander & McDonough 1991). It has the unit of m2m−2 but is most commonly given in decibels. In radar images, the backscatter coefficient is seen as levels of brightness intensity. Throughout this text, σ0 is referred to as backscatter intensity.
2.3 Synthetic Aperture Radar
As a satellite passes targets on the ground, each scatterer is recorded by many consecutive radar pulses. The synthetic aperture radar (SAR) uses the informa- tion contained in a series of radar echoes to generate a single image with a highly improved spatial resolution in the along-flight direction.
This section starts with a description of how the information from many con- secutive radar echoes can synthesise a high resolution radar image in a process called SAR focusing. The phenomena of speckle, which arise as a natural re- sponse to the coherent nature of SAR measurements, is then introduced. Lastly, the process of referencing the radar image to a geographical coordinate system is explained.
2.3.1 SAR focusing
A SAR transmits mutually coherent pulses, in contrast to a real aperture radar (ordinary radars), where the phase of each transmitted pulse is random. As a result, both range and phase of the backscattered echo are implicitly carried in its time delay. Between pulses, the satellite moves a short distance, which is seen as a phase shift in the return signal. Exploiting the phase information, the contributions from many consecutive pulses can be focused to generate a single image with a high resolution in the azimuth direction. The resulting “synthetic”
aperture is equivalent to a very long (typically several kilometres) real aperture radar. The process of combining the information from many consecutive radar echoes into a single response with highly improved spatial resolution is called SAR focusing (Henderson & Lewis 1998).
The backscatter of every transmitted pulse is recorded and stored in a memory segment in the range-azimuth memory space, as illustrated in figure 2.3. The range-azimuth memory space can be seen simply as a two-dimensional memory array, with each row corresponding to a different azimuthal satellite position and the cells of each row containing the detected backscatter at a certain range dis- tance. Since the antenna beamwidth (β) is wider than zero degrees, echoes from scatterers not orthogonal to the satellite line-of-flight give contributions to each range gate. Therefore, the wave fronts of the echoes from any particular scatterer are recorded in many consecutive receptions, while the scatterer is within the an- tenna beamwidth during the passage of the satellite (Henderson & Lewis 1998).
The variation of range to any particular point scatterer with azimuth is described by a hyperbola in the range-azimuth memory space. The position of the scatterer is given by the focus point of the hyperbola. Therefore, the position of a scatterer can be estimated through analysis of the curvature of the corresponding signal
Figure 2.3: Relation between observation space and memory space in a SAR system. The recorded signal from three scatterers at a particular satellite position is represented by three thick lines in the memory space.
locus in the range-azimuth memory space. The hyperbolic shape of the recorded signal changes with range as illustrated by figure 2.4. Before the actual focusing can be carried out, the conic shape of the recorded echoes in the memory space must be corrected. This is done separately for each position in the memory space, by re-mapping the contributions spread out over the range, as seen in figure 2.4, into a single range corresponding to the smallest antenna-scatterer distance, i.e., at the point of boresight, where the hyperbola has its minimum (Henderson &
Lewis 1998). The actual azimuthal compression is normally carried out in the fre- quency domain. The ideal backscatter response of a point target can be written as
sa(t) = E0eiφ(t), (2.12)
where E0 is the complex amplitude of the backscatter and φ (t) is the recorded
Figure 2.4: Shape of signal loci as recorded in the range-azimuth memory space. y1, y2and y3
represent the azimuth boresight position of the corresponding scatterers, as seen in figure 2.3.
The focus point of a locus give the position of the scatterer.
phase shift as a function of time (or azimuth, equivalently). The Fourier trans- form of sa(t) is then calculated and is here denoted Sa(ω). First, all phases of Sa must be aligned so that the subsequent summation results in a construc- tive interference. The phase alignment is given by a correlation of Sa(ω) with a matched filter (Curlander & McDonough 1991). The matched filter, m (t), is constructed in such a way that every point has exactly the opposite phase of the ideal impulse response, so that its Fourier transform (denoted by F {}) is given by
M (ω) = Fw (t) e−iφ(t) , (2.13)
where w (t) is the function that limits the compression of the (azimuthal) time history to the period of target illumination. In the simplest case, w (t) can be a box-like function given by
w (t) =
1 ,−tmax2 < t < tmax2 0 , otherwise.
(2.14)
Commonly it is a Hamming weighting function. The correlation, i.e., alignment and summation, is given by a product of the response with the matched filter:
V (ω) = Sa(ω) M (ω) . (2.15)
In the last step, the inverse Fourier transform of V (ω), given by
v (t) = F−1{V (ω)} , (2.16)
gives the the focused time-domain response v (t), which is centred in azimuth at the boresight pointing angle (corresponding to y1, y2 and y3 in figure 2.3).
For a point target the response has the form of a sinc-function (Curlander &
McDonough 1991). It should be noted that the compression results in improved signal-to-noise ratio, since noise, which is non-coherent in nature, is defocused during the compression and thus does not get amplified.
Azimuthal spatial resolution
The angle between the two directions where the antenna power is half that of the maximum is called half-power beamwidth. The angle, designated β, can be approximated by the expression
β = λ
l, (2.17)
where λ is the wavelength and l is the antenna length in azimuth direction (Henderson & Lewis 1998). Two targets on the ground can be separated if they are not within the antenna beam at the same time. In the case of a real aperture radar, the spatial resolution in azimuth is therefore given by
ra,real = Rλ
l, (2.18)
where R is the range distance. In the case of the synthetic aperture radar, how- ever, the effective beamwidth is produced by the synthesised antenna. The length (L) of the synthesised antenna equals the length of the full beamwidth (which is approximately twice that of the half-power beamwidth) footprint from the real antenna, so that
L = 2Rλ
l, (2.19)
The synthetic half-power beamwidth is then, in analogue with equation 2.17, given by
βs= λ L = l
2R, (2.20)
The resulting azimuthal resolution for the SAR is then
ra,synthetic = Rβs = l
2, (2.21)
so that the theoretical azimuthal spatial resolution of a SAR system equals half its antenna length. The independence of the distance between the satellite and the target can be explained by the increase in antenna beamwidth with range distance, which means that each target will be illuminated by a greater number of pulses, counterbalancing the increase in distance. Furthermore, the azimuthal spatial resolution is improved by a decrease in antenna size, since the resulting synthetic beamwidth is reduced. In the case of a real aperture radar, a smaller antenna, and thus a widened beamwidth, would result in a decreased spatial resolution.
2.3.2 Speckle
A result of the coherent nature of a SAR system is a wave interference phe- nomenon called speckle (Henderson & Lewis 1998), which can be seen as noise-like intensity variations in the radar image. It is, however, not noise. The intensity of the backscatter in a pixel in a radar image is formed as a sum of the contri- butions of many scatterers within the underlying area of that pixel. Since the transmitted radar pulse is coherent and also subsequently coherently measured (i.e., phase information is not destroyed), the resulting wave as measured by the antenna can be written as the complex sum of each sub-pixel contribution as
E0eiφ=
Ns
X
k=1
E0keiφk, (2.22)
where E is amplitude, φ phase and Ns the number of scatterers within the pixel.
The phase can be considered to be uniformly distributed in the interval [−π, π].
The summation in equation (2.22) can be visualised in the complex plane, as shown in figure 2.5. The resulting wave (left-hand side of equation (2.22)) is highlighted in the figure. Note that if any of the amplitudes E0k is much larger
Figure 2.5: Summation of complex wave contributions from individual scatterer elements.
Figure 2.6: Summation of complex wave contributions with one dominant scatterer element.
than those of the other scatterer elements, that scatterer will dominate the re- sulting measurement, as illustrated in figure 2.6. The existence of a dominating scatterer within a pixel can thus result in a measurement that basically repre- sents only that scatterer, instead of the collective properties of all underlying
contributions. While speckle contain information, it is not of interest for the purposes of this study. The influence of speckle can be reduced by averaging over measurements, at the cost of spatial resolution.
2.3.3 SAR geocoding
The SAR data in this report were geocoded by means of the Range Doppler ap- proach. It is a reverse method, which starts with geographically defined points on the ground from which the coordinates in the slant range and azimuth space are calculated. Along a mathematical model of the satellite orbit, points are picked out iteratively. It is assumed that the ground point has the known geographical coordinate ~P = {Px, Py, Pz} and is moving with the rotation of the Earth as defined by ~vp. The Pz coordinate is given by a digital elevation model. The satel- lite position and velocity are described by, respectively, ~S and ~vs. The Doppler frequency shift, fd, is given by the Doppler equation:
2 λ
~P − ~S
( ~vp − ~vs)
P − ~~ S
+ fd = 0. (2.23)
The Doppler frequency shift inferred on the coherently measured signal is given by
fd= 2v sin Θs
λ , (2.24)
where Θs is the angle between the satellite flight direction and the sensor look direction. Equations (2.23) and (2.24) is used to calculate the azimuth index of the ground point relative the nadir line of the satellite. Finally, the slant range index relative the nadir line is calculated with the range equation
P − ~~ S
R~
(2.25)
Now the azimuth and slant range index in the SAR image can be referenced to the known geographical coordinate of the ground point. After all pixels in the image have been geocoded, the image is geometrically rectified to the map projection, and the values of the projected pixels are interpolated bilinearly.
2.3.4 The ScanSAR technique
The swath width (swath is the area on the ground covered in a single satellite pass) of a SAR system is limited by the vertical beamwidth of the antenna.
ScanSAR is a technique in which rapid electronic beamsteering is used so that several sub-swaths can be scanned in sequence, allowing an expansion of the total swath width. Such a system transmits pulses to, and receives echoes from, a sub- swath for a period long enough to synthesise a radar image of the area within the beam footprint at the required resolution. It then switches beams to illuminate the next sub-swath and continues in this manner until the full-wide swath is covered, at which point it returns to the first sub-swath and the scanning cycle is repeated. An illustration of the technique, taken for ENVISAT’s SAR antenna (which will be described in detail in subsection 2.4.1), is given in figure 2.7. The ScanSAR technique is vital for monitoring dynamic processes, since the extended swath reduces the revisit time and also increases the frequency of global coverage.
Figure 2.7: Illustration of the ScanSAR scanning geometry. Credit: European Space Agency
— ESA.
2.4 ENVISAT’s Advanced Synthetic Aperture Radar
2.4.1 Platform
ENVISAT1 is an Earth observation satellite launched in March 2002 by the Eu- ropean Space Agency (ESA). The ENVISAT mission was designed to measure the Earth’s atmosphere and surface, with an emphasis on instruments that com- plement each other to support synergies between their measurements. It is the largest European satellite to date, carrying three spectrometers, two radiome- ters, two ranging instruments, a high resolution interferometer and two radars.
Figure 2.8: ENVISAT during preparation at the European Space Research and Technology Centre (ESTEC) in Noordwijk, The Netherlands. The black, rectangular ASAR antenna can be seen deployed. Credit: European Space Agency — ESA, A. Van Der Geest.
A photo of ENVISAT, with the large, rectangular SAR antenna fully deployed, can be seen in figure 2.8. Physical characteristics and system properties of the satellite are presented in tables 2.1 and 2.2.
Physical characteristics
Total mass 8211 kg
Launch dimensions length 10.5 m, envelope diam. 4.57 m Orbit dimensions 26 m x 10 m x 5 m
Table 2.1: Physical characteristics of ENVISAT.
The orbit is sun-synchronous and polar, which means that it is nearly circular,
System properties
Storage capacity 2x70 Gb + 30 Gb backup Downlink 50 Mbps or 100 Mbps Power supply (sunlight/eclipse) 3847 W / 3291 W
Table 2.2: Properties of ENVISAT’s on-board satellite system.
passing approximately over the poles and that it passes the equator at constant local times (10.00 AM for passes from the northern to the southern hemisphere).
Further orbital data are given in table 2.3.
Orbital parameters
Mean altitude 800 km Altitude range 780 – 820 km
Inclination 98.55◦ Orbit period 101 minutes Repeat cycle 35 days
Table 2.3: Orbital parameters of ENVISAT.
The repeat cycle, as seen in table 2.3, is 35 days. For sensors with wide swaths, this provides a complete coverage of the globe within seven days.
2.4.2 Sensor
The Advanced Synthetic Aperture Radar (ASAR) on-board ENVISAT is a so- phisticated instrument capable of monitoring Earth’s surface in several different modes, polarisations and incidence angles. Operating in C-band at 5.333 GHz, it ensures continuation of the data acquired with the ERS-1/2 scatterometers.
ERS-1 and ERS-2 were launched in, respectively, 1991 and 1995. ERS-1 was retired in 2000, while ERS-2 is still active. ASAR offers five mutually exclusive modes of operation; Image Mode, Alternating Polarisation Mode, Wide Swath Mode, Wave Mode and Global Monitoring Mode (GM). The first three modes
offer resolutions between 30 m and 150 m at different polarisations and incidence angles, with a maximum operation time of 30 minutes per orbit. The Wave Mode is used to retrieve information about ocean waves and winds. The data used in this study are acquired in Global Monitoring Mode, which has a resolution of 1 km and an orbit duty cycle of potentially a 100%. Further specifications are given in table 2.4.
The ASAR antenna consists of 320 transmit/receive modules which can have their phase and gain set individually, so that the beam can be steered and shaped.
This allows illumination at different look angles and therefore a widening of the swath (the ScanSAR technique, see subsection 2.3.4). High resolution in azimuth is achieved by the synthetic aperture radar technique. A high resolution in range is achieved by use of an up-chirp pulse modulation (subsection 2.2.2).
Global Monitoring Mode specifications
Spatial resolution 1000 m Swath 405 km Look angle range 15◦ – 45◦ Radiometric accuracy 1.39 dB
Polarisations HH and VV
Average revisit time 7 days at equator, 2 days at 70◦ latitude Table 2.4: Specifications for ASAR Global Monitoring Mode.
The use of the ScanSAR technique makes possible the 405 km wide swath and a high revisit frequency. The use of this technique also results in measurements being acquired at varying incidence angles, which have implications for the data processing and interpretation (Harris 1998).
2.5 Scales and scaling
2.5.1 Definition of scale
In this report, scale refers either to a characteristic length or time on which a certain process occurs. The purpose of using the term scale is not to describe an exact quantity. It is rather an indication of the typical dimension. Often, scale is referred to in terms of order of magnitude. As an example, a low-pressure system could have a spatial extent of the order of 1000 km and a lifetime of a few days.
These two quantities describe the spatial and temporal scale of the low.
There are, in particular, two scales that will occur throughout this text, namely the local and the regional ones. In this study, the local scale is the resolution of the ASAR GM data, i.e., 1 km, while the regional scale corresponds to the spatial resolution of the scatterometer data.
2.5.2 Process, measurement and model scale
When monitoring natural parameters that vary in space and time, it is essential to establish the spatial and temporal resolution at which measurements should be taken. This is called the measurement scale. The optimal situation is to have a measurement scale that is the same as the scale of the process itself, but this is seldom the case due to technological, logistical and financial limitations.
Furthermore, the scale at which the parameter is modelled must comply with the process and measurement scale. Recall the example of the low-pressure system.
The process scale is of the order of 1000 km and 1 day, in the space and time domains. Assume that there exist hourly measurements of ground air pressure from a network of meteorological stations. The measurements are taken at a point scale. Assume also that the objective is to model the pressure at a 10 km resolution. The temporal resolution in that case is suitable for the modelling of
the process, but there are differences in the spatial scales that need to be assessed.
As indicated in this example, total agreement between the process, measurement and model scales is rarely the case. The task of scaling information to overcome spatial or temporal gaps often becomes a necessity.
2.5.3 Definition of scaling
Scaling is the task of transferring information from one scale to another. Scaling can be complicated due to the spatial and temporal variability of the parameter at hand and by the fact that different processes can dominate the measurements at different scales. Upscaling is the process of transferring information from a smaller scale to a larger scale. An example is the process of performing some kind of averaging of the point measurements to the 10 km modelling scale, as required in the example from subsection 2.5.2. This report evolves around the concept of downscaling, in which information is to be transferred from a large scale to a smaller scale. There is a range of possible downscaling methods for soil moisture. An example is the method of Crow & Wood (2002), which estimates the sub-resolution variation of soil moisture by assuming that the statistical mo- ments of the soil moisture distribution follow a power-law relationship with scale.
This is mathematically expressed as
h ξλqi = λ λ0
K(q)
ξλq0 , (2.26)
where ξ is soil moisture, q the order of the statistical moment, K the scaling pa- rameter, λ and λ0 the target and nominal scale, and brackets indicating spatially averaged parameters. Evaluated over the Southern Great Plain of the United States this downscaling approach provided a result that proved simplified and sometimes inaccurate. Despite the imperfect results, when incorporated into a
land-atmosphere transfer scheme (TOPLATS) the downscaled information re- duced the error due to unexplained soil moisture variation by roughly 50%.
2.6 Soil moisture
2.6.1 Process
Soil moisture is the fraction of water and water vapour contained in the pores of the soil. The term is commonly used for the amount of soil moisture above the water table, i.e., above the level of abundant groundwater. Surface soil moisture in this report refers to the soil moisture of the top 0.5 – 2.0 cm of the soil.
In the space domain, the soil moisture process transcends many scales, from soil matrix flows and preferential flows (1 cm – 100 m) to large scale variation on the catchment or synoptic scale (10 km – 1000 km) (Bl¨oschl 1996). The process also exhibits high temporal variability, ranging from minutes to days or more, depending on spatial scale and soil depth.
2.6.2 In-situ measurement methods
Soil moisture has historically been measured at a point basis. The thermogravi- metric method consists of measuring the weight of a soil sample, which subse- quently is dried in a stove at 105 ◦C for 24 hours. Measuring the weight of the dry soil gives the mass of the evaporated water. With the time-domain reflectory (TDR) method, the moisture content is estimated from a measurement of the velocity of electromagnetic pulses through the soil. The pulse velocity is given by the dielectric properties of the soil, which in turn is determined to a large extent by the water content (Walker et al. 2004).
2.6.3 Observations from space
Soil moisture information is extracted by means of active and passive microwave remote sensing as well as from models incorporating different types of data, e.g., from thermal infrared sensors Wagner et al. (2006). As mentioned in the in- troduction, soil moisture information has been derived from scatterometer data acquired by ERS-1 and ERS-2 Active Microwave Instrument (AMI). From long time series of the data, backscatter intensities corresponding to completely wet (high backscatter) and completely dry (low backscatter) are determined. A Sur- face Water Index (SWI), is then derived from a scaling of the strength of the backscatter between the wet and dry reference. Seasonal variations in vegetation cover is taken into account with the help of the multi-incidence angle capabilities of the AMI. The SWI takes on values between zero and unity, where zero means completely dry and unity means moisture saturation. It can be interpreted as the fraction of water saturation in the pores of the upper soil layer. Soil moisture profile information is modelled from time series of SWI . The launch of ESA’s Meteorological Operational Satellite (METOP) with its Advanced Scatterometer (ASCAT) instrument is planned for 2006 and will continue the ERS-1/2 scattero- meter dataset. It will provide observations with a resolution of 25 km.
The potential of passive radiometer measurements has been well evaluated, e.g., by Jackson et al. (1999), Laymon et al. (1999) and Archer et al. (2003). NASA’s Aqua satellite with its Advanced Microwave Scanning Radiometer (AMSR-E) was launched in 2002 and has as one of its tasks to detect soil moisture. ESA’s Soil Moisture and Ocean Salinity (SMOS) mission is planned for launch in 2007.
SMOS will deliver global coverage every third day with a resolution of 35 km (Kerr et al. 2001). Attempts to model soil moisture are also conducted, e.g., by incorporation of precipitation data and thermal infrared imagery. Radiometers
have proved to give good estimations of top-soil moisture content. The spatial resolution of radiometer technology until present day is of the order of a few tens of kilometres. This is a result of the fact that radiometers rely on the weak natural emission of microwaves from the Earth’s surface. Furthermore, the synthetic processing for radiometers does not use phase information to synthesise high resolution.
3.1 Study area
3.1.1 Extent
The region chosen as the study area for this study is the Southern African Development Community (SADC) encompassing southern and central parts of the African continent, including South Africa, Lesotho, Swaziland, Namibia, Botswana, Zimbabwe, Mozambique, Zambia, Angola, Malawi, Tanzania, and the Democratic Republic of the Congo. Seychelles, Comoros, Madagascar and Mau- ritius, also part of the SADC, were not included due to a lack of data. Also the Republic of the Congo, Gabon, Equatorial Guinea, Uganda, Kenya and the southern parts of Cameroon, Central African Republic, Sudan and Ethiopia up to 6◦N were included. In total, the region covers an area of approximately 11 mil- lion km2.
3.1.2 Land cover distribution
The spatial distribution of land cover in the study area is presented in figure 3.1.
Forest, agricultural land and range land dominate the area. A small part of the region is covered by water bodies and wetlands. Since the radar backscattering over water is not used for soil moisture content retrieval, these two classes were not taken into account in the analysis. The distribution of the land cover categories used in the analysis is shown in figure 3.2. A description of the categories will be given in section 3.3.
Figure 3.1: Spatial distribution of land cover over the study area.
3.2 ASAR Global Monitoring Mode
3.2.1 General description
The foundation of this study is synthetic aperture radar data from the ASAR antenna on-board the ENVISAT satellite (see subsections 2.4.1 and 2.4.2). The data were acquired in Global Monitoring Mode (GM) with HH polarisation. The
Figure 3.2: Fractional distribution of land cover over the study area. Only categories used in the analysis are included.
geocoded data were stored as image pairs of the radar crossection (σ0) and the local incidence angle. The images had a swath of about 400 km and a north- south extent between 1000 km and 8000 km. Header information containing orbit and image parameters was available. As stated in table 2.4 in subsection 2.4.2, the radiometric and spatial resolution of the data were, respectively, 1.39 dB and 1 km. The pixel size of the images was 15 arcseconds corresponding to approximately 500 m at the equator. The ASAR GM data used in this report were acquired between December 2004 and May 2006, resulting in a data density of between 30 and 120 measurements for each point of the study area (typically
Figure 3.3: Left: detail of radar image with feature 1 (an industrial facility site in Saudi Arabia) encircled. Right: feature 2; branching of the river Kwa in Democratic Republic of the Congo.
biased toward the higher value). The images were not sampled to a fixed grid, a fact that required pre-processing before commencing the analysis. This issue will be addressed in chapter 4.
3.2.2 Geocoding quality
For the pre-processing and subsequently the analysis of the radar data to be done on a pixel basis, the geocoding accuracy had to be better than the pixel size of 500 m. Therefore, a quality analysis was carried out. Two features, with per- sistent backscatter strength in many images, were chosen. The first feature is an industrial facility site in Saudi Arabia, surrounded by desert. It is located at 24◦42′2′′N 44◦58′26′′E and has an extent of approximately 1 km times 0.5 km.
This feature showed persistently higher backscatter than its surroundings in many images, as exemplified in the left part of figure 3.3. The second feature is the branching of the river Kwa in Democratic Republic of the Congo, located at 3◦3′0′′S 16◦51′29′′E (see right part of figure 3.3). The river showed persistently lower backscatter than its surroundings. The geographical coordinates for fea- ture 1 and 2 were manually recorded in, respectively, 28 and 30 images. Despite the fact that the measuring accuracy of the analysis was estimated to ±1 pixel,
a comparison between the real coordinates of the features and averages of the recorded latitudes and longitudes showed that the precision of the geocoding was sub-pixel. Therefore, it was concluded that the processing and analysis could be done at the 500 m scale.
3.3 USGS Global Land Cover Characterization
For the analysis of land cover influence on the correlation between the local and regional backscatter, the U.S. Geological Survey’s Global Land Cover Charac- terization1 (USGS GLCC), namely the Africa Land Cover Characteristics Data Base version 2.0 with Lambert azimuthal equal area projection was used. The USGS GLCC is based on data from the NOAA satellite’s Advanced Very High Resolution Radiometer (AVHRR) data spanning April 1992 through March 1993 and have a resolution of 1 km. The classification is based on the Anderson scheme (Anderson et al. 1976). For the purpose of this report, the data were reclassified to the nine basic categories of the Anderson scheme. In the study area, five of the basic categories of interest were represented, namely urban or built-up land, agricultural land, rangeland, forest land and barren land. These are described below.
Urban or Built-up Land
Urban or built-up land is comprised of areas of intensive use with much of the land covered by structures. Included in this category are cities, towns, villages, strip developments along highways, transportation, power, and communications facilities, and areas such as those occupied by mills, shopping centres, indus- trial and commercial complexes, and institutions that may, in some instances, be isolated from urban areas.
1The data are available on the internet at http://edcsns17.cr.usgs.gov/glcc/.
Agricultural Land
Agricultural land may be defined broadly as land used primarily for production of food and fiber. This category contain, e.g., cropland and pasture, orchards, groves and vineyards.
Rangeland
Rangeland is defined as land where the potential natural vegetation is predom- inantly grasses, grasslike plants, forbs, or shrubs and where natural herbivory has been an important influence in its pre-civilization state. Some rangelands can have been seeded to introduced or domesticated plant species. This category includes tall grass (or true prairie), short grass, bunch grass and desert grass regions. Savanna is a subcategory of range land. Also included are former crop- lands or pasture lands (cleared from original forest land), which have grown up in brush in transition back to forest land to the extent that they are no longer identifiable as cropland or pasture from remote sensor imagery. Many of these brushlands are grazed by livestock and provide wildlife habitat, but not in an intensive manner.
Forest Land
Forest lands have a tree-crown areal density (crown closure percentage) of 10% or more, are stocked with trees capable of producing timber or other wood products, and exert an influence on the climate or water regime. Also included are lands from which trees have been removed to less than 10% crown closure but which have not been developed for other uses.
Barren Land
Barren land is land of limited ability to support life and in which less than one-third of the area has vegetation or other cover. In general, it is an area
of thin soil, sand, or rocks. Vegetation, if present, is more widely spaced and scrubby than that in the rangeland category. Unusual conditions, such as a heavy rainfall, occasionally result in growth of a short-lived, more luxuriant plant cover. Examples of barren land are deserts, dry salt flats, beaches and areas of bare exposed rock.
3.4 Global Lakes and Wetlands Database
Since scaling of soil moisture information over water bodies does not make any sense, large lakes, rivers and water reservoirs needed to be masked out before the land cover analysis was performed. The masking of these features was done with the Global Lakes and Wetlands Database2(GLWD) dataset (Lehner & D¨oll 2004).
This dataset was designed for purposes of hydrological modelling. Therefore, some reservoirs were represented by large circular disks, corresponding to the hydrological yield of the reservoir. For the purposes of the analysis in this work, reservoirs of this type were discarded. The classification of the data is presented in table 3.1. Only class 1, 2 and 3 were used in the masking. The resolution of the data is 30 seconds, corresponding to approximately 1 km at the equator.
2The data are available on the internet at
Cell value Lake or wetland type 1 Lake
2 Reservoir 3 River
4 Freshwater Marsh, Floodplain 5 Swamp Forest, Flooded Forest
6 Coastal Wetland (incl. Mangrove, Estuary, Delta, Lagoon) 7 Pan, Brackish/Saline Wetland
8 Bog, Fen, Mire (Peatland) 9 Intermittent Wetland/Lake 10 50 – 100% Wetland
11 25 – 50% Wetland
12 Wetland Complex (0 – 25% Wetland)
Table 3.1: Global Land and Wetlands Database land cover classification.
4.1 Purpose
The geocoded radar data were not sampled to a fixed grid, i.e., the geocoded latitude and longitude of the pixels in the images did not belong to multiples of a fixed geographical distance. For the essential operation of this report — the correlation of local radar backscatter intensity with regionally averaged back- scatter intensity — it was necessary to have the data supplied with a fixed grid, since data for specific positions in the study area were to be extracted from many radar images. Furthermore, the analysis required efficient extraction of all avail- able data for a specific geographical coordinate or specified region, a condition not fulfilled in the original database. Therefore, the radar data needed to be resampled to a fixed grid, with a data structure that facilitated, in an efficient manner, both spatial and temporal analysis.
4.2 Choice of coordinate system
The coordinate system for the resampling was chosen as a coarse grid with a spacing of 30 arcminutes in longitudinal and latitudinal directions with the origin set to 180◦W 90◦S. Within each 30 arcminute2 element (in this report called a
“gridbox”), a subgrid with a spacing of 15 arcseconds was defined. The pixel size correspond to about 500 m at the equator. This pixel size was chosen to satisfy the Nyquist sampling criterion of the radar data, which have a resolution of 1 km.
4.3 Data structure and implementation
The resampling procedures were implemented in the Interactive Data Language (IDL). For any region of interest, the software identifies all radar images overlap-
ping the region, resamples all information contained therein by means of bilinear interpolation and creates a new database. Each of the structures corresponds to the data from a specific radar image and contains the resampled backscatter intensities and local incidence angles. Also contained is a reference to the resam- pled image, the date and time of image acquisition and an indicator for the node (ascending or descending) of the satellite orbit at acquisition.
In order to investigate the degree of interrelationship between the backscatter signal on the regional and local scales an automated analysis was performed for each position in the 1 km resolution grid. At each location in the study area, a time-series of local backscatter values from all available ASAR GM data was extracted. The simultaneous regional backscatter values were generated as an average of the backscatter of the ASAR GM data over a rectangular window, centred around the position where the local backscatter was extracted. The rela- tionship between the two data sets was then analysed in a processing procedure.
The main result is the correlation coefficient, which describes the linear temporal consistency between the local and regional backscatter signal. To quantify the shift of the measurements from the one-to-one line (which has zero intercept and unit slope) the mean linear residual was calculated. The slope and intercept of the line equation fitted to the data sets were calculated together with the mean absolute error. To get a measure of any nonlinear interrelationship, in case of a weak linear correlation, the rank correlation coefficient was produced. Several quality indicators were produced, e.g., the number of measurements used and the fractional data coverage in the regional window, to give a measure of the quality of the resulting information. In the second step of the analysis, the influence of land cover on the correlation between the scales was analysed.
5.1 Data setup
The analysis was carried out on a 1 km and 25 km basis, for every 0.5 km pixel in the study area. The local backscatter time series for every position could be extracted directly from the resampled data, resulting in data-sets on the form
σloc0 =σ0loc(t1) , σloc0 (t2) , . . . , σloc0 (tN −1) , σloc0 (tN) , (5.1)
where ti is the acquisition time of the i:th image from where the data were ex- tracted. The corresponding regional mean time series was generated from the arithmetic mean of the local backscatter in a 25 km by 25 km rectangular win- dow in all images according to
σreg0 (ti) = dB
"
1 Np
Np
X
k=1
σk0(ti)
#
(5.2)
centred at the extraction point of local backscatter. Note that σk0 is inserted into equation (5.2) in linear values and that the resulting mean is converted to decibels. Np is the number of pixels (i.e., local backscatter measurements) in the regional window of each image. This results in the regional mean data-set
σreg0 =σreg0 (t1) , σ0reg(t2) , . . . , σreg0 (tN −1) , σ0reg(tN) . (5.3)
Since the coordinate system was defined with a pixel size in fixed degrees, in contrast to a fixed distance, the number of pixels in east–west direction of the 25 km wide regional window depends on the latitude. Assuming a spherical da- tum, the window stretches out in east–west direction as described by
new = n0· 1
cos ϕ, (5.4)
where n0 is the number of pixels at the equator and newis the number of pixels at latitude ϕ. Thus, the total number of pixels used to generate the regional value was
Np = n20
cos ϕ. (5.5)
To improve processing performance, the same regional mean was used for each group of nine adjacent pixels in a three by three pixel window. The regional mean was only calculated for the central pixel and was used as a close approximation for the eight surrounding pixels.
A data coverage filter calculated the spatial coverage of each regional window and removed regional means that were generated on a statistically weak data basis. A minimum coverage of 75% was postulated and used in the processing.
5.2 Correlation
The Pearson’s coefficient of correlation (R) was calculated according to
R =
1 N
PN
i=1(σloc0 (ti) − hσloc0 i) σ0reg(ti) − hσreg0 i q1
N
PN
i=1(σloc0 (ti) − hσloc0 i)2 q1
N
PN
i=1 σ0reg(ti) − hσreg0 i2, (5.6) with σ0loc and σreg0 denoting the sets of regional and local backscatter given by, respectively, equations (5.1) and (5.3). The correlation was carried out with σloc0 and σ0reg in decibels. The numerator in equation (5.6) is the covariance, which gives a measure of how much the variables vary in synchrony, i.e., the absolute value of the covariance increases if the two data sets vary from their respective means in a coupled manner. The denominator is the product of the standard deviations of the local and regional backscatter. The Pearson correlation coefficient gives a measure of the linear relationship between the two variables and can take on values between −1 and +1. A value of +1 indicates a perfectly linear relationship, i.e., the data fit perfectly to a line with positive slope, and a
close to zero indicates that no linear relationship exists between the parameters.
In a scatterplot, the measurements of such a data set are usually seen as a cloud of which the spread increases, as |R| approaches zero. To quantify the strength of the linear relationship, the coefficient of determination, here denoted R2, was formed as the square of the Pearson correlation coefficient. R2 can be directly interpreted as the fraction of the variation in one variable that is explained by the variation of the other variable (Isaaks & Srivastava 1989). As an example, a value of R2 equal to 0.6 means that 60% of the variation of one variable is explained by the variation of the other. An example of a scatterplot showing R and R2 is given in figure 5.1.
5.3 Radiometric normalisation
The backscatter intensity depends on the angle of local incidence of the radar beam at the Earth’s surface. This means that keeping all other variables (like surface roughness, vegetation cover, soil moisture, etc.) constant, a higher back- scatter is typically received from an area closer to the satellite than from an area in the far end of the swath. The correlation coefficient is sensitive to outliers, i.e., pairs of local and regional backscatter values that lie far away from the bulk of measurements in the scatterplot (see figure 5.1). Outliers can affect the result- ing correlation coefficient strongly, due to the relatively high covariance, and if the outlier is close to the 1:1 line, the correlation can be increased significantly (Isaaks & Srivastava 1989). Due to the nature of the satellite orbit, the swaths of separate images coincide. When calculating the correlation coefficient for a position, it is possible that the majority of measurements comes from a point on the far end of one “swath track”, while only few or no measurements come from a point on a swath that was close to the satellite at the time of acquisition. This can result in an uneven distribution of measurements over the incidence angle
Figure 5.1: Example of a scatterplot with regional backscatter intensity plotted against local backscatter intensity. The resulting correlation coefficient (R) and coefficient of determination (R2) indicate a strong linear relationship between the parameters.
range, with a few outliers with significantly higher backscatter. Furthermore, the incidence angle dependency expands the range of local and regional backscatter, so that the group of measurements in the scatterplot is spread out, altering the value of the correlation coefficient. An example of induced artifacts can be seen in figure 5.2. The solution was to perform a radiometric normalisation to a refer- ence angle of incidence prior to the correlation analysis. This was done for each position in the grid through a fit of the line equation to the set of local back- scatter and local incidence angles followed by a calculation of the corresponding
Figure 5.2: Artificially induced striping in the correlation coefficient can be seen as two down- ward pointing triangles from the upper part of the image toward the centre. Bright colours represent high correlation.
backscatter at a reference angle of 30◦. Linear models for incidence angle depen- dency on radar backscatter have been used, e.g., by Frison & Mougin (1996) for the ERS scatterometer, and by M¨akynen et al. (2002) for RADARSAT ScanSAR data. Furthermore, the σ0-model of ESA’s ASAR product handbook (ESA 2006) shows a closely linear relationship for incidence angles between 15◦ and 45◦, which is the angular range of the data used in this study (see table 2.4). An example of typical backscatter dependency on incidence angle is given in figure 5.3. The normalisation is mathematically written as
σ0(30) = σ0(θ) − b (θ − 30) , (5.7)
where σ0(30) is the normalised backscatter, σ0 the original backscatter, b the slope of the fitted line (shown in figure 5.3) and θ the incidence angle at which the measurement was acquired. After the normalisation, the linear dependency
