[GRASS-user] Negative raster values from i.landsat.dehaze

Wed May 2 17:10:39 EDT 2012

Hi,

I have been trying to get some results from GRASS for a couple of weeks now but am hitting a wall.  What I am trying to do is to take an area of land from a series of Landsat 7 ETM+ images and calculate NDVI over the area for the purpose of comparison over a time-period.  This obviously requires atmospheric correction, so the command sequence I am using after importing the images for a given mapset is:

# Generate radiance values
i.landsat.toar input_prefix="B" output_prefix="R" metfile="./L71174038_03820080324/L71174038_03820080324_MTL.txt" sensor="tm7"

# Atmospheric correction
i.tasscap -7 band1=R1 band2=R2 band3=R3 band4=R4 band5=R5 band7=R7 outprefix=TC
i.landsat.dehaze band1=R1 band2=R2 band3=R3 band4=R4 band5=R5 band7=R7 tasscap4=TC.4 outprefix=D

# Generate NDVI map
r.mapcalc 'ndvi=if((D.3 == 0 && D.4 == 0), -1, (D.4 - D.3) / (D.4 + D.3))'

Then I mask the NDVI map with an area raster which I already created.  There is probably a better way of doing this.

# Mask NDVI map with individual regions
r.mapcalc 'area_ndvi = ndvi * area_mask at PERMANENT'

Finally I want an average value of the NDVI over the area.  There is definitely a better way of doing this but I couldn't find a grass command that spits out a single value.  Instead I used r.stats to give the area of each group of values and then wrote a Perl script to average these.

# Process statistics
r.stats -a --q area_ndvi | ./stats.pl

Ok, so here come the questions:

1) Is there a better way to simply get an average of all values over a given area to replace my stats.pl script?

2) and this is the big one.  When I perform i.landsat.dehaze I get rasters with primarily negative values.  This renders the NDVI calculation meaningless as they are ratios of +ve and -ve values.  I am presuming this is a mistake as negative reflectance has little meaning anyway - did I use the wrong sequence of commands?  Below are info outputs for the red-band after the TOAR and dehaze operations respectively:

Many thanks for any help on this,
Jamie

 +----------------------------------------------------------------------------+
 | Layer:    R3                             Date: Wed May  2 20:15:30 2012    |
 | Mapset:   L71174038_03820080324          Login of Creator: jamie           |
 | Location: area                                                             |
 | DataBase: /home/jamie/grassdata                                            |
 | Title:     ( R3 )                                                          |
 | Timestamp: none                                                            |
 |----------------------------------------------------------------------------|
 |                                                                            |
 |   Type of Map:  raster               Number of Categories: 255             |
 |   Data Type:    DCELL                                                      |
 |   Rows:         6961                                                       |
 |   Columns:      7911                                                       |
 |   Total Cells:  55068471                                                   |
 |        Projection: UTM (zone 36)                                           |
 |            N:    3619200    S:    3410370   Res:    30                     |
 |            E:     843630    W:     606300   Res:    30                     |
 |   Range of data:    min = -0.0129146013370688  max = 0.605436510681785     |
 |                                                                            |
 |   Data Description:                                                        |
 |    generated by i.landsat.toar                                             |
 |                                                                            |
 |   Comments:                                                                |
 |     Reflectance of Landsat-7 ETM+ (method uncorrected)                     |
 |    ----------------------------------------------------------------        |
 |     Acquisition date ...................... 2008-03-24                     |
 |     Production date ....................... 2012-02-23                     |
 |                                                                            |
 |     Earth-sun distance (d) ................ 0.99699779                     |
 |     Digital number (DN) range ............. 1 to 255                       |
 |     Calibration constants (Lmin to Lmax) .. -5.000 to +234.400             |
 |     DN to Radiance (gain and bias) ........ +0.94252 and -5.94252          |
 |     Mean solar irradiance (ESUN) .......... 1551.000                       |
 |     Reflectance = Radiance divided by ..... 387.15868                      |
 |    -----------------------------------------------------------------       |
 |                                                                            |
 |    i.landsat.toar input_prefix="B" output_prefix="R" metfile="./L711740\   |
 |    38_03820080324/L71174038_03820080324_MTL.txt" sensor="tm7" method="u\   |
 |    ncorrected" percent=0.01 pixel=1000 sat_zenith=8.2000 rayleigh=0.0      |
 |                                                                            |
 +----------------------------------------------------------------------------+
 +----------------------------------------------------------------------------+
 | Layer:    D.3                            Date: Wed May  2 20:22:21 2012    |
 | Mapset:   L71174038_03820080324          Login of Creator: jamie           |
 | Location: area                                                             |
 | DataBase: /home/jamie/grassdata                                            |
 | Title:     ( D.3 )                                                         |
 | Timestamp: none                                                            |
 |----------------------------------------------------------------------------|
 |                                                                            |
 |   Type of Map:  raster               Number of Categories: 255             |
 |   Data Type:    DCELL                                                      |
 |   Rows:         334                                                        |
 |   Columns:      166                                                        |
 |   Total Cells:  55444                                                      |
 |        Projection: UTM (zone 36)                                           |
 |            N:    3493000    S:    3483000   Res: 29.94011976               |
 |            E:     649000    W:     644000   Res: 30.12048193               |
 |   Range of data:    min = -0.463376848126162  max = 0.19725822381175       |
 |                                                                            |
 |   Data Description:                                                        |
 |    generated by r.mapcalc                                                  |
 |                                                                            |
 |   Comments:                                                                |
 |    R3 - (TC.4 - 0.134975) * -3.08586                                       |
 |                                                                            |
 +----------------------------------------------------------------------------+
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/grass-user/attachments/20120502/a6f98f24/attachment-0001.html