Atmospheric Correction of Landsat Data using ARCSI

The Atmospheric and Radiometric Correction of Satellite Imagery (ARCSI) software was developed at Aberystwyth University, under contract from the Norwegian Space Centre and allows automatic correction of Landsat data to Top Of Atmosphere (TOA) or Surface Reflectance (SREF). The conversion to Surface Reflectance uses the 6S codes, accessed through the Py6S Python interface.

Getting ARCSI

If you’re running OS X or Linux, you can get ARCSI through conda. For Windows users, unless you’re very familiar with compiling software from source it is recommended you set up a virtual machine with Linux to run ARCSI.

For the latest installation instructions see the software page.

Download and untar data

This tutorial uses Landsat 8 data (although ARCSI supports other sensors). If you don’t already have a Landsat 8 scene, you can download one from http://earthexplorer.usgs.gov/.

The data is downloaded as a tarred and compressed file. After you have downloaded a scene, create a directory for the scene and copy the tar.gz file to it. Untar the file using:

tar -xf LC80430342013291LGN00.tar.gz

There are a number of tiff files (one for each each band) and associated header information in the MTL file.

After you have downloaded and untared the data create a new directory to store the output files (I’m using one called ‘OutputImages’).

Conversion to Radiance

Landsat data are provided as digital number (DN) images, the first step is to calibrate to radiance, which is the amount of energy received by the sensor per second (W) per steradian (sr) per square metre (m^{2}). This is written as W sr^{-1}m^{-2} or W sr^{-1}m^{-2}nm^{-1} for a given wavelength (in nm). This can also be defined as “The power passing through a unit area in a unit solid angle about the normal to the area (per unit spectral interval)”.

For Landsat 8, a linear equation is used to convert from DN to radiance (L):

L = (\text{DN} \times G) + I.

The gain (G) and intercept (I) are different for each band and are provided in the Landsat header file.

To convert each band to radiance in ARCSI, and stack all bands together, the following command is used:

arcsi.py -s ls8 -f KEA --stats -p RAD -o ./OutputImages \
  -i LC80430342013291LGN00/LC80430342013291LGN00_MTL.txt

Where ‘-o’ is the output directory and ‘-i’ the input header (MTL) file. Statistics and overviews (‘–stats’) are calculated for fast display.

Once this has run open the output file (ending in ‘_rad.kea’) in TuiView.

Top of Atmosphere Reflectance

When spectral radiance values are compared within and between sensors, variations frequently occur because of different sun-target-sensor configurations. Hence, these data are often converted to Top of Atmosphere (TOA) or planetary reflectance or at sensor reflectance as this eliminates the effect of different solar zenith angles at the time of image acquisition, minimises the solar irradiance band differences occurring between different sensors accounts for differences in the Earth-Sun distance. Conversion to TOA reflectance is achieved using:

\rho_{\lambda} = \frac{\pi L_{\lambda}d^2}{\text{ESUN}_{\lambda}cos(\theta_s)},

where \lambda is the wavelength, \rho_{\lambda} is the spectral (planetary or TOA) reflectance for wavelength \lambda, L_{\lambda} is the spectral radiance, d is the Earth-Sun distance (in astronomical units), ESUN is the mean solar exoatmospheric irradiance and \theta_s is the solar zenith angle.

TOA reflectance does not take into account any atmospheric effects and is simply a ratio of the at-sensor radiance with the incoming radiance from the sun, where the distance and angle of the sensor to the target are considered. The approach is easy to apply, as no external information is required.

To calibrate Landsat 8 data to TOA reflectance in ARCSI the following command is used:

arcsi.py -s ls8 -f KEA --stats -p RAD TOA -o ./OutputImages \
   -i LC80430342013291LGN00/LC80430342013291LGN00_MTL.txt

Once this has run, open the output file (ending in ‘_toa.kea’) in TuiView and compare with the radiance image (in a separate window). Note the output reflectance, which has values 0 – 1, is multiplied by a scaling factor of 1000 so the results can be stored as a 16 bit integer.

Surface Reflectance

To get to surface reflectance, an atmospheric model can be used to account for the effects of scattering and absorption in the atmosphere. Within ARCSI the 6S radiative transfer model is used. A number of parameters are required for 6S:

  • Atmospheric Profile
  • Aerosol Profile + Aerosol Optical Depth (AOD)
  • Ground Reflectance
  • Geometry
  • Altitude
  • Wavelength

Many of these parameters are calculated within ARCSI or derived from the header file (e.g., geometry, altitude and wavelength). There are a number of preset options for atmospheric profile, aerosol profile and and ground reflectance in ARCSI, to show all available options type:

arcsi.py -h

As an example we’ll run using the standard input parameters:

arcsi.py -s ls8 -f KEA --stats -p RAD SREF --aeropro Maritime \
  --atmospro MidlatitudeSummer --aot 0.25  -o ./OutputImages \
  -i LC80430342013291LGN00/LC80430342013291LGN00_MTL.txt

This will take a bit longer to run through. Once it’s finished open the top of atmosphere and surface reflectance images in linked TuiView displays and compare the profiles.

Comparison of Landsat data corrected to top of atmosphere (left) and surface reflectance (right).

Additional tasks

  1. Try changing the values of AOT (–aot) between 0 and 1, remembering to create separate output directories and observe the changes in reflectance over the range of wavelengths.
  2. You can investigate the different settings available in 6S, without applying to an image through Py6S, the quick start guide is a good place to start.

Getting help with ARCSI

The documentation for ARCSI is available at http://rsgislib.org/arcsi/

If you have any questions about ARCSI which are not covered by the documentation the RSGISLib mailing list on Google Groups is used for ARCSI support (rsgislib-support@googlegroups.com). Note you need to be a member of the group to post.

This post was adapted from notes for a module taught as part of the Masters course in GIS and Remote Sensing at Aberystwyth University, you can read more about the program here.

15 thoughts on “Atmospheric Correction of Landsat Data using ARCSI

  1. Zachary Bearss

    First off, I would like to thank you for such a wonderful write-up. It has saved me a lot of time with atmospheric correction on my imagery. When running the SREFSTDMDL command with Landsat 8 imagery, the quality of my RGB bands is improved drastically. It appears, however, that my NIR band gets extremely washed out which has a negative effect on my NDVI processing. I was wondering if you had experienced this before or had any suggestions on what I could try to change to remedy the problem. (In case you were wondering, the tile I am using is located in the Landsat 8 LC80140352014139LGN00 zip)

    Thank you for your time!

    Reply
  2. Janine Cole

    Please forgive me if this is a silly question, but do you apply a scale factor to the surface reflectance? In the images you show the values reach 500, but shouldn’t it be between 0 and 100% ?

    Thank you for making this software available

    Reply
  3. Juan

    Hi, thanks for this tutorial. I’m new at Linux and it was very helpfull.
    I have one question, How do I open an image in Tuiview?
    Can you tell me the command?

    Thanks!

    Reply
    1. Juan

      I dig a little more and now I understand how open images in Linux and its language.
      I have a problem now with the DEM (in kea and tif format I have tried) when I tried to Build the ARCSI command to process a batch of Landsat 8’s files. The DEM is not recognized. Any idea that what could be happening?

      Reply
  4. Carlos Riano

    Hi,

    I am using the CLOUDS option for the LT51710582010012MLK00. The output is fine, but I have some “fake” clouds. Can you help me to understand which variable is affecting my result.

    Thanks,

    Carlos Riano

    Reply
  5. Silke

    Hi,

    thank you for providing the ARCSI installation package! When I first installed the software (on ubuntu 14.04) via conda and using the command given above I ran into problems when using the newest miniconda version (miniconda3, based on python 3.5). Miniconda3 automatically installs python3.5, which conflicts with the tuiview installation. Furthermore, when running arcsi.py some library files have wrong versions. These problems vanish when the “old” niniconda2 (python 2.7) installer is used. Are these problems known? Are there any plans to update the package?

    Best wishes
    Silke

    Reply
  6. Silke

    Hi,

    one more question concerning the units of the output (… given on this website so I hope it’s correct to post it here…(?)). I assume that the radiance produced by arcsi is given in W/(sr*m^2*nm) (in the formulas in the section “conversion to radiance” it says per m^3, I guess that’s a typo?). Furthermore, I understand that TOA and SREF BOTH give the fraction of the radiation that is reflected, i.e. are values between 0 and 1, scaled by a factor 1000. Correct? Thanks!

    Silke

    Reply

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s