[GRASS-SVN] r31857 - in grass-addons/gipe: . i.dn2full.l7 i.eb.eta i.water

svn_grass at osgeo.org svn_grass at osgeo.org
Thu Jun 26 21:43:43 EDT 2008


Author: ychemin
Date: 2008-06-26 21:43:43 -0400 (Thu, 26 Jun 2008)
New Revision: 31857

Added:
   grass-addons/gipe/i.water/
   grass-addons/gipe/i.water/Makefile
   grass-addons/gipe/i.water/description.html
   grass-addons/gipe/i.water/main.c
   grass-addons/gipe/i.water/water.c
   grass-addons/gipe/i.water/water_modis.c
Modified:
   grass-addons/gipe/i.dn2full.l7/main.c
   grass-addons/gipe/i.eb.eta/main.c
Log:
water detection module

Modified: grass-addons/gipe/i.dn2full.l7/main.c
===================================================================
--- grass-addons/gipe/i.dn2full.l7/main.c	2008-06-27 01:00:59 UTC (rev 31856)
+++ grass-addons/gipe/i.dn2full.l7/main.c	2008-06-27 01:43:43 UTC (rev 31857)
@@ -48,13 +48,10 @@
 	int nrows, ncols;
 	int row,col;
 
-	int verbose=1;
 	struct GModule *module;
 	struct Option *input,*output;
 	
-	struct Flag *flag1;
 	struct History history; //metadata
-
 	/************************************/
 	/* FMEO Declarations*****************/
 	char history_buf[200];
@@ -74,7 +71,7 @@
 	int i=0,j=0;
 	
 	void *inrast[MAXFILES];
-	unsigned char *outrast[MAXFILES];
+	DCELL *outrast[MAXFILES];
 	RASTER_MAP_TYPE in_data_type[MAXFILES];
 	RASTER_MAP_TYPE out_data_type = DCELL_TYPE; /* 0=numbers  1=text */
 
@@ -110,19 +107,10 @@
 	input->gisprompt  = _("old_file,file,file");
 	input->description= _("Landsat 7ETM+ Header File (.met)");
 
-	output = G_define_option() ;
+	output = G_define_standard_option(G_OPT_R_OUTPUT) ;
 	output->key        = _("output");
-	output->type       = TYPE_STRING;
-	output->required   = YES;
-	output->gisprompt  = _("new,cell,raster");
 	output->description= _("Base name of the output layers (will add .x)");
 
-	/* Define the different flags */
-
-	flag1 = G_define_flag() ;
-	flag1->key         =_('q');
-	flag1->description =_("Quiet");
-
 	/********************/
 	if (G_parser(argc, argv))
 		exit (-1);
@@ -130,7 +118,6 @@
 	metfName	= input->answer;
 	result		= output->answer;
 
-	verbose = (!flag1->answer);
 	//******************************************
 	//Fetch parameters for DN2Rad2Ref correction
 	l7_in_read(metfName,b1,b2,b3,b4,b5,b61,b62,b7,b8,lmin,lmax,qcalmin,qcalmax,&sun_elevation,&sun_azimuth,&day,&month,&year);
@@ -410,9 +397,7 @@
 	DCELL d[MAXFILES];
 	for (row = 0; row < nrows; row++)
 	{
-		if (verbose){
-			G_percent (row, nrows, 2);
-		}
+		G_percent (row, nrows, 2);
 		/* read input map */
 		for (i=0;i<MAXFILES;i++)
 		{
@@ -429,14 +414,14 @@
 				dout[i]=dn2rad_landsat7(lmin[i],lmax[i],qcalmax[i],qcalmin[i],d[i]);
 				if(i==5||i==6){//if band 61/62, process brightness temperature
 					if(dout[i]<=0.0){
-						dout[i]=-999.990;
+						G_set_d_null_value(&outrast[i][col],1);
 					}else{
 						dout[i]=tempk_landsat7(dout[i]);
 					}
 				}else{//process reflectance
 					dout[i]=rad2ref_landsat7(dout[i],doy,sun_elevation,kexo[i]);
 			}
-			((DCELL *) outrast[i])[col] = dout[i];
+			outrast[i][col] = dout[i];
  			}
 		}
 		for(i=0;i<MAXFILES;i++){

Modified: grass-addons/gipe/i.eb.eta/main.c
===================================================================
--- grass-addons/gipe/i.eb.eta/main.c	2008-06-27 01:00:59 UTC (rev 31856)
+++ grass-addons/gipe/i.eb.eta/main.c	2008-06-27 01:43:43 UTC (rev 31857)
@@ -47,7 +47,7 @@
 	int i=0,j=0;
 	
 	void *inrast_rnetday, *inrast_evapfr, *inrast_tempk;
-	unsigned char *outrast1;
+	DCELL *outrast1;
 	RASTER_MAP_TYPE data_type_output=DCELL_TYPE;
 	RASTER_MAP_TYPE data_type_rnetday;
 	RASTER_MAP_TYPE data_type_evapfr;
@@ -180,17 +180,15 @@
 					d_tempk = ((DCELL *) inrast_tempk)[col];
 					break;
 			}
-			if(G_is_d_null_value(&d_rnetday)){
-				((DCELL *) outrast1)[col] = -999.99;
-			}else if(G_is_d_null_value(&d_evapfr)){
-				((DCELL *) outrast1)[col] = -999.99;
-			}else if(G_is_d_null_value(&d_tempk)){
-				((DCELL *) outrast1)[col] = -999.99;
+			if(G_is_d_null_value(&d_rnetday)||
+			G_is_d_null_value(&d_evapfr)||
+			G_is_d_null_value(&d_tempk)){
+				G_set_d_null_value(&outrast1[col],1);
 			}else {
 				/************************************/
 				/* calculate soil heat flux	    */
 				d = et_a(d_rnetday,d_evapfr,d_tempk);
-				((DCELL *) outrast1)[col] = d;
+				outrast1[col] = d;
 			}
 		}
 		if (G_put_raster_row (outfd1, outrast1, data_type_output) < 0)

Added: grass-addons/gipe/i.water/Makefile
===================================================================
--- grass-addons/gipe/i.water/Makefile	                        (rev 0)
+++ grass-addons/gipe/i.water/Makefile	2008-06-27 01:43:43 UTC (rev 31857)
@@ -0,0 +1,14 @@
+MODULE_TOPDIR = ../..
+
+PGM = i.water
+
+LIBES = $(GISLIB)
+DEPENDENCIES = $(GISDEP)
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+ifneq ($(USE_LARGEFILES),)
+	EXTRA_CFLAGS = -D_FILE_OFFSET_BITS=64
+endif
+
+default: cmd

Added: grass-addons/gipe/i.water/description.html
===================================================================
--- grass-addons/gipe/i.water/description.html	                        (rev 0)
+++ grass-addons/gipe/i.water/description.html	2008-06-27 01:43:43 UTC (rev 31857)
@@ -0,0 +1,32 @@
+<H2>DESCRIPTION</H2>
+
+<EM>i.water</EM> calculates a water layer after Xiao et al (2005) and Roy et al (2005).
+It takes input of Diurnal Net Radiation (see r.sun), evaporative fraction (see r.eb.evapfr) and surface skin temprature. 
+
+<H2>REFERENCES</H2>
+
+Xiao X., Boles S., Liu J., Zhuang D., Frokling S., Li C., Salas W., Moore III B. (2005). Mapping paddy rice agriculture in southern China using multi-temporal MODIS images. Remote Sensing of Environment 95:480-492.
+
+Roy D.P., Jin Y., Lewis P.E., Justice C.O. (2005). Prototyping a global algorithm for systematic fire-affected area mapping using MODIS time series data. Remote Sensing of Environment 97:137-162.
+
+<H2>NOTES</H2>
+
+<H2>TODO</H2>
+
+
+<H2>SEE ALSO</H2>
+
+<em>
+<A HREF="r.sun.html">r.sun</A><br>
+<A HREF="i.albedo.html">i.albedo</A><br>
+<A HREF="i.vi.html">i.vi</A><br>
+</em>
+
+
+<H2>AUTHORS</H2>
+
+Yann Chemin, International Rice Research Insitute, The Philippines<BR>
+
+
+<p>
+<i>Last changed: $Date: 2008/06/27 10:00:00 $</i>

Added: grass-addons/gipe/i.water/main.c
===================================================================
--- grass-addons/gipe/i.water/main.c	                        (rev 0)
+++ grass-addons/gipe/i.water/main.c	2008-06-27 01:43:43 UTC (rev 31857)
@@ -0,0 +1,242 @@
+/****************************************************************************
+ *
+ * MODULE:       i.water
+ * AUTHOR(S):    Yann Chemin - yann.chemin at gmail.com
+ * PURPOSE:      Calculates if water is there (value=1)
+ * 		 two versions, 1) generic (albedo,ndvi)
+ * 		 2) Modis (surf_refl_7,ndvi)
+ *
+ * COPYRIGHT:    (C) 2008 by the GRASS Development Team
+ *
+ *               This program is free software under the GNU General Public
+ *   	    	 License (>=v2). Read the file COPYING that comes with GRASS
+ *   	    	 for details.
+ *
+ *****************************************************************************/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <grass/gis.h>
+#include <grass/glocale.h>
+
+
+double water(double albedo, double ndvi);
+double water_modis(double surf_ref_7, double ndvi);
+
+int main(int argc, char *argv[])
+{
+	struct Cell_head cellhd; //region+header info
+	char *mapset; // mapset name
+	int nrows, ncols;
+	int row,col;
+
+	struct GModule *module;
+	struct Option *input1, *input2, *input3, *output1;
+	
+	struct Flag *flag1;	
+	struct History history; //metadata
+	
+	/************************************/
+	/* FMEO Declarations*****************/
+	char *name;   // input raster name
+	char *result1; //output raster name
+	//File Descriptors
+	int infd_ndvi, infd_albedo, infd_ref7;
+	int outfd1;
+	
+	char *ndvi,*albedo,*ref7;
+	int i=0,j=0;
+	
+	void *inrast_ndvi, *inrast_albedo, *inrast_ref7;
+	CELL *outrast1;
+	RASTER_MAP_TYPE data_type_output=CELL_TYPE;
+	RASTER_MAP_TYPE data_type_ndvi;
+	RASTER_MAP_TYPE data_type_albedo;
+	RASTER_MAP_TYPE data_type_ref7;
+	/************************************/
+	G_gisinit(argv[0]);
+
+	module = G_define_module();
+	module->keywords = _("water, detection");
+	module->description = _("Water detection, 1 if found, 0 if not");
+
+	/* Define the different options */
+	input1 = G_define_standard_option(G_OPT_R_INPUT) ;
+	input1->key	   = _("ndvi");
+	input1->description=_("Name of the NDVI map [-]");
+	input1->answer     =_("ndvi");
+
+	input2 = G_define_standard_option(G_OPT_R_INPUT) ;
+	input2->key        =_("albedo");
+	input2->required   = NO ;
+	input2->description=_("Name of the Albedo map [-]");
+
+	input3 = G_define_standard_option(G_OPT_R_INPUT) ;
+	input3->key        =_("Modref7");
+	input3->required   = NO ;
+	input3->description=_("Name of the Modis surface reflectance band 7 [-]");
+
+	output1 = G_define_standard_option(G_OPT_R_OUTPUT) ;
+	output1->key        =_("water");
+	output1->description=_("Name of the output water layer [0/1]");
+	output1->answer     =_("water");
+	/********************/
+	if (G_parser(argc, argv))
+		exit (EXIT_FAILURE);
+
+	ndvi	 	= input1->answer;
+	if(input2->answer)
+	albedo	 	= input2->answer;
+	if(input3->answer)
+	ref7		= input3->answer;
+	
+
+	if(!input2->answer&&!input3->answer){
+		G_fatal_error(_("ERROR: needs either Albedo or Modis surface reflectance in band 7, bailing out."));
+	}
+	
+	result1  = output1->answer;
+	/***************************************************/
+	mapset = G_find_cell2(ndvi, "");
+	if (mapset == NULL) {
+		G_fatal_error(_("cell file [%s] not found"), ndvi);
+	}
+	data_type_ndvi = G_raster_map_type(ndvi,mapset);
+	if ( (infd_ndvi = G_open_cell_old (ndvi,mapset)) < 0)
+		G_fatal_error (_("Cannot open cell file [%s]"), ndvi);
+	if (G_get_cellhd (ndvi, mapset, &cellhd) < 0)
+		G_fatal_error (_("Cannot read file header of [%s])"), ndvi);
+	inrast_ndvi = G_allocate_raster_buf(data_type_ndvi);
+	/***************************************************/
+	if(input2->answer){
+	mapset = G_find_cell2 (albedo, "");
+	if (mapset == NULL) {
+		G_fatal_error(_("cell file [%s] not found"),albedo);
+	}
+	data_type_albedo = G_raster_map_type(albedo,mapset);
+	if ( (infd_albedo = G_open_cell_old (albedo,mapset)) < 0)
+		G_fatal_error(_("Cannot open cell file [%s]"), albedo);
+	if (G_get_cellhd (albedo, mapset, &cellhd) < 0)
+		G_fatal_error(_("Cannot read file header of [%s]"), albedo);
+	inrast_albedo = G_allocate_raster_buf(data_type_albedo);
+	}
+	/***************************************************/
+	if(input3->answer){
+	mapset = G_find_cell2 (ref7, "");
+	if (mapset == NULL) {
+		G_fatal_error(_("Cell file [%s] not found"), ref7);
+	}
+	data_type_ref7 = G_raster_map_type(ref7,mapset);
+	if ( (infd_ref7 = G_open_cell_old (ref7,mapset)) < 0)
+		G_fatal_error(_("Cannot open cell file [%s]"), ref7);
+	if (G_get_cellhd (ref7, mapset, &cellhd) < 0)
+		G_fatal_error(_("Cannot read file header of [%s]"), ref7);
+	inrast_ref7 = G_allocate_raster_buf(data_type_ref7);
+	}
+	/***************************************************/
+	G_debug(3, "number of rows %d",cellhd.rows);
+	nrows = G_window_rows();
+	ncols = G_window_cols();
+	outrast1 = G_allocate_raster_buf(data_type_output);
+	/* Create New raster files */
+	if ( (outfd1 = G_open_raster_new (result1,data_type_output)) < 0)
+		G_fatal_error(_("Could not open <%s>"),result1);
+	/* Process pixels */
+	for (row = 0; row < nrows; row++)
+	{
+		CELL d;
+		DCELL d_ndvi;
+		DCELL d_albedo;
+		DCELL d_ref7;
+		G_percent(row,nrows,2);
+		/* read input maps */	
+		if(G_get_raster_row(infd_ndvi,inrast_ndvi,row,data_type_ndvi)<0)
+			G_fatal_error(_("Could not read from <%s>"),ndvi);
+		if(input2->answer){
+		if(G_get_raster_row(infd_albedo,inrast_albedo,row,data_type_albedo)<0)
+			G_fatal_error(_("Could not read from <%s>"),albedo);
+		}
+		if(input3->answer){
+		if(G_get_raster_row(infd_ref7,inrast_ref7,row,data_type_ref7)<0)
+			G_fatal_error(_("Could not read from <%s>"),ref7);
+		}
+		/*process the data */
+		for (col=0; col < ncols; col++)
+		{
+			switch(data_type_ndvi){
+				case CELL_TYPE:
+					d_ndvi = (double) ((CELL *) inrast_ndvi)[col];
+					break;
+				case FCELL_TYPE:
+					d_ndvi = (double) ((FCELL *) inrast_ndvi)[col];
+					break;
+				case DCELL_TYPE:
+					d_ndvi = ((DCELL *) inrast_ndvi)[col];
+					break;
+			}
+			if(input2->answer){
+			switch(data_type_albedo){
+				case CELL_TYPE:
+					d_albedo = (double) ((CELL *) inrast_albedo)[col];
+					break;
+				case FCELL_TYPE:
+					d_albedo = (double) ((FCELL *) inrast_albedo)[col];
+					break;
+				case DCELL_TYPE:
+					d_albedo = ((DCELL *) inrast_albedo)[col];
+					break;
+			}
+			}
+			if(input3->answer){
+			switch(data_type_ref7){
+				case CELL_TYPE:
+					d_ref7 = (double) ((CELL *) inrast_ref7)[col];
+					break;
+				case FCELL_TYPE:
+					d_ref7 = (double) ((FCELL *) inrast_ref7)[col];
+					break;
+				case DCELL_TYPE:
+					d_ref7 = ((DCELL *) inrast_ref7)[col];
+					break;
+			}
+			}
+			if(G_is_d_null_value(&d_ndvi)||
+			(input2->answer&&G_is_d_null_value(&d_albedo))||
+			(input3->answer&&G_is_d_null_value(&d_ref7))){
+				G_set_c_null_value(&outrast1[col],1);
+			}else {
+				/************************************/
+				/* calculate water detection	    */
+				if(input2->answer){
+					d = water(d_albedo,d_ndvi);
+				} else if (input3->answer){
+					d = water_modis(d_ref7,d_ndvi);
+				}
+				outrast1[col] = d;
+			}
+		}
+		if (G_put_raster_row (outfd1, outrast1, data_type_output) < 0)
+			G_fatal_error(_("Cannot write to output raster file"));
+	}
+
+	G_free (inrast_ndvi);
+	G_close_cell (infd_ndvi);
+	if(input2->answer){
+		G_free (inrast_albedo);
+		G_close_cell (infd_albedo);
+	}
+	if(input3->answer){
+		G_free (inrast_ref7);
+		G_close_cell (infd_ref7);
+	}
+	G_free (outrast1);
+	G_close_cell (outfd1);
+
+	G_short_history(result1, "raster", &history);
+	G_command_history(&history);
+	G_write_history(result1,&history);
+
+	exit(EXIT_SUCCESS);
+}
+

Added: grass-addons/gipe/i.water/water.c
===================================================================
--- grass-addons/gipe/i.water/water.c	                        (rev 0)
+++ grass-addons/gipe/i.water/water.c	2008-06-27 01:43:43 UTC (rev 31857)
@@ -0,0 +1,13 @@
+#include<stdio.h>
+
+int water(double albedo, double ndvi)
+{
+	double result;
+	if (albedo<0.1&&ndvi<0.1){
+		result = 1 ;
+	} else {
+		result = 0 ;
+	}
+	return result;
+}
+

Added: grass-addons/gipe/i.water/water_modis.c
===================================================================
--- grass-addons/gipe/i.water/water_modis.c	                        (rev 0)
+++ grass-addons/gipe/i.water/water_modis.c	2008-06-27 01:43:43 UTC (rev 31857)
@@ -0,0 +1,13 @@
+#include<stdio.h>
+
+int water_modis(double surf_ref_7, double ndvi)
+{
+	double result;
+	if (surf_ref_7<0.04&&ndvi<0.1){
+		result = 1 ;
+	} else {
+		result = 0 ;
+	}
+	return result;
+}
+



More information about the grass-commit mailing list