[GRASS-user] Estimating Albedo from Landsat

Tue Sep 14 07:30:15 EDT 2010

Hi Hamish!

I was about to report things related to what you describe below but did not do 
so when I realised that I don't really have the "correct" metadata file for 
the scene I was given. The "meta" file available at Glovis is not really the 
complete set of meta-information (is my guess -- missing "sat_zenith=" value 
for example).

I ordered the scene again just to get the metadata stuff (since I can't 
contact at the moment the person who provided the Landsat5 TM image) and wait 
for it.

Hamish wrote:

> I've just made some small changes to i.landsat.acca and i.landsat.toar in
> svn which shouldn't change the behaviour of the module at all, just do
> some more error checks...

The "long and never ending process" was due to missing parameters I guess. 
Adding some arbitrary "sat_zenith=" value (actually I used the sun_azimuth 
given in the ".meta" file :-)), along with "product_date=2007-09-05 
date=2007-09-05 solar_elevation=52.48237006", did the trick and completed the 
task within a short time.

FWIW,

the output values (even if the sat_zenith value I used is non-sense) seem to 
be as expected for radiance/reflectance values (e.g. band 1, 
min=-0.00308741863579021, max=0.39202091888652) .

> ...and make better use of libgis (writes out more
> metadata etc. e.g. added units=Kelvin to band 6 metadata, not sure about
> W/m^2 for others or how safe it is to include non-ascii chars in the units
> field? (micron "u", etc.)).

Is the following the problem (you expect-ed)? I see strange characters, e.g.:

-------------------
 BAND 1  (code 1)
   calibrated digital number (DN): 1.0 to 255.0
   calibration constants (L): -1.520 to 193.000
   at-sensor radiance = 0.76583 � DN + -2.28583
   mean solar exoatmospheric irradiance (ESUN): 1983.000
   at-sensor reflectance = radiance / 492.32067
-------------------
 BAND 2  (code 2)
   calibrated digital number (DN): 1.0 to 255.0
   calibration constants (L): -2.840 to 365.000
   at-sensor radiance = 1.44819 � DN + -4.28819
   mean solar exoatmospheric irradiance (ESUN): 1796.000
   at-sensor reflectance = radiance / 445.89406

> there were some troubles if e.g. production date was not set, it would
> continue anyway with a very old date. (it would be nice to have a check
> that production date is not before acq. date as well, but that's not
> in it yet.)
> also if the date format was not right it would continue anyway with a
> bogus value and crazy value for solar distance, etc.
> 
> I've got an old landsat-5 TM scene here from 2004 which seemed to run ok.
> i.landsat.toar didn't recognize the header.dat metadata file, so I had
> to enter some values by hand.
> 
> band 6 temperatures do not look very plausible (30 degC ocean and -6 degC
> land) on a sunny autumn day.

> i.landsat.acca then does a really bad job finding the clouds, but I think
> that's just because the toar.6 band is so badly wrong.

(
Are the "raw" values given in band 6 in Kelvin or are they the result of some 
(re-)scaling to 0-255 grey values?
)

> ISTR doing the
> temperature calc with r.mapcalc from Markus and Helena's book some years
> ago, I don't remember it being so badly wrong then. there's some commented
> out landsat-5 code in i.landsat.acca though, maybe it is as yet unfinished?
> 
> 
> It would be good to compare with i.atcorr, but I haven't figured out the
> exact recipe for that module yet.
> 
> Probably we should concentrate on the North Carolina dataset for common
> tests... its landsat mapset has:
> 
> lsat5_1987_10
> lsat5_1987_20
> lsat5_1987_30
> lsat5_1987_40
> lsat5_1987_50
> lsat5_1987_60
> lsat5_1987_70
> 
> not enough metadata for toar, might have to redownload the scene from
> GloVis:
> IMAGE_ID=P016R35_5T871014,PATH=16,ROW=35,DATE=10/14/87,PLATFORM=LANDSAT5
> 
> and
> lsat7_2000_10
> lsat7_2000_20
> lsat7_2000_30
> lsat7_2000_40
> lsat7_2000_50
> lsat7_2000_61
> lsat7_2000_70
> lsat7_2000_80
> 
> metadata from hist/ file:
> SPACECRAFT_ID=Landsat7,SENSOR_ID=ETM+,ACQUISITION_DATE=2000-03-31,WRS_P
> ATH=16,CPF_FILE_NAME=L7CPF20000101_20000331_12,LMAX_BAND1=191.600,LMIN_
> BAND1=-6.200,LMAX_BAND2=196.500,LMIN_BAND2=-6.400,LMAX_BAND3=152.900,LM
> IN_BAND3=-5.000,LMAX_BAND4=241.100,LMIN_BAND4=-5.100,LMAX_BAND5=31.060,
> LMIN_BAND5=-1.000,LMAX_BAND61=17.040,LMIN_BAND61=0.000,LMAX_BAND62=12.6
> 50,LMIN_BAND62=3.200,LMAX_BAND7=10.800,LMIN_BAND7=-0.350,LMAX_BAND8=243
> .100,LMIN_BAND8=-4.700,QCALMAX_BAND1=255.0,QCALMIN_BAND1=1.0,QCALMAX_BA
> ND2=255.0,QCALMIN_BAND2=1.0,QCALMAX_BAND3=255.0,QCALMIN_BAND3=1.0,QCALM
> AX_BAND4=255.0,QCALMIN_BAND4=1.0,QCALMAX_BAND5=255.0,QCALMIN_BAND5=1.0,
> QCALMAX_BAND61=255.0,QCALMIN_BAND61=1.0,QCALMAX_BAND62=255.0,QCALMIN_BA
> ND62=1.0,QCALMAX_BAND7=255.0,QCALMIN_BAND7=1.0,QCALMAX_BAND8=255.0,QCAL
> MIN_BAND8=1.0,SUN_AZIMUTH=139.6033279,SUN_ELEVATION=51.524652

> should we be calculating the sat_zenith= param from SUN_AZIMUTH in the
> metadata file, or is that something else? I just used the default 8.2 deg.

Good question! As described above, for the sake of testing I just fed the 
"sat_zenith" parameter the "sun_azimuth" value and the module run without 
complaints. My guess is that this parameter is important and will make a 
difference depending on the date and time of the acquisition.

> > > > 1. the "LT51830332007248MOR00.meta" isn't read correctly or
> > > > the module is  looking for a string other than the "date_acquired =
> > > > 20070905" string which is the only "date-string" included in the
> > > > ".meta" file
> 
> the program wants "yyyy-mm-dd". (the test for that could still be improved)

yep, finally fed it manually

> > > > 2. does not recognize the solar elevation value
> 
> for me that worked (at least it was reported as given on the command line
> when the -v show settings flag was used)
> 
> > > > 3. after manually providing the "solar_elevation" and the
> > > > "product_date" and using the "-v" flag, some random (?) segfaults
> > > > occur.
> 
> do you still get those with the new svn version?

In my last test(s) no segfault occured... with or without "-v" !??
>
> > > > 4. without the "-v" flag it runs but seems to take too much
> > > > time? OK, the image set is large ( 7654 cols x 7127 rows x 7 bands).
> 
> were any other values bogus? e.g. earth-sun distance not near 1AU?

no I don't think so. I can try to repeat it but does it worth the time 
spending? I mean, after feeding the module with the "expected" values it runs 
fine.

> probably a memory bug if the time variable was in the millions, etc.
> 
> > Oh man, I guess I am doing something wrong [ again :-( ].
> > It takes too much  time here. I broke the process after >30 min
> > and set up a new (smaller) computational region (of my interest).
> 
> AFAIU at least i.landsat.acca ignores the region settings (why?), and
> I'm not sure if i.landsat.toar does or not, I suspect it ignores it as
> well and works on the full image.

confirmed, the "toar" module ignores the region settings. Working on

rows=4906
cols=3707

gives

Rows:            7127
Columns:      7654 

> > Thanks for all of your work Hamish, Nikos

> really I've just made some minor tests and tidying of the code, the credit
> belongs to others. I'm still learning these modules as much as anybody
> else ..

Well, to me your work (as well as the work of all others) is invaluable. Even 
though "nice words" are inexpensive, I still think they are important (when 
they are honest). From my side, I just need some timing to contribute back - 
all seems to go wrong currently :-p

> ps- the grass7 patches will be badly out of date now, as they are a pain
> to update for each rev, I'll wait at least until the libgis/libraster
> G_fatal_error() and other messages are sync'd with the standardized
> versions before cutting a new one of those. or at least until a little
> dust settles.