GIS Ag Maps.org

Landsat 7 & 5 Atmospheric Correction / Conversion to COST (Chavez, 1996) Surface Reflectance & Relative Scatter Lookup Table (for ArcGIS & QGIS)

Landsat surface reflectance (SR) can be ordered for free from the USGS; this tutorial shows how to convert to SR independently - please read: Background information about USGS & GIS Ag Maps surface reflectance.

Relevant Pages: Landsat Tiers & RMSE (offsite page)    Landsat 8 / Sentinel-2 Rare Comparison & Data

The steps below guide you through the process of converting Landsat 7 and 5 to COST surface reflectance (Chavez, 1996). Continuous relative scatter is applied here (suggested by Chavez [1988]). The steps below calculate COST SR in a slightly different manner than the equation (as written in previous link and the next link [original article), but the SR values calculated are the same (following the steps below is another way to calculate COST surface reflectance). IMPORTANT: COST SR should only be applied when solar elevation is > 45⁰. Chavez (1996) stated further research for the COST model is needed for solar elevations less than 35⁰ (solar zenith angles greater than 55⁰). (Access the Chavez [1996] COST article here for methods and results; PDF on this website.)

* IMPORTANT: If converting to SR for visible bands, try to use imagery with a solar elevation > 45⁰ (for mid-latitudes in the northern hemisphere, this is from about the middle of March to the middle of September for Landsat and Sentinel-2 [opposite in southern hemisphere). NIR and SWIR can be converted to SR at any solar elevation.Retrieved Landsat 8 SR with the methods described on this website have only been tested here for accuracy for solar elevations greater than 50⁰ and less than 30⁰; we have found that visible band, particularly the blue band, SR values become much less accurate (too high) at solar elevations less than 30⁰ (because TOA becomes disproportionatly higher for shorter wavelengths at low solar elevations). We have not completed research for solar elevations from 30 to 50⁰. It is very likely that visible band retrieved SR, particularly the blue band, becomes increasingly less accurate (too high) as solar elevations decrease from about 45⁰ - NIR AND SWIR BANDS REMAIN ACCURATE AT VERY LOW SOLAR ELEVATIONS. Chavez (1996) stated further research for the COST model is needed for solar elevations less than 35⁰ (solar zenith angles greater than 55⁰). CLICK HERE FOR LOW SOLAR ELEVATION IMAGERY SR COMPARISONS & DOWNLOADS.

TUTORIAL STEPS

1) Download Landsat 7 or 5 imagery from GloVis Next, which can be accessed on top menu or in the Data Sources folder to the left. You need account to download data (easy and free process). Convert the imagery to COST (Chavez, 1996) Top of Atmosphere (TOA) reflectance - for COST TOA, the denominator must have the squared cosine of the solar zenith (cosθs²). To do this access the MTL text file that is in the download (may be easier to copy the information in the text file and paste in into a Word document) and find the following information for each band you want to convert (in this example, which is Landsat 7, bands 1 - 4 will be converted):

● RADIANCE_MULT_BAND_1 = 7.7874E-01 (equals .77874); RADIANCE_MULT_BAND_2 = 7.9882E-01 (equals .79882); RADIANCE_MULT_BAND_3 = 6.2165E-01 (equals .62165); RADIANCE_MULT_BAND_4 = 9.6929E-01 equals .96929). These are Landsat 7 values and may be the same or vary for other Landsat 7 scenes, so it is recommended to do this on an image by image basis. Landsat 5 values are different than Landsat 7.

● RADIANCE_ADD_BAND_1 = -6.97874; RADIANCE_ADD_BAND_2 = -7.19882; RADIANCE_ADD_BAND_3 = -5.62165; RADIANCE_ADD_BAND_4 = -6.06929. These are Landsat 7 values and may be the same or vary for other Landsat 7 scenes, so it is recommended to do this on an image by image basis. Landsat 5 values are different than Landsat 7.

● EARTH_SUN_DISTANCE = 1.0162916 (older Landsat 5 or 7 may  not have this value in the MTL file; if not access this document [PDF on this website] to find the earth-sun distance based on imagery day of year).

● SUN_ELEVATION = 65.57069113

● You will need to use the value of pi (π) = 3.14159265

● You need to use the ESUN value for each band, which can be found in this document (PDF on this website; same as previous link). For example, the value for Landsat 7 band 3 is 1,533.

2) Calculating the Numerator in the COST TOA equation: To start converting to COST TOA reflectance, convert the Landsat 5 or 7 band raster (which should be in digital number [DN] format) to radiance by multiplying it with a Raster Calculator in GIS software by the proper "RADIANCE_MULT_BAND" value above and adding the proper "RADIANCE_ADD_BAND" value (which are all negative numbers, so you are essentially subtracting here). The equation to convert to radiance is, (DN x RADIANCE_MULT_BAND) + RADIANCE_ADD_BAND.

Multiply the radiance by pi (3.14159265) and by the square of the earth-sun distance (1.0328486), which is radiance x pi x earth-sun distance². You now have calculated the numerator in the COST TOA equation (scatter will be deducted during a later step in order to convert to surface reflectance).

3) Calculating the Denominator in the COST TOA equation: First acquire the proper band ESUN value for the correct satellite by viewing the document that can be accessed through the above link. If calculating for Landsat 7 band 3, the ESUN value is 1,533. You now need to convert the "SUN_ELEVATION" value from above to the solar zenith and then square the solar zenith. To do this you can use the following calculator or access a website below it:

Cosine of Solar Zenith (Embedded Excel) Calculator (Cloud-Based; MAY TAKE A FEW MOMENTS TO LOAD):

The cosine of the solar zenith can also be calculated at the following website (solar zenith = [90 - solar elevation]): http://www.rapidtables.com/calc/math/Cos_Calculator.htm; https://web2.0calc.com/

The value of cosine of the sun (solar) elevation from above (65.57069113) is 0.91047222. This value needs to be squared. The square of 0.91047222 = 0.82895966.

To calculate the denominator of the COST TOA equation, multiply the ESUN value by the square of the cosine. For band 3 in this example, the denominator = 1,533 x 0.8289566.

4) CALCULATE COST TOA REFLECTANCE by dividing the numerator from Step 2 by the denominator from Step 3.

IMPORTANT: We highly recommend the red band as the starting scatter band to base relative scatter on (especially if you are using the red band in an index with NIR), but it is not a requirement (another visible band can be used). Consistency in the method used is important, however.

5) Band 3 (red) is used for starting scatter in this tutorial. Select the band 3 scatter/haze digital number (DN) from the histogram, keeping in mind there are erroneous low, nonzero value along and near image edges (contact us if you have questions about this, as this is a vital step). You need to make sure that the scatter DN value selected is valid; do not select a value within the tail of the histogram (value should be near the base). View the the Chavez (1988) selection method here for guidance (the selected value is not necessarily the lowest value); if band 3 low values seems too incoherent to select a value from that represents the base of the histogram, choose a scatter value from another visible band that has a more distinct increase in values at the low end of the histogram.

6) Convert the starting haze DN value from Step 2 to COST TOA by the same method as in Step 1 - Step 4 (must also use the squared cosine of the solar zenith). This can be calculated in a spreadsheet or calculator then deducted from the image with a Raster Calculator in GIS software.

7) Deduct .01 from the scatter value from Step 3 (per Chavez [1996]); this will establish the low value having a surface reflectance of .01, because of the "fact that very few targets on Earth's surface are absolute black, so an assumed one-percent minimum reflectance is better than zero percent" (Chavez, 1996).

8) Find the closest value to one-percent reflectance from Step 4 and use the Landsat 5 & 7 Relative Scatter Table below to determine scatter for the other bands. Deduct the scatter from COST TOA using a raster calculator for the corresponding bands to calculate surface reflectance (use the precise starting band starting value and the relative values closest to it).

Landsat 7 & 5 Relative Scatter Lookup Table

L5 is Landsat 5; L7 is Landsat 7; B is Band: B1 is Blue, B2 is Green, B3 is Red, B4 is NIR; Band center wavelengths in micrometers are included. Landsat 5 & 7 have the same band center wavelengths for bands 1 - 3; band 4 is slightly different. The power correlation between band center wavelength and scatter value is R² = 1.0 in all cases (as it theoretically should be) and the power exponents are appropriate. When merged together, the power correlation between Landsat 5, 7, 8, and Sentinel-2 band center wavelength Relative Scatter Lookup Table values is R² = 1.0 in all cases (as it theoretically should be). The starting band 4 value of .01000 on the table below represents an unrealistically clear atmosphere; starting haze values will likely be from about .02 or more (after the .01 deduction described above). Landsat 5 & 7 have the same band center wave

 References

Chavez, P.S., Jr. 1996. Image-based atmospheric corrections–revisited and improved. Photogrammetric Engineering and Remote Sensing 62(9): pp.1025-1036.

Chavez, P.S., Jr. 1988. An improved dark-object subtraction technique for atmospheric scattering correction of multispectral data. Remote Sensing of Environment 24: pp.459-479.