[GRASS-SVN] r38036 - in grass/trunk: gui/wxpython/xml imagery
imagery/i.biomass
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Jun 23 09:50:16 EDT 2009
Author: ychemin
Date: 2009-06-23 09:50:15 -0400 (Tue, 23 Jun 2009)
New Revision: 38036
Added:
grass/trunk/imagery/i.biomass/
grass/trunk/imagery/i.biomass/Makefile
grass/trunk/imagery/i.biomass/biomass.c
grass/trunk/imagery/i.biomass/i.biomass.html
grass/trunk/imagery/i.biomass/main.c
grass/trunk/imagery/i.biomass/solar_day.c
Modified:
grass/trunk/gui/wxpython/xml/menudata.xml
grass/trunk/imagery/Makefile
Log:
Added biomass growth module
Modified: grass/trunk/gui/wxpython/xml/menudata.xml
===================================================================
--- grass/trunk/gui/wxpython/xml/menudata.xml 2009-06-23 11:46:47 UTC (rev 38035)
+++ grass/trunk/gui/wxpython/xml/menudata.xml 2009-06-23 13:50:15 UTC (rev 38036)
@@ -2258,6 +2258,13 @@
<handler>self.OnMenuCmd</handler>
<command>i.vi</command>
</menuitem>
+ <menuitem>
+ <label>Biomass Growth</label>
+ <help>Calculates the growth of biomass for a day.</help>
+ <keywords>imagery,biomass growth,vegetation index,biophysical parameters</keywords>
+ <handler>self.OnMenuCmd</handler>
+ <command>i.biomass</command>
+ </menuitem>
<separator />
<menuitem>
<label>Emissivity</label>
Modified: grass/trunk/imagery/Makefile
===================================================================
--- grass/trunk/imagery/Makefile 2009-06-23 11:46:47 UTC (rev 38035)
+++ grass/trunk/imagery/Makefile 2009-06-23 13:50:15 UTC (rev 38036)
@@ -3,6 +3,7 @@
SUBDIRS = \
i.albedo \
i.atcorr \
+ i.biomass \
i.cca \
i.cluster \
i.eb.evapfr \
Added: grass/trunk/imagery/i.biomass/Makefile
===================================================================
--- grass/trunk/imagery/i.biomass/Makefile (rev 0)
+++ grass/trunk/imagery/i.biomass/Makefile 2009-06-23 13:50:15 UTC (rev 38036)
@@ -0,0 +1,14 @@
+MODULE_TOPDIR = ../..
+
+PGM = i.biomass
+
+LIBES = $(GISLIB) $(RASTERLIB)
+DEPENDENCIES = $(GISDEP) $(RASTERDEP)
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+ifneq ($(USE_LARGEFILES),)
+ EXTRA_CFLAGS = -D_FILE_OFFSET_BITS=64
+endif
+
+default: cmd
Added: grass/trunk/imagery/i.biomass/biomass.c
===================================================================
--- grass/trunk/imagery/i.biomass/biomass.c (rev 0)
+++ grass/trunk/imagery/i.biomass/biomass.c 2009-06-23 13:50:15 UTC (rev 38036)
@@ -0,0 +1,21 @@
+#include<stdio.h>
+#include<math.h>
+#include<stdlib.h>
+ double biomass(double fpar, double solar_day, double evap_fr,
+ double light_use_ef)
+{
+ double result, apar, conversion_coefs;
+
+
+ /*fPAR in Bastiaassen and Ali = a (= 1.257) * NDVI + b (= -0.161); */
+ /*light use efficiency is variable for each crop! */
+ /* Cotton light_use_ef in Uzbekistan is 1.9 */
+
+ /* apar = (b+a*ndvi); */
+ conversion_coefs = 0.0864 * 10;
+ apar = fpar * (0.48 * solar_day);
+ result = apar * light_use_ef * evap_fr * conversion_coefs;
+ return result;
+}
+
+
Added: grass/trunk/imagery/i.biomass/i.biomass.html
===================================================================
--- grass/trunk/imagery/i.biomass/i.biomass.html (rev 0)
+++ grass/trunk/imagery/i.biomass/i.biomass.html 2009-06-23 13:50:15 UTC (rev 38036)
@@ -0,0 +1,40 @@
+<H2>DESCRIPTION</H2>
+
+<EM>i.biomass</EM> Calculates the biomass growth for a day after [1][2].
+
+Input:
+<ul>
+ <li>fPAR, the modified Photosynthetic Active Radiation for crops.
+ <li>Light Use Efficiency [0.0-1.0], in Uzbekistan cotton is at 1.9 most of the time.
+ <li>Latitude [0.0-90.0], from i.latlong.
+ <li>DOY [1-366].
+ <li>Transmissivity of the Atmosphere single-way [0.0-1.0], mostly around 0.7+ in clear sky.
+ <li>Water availability [0.0-1.0], possibly using direct output from i.eb.evapfr.
+</ul>
+
+<H2>NOTES</H2>
+It can use the output of i.eb.evapfr directly as water availability input.
+
+<H2>TODO</H2>
+remove Latitude, DOY and Tsw from input and replace with a raster input compatible with r.sun output.
+
+<H2>SEE ALSO</H2>
+
+<em>
+<A HREF="i.eb.evapfr">i.eb.evapfr</A><br>
+</em>
+
+<H2>REFERENCES</H2>
+
+ <p>[1] Bastiaanssen, W.G.M., Ali, S., 2002. A new crop yield forecasting model based on satellite measurements applied across the Indus Basin, Pakistan. Agriculture, Ecosystems and Environment, 94(3):321-340.
+
+ <p>[2] Chemin, Y., Platonov, A., Abdullaev, I., Ul-Hassan, M. 2005. Supplementing farm level water productivity assessment by remote sensing in transition economies. Water International. 30(4):513-521.
+
+
+<H2>AUTHORS</H2>
+
+Yann Chemin, Bec de Mortagne, France<BR>
+
+
+<p>
+<i>Last changed: $Date:$</i>
Added: grass/trunk/imagery/i.biomass/main.c
===================================================================
--- grass/trunk/imagery/i.biomass/main.c (rev 0)
+++ grass/trunk/imagery/i.biomass/main.c 2009-06-23 13:50:15 UTC (rev 38036)
@@ -0,0 +1,203 @@
+
+/****************************************************************************
+ *
+ * MODULE: i.biomass
+ * AUTHOR(S): Yann Chemin - yann.chemin at gmail.com
+ * PURPOSE: Creates a map of biomass growth
+ *
+ *
+ * COPYRIGHT: (C) 2002-2009 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/raster.h>
+#include <grass/glocale.h>
+
+double solar_day(double lat, double doy, double tsw);
+double biomass(double apar, double solar_day, double evap_fr,
+ double light_use_ef);
+
+int main(int argc, char *argv[])
+{
+ int nrows, ncols;
+ int row, col;
+ struct GModule *module;
+ struct Option *input1, *input2, *input3, *input4, *input5, *input6;
+ struct Option *output1;
+ struct History history; /*metadata */
+ struct Colors colors;
+ char *result1; /*output raster name */
+ int infd_fpar, infd_luf, infd_lat;
+ int infd_doy, infd_tsw, infd_wa;
+ int outfd1;
+ char *fpar, *luf, *lat, *doy, *tsw, *wa;
+ void *inrast_fpar, *inrast_luf, *inrast_lat;
+ void *inrast_doy, *inrast_tsw, *inrast_wa;
+ DCELL * outrast1;
+ /************************************/
+ G_gisinit(argv[0]);
+ module = G_define_module();
+ module->keywords = _("biomass, fpar, yield");
+ module->description =
+ _("biomass growth calculation, precursor of crop yield calculation");
+
+ /* Define the different options */
+ input1 = G_define_standard_option(G_OPT_R_INPUT);
+ input1->key = "fpar";
+ input1->description = _("Name of the fPAR map");
+
+ input2 = G_define_standard_option(G_OPT_R_INPUT);
+ input2->key = "luf";
+ input2->description =
+ _("Name of the light use efficiency map (UZB:cotton=1.9)");
+
+ input3 = G_define_standard_option(G_OPT_R_INPUT);
+ input3->key = "lat";
+ input3->description = _("Name of the degree latitude map [dd.ddd]");
+
+ input4 = G_define_standard_option(G_OPT_R_INPUT);
+ input4->key = "doy";
+ input4->description = _("Name of the Day of Year map [1-366]");
+
+ input5 = G_define_standard_option(G_OPT_R_INPUT);
+ input5->key = "tsw";
+ input5->description =
+ _("Name of the single-way transmissivity map [0.0-1.0]");
+
+ input6 = G_define_standard_option(G_OPT_R_INPUT);
+ input6->key = "wa";
+ input6->description =
+ _("Value of the water availability map [0.0-1.0]");
+
+ output1 = G_define_standard_option(G_OPT_R_INPUT);
+ output1->description =
+ _("Name of the output daily biomass growth layer [kg/ha/d]");
+
+ /********************/
+ if (G_parser(argc, argv))
+ exit(EXIT_FAILURE);
+ fpar = input1->answer;
+ luf = input2->answer;
+ lat = input3->answer;
+ doy = input4->answer;
+ tsw = input5->answer;
+ wa = input6->answer;
+ result1 = output1->answer;
+
+ /***************************************************/
+ if ((infd_fpar = Rast_open_cell_old(fpar, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), fpar);
+ inrast_fpar = Rast_allocate_d_raster_buf();
+
+ if ((infd_luf = Rast_open_cell_old(luf, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), luf);
+ inrast_luf = Rast_allocate_d_raster_buf();
+
+ if ((infd_lat = Rast_open_cell_old(lat, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), lat);
+ inrast_lat = Rast_allocate_d_raster_buf();
+
+ if ((infd_doy = Rast_open_cell_old(doy, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), doy);
+ inrast_doy = Rast_allocate_d_raster_buf();
+
+ if ((infd_tsw = Rast_open_cell_old(tsw, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), tsw);
+ inrast_tsw = Rast_allocate_d_raster_buf();
+
+ if ((infd_wa = Rast_open_cell_old(wa, "")) < 0)
+ G_fatal_error(_("Unable to open raster map <%s>"), wa);
+ inrast_wa = Rast_allocate_d_raster_buf();
+
+ nrows = G_window_rows();
+ ncols = G_window_cols();
+ outrast1 = Rast_allocate_d_raster_buf();
+
+ /* Create New raster files */
+ if ((outfd1 = Rast_open_raster_new(result1, DCELL_TYPE)) < 0)
+ G_fatal_error(_("Unable to create raster map <%s>"), result1);
+
+ /* Process pixels */
+ for (row = 0; row < nrows; row++)
+ {
+ DCELL d;
+ DCELL d_fpar;
+ DCELL d_luf;
+ DCELL d_lat;
+ DCELL d_doy;
+ DCELL d_tsw;
+ DCELL d_solar;
+ DCELL d_wa;
+ G_percent(row, nrows, 2);
+
+ /* read input maps */
+ if (Rast_get_d_raster_row(infd_fpar, inrast_fpar, row)< 0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),fpar,row);
+ if (Rast_get_d_raster_row(infd_luf,inrast_luf,row)<0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),luf,row);
+ if (Rast_get_d_raster_row(infd_lat, inrast_lat, row) < 0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),lat,row);
+ if (Rast_get_d_raster_row(infd_doy, inrast_doy, row) < 0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),doy,row);
+ if (Rast_get_d_raster_row(infd_tsw, inrast_tsw, row) < 0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),tsw,row);
+ if (Rast_get_d_raster_row(infd_wa,inrast_wa,row)<0)
+ G_fatal_error(_("Unable to read from <%s> row %d"),wa,row);
+ /*process the data */
+ for (col = 0; col < ncols; col++)
+ {
+ d_fpar = ((DCELL *) inrast_fpar)[col];
+ d_luf = ((DCELL *) inrast_luf)[col];
+ d_lat = ((DCELL *) inrast_lat)[col];
+ d_doy = ((DCELL *) inrast_doy)[col];
+ d_tsw = ((DCELL *) inrast_tsw)[col];
+ d_wa = ((DCELL *) inrast_wa)[col];
+ if (Rast_is_d_null_value(&d_fpar) ||
+ Rast_is_d_null_value(&d_luf) ||
+ Rast_is_d_null_value(&d_lat) ||
+ Rast_is_d_null_value(&d_doy) ||
+ Rast_is_d_null_value(&d_tsw) ||
+ Rast_is_d_null_value(&d_wa))
+ Rast_set_d_null_value(&outrast1[col], 1);
+ else {
+ d_solar = solar_day(d_lat, d_doy, d_tsw);
+ d = biomass(d_fpar,d_solar,d_wa,d_luf);
+ outrast1[col] = d;
+ }
+ }
+ if (Rast_put_d_raster_row(outfd1, outrast1) < 0)
+ G_fatal_error(_("Failed writing raster map <%s>"), result1);
+ }
+
+ /* Color table for biomass */
+ Rast_init_colors(&colors);
+ Rast_add_color_rule(0, 0, 0, 0, 1, 255, 255, 255, &colors);
+ G_free(inrast_fpar);
+ G_free(inrast_luf);
+ G_free(inrast_lat);
+ G_free(inrast_doy);
+ G_free(inrast_tsw);
+ G_free(inrast_wa);
+ G_free(outrast1);
+ Rast_close_cell(infd_fpar);
+ Rast_close_cell(infd_luf);
+ Rast_close_cell(infd_lat);
+ Rast_close_cell(infd_doy);
+ Rast_close_cell(infd_tsw);
+ Rast_close_cell(infd_wa);
+ Rast_close_cell(outfd1);
+ Rast_short_history(result1, "raster", &history);
+ Rast_command_history(&history);
+ Rast_write_history(result1, &history);
+ exit(EXIT_SUCCESS);
+}
+
+
Added: grass/trunk/imagery/i.biomass/solar_day.c
===================================================================
--- grass/trunk/imagery/i.biomass/solar_day.c (rev 0)
+++ grass/trunk/imagery/i.biomass/solar_day.c 2009-06-23 13:50:15 UTC (rev 38036)
@@ -0,0 +1,25 @@
+#include<stdio.h>
+#include<math.h>
+#include<stdlib.h>
+
+ /*
+ * Average Solar Diurnal Radiation after Bastiaanssen (1995)
+ */
+
+#define PI 3.1415927
+double solar_day(double lat, double doy, double tsw)
+{
+ double ws, cosun, n10_temp, delta, ds, result;
+
+ ds = 1.0 + 0.01672 * sin(2 * PI * (doy - 93.5) / 365.0);
+ delta = 0.4093 * sin((2 * PI * doy / 365) - 1.39);
+ n10_temp = lat * PI / 180.0;
+ ws = acos(-tan(n10_temp) * tan(delta * PI / 180.0));
+ cosun =
+ ws * sin(delta * PI / 180.0) * sin(n10_temp) +
+ cos(delta * PI / 180.0) * cos(n10_temp) * sin(ws);
+ result = (cosun * 1367 * tsw) / (PI * ds * ds);
+ return result;
+}
+
+
More information about the grass-commit
mailing list