[GRASS-SVN] r30177 - grass/trunk/raster/simwe/r.sim.water
svn_grass at osgeo.org
svn_grass at osgeo.org
Sat Feb 16 00:34:17 EST 2008
Author: helena
Date: 2008-02-16 00:34:16 -0500 (Sat, 16 Feb 2008)
New Revision: 30177
Modified:
grass/trunk/raster/simwe/r.sim.water/main.c
Log:
improved description, option for rain and mannings as values added, non-functioning flag -m removed
Modified: grass/trunk/raster/simwe/r.sim.water/main.c
===================================================================
--- grass/trunk/raster/simwe/r.sim.water/main.c 2008-02-15 22:42:03 UTC (rev 30176)
+++ grass/trunk/raster/simwe/r.sim.water/main.c 2008-02-16 05:34:16 UTC (rev 30177)
@@ -1,6 +1,6 @@
/****************************************************************************
*
- * MODULE: r.simwe.hydro: main program for hydrologic and sediment transport
+ * MODULE: r.sim.water: main program for hydrologic and sediment transport
* simulation (SIMWE)
*
* AUTHOR(S): L. Mitas, H. Mitasova, J. Hofierka
@@ -15,7 +15,7 @@
*****************************************************************************/
/*-
- * r.simwe.hydro: main program for hydrologic and sediment transport
+ * r.sim.water: main program for hydrologic and sediment transport
* simulation (SIMWE)
*
* Original program (2002) and various modifications:
@@ -46,15 +46,23 @@
*
*/
+/********************************/
+/* DEFINE GLOB VAR */
+/********************************/
#define NWALK "2000000"
#define DIFFC "0.8"
-#define HMAX "0.4"
+#define HMAX "0.3"
#define HALPHA "4.0"
#define HBETA "0.5"
-#define NITER "1200"
-#define ITEROUT "300"
+#define NITER "10"
+#define ITEROUT "2"
#define DENSITY "200"
+#define RAINVAL "30"
+#define MANINVAL "0.1"
+/********************************/
+/* INCLUDES */
+/********************************/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
@@ -70,6 +78,9 @@
#include <grass/site.h>
#include <grass/glocale.h>
+/********************************/
+/* Specific stuff */
+/********************************/
#define MAIN
#include <grass/waterglobs.h>
@@ -83,352 +94,379 @@
char msg[1024];
+/****************************************/
+/* MAIN */
+/****************************************/
int main ( int argc, char *argv[])
{
- int i,ii,j,l;
- int ret_val;
- double x_orig, y_orig;
- static int rand1 = 12345;
- static int rand2 = 67891;
+ int i,ii,j,l;
+ int ret_val;
+ double x_orig, y_orig;
+ static int rand1 = 12345;
+ static int rand2 = 67891;
+
+ G_gisinit (argv[0]);
+
+ module = G_define_module();
+ module->keywords = _("raster");
+ module->description =
+ _("Overland flow hydrologic simulation using "
+ "path sampling method (SIMWE)");
+
+ if (G_get_set_window (&cellhd) == -1)
+ exit (EXIT_FAILURE);
- G_gisinit (argv[0]);
+ conv = G_database_units_to_meters_factor();
+
+ mixx = conv * cellhd.west;
+ maxx = conv * cellhd.east;
+ miyy = conv * cellhd.south;
+ mayy = conv * cellhd.north;
+
+ stepx = cellhd.ew_res * conv;
+ stepy = cellhd.ns_res * conv;
+ /* step = amin1(stepx,stepy);*/
+ step = (stepx + stepy) / 2.;
+ mx = cellhd.cols;
+ my = cellhd.rows;
+ x_orig = cellhd.west * conv;
+ y_orig = cellhd.south * conv;/* do we need this? */
+ xmin = 0.;
+ ymin = 0.;
+ xp0 = xmin + stepx / 2.;
+ yp0 = ymin + stepy / 2.;
+ xmax = xmin + stepx * (float) mx;
+ ymax = ymin + stepy * (float) my;
+ /*
+ bxmi=2093113. * conv;
+ bymi=731331. * conv;
+ bxma=2093461. * conv;
+ byma=731529. * conv;
+ bresx=2. * conv;
+ bresy=2. * conv;
+ maxwab=100000;
+
+ mx2o= (int)((bxma-bxmi)/bresx);
+ my2o= (int)((byma-bymi)/bresy);
+ */
+ /* relative small box coordinates: leave 1 grid layer for overlap */
+ /*
+ bxmi = bxmi - mixx + stepx;
+ bymi = bymi - miyy + stepy;
+ bxma = bxma - mixx - stepx;
+ byma = byma - miyy - stepy;
+ mx2 = mx2o - 2*((int) (stepx / bresx));
+ my2 = my2o - 2*((int) (stepy / bresy));
+ */
+ parm.elevin = G_define_option();
+ parm.elevin->key = "elevin";
+ parm.elevin->type = TYPE_STRING;
+ parm.elevin->required = YES;
+ parm.elevin->gisprompt = "old,cell,raster";
+ parm.elevin->description = _("Name of the elevation raster map, same units as N,E: [m] or [ft] ");
+ parm.elevin->guisection = _("Input_options");
+
+ parm.dxin = G_define_option();
+ parm.dxin->key = "dxin";
+ parm.dxin->type = TYPE_STRING;
+ parm.dxin->required = YES;
+ parm.dxin->gisprompt = "old,cell,raster";
+ parm.dxin->description = _("Name of the x-derivatives raster map [m/m] or [ft/ft]");
+ parm.dxin->guisection = _("Input_options");
+
+ parm.dyin = G_define_option();
+ parm.dyin->key = "dyin";
+ parm.dyin->type = TYPE_STRING;
+ parm.dyin->required = YES;
+ parm.dyin->gisprompt = "old,cell,raster";
+ parm.dyin->description = _("Name of the y-derivatives raster map [m/m] or [ft/ft]");
+ parm.dyin->guisection = _("Input_options");
+
+ parm.rain = G_define_option();
+ parm.rain->key = "rain";
+ parm.rain->type = TYPE_STRING;
+ parm.rain->required = NO;
+ parm.rain->gisprompt = "old,cell,raster";
+ parm.rain->description = _("Name of the rainfall excess rate (rain-infilt) raster map [mm/hr]");
+ parm.rain->guisection = _("Input_options");
+
+ parm.rainval = G_define_option();
+ parm.rainval->key = "rain_val";
+ parm.rainval->type = TYPE_DOUBLE;
+ parm.rainval->answer = RAINVAL;
+ parm.rainval->required = NO;
+ parm.rainval->description = _("Rainfall excess rate unique value [mm/hr]");
+ parm.rainval->guisection = _("Input_options");
+
+ parm.infil = G_define_option();
+ parm.infil->key = "infil";
+ parm.infil->type = TYPE_STRING;
+ parm.infil->required = YES;
+ parm.infil->gisprompt = "old,cell,raster";
+ parm.infil->description = _("Name of the runoff infiltration rate raster map [mm/hr]");
+ parm.infil->guisection = _("Input_options");
+
+ parm.traps = G_define_option();
+ parm.traps->key = "traps";
+ parm.traps->type = TYPE_STRING;
+ parm.traps->required = NO;
+ parm.traps->gisprompt = "old,cell,raster";
+ parm.traps->description = _("Name of the flow controls raster map (permeability ratio 0-1)");
+ parm.traps->guisection = _("Input_options");
+
+ parm.manin = G_define_option();
+ parm.manin->key = "manin";
+ parm.manin->type = TYPE_STRING;
+ parm.manin->required = NO;
+ parm.manin->gisprompt = "old,cell,raster";
+ parm.manin->description = _("Name of the Mannings n raster map");
+ parm.manin->guisection = _("Input_options");
+
+ parm.maninval = G_define_option();
+ parm.maninval->key = "manin_val";
+ parm.maninval->type = TYPE_DOUBLE;
+ parm.maninval->answer = MANINVAL;
+ parm.maninval->required = NO;
+ parm.maninval->description = _("Mannings n unique value");
+ parm.maninval->guisection = _("Input_options");
+
+ /* needs to be updated to GRASS 6 vector format !! */
+ parm.sfile = G_define_option ();
+ parm.sfile->key = "vector";
+ parm.sfile->type = TYPE_STRING;
+ parm.sfile->required = NO;
+ parm.sfile->gisprompt = "old,vector,vector";
+ parm.sfile->description = _("Name of the sampling locations vector points map");
+ parm.sfile->guisection = _("Input_options");
+
+ parm.depth = G_define_option();
+ parm.depth->key = "depth";
+ parm.depth->type = TYPE_STRING;
+ parm.depth->required = NO;
+ parm.depth->gisprompt = "new,cell,raster";
+ parm.depth->description = _("Output water depth raster map [m]");
+ parm.depth->guisection = _("Output_options");
+
+ parm.disch = G_define_option();
+ parm.disch->key = "disch";
+ parm.disch->type = TYPE_STRING;
+ parm.disch->required = NO;
+ parm.disch->gisprompt = "new,cell,raster";
+ parm.disch->description = _("Output water discharge raster map [m3/s]");
+ parm.disch->guisection = _("Output_options");
+
+ parm.err = G_define_option();
+ parm.err->key = "err";
+ parm.err->type = TYPE_STRING;
+ parm.err->required = NO;
+ parm.err->gisprompt = "new,cell,raster";
+ parm.err->description = _("Output simulation error raster map [m]");
+ parm.err->guisection = _("Output_options");
+
+ parm.outwalk = G_define_option ();
+ parm.outwalk->key = "outwalk";
+ parm.outwalk->type = TYPE_STRING;
+ parm.outwalk->required = NO;
+ parm.outwalk->gisprompt = "new,vector,vector";
+ parm.outwalk->description = _("Name of the output walkers vector points map");
+ parm.outwalk->guisection = _("Output_options");
+
+ parm.nwalk = G_define_option();
+ parm.nwalk->key = "nwalk";
+ parm.nwalk->type = TYPE_INTEGER;
+ parm.nwalk->answer = NWALK;
+ parm.nwalk->required = NO;
+ parm.nwalk->description = _("Number of walkers");
+ parm.nwalk->guisection = _("Parameters");
+
+ parm.niter = G_define_option();
+ parm.niter->key = "niter";
+ parm.niter->type = TYPE_INTEGER;
+ parm.niter->answer = NITER;
+ parm.niter->required = NO;
+ parm.niter->description = _("Time used for iterations [minutes]");
+ parm.niter->guisection = _("Parameters");
+
+ parm.outiter = G_define_option();
+ parm.outiter->key = "outiter";
+ parm.outiter->type = TYPE_INTEGER;
+ parm.outiter->answer = ITEROUT;
+ parm.outiter->required = NO;
+ parm.outiter->description = _("Time interval for creating output maps [minutes]");
+ parm.outiter->guisection = _("Parameters");
+
+ parm.density = G_define_option();
+ parm.density->key = "density";
+ parm.density->type = TYPE_INTEGER;
+ parm.density->answer = DENSITY;
+ parm.density->required = NO;
+ parm.density->description = _("Density of output walkers");
+ parm.density->guisection = _("Parameters");
+
+ parm.diffc = G_define_option();
+ parm.diffc->key = "diffc";
+ parm.diffc->type = TYPE_DOUBLE;
+ parm.diffc->answer = DIFFC;
+ parm.diffc->required = NO;
+ parm.diffc->description = _("Water diffusion constant");
+ parm.diffc->guisection = _("Parameters");
+
+ parm.hmax = G_define_option();
+ parm.hmax->key = "hmax";
+ parm.hmax->type = TYPE_DOUBLE;
+ parm.hmax->answer = HMAX;
+ parm.hmax->required = NO;
+ parm.hmax->description = _("Threshold water depth [m] (diffusion increases after this water depth is reached)");
+ parm.hmax->guisection = _("Parameters");
+
+ parm.halpha = G_define_option();
+ parm.halpha->key = "halpha";
+ parm.halpha->type = TYPE_DOUBLE;
+ parm.halpha->answer = HALPHA;
+ parm.halpha->required = NO;
+ parm.halpha->description = _("Diffusion increase constant");
+ parm.halpha->guisection = _("Parameters");
+
+ parm.hbeta = G_define_option();
+ parm.hbeta->key = "hbeta";
+ parm.hbeta->type = TYPE_DOUBLE;
+ parm.hbeta->answer = HBETA;
+ parm.hbeta->required = NO;
+ parm.hbeta->description = _("Weighting factor for water flow velocity vector");
+ parm.hbeta->guisection = _("Parameters");
+ /*
+ flag.mscale = G_define_flag ();
+ flag.mscale->key = 'm';
+ flag.mscale->description = _("Multiscale simulation");
+ */
+ flag.tserie = G_define_flag ();
+ flag.tserie->key = 't';
+ flag.tserie->description = _("Time-series output");
- module = G_define_module();
- module->keywords = _("raster");
- module->description =
- _("Overland flow hydrologic model based on duality "
- "particle-field concept (SIMWE)");
-
- if (G_get_set_window (&cellhd) == -1)
- exit (EXIT_FAILURE);
- conv = G_database_units_to_meters_factor();
+ if (G_parser (argc, argv))
+ exit (EXIT_FAILURE);
- mixx = conv * cellhd.west;
- maxx = conv * cellhd.east;
- miyy = conv * cellhd.south;
- mayy = conv * cellhd.north;
-
- stepx = cellhd.ew_res * conv;
- stepy = cellhd.ns_res * conv;
-/* step = amin1(stepx,stepy);*/
- step = (stepx + stepy) / 2.;
- mx = cellhd.cols;
- my = cellhd.rows;
- x_orig = cellhd.west * conv;
- y_orig = cellhd.south * conv;/* do we need this? */
- xmin = 0.;
- ymin = 0.;
- xp0 = xmin + stepx / 2.;
- yp0 = ymin + stepy / 2.;
- xmax = xmin + stepx * (float) mx;
- ymax = ymin + stepy * (float) my;
+ /* mscale=flag.mscale->answer;*/
+ ts=flag.tserie->answer;
- bxmi=2093113. * conv;
- bymi=731331. * conv;
- bxma=2093461. * conv;
- byma=731529. * conv;
- bresx=2. * conv;
- bresy=2. * conv;
- maxwab=100000;
+ elevin = parm.elevin->answer;
+ dxin = parm.dxin->answer;
+ dyin = parm.dyin->answer;
+ /* maskmap = parm.maskmap->answer;*/
+ rain = parm.rain->answer;
+ infil = parm.infil->answer;
+ traps = parm.traps->answer;
+ manin = parm.manin->answer;
+ depth = parm.depth->answer;
+ disch = parm.disch->answer;
+ err = parm.err->answer;
+ outwalk = parm.outwalk->answer;
+ sfile = parm.sfile->answer;
+
+ sscanf(parm.nwalk->answer, "%d", &maxwa);
+ sscanf(parm.niter->answer, "%d", ×ec);
+ sscanf(parm.outiter->answer, "%d", &iterout);
+ sscanf(parm.density->answer, "%d", &ldemo);
+ sscanf(parm.diffc->answer, "%lf", &frac);
+ sscanf(parm.hmax->answer, "%lf", &hhmax);
+ sscanf(parm.halpha->answer, "%lf", &halpha);
+ sscanf(parm.hbeta->answer, "%lf", &hbeta);
+ /*Check for Rain Unique Value Input*/
+ if(parm.rainval->answer==NULL){
+ rain_val = -999.99;
+ printf("rain_val is set to: %f\n",rain_val);
+ } else {
+ sscanf(parm.rainval->answer, "%lf", &rain_val);
+ }
+ /*Check for Manin Unique Value Input*/
+ if(parm.maninval->answer==NULL){
+ manin_val = -999.99;
+ printf("manin_val is set to: %f\n",manin_val);
+ } else {
+ sscanf(parm.maninval->answer, "%lf", &manin_val);
+ }
- mx2o= (int)((bxma-bxmi)/bresx);
- my2o= (int)((byma-bymi)/bresy);
+ /* Recompute timesec from user input in minutes
+ * to real timesec in seconds */
+ timesec = timesec * 60.0;
+ iterout = iterout * 60.0;
+ if ((timesec/iterout) > 100.0)
+ G_message(_("More than 100 files are going to be created !!!!!"));
+ /* compute how big the raster is and set this to appr 2 walkers per cell */
+ rwalk = (double) maxwa;
-/* relative small box coordinates: leave 1 grid layer for overlap */
+ if (conv != 1.0)
+ G_message(_("Using metric conversion factor %f, step=%f"), conv, step);
- bxmi = bxmi - mixx + stepx;
- bymi = bymi - miyy + stepy;
- bxma = bxma - mixx - stepx;
- byma = byma - miyy - stepy;
- mx2 = mx2o - 2*((int) (stepx / bresx));
- my2 = my2o - 2*((int) (stepy / bresy));
+ /*
+ * G_set_embedded_null_value_mode(1);
+ */
- parm.elevin = G_define_option();
- parm.elevin->key = "elevin";
- parm.elevin->type = TYPE_STRING;
- parm.elevin->required = YES;
- parm.elevin->gisprompt = "old,cell,raster";
- parm.elevin->description = _("Name of the elevation raster map, same units as N,E: [m] or [ft] ");
- parm.elevin->guisection = _("Input_options");
+ if ((depth == NULL) && (disch == NULL) && (err == NULL) && (outwalk == NULL))
+ G_warning(_("You are not outputting any raster or vector points maps"));
+ ret_val = input_data();
+ if (ret_val != 1)
+ G_fatal_error(_("Input failed"));
- parm.dxin = G_define_option();
- parm.dxin->key = "dxin";
- parm.dxin->type = TYPE_STRING;
- parm.dxin->required = YES;
- parm.dxin->gisprompt = "old,cell,raster";
- parm.dxin->description = _("Name of the x-derivatives raster map [m/m] or [ft/ft]");
- parm.dxin->guisection = _("Input_options");
+ /* memory allocation for output grids */
- parm.dyin = G_define_option();
- parm.dyin->key = "dyin";
- parm.dyin->type = TYPE_STRING;
- parm.dyin->required = YES;
- parm.dyin->gisprompt = "old,cell,raster";
- parm.dyin->description = _("Name of the y-derivatives raster map [m/m] or [ft/ft]");
- parm.dyin->guisection = _("Input_options");
+ gama = (double **)G_malloc (sizeof(double *)*(my));
+ for(l=0;l<my;l++)
+ {
+ gama[l] = (double*)G_malloc (sizeof(double)*(mx));
+ }
+ for (j = 0; j < my; j++)
+ {
+ for (i = 0; i < mx; i++)
+ gama[j][i] = 0.;
+ }
- parm.rain = G_define_option();
- parm.rain->key = "rain";
- parm.rain->type = TYPE_STRING;
- parm.rain->required = YES;
- parm.rain->gisprompt = "old,cell,raster";
- parm.rain->description = _("Name of the rainfall excess rate (rain-infilt) raster map [mm/hr]");
- parm.rain->guisection = _("Input_options");
+ if (err != NULL)
+ {
+ gammas = (double **)G_malloc (sizeof(double *)*(my));
+ for(l=0;l<my;l++)
+ {
+ gammas[l] = (double*)G_malloc (sizeof(double)*(mx));
+ }
+ for (j = 0; j < my; j++)
+ {
+ for (i = 0; i < mx; i++)
+ gammas[j][i] = 0.;
+ }
+ }
- parm.infil = G_define_option();
- parm.infil->key = "infil";
- parm.infil->type = TYPE_STRING;
- parm.infil->required = YES;
- parm.infil->gisprompt = "old,cell,raster";
- parm.infil->description = _("Name of the runoff infiltration rate raster map [mm/hr]");
- parm.infil->guisection = _("Input_options");
+ dif = (float **)G_malloc (sizeof(float *)*(my));
+ for(l=0;l<my;l++)
+ {
+ dif[l] = (float*)G_malloc (sizeof(float)*(mx));
+ }
+ for (j = 0; j < my; j++)
+ {
+ for (i = 0; i < mx; i++)
+ dif[j][i] = 0.;
+ }
+
+ /* if (maskmap != NULL)
+ bitmask = BM_create (cols, rows);
+ IL_create_bitmask (¶ms, bitmask);
+ */
+ seeds(rand1,rand2);
+ grad_check();
+ main_loop();
- parm.traps = G_define_option();
- parm.traps->key = "traps";
- parm.traps->type = TYPE_STRING;
- parm.traps->required = NO;
- parm.traps->gisprompt = "old,cell,raster";
- parm.traps->description = _("Name of the flow controls raster map (permeability ratio 0-1)");
- parm.traps->guisection = _("Input_options");
+ if (ts == 0) {
+ ii=output_data (0,1.);
+ if (ii != 1)
+ G_fatal_error(_("Cannot write raster maps"));
+ }
- parm.manin = G_define_option();
- parm.manin->key = "manin";
- parm.manin->type = TYPE_STRING;
- parm.manin->required = YES;
- parm.manin->gisprompt = "old,cell,raster";
- parm.manin->description = _("Name of the Mannings n raster map");
- parm.manin->guisection = _("Input_options");
+ if (fdwalkers != NULL)
+ fclose (fdwalkers);
-/* needs to be updated to GRASS 6 vector format !! */
- parm.sfile = G_define_option ();
- parm.sfile->key = "vector";
- parm.sfile->type = TYPE_STRING;
- parm.sfile->required = NO;
- parm.sfile->gisprompt = "old,vector,vector";
- parm.sfile->description = _("Name of the vector points map with x,y locations");
- parm.sfile->guisection = _("Input_options");
-
- parm.depth = G_define_option();
- parm.depth->key = "depth";
- parm.depth->type = TYPE_STRING;
- parm.depth->required = NO;
- parm.depth->gisprompt = "new,cell,raster";
- parm.depth->description = _("Output water depth raster map");
- parm.depth->guisection = _("Output_options");
-
- parm.disch = G_define_option();
- parm.disch->key = "disch";
- parm.disch->type = TYPE_STRING;
- parm.disch->required = NO;
- parm.disch->gisprompt = "new,cell,raster";
- parm.disch->description = _("Output water discharge raster map");
- parm.disch->guisection = _("Output_options");
-
- parm.err = G_define_option();
- parm.err->key = "err";
- parm.err->type = TYPE_STRING;
- parm.err->required = NO;
- parm.err->gisprompt = "new,cell,raster";
- parm.err->description = _("Output simulation error raster map");
- parm.err->guisection = _("Output_options");
-
- parm.outwalk = G_define_option ();
- parm.outwalk->key = "outwalk";
- parm.outwalk->type = TYPE_STRING;
- parm.outwalk->required = NO;
- parm.outwalk->gisprompt = "new,vector,vector";
- parm.outwalk->description = _("Name of the output walkers vector points map");
- parm.outwalk->guisection = _("Output_options");
-
- parm.nwalk = G_define_option();
- parm.nwalk->key = "nwalk";
- parm.nwalk->type = TYPE_INTEGER;
- parm.nwalk->answer = NWALK;
- parm.nwalk->required = NO;
- parm.nwalk->description = _("Number of walkers");
- parm.nwalk->guisection = _("Parameters");
-
- parm.niter = G_define_option();
- parm.niter->key = "niter";
- parm.niter->type = TYPE_INTEGER;
- parm.niter->answer = NITER;
- parm.niter->required = NO;
- parm.niter->description = _("Number of time iterations (sec.)");
- parm.niter->guisection = _("Parameters");
-
- parm.outiter = G_define_option();
- parm.outiter->key = "outiter";
- parm.outiter->type = TYPE_INTEGER;
- parm.outiter->answer = ITEROUT;
- parm.outiter->required = NO;
- parm.outiter->description = _("Time step for saving output maps (sec.)");
- parm.outiter->guisection = _("Parameters");
-
- parm.density = G_define_option();
- parm.density->key = "density";
- parm.density->type = TYPE_INTEGER;
- parm.density->answer = DENSITY;
- parm.density->required = NO;
- parm.density->description = _("Density of output walkers");
- parm.density->guisection = _("Parameters");
-
- parm.diffc = G_define_option();
- parm.diffc->key = "diffc";
- parm.diffc->type = TYPE_DOUBLE;
- parm.diffc->answer = DIFFC;
- parm.diffc->required = NO;
- parm.diffc->description = _("Water diffusion constant");
- parm.diffc->guisection = _("Parameters");
-
- parm.hmax = G_define_option();
- parm.hmax->key = "hmax";
- parm.hmax->type = TYPE_DOUBLE;
- parm.hmax->answer = HMAX;
- parm.hmax->required = NO;
- parm.hmax->description = _("Threshold water depth (diffusion increases after this water depth is reached)");
- parm.hmax->guisection = _("Parameters");
-
- parm.halpha = G_define_option();
- parm.halpha->key = "halpha";
- parm.halpha->type = TYPE_DOUBLE;
- parm.halpha->answer = HALPHA;
- parm.halpha->required = NO;
- parm.halpha->description = _("Diffusion increase constant");
- parm.halpha->guisection = _("Parameters");
-
- parm.hbeta = G_define_option();
- parm.hbeta->key = "hbeta";
- parm.hbeta->type = TYPE_DOUBLE;
- parm.hbeta->answer = HBETA;
- parm.hbeta->required = NO;
- parm.hbeta->description = _("Weighting factor for water flow velocity vector");
- parm.hbeta->guisection = _("Parameters");
-
- flag.mscale = G_define_flag ();
- flag.mscale->key = 'm';
- flag.mscale->description = _("Multiscale simulation");
-
- flag.tserie = G_define_flag ();
- flag.tserie->key = 't';
- flag.tserie->description = _("Time-series (dynamic) output");
-
-
- if (G_parser (argc, argv))
- exit (EXIT_FAILURE);
-
- mscale=flag.mscale->answer;
- ts=flag.tserie->answer;
-
- elevin = parm.elevin->answer;
- dxin = parm.dxin->answer;
- dyin = parm.dyin->answer;
-/* maskmap = parm.maskmap->answer;*/
- rain = parm.rain->answer;
- infil = parm.infil->answer;
- traps = parm.traps->answer;
- manin = parm.manin->answer;
- depth = parm.depth->answer;
- disch = parm.disch->answer;
- err = parm.err->answer;
- outwalk = parm.outwalk->answer;
- sfile = parm.sfile->answer;
-
- sscanf(parm.nwalk->answer, "%d", &maxwa);
- sscanf(parm.niter->answer, "%d", ×ec);
- sscanf(parm.outiter->answer, "%d", &iterout);
- sscanf(parm.density->answer, "%d", &ldemo);
- sscanf(parm.diffc->answer, "%lf", &frac);
- sscanf(parm.hmax->answer, "%lf", &hhmax);
- sscanf(parm.halpha->answer, "%lf", &halpha);
- sscanf(parm.hbeta->answer, "%lf", &hbeta);
-
- /* compute how big the raster is and set this to appr 2 walkers per cell */
- rwalk = (double) maxwa;
-
- if (conv != 1.0)
- G_message(_("Using metric conversion factor %f, step=%f"), conv, step);
-
-
- /*
- * G_set_embedded_null_value_mode(1);
- */
-
- if ((depth == NULL) && (disch == NULL) && (err == NULL) && (outwalk == NULL))
- G_warning(_("You are not outputting any raster or vector points maps"));
-
- ret_val = input_data();
- if (ret_val != 1)
- G_fatal_error(_("Input failed"));
-
-/* memory allocation for output grids */
-
- gama = (double **)G_malloc (sizeof(double *)*(my));
- for(l=0;l<my;l++)
- {
- gama[l] = (double*)G_malloc (sizeof(double)*(mx));
- }
- for (j = 0; j < my; j++)
- {
- for (i = 0; i < mx; i++)
- gama[j][i] = 0.;
- }
-
- if (err != NULL)
- {
- gammas = (double **)G_malloc (sizeof(double *)*(my));
- for(l=0;l<my;l++)
- {
- gammas[l] = (double*)G_malloc (sizeof(double)*(mx));
- }
- for (j = 0; j < my; j++)
- {
- for (i = 0; i < mx; i++)
- gammas[j][i] = 0.;
- }
- }
-
- dif = (float **)G_malloc (sizeof(float *)*(my));
- for(l=0;l<my;l++)
- {
- dif[l] = (float*)G_malloc (sizeof(float)*(mx));
- }
- for (j = 0; j < my; j++)
- {
- for (i = 0; i < mx; i++)
- dif[j][i] = 0.;
- }
-
-
-
-
-/* if (maskmap != NULL)
- bitmask = BM_create (cols, rows);
-
- IL_create_bitmask (¶ms, bitmask);
-*/
-
-
- seeds(rand1,rand2);
-
- grad_check();
-
-/* fprintf (stdout, "\n Percent complete: ");*/
-
- main_loop();
-
- if (ts == 0) {
- ii=output_data (0,1.);
- if (ii != 1)
- G_fatal_error(_("Cannot write raster maps"));
- }
-
- if (fdwalkers != NULL)
- fclose (fdwalkers);
-
- if (sfile != NULL)
- fclose(fw);
-
- exit (EXIT_SUCCESS);
+ if (sfile != NULL)
+ fclose(fw);
+ /* Exit with Success */
+ exit (EXIT_SUCCESS);
}
More information about the grass-commit
mailing list