[GRASS-SVN] r60180 - grass-addons/grass7/vector/v.habitat.dem

svn_grass at osgeo.org svn_grass at osgeo.org
Sun May 11 09:01:50 PDT 2014


Author: hellik
Date: 2014-05-11 09:01:50 -0700 (Sun, 11 May 2014)
New Revision: 60180

Modified:
   grass-addons/grass7/vector/v.habitat.dem/Makefile
   grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.html
   grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.py
Log:
svn propset

Modified: grass-addons/grass7/vector/v.habitat.dem/Makefile
===================================================================
--- grass-addons/grass7/vector/v.habitat.dem/Makefile	2014-05-11 15:52:33 UTC (rev 60179)
+++ grass-addons/grass7/vector/v.habitat.dem/Makefile	2014-05-11 16:01:50 UTC (rev 60180)
@@ -1,7 +1,7 @@
-MODULE_TOPDIR = ../..
-
-PGM = v.habitat.dem
-
-include $(MODULE_TOPDIR)/include/Make/Script.make
-
-default: script
+MODULE_TOPDIR = ../..
+
+PGM = v.habitat.dem
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script


Property changes on: grass-addons/grass7/vector/v.habitat.dem/Makefile
___________________________________________________________________
Added: svn:mime-type
   + text/x-makefile
Added: svn:eol-style
   + native

Modified: grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.html
===================================================================
--- grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.html	2014-05-11 15:52:33 UTC (rev 60179)
+++ grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.html	2014-05-11 16:01:50 UTC (rev 60180)
@@ -1,210 +1,210 @@
-<h2>DESCRIPTION</h2>
-
-<em>v.habitat.dem</em> calculates DEM and solar derived characteristics of habitat 
-vector polygons. The user must specify the input <b>elevation raster</b> map,
-a <b>habitat vector</b> map with a <b>numeric unique ID</b> column and a <b>prefix</b>
-used for all results.
-
-<p>
-A preliminary visual delineation of habitats based upon digital orthophotos is a 
-common task for an ecologist before fieldwork. Ecological site conditions of habitats 
-are often influenced amongst others by terrain forms, solar irradiance and irradiation.
-<em>v.habitat.dem</em> gives some DEM derived characteristics for a quick validation
-of the preliminary visual habitat delineation.
-</p>
-
-<h2>NOTES</h2>
-
-The location has to be in a projected coordination system. Before running <em>v.habitat.dem</em> 
-the region has to be aligned to the <b>elevation raster</b> map and the <b>habitat vector</b> map by 
-<em>g.region</em>. During calculations, especially for solar characteristics, the region will be 
-extended by a user input (default 5.000m). The results are as good as the DEM quality and resolution is.
-
-
-<h3>Terrain characteristics</h3>
-
-<b>Slope</b> and <b>aspect</b> are calculated by <em>r.slope.aspect</em>.
-
-<p>
-The <b>slope</b> output raster map contains slope values, stated in degrees of
-inclination from the horizontal.
-<br>
-The <b>aspect</b> output raster map indicates the direction that slopes are
-facing. The aspect categories represent the number degrees of east.
-</p>
-
-<b>Accumulation</b>, <b>drainage direction</b> and <b>topographic index</b> are 
-calculated by <em>r.watershed</em>. The flag <b>-a</b> (use positive flow accumulation 
-even for likely underestimates) is used as default.
-
-<p>
-The <b>accumulation</b> map contains the absolute value of each cell in this output 
-map and is the amount of overland flow that traverses the cell. This value will be 
-the number of upland cells plus one if no overland flow map is given.
-<br>
-The <b>drainage direction</b> map contains drainage direction. Provides the "aspect" 
-for each cell measured CCW from East.
-<br>
-The <b>topographic index</b> raster map contains topographic index TCI and is
-computed as <tt>ln(α / tan(β))</tt> where α a is the cumulative
-upslope area draining through a point per unit contour length and <tt>tan(β)</tt> 
-is the local slope angle. The TCI reflects the tendency of water to accumulate at 
-any point in the catchment and the tendency for gravitaional forces to move that 
-water downslope.  This value will be negative if <tt>α / tan(β) < 1</tt>.
-</p>
-
-Terrain forms are calculated by <em>r.geomorphon</em>.
-
-<p>
-Geomorphon is a new concept of presentation and analysis of terrain forms using machine vision 
-approach. This concept utilises 8-tuple pattern of the visibility neighbourhood and breaks 
-well known limitation of standard calculus approach where all terrain forms cannot 
-be detected in a single window size. The pattern arises from a comparison of a focus 
-pixel with its eight neighbours starting from the one located to the east and continuing
-counterclockwise producing a ternary operator. All options in the <em>r.geomorphon</em>-calculation 
-are set to default (<b>skip = 0</b>, <b>search = 3</b>, <b>flat = 1</b>, <b>dist = 0</b>) where 
-<b>search</b> determines the length on the geodesic distances in all eight directions where 
-line-of-sight is calculated, <b>skip</b> determines length on the geodesic distances at the 
-beginning of calculation all eight directions where line-of-sight is yet calculated, <b>flat</b> defines 
-the difference (in degrees) between zenith and nadir line-of-sight which indicate flat direction and <b>dist</b>
-determines > flat distance..
-</p>
-
- The most common terrain forms calculated by <em>r.geomorphon</em> are:
-
-<ul>
- <li>flat</li>
- <li>summit</li>
- <li>ridge</li>
- <li>shoulder</li>
- <li>spur</li>
- <li>slope</li>
- <li>hollow</li>
- <li>footslope</li>
- <li>valley</li>
- <li>depression</li>
-</ul>
-
-The LS factor
-
-<p>
-The LS is the slope length-gradient factor. The LS factor represents a ratio of soil loss under given conditions to 
-that at a site with the "standard" slope steepness of 9% and slope length of 22.13m. The steeper and longer the slope, 
-the higher the risk for erosion.
-</p>
-
-The LS factor is calculated accordingly Neteler & Mitasova 2008 in <em>r.mapcalc</em> with flow accumulation of <em>r.flow</em>
-and slope of <em>r.slope.aspect</em>
-
-<div class="code">
- <pre>
-  1.4 * exp(flow_accumulation / 22.1, 0.4) * exp(sin(slope) 0.09, 1.2)
- </pre>
-</div>
-
-The colors of the LS factor map are set to:
-
-<ul>
- <li>0 white</li>
- <li>3 yellow</li>
- <li>6 orange</li>
- <li>10 red</li>
- <li>50 magenta</li>
- <li>100 violet</li>
-</ul>
-
-<h3>Terrain characteristics uploaded to the habitat vector attribute table per polygon</h3>
-
-<ul>
- <li><b>DEM altitude</b>: minimum, maximum, range, average, median</li>
- <li><b>slope</b>: minimum, maximum, range, average, median</li>
- <li><b>aspect</b>: minimum, maximum, range, average, median</li>
- <li><b>geomorphons</b>: absolute area of flat, summit, ridge, shoulder, spur, slope, hollow, footslope, valley, depression</li> 
-</ul>
-
-Additionally the mutual occurrence by <em>r.coin</em> of unique habitat ID and geomorphons in percent of the row is printed to the output.
-
-<h3>Simple check of terrain characteristics</h3>
-
-Simple checks regarding aspect and slope per unique habitat ID are evaluated and marked in the attribute table as follow:
-
-<ul>
-<li><b>simple check regarding aspect range:</b></li>
-<li>aspect range 100-200 *</li>
-<li>aspect range 201-300 **</li>
-<li>aspect range ≥ 300 ***</li>
-<li><b>simple checks regarding aspect range and slope median:</b></li>
-<li>aspect range 100-200 and median slope < 5 *</li>
-<li>aspect range 201-300 and median slope < 5 **</li>
-<li>aspect range ≥ 300 and median slope < 5 ***</li>
-</ul>
-
-These simple checks may indicate reconsidering of some preliminary visual habitat delineations.
-
-<h3>Solar characteristics</h3>
-
-<p>
-The solar characterstics (direct sunlight / shadows caused by terrain for a certain day in the year) are 
-calculated by <em>r.sun.hourly</em> based upon <em>r.sun</em>. The <b>-b</b>-flag is used to create binary rasters 
-instead of irradiation rasters. The user can define start time of interval, end time of interval, time step for 
-running <em>r.sun</em>, number of day of the year and the year. As default is set summer solstice (21st June 2014, 8:00-18:00, 1 hour time step).
-</p>
-
-The results of the <em>r.sun.hourly</em>-analysis are automatically registered into a temporal database. The space time raster dataset 
-can be easily animated in the <em>g.gui.animation</em>-tool.
-
-<h2>EXAMPLE</h2>
-
-<div class="code">
- <pre>
-  # align region to DEM and habitat vector
-  g.region -a rast=DEM vect=myhabitats align=DEM
-
-  # run <em>v.habitat.dem</em>
-  v.habitat.dem elevation=DEM vector=myhabitats column=Id prefix=a dir=C:\wd
-
-  # do <em>r.null</em> to the <em>r.sun.hourly</em> output to get maps without direct beam
-  r.null map=a_beam_rad_08.00 setnull=1
-  [...]
-  r.null map=a_beam_rad_18.00 setnull=1
-
-  # animate the <em>r.sun.hourly</em> output by the <em>g.gui.animation</em>-tool
-  g.gui.animation strds=a_beam_rad
- </pre>
-</div>
-
-<h2>DEPENDENCIES</h2>
- 
-<ul>
-<li>r.geomorphon (addon)</li>
-<li>r.sun.hourly (addon)</li>
-</ul>
-
-<h2>SEE ALSO</h2>
-
-<em>
-<a href="g.gui.animation.html">g.gui.animation</a>
-<a href="g.region.html">g.region</a>,
-<a href="r.coin.html">r.coin</a>,
-<a href="r.geomorphon.html">r.geomorphon</a>,
-<a href="r.mapcalc.html">r.mapcalc</a>,
-<a href="r.slope.aspect.html">r.slope.aspect</a>,
-<a href="r.sun.html">r.sun</a>,
-<a href="r.sun.hourly.html">r.sun.hourly</a>,
-<a href="r.stats.html">r.stats</a>,
-<a href="r.watershed.html">r.watershed</a>,
-<a href="v.rast.stats.html">v.rast.stats</a>,
-<a href="v.to.rast.html">v.to.rast</a>
-</em>
-
-<h2>REFERENCES</h2>
-
-Neteler, M. & Mitasova, H. 2008. Open Source GIS. A GRASS GIS Aproach. Third Edition. Springer.
-
-<h2>AUTHOR</h2>
-
-Helmut Kudrnovsky
-
-<p>
-<i>Last changed: $Date: 2014-05-11 12:00:00 +0200 (So, 11 May 2014) $</i>
+<h2>DESCRIPTION</h2>
+
+<em>v.habitat.dem</em> calculates DEM and solar derived characteristics of habitat 
+vector polygons. The user must specify the input <b>elevation raster</b> map,
+a <b>habitat vector</b> map with a <b>numeric unique ID</b> column and a <b>prefix</b>
+used for all results.
+
+<p>
+A preliminary visual delineation of habitats based upon digital orthophotos is a 
+common task for an ecologist before fieldwork. Ecological site conditions of habitats 
+are often influenced amongst others by terrain forms, solar irradiance and irradiation.
+<em>v.habitat.dem</em> gives some DEM derived characteristics for a quick validation
+of the preliminary visual habitat delineation.
+</p>
+
+<h2>NOTES</h2>
+
+The location has to be in a projected coordination system. Before running <em>v.habitat.dem</em> 
+the region has to be aligned to the <b>elevation raster</b> map and the <b>habitat vector</b> map by 
+<em>g.region</em>. During calculations, especially for solar characteristics, the region will be 
+extended by a user input (default 5.000m). The results are as good as the DEM quality and resolution is.
+
+
+<h3>Terrain characteristics</h3>
+
+<b>Slope</b> and <b>aspect</b> are calculated by <em>r.slope.aspect</em>.
+
+<p>
+The <b>slope</b> output raster map contains slope values, stated in degrees of
+inclination from the horizontal.
+<br>
+The <b>aspect</b> output raster map indicates the direction that slopes are
+facing. The aspect categories represent the number degrees of east.
+</p>
+
+<b>Accumulation</b>, <b>drainage direction</b> and <b>topographic index</b> are 
+calculated by <em>r.watershed</em>. The flag <b>-a</b> (use positive flow accumulation 
+even for likely underestimates) is used as default.
+
+<p>
+The <b>accumulation</b> map contains the absolute value of each cell in this output 
+map and is the amount of overland flow that traverses the cell. This value will be 
+the number of upland cells plus one if no overland flow map is given.
+<br>
+The <b>drainage direction</b> map contains drainage direction. Provides the "aspect" 
+for each cell measured CCW from East.
+<br>
+The <b>topographic index</b> raster map contains topographic index TCI and is
+computed as <tt>ln(α / tan(β))</tt> where α a is the cumulative
+upslope area draining through a point per unit contour length and <tt>tan(β)</tt> 
+is the local slope angle. The TCI reflects the tendency of water to accumulate at 
+any point in the catchment and the tendency for gravitaional forces to move that 
+water downslope.  This value will be negative if <tt>α / tan(β) < 1</tt>.
+</p>
+
+Terrain forms are calculated by <em>r.geomorphon</em>.
+
+<p>
+Geomorphon is a new concept of presentation and analysis of terrain forms using machine vision 
+approach. This concept utilises 8-tuple pattern of the visibility neighbourhood and breaks 
+well known limitation of standard calculus approach where all terrain forms cannot 
+be detected in a single window size. The pattern arises from a comparison of a focus 
+pixel with its eight neighbours starting from the one located to the east and continuing
+counterclockwise producing a ternary operator. All options in the <em>r.geomorphon</em>-calculation 
+are set to default (<b>skip = 0</b>, <b>search = 3</b>, <b>flat = 1</b>, <b>dist = 0</b>) where 
+<b>search</b> determines the length on the geodesic distances in all eight directions where 
+line-of-sight is calculated, <b>skip</b> determines length on the geodesic distances at the 
+beginning of calculation all eight directions where line-of-sight is yet calculated, <b>flat</b> defines 
+the difference (in degrees) between zenith and nadir line-of-sight which indicate flat direction and <b>dist</b>
+determines > flat distance..
+</p>
+
+ The most common terrain forms calculated by <em>r.geomorphon</em> are:
+
+<ul>
+ <li>flat</li>
+ <li>summit</li>
+ <li>ridge</li>
+ <li>shoulder</li>
+ <li>spur</li>
+ <li>slope</li>
+ <li>hollow</li>
+ <li>footslope</li>
+ <li>valley</li>
+ <li>depression</li>
+</ul>
+
+The LS factor
+
+<p>
+The LS is the slope length-gradient factor. The LS factor represents a ratio of soil loss under given conditions to 
+that at a site with the "standard" slope steepness of 9% and slope length of 22.13m. The steeper and longer the slope, 
+the higher the risk for erosion.
+</p>
+
+The LS factor is calculated accordingly Neteler & Mitasova 2008 in <em>r.mapcalc</em> with flow accumulation of <em>r.flow</em>
+and slope of <em>r.slope.aspect</em>
+
+<div class="code">
+ <pre>
+  1.4 * exp(flow_accumulation / 22.1, 0.4) * exp(sin(slope) 0.09, 1.2)
+ </pre>
+</div>
+
+The colors of the LS factor map are set to:
+
+<ul>
+ <li>0 white</li>
+ <li>3 yellow</li>
+ <li>6 orange</li>
+ <li>10 red</li>
+ <li>50 magenta</li>
+ <li>100 violet</li>
+</ul>
+
+<h3>Terrain characteristics uploaded to the habitat vector attribute table per polygon</h3>
+
+<ul>
+ <li><b>DEM altitude</b>: minimum, maximum, range, average, median</li>
+ <li><b>slope</b>: minimum, maximum, range, average, median</li>
+ <li><b>aspect</b>: minimum, maximum, range, average, median</li>
+ <li><b>geomorphons</b>: absolute area of flat, summit, ridge, shoulder, spur, slope, hollow, footslope, valley, depression</li> 
+</ul>
+
+Additionally the mutual occurrence by <em>r.coin</em> of unique habitat ID and geomorphons in percent of the row is printed to the output.
+
+<h3>Simple check of terrain characteristics</h3>
+
+Simple checks regarding aspect and slope per unique habitat ID are evaluated and marked in the attribute table as follow:
+
+<ul>
+<li><b>simple check regarding aspect range:</b></li>
+<li>aspect range 100-200 *</li>
+<li>aspect range 201-300 **</li>
+<li>aspect range ≥ 300 ***</li>
+<li><b>simple checks regarding aspect range and slope median:</b></li>
+<li>aspect range 100-200 and median slope < 5 *</li>
+<li>aspect range 201-300 and median slope < 5 **</li>
+<li>aspect range ≥ 300 and median slope < 5 ***</li>
+</ul>
+
+These simple checks may indicate reconsidering of some preliminary visual habitat delineations.
+
+<h3>Solar characteristics</h3>
+
+<p>
+The solar characterstics (direct sunlight / shadows caused by terrain for a certain day in the year) are 
+calculated by <em>r.sun.hourly</em> based upon <em>r.sun</em>. The <b>-b</b>-flag is used to create binary rasters 
+instead of irradiation rasters. The user can define start time of interval, end time of interval, time step for 
+running <em>r.sun</em>, number of day of the year and the year. As default is set summer solstice (21st June 2014, 8:00-18:00, 1 hour time step).
+</p>
+
+The results of the <em>r.sun.hourly</em>-analysis are automatically registered into a temporal database. The space time raster dataset 
+can be easily animated in the <em>g.gui.animation</em>-tool.
+
+<h2>EXAMPLE</h2>
+
+<div class="code">
+ <pre>
+  # align region to DEM and habitat vector
+  g.region -a rast=DEM vect=myhabitats align=DEM
+
+  # run <em>v.habitat.dem</em>
+  v.habitat.dem elevation=DEM vector=myhabitats column=Id prefix=a dir=C:\wd
+
+  # do <em>r.null</em> to the <em>r.sun.hourly</em> output to get maps without direct beam
+  r.null map=a_beam_rad_08.00 setnull=1
+  [...]
+  r.null map=a_beam_rad_18.00 setnull=1
+
+  # animate the <em>r.sun.hourly</em> output by the <em>g.gui.animation</em>-tool
+  g.gui.animation strds=a_beam_rad
+ </pre>
+</div>
+
+<h2>DEPENDENCIES</h2>
+ 
+<ul>
+<li>r.geomorphon (addon)</li>
+<li>r.sun.hourly (addon)</li>
+</ul>
+
+<h2>SEE ALSO</h2>
+
+<em>
+<a href="g.gui.animation.html">g.gui.animation</a>
+<a href="g.region.html">g.region</a>,
+<a href="r.coin.html">r.coin</a>,
+<a href="r.geomorphon.html">r.geomorphon</a>,
+<a href="r.mapcalc.html">r.mapcalc</a>,
+<a href="r.slope.aspect.html">r.slope.aspect</a>,
+<a href="r.sun.html">r.sun</a>,
+<a href="r.sun.hourly.html">r.sun.hourly</a>,
+<a href="r.stats.html">r.stats</a>,
+<a href="r.watershed.html">r.watershed</a>,
+<a href="v.rast.stats.html">v.rast.stats</a>,
+<a href="v.to.rast.html">v.to.rast</a>
+</em>
+
+<h2>REFERENCES</h2>
+
+Neteler, M. & Mitasova, H. 2008. Open Source GIS. A GRASS GIS Aproach. Third Edition. Springer.
+
+<h2>AUTHOR</h2>
+
+Helmut Kudrnovsky
+
+<p>
+<i>Last changed: $Date$</i>
 </p>
\ No newline at end of file


Property changes on: grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.html
___________________________________________________________________
Added: svn:mime-type
   + text/html
Added: svn:keywords
   + Author Date Id
Added: svn:eol-style
   + native

Modified: grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.py
===================================================================
--- grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.py	2014-05-11 15:52:33 UTC (rev 60179)
+++ grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.py	2014-05-11 16:01:50 UTC (rev 60180)
@@ -1,519 +1,519 @@
-#!/usr/bin/env python
-
-"""
-MODULE:    v.habitat.dem
-
-AUTHOR(S): Helmut Kudrnovsky <alectoria AT gmx at>
-
-PURPOSE:   DEM and solar derived characteristics of habitats
-
-COPYRIGHT: (C) 2014 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.
-"""
-
-#%module
-#% description: DEM derived characteristics of habitats
-#% keywords: vector
-#% keywords: raster
-#% keywords: terrain
-#% keywords: statistics
-#% keywords: sun
-#% keywords: zonal statistics
-#%end
-
-#%option G_OPT_R_ELEV
-#% key: elevation
-#% description: Name of elevation raster map 
-#% required: yes
-#%end
-
-#%option G_OPT_V_INPUT
-#% key: vector
-#% description: Name of habitat vector map 
-#% required: yes
-#%end
-
-#%option G_OPT_DB_COLUMN
-#% description: Name of attribute column with a unique habitat ID (must be numeric)
-#% required: yes
-#%end
-
-#%option
-#% key: prefix
-#% type: string
-#% key_desc: prefix
-#% description: output prefix (must start with a letter)
-#% required: yes
-#%end
-
-#%option G_OPT_M_DIR
-#% key: dir
-#% description: Directory where the output will be found
-#% required : yes
-#%end
-
-#%option
-#% key: region_extension
-#% type: double
-#% key_desc: float
-#% description: region extension
-#% required : no
-#% answer: 5000
-#%end
-
-#%option
-#% key: start_time
-#% type: double
-#% label: Start time of interval
-#% description: Use up to 2 decimal places
-#% options: 0-24
-#% answer: 8
-#%end
-
-#%option
-#% key: end_time
-#% type: double
-#% label: End time of interval
-#% description: Use up to 2 decimal places
-#% options: 0-24
-#% answer: 18
-#%end
-
-#%option
-#% key: time_step
-#% type: double
-#% label: Time step for running r.sun [decimal hours]
-#% description: Use up to 2 decimal places
-#% options: 0-24
-#% answer: 1
-#%end
-
-#%option
-#% key: day
-#% type: integer
-#% description: No. of day of the year
-#% options: 1-365
-#% answer: 172
-#%end
-
-#%option
-#% key: year
-#% type: integer
-#% label: Year used for map registration into temporal dataset or r.timestamp
-#% description: This value is not used in r.sun calculations
-#% options: 1900-9999
-#% answer: 2014
-#%end
-
-import sys
-import os
-import grass.script as grass
-import math
-from numpy import array
-from numpy import zeros
-import csv
-
-if not os.environ.has_key("GISBASE"):
-    grass.message( "You must be in GRASS GIS to run this program." )
-    sys.exit(1)
-
-def main():
-    r_elevation = options['elevation'].split('@')[0] 
-    v_habitat = options['vector'].split('@')[0] 
-    v_column = options['column']
-    regext = options['region_extension']
-    directory = options['dir']
-    prefix = options['prefix']
-    r_accumulation = prefix+'_accumulation'
-    r_drainage = prefix+'_drainage'
-    r_tci = prefix+'_topographicindex'
-    r_slope = prefix+'_slope'
-    r_aspect = prefix+'_aspect'
-    r_flow_accum = prefix+'_flow_accumulation'
-    r_LS = prefix+'_LS_factor'
-    r_habitat = prefix+'_'+v_habitat+'_rast'
-    r_geomorphon = prefix+'_'+r_elevation+'_geomorphon'
-    r_beam_rad = prefix+'_beam_rad'
-    f_habitat_small_areas = prefix+'_small_areas.csv'
-    f_habitat_geomorphons = prefix+'_habitat_geomorphons.csv'
-    f_habitat_geomorphons_pivot = prefix+'_habitat_geomorphons_pivot.csv'	
-    f_LS_colorrules = prefix+'_LS_color.txt'
-    t_habitat_geomorphons = prefix+'_habgeo_tab'
-    t_habitat_geomorphons_pivot = prefix+'_habgeo_tab_pivot'
-    d_start_time = options['start_time']
-    d_end_time = options['end_time']
-    d_time_step = options['time_step']
-    d_day = options['day']
-    d_year = options['year']
-    saved_region = prefix+'_region_ori'
-    current_region = grass.region()
-    X = float(regext)
-    N = current_region["n"]
-    S = current_region["s"]
-    E = current_region["e"]
-    W = current_region["w"]
-    Yres = current_region['nsres']
-    Xres = current_region['ewres']	
-    PixelArea = Xres * Yres	
-    global tmp
- 
-    ## check if r.geomorphon addon is installed
-    if not grass.find_program('r.geomorphon', '--help'):
-            grass.fatal(_("The 'r.geomorphon' addon was not found, install it first:") +
-                        "\n" +
-                        "g.extension r.geomorphon")	
-
-    ## check if r.geomorphon addon is installed
-    if not grass.find_program('r.sun.hourly', '--help'):
-            grass.fatal(_("The 'r.sun.hourly' addon was not found, install it first:") +
-                        "\n" +
-                        "g.extension r.sun.hourly")	
-						
-    # Region settings
-    grass.message( "Current region will be saved and extended." )	
-    grass.message( "----" )
-	
-    # Print and save current region
-    grass.run_command('g.region', flags = 'p', save = saved_region, overwrite = True)
-    grass.message( "Current region saved." )	
-    grass.message( "----" )
-	
-    # Align region to elevation raster and habitat vector
-    # grass.message( "Align region to elevation raster and habitat vector ..." )		
-    # grass.run_command('g.region', flags = 'ap',
-    #                                 rast = r_elevation,
-    #                                 vect = v_habitat,
-    #                                 align = r_elevation)
-    #grass.message( "Alignment done." )	
-    #grass.message( "----" )									 
-									 
-    # Extend region
-    grass.message( "Extend region ..." )
-    grass.message( "n, s, e, w" )
-    grass.message( [current_region[key] for key in "nsew"] )
-    grass.message( "by" )
-    grass.message( regext )
-    grass.message( "in all directions" )
-	
-    grass.run_command('g.region', n = N+X,
-                                     s = S-X,
-                                     e = E+X,
-                                     w = W-X)
-    grass.message( "Region extension done." )
-	
-    grass.message( "Extended region:" )
-    grass.run_command('g.region', flags = 'p')
-    grass.message( "----" )
-	
-    # Watershed calculation: accumulation, drainage direction, topographic index
-    grass.message( "Calculation of accumulation, drainage direction, topographic index by r.watershed ..." )
-    grass.run_command('r.watershed', elevation = r_elevation, 
-                                     accumulation = r_accumulation, 
-                                     drainage = r_drainage,
-                                     tci = r_tci, 									 
-                                     convergence = 5, 
-                                     flags = 'am',
-                                     overwrite = True)
-    grass.message( "Calculation of accumulation, drainage direction, topographic done." )
-    grass.message( "----" )
-									 
-    # Calculation of slope and aspect maps
-    grass.message( "Calculation of slope and aspect by r.slope.aspect ..." )
-    grass.run_command('r.slope.aspect', elevation = r_elevation,
-                                     slope = r_slope, 
-                                     aspect = r_aspect,
-                                     overwrite = True)
-    grass.message( "Calculation of slope and aspect done." )
-    grass.message( "----" )
- 
-    # Calculate pixel area by nsres x ewres
-    grass.message( "Pixel area:" )	
-    grass.message( PixelArea )	
-    grass.message( "----" )
-
-    # Calculate habitat area and populate it to the attribute table
-    grass.message( "Calculate habitat's areas and populate it to the attribute table ..." )	
-    grass.run_command("v.db.addcolumn", map = v_habitat,
-                                     layer = 1, 
-                                     columns = "habarea double")		
-	
-    grass.run_command("v.to.db", map = v_habitat,
-                                     option = 'area',
-                                     layer = 1, 
-                                     columns = 'habarea',
-                                     overwrite = True)
-
-    grass.message( "Calculate habitat's areas done." )	
-    grass.message( "----" )	 	
-
-    # Show habitat areas smaller than Pixel Area
-    grass.message( "Habitat areas smaller than pixel area." )	
-    grass.run_command("v.db.select", map = v_habitat,
-                                    flags = 'v',
-                                    layer = 1,
-                                    columns = v_column,
-                                    where = "habarea < %s" % (PixelArea))
-
-    smallareacsv = os.path.join( directory, f_habitat_small_areas )								
-									
-    grass.run_command("v.db.select", map = v_habitat,
-                                    flags = 'v',
-                                    layer = 1,
-                                    columns = v_column,
-                                    file = smallareacsv,
-                                    where = "habarea < %s" % (PixelArea))
-
-    grass.message( "A list of habitat areas smaller than pixel area can be found in: " )
-    grass.message( smallareacsv )	
-    grass.message( "----" )	    
-
-    # Mark habitats smaller than pixel area in attribute table
-    grass.message( "Mark habitats smaller than pixel area in attribute table ..." )	
-    grass.run_command("v.db.addcolumn", map = v_habitat,
-                                     layer = 1, 
-                                     columns = "%s_smallarea varchar(1)" % prefix )
-    grass.run_command("v.db.update", map = v_habitat, layer = 1, column = '%s_smallarea' % (prefix), value = '*', where = 'habarea < %s' % (PixelArea))
-    grass.message( "See column" )
-    grass.message( '%s_smallarea' % prefix )
-    grass.message( "marked by *." )
-    grass.message( "----" )	  
-
-    # Upload DEM zonal statistics to the attribute table
-    grass.message( "Upload DEM zonal statistics to the attribute table ..." )
-    grass.run_command("v.rast.stats", map = v_habitat,
-                                     flags = 'c',
-                                     layer = 1,
-                                     raster = r_elevation,
-                                     column_prefix = prefix+'_dem',
-                                     method = 'minimum,maximum,range,average,median')		
-    grass.message( "Upload DEM zonal statistics done." )		
-    grass.message( "----" )
-	
-    # Upload slope zonal statistics to the attribute table
-    grass.message( "Upload slope zonal statistics to the attribute table ..." )		
-    grass.run_command("v.rast.stats", map = v_habitat,
-                                     flags = 'c',
-                                     layer = 1,
-                                     raster = r_slope,
-                                     column_prefix = prefix+'_slope',
-                                     method = 'minimum,maximum,range,average,median')		
-    grass.message( "Upload slope zonal statistics done." )
-    grass.message( "----" )
-
-    # Upload slope zonal statistics to the attribute table
-    grass.message( "Upload aspect zonal statistics to the attribute table ..." )		
-    grass.run_command("v.rast.stats", map = v_habitat,
-                                     flags = 'c',
-                                     layer = 1,
-                                     raster = r_aspect,
-                                     column_prefix = prefix+'_aspect',
-                                     method = 'minimum,maximum,range,average,median')		
-    grass.message( "Upload aspect zonal statistics done." )
-    grass.message( "----" )
-
-    # Do some simple checks	regarding aspect range
-    grass.message( "Do some simple checks regarding aspect range and populate it to the attribute table..." )	
-    grass.message( "aspect range 100-200 *" )
-    grass.message( "aspect range 201-300 **" )
-    grass.message( "aspect range >= 300 ***" )
-    grass.run_command("v.db.addcolumn", map = v_habitat,
-                                     layer = 1, 
-                                     columns = "%s varchar(3)" % (prefix+'_check_aspect_range'))
-
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='*' WHERE %s < 200 AND %s >= 100" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range', prefix+'_aspect_range'))									 
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='**' WHERE %s < 300 AND %s >= 200" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range', prefix+'_aspect_range'))									 
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='***' WHERE %s >= 300" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range'))								 
-
-									 
-    grass.message( "Simple checks regarding aspect range done." )
-    grass.message( "----" )
-
-    # Do some simple checks	regarding aspect and and slope
-    grass.message( "Do some simple checks regarding aspect range and slope median and populate it to the attribute table..." )		
-    grass.message( "aspect range 100-200 and median slope < 5 *" )
-    grass.message( "aspect range 201-300 and median slope < 5 **" )
-    grass.message( "aspect range >= 300 and median slope < 5 ***" )
-    grass.run_command("v.db.addcolumn", map = v_habitat,
-                                     layer = 1, 
-                                     columns = "%s varchar(3)" % (prefix+'_check_aspect_slope'))	
-
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='*' WHERE (%s < 200 AND %s >= 100) AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_aspect_range', prefix+'_slope_median'))									 
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='**' WHERE (%s < 300 AND %s >= 200) AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_aspect_range', prefix+'_slope_median'))									 
-    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='***' WHERE %s >= 300 AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_slope_median'))
-
-    grass.message( "Simple checks regarding aspect range and median slope done." )
-    grass.message( "----" )
-
-    # Do some simple checks	regarding aspect and and slope
-    grass.message( "Convert habitat vector to raster ..." )
-    grass.run_command("v.to.rast", input = v_habitat,
-                                     layer = 1,
-                                     type = 'area',
-                                     use = 'attr',
-                                     attrcolumn = v_column,
-                                     output = r_habitat)
-    grass.message( "Conversion done." )
-    grass.message( "----" )
-
-    # Calculate the most common geomorphons
-    grass.message( "Calculate the most common geomorphons ..." )
-    grass.run_command("r.geomorphon", dem = r_elevation,
-                                     skip = 0,
-                                     search = 3,
-                                     flat = 1,
-                                     dist = 0,
-                                     forms = r_geomorphon)
-    grass.message( "Geomorphon calculations done." )
-    grass.message( "----" )
-
-    # Mutual occurrence (coincidence) of categories of habitats and geomorphons
-    grass.message( "Calculate mutual occurrences of habitats and geomorphons ..." )
-    grass.message( "1 - flat" )
-    grass.message( "2 - summit" )
-    grass.message( "3 - ridge" )
-    grass.message( "4 - shoulder" )
-    grass.message( "5 - spur" )
-    grass.message( "6 - slope" )
-    grass.message( "7 - hollow" )
-    grass.message( "8 - footslope" )
-    grass.message( "9 - valley" )
-    grass.message( "10 - depression" )
-    grass.message( " " )
-    grass.message( "Mutual occurrence in percent of the row" )	
-    grass.run_command("r.coin", map1 = r_habitat,
-                                     map2 = r_geomorphon,
-                                     flags = 'w',
-                                     units = 'y')
-    grass.message( "Calculations of mutual occurrences done." )
-    grass.message( "----" )
-
-    # Join geomorphons to habitat attribute table	
-    grass.message( "Join geomorphon information to habitat attribute table ...." ) 
-    
-    habgeocsv = os.path.join( directory, f_habitat_geomorphons)
-	
-    grass.run_command("r.stats", input = [r_habitat, r_geomorphon],
-                                     flags = 'aln',
-                                     separator = ';',
-                                     output = habgeocsv)	
-
-    grass.run_command("db.in.ogr", dsn = habgeocsv,
-                                     output = t_habitat_geomorphons)
-
-    grass.run_command("db.dropcolumn", table = t_habitat_geomorphons,
-                                     column = 'field_2',
-                                     flags = 'f')									 
-									 
-    grass.run_command("db.dropcolumn", table = t_habitat_geomorphons,
-                                     column = 'field_3',
-                                     flags = 'f')									 
-    								 
-    habgeocsv_pivot = os.path.join( directory, f_habitat_geomorphons_pivot)
-
-    grass.run_command("db.select", separator = ';',
-                                     output = habgeocsv_pivot,	
-                                     sql = "SELECT field_1, sum(case when field_4 = 'flat' then field_5 end) as flat, sum(case when field_4 = 'summit' then field_5 end) as summit, sum(case when field_4 = 'ridge' then field_5 end) as ridge, sum(case when field_4 = 'shoulder' then field_5 end) as shoulder, sum(case when field_4 = 'spur' then field_5 end) as spur, sum(case when field_4 = 'slope' then field_5 end) as slope, sum(case when field_4 = 'hollow' then field_5 end) as hollow, sum(case when field_4 = 'footslope' then field_5 end) as footslope, sum(case when field_4 = 'valley' then field_5 end) as valley, sum(case when field_4 = 'depression' then field_5 end) as depression , sum(field_5) as SubTotal FROM %s GROUP BY field_1" % t_habitat_geomorphons)
-
-    grass.run_command("db.in.ogr", dsn = habgeocsv_pivot,
-                                     output = t_habitat_geomorphons_pivot)
-
-    grass.run_command("v.db.join", map = v_habitat,
-                                     column = v_column,
-                                     otable = t_habitat_geomorphons_pivot,
-                                     ocolumn = 'field_1')
-    grass.message( "Geomorphon information joint to habitat attribute table." )
-    grass.message( "----" )
-
-    # Give information where output files are									 
-    grass.message( "Geomorphon information:" )
-    grass.message( habgeocsv )
-    grass.message( "Geomorphon information in pivot format:" )
-    grass.message( habgeocsv_pivot )
-    grass.message( "----" )
-
-    # Calculate LS factor see Neteler & Mitasova 2008. Open Source GIS - A GRASS GIS Approach									 
-    grass.message( "Calculate LS factor ..." )
-    grass.run_command("r.flow", elevation = r_elevation, 
-                                     aspect = r_aspect, 
-                                     flowaccumulation = r_flow_accum)
-
-    grass.message( "..." )									 
-    grass.mapcalc("$outmap = 1.4 * exp($flowacc / 22.1, 0.4) * exp(sin($slope) / 0.09, 1.2)",
-                                     outmap = r_LS, 
-                                     flowacc = r_flow_accum, 
-                                     slope = r_slope)
-
-    # create and define color rules file for LS factor map
-    ls_color_rules_out = os.path.join( directory, f_LS_colorrules )
-    with open(ls_color_rules_out, 'w') as f:
-         writer = csv.writer(f)
-         writer.writerow(['0 white'])
-         writer.writerow(['3 yellow'])
-         writer.writerow(['6 orange'])
-         writer.writerow(['10 red'])
-         writer.writerow(['50 magenta'])
-         writer.writerow(['100 violet'])	
-
-    grass.run_command("r.colors", map = r_LS, 
-                                     rules = ls_color_rules_out)		 
-
-    grass.message( "Calculation LS factor done." )
-    grass.message( "----" )
-	
-    # Run r.sun.hourly in binary mode for light/shadow
-    grass.message( "Run r.sun.hourly in binary mode for light/shadow for a certain day in the year ..." )	
-    grass.run_command("r.sun.hourly", elev_in = r_elevation,
-                                     flags = 'tb',
-                                     asp_in = r_aspect,
-                                     slope_in = r_slope,
-                                     start_time = d_start_time,
-                                     end_time = d_end_time,
-                                     day = d_day,
-                                     year = d_year,
-                                     beam_rad_basename = r_beam_rad)
-
-    grass.message( "----" )									 
-    grass.message( "Light/shadow conditions calculated for year" )
-    grass.message( d_year )	
-    grass.message( "and day" )
-    grass.message( d_day )
-    grass.message( 'from' )	
-    grass.message( d_start_time )	
-    grass.message( 'to' )
-    grass.message( d_end_time )	
-    grass.message( 'done.' )
-    grass.message( "----" )								 
-    grass.run_command("t.info", flags = 'h',
-                                     input = r_beam_rad)
-    grass.message( "----" )
-	
-    # Set region to original 
-    grass.message( "Restore original region settings:" )	
-    grass.run_command("g.region", flags = 'p', 
-                                     region = saved_region)
-    grass.message( "----" )	
-	
-    # clean up some temporay files and maps
-    grass.message( "Some clean up ..." )	
-    grass.run_command("g.remove", region = saved_region)
-    grass.run_command("g.remove", rast = r_flow_accum)
-    grass.run_command("g.remove", rast = r_habitat)	
-    grass.run_command("db.droptable", flags = 'f', 
-                                     table = t_habitat_geomorphons)
-    grass.run_command("db.droptable", flags = 'f', 
-                                     table = t_habitat_geomorphons_pivot)
-    grass.run_command("db.dropcolumn", flags = 'f',
-                                     table = v_habitat,	
-                                     column = 'field_1')
-    grass.message( "Clean up done." )
-    grass.message( "----" )	
-
-    # v.habitat.dem done!	
-    grass.message( "v.habitat.dem done!" )	
-
-if __name__ == "__main__":
-    options, flags = grass.parser()
-    sys.exit(main())
+#!/usr/bin/env python
+
+"""
+MODULE:    v.habitat.dem
+
+AUTHOR(S): Helmut Kudrnovsky <alectoria AT gmx at>
+
+PURPOSE:   DEM and solar derived characteristics of habitats
+
+COPYRIGHT: (C) 2014 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.
+"""
+
+#%module
+#% description: DEM derived characteristics of habitats
+#% keywords: vector
+#% keywords: raster
+#% keywords: terrain
+#% keywords: statistics
+#% keywords: sun
+#% keywords: zonal statistics
+#%end
+
+#%option G_OPT_R_ELEV
+#% key: elevation
+#% description: Name of elevation raster map 
+#% required: yes
+#%end
+
+#%option G_OPT_V_INPUT
+#% key: vector
+#% description: Name of habitat vector map 
+#% required: yes
+#%end
+
+#%option G_OPT_DB_COLUMN
+#% description: Name of attribute column with a unique habitat ID (must be numeric)
+#% required: yes
+#%end
+
+#%option
+#% key: prefix
+#% type: string
+#% key_desc: prefix
+#% description: output prefix (must start with a letter)
+#% required: yes
+#%end
+
+#%option G_OPT_M_DIR
+#% key: dir
+#% description: Directory where the output will be found
+#% required : yes
+#%end
+
+#%option
+#% key: region_extension
+#% type: double
+#% key_desc: float
+#% description: region extension
+#% required : no
+#% answer: 5000
+#%end
+
+#%option
+#% key: start_time
+#% type: double
+#% label: Start time of interval
+#% description: Use up to 2 decimal places
+#% options: 0-24
+#% answer: 8
+#%end
+
+#%option
+#% key: end_time
+#% type: double
+#% label: End time of interval
+#% description: Use up to 2 decimal places
+#% options: 0-24
+#% answer: 18
+#%end
+
+#%option
+#% key: time_step
+#% type: double
+#% label: Time step for running r.sun [decimal hours]
+#% description: Use up to 2 decimal places
+#% options: 0-24
+#% answer: 1
+#%end
+
+#%option
+#% key: day
+#% type: integer
+#% description: No. of day of the year
+#% options: 1-365
+#% answer: 172
+#%end
+
+#%option
+#% key: year
+#% type: integer
+#% label: Year used for map registration into temporal dataset or r.timestamp
+#% description: This value is not used in r.sun calculations
+#% options: 1900-9999
+#% answer: 2014
+#%end
+
+import sys
+import os
+import grass.script as grass
+import math
+from numpy import array
+from numpy import zeros
+import csv
+
+if not os.environ.has_key("GISBASE"):
+    grass.message( "You must be in GRASS GIS to run this program." )
+    sys.exit(1)
+
+def main():
+    r_elevation = options['elevation'].split('@')[0] 
+    v_habitat = options['vector'].split('@')[0] 
+    v_column = options['column']
+    regext = options['region_extension']
+    directory = options['dir']
+    prefix = options['prefix']
+    r_accumulation = prefix+'_accumulation'
+    r_drainage = prefix+'_drainage'
+    r_tci = prefix+'_topographicindex'
+    r_slope = prefix+'_slope'
+    r_aspect = prefix+'_aspect'
+    r_flow_accum = prefix+'_flow_accumulation'
+    r_LS = prefix+'_LS_factor'
+    r_habitat = prefix+'_'+v_habitat+'_rast'
+    r_geomorphon = prefix+'_'+r_elevation+'_geomorphon'
+    r_beam_rad = prefix+'_beam_rad'
+    f_habitat_small_areas = prefix+'_small_areas.csv'
+    f_habitat_geomorphons = prefix+'_habitat_geomorphons.csv'
+    f_habitat_geomorphons_pivot = prefix+'_habitat_geomorphons_pivot.csv'	
+    f_LS_colorrules = prefix+'_LS_color.txt'
+    t_habitat_geomorphons = prefix+'_habgeo_tab'
+    t_habitat_geomorphons_pivot = prefix+'_habgeo_tab_pivot'
+    d_start_time = options['start_time']
+    d_end_time = options['end_time']
+    d_time_step = options['time_step']
+    d_day = options['day']
+    d_year = options['year']
+    saved_region = prefix+'_region_ori'
+    current_region = grass.region()
+    X = float(regext)
+    N = current_region["n"]
+    S = current_region["s"]
+    E = current_region["e"]
+    W = current_region["w"]
+    Yres = current_region['nsres']
+    Xres = current_region['ewres']	
+    PixelArea = Xres * Yres	
+    global tmp
+ 
+    ## check if r.geomorphon addon is installed
+    if not grass.find_program('r.geomorphon', '--help'):
+            grass.fatal(_("The 'r.geomorphon' addon was not found, install it first:") +
+                        "\n" +
+                        "g.extension r.geomorphon")	
+
+    ## check if r.geomorphon addon is installed
+    if not grass.find_program('r.sun.hourly', '--help'):
+            grass.fatal(_("The 'r.sun.hourly' addon was not found, install it first:") +
+                        "\n" +
+                        "g.extension r.sun.hourly")	
+						
+    # Region settings
+    grass.message( "Current region will be saved and extended." )	
+    grass.message( "----" )
+	
+    # Print and save current region
+    grass.run_command('g.region', flags = 'p', save = saved_region, overwrite = True)
+    grass.message( "Current region saved." )	
+    grass.message( "----" )
+	
+    # Align region to elevation raster and habitat vector
+    # grass.message( "Align region to elevation raster and habitat vector ..." )		
+    # grass.run_command('g.region', flags = 'ap',
+    #                                 rast = r_elevation,
+    #                                 vect = v_habitat,
+    #                                 align = r_elevation)
+    #grass.message( "Alignment done." )	
+    #grass.message( "----" )									 
+									 
+    # Extend region
+    grass.message( "Extend region ..." )
+    grass.message( "n, s, e, w" )
+    grass.message( [current_region[key] for key in "nsew"] )
+    grass.message( "by" )
+    grass.message( regext )
+    grass.message( "in all directions" )
+	
+    grass.run_command('g.region', n = N+X,
+                                     s = S-X,
+                                     e = E+X,
+                                     w = W-X)
+    grass.message( "Region extension done." )
+	
+    grass.message( "Extended region:" )
+    grass.run_command('g.region', flags = 'p')
+    grass.message( "----" )
+	
+    # Watershed calculation: accumulation, drainage direction, topographic index
+    grass.message( "Calculation of accumulation, drainage direction, topographic index by r.watershed ..." )
+    grass.run_command('r.watershed', elevation = r_elevation, 
+                                     accumulation = r_accumulation, 
+                                     drainage = r_drainage,
+                                     tci = r_tci, 									 
+                                     convergence = 5, 
+                                     flags = 'am',
+                                     overwrite = True)
+    grass.message( "Calculation of accumulation, drainage direction, topographic done." )
+    grass.message( "----" )
+									 
+    # Calculation of slope and aspect maps
+    grass.message( "Calculation of slope and aspect by r.slope.aspect ..." )
+    grass.run_command('r.slope.aspect', elevation = r_elevation,
+                                     slope = r_slope, 
+                                     aspect = r_aspect,
+                                     overwrite = True)
+    grass.message( "Calculation of slope and aspect done." )
+    grass.message( "----" )
+ 
+    # Calculate pixel area by nsres x ewres
+    grass.message( "Pixel area:" )	
+    grass.message( PixelArea )	
+    grass.message( "----" )
+
+    # Calculate habitat area and populate it to the attribute table
+    grass.message( "Calculate habitat's areas and populate it to the attribute table ..." )	
+    grass.run_command("v.db.addcolumn", map = v_habitat,
+                                     layer = 1, 
+                                     columns = "habarea double")		
+	
+    grass.run_command("v.to.db", map = v_habitat,
+                                     option = 'area',
+                                     layer = 1, 
+                                     columns = 'habarea',
+                                     overwrite = True)
+
+    grass.message( "Calculate habitat's areas done." )	
+    grass.message( "----" )	 	
+
+    # Show habitat areas smaller than Pixel Area
+    grass.message( "Habitat areas smaller than pixel area." )	
+    grass.run_command("v.db.select", map = v_habitat,
+                                    flags = 'v',
+                                    layer = 1,
+                                    columns = v_column,
+                                    where = "habarea < %s" % (PixelArea))
+
+    smallareacsv = os.path.join( directory, f_habitat_small_areas )								
+									
+    grass.run_command("v.db.select", map = v_habitat,
+                                    flags = 'v',
+                                    layer = 1,
+                                    columns = v_column,
+                                    file = smallareacsv,
+                                    where = "habarea < %s" % (PixelArea))
+
+    grass.message( "A list of habitat areas smaller than pixel area can be found in: " )
+    grass.message( smallareacsv )	
+    grass.message( "----" )	    
+
+    # Mark habitats smaller than pixel area in attribute table
+    grass.message( "Mark habitats smaller than pixel area in attribute table ..." )	
+    grass.run_command("v.db.addcolumn", map = v_habitat,
+                                     layer = 1, 
+                                     columns = "%s_smallarea varchar(1)" % prefix )
+    grass.run_command("v.db.update", map = v_habitat, layer = 1, column = '%s_smallarea' % (prefix), value = '*', where = 'habarea < %s' % (PixelArea))
+    grass.message( "See column" )
+    grass.message( '%s_smallarea' % prefix )
+    grass.message( "marked by *." )
+    grass.message( "----" )	  
+
+    # Upload DEM zonal statistics to the attribute table
+    grass.message( "Upload DEM zonal statistics to the attribute table ..." )
+    grass.run_command("v.rast.stats", map = v_habitat,
+                                     flags = 'c',
+                                     layer = 1,
+                                     raster = r_elevation,
+                                     column_prefix = prefix+'_dem',
+                                     method = 'minimum,maximum,range,average,median')		
+    grass.message( "Upload DEM zonal statistics done." )		
+    grass.message( "----" )
+	
+    # Upload slope zonal statistics to the attribute table
+    grass.message( "Upload slope zonal statistics to the attribute table ..." )		
+    grass.run_command("v.rast.stats", map = v_habitat,
+                                     flags = 'c',
+                                     layer = 1,
+                                     raster = r_slope,
+                                     column_prefix = prefix+'_slope',
+                                     method = 'minimum,maximum,range,average,median')		
+    grass.message( "Upload slope zonal statistics done." )
+    grass.message( "----" )
+
+    # Upload slope zonal statistics to the attribute table
+    grass.message( "Upload aspect zonal statistics to the attribute table ..." )		
+    grass.run_command("v.rast.stats", map = v_habitat,
+                                     flags = 'c',
+                                     layer = 1,
+                                     raster = r_aspect,
+                                     column_prefix = prefix+'_aspect',
+                                     method = 'minimum,maximum,range,average,median')		
+    grass.message( "Upload aspect zonal statistics done." )
+    grass.message( "----" )
+
+    # Do some simple checks	regarding aspect range
+    grass.message( "Do some simple checks regarding aspect range and populate it to the attribute table..." )	
+    grass.message( "aspect range 100-200 *" )
+    grass.message( "aspect range 201-300 **" )
+    grass.message( "aspect range >= 300 ***" )
+    grass.run_command("v.db.addcolumn", map = v_habitat,
+                                     layer = 1, 
+                                     columns = "%s varchar(3)" % (prefix+'_check_aspect_range'))
+
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='*' WHERE %s < 200 AND %s >= 100" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range', prefix+'_aspect_range'))									 
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='**' WHERE %s < 300 AND %s >= 200" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range', prefix+'_aspect_range'))									 
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='***' WHERE %s >= 300" % (v_habitat, prefix+'_check_aspect_range', prefix+'_aspect_range'))								 
+
+									 
+    grass.message( "Simple checks regarding aspect range done." )
+    grass.message( "----" )
+
+    # Do some simple checks	regarding aspect and and slope
+    grass.message( "Do some simple checks regarding aspect range and slope median and populate it to the attribute table..." )		
+    grass.message( "aspect range 100-200 and median slope < 5 *" )
+    grass.message( "aspect range 201-300 and median slope < 5 **" )
+    grass.message( "aspect range >= 300 and median slope < 5 ***" )
+    grass.run_command("v.db.addcolumn", map = v_habitat,
+                                     layer = 1, 
+                                     columns = "%s varchar(3)" % (prefix+'_check_aspect_slope'))	
+
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='*' WHERE (%s < 200 AND %s >= 100) AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_aspect_range', prefix+'_slope_median'))									 
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='**' WHERE (%s < 300 AND %s >= 200) AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_aspect_range', prefix+'_slope_median'))									 
+    grass.run_command("db.execute", sql = "UPDATE %s SET %s ='***' WHERE %s >= 300 AND %s < 5" % (v_habitat, prefix+'_check_aspect_slope', prefix+'_aspect_range', prefix+'_slope_median'))
+
+    grass.message( "Simple checks regarding aspect range and median slope done." )
+    grass.message( "----" )
+
+    # Do some simple checks	regarding aspect and and slope
+    grass.message( "Convert habitat vector to raster ..." )
+    grass.run_command("v.to.rast", input = v_habitat,
+                                     layer = 1,
+                                     type = 'area',
+                                     use = 'attr',
+                                     attrcolumn = v_column,
+                                     output = r_habitat)
+    grass.message( "Conversion done." )
+    grass.message( "----" )
+
+    # Calculate the most common geomorphons
+    grass.message( "Calculate the most common geomorphons ..." )
+    grass.run_command("r.geomorphon", dem = r_elevation,
+                                     skip = 0,
+                                     search = 3,
+                                     flat = 1,
+                                     dist = 0,
+                                     forms = r_geomorphon)
+    grass.message( "Geomorphon calculations done." )
+    grass.message( "----" )
+
+    # Mutual occurrence (coincidence) of categories of habitats and geomorphons
+    grass.message( "Calculate mutual occurrences of habitats and geomorphons ..." )
+    grass.message( "1 - flat" )
+    grass.message( "2 - summit" )
+    grass.message( "3 - ridge" )
+    grass.message( "4 - shoulder" )
+    grass.message( "5 - spur" )
+    grass.message( "6 - slope" )
+    grass.message( "7 - hollow" )
+    grass.message( "8 - footslope" )
+    grass.message( "9 - valley" )
+    grass.message( "10 - depression" )
+    grass.message( " " )
+    grass.message( "Mutual occurrence in percent of the row" )	
+    grass.run_command("r.coin", map1 = r_habitat,
+                                     map2 = r_geomorphon,
+                                     flags = 'w',
+                                     units = 'y')
+    grass.message( "Calculations of mutual occurrences done." )
+    grass.message( "----" )
+
+    # Join geomorphons to habitat attribute table	
+    grass.message( "Join geomorphon information to habitat attribute table ...." ) 
+    
+    habgeocsv = os.path.join( directory, f_habitat_geomorphons)
+	
+    grass.run_command("r.stats", input = [r_habitat, r_geomorphon],
+                                     flags = 'aln',
+                                     separator = ';',
+                                     output = habgeocsv)	
+
+    grass.run_command("db.in.ogr", dsn = habgeocsv,
+                                     output = t_habitat_geomorphons)
+
+    grass.run_command("db.dropcolumn", table = t_habitat_geomorphons,
+                                     column = 'field_2',
+                                     flags = 'f')									 
+									 
+    grass.run_command("db.dropcolumn", table = t_habitat_geomorphons,
+                                     column = 'field_3',
+                                     flags = 'f')									 
+    								 
+    habgeocsv_pivot = os.path.join( directory, f_habitat_geomorphons_pivot)
+
+    grass.run_command("db.select", separator = ';',
+                                     output = habgeocsv_pivot,	
+                                     sql = "SELECT field_1, sum(case when field_4 = 'flat' then field_5 end) as flat, sum(case when field_4 = 'summit' then field_5 end) as summit, sum(case when field_4 = 'ridge' then field_5 end) as ridge, sum(case when field_4 = 'shoulder' then field_5 end) as shoulder, sum(case when field_4 = 'spur' then field_5 end) as spur, sum(case when field_4 = 'slope' then field_5 end) as slope, sum(case when field_4 = 'hollow' then field_5 end) as hollow, sum(case when field_4 = 'footslope' then field_5 end) as footslope, sum(case when field_4 = 'valley' then field_5 end) as valley, sum(case when field_4 = 'depression' then field_5 end) as depression , sum(field_5) as SubTotal FROM %s GROUP BY field_1" % t_habitat_geomorphons)
+
+    grass.run_command("db.in.ogr", dsn = habgeocsv_pivot,
+                                     output = t_habitat_geomorphons_pivot)
+
+    grass.run_command("v.db.join", map = v_habitat,
+                                     column = v_column,
+                                     otable = t_habitat_geomorphons_pivot,
+                                     ocolumn = 'field_1')
+    grass.message( "Geomorphon information joint to habitat attribute table." )
+    grass.message( "----" )
+
+    # Give information where output files are									 
+    grass.message( "Geomorphon information:" )
+    grass.message( habgeocsv )
+    grass.message( "Geomorphon information in pivot format:" )
+    grass.message( habgeocsv_pivot )
+    grass.message( "----" )
+
+    # Calculate LS factor see Neteler & Mitasova 2008. Open Source GIS - A GRASS GIS Approach									 
+    grass.message( "Calculate LS factor ..." )
+    grass.run_command("r.flow", elevation = r_elevation, 
+                                     aspect = r_aspect, 
+                                     flowaccumulation = r_flow_accum)
+
+    grass.message( "..." )									 
+    grass.mapcalc("$outmap = 1.4 * exp($flowacc / 22.1, 0.4) * exp(sin($slope) / 0.09, 1.2)",
+                                     outmap = r_LS, 
+                                     flowacc = r_flow_accum, 
+                                     slope = r_slope)
+
+    # create and define color rules file for LS factor map
+    ls_color_rules_out = os.path.join( directory, f_LS_colorrules )
+    with open(ls_color_rules_out, 'w') as f:
+         writer = csv.writer(f)
+         writer.writerow(['0 white'])
+         writer.writerow(['3 yellow'])
+         writer.writerow(['6 orange'])
+         writer.writerow(['10 red'])
+         writer.writerow(['50 magenta'])
+         writer.writerow(['100 violet'])	
+
+    grass.run_command("r.colors", map = r_LS, 
+                                     rules = ls_color_rules_out)		 
+
+    grass.message( "Calculation LS factor done." )
+    grass.message( "----" )
+	
+    # Run r.sun.hourly in binary mode for light/shadow
+    grass.message( "Run r.sun.hourly in binary mode for light/shadow for a certain day in the year ..." )	
+    grass.run_command("r.sun.hourly", elev_in = r_elevation,
+                                     flags = 'tb',
+                                     asp_in = r_aspect,
+                                     slope_in = r_slope,
+                                     start_time = d_start_time,
+                                     end_time = d_end_time,
+                                     day = d_day,
+                                     year = d_year,
+                                     beam_rad_basename = r_beam_rad)
+
+    grass.message( "----" )									 
+    grass.message( "Light/shadow conditions calculated for year" )
+    grass.message( d_year )	
+    grass.message( "and day" )
+    grass.message( d_day )
+    grass.message( 'from' )	
+    grass.message( d_start_time )	
+    grass.message( 'to' )
+    grass.message( d_end_time )	
+    grass.message( 'done.' )
+    grass.message( "----" )								 
+    grass.run_command("t.info", flags = 'h',
+                                     input = r_beam_rad)
+    grass.message( "----" )
+	
+    # Set region to original 
+    grass.message( "Restore original region settings:" )	
+    grass.run_command("g.region", flags = 'p', 
+                                     region = saved_region)
+    grass.message( "----" )	
+	
+    # clean up some temporay files and maps
+    grass.message( "Some clean up ..." )	
+    grass.run_command("g.remove", region = saved_region)
+    grass.run_command("g.remove", rast = r_flow_accum)
+    grass.run_command("g.remove", rast = r_habitat)	
+    grass.run_command("db.droptable", flags = 'f', 
+                                     table = t_habitat_geomorphons)
+    grass.run_command("db.droptable", flags = 'f', 
+                                     table = t_habitat_geomorphons_pivot)
+    grass.run_command("db.dropcolumn", flags = 'f',
+                                     table = v_habitat,	
+                                     column = 'field_1')
+    grass.message( "Clean up done." )
+    grass.message( "----" )	
+
+    # v.habitat.dem done!	
+    grass.message( "v.habitat.dem done!" )	
+
+if __name__ == "__main__":
+    options, flags = grass.parser()
+    sys.exit(main())


Property changes on: grass-addons/grass7/vector/v.habitat.dem/v.habitat.dem.py
___________________________________________________________________
Added: svn:mime-type
   + text/x-python
Added: svn:eol-style
   + native



More information about the grass-commit mailing list