(adapted from "Atmospheric and radiometric correction of Landsat TM 5 Data using the COST model of Chavez, 1996", S. M. Skirvin)
Introduction
The image-based COST method for atmospheric correction has been combined with radiometric calibration in a model which streamlines this time-consuming but crucial processing step. Although some initial scene-dependent calculations must be performed to obtain necessary model constants, subsequent processing of raw DN values to atmospherically corrected reflectance values is performed in a single process. The model implements the Chavez (1996) improved dark-object atmospheric correction for Landsat TM5 multispectral data (bands 1-5 and 7). The input image file is assumed to have only these 6 bands.
Please read through the following steps carefully prior to beginning this procedure with your image, just to see what is coming up and to make sure that you have all of the necessary components in your working directory. Plan on at least two hours to complete this procedure and run the models (more depending on how busy the network is and/or what else your machine is doing). You will need a few GB of disk space for the input and output images involved if you are processing an entire Landsat scene.
Background
The inputs to the model are the Earth-Sun Distance, sun elevation angle, and minimum DN values for each band.
The model first converts each minimum DN value to an at-satellite minimum spectral radiance value:
For each band, the theoretical radiance of a dark object (assumed to have a reflectance of 1%
by Chavez 1996 and Moran et al. 1992) is computed:
spectral irradiance from Table 4 of Markham and Barker (1986), d is the sun-earth distance,
and theta is the solar zenith angle (90-sun elev).A haze correction is computed using the computed dark object values (Chavez 1996):
.
The fundamental radiance to reflectance (rho) equation (eq. 2 of Chavez 1996) is:
.
The output from the model is in reflectance units, which should range from 0 to 1. Occasionally, values greater than 1 are obtained, and these generally correspond to bright objects (e.g., clouds, snow, playa) or sometimes noise (randomly or systematically scattered saturated pixels). The output data format is float single. If you plan to compute indices or classify using reflectance data, you should keep the images in this float-single format. If you have some other use that does not require this radiometric resolution, you may want to rescale the images to 8-bit format to save disk space.
Running the Model
First, obtain the relevant values for your scene:
The values that you will need to modify in the graphical model are Lhaze, sun elevation angle, and Earth-Sun distance. The initial inputs for the model are entered into an Excel spreadsheet, in which the Lhaze value will be computed for you once the other values are entered. NOTE: If the computed Lhaze is negative, it is probably due to a minimum DN that is too small - in this case, enter '0' (zero) in the model equation for this parameter rather than the negative value.
COST_TM5-input.xls
(right-click to download file)
The COST model for TM5 images runs from an Imagine graphical model. Download this to your working directory and open it in Imagine Modeler> Model Maker.
COST_TM5.gmd
(right-click to download model)
The graphical model looks like the image shown below. Modify the top and bottom icons to represent the input and output images, respectively. Modify the equation in each of the icons (circles) in the second row using the values from the spreadsheet for bands 1,2,3,4,5,7 in that order, left to right. Do not modify anything else. This example equation is also shown in the Excel file. For each band, modify the colored/underlined portions only, using values from the corresponding columns in the spreadsheet.
MODEL = (( -2.8890805 + (0.0602353 * $n1_gb_tm_or(1) - 0.15)) * PI * 0.9932554 ** 2) / (195.7 * COS (PI/180 * (90 - 52.21)) ** 2)

Save the spreadsheet and the modified graphical model with your image. Then run the model by clicking on the red 'lightning bolt' icon in the modeler window. Wait - it could take up to an hour for an entire scene.
Zero the Background in the Output Image
This process fixes the phenomenon that band 7 is not completely zeroed in the COST model, resulting in a non-transparent background. You can check this by opening the image and using the inquire cursor to see what the values of pixels outside of the image area. If necessary, the atmospherically corrected image will need to be run through the following model:
zero_cost_output.gmd
(right-click to download model)
As usual, open the model in Modeler> Model Maker and manually
enter the input and output images by double-clicking on the accordion-shaped
icons.
OBTAINING VALUES FOR ATMOSPHERIC/RADIOMETRIC CORRECTION INPUT VARIABLES
Sun Elevation Angle
For discussion regarding the influence of the sun elevation angle term in the COST model, refer to "Notes on the Sun Elevation Angle and 1/cos(theta) term in the COST Model". This is of particular relevance if the sun elevation angle for your scene is below about 40 degrees.
Sources of the Sun Elevation Angle: Check to see if the image has a header available that came with your order or that is on ARIA. If there is no header available, you should be able to retrieve this information online using one of the following sources:
1) Search for the image by path/row and date on the EOS-DIS Data Gateway:
http://edcimswww.cr.usgs.gov/pub/imswelcome/
Navigate to the Data Granule Attributes page, then scroll down toward the bottom of the table and take the Sun_Elevation_Angle.
2) Search for the image by path/row and date in USGS Earth Explorer:
http://edcsns17.cr.usgs.gov/EarthExplorer/
In the results table, choose "Show All Fields" next to your image and scroll down to find Sun Elevation.
** If the Sun Elevation Angle for your scene is not available by any of these means, the approximate value may be obtained from the Ephemeris Generator described below for obtaining Earth-Sun Distance. In the table that is generated, the relevant value is under the heading of Sun_Elev.
The overpass time that you enter into the Ephemeris Generator can make an important difference to the sun elevation angle and, thus, to the reflectance values. Don't assume that the hypothetical 9:30 local time for a Landsat 5 overpass will apply to your scene - it seldom actually does.
If overpass time (Scene Start Time) is not included in the image header AND your scene is relatively recent (2001-), the overpass time for your scene can be estimated by using the NASA Earth Observatory Satellite Overpass Predictor (http://earthobservatory.nasa.gov/MissionControl/overpass.html):
UTC Zone 0 is Greenwich Mean Time
Arizona is UTC-7 (7 hours earlier than GMT)
Obtaining the Earth-Sun distance
Generate an ephemeris (http://ssd.jpl.nasa.gov/cgi-bin/eph/) for the overpass time indicated in EOS-DIS (Arizona time is the UT -7 time zone, but this may differ if you are using a scene from elsewhere) on the acquisition date for the center of the scene. The value of d (delta in the ephemeris output) is the Earth-Sun distance in astronomical units (AU) to enter into the model. The value will be between 0.9 and 1.1 with several decimal places.
References
Chavez, P. S., jr. 1996. Image-based atmospheric corrections
- Revisited and Improved. Photogrammetric Engineering and Remote
Sensing 62 (9): 1025-1036.
Markham, B.L., and J.L. Barker. 1986. Landsat MSS and TM post-calibration
dynamic ranges, exoatmospheric reflectances and at-satellite temperatures.
EOSAT Technical Notes, August 1986.
Moran, M.S., R.D. Jackson, P.N. Slater, and P.M. Teillet. 1992.
Evaluation of simplified procedures for retrieval of land surface
reflectance factors from satellite sensor output. Remote Sensing
of Environment 41:169-184.
.
|
Arizona Remote Sensing Center, Office of
Arid Lands Studies, Last updated: 30 December 2002 |