[GRASS-SVN] r53209 - grass-addons/grass6/imagery/i.landsat.toar
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Sep 18 09:29:18 PDT 2012
Author: ejtizado
Date: 2012-09-18 09:29:17 -0700 (Tue, 18 Sep 2012)
New Revision: 53209
Modified:
grass-addons/grass6/imagery/i.landsat.toar/description.html
grass-addons/grass6/imagery/i.landsat.toar/landsat_met.c
grass-addons/grass6/imagery/i.landsat.toar/landsat_set.c
grass-addons/grass6/imagery/i.landsat.toar/local_proto.h
grass-addons/grass6/imagery/i.landsat.toar/main.c
Log:
Update to read metadata from .met and MTL.txt for all landsats
Modified: grass-addons/grass6/imagery/i.landsat.toar/description.html
===================================================================
--- grass-addons/grass6/imagery/i.landsat.toar/description.html 2012-09-18 15:09:24 UTC (rev 53208)
+++ grass-addons/grass6/imagery/i.landsat.toar/description.html 2012-09-18 16:29:17 UTC (rev 53209)
@@ -12,10 +12,10 @@
needed the gain (high or low) of the nine respective bands.
<p>
-Optionally, the data can be read from header file (.met) for all
+Optionally (recommended), the data can be read from metadata file (.met or MTL.txt) for all
Landsat MSS, TM and ETM+. However, if the solar elevation or the
-product creation date are given the values of the metfile are
-overwriten. This is necessary when the data in the metfile is
+product creation date are given the values of the metadata file are
+overwriten (only with .met files). This is necessary when the data in the .met file is
incorrect or not accurate.
<p>
@@ -127,12 +127,6 @@
satellite data and the parameters used in the transformations.
<p>
-In L5_MTL mode (flag <b>-t</b>), the Landsat 5TM imagery that has a
-_MTL.txt metadata file can be processed. Landsat 7 ETM+ does not need
-a flag since .met and _MTL.txt are sufficient compatible for this
-sensor.
-
-<p>
Production date is not an exact value but it is necessary to apply
correct calibration constants, which were changed in the dates:
<ul>
@@ -155,23 +149,23 @@
<div class="code"><pre>
i.landsat.toar input_prefix=203_30. output_prefix=_toar \
- sensor=tm7 met=p203r030_7x20010620.met
+ metfile=p203r030_7x20010620.met
</pre></div>
or
<div class="code"><pre>
-i.landsat.toar input_prefix=203_30. output_prefix=_toar \
- sensor=tm7 product_date=2004-06-07 date=2001-06-20 \
- solar_elevation=64.3242970 gain="HHHLHLHHL"
+i.landsat.toar input_prefix=L5121060_06020060714. \
+ output_prefix=L5121060_06020060714_toar \
+ metfile=L5121060_06020060714_MTL.txt
</pre></div>
or
<div class="code"><pre>
-i.landsat.toar input_prefix=L5121060_06020060714. \
- output_prefix=L5121060_06020060714_toar sensor=tm5 \
- metfile=L5121060_06020060714_MTL.txt -t
+i.landsat.toar input_prefix=203_30. output_prefix=_toar \
+ sensor=tm7 product_date=2004-06-07 date=2001-06-20 \
+ sun_elevation=64.3242970 gain="HHHLHLHHL"
</pre></div>
<h2>REFERENCES</h2>
Modified: grass-addons/grass6/imagery/i.landsat.toar/landsat_met.c
===================================================================
--- grass-addons/grass6/imagery/i.landsat.toar/landsat_met.c 2012-09-18 15:09:24 UTC (rev 53208)
+++ grass-addons/grass6/imagery/i.landsat.toar/landsat_met.c 2012-09-18 16:29:17 UTC (rev 53209)
@@ -1,5 +1,6 @@
#include<stdio.h>
#include<stdlib.h>
+#include<string.h>
#include <math.h>
#include <grass/gis.h>
@@ -8,199 +9,239 @@
#include "local_proto.h"
#include "earth_sun.h"
-#define ETM_MET_SIZE 5600 /* .met file size 5516 bytes */
-#define TM5_MET_SIZE 28700 /* .met file size 28686 bytes */
-#define MAX_STR 127
+#define MAX_STR 127
+#define METADATA_SIZE 65535 /* MTL.txt file size 65535 bytes */
+#define TM5_MET_SIZE 28700 /* .met file size 28686 bytes */
-inline void chrncpy(char *dest, char *src, int n)
+
+inline void chrncpy(char *dest, char src[], int n)
{
if (src == NULL)
- n = 1;
+ {
+ dest[0] = '\0';
+ }
else
- strncpy(dest, src, n);
- dest[n - 1] = '\0';
+ {
+ int i;
+ for (i = 0; i < n && src[i] != '\0' && src[i] != '\"'; i++) dest[i] = src[i];
+ dest[i] = '\0';
+ }
}
+
/****************************************************************************
- * PURPOSE: Read values of Landsat-7 ETM+ from header (.met) file
+ * PURPOSE: Read values of Landsat MSS/TM from header (.met) file
*****************************************************************************/
-void get_value_met7(const char mettext[], char *text, char value[])
+void get_metdata(const char mettext[], char *text, char value[])
{
char *ptr;
- value[0] = '\0';
+
ptr = strstr(mettext, text);
if (ptr == NULL)
- return;
+ {
+ value[0] = 0;
+ return;
+ }
- while (*ptr++ != '=') ;
- sscanf(ptr, "%s", value);
+ ptr = strstr(ptr, " VALUE ");
+ if (ptr == NULL) return;
+ while (*ptr++ != '\"') ;
+ int i = 0;
+ while (*ptr != '\"' && i < MAX_STR) value[i++] = *ptr++;
+ value[i] = '\0';
+
return;
}
-void met_ETM(char *metfile, lsat_data * lsat)
+void lsat_metdata(char *metfile, lsat_data * lsat)
{
FILE *f;
- char mettext[ETM_MET_SIZE];
+ char mettext[TM5_MET_SIZE];
char name[MAX_STR], value[MAX_STR];
- int i;
- static int band[] = { 1, 2, 3, 4, 5, 6, 6, 7, 8 };
- static int code[] = { 1, 2, 3, 4, 5, 61, 62, 7, 8 };
+ /* char metdate[MAX_STR]; */
- static double esun[] =
- { 1969., 1840., 1551., 1044., 225.7, 0., 82.07, 1368. };
-
if ((f = fopen(metfile, "r")) == NULL)
- G_fatal_error(_("Metadata file <%s> not found"), metfile);
+ G_fatal_error(_("Metadata file <%s> not found"), metfile);
- fread(mettext, 1, ETM_MET_SIZE, f);
+ fread(mettext, TM5_MET_SIZE, 1, f);
/* --------------------------------------- */
- lsat->number = 7;
+ get_metdata(mettext, "PLATFORMSHORTNAME", value);
+ chrncpy(name, value + 8, 1);
+ lsat->number = atoi(name);
- get_value_met7(mettext, "SENSOR_ID", value);
- chrncpy(lsat->sensor, value + 1, 5);
+ get_metdata(mettext, "SENSORSHORTNAME", value);
+ chrncpy(lsat->sensor, value + 1, 4);
- if (lsat->creation[0] == 0) {
- get_value_met7(mettext, "CREATION_TIME", value);
- chrncpy(lsat->creation, value, 11);
+ get_metdata(mettext, "CALENDARDATE", value);
+ chrncpy(lsat->date, value, 10);
+
+ if (lsat->creation[0] == 0)
+ {
+ get_metdata(mettext, "PRODUCTIONDATETIME", value);
+ if (!value[0])
+ G_fatal_error(_("Product creation date not in metadata file <%s>, input this data in the command line parameters"), metfile);
+ chrncpy(lsat->creation, value, 10);
}
- get_value_met7(mettext, "ACQUISITION_DATE", value);
- chrncpy(lsat->date, value, 11);
- lsat->dist_es = earth_sun(lsat->date);
+ if (lsat->sun_elev == 0)
+ {
+ get_metdata(mettext, "SolarElevation", value);
+ if (!value[0])
+ G_fatal_error(_("Unable to read solar elevation from metadata file in metadata file <%s>, input this data in the command line parameters"), metfile);
+ lsat->sun_elev = atof(value);
+ }
- get_value_met7(mettext, "SUN_ELEVATION", value);
- lsat->sun_elev = atof(value);
-
- lsat->bands = 9;
- for (i = 0; i < lsat->bands; i++) {
- lsat->band[i].number = *(band + i);
- lsat->band[i].code = *(code + i);
- lsat->band[i].esun = *(esun + lsat->band[i].number - 1);
- snprintf(name, MAX_STR, "LMAX_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].lmax = atof(value);
- snprintf(name, MAX_STR, "LMIN_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].lmin = atof(value);
- snprintf(name, MAX_STR, "QCALMAX_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].qcalmax = atof(value);
- snprintf(name, MAX_STR, "QCALMIN_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].qcalmin = atof(value);
- if (lsat->band[i].number == 6) {
- lsat->band[i].thermal = 1;
- lsat->band[i].K1 = 666.09;
- lsat->band[i].K2 = 1282.71;
- }
- else {
- lsat->band[i].thermal = 0;
- }
+ // Fill data with the sensor_XXX functions
+ switch(lsat->number)
+ {
+ case 1:
+ set_MSS1(lsat);
+ break;
+ case 2:
+ set_MSS2(lsat);
+ break;
+ case 3:
+ set_MSS3(lsat);
+ break;
+ case 4:
+ if (lsat->sensor[0] == 'M')
+ set_MSS4(lsat);
+ else
+ set_TM4(lsat);
+ break;
+ case 5:
+ if (lsat->sensor[0] == 'M')
+ set_MSS5(lsat);
+ else
+ set_TM5(lsat);
+ break;
+ default:
+ G_warning("Unable to recognize satellite platform [%d]", lsat->number);
+ break;
}
+ /* --------------------------------------- */
(void)fclose(f);
return;
}
+
/****************************************************************************
- * PURPOSE: Read values of Landsat MSS/TM from header (.met) file
+ * PURPOSE: Read values of Landsat from MTL metadata (MTL.txt) file
*****************************************************************************/
-void get_value_met(const char mettext[], char *text, char value[])
+
+void get_mtldata(const char mtltext[], char *text, char value[])
{
char *ptr;
- int i;
- value[0] = 0;
-
- ptr = strstr(mettext, text);
+ ptr = strstr(mtltext, text);
if (ptr == NULL)
- return;
+ {
+ value[0] = '\0';
+ return;
+ }
- ptr = strstr(ptr, " VALUE ");
- if (ptr == NULL)
- return;
-
- i = 0;
- while (*ptr++ != '\"') ;
- while (*ptr != '\"' && i < MAX_STR)
- value[i++] = *ptr++;
+ while (*ptr++ != '=') ;
+ while (*ptr <= ' ' || *ptr == '\"') *ptr++;
+ int i = 0;
+ while (i < MAX_STR && *ptr != '\"' && *ptr > ' ') value[i++] = *ptr++;
value[i] = '\0';
return;
}
-void met_TM5(char *metfile, lsat_data * lsat)
+void lsat_mtldata(char *mtlfile, lsat_data * lsat)
{
FILE *f;
- char mettext[TM5_MET_SIZE];
- char value[MAX_STR];
+ char mtldata[METADATA_SIZE];
+ char name[MAX_STR], value[MAX_STR];
+ int i;
- /* char metdate[MAX_STR]; */
+ if ((f = fopen(mtlfile, "r")) == NULL)
+ G_fatal_error(_("Metadata file <%s> not found"), mtlfile);
- if ((f = fopen(metfile, "r")) == NULL)
- G_fatal_error(_("Metadata file <%s> not found"), metfile);
+ fread(mtldata, METADATA_SIZE, 1, f);
- fread(mettext, 1, TM5_MET_SIZE, f);
-
/* --------------------------------------- */
- get_value_met(mettext, "CALENDARDATE", value);
- chrncpy(lsat->date, value, 11);
+ get_mtldata(mtldata, "SPACECRAFT_ID", value);
+ chrncpy(name, value + 7, 1);
+ lsat->number = atoi(name);
- if (lsat->creation[0] == 0) {
- get_value_met(mettext, "PRODUCTIONDATETIME", value);
- chrncpy(lsat->creation, value, 11);
- }
+ get_mtldata(mtldata, "SENSOR_ID", value);
+ chrncpy(lsat->sensor, value, 4);
- if (lsat->creation[0] == 0)
- G_fatal_error(_("Product creation date not in metadata file <%s>"),
- metfile);
- G_debug(1, "met_TM5: Product creation date = [%s]", lsat->creation);
+ get_mtldata(mtldata, "ACQUISITION_DATE", value);
+ chrncpy(lsat->date, value, 10);
+ get_mtldata(mtldata, "CREATION_TIME", value);
+ chrncpy(lsat->creation, value, 10);
- get_value_met(mettext, "SolarElevation", value);
- if (!value[0])
- G_warning("Unable to read solar elevation from metadata file");
- else
- lsat->sun_elev = atof(value);
- G_debug(1, "met_TM5: value=[%s], SolarElevation = %.2f", value,
- lsat->sun_elev);
+ get_mtldata(mtldata, "SUN_ELEVATION", value);
+ lsat->sun_elev = atof(value);
+ // Fill data with the sensor_XXX functions
+ switch(lsat->number)
+ {
+ case 1:
+ set_MSS1(lsat);
+ break;
+ case 2:
+ set_MSS2(lsat);
+ break;
+ case 3:
+ set_MSS3(lsat);
+ break;
+ case 4:
+ if (lsat->sensor[0] == 'M')
+ set_MSS4(lsat);
+ else
+ set_TM4(lsat);
+ break;
+ case 5:
+ if (lsat->sensor[0] == 'M')
+ set_MSS5(lsat);
+ else
+ set_TM5(lsat);
+ break;
+ case 7:
+ get_mtldata(mtldata, "BAND1_GAIN", value);
+ get_mtldata(mtldata, "BAND2_GAIN", value + 1);
+ get_mtldata(mtldata, "BAND3_GAIN", value + 2);
+ get_mtldata(mtldata, "BAND4_GAIN", value + 3);
+ get_mtldata(mtldata, "BAND5_GAIN", value + 4);
+ get_mtldata(mtldata, "BAND6_GAIN1", value + 5);
+ get_mtldata(mtldata, "BAND6_GAIN2", value + 6);
+ get_mtldata(mtldata, "BAND7_GAIN", value + 7);
+ get_mtldata(mtldata, "BAND8_GAIN", value + 8);
+ value[9] = '\0';
+ set_ETM(lsat, value);
+ break;
+ default:
+ G_warning("Unable to recognize satellite platform [%d]", lsat->number);
+ break;
+ }
- get_value_met(mettext, "PLATFORMSHORTNAME", value);
- G_debug(1, "met_TM5: PLATFORMSHORTNAME=[%s]", value);
- switch (value[8]) {
- case '1':
- set_MSS1(lsat);
- break;
- case '2':
- set_MSS2(lsat);
- break;
- case '3':
- set_MSS3(lsat);
- break;
- case '4':
- get_value_met(mettext, "SENSORSHORTNAME", value);
- if (value[0] == 'M')
- set_MSS4(lsat);
- else
- set_TM4(lsat);
- break;
- case '5':
- get_value_met(mettext, "SENSORSHORTNAME", value);
- if (value[0] == 'M')
- set_MSS5(lsat);
- else
- set_TM5(lsat);
- break;
- default:
- G_warning("Unable to recognize satellite platform [%s]", value);
- break;
+ // Update the information from metadata file
+ for (i = 0; i < lsat->bands; i++) {
+ snprintf(name, MAX_STR, "LMAX_BAND%d", lsat->band[i].code);
+ get_mtldata(mtldata, name, value);
+ lsat->band[i].lmax = atof(value);
+ snprintf(name, MAX_STR, "LMIN_BAND%d", lsat->band[i].code);
+ get_mtldata(mtldata, name, value);
+ lsat->band[i].lmin = atof(value);
+ snprintf(name, MAX_STR, "QCALMAX_BAND%d", lsat->band[i].code);
+ get_mtldata(mtldata, name, value);
+ lsat->band[i].qcalmax = atof(value);
+ snprintf(name, MAX_STR, "QCALMIN_BAND%d", lsat->band[i].code);
+ get_mtldata(mtldata, name, value);
+ lsat->band[i].qcalmin = atof(value);
}
+ /* --------------------------------------- */
(void)fclose(f);
return;
@@ -208,78 +249,3 @@
-/****************************************************************************
- * PURPOSE: Read values of Landsat TM5 from header (.mtl) file
- *****************************************************************************/
-
-/****************************************************************************
- * EXPLANATION: This module is a modification of the met_ETM() found before
- * to allow TM5 from GLOVIS to use .MTL extension that responds
- * near to perfectly to the .met parser. While L7 files using
- * .MTL from GLOVIS can be processed as if having .met files
- * seemlessly, TM5 using .MTL need to read basic info and
- * additionally the LMIN, LMAX, QCALMIN, QCALMAX being explicitely
- * provided in the .MTL as if in a .met file.
- *****************************************************************************/
-void mtl_TM5(char *metfile, lsat_data * lsat)
-{
- FILE *f;
- char mettext[ETM_MET_SIZE];
- char name[MAX_STR], value[MAX_STR];
- int i;
-
- static int code[] = { 1, 2, 3, 4, 5, 6, 7 };
-
- if ((f = fopen(metfile, "r")) == NULL)
- G_fatal_error(_("Metadata file <%s> not found"), metfile);
-
- fread(mettext, 1, ETM_MET_SIZE, f);
-
- /* --------------------------------------- */
- get_value_met7(mettext, "SENSOR_ID", value);
- if (value[1] == 'M')
- chrncpy(lsat->sensor, value + 1, 4);
- else
- chrncpy(lsat->sensor, value + 1, 3);
-
- if (lsat->creation[0] == 0) {
- get_value_met7(mettext, "PRODUCT_CREATION_TIME", value);
- chrncpy(lsat->creation, value, 11);
- }
-
- get_value_met7(mettext, "ACQUISITION_DATE", value);
- chrncpy(lsat->date, value, 11);
- lsat->dist_es = earth_sun(lsat->date);
-
- get_value_met7(mettext, "SUN_ELEVATION", value);
- lsat->sun_elev = atof(value);
-
- /* We still have to initialize most of the info */
- /* So instead of rewriting a new function, we use set_TM5()... */
- set_TM5(lsat);
- /* ... and we rewrite the necessary 'a la Landsat 7' */
-
- if (strcmp(lsat->sensor, "MSS") == 0)
- lsat->bands = 4;
- else
- lsat->bands = 7;
- for (i = 0; i < lsat->bands; i++) {
- lsat->band[i].code = *(code + i);
- snprintf(name, MAX_STR, "LMAX_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].lmax = atof(value);
- snprintf(name, MAX_STR, "LMIN_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].lmin = atof(value);
- snprintf(name, MAX_STR, "QCALMAX_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].qcalmax = atof(value);
- snprintf(name, MAX_STR, "QCALMIN_BAND%d", lsat->band[i].code);
- get_value_met7(mettext, name, value);
- lsat->band[i].qcalmin = atof(value);
- if (lsat->band[i].number == 6)
- lsat->band[i].thermal = 1;
- }
- (void)fclose(f);
- return;
-}
Modified: grass-addons/grass6/imagery/i.landsat.toar/landsat_set.c
===================================================================
--- grass-addons/grass6/imagery/i.landsat.toar/landsat_set.c 2012-09-18 15:09:24 UTC (rev 53208)
+++ grass-addons/grass6/imagery/i.landsat.toar/landsat_set.c 2012-09-18 16:29:17 UTC (rev 53209)
@@ -18,7 +18,7 @@
int code[] = { 4, 5, 6, 7 };
double wmax[] = { 0.6, 0.7, 0.8, 1.1 };
double wmin[] = { 0.5, 0.6, 0.7, 0.8 };
- /* 79-82, 79-82, 79-82, 79-82 */
+ /* 68x83, 68x83, 68x83, 68x83 */
strcpy(lsat->sensor, "MSS");
@@ -89,7 +89,7 @@
/** **********************************************
- ** Before access to this function ...
+ ** Before access to these functions ...
** store previously
** >>> adquisition date,
** >>> creation date, and
@@ -98,21 +98,20 @@
/****************************************************************************
* PURPOSE: Store values of Landsat-1 MSS
- * July 23, 1972 to January 7, 1978
+ * July 23, 1972 to January 6, 1978
*****************************************************************************/
void set_MSS1(lsat_data * lsat)
{
int i, j;
- /** Markham and Barker. EOSAT Landsat Technical Notes, No. 1, 1986;
- Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009)*/
+ /** USGS Calibration Parameter Files 2012 */
/* Spectral radiances at detector */
double lmax[] = { 248., 200., 176., 153. };
double lmin[] = { 0., 0., 0., 0. };
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1823., 1559., 1276., 880.1 };
+ double esun[] = { 1824., 1570., 1249., 853.4 };
lsat->number = 1;
sensor_MSS(lsat);
@@ -125,33 +124,33 @@
lsat->band[i].lmax = *(lmax + j);
lsat->band[i].lmin = *(lmin + j);
}
+ G_debug(1, "Landsat-1 MSS");
return;
}
/****************************************************************************
* PURPOSE: Store values of Landsat-2 MSS
- * January 22, 1975 to February 25, 1982
+ * January 22, 1975 to February 5, 1982
*****************************************************************************/
void set_MSS2(lsat_data * lsat)
{
int i, j;
double julian, *lmax, *lmin;
- /** Markham and Barker. EOSAT Landsat Technical Notes, No. 1, 1986;
- Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009)*/
+ /** USGS Calibration Parameter Files 2012 */
/* Spectral radiances at detector */
double Lmax[][4] = {
{210., 156., 140., 138.}, /* before July 16, 1975 */
- {263., 176., 152., 130.333} /* on or after July 16, 1975 */
+ {263., 176., 152., 130.} /* on or after July 16, 1975 */
};
double Lmin[][4] = {
{10., 7., 7., 5.},
- {8., 6., 6., 3.667}
+ {8., 6., 6., 4.}
};
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1829., 1539., 1268., 886.6 };
+ double esun[] = { 1824., 1570., 1249., 853.4 };
julian = julian_char(lsat->creation);
if (julian < julian_char("1975-07-16"))
@@ -172,6 +171,7 @@
lsat->band[i].lmax = *(lmax + j);
lsat->band[i].lmin = *(lmin + j);
}
+ G_debug(1, "Landsat-2 MSS");
return;
}
@@ -186,8 +186,8 @@
int i, j;
double julian, *lmax, *lmin;
- /** Markham and Barker. EOSAT Landsat Technical Notes, No. 1, 1986;
- Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009)*/
+ /** USGS Calibration Parameter Files 2012 */
+
/* Spectral radiances at detector */
double Lmax[][4] = {
{220., 175., 145., 147.}, /* before June 1, 1978 */
@@ -198,7 +198,7 @@
{4., 3., 3., 1.}
};
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1839., 1555., 1291., 887.9 };
+ double esun[] = { 1824., 1570., 1249., 853.4 };
julian = julian_char(lsat->creation);
if (julian < julian_char("1978-06-01"))
@@ -219,6 +219,7 @@
lsat->band[i].lmax = *(lmax + j);
lsat->band[i].lmin = *(lmin + j);
}
+ G_debug(1, "Landsat-3 MSS");
return;
}
@@ -231,8 +232,7 @@
int i, j;
double julian, *lmax, *lmin;
- /** Markham and Barker. EOSAT Landsat Technical Notes, No. 1, 1986;
- Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009)*/
+ /** USGS Calibration Parameter Files 2012 */
/* Spectral radiances at detector */
double Lmax[][4] = {
@@ -247,15 +247,15 @@
};
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1827., 1569., 1260., 866.4 };
+ double esun[] = { 1824., 1570., 1249., 853.4 };
julian = julian_char(lsat->creation);
if (julian < julian_char("1982-08-26"))
- i = 0;
+ i = 0;
else if (julian < julian_char("1983-03-31"))
- i = 1;
+ i = 1;
else
- i = 2;
+ i = 2;
lmax = Lmax[i];
lmin = Lmin[i];
@@ -270,6 +270,7 @@
lsat->band[i].lmax = *(lmax + j);
lsat->band[i].lmin = *(lmin + j);
}
+ G_debug(1, "Landsat-4 MSS");
return;
}
@@ -278,34 +279,32 @@
int i, j;
double julian, *lmax, *lmin;
- /** Brian L. Markham and John L. Barker.
- EOSAT Landsat Technical Notes, No. 1, 1986 */
+ /** USGS Calibration Parameter Files 2012 */
+
/* Spectral radiances at detector */
double Lmax[][7] = {
{158.42, 308.17, 234.63, 224.32, 32.42, 15.64, 17.00}, /* before August 1983 */
{142.86, 291.25, 225.00, 214.29, 30.00, 12.40, 15.93}, /* before January 15, 1984 */
- {152.10, 296.81, 204.30, 206.20, 27.19, 15.3032, 14.38} /* after Jaunary 15, 1984 */
+ {171.00, 336.00, 254.00, 221.00, 31.40, 15.303, 16.60} /* after Jaunary 15, 1984 */
};
double Lmin[][7] = {
{-1.52, -2.84, -1.17, -1.51, -0.37, 2.00, -0.15},
{0.00, 0.00, 0.00, 0.00, 0.00, 4.84, 0.00},
- {-1.50, -2.80, -1.20, -1.50, -0.37, 1.2378, -0.15}
+ {-1.52, -2.84, -1.17, -1.51, -0.37, 1.2378, -0.15}
};
- /** Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009) */
-
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1983., 1795., 1539., 1028., 219.8, 0., 83.49 };
+ double esun[] = { 1957., 1825., 1557., 1033., 214.9, 0., 80.72 };
/* Thermal band calibration constants: K1 = 671.62 K2 = 1284.30 */
julian = julian_char(lsat->creation);
if (julian < julian_char("1983-08-01"))
- i = 0;
+ i = 0;
else if (julian < julian_char("1984-01-15"))
- i = 1;
+ i = 1;
else
- i = 2;
+ i = 2;
lmax = Lmax[i];
lmin = Lmin[i];
@@ -324,6 +323,7 @@
lsat->band[i].K2 = 1284.30;
}
}
+ G_debug(1, "Landsat-4 TM");
return;
}
@@ -337,8 +337,8 @@
int i, j;
double julian, *lmax, *lmin;
- /** Brian L. Markham and John L. Barker.
- EOSAT Landsat Technical Notes, No. 1, 1986 */
+ /** USGS Calibration Parameter Files 2012 */
+
/* Spectral radiances at detector */
double Lmax[][4] = {
{240., 170., 150., 127.}, /* before April 6, 1984 */
@@ -357,7 +357,7 @@
julian = julian_char(lsat->creation);
if (julian < julian_char("1984-04-06"))
i = 0;
- else if (julian < julian_char("1984-11-09"))
+ else if (julian < julian_char("1984-11-08"))
i = 1;
else
i = 2;
@@ -375,6 +375,7 @@
lsat->band[i].lmax = *(lmax + j);
lsat->band[i].lmin = *(lmin + j);
}
+ G_debug(1, "Landsat-5 MSS");
return;
}
@@ -383,8 +384,7 @@
int i, j;
double julian, *lmax, *lmin, jbuf;
- /** Gyanesh Chander and Brian Markham.
- IEEE Transactions On Geoscience And Remote Sensing, Vol. 41, No. 11, November 2003 */
+ /** USGS Calibration Parameter Files 2012 */
/* Spectral radiances at detector */
double Lmax[][7] = {
@@ -398,20 +398,18 @@
{-1.52, -2.84, -1.17, -1.51, -0.37, 1.2378, -0.15}
};
- /** Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009) */
-
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1983., 1796., 1536., 1031., 220.0, 0., 83.44 };
+ double esun[] = { 1957., 1826., 1554., 1036., 215.0, 0., 80.67 };
/* Thermal band calibration constants: K1 = 607.76 K2 = 1260.56 */
julian = julian_char(lsat->creation);
if (julian < julian_char("2003-05-04"))
- i = 0;
+ i = 0;
else if (julian < julian_char("2007-04-02"))
- i = 1;
+ i = 1;
else
- i = 2;
+ i = 2;
lmax = Lmax[i];
lmin = Lmin[i];
if (i == 2) { /* in Chander, Markham and Barsi 2007 */
@@ -445,21 +443,21 @@
lsat->band[i].K2 = 1260.56;
}
}
+ G_debug(1, "Landsat-5 TM");
return;
}
/****************************************************************************
* PURPOSE: Store values of Landsat-7 ETM+
- * April 15, 1999 to today
+ * April 15, 1999 to May 31, 2003 (SLC failure)
*****************************************************************************/
void set_ETM(lsat_data * lsat, char gain[])
{
int i, k, j;
double julian, *lmax, *lmin;
- /** Richard Irish.
- Landsat 7. Science Data Users Handbook. Last update: February 17, 2007 */
+ /** USGS Calibration Parameter Files 2012 */
/* Spectral radiances at detector */
/* - LOW GAIN - */
@@ -481,10 +479,8 @@
{-6.2, -6.4, -5.0, -5.1, -1.0, 3.2, -0.35, -4.7}
};
- /** Chander, Markham and Helder. Remote Sensing of Environment, 113 (2009) */
-
/* Solar exoatmospheric spectral irradiances */
- double esun[] = { 1997., 1812., 1533., 1039., 230.8, 0., 84.90, 1362. };
+ double esun[] = { 1969., 1840., 1551., 1044., 225.7, 0., 82.07, 1368. };
/* Thermal band calibration constants: K1 = 666.09 K2 = 1282.71 */
@@ -517,5 +513,6 @@
lsat->band[i].K2 = 1282.71;
}
}
+ G_debug(1, "Landsat-7 ETM+");
return;
}
Modified: grass-addons/grass6/imagery/i.landsat.toar/local_proto.h
===================================================================
--- grass-addons/grass6/imagery/i.landsat.toar/local_proto.h 2012-09-18 15:09:24 UTC (rev 53208)
+++ grass-addons/grass6/imagery/i.landsat.toar/local_proto.h 2012-09-18 16:29:17 UTC (rev 53209)
@@ -4,9 +4,8 @@
#include <string.h>
#include "landsat.h"
-void met_ETM(char *, lsat_data *);
-void met_TM5(char *, lsat_data *);
-void mtl_TM5(char *, lsat_data *);
+void lsat_mtldata(char *, lsat_data *);
+void lsat_metdata(char *, lsat_data *);
void set_MSS1(lsat_data *);
void set_MSS2(lsat_data *);
@@ -16,6 +15,7 @@
void set_TM4(lsat_data *);
void set_TM5(lsat_data *);
+
void set_ETM(lsat_data *, char[]);
#endif
Modified: grass-addons/grass6/imagery/i.landsat.toar/main.c
===================================================================
--- grass-addons/grass6/imagery/i.landsat.toar/main.c 2012-09-18 15:09:24 UTC (rev 53208)
+++ grass-addons/grass6/imagery/i.landsat.toar/main.c 2012-09-18 16:29:17 UTC (rev 53209)
@@ -4,13 +4,13 @@
* MODULE: i.landsat.toar
*
* AUTHOR(S): E. Jorge Tizado - ej.tizado at unileon.es
- * Hamish Bowman (small grassification cleanups)
- * Yann Chemin (v7 + L5TM _MTL.txt support)
+ * Hamish Bowman (small grassification cleanups)
+ * Yann Chemin (v7 + L5TM _MTL.txt support) [removed after update]
*
* PURPOSE: Calculate TOA Radiance or Reflectance and Kinetic Temperature
* for Landsat 1/2/3/4/5 MS, 4/5 TM or 7 ETM+
*
- * COPYRIGHT: (C) 2002, 2005, 2008, 2010 by the GRASS Development Team
+ * COPYRIGHT: (C) 2012 by the GRASS Development Team
*
* This program is free software under the GNU General
* Public License (>=v2). Read the file COPYING that
@@ -43,7 +43,7 @@
struct Option *input_prefix, *output_prefix, *metfn, *sensor, *adate, *pdate, *elev,
*bgain, *metho, *perc, *dark, *atmo;
char *inputname, *met, *outputname, *sensorname;
- struct Flag *msss, *frad, *l5_mtl;
+ struct Flag *frad;
lsat_data lsat;
char band_in[GNAME_MAX], band_out[GNAME_MAX];
@@ -82,35 +82,37 @@
metfn = G_define_standard_option(G_OPT_F_INPUT);
metfn->key = "metfile";
metfn->required = NO;
- metfn->description = _("Name of Landsat ETM+ or TM5 header file (.met/MTL.txt)");
+ metfn->description = _("Name of Landsat metadata file (.met or MTL.txt)");
metfn->guisection = _("Metadata");
+ metho = G_define_option();
+ metho->key = "method";
+ metho->type = TYPE_STRING;
+ metho->required = NO;
+ metho->options = "uncorrected,corrected,dos1,dos2,dos2b,dos3,dos4";
+ metho->label = _("Atmospheric correction method");
+ metho->description = _("Atmospheric correction method");
+ metho->answer = "uncorrected";
+ metho->guisection = _("Metadata");
+
sensor = G_define_option();
sensor->key = "sensor";
sensor->type = TYPE_STRING;
sensor->label = _("Spacecraft sensor");
sensor->description = _("Required only if 'metfile' not given");
- sensor->options = "mss1,mss2,mss3,tm4,tm5,tm7";
+ sensor->options = "mss1,mss2,mss3,mss4,mss5,tm4,tm5,tm7";
sensor->descriptions =
_("mss1;Landsat-1 MSS;"
"mss2;Landsat-2 MSS;"
- "mss3;Landsat-3 MSS;"
- "tm4;Landsat-4 TM;"
- "tm5;Landsat-5 TM;"
+ "mss3;Landsat-3 MSS;"
+ "mss4;Landsat-4 MSS;"
+ "mss5;Landsat-5 MSS;"
+ "tm4;Landsat-4 TM;"
+ "tm5;Landsat-5 TM;"
"tm7;Landsat-7 ETM+");
sensor->required = NO;
sensor->guisection = _("Metadata");
- metho = G_define_option();
- metho->key = "method";
- metho->type = TYPE_STRING;
- metho->required = NO;
- metho->options = "uncorrected,corrected,dos1,dos2,dos2b,dos3,dos4";
- metho->label = _("Atmospheric correction method");
- metho->description = _("Required only if 'metfile' not given");
- metho->answer = "uncorrected";
- metho->guisection = _("Metadata");
-
adate = G_define_option();
adate->key = "date";
adate->type = TYPE_STRING;
@@ -121,10 +123,10 @@
adate->guisection = _("Metadata");
elev = G_define_option();
- elev->key = "solar_elevation";
+ elev->key = "sun_elevation";
elev->type = TYPE_DOUBLE;
elev->required = NO;
- elev->label = _("Solar elevation in degrees");
+ elev->label = _("Sun elevation in degrees");
elev->description = _("Required only if 'metfile' not given");
elev->guisection = _("Metadata");
@@ -141,8 +143,8 @@
bgain->key = "gain";
bgain->type = TYPE_STRING;
bgain->required = NO;
- bgain->label =
- _("Gain (H/L) of all Landsat ETM+ bands (1-5,61,62,7,8)");
+ bgain->label = _("Gain (H/L) of all Landsat ETM+ bands (1-5,61,62,7,8)");
+ bgain->description = _("Required only if 'metfile' not given");
bgain->guisection = _("Settings");
perc = G_define_option();
@@ -166,25 +168,16 @@
atmo->key = "rayleigh";
atmo->type = TYPE_DOUBLE;
atmo->required = NO;
- atmo->description = _("Rayleigh atmosphere"); /* scattering coefficient? */
+ atmo->description = _("Rayleigh atmosphere (diffuse sky irradiance)"); /* scattering coefficient? */
atmo->answer = "0.0";
atmo->guisection = _("Settings");
+
/* define the different flags */
frad = G_define_flag();
frad->key = 'r';
frad->description = _("Output at-sensor radiance for all bands");
- msss = G_define_flag();
- msss->key = 's';
- msss->description = _("Set sensor of Landsat TM4/5 to MSS");
- msss->guisection = _("Settings");
-
- l5_mtl = G_define_flag();
- l5_mtl->key = 't';
- l5_mtl->description = _("Landsat ETM+/TM5 has a MTL.txt file instead of .met");
- l5_mtl->guisection = _("Metadata");
-
/* options and afters parser */
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
@@ -197,7 +190,7 @@
met = metfn->answer;
inputname = input_prefix->answer;
outputname = output_prefix->answer;
- sensorname = sensor -> answer ? sensor->answer: "";
+ sensorname = sensor->answer ? sensor->answer: "";
G_zero(&lsat, sizeof(lsat));
@@ -226,24 +219,30 @@
pixel = atoi(dark->answer);
rayleigh = atof(atmo->answer);
- /* Data from MET file: only Landsat-7 ETM+ and Landsat-5 TM */
- if (met != NULL) {
- if (strcmp(sensorname, "tm7") == 0)
- met_ETM(met, &lsat);
- else if (l5_mtl->answer)
- mtl_TM5(met, &lsat);
- else
- met_TM5(met, &lsat);
+ /* Data from metadata file */
+ if (met != NULL)
+ {
+ i = strlen(met);
+ if (strcmp(met + i - 7, "MTL.txt") == 0)
+ {
+ lsat_mtldata(met, &lsat);
+ }
+ else if (strcmp(met + i - 4, ".met") == 0)
+ {
+ if (strcmp(sensorname, "tm7") == 0)
+ lsat_mtldata(met, &lsat); /* .met of Landsat-7 = new MTL file */
+ else
+ lsat_metdata(met, &lsat);
+ }
+ G_debug(1, "lsat.number = %d, lsat.sensor = [%s]", lsat.number, lsat.sensor);
+ if (!lsat.sensor || lsat.number > 7 || lsat.number < 1)
+ G_fatal_error(_("Failed to identify satellite"));
- G_debug(1, "lsat.number = %d, lsat.sensor = [%s]", lsat.number,
- lsat.sensor);
- if (!lsat.sensor || lsat.number > 7 || lsat.number < 1)
- G_fatal_error(_("Failed to identify satellite"));
+ G_debug(1, "Landsat-%d %s with data set in met file [%s]", lsat.number, lsat.sensor, met);
- G_debug(1, "Landsat-%d %s with data set in met file [%s]",
- lsat.number, lsat.sensor, met);
- if (elev->answer != NULL)
- lsat.sun_elev = atof(elev->answer); /* Overwrite solar elevation of met file */
+ /* Overwrite solar elevation of metadata file */
+ if (elev->answer != NULL)
+ lsat.sun_elev = atof(elev->answer);
}
/* Data from date and solar elevation */
else if (adate->answer == NULL || elev->answer == NULL) {
@@ -252,61 +251,48 @@
}
else {
if (strcmp(sensorname, "tm7") == 0) { /* Need gain */
- if (bgain->answer != NULL && strlen(bgain->answer) == 9) {
- set_ETM(&lsat, bgain->answer);
- G_debug(1, "Landsat 7 ETM+");
- }
- else {
- G_fatal_error(_("Landsat-7 requires band gain with 9 (H/L) data"));
- }
- }
- else { /* Not need gain */
- if (strcmp(sensorname, "tm5") == 0) {
- if (msss->answer)
- set_MSS5(&lsat);
- else
- set_TM5(&lsat);
- G_debug(1, "Landsat-5 %s", lsat.sensor);
- }
- else if (strcmp(sensorname, "tm4") == 0) {
- if (msss->answer)
- set_MSS4(&lsat);
- else
- set_TM4(&lsat);
- G_debug(1, "Landsat-4 %s", lsat.sensor);
- }
- else if (strcmp(sensorname, "mss3") == 0) {
- set_MSS3(&lsat);
- G_debug(1, "Landsat-3 MSS");
- }
- else if (strcmp(sensorname, "mss2") == 0) {
- set_MSS2(&lsat);
- G_debug(1, "Landsat-2 MSS");
- }
- else if (strcmp(sensorname, "mss1") == 0) {
- set_MSS1(&lsat);
- G_debug(1, "Landsat-1 MSS");
- }
- else {
- G_fatal_error(_("Unknown satellite type (defined by '%s')"), sensor->key);
- }
- }
+ if (bgain->answer != NULL && strlen(bgain->answer) == 9) {
+ set_ETM(&lsat, bgain->answer);
+ G_debug(1, "Landsat 7 ETM+");
+ }
+ else {
+ G_fatal_error(_("Landsat-7 requires band gain with 9 (H/L) data"));
+ }
}
+ else { /* Not need gain */
+ if (strcmp(sensorname, "tm5") == 0)
+ set_TM5(&lsat);
+ else if (strcmp(sensorname, "tm4") == 0)
+ set_TM4(&lsat);
+ else if (strcmp(sensorname, "mss5") == 0)
+ set_MSS5(&lsat);
+ else if (strcmp(sensorname, "mss4") == 0)
+ set_MSS4(&lsat);
+ else if (strcmp(sensorname, "mss3") == 0)
+ set_MSS3(&lsat);
+ else if (strcmp(sensorname, "mss2") == 0)
+ set_MSS2(&lsat);
+ else if (strcmp(sensorname, "mss1") == 0)
+ set_MSS1(&lsat);
+ else
+ G_fatal_error(_("Unknown satellite type (defined by '%s')"), sensor->key);
+ }
+ }
/*****************************************
* ------------ PREPARATION --------------
*****************************************/
- if (G_strcasecmp(metho->answer, "corrected") == 0)
+ if (strcasecmp(metho->answer, "corrected") == 0)
method = CORRECTED;
- else if (G_strcasecmp(metho->answer, "dos1") == 0)
+ else if (strcasecmp(metho->answer, "dos1") == 0)
method = DOS1;
- else if (G_strcasecmp(metho->answer, "dos2") == 0)
+ else if (strcasecmp(metho->answer, "dos2") == 0)
method = DOS2;
- else if (G_strcasecmp(metho->answer, "dos2b") == 0)
+ else if (strcasecmp(metho->answer, "dos2b") == 0)
method = DOS2b;
- else if (G_strcasecmp(metho->answer, "dos3") == 0)
+ else if (strcasecmp(metho->answer, "dos3") == 0)
method = DOS3;
- else if (G_strcasecmp(metho->answer, "dos4") == 0)
+ else if (strcasecmp(metho->answer, "dos4") == 0)
method = DOS4;
else
method = UNCORRECTED;
@@ -402,7 +388,7 @@
G_fatal_error(_("Unknown production date (defined by '%s')"), pdate->key);
if (G_verbose() > G_verbose_std()) {
- fprintf(stderr, " SENSOR: %s\n", lsat.sensor);
+ fprintf(stderr, " LANSAT: %d SENSOR: %s\n", lsat.number, lsat.sensor);
fprintf(stderr, " ACQUISITION DATE %s [production date %s]\n",
lsat.date, lsat.creation);
fprintf(stderr, " earth-sun distance = %.8lf\n", lsat.dist_es);
More information about the grass-commit
mailing list