|
@@ -7,37 +7,34 @@ available at the
|
|
|
<a href="http://modis-sr.ltdri.org/">Land Surface
|
|
|
Reflectance Science Computing Facility website</a>.
|
|
|
|
|
|
-<p><em>Important note: Current region settings are ignored!</em> The region is adjusted
|
|
|
-to cover the input raster map before the atmospheric correction is
|
|
|
-performed. The previous settings are restored afterwards.
|
|
|
+<p><em>Important: Current region settings are ignored!</em> The
|
|
|
+region is adjusted to cover the input raster map before the atmospheric
|
|
|
+correction is performed. The previous settings are restored afterwards.
|
|
|
|
|
|
-This flag tells <em>i.atcorr</em> to try and speedup calculations.
|
|
|
-However, this option will increase memory requirements.
|
|
|
-
|
|
|
-<p>If flag <b>-r</b> is used, the input raster data are treated as
|
|
|
-<em>reflectance</em>. Otherwise, the input raster data are treated
|
|
|
-as <em>radiance</em> values and are converted to reflectance at
|
|
|
+<p>If the <b>-r</b> flag is used, the input raster map is treated as
|
|
|
+<em>reflectance</em>. Otherwise, the input raster map is treated
|
|
|
+as <em>radiance</em> values and it is converted to reflectance at
|
|
|
the <em>i.atcorr</em> runtime. The output data are always reflectance.
|
|
|
|
|
|
-<p>Note that the satellite overpass time has to be specified in Greenwich
|
|
|
+<p>The satellite overpass time has to be specified in Greenwich
|
|
|
Mean Time (GMT).
|
|
|
|
|
|
-<p>An example 6S parameters:
|
|
|
+<p>An example of the 6S parameters could be:
|
|
|
|
|
|
<div class="code"><pre>
|
|
|
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 mode=tropical
|
|
|
+1 - atmospheric model=tropical
|
|
|
1 - aerosols model=continental
|
|
|
15 - visibility [km] (aerosol model concentration)
|
|
|
--0.600 - mean target elevation above sea level [km] (here 600m asl)
|
|
|
+-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
|
|
|
</pre></div>
|
|
|
|
|
|
If the position is not available in longitude-latitude (WGS84), the
|
|
|
<em><a href="m.proj.html">m.proj</a></em> conversion module can be
|
|
|
-used to reproject from a different projection.
|
|
|
+used to reproject from a different reference system.
|
|
|
|
|
|
<h2>6S CODE PARAMETER CHOICES</h2>
|
|
|
|
|
@@ -377,7 +374,7 @@ ozone content:
|
|
|
<br>c(3) = volumic % of oceanic
|
|
|
<br>c(4) = volumic % of soot
|
|
|
<br>
|
|
|
-<br>All values between 0 and 1.</td>
|
|
|
+<br>All values should be between 0 and 1.</td>
|
|
|
</tr>
|
|
|
|
|
|
<tr>
|
|
@@ -400,9 +397,7 @@ ozone content:
|
|
|
|
|
|
<tr>
|
|
|
<td>11</td>
|
|
|
-
|
|
|
<td>define your own model</td>
|
|
|
-
|
|
|
<td>Sun-photometer measurements, 50 values max, entered as:
|
|
|
<br>
|
|
|
<br>r and d V / d (logr)
|
|
@@ -421,14 +416,14 @@ refractive index.</td>
|
|
|
<h3>D. Aerosol concentration model (visibility)</h3>
|
|
|
|
|
|
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).
|
|
|
+v, enter directly the value of v [km] (the aerosol optical depth (AOD)
|
|
|
+will be computed from a standard aerosol profile).
|
|
|
<p>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:<br>
|
|
|
<div class="code"><pre>
|
|
|
0 - visibility
|
|
|
-0.112 - aerosol optical depth 550 nm
|
|
|
+0.112 - aerosol optical depth at 550 nm
|
|
|
</pre></div>
|
|
|
|
|
|
<p>NOTE: if iaer is 0, enter -1 for visibility.
|
|
@@ -439,9 +434,9 @@ Target altitude (xps, in negative [km]):
|
|
|
<blockquote>xps >= 0 means the target is at the sea level.
|
|
|
<br>otherwise xps expresses the altitude of the target (e.g., mean elevation)
|
|
|
in [km], given as negative value
|
|
|
-</blockquote>
|
|
|
+</blockquote><br>
|
|
|
|
|
|
-<p>Sensor platform (xpp, in negative [km] or -1000):
|
|
|
+Sensor platform (xpp, in negative [km] or -1000):
|
|
|
<blockquote>
|
|
|
<br>xpp = -1000 means that the sensor is on board a satellite.
|
|
|
<br>xpp = 0 means that the sensor is at the ground level.
|
|
@@ -456,7 +451,7 @@ puw,po3 (water vapor content,ozone content between the aircraft and the surface)
|
|
|
surface)
|
|
|
<p>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 2km
|
|
|
+to the values at the ground level; taerp will be computed according to a 2 km
|
|
|
exponential profile for aerosol.
|
|
|
</blockquote>
|
|
|
|
|
@@ -490,13 +485,13 @@ but step by step output will be printed.</td>
|
|
|
<tr>
|
|
|
<td>0</td>
|
|
|
<td>Enter wlinf, wlsup.
|
|
|
-<br>The filter function will be equal to 1over the whole band.</td>
|
|
|
+<br>The filter function will be equal to 1 over the whole band.</td>
|
|
|
</tr>
|
|
|
|
|
|
<tr>
|
|
|
<td>1</td>
|
|
|
-<td>Enter wlinf, wlsup and user's filter function s(lambda) by step of 0.0025
|
|
|
-micrometer.</td>
|
|
|
+<td>Enter wlinf, wlsup and user's filter function s (lambda) by step of
|
|
|
+0.0025 micrometer.</td>
|
|
|
</tr>
|
|
|
</table>
|
|
|
|
|
@@ -617,7 +612,7 @@ micrometer.</td>
|
|
|
<tr><td>88</td><td><b>RapidEye</b> Blue band (440nm - 512nm)</td></tr>
|
|
|
<tr><td>89</td><td>RapidEye Green band (515nm - 592nm)</td></tr>
|
|
|
<tr><td>90</td><td>RapidEye Red band (628nm - 687nm)</td></tr>
|
|
|
-<tr><td>91</td><td>RapidEye RedEdge band (685nm - 735nm)</td></tr>
|
|
|
+<tr><td>91</td><td>RapidEye Red edge band (685nm - 735nm)</td></tr>
|
|
|
<tr><td>92</td><td>RapidEye NIR band (750nm - 860nm)</td></tr>
|
|
|
|
|
|
<tr><td>93</td><td><b>VGT1 (SPOT4)</b> band 0 (420nm - 497nm)</td></tr>
|
|
@@ -636,7 +631,7 @@ micrometer.</td>
|
|
|
<tr><td>104</td><td>WorldView 2 Green band (503nm - 587nm)</td></tr>
|
|
|
<tr><td>105</td><td>WorldView 2 Yellow band (583nm - 632nm)</td></tr>
|
|
|
<tr><td>106</td><td>WorldView 2 Red band (623nm - 695nm)</td></tr>
|
|
|
-<tr><td>107</td><td>WorldView 2 Red Edge band (698nm - 750nm)</td></tr>
|
|
|
+<tr><td>107</td><td>WorldView 2 Red edge band (698nm - 750nm)</td></tr>
|
|
|
<tr><td>108</td><td>WorldView 2 NIR1 band (760nm - 905nm)</td></tr>
|
|
|
<tr><td>109</td><td>WorldView 2 NIR2 band (853nm - 1047nm)</td></tr>
|
|
|
|
|
@@ -646,21 +641,21 @@ micrometer.</td>
|
|
|
<tr><td>113</td><td>QuickBird Red band (560nm - 747nm)</td></tr>
|
|
|
<tr><td>114</td><td>QuickBird NIR1 band (650nm - 935nm)</td></tr>
|
|
|
|
|
|
-<tr><td>115</td><td><b>Landsat 8 </b> Coastal Aerosol Band (433nm - 455nm)</td></tr>
|
|
|
-<tr><td>116</td><td>Landsat 8 Blue Band (448nm - 515nm)</td></tr>
|
|
|
-<tr><td>117</td><td>Landsat 8 Green Band (525nm - 595nm)</td></tr>
|
|
|
-<tr><td>118</td><td>Landsat 8 Red Band (633nm - 677nm)</td></tr>
|
|
|
-<tr><td>119</td><td>Landsat 8 Panchromatic Band (498nm - 682nm)</td></tr>
|
|
|
-<tr><td>120</td><td>Landsat 8 NIR Band (845nm - 885nm)</td></tr>
|
|
|
-<tr><td>121</td><td>Landsat 8 Cirrus Band (1355nm - 1390nm)</td></tr>
|
|
|
-<tr><td>122</td><td>Landsat 8 SWIR1 Band (1540nm - 1672nm)</td></tr>
|
|
|
-<tr><td>123</td><td>Landsat 8 SWIR2 Band (2073nm - 2322nm)</td></tr>
|
|
|
+<tr><td>115</td><td><b>Landsat 8 </b> Coastal aerosol band (433nm - 455nm)</td></tr>
|
|
|
+<tr><td>116</td><td>Landsat 8 Blue band (448nm - 515nm)</td></tr>
|
|
|
+<tr><td>117</td><td>Landsat 8 Green band (525nm - 595nm)</td></tr>
|
|
|
+<tr><td>118</td><td>Landsat 8 Red band (633nm - 677nm)</td></tr>
|
|
|
+<tr><td>119</td><td>Landsat 8 Panchromatic band (498nm - 682nm)</td></tr>
|
|
|
+<tr><td>120</td><td>Landsat 8 NIR band (845nm - 885nm)</td></tr>
|
|
|
+<tr><td>121</td><td>Landsat 8 Cirrus band (1355nm - 1390nm)</td></tr>
|
|
|
+<tr><td>122</td><td>Landsat 8 SWIR1 band (1540nm - 1672nm)</td></tr>
|
|
|
+<tr><td>123</td><td>Landsat 8 SWIR2 band (2073nm - 2322nm)</td></tr>
|
|
|
|
|
|
<tr><td>115</td><td><b>GeoEye 1</b> Panchromatic band (448nm - 812nm)</td></tr>
|
|
|
-<tr><td>116</td><td>GeoEye 1 Blue Band (443nm - 525nm)</td></tr>
|
|
|
-<tr><td>117</td><td>GeoEye 1 Green Band (503nm - 587nm)</td></tr>
|
|
|
-<tr><td>118</td><td>GeoEye 1 Red Band (653nm - 697nm)</td></tr>
|
|
|
-<tr><td>120</td><td>GeoEye 1 NIR Band (770nm - 932nm)</td></tr>
|
|
|
+<tr><td>116</td><td>GeoEye 1 Blue band (443nm - 525nm)</td></tr>
|
|
|
+<tr><td>117</td><td>GeoEye 1 Green band (503nm - 587nm)</td></tr>
|
|
|
+<tr><td>118</td><td>GeoEye 1 Red band (653nm - 697nm)</td></tr>
|
|
|
+<tr><td>120</td><td>GeoEye 1 NIR band (770nm - 932nm)</td></tr>
|
|
|
|
|
|
<tr><td>129</td><td><b>Spot6</b> Blue band (440nm - 532nm)</td></tr>
|
|
|
<tr><td>130</td><td>Spot6 Green band (515nm - 600nm)</td></tr>
|
|
@@ -738,32 +733,189 @@ micrometer.</td>
|
|
|
|
|
|
<h2>EXAMPLES</h2>
|
|
|
|
|
|
-<h3>Atmospheric correction of a LANDSAT-7 channel</h3>
|
|
|
+<h3>Atmospheric correction of a Sentinel-2 band</h3>
|
|
|
+<p>This example illustrates how to perform atmospheric correction of a
|
|
|
+Sentinel-2 scene in the North Carolina location.
|
|
|
+
|
|
|
+<p>Let's assume that the Sentinel-2 L1C scene
|
|
|
+<em>S2A_OPER_PRD_MSIL1C_PDMC_20161029T092602_R054_V20161028T155402_20161028T155402</em>
|
|
|
+was downloaded and imported with region cropping
|
|
|
+(see <a href="r.import.html">r.import</a>)
|
|
|
+into the <em>PERMANENT</em> mapset of the North Carolina location. The
|
|
|
+computational region was set to the extent of the <em>elevation</em>
|
|
|
+map in the North Carolina dataset. Now, we have 13 individual bands
|
|
|
+(<em>B01-B12</em>) that we want to apply the atmospheric correction to.
|
|
|
+The following steps are applied to each band separately.
|
|
|
+
|
|
|
+<p><b>Create the parameters file for i.atcorr</b>
|
|
|
+<p>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:
|
|
|
+<ul>
|
|
|
+ <li> geometrical conditions,
|
|
|
+ <li> moth, day, decimal hours in GMT, decimal longitude and latitude of measurement,
|
|
|
+ <li> atmospheric model,
|
|
|
+ <li> aerosol model,
|
|
|
+ <li> visibility or aerosol optical depth,
|
|
|
+ <li> mean target elevation above sea level,
|
|
|
+ <li> sensor height and,
|
|
|
+ <li> sensor band.
|
|
|
+</ul>
|
|
|
+
|
|
|
+<ol>
|
|
|
+<li><em>Geometrical conditions</em>
|
|
|
+<p>For Sentinel-2A, the geometrical conditions take the value <tt>25</tt> and for
|
|
|
+Sentinel-2B, the geometrical conditions value is <tt>26</tt> (See table A).
|
|
|
+Our scene comes from the Sentinel-2A mission (the file name begins with
|
|
|
+S2A_...).
|
|
|
+<br><br>
|
|
|
+<li><em>Day, time, longitude and latitude of measurement</em>
|
|
|
+<p>Day and time of the measurement are hidden in the filename (i.e., the
|
|
|
+second datum in the file name with format <tt>YYYYMMDDTHHMMSS</tt>),
|
|
|
+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.
|
|
|
+
|
|
|
+<p>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:
|
|
|
+
|
|
|
+<div class="code"><pre>
|
|
|
+g.region -bg
|
|
|
+</pre></div>
|
|
|
+
|
|
|
+<p>The longitude and latitude of the centre are stored in <em>ll_clon</em>
|
|
|
+and <em>ll_clat</em>. In our case, <tt>ll_clon=-78.691</tt> and
|
|
|
+<tt>ll_clat=35.749</tt>.
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Atmospheric model</em>
|
|
|
+<p>We can choose between various atmospheric models as defined at the
|
|
|
+beginning of this manual. For North Carolina, we can choose <tt>2 -
|
|
|
+midlatitude summer</tt>.
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Aerosol model</em>
|
|
|
+<p>We can also choose between various aerosol models as defined at the
|
|
|
+beginning of this manual. For North Carolina, we can choose <tt>1 -
|
|
|
+continental model</tt>.
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Visibility or Aerosol Optical Depth</em>
|
|
|
+<p>For Sentinel-2 scenes, the visibility is not measured, and therefore
|
|
|
+we have to estimate the aerosol optical depth instead, e.g. from
|
|
|
+<a href="https://aeronet.gsfc.nasa.gov">AERONET</a>. 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 <em>EPA-Res_Triangle_Pk</em> station
|
|
|
+measured AOD = 0.07 (approximately).
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Mean target elevation above sea level</em>
|
|
|
+<p>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:
|
|
|
+
|
|
|
+<div class="code"><pre>
|
|
|
+r.univar -g elevation
|
|
|
+</pre></div>
|
|
|
+
|
|
|
+<p>The mean elevation is stored in <em>mean</em>. In our case,
|
|
|
+<tt>mean=110</tt>. In the 6S file it will be displayed in [-km],
|
|
|
+i.e., <tt>-0.110</tt>.
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Sensor height</em>
|
|
|
+<p>Since the sensor is on board a satellite, the sensor height will be
|
|
|
+set to <tt>-1000</tt>.
|
|
|
+
|
|
|
+<br><br>
|
|
|
+<li><em>Sensor band</em>
|
|
|
+<p>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.
|
|
|
+</ol>
|
|
|
+
|
|
|
+<p>Finally, here is what the 6S file would look like for Band 02 of our
|
|
|
+scene. In order to use it in the <em>i.atcorr</em> module, we can save
|
|
|
+it in a text file, for example <tt>params_B02.txt</tt>.
|
|
|
+<div class="code"><pre>
|
|
|
+25
|
|
|
+10 28 15.901 -78.691 35.749
|
|
|
+2
|
|
|
+1
|
|
|
+0
|
|
|
+0.07
|
|
|
+-0.110
|
|
|
+-1000
|
|
|
+167
|
|
|
+</pre></div>
|
|
|
+
|
|
|
+<p><b>Compute atmospheric correction</b>
|
|
|
+<p>In the next step we run <em>i.atcorr</em> for the selected band
|
|
|
+<em>B02</em> of our Sentinel 2 scene. We have to specify the following
|
|
|
+parameters:
|
|
|
+<ul>
|
|
|
+ <li><b>input</b> = raster band to be processed,
|
|
|
+ <li><b>parameters</b> = path to 6S file created in the previous step (we could also enter the values directly),
|
|
|
+ <li><b>output</b> = name for the output corrected raster band,
|
|
|
+ <li><b>range</b> = from 1 to the <tt>QUANTIFICATION_VALUE</tt> stored in the metadata file. It is <tt>10000</tt> for both Sentinel-2A and Sentinel-2B.
|
|
|
+ <li><b>rescale</b> = 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.
|
|
|
+</ul>
|
|
|
+<p>If the data is available, the following parameters can be specified
|
|
|
+as well:
|
|
|
+<ul>
|
|
|
+ <li><b>elevation</b> = raster of digital elevation model,
|
|
|
+ <li><b>visibility</b> = raster of visibility model.
|
|
|
+</ul>
|
|
|
+
|
|
|
+<p>Finally, this is how the command would look like to apply atmospheric
|
|
|
+correction to band <em>B02</em>:
|
|
|
+
|
|
|
+<div class="code"><pre>
|
|
|
+i.atcorr input=B02 parameters=params_B02.txt output=B02.atcorr range=1,10000 rescale=0,255 elevation=elevation
|
|
|
+</pre></div>
|
|
|
+
|
|
|
+<p>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.
|
|
|
+
|
|
|
+<div align="center" style="margin: 10px">
|
|
|
+<a href="i_atcorr_B02_atcorr.png">
|
|
|
+<img src="i_atcorr_B02_atcorr.png" width="600" height="486" alt="i.atcorr example" border="0">
|
|
|
+</a><br>
|
|
|
+<i>Figure: Sentinel-2A Band 02 with applied atmospheric correction (histogram equalization grayscale color scheme)</i>
|
|
|
+</div>
|
|
|
+
|
|
|
+<h3>Atmospheric correction of a Landsat-7 band</h3>
|
|
|
+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:
|
|
|
|
|
|
-The example is based on the North Carolina sample dataset (GMT -5 hours).
|
|
|
-First we set the computational region to the satellite map, e.g. channel 4:
|
|
|
<div class="code"><pre>
|
|
|
g.region raster=lsat7_2002_40 -p
|
|
|
</pre></div>
|
|
|
|
|
|
-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
|
|
|
+<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
|
|
|
<a href="ftp://ftp.glcf.umd.edu/glcf/Landsat/WRS2/p016/r035/p016r035_7x20020524.ETM-EarthSat-Orthorectified/p016r035_7x20020524.met">metadata</a>
|
|
|
file (<a href="http://www.grassbook.org/wp-content/uploads/ncexternal/landsat/2002/p016r035_7x20020524.met">file copy</a>; North Carolina
|
|
|
-sample dataset). In case of the North Carolina sample dataset, values
|
|
|
-have been stored for each channel and can be retrieved like this:
|
|
|
+sample dataset). In the case of the North Carolina sample dataset, these
|
|
|
+values have been stored for each channel and can be retrieved with:
|
|
|
|
|
|
<div class="code"><pre>
|
|
|
r.info lsat7_2002_40
|
|
|
</pre></div>
|
|
|
|
|
|
In this case, we have: SUN_AZIMUTH = 120.8810347, SUN_ELEVATION = 64.7730999.
|
|
|
-<p>
|
|
|
-If the sun position metadata are unavailable, we can also calculate
|
|
|
+<p>If the sun position metadata are unavailable, we can also calculate
|
|
|
them from the overpass time as follows
|
|
|
(<em><a href="r.sunmask.html">r.sunmask</a></em>
|
|
|
uses <a href="http://www.nrel.gov/midc/solpos/solpos.html">SOLPOS</a>):
|
|
|
+
|
|
|
<div class="code"><pre>
|
|
|
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
|
|
@@ -772,20 +924,22 @@ r.sunmask -s elev=elevation out=dummy year=2002 month=5 day=24 hour=10 min=42 se
|
|
|
If the overpass time is unknown, use the
|
|
|
<a href="http://cloudsgate2.larc.nasa.gov/cgi-bin/predict/predict.cgi">NASA LaRC Satellite Overpass Predictor</a>.
|
|
|
|
|
|
-<h4>Conversion of digital number (DN) to radiance at top-of-atmosphere (TOA)</h4>
|
|
|
+<h4>Convert digital numbers (DN) to radiance at top-of-atmosphere (TOA)</h4>
|
|
|
|
|
|
For Landsat and ASTER, the conversion can be conveniently done with
|
|
|
-<a href="i.landsat.toar.html">i.landsat.toar</a> or <a href="i.aster.toar.html">i.aster.toar</a>,
|
|
|
-respectively.
|
|
|
+<em><a href="i.landsat.toar.html">i.landsat.toar</a></em> or
|
|
|
+<em><a href="i.aster.toar.html">i.aster.toar</a></em>, respectively.
|
|
|
+
|
|
|
+<p>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:
|
|
|
|
|
|
-<p>
|
|
|
-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
|
|
|
<div class="code"><pre>
|
|
|
# formula depends on satellite sensor, see respective metadata
|
|
|
Lλ = ((LMAXλ - LMINλ)/(QCALMAX-QCALMIN)) * (QCAL-QCALMIN) + LMINλ
|
|
|
</pre></div>
|
|
|
-where:
|
|
|
+
|
|
|
+where,
|
|
|
<ul>
|
|
|
<li> Lλ = Spectral Radiance at the sensor's aperture in Watt/(meter squared * ster * µm), the
|
|
|
apparent radiance as seen by the satellite sensor;</li>
|
|
@@ -796,15 +950,16 @@ where:
|
|
|
<li> QCALMAX = the maximum quantized calibrated pixel value (corresponding to LMAXλ) in DN=255.</li>
|
|
|
</ul>
|
|
|
|
|
|
-LMINλ and LMAXλ are the radiances related to the minimal and
|
|
|
-maximal DN value, and are reported in the metadata file for each image, or in
|
|
|
-the table 1. High gain or low gain is also reported in the metadata file of each
|
|
|
-satellite image. For Landsat, the minimal DN value (QCALMIN) is 1 for Landsat ETM+
|
|
|
-images (see
|
|
|
-<a href="http://landsathandbook.gsfc.nasa.gov/pdfs/Landsat7_Handbook.pdf">Landsat handbook</a>, see chapter 11),
|
|
|
+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 <a href="http://landsathandbook.gsfc.nasa.gov/pdfs/Landsat7_Handbook.pdf">Landsat handbook</a>, chapter 11),
|
|
|
and the maximal DN value (QCALMAX) is 255. QCAL is the DN value for every
|
|
|
separate pixel in the Landsat image.
|
|
|
-<p>We extract the coefficients and apply them in order to obtain the radiance map:
|
|
|
+<p>We extract the coefficients and apply them in order to obtain the
|
|
|
+radiance map:
|
|
|
+
|
|
|
<div class="code"><pre>
|
|
|
CHAN=4
|
|
|
r.info lsat7_2002_${CHAN}0 -h | tr '\n' ' ' | sed 's+ ++g' | tr ':' '\n' | grep "LMIN_BAND${CHAN}\|LMAX_BAND${CHAN}"
|
|
@@ -814,36 +969,33 @@ QCALMAX_BAND4=255.0,p016r035_7x20020524.met
|
|
|
QCALMIN_BAND4=1.0,p016r035_7x20020524.met
|
|
|
</pre></div>
|
|
|
|
|
|
-Conversion to radiance (this calculation is done for band 4, for the other bands, the numbers in italics
|
|
|
-need to be replaced with their related values):
|
|
|
+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):
|
|
|
|
|
|
<div class="code"><pre>
|
|
|
r.mapcalc "lsat7_2002_40_rad = ((241.1 - (-5.1)) / (255.0 - 1.0)) * (lsat7_2002_40 - 1.0) + (-5.1)"
|
|
|
</pre></div>
|
|
|
|
|
|
-Again, the <em>r.mapcalc</em> calculation is only needed when working with satellite data
|
|
|
-other than Landsat or ASTER.
|
|
|
+Again, the <em>r.mapcalc</em> calculation is only needed when working
|
|
|
+with satellite data other than Landsat or ASTER.
|
|
|
|
|
|
-<h4>Creation of parameter file for i.atcorr</h4>
|
|
|
+<h4>Create the parameters file for i.atcorr</h4>
|
|
|
|
|
|
-The underlying 6S model is parametrized through a control file, indicated with the
|
|
|
-<em>parameter</em> option. This is a text file defining geometrical and atmospherical
|
|
|
-conditions of the satellite overpass. Below some details:
|
|
|
+The underlying 6S model is parametrized through a control file,
|
|
|
+indicated with the <b>parameters</b> option. This is a text file
|
|
|
+defining geometrical and atmospherical conditions of the satellite
|
|
|
+overpass.
|
|
|
|
|
|
-<p>
|
|
|
-<div class="code"><pre>
|
|
|
-# find mean elevation (target above sea level, used as initialization value in control file)
|
|
|
-r.univar elevation
|
|
|
-</pre></div>
|
|
|
-
|
|
|
-Create a control file 'icnd.txt' for channel 4 (NIR), based on metadata. For the overpass time,
|
|
|
-we need to define decimal hours:<br>
|
|
|
-10:42:07 NC local time = 10.70 decimal hours (decimal minutes: 42 * 100 / 60) which is 15.70 GMT:
|
|
|
+Here we create a control file <tt>icnd_lsat4.txt</tt> 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.
|
|
|
|
|
|
<div class="code"><pre>
|
|
|
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 mode=midlatitude summer
|
|
|
+2 - atmospheric model=midlatitude summer
|
|
|
1 - aerosols model=continental
|
|
|
50 - visibility [km] (aerosol model concentration)
|
|
|
-0.110 - mean target elevation above sea level [km]
|
|
@@ -851,27 +1003,30 @@ we need to define decimal hours:<br>
|
|
|
64 - 4th band of ETM+ Landsat 7
|
|
|
</pre></div>
|
|
|
|
|
|
-Finally, run the atmospheric correction (-r for reflectance input map; -a for date >July 2000):
|
|
|
+Finally, run the atmospheric correction (-r for reflectance input map;
|
|
|
+-a for date > July 2000):
|
|
|
+
|
|
|
<div class="code"><pre>
|
|
|
-i.atcorr -r -a lsat7_2002_40_rad elev=elevation parameters=icnd_lsat4.txt output=lsat7_2002_40_atcorr
|
|
|
+i.atcorr -r -a lsat7_2002_40_rad elevation=elevation parameters=icnd_lsat4.txt output=lsat7_2002_40_atcorr
|
|
|
</pre></div>
|
|
|
|
|
|
-Note that the altitude value from 'icnd_lsat4.txt' file is read at the beginning
|
|
|
-to compute the initial transform. It is necessary to give a value which could
|
|
|
-be the mean value of the elevation model. For the atmospheric correction then
|
|
|
-the raster elevation values are used from the map.
|
|
|
-<p>Note that the process is computationally intensive.<br>
|
|
|
-Note also, that <em>i.atcorr</em> reports solar elevation angle above horizon rather than solar zenith angle.
|
|
|
+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
|
|
|
+(<tt>r.univar elevation</tt>). For the atmospheric correction per se, the
|
|
|
+elevation values from the raster map are used.
|
|
|
+<p>Note that the process is computationally intensive. Note also, that
|
|
|
+<em>i.atcorr</em> reports solar elevation angle above horizon rather
|
|
|
+than solar zenith angle.
|
|
|
|
|
|
<h2>REMAINING DOCUMENTATION ISSUES</h2>
|
|
|
-1. The influence and importance of the visibility value or map should be
|
|
|
+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.
|
|
|
|
|
|
<h2>SEE ALSO</h2>
|
|
|
|
|
|
-GRASS Wiki page about
|
|
|
- <a href="http://grasswiki.osgeo.org/wiki/Atmospheric_correction">Atmospheric correction</a>
|
|
|
+GRASS Wiki page about <a href="http://grasswiki.osgeo.org/wiki/Atmospheric_correction">Atmospheric correction</a>
|
|
|
<p>
|
|
|
<em>
|
|
|
<a href="i.aster.toar.html">i.aster.toar</a>,
|
|
@@ -884,7 +1039,7 @@ GRASS Wiki page about
|
|
|
<h2>REFERENCES</h2>
|
|
|
|
|
|
<ul>
|
|
|
-<li> Vermote, E.F., Tanre, D., Deuze, J.L., Herman, M., and Morcrette, J.J., 1997,
|
|
|
+<li>Vermote, E.F., Tanre, D., Deuze, J.L., Herman, M., and Morcrette, J.J., 1997,
|
|
|
Second simulation of the satellite signal in the solar spectrum, 6S: An
|
|
|
overview., IEEE Trans. Geosc. and Remote Sens. 35(3):675-686.
|
|
|
<!-- too new:
|
|
@@ -893,14 +1048,15 @@ overview., IEEE Trans. Geosc. and Remote Sens. 35(3):675-686.
|
|
|
at the <a href="http://6s.ltdri.org">6S homepage</a> of the Land Surface Reflectance
|
|
|
Science Computing Facility
|
|
|
-->
|
|
|
-<li> 6S Manual: <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P1.pdf">PDF1</a>,
|
|
|
- <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P2.pdf">PDF2</a>,
|
|
|
- and <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P3.pdf">PDF3</a>
|
|
|
+<li>6S Manual: <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P1.pdf">PDF1</a>,
|
|
|
+ <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P2.pdf">PDF2</a>,
|
|
|
+ and <a href="http://www.rsgis.ait.ac.th/~honda/textbooks/advrs/6smanv2.0_P3.pdf">PDF3</a>
|
|
|
<li>RapidEye sensors have been provided by <a href="http://www.rapideye.com/">RapidEye AG, Germany</a>
|
|
|
-<li> Julia A. Barsi, Brian L. Markham and Jeffrey A. Pedelty "The operational land imager: spectral response and spectral uniformity", Proc. SPIE 8153, 81530G (2011); doi:10.1117/12.895438
|
|
|
+<li>Barsi, J.A., Markham, B.L. and Pedelty, J.A., 2011, The operational
|
|
|
+land imager: spectral response and spectral uniformity., Proc. SPIE 8153,
|
|
|
+81530G; doi:10.1117/12.895438
|
|
|
</ul>
|
|
|
|
|
|
-
|
|
|
<h2>AUTHORS</h2>
|
|
|
|
|
|
<p><em>Original version of the program for GRASS 5:</em>
|