Your company here — click to reach over 10,000 unique daily visitors

i.atcorr.1grass - Man Page


i.atcorr — Performs atmospheric correction using the 6S algorithm.
6S — Second Simulation of Satellite Signal in the Solar Spectrum.


imagery, atmospheric correction, radiometric conversion, radiance, reflectance, satellite


i.atcorr --help
i.atcorr [-irab] input=name  [range=min,max]   [elevation=name]   [visibility=name]  parameters=name output=name  [rescale=min,max]   [--overwrite]  [--help]  [--verbose]  [--quiet]  [--ui]



Output raster map as integer


Input raster map converted to reflectance (default is radiance)


Input from ETM+ image taken after July 1, 2000


Input from ETM+ image taken before July 1, 2000


Allow output files to overwrite existing files


Print usage summary


Verbose module output


Quiet module output


Force launching GUI dialog


input=name [required]

Name of input raster map


Input range
Default: 0,255


Name of input elevation raster map (in m)


Name of input visibility raster map (in km)

parameters=name [required]

Name of input text file with 6S parameters

output=name [required]

Name for output raster map


Rescale output raster map
Default: 0,255


i.atcorr performs atmospheric correction on the input raster map using the 6S algorithm (Second Simulation of Satellite Signal in the Solar Spectrum). A detailed algorithm description is available at the Land Surface Reflectance Science Computing Facility website.

Important: Current region settings are ignored! The region is adjusted to cover the input raster map before the atmospheric correction is performed. The previous settings are restored afterwards.

If the -r flag is used, the input raster map is treated as reflectance. Otherwise, the input raster map is treated as radiance values and it is converted to reflectance at the i.atcorr runtime. The output data are always reflectance.

The satellite overpass time has to be specified in Greenwich Mean Time (GMT).

An example of the 6S parameters could be:

8                            - geometrical conditions=Landsat ETM+
2 19 13.00 -47.410 -20.234   - month day hh.ddd longitude latitude ("hh.ddd" is in decimal hours GMT)
1                            - atmospheric model=tropical
1                            - aerosols model=continental
15                           - visibility [km] (aerosol model concentration)
-0.600                       - mean target elevation above sea level [km] (here 600 m asl)
-1000                        - sensor height (here, sensor on board a satellite)
64                           - 4th band of ETM+ Landsat 7

If the position is not available in longitude-latitude (WGS84), the m.proj conversion module can be used to reproject from a different reference system.

6S Code Parameter Choices

A. Geometrical conditions

1meteosat observationenter month,day,decimal hour (universal time-hh.ddd)                       n. of column,n. of line. (full scale 5000*2500) 
2goes east observationenter month,day,decimal hour (universal time-hh.ddd)                       n. of column,n. of line. (full scale 17000*12000)c
3goes west observationenter month,day,decimal hour (universal time-hh.ddd)                       n. of column,n. of line. (full scale 17000*12000)
4avhrr (PM noaa)enter month,day,decimal hour (universal time-hh.ddd)                       n. of column(1-2048),xlonan,hna                       give long.(xlonan) and overpass hour (hna) at                       the ascendant node at equator
5avhrr (AM noaa)enter month,day,decimal hour (universal time-hh.ddd)                       n. of column(1-2048),xlonan,hna                       give long.(xlonan) and overpass hour (hna) at                       the ascendant node at equator
6hrv (spot)enter month,day,hh.ddd,long.,lat. *
7tm (landsat)enter month,day,hh.ddd,long.,lat. *
8etm+ (landsat7)enter month,day,hh.ddd,long.,lat. *
9liss (IRS 1C)enter month,day,hh.ddd,long.,lat. *
10asterenter month,day,hh.ddd,long.,lat. *
11avnirenter month,day,hh.ddd,long.,lat. *
12ikonosenter month,day,hh.ddd,long.,lat. *
13RapidEyeenter month,day,hh.ddd,long.,lat. *
14VGT1 (SPOT4)enter month,day,hh.ddd,long.,lat. *
15VGT2 (SPOT5)enter month,day,hh.ddd,long.,lat. *
16WorldView 2enter month,day,hh.ddd,long.,lat. *
17QuickBirdenter month,day,hh.ddd,long.,lat. *
18LandSat 8enter month,day,hh.ddd,long.,lat. *
19Geoeye 1enter month,day,hh.ddd,long.,lat. *
20Spot6enter month,day,hh.ddd,long.,lat. *
21Spot7enter month,day,hh.ddd,long.,lat. *
22Pleiades1Aenter month,day,hh.ddd,long.,lat. *
23Pleiades1Benter month,day,hh.ddd,long.,lat. *
24Worldview3enter month,day,hh.ddd,long.,lat. *
25Sentinel-2Aenter month,day,hh.ddd,long.,lat. *
26Sentinel-2Benter month,day,hh.ddd,long.,lat. *
27PlanetScope 0c 0denter month,day,hh.ddd,long.,lat. *
28PlanetScope 0eenter month,day,hh.ddd,long.,lat. *
29PlanetScope 0f 10enter month,day,hh.ddd,long.,lat. *
30Worldview4enter month,day,hh.ddd,long.,lat. *

NOTE: for HRV, TM, ETM+, LISS and ASTER experiments, longitude and latitude are the coordinates of the scene center. Latitude must be > 0 for northern hemisphere and < 0 for southern. Longitude must be > 0 for eastern hemisphere and < 0 for western.

B. Atmospheric model

0no gaseous absorption
2midlatitude summer
3midlatitude winter
4subarctic summer
5subarctic winter
6us standard 62
7Define your own atmospheric model as a set of the following 5 parameters per each measurement: altitude [km] pressure [mb] temperature [k] h2o density [g/m3] o3 density [g/m3] For example: there is one radiosonde measurement for each altitude of 0-25km at a step of 1km, one measurement for each altitude of 25-50km at a step of 5km, and two single measurements for altitudes 70km and 100km. This makes 34 measurements. In that case, there are 34*5 values to input.
8Define your own atmospheric model providing values of the water vapor and ozone content: uw [g/cm2] uo3 [cm-atm] The profile is taken from us62.

C. Aerosols model

0no aerosols 
1continental model 
2maritime model 
3urban model 
4shettle model for background desert aerosol 
5biomass burning 
6stratospheric model 
7define your own modelEnter the volumic percentage of each component: c(1) = volumic % of dust-like c(2) = volumic % of water-soluble c(3) = volumic % of oceanic c(4) = volumic % of soot All values should be between 0 and 1.
8define your own modelSize distribution function: Multimodal Log Normal (up to 4 modes).
9define your own modelSize distribution function: Modified gamma.
10define your own modelSize distribution function: Junge Power-Law.
11define your own modelSun-photometer measurements, 50 values max, entered as: r and d V / d (logr) where r is the radius [micron], V is the volume, d V / d (logr) [cm3/cm2/micron]. Followed by: nr and ni for each wavelength where nr and ni are respectively the real and imaginary part of the refractive index.

D. Aerosol concentration model (visibility)

If you have an estimate of the meteorological parameter visibility v, enter directly the value of v [km] (the aerosol optical depth (AOD) will be computed from a standard aerosol profile).

If you have an estimate of aerosol optical depth, enter 0 for the visibility and in a following line enter the aerosol optical depth at 550nm (iaer means ’i’ for input and ’aer’ for aerosol), for example:

0                            - visibility
0.112                        - aerosol optical depth at 550 nm

NOTE: if iaer is 0, enter -1 for visibility.

NOTE: if a visibility map is provided, these parameters are ignored.

E. Target altitude (xps), sensor platform (xpp)

Target altitude (xps, in negative [km]): xps >= 0 means the target is at the sea level.
otherwise xps expresses the altitude of the target (e.g., mean elevation) in [km], given as negative value
Sensor platform (xpp, in negative [km] or -1000):
xpp = -1000 means that the sensor is on board a satellite.
xpp = 0 means that the sensor is at the ground level.
-100 < xpp < 0 defines the altitude of the sensor expressed in [km]; this altitude is given relative to the target altitude as negative value.

For aircraft simulations only (xpp is neither equal to 0 nor equal to -1000): puw,po3 (water vapor content,ozone content between the aircraft and the surface)
taerp (the aerosol optical thickness at 550nm between the aircraft and the surface)

If these data are not available, enter negative values for all of them. puw,po3 will then be interpolated from the us62 standard profile according to the values at the ground level; taerp will be computed according to a 2 km exponential profile for aerosol.

F. Sensor band

There are two possibilities: either define your own spectral conditions (codes -2, -1, 0, or 1) or choose a code indicating the band of one of the pre-defined satellites.

Define your own spectral conditions:

-2Enter wlinf, wlsup. The filter function will be equal to 1 over the whole band (as iwave=0) but step by step output will be printed.
-1Enter wl (monochr. cond, gaseous absorption is included).
0Enter wlinf, wlsup. The filter function will be equal to 1 over the whole band.
1Enter wlinf, wlsup and user’s filter function s (lambda) by step of 0.0025 micrometer.

Pre-defined satellite bands:

CodeBand name (peak response)
2meteosat vis band (0.350-1.110)
3goes east band vis (0.490-0.900)
4goes west band vis (0.490-0.900)
5avhrr (noaa6) band 1 (0.550-0.750)
6avhrr (noaa6) band 2 (0.690-1.120)
7avhrr (noaa7) band 1 (0.500-0.800)
8avhrr (noaa7) band 2 (0.640-1.170)
9avhrr (noaa8) band 1 (0.540-1.010)
10avhrr (noaa8) band 2 (0.680-1.120)
11avhrr (noaa9) band 1 (0.530-0.810)
12avhrr (noaa9) band 1 (0.680-1.170)
13avhrr (noaa10) band 1 (0.530-0.780)
14avhrr (noaa10) band 2 (0.600-1.190)
15avhrr (noaa11) band 1 (0.540-0.820)
16avhrr (noaa11) band 2 (0.600-1.120)
17hrv1 (spot1) band 1 (0.470-0.650)
18hrv1 (spot1) band 2 (0.600-0.720)
19hrv1 (spot1) band 3 (0.730-0.930)
20hrv1 (spot1) band pan (0.470-0.790)
21hrv2 (spot1) band 1 (0.470-0.650)
22hrv2 (spot1) band 2 (0.590-0.730)
23hrv2 (spot1) band 3 (0.740-0.940)
24hrv2 (spot1) band pan (0.470-0.790)
25tm (landsat5) band 1 (0.430-0.560)
26tm (landsat5) band 2 (0.500-0.650)
27tm (landsat5) band 3 (0.580-0.740)
28tm (landsat5) band 4 (0.730-0.950)
29tm (landsat5) band 5 (1.5025-1.890)
30tm (landsat5) band 7 (1.950-2.410)
31mss (landsat5) band 1 (0.475-0.640)
32mss (landsat5) band 2 (0.580-0.750)
33mss (landsat5) band 3 (0.655-0.855)
34mss (landsat5) band 4 (0.785-1.100)
35MAS (ER2) band 1 (0.5025-0.5875)
36MAS (ER2) band 2 (0.6075-0.7000)
37MAS (ER2) band 3 (0.8300-0.9125)
38MAS (ER2) band 4 (0.9000-0.9975)
39MAS (ER2) band 5 (1.8200-1.9575)
40MAS (ER2) band 6 (2.0950-2.1925)
41MAS (ER2) band 7 (3.5800-3.8700)
42MODIS band 1 (0.6100-0.6850)
43MODIS band 2 (0.8200-0.9025)
44MODIS band 3 (0.4500-0.4825)
45MODIS band 4 (0.5400-0.5700)
46MODIS band 5 (1.2150-1.2700)
47MODIS band 6 (1.6000-1.6650)
48MODIS band 7 (2.0575-2.1825)
49avhrr (noaa12) band 1 (0.500-1.000)
50avhrr (noaa12) band 2 (0.650-1.120)
51avhrr (noaa14) band 1 (0.500-1.110)
52avhrr (noaa14) band 2 (0.680-1.100)
53POLDER band 1 (0.4125-0.4775)
54POLDER band 2 (non polar) (0.4100-0.5225)
55POLDER band 3 (non polar) (0.5325-0.5950)
56POLDER band 4 P1 (0.6300-0.7025)
57POLDER band 5 (non polar) (0.7450-0.7800)
58POLDER band 6 (non polar) (0.7000-0.8300)
59POLDER band 7 P1 (0.8100-0.9200)
60POLDER band 8 (non polar) (0.8650-0.9400)
61etm+ (landsat7) band 1 blue (435nm - 517nm)
62etm+ (landsat7) band 2 green (508nm - 617nm)
63etm+ (landsat7) band 3 red (625nm - 702nm)
64etm+ (landsat7) band 4 NIR (753nm - 910nm)
65etm+ (landsat7) band 5 SWIR (1520nm - 1785nm)
66etm+ (landsat7) band 7 SWIR (2028nm - 2375nm)
67etm+ (landsat7) band 8 PAN (505nm - 917nm)
68liss (IRC 1C) band 2 (0.502-0.620)
69liss (IRC 1C) band 3 (0.612-0.700)
70liss (IRC 1C) band 4 (0.752-0.880)
71liss (IRC 1C) band 5 (1.452-1.760)
72aster band 1 (0.480-0.645)
73aster band 2 (0.588-0.733)
74aster band 3N (0.723-0.913)
75aster band 4 (1.530-1.750)
76aster band 5 (2.103-2.285)
77aster band 6 (2.105-2.298)
78aster band 7 (2.200-2.393)
79aster band 8 (2.248-2.475)
80aster band 9 (2.295-2.538)
81avnir band 1 (408nm - 517nm)
82avnir band 2 (503nm - 612nm)
83avnir band 3 (583nm - 717nm)
84avnir band 4 (735nm - 922nm)
85Ikonos Green band (408nm - 642nm)
86Ikonos Red band (448nm - 715nm)
87Ikonos NIR band (575nm - 787nm)
88RapidEye Blue band (440nm - 512nm)
89RapidEye Green band (515nm - 592nm)
90RapidEye Red band (628nm - 687nm)
91RapidEye Red edge band (685nm - 735nm)
92RapidEye NIR band (750nm - 860nm)
93VGT1 (SPOT4) band 0 (420nm - 497nm)
94VGT1 (SPOT4) band 2 (603nm - 747nm)
95VGT1 (SPOT4) band 3 (740nm - 942nm)
96VGT1 (SPOT4) MIR band (1540nm - 1777nm)
97VGT2 (SPOT5) band 0 (423nm - 492nm)
98VGT2 (SPOT5) band 2 (600nm - 737nm)
99VGT2 (SPOT5) band 3 (745nm - 945nm)
100VGT2 (SPOT5) MIR band (1523nm - 1757nm)
101WorldView2 Panchromatic band (448nm - 812nm)
102WorldView2 Coastal Blue band (395nm - 457nm)
103WorldView2 Blue band (440nm - 517nm)
104WorldView2 Green band (503nm - 587nm)
105WorldView2 Yellow band (583nm - 632nm)
106WorldView2 Red band (623nm - 695nm)
107WorldView2 Red edge band (698nm - 750nm)
108WorldView2 NIR1 band (760nm - 905nm)
109WorldView2 NIR2 band (853nm - 1047nm)
110QuickBird Panchromatic band (385nm - 1060nm)
111QuickBird Blue band (420nm - 585nm)
112QuickBird Green band (448nm - 682nm)
113QuickBird Red band (560nm - 747nm)
114QuickBird NIR1 band (650nm - 935nm)
115Landsat 8 Coastal aerosol band (433nm - 455nm)
116Landsat 8 Blue band (448nm - 515nm)
117Landsat 8 Green band (525nm - 595nm)
118Landsat 8 Red band (633nm - 677nm)
119Landsat 8 Panchromatic band (498nm - 682nm)
120Landsat 8 NIR band (845nm - 885nm)
121Landsat 8 Cirrus band (1355nm - 1390nm)
122Landsat 8 SWIR1 band (1540nm - 1672nm)
123Landsat 8 SWIR2 band (2073nm - 2322nm)
124GeoEye 1 Panchromatic band (448nm - 812nm)
125GeoEye 1 Blue band (443nm - 525nm)
126GeoEye 1 Green band (503nm - 587nm)
127GeoEye 1 Red band (653nm - 697nm)
128GeoEye 1 NIR band (770nm - 932nm)
129Spot6 Blue band (440nm - 532nm)
130Spot6 Green band (515nm - 600nm)
131Spot6 Red band (610nm - 710nm)
132Spot6 NIR band (738nm - 897nm)
133Spot6 Pan band (438nm - 760nm)
134Spot7 Blue band (445nm - 532nm)
135Spot7 Green band (525nm - 607nm)
136Spot7 Red band (610nm - 727nm)
137Spot7 NIR band (745nm - 902nm)
138Spot7 Pan band (443nm - 760nm)
139Pleiades1A Blue band (433nm - 560nm)
140Pleiades1A Green band (500nm - 617nm)
141Pleiades1A Red band (590nm - 722nm)
142Pleiades1A NIR band (740nm - 945nm)
143Pleiades1A Pan band (460nm - 845nm)
144Pleiades1B Blue band 438nm - 560nm)
145Pleiades1B Green band (498nm - 615nm)
146Pleiades1B Red band (608nm - 727nm)
147Pleiades1B NIR band (750nm - 945nm)
148Pleiades1B Pan band (460nm - 845nm)
149Worldview3 Pan band (445nm - 812nm)
150Worldview3 Coastal blue band (395nm - 455nm)
151Worldview3 Blue band (443nm - 517nm)
152Worldview3 Green band (508nm - 587nm)
153Worldview3 Yellow band (580nm - 630nm)
154Worldview3 Red band (625nm - 697nm)
155Worldview3 Red edge band (698nm - 752nm)
156Worldview3 NIR1 band (760nm - 902nm)
157Worldview3 NIR2 band (855nm - 1042nm)
158Worldview3 SWIR1 band (1178nm - 1242nm)
159Worldview3 SWIR2 band (1545nm - 1600nm)
160Worldview3 SWIR3 band (1633nm - 1687nm)
161Worldview3 SWIR4 band (1698nm - 1762nm)
162Worldview3 SWIR5 band (2133nm - 2195nm)
163Worldview3 SWIR6 band (2170nm - 2235nm)
164Worldview3 SWIR7 band (2225nm - 2295nm)
165Worldview3 SWIR8 band (2283nm - 2377nm)
166Sentinel2A Coastal blue band B1 (430nm - 455nm)
167Sentinel2A Blue band B2 (440nm - 530nm)
168Sentinel2A Green band B3 (540nm - 580nm)
169Sentinel2A Red band B4 (648nm - 682nm)
170Sentinel2A Red edge band B5 (695nm - 712nm)
171Sentinel2A Red edge band B6 (733nm - 747nm)
172Sentinel2A Red edge band B7 (770nm - 795nm)
173Sentinel2A NIR band B8 (775nm - 905nm)
174Sentinel2A Red edge band B8A (850nm - 880nm)
175Sentinel2A Water vapour band B9 (933nm - 957nm)
176Sentinel2A SWIR Cirrus band B10 (1355nm - 1392nm)
177Sentinel2A SWIR band B11 (1558nm - 1667nm)
178Sentinel2A SWIR band B12 (2088nm - 2315nm)
179Sentinel2B Coastal blue band B1 (430nm - 455nm)
180Sentinel2B Blue band B2 (440nm - 530nm)
181Sentinel2B Green band B3 (538nm - 580nm)
182Sentinel2B Red band B4 (648nm - 682nm)
183Sentinel2B Red edge band B5 (695nm - 712nm)
184Sentinel2B Red edge band B6 (730nm - 747nm)
185Sentinel2B Red edge band B7 (768nm - 792nm)
186Sentinel2B NIR band B8 (778nm - 905nm)
187Sentinel2B Red edge band B8A (850nm - 877nm)
188Sentinel2B Water vapour band B9 (930nm - 955nm)
189Sentinel2B SWIR Cirrus band B10 (1358nm - 1397nm)
190Sentinel2B SWIR band B11 (1555nm - 1667nm)
191Sentinel2B SWIR band B12 (2075nm - 2300nm)
192PlanetScope 0c 0d Blue band B1 (440nm - 570nm)
193PlanetScope 0c 0d Green band B2 (450nm - 690nm)
194PlanetScope 0c 0d Red band B3 (460nm - 700nm)
195PlanetScope 0c 0d NIR band B4 (770nm - 880nm)
196PlanetScope 0e Blue band B1 (430nm - 700nm)
197PlanetScope 0e Green band B2 (450nm - 700nm)
198PlanetScope 0e Red band B3 (460nm - 700nm)
199PlanetScope 0e NIR band B4 (760nm - 880nm)
200PlanetScope 0f 10 Blue band B1 (450nm - 680nm)
201PlanetScope 0f 10 Green band B2 (450nm - 680nm)
202PlanetScope 0f 10 Red band B3 (450nm - 680nm)
203PlanetScope 0f 10 NIR band B4 (760nm - 870nm)
204Worldview4 Pan band (424nm - 842nm)
205Worldview4 Blue band (416nm - 567nm)
206Worldview4 Green band (488nm - 626nm)
207Worldview4 Red band (639nm - 711nm)
208Worldview4 NIR1 band (732nm - 962nm)


Atmospheric correction of a Sentinel-2 band

This example illustrates how to perform atmospheric correction of a Sentinel-2 scene in the North Carolina location.

Let’s assume that the Sentinel-2 L1C scene S2A_OPER_PRD_MSIL1C_PDMC_20161029T092602_R054_V20161028T155402_20161028T155402 was downloaded and imported with region cropping (see r.import) into the PERMANENT mapset of the North Carolina location. The computational region was set to the extent of the elevation map in the North Carolina dataset. Now, we have 13 individual bands (B01-B12) that we want to apply the atmospheric correction to. The following steps are applied to each band separately.

Create the parameters file for i.atcorr

In the first step we create a file containing the 6S parameters for a particular scene and band. To create a 6S file, we need to obtain the following information:

  • geometrical conditions,
  • moth, day, decimal hours in GMT, decimal longitude and latitude of measurement,
  • atmospheric model,
  • aerosol model,
  • visibility or aerosol optical depth,
  • mean target elevation above sea level,
  • sensor height and,
  • sensor band.

Geometrical conditions

For Sentinel-2A, the geometrical conditions take the value 25 and for Sentinel-2B, the geometrical conditions value is 26 (See table A). Our scene comes from the Sentinel-2A mission (the file name begins with S2A_...).


Day, time, longitude and latitude of measurement

Day and time of the measurement are hidden in the filename (i.e., the second datum in the file name with format YYYYMMDDTHHMMSS), and are also noted in the metadata file, which is included in the downloaded scene (file with .xml extension). Our sample scene was taken on October 28th (20161028) at 15:54:02 (155402). Note that the time has to be specified in decimal hours in Greenwich Mean Time (GMT). Luckily, the time in the scene name is in GMT and we can convert it to decimal hours as follows: 15 + 54/60 + 2/3600 = 15.901.

Longitude and latitude refer to the centre of the computational region (which can be smaller than the scene), and must be in WGS84 decimal coordinates. To obtain the coordinates of the centre, we can run:

g.region -bg

The longitude and latitude of the centre are stored in ll_clon and ll_clat. In our case, ll_clon=-78.691 and ll_clat=35.749.


Atmospheric model

We can choose between various atmospheric models as defined at the beginning of this manual. For North Carolina, we can choose 2 - midlatitude summer.


Aerosol model

We can also choose between various aerosol models as defined at the beginning of this manual. For North Carolina, we can choose 1 - continental model.


Visibility or Aerosol Optical Depth

For Sentinel-2 scenes, the visibility is not measured, and therefore we have to estimate the aerosol optical depth instead, e.g. from AERONET. With a bit of luck, you can find a station nearby your location, which measured the Aerosol Optical Depth at 500 nm at the same time as the scene was taken. In our case, on 28th October 2016, the EPA-Res_Triangle_Pk station measured AOD = 0.07 (approximately).


Mean target elevation above sea level

Mean target elevation above sea level refers to the mean elevation of the computational region. You can estimate it from the digital elevation model, e.g. by running:

r.univar -g elevation

The mean elevation is stored in mean. In our case, mean=110. In the 6S file it will be displayed in [-km], i.e., -0.110.


Sensor height

Since the sensor is on board a satellite, the sensor height will be set to -1000.


Sensor band

The overview of satellite bands can be found in table F (see above). For Sentinel-2A, the band numbers span from 166 to 178, and for Sentinel-2B, from 179 to 191.

Finally, here is what the 6S file would look like for Band 02 of our scene. In order to use it in the i.atcorr module, we can save it in a text file, for example params_B02.txt.

10 28 15.901 -78.691 35.749

Compute atmospheric correction

In the next step we run i.atcorr for the selected band B02 of our Sentinel 2 scene. We have to specify the following parameters:

  • input = raster band to be processed,
  • parameters = path to 6S file created in the previous step (we could also enter the values directly),
  • output = name for the output corrected raster band,
  • range = from 1 to the QUANTIFICATION_VALUE stored in the metadata file. It is 10000 for both Sentinel-2A and Sentinel-2B.
  • rescale = the output range of values for the corrected bands. This is up to the user to choose, for example: 0-255, 0-1, 1-10000.

If the data is available, the following parameters can be specified as well:

  • elevation = raster of digital elevation model,
  • visibility = raster of visibility model.

Finally, this is how the command would look like to apply atmospheric correction to band B02:

i.atcorr input=B02 parameters=params_B02.txt output=B02.atcorr range=1,10000 rescale=0,255 elevation=elevation

To apply atmospheric correction to the remaining bands, only the last line in the 6S parameters file (i.e., the sensor band) needs to be changed. The other parameters will remain the same.
Figure: Sentinel-2A Band 02 with applied atmospheric correction (histogram equalization grayscale color scheme)

Atmospheric correction of a Landsat-7 band

This example is also based on the North Carolina sample dataset (GMT -5 hours). First we set the computational region to the satellite map, e.g. band 4:

g.region raster=lsat7_2002_40 -p

It is important to verify the available metadata for the sun position which has to be defined for the atmospheric correction. An option is to check the satellite overpass time with sun position as reported in the metadata file (file copy; North Carolina sample dataset). In the case of the North Carolina sample dataset, these values have been stored for each channel and can be retrieved with:

r.info lsat7_2002_40

In this case, we have: SUN_AZIMUTH = 120.8810347, SUN_ELEVATION = 64.7730999.

If the sun position metadata are unavailable, we can also calculate them from the overpass time as follows (r.sunmask uses SOLPOS):

r.sunmask -s elev=elevation out=dummy year=2002 month=5 day=24 hour=10 min=42 sec=7 timezone=-5
# .. reports: sun azimuth: 121.342461, sun angle above horz.(refraction corrected): 65.396652

If the overpass time is unknown, use the NASA LaRC Satellite Overpass Predictor.

Convert digital numbers (DN) to radiance at top-of-atmosphere (TOA)

For Landsat and ASTER, the conversion can be conveniently done with i.landsat.toar or i.aster.toar, respectively.

In case of different satellites, the conversion of DN (digital number = pixel values) to radiance at top-of-atmosphere (TOA) can also be done manually, using e.g. the formula:

# formula depends on satellite sensor, see respective metadata
  • Lλ = Spectral Radiance at the sensor’s aperture in Watt/(meter squared * ster * µm), the apparent radiance as seen by the satellite sensor;
  • QCAL = the quantized calibrated pixel value in DN;
  • LMINλ = the spectral radiance that is scaled to QCALMIN in watts/(meter squared * ster * µm);
  • LMAXλ = the spectral radiance that is scaled to QCALMAX in watts/(meter squared * ster * µm);
  • QCALMIN = the minimum quantized calibrated pixel value (corresponding to LMINλ) in DN;
  • QCALMAX = the maximum quantized calibrated pixel value (corresponding to LMAXλ) in DN=255.

LMINλ and LMAXλ are the radiances related to the minimal and maximal DN value, and they are reported in the metadata file of each image. High gain or low gain is also reported in the metadata file of each satellite image. For Landsat ETM+, the minimal DN value (QCALMIN) is 1 (see Landsat handbook, chapter 11), and the maximal DN value (QCALMAX) is 255. QCAL is the DN value for every separate pixel in the Landsat image.

We extract the coefficients and apply them in order to obtain the radiance map:

r.info lsat7_2002_${CHAN}0 -h | tr ’\n’ ’ ’ | sed ’s+ ++g’ | tr ’:’ ’\n’ | grep "LMIN_BAND${CHAN}\|LMAX_BAND${CHAN}"

Conversion to radiance (this calculation is done for band 4, for the other bands, the numbers will need to be replaced with their related values):

r.mapcalc "lsat7_2002_40_rad = ((241.1 - (-5.1)) / (255.0 - 1.0)) * (lsat7_2002_40 - 1.0) + (-5.1)"

Again, the r.mapcalc calculation is only needed when working with satellite data other than Landsat or ASTER.

Create the parameters file for i.atcorr

The underlying 6S model is parametrized through a control file, indicated with the parameters option. This is a text file defining geometrical and atmospherical conditions of the satellite overpass. Here we create a control file icnd_lsat4.txt for band 4 (NIR), based on metadata. For the overpass time, we need to define decimal hours: 10:42:07 NC local time = 10.70 decimal hours (decimal minutes: 42 * 100 / 60) which is 15.70 GMT.

8                            - geometrical conditions=Landsat ETM+
5 24 15.70 -78.691 35.749    - month day hh.ddd longitude latitude ("hh.ddd" is in GMT decimal hours)
2                            - atmospheric model=midlatitude summer
1                            - aerosols model=continental
50                           - visibility [km] (aerosol model concentration)
-0.110                       - mean target elevation above sea level [km]
-1000                        - sensor on board a satellite
64                           - 4th band of ETM+ Landsat 7

Finally, run the atmospheric correction (-r for reflectance input map; -a for date > July 2000):

i.atcorr -r -a lsat7_2002_40_rad elevation=elevation parameters=icnd_lsat4.txt output=lsat7_2002_40_atcorr

Note that the altitude value from ’icnd_lsat4.txt’ file is read at the beginning to compute the initial transform. Therefore, it is necessary to provide a value that might be the mean value of the elevation model (r.univar elevation). For the atmospheric correction per se, the elevation values from the raster map are used.

Note that the process is computationally intensive. Note also, that i.atcorr reports solar elevation angle above horizon rather than solar zenith angle.

Remaining Documentation Issues

The influence and importance of the visibility value or map should be explained, also how to obtain an estimate for either visibility or aerosol optical depth at 550nm.


See Also

GRASS Wiki page about Atmospheric correction

i.aster.toar, i.colors.enhance, i.landsat.toar, r.info, r.mapcalc, r.univar


Original version of the program for GRASS 5:
Christo Zietsman, 13422863(at)sun.ac.za

Code clean-up and port to GRASS 6.3, 15.12.2006:
Yann Chemin, ychemin(at)gmail.com

Documentation clean-up + IRS LISS sensor addition 5/2009:
Markus Neteler, FEM, Italy

ASTER sensor addition 7/2009:
Michael Perdue, Canada

AVNIR, IKONOS sensors addition 7/2010:
Daniel Victoria, Anne Ghisla

RapidEye sensors addition 11/2010:
Peter Löwe, Anne Ghisla

VGT1 and VGT2 sensors addition from 6SV-1.1 sources, addition 07/2011:
Alfredo Alessandrini, Anne Ghisla

Added Landsat 8 from NASA sources, addition 05/2014:
Nikolaos Ves

Geoeye1 addition 7/2015:
Marco Vizzari

Worldview3 addition 8/2016:
Markus Neteler, mundialis.de, Germany

Sentinel-2A addition 12/2016:
Markus Neteler, mundialis.de, Germany

Sentinel-2B addition 1/2018:
Stefan Blumentrath, Zofie Cimburova, Norwegian Institute for Nature Research, NINA, Oslo, Norway

Worldview4 addition 12/2018:
Markus Neteler, mundialis.de, Germany

Source Code

Available at: i.atcorr source code (history)

Accessed: Tuesday May 14 13:41:55 2024

Main index | Imagery index | Topics index | Keywords index | Graphical index | Full index

© 2003-2024 GRASS Development Team, GRASS GIS 8.3.2 Reference Manual


GRASS 8.3.2 GRASS GIS User's Manual