ASSESSMENT OF RADIOMETRIC CORRECTION TECHNIQUES IN ANALYZING VEGETATION VARIABILITY AND CHANGE USING TIME SERIES OF LANDSAT IMAGES
Sergio M. Vicente-Serrano1*, Fernando Pérez-Cabello2 and Teodoro Lasanta1
1 Instituto Pirenaico de Ecología, CSIC (Spanish Research Council), Campus de Aula Dei, P.O. Box 202, Zaragoza 50080, Spain
2 Departamento de Geografía. Universidad de Zaragoza. C/ Pedro Cerbuna 12. 50009. Zaragoza. Spain.
*
Abstract. The homogeneity of time series of satellite images is crucial when studying abrupt or gradual changes in vegetation cover via remote sensing data. Various sources of noise affect the information received by satellites, making it difficult to differentiate the surface signal from noise and complicates attempts to obtain homogeneous time series. We compare different procedures developed to create homogeneous time series of Landsat images, including sensor calibration, atmospheric and topographic correction, and radiometric normalization. Two seasonal time series of Landsat images were created for the middle EbroValley (NE Spain) covering the period 1984–2007. Different processing steps were tested and the best option selected according to quantitative statistics obtained from invariant areas, simultaneous medium-resolution images, and field measurements. The optimum procedure includes cross-calibration between Landsat sensors, atmospheric correction using complex radiative transfer models, a non-lambertian topographic correction, and a relative radiometric normalization using an automatic procedure. Finally, three case studies are presented to illustrate the role of the different radiometric correction procedures when analyzing and explaining gradual and abrupt temporal changes in vegetation cover, as well as temporal variability.We have shown that to analysedifferent vegetation processes with Landsat data, it is necessary to accurately ensure the homogeneity of the multitemporal datasets by means of complex radiometric correction procedures. Failure to follow such a procedure may mean that the analyzed processes are non-recognizable and that the obtained results are invalid.
Key-words. Landsat time series, TM-ETM+ cross-calibration, atmospheric correction, relative normalization, vegetation change, Ebro Valley, Spain
1. Introduction
The Landsat program for Earth observation has provided invaluable information on the Earth’s surface characteristics over the past four decades (Cohen and Goward, 2004). Landsat images have been widely used for land cover mapping and the creation of vegetation inventories at different spatial scales (e.g., Bossard et al., 2000). Moreover, the systematic archiving of Landsat data makes this information highly valuable for retrospective analyses of land surface characteristics. Together with other types of satellite images such as NOAA-AVHRR (National Oceanic and Atmospheric Administration-Advanced Very High Resolution Radiometer) and MODIS (Moderate Resolution Imaging Spectroradiometer), multitemporal Landsat data have been widely used in recent decades to identify changes in land cover (e.g., Lenney et al., 1996). The spatial and spectral resolution of Landsat images makes these data highly suitable in analyzing both abrupt and gradual changes in vegetation cover and monitoring environmental processes such as degradation and desertification (Almeida-Filho and Shimabukuro, 2002), deforestation (Cohen et al., 2002; Huang et al., 2007), habitat fragmentation (Millington et al., 2003), forest succession (Song and Woodcock, 2003; Song et al., 2007), overgrazing (Jano et al., 1998; Pickup and Chewing, 1994), rangeland monitoring (Hostert et al., 2003), and vegetation recovery after natural disturbance such as volcanism (Lawrence and Ripple, 1999) and forest fires (Viedma et al., 1997; Lozano et al., 2007).
Nevertheless, there exist constraints in using Landsat data for multitemporal studies because of problems in obtaining homogeneous time series, not affected by non-surface noise, and in which the images are comparable between different dates since the data only report on surface conditions. Multitemporal satellite-image datasets are affected by different sources of noise related to the stability of sensors, changes in satellite responsivity, changes in illumination, atmospheric effects, etc. These problems are not unique to Landsat images, but they are more difficult to overcome in Landsat images compared with images compiled by other satellites because the low temporal resolution of Landsat images makes it impossible to apply simple homogenization procedures such as composite creation (Holben, 1986) and temporal filtering (Van Dijk et al., 1987), which markedly reduce the non-surface noise.
Notable efforts have been made to reduce non-surface noise in Landsat images, including attempts to calibrate the sensor to correct lifetime radiometric trends (Teillet et al., 2004; de Vries et al., 2007), cross-calibrate the images obtained from different sensors (Teillet et al., 2006; Röder et al., 2005), minimize atmospheric noise (Chavez, 1989; Ouaidrari and Vermote, 1999; Liang et al., 2001), and reduce the influence of topography (Civco, 1989; Gu and Gillespie, 1998; Pons and Solé, 1994); however, other studies have shown that the application of accurate sensor calibrations and complex atmospheric corrections does not guarantee the multitemporal homogeneity of Landsat datasets (Schroeder et al., 2006) because complete atmospheric properties are difficult to quantify and simplifications are commonly assumed.This problem has led to the development of relative radiometric normalization techniques based on adjustments to the radiometric properties of an image time series to match that of a single reference image. In recent years, efforts have been made to develop methods to select invariant pixels for a reliable application of relative normalization techniques(Du et al., 2002; Chen et al., 2005; Canty et al., 2004).
Protocols have been proposed in processing multitemporal Landsat datasets (e.g., Hall et al., 1991; Hill and Strum, 1991; Han et al., 2007), comprising the following steps: i) geometric correction, ii) calibration of the satellite signal to obtain Top of the Atmosphere Radiance, iii) atmospheric correction to estimatesurface reflectance, iv) topographic correction, and v) relative radiometric normalization between images obtained on different dates. Radiometric processing is recommended to be done prior to geometric processing since this resampling step generally smoothness the data set. Nevertheless, Landsat users in Europe commonly do not have access to data that have not already been geometrically corrected (the Level 1 System Corrected -1G- from the European Space Agency (ESA) is considered the standard product for most users since previous levels require extensive processing). When these protocols are applied, several decisions must be made due to the different procedures available at each step.
The high degree of interest in using multitemporal Landsat datasets is in contrast with the small number of studies that have tested different procedures of radiometric correction with the aim of obtaining temporally homogeneous images. The majority of studies have tested the influence of atmospheric correction (Song and Woodcock, 2003; Norjamäki and Tokola, 2007) or relative normalization (Yuan and Elvidge, 1996; Tokola et al., 1998; Heo and Fitzhugh, 2000; Olthof et al., 2005) on the temporal stability of time series of Landsat images. In contrast, few studies have tested the influence of each of the above steps in obtaining temporally stable time series of images and few have assessed the relative importance of using simple or complex techniques at each step. Existing studies have shown that relative normalization is the most critical step in obtaining temporal stability in a series of images (Schroeder et al., 2006; Janzen et al., 2006), although the calibration procedure also has a noticeable effect on the final results (Paolini et al., 2006) and atmospheric correction is important in obtaining accurately estimate surface reflectance's and accurate magnitudes of vegetation indices over time (Song and Woodcock, 2003).
Few studies have analyzed the role of complete radiometric correction protocols in processing multitemporal Landsat data when a number of different vegetation processes are of interest. Nevertheless, this is a crucial topic when analysing vegetation multitemporal dynamic using Landsat data. Thus different results have been found for land classification (Song et al., 2001; Paolini et al., 2006; Norjamäki and Tokola, 2007) and forest succession (Song and Woodcock, 2003; Schroeder et al., 2006) as a function of the radiometric correction applied.
The present study involves a complete evaluation of various radiometric correction processes required to obtain accurate time series of Landsat imagery, including calibration, atmospheric correction, topographic correction, and relative normalization. The objective is to identify the influence of different processing steps on multitemporal spectral reflectance trajectories developed with Landsat data and also to detect their performance to analyse different vegetation processes.Three case studies related to different processes of vegetation change and variability are provided to illustrate how different results can be obtained as a function of the employed radiometric correction protocol and the ecological process of interest: land cover change, forest regeneration after fire and ecosystem response to climate variability. The case studies were carried out in a complex ecological region in which several natural and human-induced processes of changes in vegetation cover have been recorded in recent decades.
2. Study area
The study area is the central Ebro Valley, Spain (Landsat Path 199 Row 31), located in the northernmost semiarid region in Europe. Figure 1 shows the location of the study area, including a detailed land cover map. Surrounded by mountain chains, the EbroValley has a Mediterranean climate with continental characteristics, with marked spatial and seasonal variations in precipitation; the dry season occurs during the summer months. Theelevation ranges from less than 300 mabove sea level in the middle of the Ebro Depression to more than 2000 m a.s.l. in the Prepyrenees (North) and 2000 m a.s.l. in the highest peaks of the IberianRange (South). In the central part of the valley, mean annual precipitation is 326 mm, with marked seasonality. The study area shows a strongly negative water balance (precipitation minus potential evapotranspiration), greater than 900 mm. The dominant land covers include steppes and herbaceous cultivation in dry farming areas. Coniferous forests (mainly Pinus halepensis) cover the few slopes of the region. Vegetation distribution is strongly controlled by aridity (Vicente-Serrano et al., 2006), and drought conditions have a marked influence on vegetation cover and activity (Vicente-Serrano, 2007).The lithology of the area ischaracterised by limestones and gypsums (Peña et al., 2002) that contributeto its aridity, since the soils are unable to retain the water as a consequenceof the high hydraulic conductivity.
3. Methodology
3.1. Satellite imagery database
To determine the patterns of vegetation variability and change, it is important to take into account seasonal variations in the distributions of different vegetation types. In the middle EbroValley, herbaceous species, shrubs, and forests show contrasting seasonal variations in vegetation activity (Vicente-Serrano et al., 2006). The peak activity for herbaceous species and shrubs occurs in spring; for forests in summer. These seasonal fluctuations make it difficult to capture the full range of vegetation processes with just a single season of imagery.
Abrupt changes in vegetation activity commonly occur as a consequence of climate seasonality; consequently, it is important that the capture dates of images are similar in different years. It is not possible to combine images taken in different months, as this would have a strong influence on the temporal homogeneity of the dataset.
We reviewed all of the available Landsat-Thematic Mapper (TM) and -Enhanced Thematic Mapper (ETM+) images in the archives of the European Space Agency (ESA). Since most studies based on Landsat data consider ETM+ and TM radiometry to be comparable,We used TM data from 1984 and, then, switched to the ETM+ data when it became available in 1999 due to the better calibration of this sensor (Teillet et al., 2001); we subsequently switched back to using TM date due to the ETM+ SLC failure in 2003. Frequent cloud cover in spring means that only the month of March has a sufficient number of images to enable an analysis of vegetation activity. Months in summer pose fewer problems in terms of obtaining reliable clouds-free images; therefore, we selected the month of August to create a second time series. A total of 28 Landsat-TM and -ETM+ images taken between 1984 and 2007 were acquired from ESA, 16 corresponding to the summer season and 12 to spring. Table 1 lists the dates and types (TM or ETM+) of the selected images.
3.2. Geometric correction
The images were orthorectified using control points and a 1-m digital elevation model obtained from stereo pairs of aerial photographs, and resampled to 30 m to match the TM and ETM+ resolutions. The image taken on 8 August 2000 with good visibility and being free of clouds, was registered using orthorectified digital aerial photographs as a reference. The rest of the images were co-registered to this image using control points. The images were orthorectified following the method of Palà and Pons (1995), which includes elevation data in performing a polynomial geometric correction. The X and Y root mean square error wasless than 15 m (0.5 pixels) for all images, guaranteeing a precisegeometric match among them. After geometric correction, cloud-covered and cloud shadows were manually digitized and eliminated.
3.3. Calibration
A precise calibration is required to convert DN (Digital Numbers) to satellite radiances in W m–2 sr–1m–1. This is highly problematic for Landsat5 imagery because calibration coefficients are not constant over time. A large variation, on the order of 20%, has been observed between the prelaunch gain coefficients and postlaunch cross-calibration (Chander et al., 2004; Chander et al., 2007). Moreover, the different spectral responsesof TM and ETM+ sensors introduces comparability problems between Landsat5 and Landsat 7 images(Teillet et al., 2001).
For Landsat 5 imagery, the European Spatial Agency (ESA) has used constant calibration dynamic ranges since 1984 ( embedded within the product format. Chander et al. (2004) and Teillet et al. (2004) demonstrated an exponential decay of the solar reflective bands of Landsat 5-TM since 1984, with some differences between bands. We corrected the calibration coefficients embedded in the ESA products with reference to the time after launch of Landsat 5 in 1984, applying the equations formulated by Teillet et al. (2004). This procedure is recommended by the ESA in recalibrating Landsat products (Saunier and Rodríguez, 2006). The coefficients of calibration for Landsat 7-ETM+ were obtained according to upper and lower at-satellite radiance indicated in the Landsat-7 Science Data Users Handbook ( using the values corresponding to the date (before or after July 1, 2000) and the type of gain (high or low, as embedded within the product format).
Although the Landsat 5-TM and Landsat 7-ETM+ bands are commonly considered to be comparable, previous studies have demonstrated important differences between the two that may affect comparability, suggesting the need to apply cross-calibration procedures to both sensors because for the majority of bands (2, 3, 4, and 7) the TM sensor underestimates the radiance values regarding ETM+ (Teillet et al., 2001; Vogelmann et al., 2001).
Vogelmann et al. (2001) developed a simple procedure to cross-calibrate the Landsat-TM and ETM+ images, applicable to Level 1G formats. These authors proposed empirically derived slope and intercept values to convert Landsat 7 ETM+ DNs to Landsat 5 TM DNs based on two simultaneous images taken on June 2, 1999. We used these slope and intercept values in the present study to adjust the Landsat 7-ETM+ DNs to Landsat 5-TM DNs.
Satellite radiance was obtained for Landsat 5-TM quantized calibrated pixel values inDNs, Landsat 7-ETM+ quantized calibrated pixel values in DNs, and the cross-calibrated Landsat 7 ETM+ quantized calibrated pixel values in DNs to Landsat 5-TM quantized calibrated pixel values in DNs according to:
(1)
where Lsat is the satellite radiance in W m–2 sr–1 μm–1 for band . Grescale and Brescale are the calculated band-specific rescaling factors.Satellite radiances were converted to Top Of the Atmosphere (TOA) reflectances according to
(2)
where is the TOA reflectance for band , d is the earth–sun distance in astronomical units, ESUN is the mean solar exoatmospheric irradiance for band , and s is the solar zenith angle in degrees. ESUN values were obtained from Chander and Markham (2003) for TM images and from the Landsat-7 Science Data User Handbook for ETM+ images (
3.4. Atmospheric correction
Although TOA reflectance values are widely used in inventory and ecosystem studies, they do not take into account signal attenuation by the atmosphere, which strongly affects the intercomparability of satellite images taken on different dates. Upward and downward irradiance is modified by two atmospheric processes: absorption by gases and scattering by aerosols and water molecules (Vermote et al., 1997a; Lu et al., 2002). There are noticeable daily variations in atmospheric gas concentrations and aerosol volumes; these must be taken into account in multitemporal studies. Moreover, the effect of atmospheric processes is spectral-dependent (Teillet and Holben, 1993; Cachorro et al., 2000), affecting the magnitude of band ratio operations such as vegetation indices (Myneni and Asrar, 1994).
Different models have been developed to minimize the noise introduced by atmospheric processes on the signal received by the satellite, ranging from simple methods based on information contained in the image [e.g., Dark Object Subtraction (DOS)-based methods (Chavez, 1996)] to complex radiative transference modelssuch as SMAC (Rahman and Dedieu, 1994), 6S (Vermote et al., 1997a), MODTRAN (Berk et al., 1999), and ATCOR (Richter, 1996) that simulate the atmosphere/light interactions between the sun surface and surface-sensor trajectories.
In this paper we tested two methods developed to minimize the atmospheric effects on Landsat time series: i) one based on a modification of the DOS method that includes some atmospheric information (Dark Target Approach -DTA-), and ii) one based on a complex radiative transfer code that includes some atmospheric information commonly available in global climate databases and some parameters obtained from the images themselves.
3.4.1. DTA method
Dark Target Approach (DTA)-based methods assume that some areas of the images have near-zero reflectance. Although there are various DTA-based methods (Chavez, 1996), Song et al. (2001) reported better results when including in the method atmospheric transmittance and Rayleigh scattering. We therefore followed this approach in the present study. Surface reflectance can be obtained as follows (Song et al., 2001):