Import and atmospherically/topographically correct North Carolina Dataset

import_nc.sh

# Open the North Carolina dataset and create a new mapset.
# You will also need to download the *.met file from here:
#    http://www.grassbook.org/wp-content/uploads/ncexternal/landsat/2002/p016r035_7x20020524.met.gz

# toar variables
lsat_prefix=lsat7_2002_
toar_prefix=R
metadata_path=/home/grassuser/repo/lsat7_2000.met   # download it from here:
toar_method=dos3

# topo.corr variables
# get numeric values from *.met file
# by the way "sun zenith" is derived from the following equation:
#       $sun_zenith = 90 - $sun_elevation
elevation='elevation@PERMANENT'
resampled_elevation=resampled_elevation
illumination=illumination
topo_method=c-factor
topo_prefix=T
resampling_method=average
resampling_method=mode
sun_azimuth=120.8810347
sun_elevation=64.7730999
sun_zenith=25.2269001

# setup region
g.region raster=${lsat_prefix}10 align=${lsat_prefix}10 -p

# copy landsat7 images from PERMANENT to the current mapset.
# We also rename the files so that i.landsat.toar can run without problems
g.copy raster=lsat7_2002_10@PERMANENT,lsat7_2002_1
g.copy raster=lsat7_2002_20@PERMANENT,lsat7_2002_2
g.copy raster=lsat7_2002_30@PERMANENT,lsat7_2002_3
g.copy raster=lsat7_2002_40@PERMANENT,lsat7_2002_4
g.copy raster=lsat7_2002_50@PERMANENT,lsat7_2002_5
g.copy raster=lsat7_2002_61@PERMANENT,lsat7_2002_61
g.copy raster=lsat7_2002_62@PERMANENT,lsat7_2002_62
g.copy raster=lsat7_2002_70@PERMANENT,lsat7_2002_7
g.copy raster=lsat7_2002_80@PERMANENT,lsat7_2002_8

# Calculate TOA reflectance.
i.landsat.toar \
    --verbose \
    --overwrite \
    input=$lsat_prefix \
    output=$toar_prefix \
    metfile=$metadata_path \
    method=$toar_method

# remove the DNs from this mapset
g.remove -f raster pattern=lsat7_2002_*

# It makes sense to follow the same naming scheme as the original data.
g.rename raster=R1,R10
g.rename raster=R2,R20
g.rename raster=R3,R30
g.rename raster=R4,R40
g.rename raster=R5,R50
g.rename raster=R6,R61
g.rename raster=R6,R62
g.rename raster=R7,R70
g.rename raster=R8,R80

# We can now move on to the topographic correction.
# We must first make sure that DEM has the same resolution as the satellite imagery.
# We can check the respective resolutions using r.info |
#       r.info -gre <map_name> | grep "..res"
r.info -gre R10 | grep "..res"
r.info -gre elevation | grep "..res"
# In our case the respective resolutions are equal to:
# * R10 -> 28.5m
# * elevation -> 10
# Since we are resampling to a coarser resolution we are going to use r.resamp.stats
# If we wanted to resample to a finer resolution we would have to choose a different
# module. For more info check https://grasswiki.osgeo.org/wiki/Interpolation#Resampling_of_raster_maps_to_coarser_resolution
r.resamp.stats \
    -w \
    --verbose \
    --overwrite \
    input=$elevation \
    output=$resampled_elevation \
    method=$resampling_method

# Create illumination model
i.topo.corr \
    -i \
    --verbose \
    --overwrite \
    basemap=$resampled_elevation \
    zenith=$sun_zenith \
    azimuth=$sun_azimuth \
    output=$illumination \
    method=$topo_method

i.topo.corr \
    -s \
    --verbose \
    --overwrite \
   in=R10,R20,R30,R40,R50,R61,R62,R70,R80 \
   out=$topo_prefix \
   zenith=$sun_zenith \
   basemap=$illumination

univars.sh

GRASS 7.0.1 (nc_spm_08_grass7):~/repo > r.univar illumination -ge
n=247275
null_cells=3050
cells=250325
min=0.717182085952387
max=0.996726780887321
range=0.279544694934934
mean=0.903771220913367
mean_of_abs=0.903771220913367
stddev=0.0214276127990178
variance=0.000459142590264632
coeff_var=2.37091116680643
sum=223480.028651353
first_quartile=0.891549
median=0.905526
third_quartile=0.917649
percentile_90=0.928806

GRASS 7.0.1 (nc_spm_08_grass7):~/repo > r.univar resampled_elevation -ge
n=250325
null_cells=0
cells=250325
min=55.5787925720215
max=156.201751708984
range=100.622959136963
mean=109.938668436965
mean_of_abs=109.938668436965
stddev=20.3289250189747
variance=413.265192427097
coeff_var=18.4911508461927
sum=27520397.1764832
first_quartile=94.3594
median=108.412
third_quartile=126.33
percentile_90=138.285