r.horizon + r.sun: parameter changes (backported trunk https://trac.osgeo.org/grass/changeset/61096, https://trac.osgeo.org/grass/changeset/61098, https://trac.osgeo.org/grass/changeset/61099)

git-svn-id: https://svn.osgeo.org/grass/grass/branches/releasebranch_7_0@62430 15284696-431f-4ddb-bdfa-cd5b030d7da7

raster/r.horizon/r.horizon.html

@@ -35,24 +35,33 @@ not applied.
 
 
 <h3>Input parameters:</h3>
-<p>The <i>elev_in</i> parameter is an input elevation raster map. If
+<p>The <i>elevation</i> parameter is an input elevation raster map. If
 the buffer options are used (see below), this raster should extend
 over the area that accommodate the presently defined region plus
 defined buffer zones. 
 
-<p>The <i>horizon_step</i> parameter gives the angle step (in degrees)
+<p>The <i>step</i> parameter gives the angle step (in degrees)
 between successive azimuthal directions for the calculation of the
-horizon. Thus, a value of 5 for the <i>horizon_step</i> will give a total of
+horizon. Thus, a value of 5 for the <i>step</i> will give a total of
 360/5=72 directions (72 raster maps if used in the raster map mode). 
 
+<p>The <i>start</i> parameter gives the angle start (in degrees)
+for the calculation of the horizon. The default value is 0.
+
+<p>The <i>end</i> parameter gives the angle end (in degrees)
+for the calculation of the horizon. The end point is omitted!
+So for example if we run r.horizon with step=10, start=30 and end=70
+the raster maps generated by r.horizon will be only for angles:
+30, 40, 50, 60. The default value is 360.
+
 <p>The <i>direction</i> parameter gives the initial direction of the
 first output. This parameter acts as an direction angle offset. For
 example, if you want to get horizon angles for directions 45 and 225
-degrees, the <i>direction</i> should be set to 45 and <i>horizon_step</i> to
+degrees, the <i>direction</i> should be set to 45 and <i>step</i> to
 180. If you only want one single direction, use this parameter to
-specify desired direction of horizon angle, and set the <i>horizon_step</i>
+specify desired direction of horizon angle, and set the <i>step</i>
 size to 0 degrees. Otherwise all angles for a given starting <i>direction</i>
-with step of <i>horizon_step</i> are calculated.
+with step of <i>step</i> are calculated.
 
 <p>The <i>distance</i> controls the sampling distance step size for the
 search for horizon along the line of sight. The default value is 1.0
@@ -87,10 +96,11 @@ a buffer size individually for each of the four directions. The
 buffer parameters influence only size of the read elevation map,
 while the analysis in the raster mode will be done only for the
 area specified by the current region definition.
-<p>The <i>horizon </i>parameter gives the prefix of the output
+
+<p>The <i>basename</i> parameter gives the basename of the output
 horizon raster maps. The raster name of each horizon direction
-raster will be constructed as <i>horizon_</i>NNN , where NNN counts
-upwards from 0 to total number of directions. If you use <b>r.horizon</b>
+raster will be constructed as <i>basename_</i>ANGLE, where ANGLE
+is the angle in degrees with the direction. If you use <b>r.horizon</b>
 in the single point mode this option will be ignored. 
 
 <p>The <i>output</i> parameter allows to save the resulting horizon
@@ -116,7 +126,7 @@ higher than any point in the region or until it reaches the border of
 the region (see also the <i>bufferzone,e_buff</i>, <i>n_buff</i>,
 <i>s_buff</i>, and <i>w_buff</i>). The the number of lines of sight (azimuth 
 directions) is determined from the <i>direction</i> and
-<i>horizon_step</i> parameters. The method takes into account the curvature
+<i>step</i> parameters. The method takes into account the curvature
 of the Earth whereby remote features will seem to be lower than they
 actually are. It also accounts for the changes of angles towards
 cardinal directions caused by the projection (see above). 
@@ -128,48 +138,52 @@ The examples are intended for the North Carolina sample dataset.
 
 <h3>Single point mode</h3>
 
-Single point mode (output of horizon angles CCW from East):
+<b>Example 1</b>: determine horizon angle in 225 degree direction (output
+of horizon angles CCW from East):
+
 <div class="code"><pre>
-# determine horizon angle in 225 degree direction:
-r.horizon elev_in=elevation direction=215 horizon_step=0 bufferzone=200 \
+g.region rast=elevation -p
+r.horizon elevation=elevation direction=215 step=0 bufferzone=200 \
     coordinate=638871.6,223384.4 maxdistance=5000
+</pre></div>
+<p>
+<b>Example 2</b>: determine horizon values starting at 90 deg (North), 
+step size of 5 deg, saving result as CSV file:
 
-# determine horizon values starting at 90 deg (North), step size of 5 deg:
-# save result in CSV file
-r.horizon elev_in=elevation direction=90 horizon_step=5 bufferzone=200 \
+<div class="code"><pre>
+r.horizon elevation=elevation direction=90 step=5 bufferzone=200 \
     coordinate=638871.6,223384.4 maxdistance=5000 output=horizon.csv
+</pre></div>
+
+<p>
+<b>Example 3</b>: test point near highway intersection, saving result
+as CSV file for plotting the horizon around the highway intersection:
 
-# test point near high way intersection
+<div class="code"><pre>
 g.region n=223540 s=220820 w=634650 e=638780 res=10 -p
-r.horizon elev_in=elevation direction=0 horizon_step=5 bufferzone=200 \
+r.horizon elevation=elevation direction=0 step=5 bufferzone=200 \
     coordinate=636483.54,222176.25 maxdistance=5000 -d output=horizon.csv
 </pre></div>
 
+<center>
 <img src="rhorizon_shaded_dem_point.png"><br>
 Test point near high way intersection (North Carolina sample dataset)
 <p>
 <img src="rhorizon_singlepoint_plot.png"><br>
 Horizon angles for test point (CCW from East)
-<p>
+</center>
 
 <h3>Raster map mode</h3>
 
 Raster map mode (output maps "horangle*" become input for <em>r.sun</em>):
 <div class="code"><pre>
 # we put a bufferzone of 10% of maxdistance around the study area
-r.horizon elev_in=elevation horizon_step=30 bufferzone=200 horizon=horangle \
-    maxdistance=5000
+# compute only direction between 90 and 270 degrees
+r.horizon elevation=elevation step=30 start=90 end=300 \
+    bufferzone=200 basename=horangle maxdistance=5000
 </pre></div>
 
 
-<h2>SEE ALSO</h2>
-
-<em>
-<a href="r.sun.html">r.sun</a>,
-<a href="r.sunmask.html">r.sunmask</a>,
-<a href="r.viewshed.html">r.viewshed</a></em>
-
-
 <h2>REFERENCES</h2>
 <p>Hofierka J., 1997. Direct solar radiation modelling within an
 open GIS environment. Proceedings of JEC-GI'97 conference in Vienna,
@@ -190,6 +204,13 @@ A New GIS-based Solar Radiation Model and Its Application for
 Photovoltaic Assessments. <a href="http://www.blackwellpublishing.com/toc.asp?ref=1361-1682">Transactions
 in GIS</a>, 8(2), 175-190
 
+<h2>SEE ALSO</h2>
+
+<em>
+<a href="r.sun.html">r.sun</a>,
+<a href="r.sunmask.html">r.sunmask</a>,
+<a href="r.viewshed.html">r.viewshed</a>
+</em>
 
 <h2>AUTHORS</h2>
 <p>Thomas Huld, Joint Research Centre of

raster/r.sun/main.c

@@ -98,6 +98,10 @@ const char *glob_rad = NULL;
 const char *mapset = NULL;
 const char *per;
 
+size_t decimals;
+char *str_step;
+
+
 struct Cell_head cellhd;
 struct pj_info iproj;
 struct pj_info oproj;
@@ -202,8 +206,6 @@ double distance(double x1, double x2, double y1, double y2)
 }
 
 
-
-
 int main(int argc, char *argv[])
 {
     double singleSlope;
@@ -249,7 +251,7 @@ int main(int argc, char *argv[])
 	 "the topography is optionally incorporated.");
 
     parm.elevin = G_define_option();
-    parm.elevin->key = "elev_in";
+    parm.elevin->key = "elevation";
     parm.elevin->type = TYPE_STRING;
     parm.elevin->required = YES;
     parm.elevin->gisprompt = "old,cell,raster";
@@ -258,7 +260,7 @@ int main(int argc, char *argv[])
     parm.elevin->guisection = _("Input options");
 
     parm.aspin = G_define_option();
-    parm.aspin->key = "asp_in";
+    parm.aspin->key = "aspect";
     parm.aspin->type = TYPE_STRING;
     parm.aspin->required = NO;
     parm.aspin->gisprompt = "old,cell,raster";
@@ -267,7 +269,7 @@ int main(int argc, char *argv[])
     parm.aspin->guisection = _("Input options");
 
     parm.aspect = G_define_option();
-    parm.aspect->key = "aspect";
+    parm.aspect->key = "aspect_value";
     parm.aspect->type = TYPE_DOUBLE;
     parm.aspect->answer = ASPECT;
     parm.aspect->required = NO;
@@ -276,7 +278,7 @@ int main(int argc, char *argv[])
     parm.aspect->guisection = _("Input options");
 
     parm.slopein = G_define_option();
-    parm.slopein->key = "slope_in";
+    parm.slopein->key = "slope";
     parm.slopein->type = TYPE_STRING;
     parm.slopein->required = NO;
     parm.slopein->gisprompt = "old,cell,raster";
@@ -285,7 +287,7 @@ int main(int argc, char *argv[])
     parm.slopein->guisection = _("Input options");
 
     parm.slope = G_define_option();
-    parm.slope->key = "slope";
+    parm.slope->key = "slope_value";
     parm.slope->type = TYPE_DOUBLE;
     parm.slope->answer = SLOPE;
     parm.slope->required = NO;
@@ -293,7 +295,7 @@ int main(int argc, char *argv[])
     parm.slope->guisection = _("Input options");
 
     parm.linkein = G_define_option();
-    parm.linkein->key = "linke_in";
+    parm.linkein->key = "linke";
     parm.linkein->type = TYPE_STRING;
     parm.linkein->required = NO;
     parm.linkein->gisprompt = "old,cell,raster";
@@ -333,7 +335,7 @@ int main(int argc, char *argv[])
     }
 
     parm.latin = G_define_option();
-    parm.latin->key = "lat_in";
+    parm.latin->key = "lat";
     parm.latin->type = TYPE_STRING;
     parm.latin->required = NO;
     parm.latin->gisprompt = "old,cell,raster";
@@ -342,7 +344,7 @@ int main(int argc, char *argv[])
     parm.latin->guisection = _("Input options");
 
     parm.longin = G_define_option();
-    parm.longin->key = "long_in";
+    parm.longin->key = "long";
     parm.longin->type = TYPE_STRING;
     parm.longin->required = NO;
     parm.longin->gisprompt = "old,cell,raster";
@@ -368,12 +370,11 @@ int main(int argc, char *argv[])
 	_("Name of real-sky diffuse radiation coefficient (haze) input raster map [0-1]");
     parm.coefdh->guisection = _("Input options");
 
-    parm.horizon = G_define_option();
-    parm.horizon->key = "horizon";
-    parm.horizon->type = TYPE_STRING;
+    parm.horizon = G_define_standard_option(G_OPT_R_BASENAME_INPUT);
+    parm.horizon->key = "horizon_basename";
     parm.horizon->required = NO;
     parm.horizon->gisprompt = "old,cell,raster";
-    parm.horizon->description = _("The horizon information input map prefix");
+    parm.horizon->description = _("The horizon information input map basename");
     parm.horizon->guisection = _("Input options");
 
     parm.horizonstep = G_define_option();
@@ -599,6 +600,7 @@ int main(int argc, char *argv[])
     if (parm.horizonstep->answer != NULL) {
 	if (sscanf(parm.horizonstep->answer, "%lf", &horizonStep) != 1)
 	    G_fatal_error(_("Error reading horizon step size"));
+        str_step = parm.horizonstep->answer;
 	if (horizonStep > 0.)
 	    setHorizonInterval(deg2rad * horizonStep);
 	else
@@ -784,8 +786,9 @@ int INPUT_part(int offset, double *zmax)
     static int *fd_shad;
     int fr1 = -1, fr2 = -1;
     int l, i, j;
-    char shad_filename[256];
+    char *shad_filename;
     char formatString[10];
+    double angle_deg = 0.;
 
     finalRow = m - offset - m / numPartitions;
     if (finalRow < 0) {
@@ -917,13 +920,24 @@ int INPUT_part(int offset, double *zmax)
 	 * else
 	 * {
 	 */
-	numDigits = (int)(log10(1. * arrayNumInt)) + 1;
+        decimals = G_get_num_decimals(str_step);
+        angle_deg = 0;
+        for (i = 0; i < arrayNumInt; i++) {
+            horizonbuf[i] = Rast_allocate_f_buf();
+            shad_filename = G_generate_basename(horizon, angle_deg, 
+                                                3, decimals);
+            fd_shad[i] = Rast_open_old(shad_filename, "");
+            angle_deg += horizonStep;
+            G_free(shad_filename);
+        }
+        /*
+        numDigits = (int)(log10(1. * arrayNumInt)) + 1;
 	sprintf(formatString, "%%s_%%0%dd", numDigits);
 	for (i = 0; i < arrayNumInt; i++) {
 	    horizonbuf[i] = Rast_allocate_f_buf();
 	    sprintf(shad_filename, formatString, horizon, i);
 	    fd_shad[i] = Rast_open_old(shad_filename, "");
-	}
+	} */
     }
     /*
      * }

raster/r.sun/r.sun.html

@@ -45,8 +45,8 @@ with the <a href="r.info.html">r.info</a> command.
 
 <p>
 The solar incidence angle raster map <i>incidout</i> is computed specifying 
-elevation raster map <i>elev_in</i>, aspect raster map <i>asp_in</i>, slope 
-steepness raster map <i>slope_in,</i> given the day <i>day</i> and local time
+elevation raster map <i>elevation</i>, aspect raster map <i>aspect</i>, slope 
+steepness raster map <i>slope,</i> given the day <i>day</i> and local time
 <i>time</i>. There is no need to define latitude for locations with known
 and defined projection/coordinate system (check it with the <a href="g.proj.html">
 g.proj</a>
@@ -70,7 +70,7 @@ of geographical latitude can be set as a constant for the whole computed
 region or, as an option, a grid representing spatially distributed values
 over a large region. The geographical latitude must be also in decimal system
 with positive values for northern hemisphere and negative for southern one.
-In similar principle the Linke turbidity factor (<i>linke_in</i>, <i>lin</i>
+In similar principle the Linke turbidity factor (<i>linke</i>, <i>lin</i>
 ) and ground albedo (<i>albedo</i>, <i>alb</i>) can be set. 
 <p>Besides clear-sky radiations, the user can compute a real-sky radiation (beam,
 diffuse) using <i>coef_bh</i> and <i>coef_dh</i> input raster maps defining
@@ -83,7 +83,7 @@ for beam and diffuse radiation (see Suri and Hofierka, 2004, section 3.2).
 <p>
 The solar irradiation or irradiance raster maps <i>beam_rad</i>, <i>diff_rad</i>,
 <i>refl_rad</i> are computed for a given day <i>day,</i> latitude <i>lat_in</i>,
-elevation <i>elev_in</i>, slope <i>slope_in</i> and aspect <i>asp_in</i> raster maps.
+elevation <i>elevation</i>, slope <i>slope</i> and aspect <i>aspect</i> raster maps.
 For convenience, the output raster given as <i>glob_rad</i>
 will output the sum of the three radiation components. The program uses 
 the Linke atmosphere turbidity factor and ground albedo coefficient. 
@@ -258,19 +258,20 @@ North Carolina example (considering also cast shadows):
 g.region rast=elevation -p
 
 # calculate horizon angles (to speed up the subsequent r.sun calculation)
-r.horizon elev_in=elevation horizon_step=30 bufferzone=200 horizon=horangle \
+r.horizon elevation=elevation step=30 bufferzone=200 basename=horangle \
     maxdistance=5000
 
 # slope + aspect
 r.slope.aspect elevation=elevation aspect=aspect.dem slope=slope.dem
 
-# calculate global radiation for day 180 at 14:00hrs, using r.horizon output
-r.sun elev_in=elevation horizon=horangle horizon_step=30 asp_in=aspect.dem \
-      slope_in=slope.dem glob_rad=global_rad day=180 time=14
+# calculate global radiation for day 180 at 2p.m., using r.horizon output
+r.sun elevation=elevation horizon_basename=horangle horizon_step=30 \
+      aspect=aspect.dem slope=slope.dem glob_rad=global_rad day=180 time=14
+# result: output global (total) irradiance/irradiation [W.m-2] for given day/time
+r.univar global_rad
 </pre></div>
 
 <p>
-
 Calculation of the integrated daily irradiation for a region in North-Carolina
 for a given day of the year at 30m resolution. Here day 172 (i.e., 21 June
 in non-leap years):
@@ -279,7 +280,7 @@ in non-leap years):
 g.region rast=elev_ned_30m -p
 
 # considering cast shadows
-r.sun elev_in=elev_ned_30m lin=2.5 alb=0.2 day=172 \
+r.sun elevation=elev_ned_30m lin=2.5 alb=0.2 day=172 \
       beam_rad=b172 diff_rad=d172 \
       refl_rad=r172 insol_time=it172