[GRASS-SVN] r69855 - in grass-addons/grass7/raster3d: . r3.forestfrag r3.forestfrag/testsuite

svn_grass at osgeo.org svn_grass at osgeo.org
Sun Nov 20 14:11:34 PST 2016


Author: wenzeslaus
Date: 2016-11-20 14:11:34 -0800 (Sun, 20 Nov 2016)
New Revision: 69855

Added:
   grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r3_forestfrag_trivial.py
Removed:
   grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_ncspm.py
   grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_trivial.py
Modified:
   grass-addons/grass7/raster3d/Makefile
   grass-addons/grass7/raster3d/r3.forestfrag/Makefile
   grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.html
   grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.py
   grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_xy.py
Log:
r3.forestfrag: 3D version of Riitters et al. 2000 forest fragmentation index

Modified: grass-addons/grass7/raster3d/Makefile
===================================================================
--- grass-addons/grass7/raster3d/Makefile	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/Makefile	2016-11-20 22:11:34 UTC (rev 69855)
@@ -2,6 +2,7 @@
 
 SUBDIRS = \
 	r3.count.categories \
+	r3.r3.forestfrag \
 	r3.profile \
 	r3.scatterplot \
 	r3.what

Modified: grass-addons/grass7/raster3d/r3.forestfrag/Makefile
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/Makefile	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/Makefile	2016-11-20 22:11:34 UTC (rev 69855)
@@ -1,6 +1,6 @@
 MODULE_TOPDIR = ../..
 
-PGM = r.forestfrag
+PGM = r3.forestfrag
 
 include $(MODULE_TOPDIR)/include/Make/Script.make
 

Modified: grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.html
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.html	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.html	2016-11-20 22:11:34 UTC (rev 69855)
@@ -1,7 +1,9 @@
 <h2>DESCRIPTION</h2>
 
+TODO: Create a description specific for 3D version
+
 <em>r.forestfrag</em> Computes the forest fragmentation following 
-the methodology proposed by Riitters et al. (2000). See 
+the methodology proposed by Riitters et. al (2000). See 
 <a href="http://www.consecol.org/vol4/iss2/art3/">this article</a> for a 
 detailed explanation.
 
@@ -30,10 +32,10 @@
 
 <div class="code"><pre>
 interior:       Pf = 1.0
-patch:          Pf < 0.4
-transitional:   0.4 ≤ Pf < 0.6
-edge:           Pf ≥ 0.6 and Pf - Pff < 0
-perforated:     Pf ≥ 0.6 and Pf - Pff > 0
+patch:          Pf < 0.4
+transitional:   0.4 ≤ Pf < 0.6
+edge:           Pf ≥ 0.6 and Pf - Pff < 0
+perforated:     Pf ≥ 0.6 and Pf - Pff > 0
 undetermined:   Pf ≥ 0.6 and Pf = Pff
 </pre></div>
 
@@ -52,61 +54,19 @@
 to set the region to match the input layer.
 </ul>
 
-
-<h2>EXAMPLE</h2>
-
-In the North Carolina sample Location, set the computational region
-to match the land classification raster map:
-
-<div class="code"><pre>
-g.region raster=landclass96
-</pre></div>
-
-Then mark all cells which are forest as 1 and everything else as zero:
-
-<div class="code"><pre>
-r.mapcalc "forest = if(landclass96 == 5, 1, 0)"
-</pre></div>
-
-Use the new forest presence raster map to compute the forest
-fragmentation index with window size 7:
-
-<div class="code"><pre>
-r.forestfrag input=forest output=fragmentation window=7
-</pre></div>
-
-<center>
-<img src="r_forestfrag_window_7.png" alt="r.forestfrag output with window size 7">
-<img src="r_forestfrag_window_11.png" alt="r.forestfrag output with window size 11">
-<p><em>
-Two forest fragmentation indices with window size 7 (left) and
-11 (right) show how increasing window size increases the amount of
-edges.
-</em></p>
-</center>
-
-
 <h2>SEE ALSO</h2>
-
-<em>
-<a href="r.mapcalc.html">r.mapcalc</a>,
-<a href="r.li.html">r.li</a>
-</em>
-
-<p>
 The addon is based on the
-<a href="https://grasswiki.osgeo.org/wiki/AddOns/GRASS_6#r.forestfrag">
+<a href="http://grasswiki.osgeo.org/wiki/AddOns/GRASS_6#r.forestfrag>">
 r.forestfrag.sh</a> script, with as extra options user-defined
 moving window size, option to trim the region (by default it
 respects the region) and a better handling of no-data cells.
 
-
 <h2>AUTHORS</h2>
-<em>Paulo van Breugel</em>: Rewritten in Python, added option to set
-moving window size, and improved handling no-data cells<br> 
+<em>Paulo van Breugel<em>: Rewritten in Python, added option to set
+moving window size, and improved handling no-data cells 
 <em>Stefan Sylla</em>: Revised / updated the r.forestfrag.sh for
 GRASS GIS 6<br>
-<em>Emmanuel Sambale</em>: Author original shell script
+<em>Emmanuel Sambale</em>: Author original shell script<br>
 
 <h2>REFERENCES</h2> 
 Riitters, K., J. Wickham, R. O'Neill, B. Jones, 

Modified: grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.py
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.py	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/r3.forestfrag.py	2016-11-20 22:11:34 UTC (rev 69855)
@@ -3,12 +3,12 @@
 
 ############################################################################
 #
-# MODULE:    r.forestfrag
+# MODULE:    r3.forestfrag
 #
-# AUTHOR(S): Emmanuel Sambale (original shell version)
+# AUTHOR(S): Vaclav Petras (based on r.forestfrag, wenzeslaus gmail com)
+#            Emmanuel Sambale (original shell version)
 #            Stefan Sylla (original shell version)
 #            Paulo van Breugel (Python version, paulo at ecodiv.org)
-#            Vaclav Petras (major code clean up, wenzeslaus gmail com)
 #
 # PURPOSE:   Creates forest fragmentation index map from a
 #            forest-non-forest raster; The index map is based on
@@ -27,19 +27,20 @@
 
 #%module
 #% description: Computes the forest fragmentation index (Riitters et al. 2000)
-#% keyword: raster
+#% keyword: raster3d
 #% keyword: landscape structure analysis
+#% keyword: vegetation structure analysis
 #% keyword: forest
 #% keyword: fragmentation index
 #% keyword: Riitters
 #%end
 
-#%option G_OPT_R_INPUT
+#%option G_OPT_R3_INPUT
 #% description: Name of forest raster map (where forest=1, non-forest=0)
 #% required: yes
 #%end
 
-#%option G_OPT_R_OUTPUT
+#%option G_OPT_R3_OUTPUT
 #% required: yes
 #%end
 
@@ -50,7 +51,7 @@
 #% key_desc: number
 #% options: 3-
 #% answer : 3
-#% required: no
+#% required: yes
 #%end
 
 #%option G_OPT_R_OUTPUT
@@ -67,6 +68,39 @@
 #% required: no
 #%end
 
+#%option
+#% key: transitional_limit
+#% type: double
+#% description: transitional_limit
+#% answer: 0.6
+#% required: yes
+#%end
+
+#%option
+#% key: patch_limit
+#% type: double
+#% description: patch_limit
+#% answer: 0.4
+#% required: yes
+#%end
+
+#%option
+#% key: interior_limit
+#% type: double
+#% description: interior_limit
+#%end
+
+#%option
+#% key: color
+#% type: string
+#% label: Source raster for colorization
+#% description: Input and color_input are taken from input and color_input options respectively. The rest is computed using r.slope.aspect
+#% required: no
+#% options: sambale,riitters,perceptual
+#% descriptions: sambale;Sambale, Stefan Sylla;riitters; Riitters et. al 2000;perceptual;Perceptually uniform
+#% answer: sambale
+#%end
+
 #%flag
 #% key: r
 #% description: Set computational region to input raster map
@@ -87,13 +121,6 @@
 #% description: Trim the output map to avoid border effects
 #%end
 
-#%option
-#% key: window
-#% type: integer
-#% label: This option is deprecated, use the option size instead
-#% options: 3-
-#% required: no
-#%end
 
 import os
 import sys
@@ -102,20 +129,10 @@
 import tempfile
 import string
 import grass.script as gs
-# neutral naming for better compatibility between 2D and 3D version
-from grass.script.raster import mapcalc
+# neutral naming for better compatibility with 2D version
+from grass.script.raster3d import mapcalc3d as mapcalc
 
 
-LABELS = """\
-0 exterior
-1 patch
-2 transitional
-3 edge
-4 perforated
-5 interior
-6 undetermined
-"""
-
 COLORS_SAMBALE = """\
 0 255:255:0
 1 215:48:39
@@ -126,6 +143,37 @@
 6 145:207:96
 """
 
+COLORS_RIITTERS = """\
+0 255:255:255
+1 34:34:220
+2 153:204:255
+3 255:153:102
+4 255:255:102
+5 34:255:34
+6 230:255:0
+"""
+
+COLORS_PERCEPTUAL = """\
+0 245:244:68
+1 35:60:37
+2 172:92:80
+3 192:126:73
+4 157:182:90
+5 107:214:72
+6 195:233:82
+"""
+
+COLORS_PERCEPTUAL_WHITE = """\
+0 255:255:255
+1 172:92:80
+2 192:126:73
+3 157:182:90
+4 245:244:68
+5 107:214:72
+6 195:233:82
+"""
+
+
 # create set to store names of temporary maps to be deleted upon exit
 CLEAN_RAST = []
 
@@ -134,13 +182,13 @@
     """Remove temporary maps specified in the global list"""
     cleanrast = list(reversed(CLEAN_RAST))
     for rast in cleanrast:
-        gs.run_command("g.remove", flags="f", type="raster", name=rast,
+        gs.run_command("g.remove", flags="f", type="raster3d", name=rast,
                        quiet=True)
 
 
 def raster_exists(name):
     """Check if the raster map exists, call GRASS fatal otherwise"""
-    ffile = gs.find_file(name, element='cell')
+    ffile = gs.find_file(name, element='grid3')
     if not ffile['fullname']:
         gs.fatal(_("Raster map <%s> not found") % name)
 
@@ -167,16 +215,23 @@
     :param aggregate_op: operator used to combine all pairs together
     """
     s = max_index  # just to be brief
-    base_expr = "({m}[{a},{b}] {o} {m}[{c},{d}])"
+    base_expr = "(int({m}[{a},{b},{c}]) {o} int({m}[{d},{e},{f}]))"
     expr = []
     for j in range(-s, s + 1):
-        for i in range(-s, s):
-            expr.append(base_expr.format(
-                m=map_name, o=combine_op, a=i, b=j, c=i + 1, d=j))
+        for k in range(-s, s + 1):
+            for i in range(-s, s):
+                expr.append(base_expr.format(
+                    m=map_name, o=combine_op, a=i, b=j, c=k, d=i + 1, e=j, f=k))
+    for k in range(-s, s + 1):
+        for i in range(-s, s + 1):
+            for j in range(-s, s):
+                expr.append(base_expr.format(
+                    m=map_name, o=combine_op, a=i, b=j, c=k, d=i, e=j + 1, f=k))
     for i in range(-s, s + 1):
-        for j in range(-s, s):
-            expr.append(base_expr.format(
-                m=map_name, o=combine_op, a=i, b=j, c=i, d=j + 1))
+        for j in range(-s, s + 1):
+            for k in range(-s, s):
+                expr.append(base_expr.format(
+                    m=map_name, o=combine_op, a=i, b=j, c=k, d=i, e=j, f=k + 1))
     return aggregate_op.join(expr)
 
 
@@ -185,34 +240,29 @@
     ipl = options['input']
     raster_exists(ipl)
     opl = options['output']
-    # size option backwards compatibility with window
-    if not options['size'] and not options['window']:
-        gs.fatal(_("Required parameter <%s> not set") % 'size')
-    if options['size']:
-        wz = int(options['size'])
-    if options['window']:
-        gs.warning(_("The window option is deprecated, use the option"
-                     " size instead"))
-        wz = int(options['window'])
-    if options['size'] and options['size'] != '3' and options['window']:
-        gs.warning(_("When the obsolete window option is used, the"
-                     " new size option is ignored"))
-    if wz % 2 == 0:
-        gs.fatal(_("Please provide an odd number for the moving"
-                   " window size, not %d") % wz)
+    size  = int(options['size'])
+    if size % 2 == 0:
+        gs.fatal("Please provide an odd number for the moving window")
+
+    transitional_limit = float(options['transitional_limit'])
+    patch_limit = float(options['patch_limit'])
+    if patch_limit >= transitional_limit:
+        gs.fatal("Patch limit must be lower than transitional limit")
+    if options['interior_limit']:
+        interior_limit = float(options['interior_limit'])
+    else:
+        # we choose high value just to be sure
+        float_epsilon = 2.0e-06
+        interior_limit = float_epsilon
+    if 1 - transitional_limit <= interior_limit:
+        gs.fatal("Interior tolerance is too high or transitional limit is too high")
+    # TODO: set those using flags
+    interior_upper_test = False
+    interior_circle_test = True
+    interior_limit = "{:.7f}".format(interior_limit)
     # user wants pf or pff
     user_pf = options['pf']
     user_pff = options['pff']
-    # backwards compatibility
-    if flags['t']:
-        gs.warning(_("The -t flag is deprecated, use pf and pff options"
-                     " instead"))
-    if not user_pf and not user_pff and flags['t']:
-        user_pf = opl + '_pf'
-        user_pff = opl + '_pff'
-    elif flags['t']:
-        gs.warning(_("When pf or pff option is used, the -t flag"
-                     " is ignored"))
     flag_r = flags['r']
     flag_s = flags['s']
     clip_output = flags['a']
@@ -223,20 +273,9 @@
     # it makes sense to use it from now on (but we should reconsider)
     if flag_r:
         gs.message(_("Setting region to input map..."))
-        gs.run_command('g.region', raster=ipl, quiet=True)
+        gs.run_command('g.region', raster3d=ipl, quiet=True)
 
-    # check if map values are limited to 1 and 0
-    input_info = gs.raster_info(ipl)
-    # we know what we are doing only when input is integer
-    if input_info['datatype'] != 'CELL':
-        gs.fatal(_("The input raster map must have type CELL"
-                   " (integer)"))
-    # for integer, we just need to text min and max
-    if input_info['min'] != 0 or input_info['max'] != 1:
-        gs.fatal(_("The input raster map must be a binary raster,"
-                   " i.e. it should contain only values 0 and 1"
-                   " (now the minimum is %d and maximum is %d)")
-                 % (input_info['min'], input_info['max']))
+    # TODO: check if map values are limited to 1 and 0
 
     # computing pf values
     # let forested pixels be x and number of all pixels in moving window
@@ -244,20 +283,13 @@
 
     gs.info(_("Step 1: Computing Pf values..."))
 
-    # generate grid with pixel-value=number of forest-pixels in window
-    # generate grid with pixel-value=number of pixels in moving window:
-    tmpA2 = tmpname('tmpA01_')
-    tmpC3 = tmpname('tmpA02_')
-    gs.run_command("r.neighbors", quiet=True, input=ipl,
-                   output=[tmpA2, tmpC3], method=["sum", "count"], size=wz)
-
     # create pf map
     if user_pf:
         pf = user_pf
     else:
         pf = tmpname('tmpA03_')
-    mapcalc("$pf = if( $ipl >=0, float($tmpA2) / float($tmpC3))",
-            ipl=ipl, pf=pf, tmpA2=tmpA2, tmpC3=tmpC3)
+    gs.run_command("r3.neighbors", input=ipl, output=pf,
+                      method="average", window=(size, size, size))
 
     # computing pff values
     # Considering pairs of pixels in cardinal directions in
@@ -269,11 +301,11 @@
 
     # create copy of forest map and convert NULL to 0 (if any)
     tmpC4 = tmpname('tmpA04_')
-    gs.run_command("g.copy", raster=[ipl, tmpC4], quiet=True)
-    gs.run_command("r.null", map=tmpC4, null=0, quiet=True)
+    gs.run_command("g.copy", raster3d=[ipl, tmpC4], quiet=True)
+    gs.run_command("r3.null", map=tmpC4, null=0, quiet=True)
 
     # window dimensions
-    max_index = int((wz - 1) / 2)
+    max_index = int((size - 1) / 2)
     # number of 'forest-forest' pairs
     expr1 = pairs_expression(map_name=tmpC4, max_index=max_index,
                              combine_op='&')
@@ -285,9 +317,8 @@
         pff = user_pff
     else:
         pff = tmpname('tmpA07_')
-    # potentially this can be split and parallelized
-    mapcalc("$pff = if($ipl >= 0, float($tmpl4) / float($tmpl5))",
-            ipl=ipl, tmpl4=expr1, tmpl5=expr2, pff=pff)
+    mapcalc("$pff = float($e1) / float($e2)",
+            pff=pff, e1=expr1, e2=expr2)
 
     # computing fragmentation index
     # (a b) name, condition
@@ -303,20 +334,30 @@
 
     gs.info(_("Step 3: Computing fragmentation index..."))
 
+    # equation for the interior
+    if interior_upper_test:
+        # using abs just be be sure in case floating pf goes little bit over 1
+        interior_eq = "abs($pf - 1) < $in_limit"
+    elif interior_circle_test:
+        interior_eq = "($pff - 1)^2 + ($pf - 1)^2 < $in_limit^2"
+    else:
+        interior_eq = "$pf == 1"
+
     if clip_output:
         indexfin2 = tmpname('tmpA16_')
     else:
         indexfin2 = opl
-    mapcalc(
+    expression = (
         "eval("
         "dpf = $pf - $pff,"
+        "inter = " + interior_eq + ","  # using plus to avoid formating twice
         # individual classes
-        "patch = if($pf < 0.4, 1, 0),"
-        "transitional = if($pf >= 0.4 && $pf < 0.6, 2, 0),"
-        "edge = if($pf >= 0.6 && dpf<0,3,0),"
-        "perforated = if($pf > 0.6 && $pf < 1 && dpf > 0, 4, 0),"
-        "interior = if($pf == 1, 5, 0),"
-        "undetermined = if($pf > 0.6 && $pf < 1 && dpf == 0, 6, 0),"
+        "patch = if($pf < $pa_limit, 1, 0),"
+        "transitional = if($pf >= $pa_limit && $pf < $tr_limit, 2, 0),"
+        "edge = if($pf >= $tr_limit && not(inter) && dpf < 0, 3, 0),"
+        "perforated = if($pf > $tr_limit && not(inter) && dpf > 0, 4, 0),"
+        "interior = if(inter, 5, 0),"  # TODO: we could skip this
+        "undetermined = if($pf > $tr_limit && not(inter) && dpf == 0, 6, 0),"
         # null is considered as non-forest and we need to do it before +
         "patch = if(isnull(patch), 0, patch),"
         "transitional = if(isnull(transitional), 0, transitional),"
@@ -331,82 +372,25 @@
         ")\n"
         # mask result by non-forest (according to the input)
         # removes the nonsense data created in the non-forested areas
-        "$out = all * $binary_forest",
-        out=indexfin2, binary_forest=ipl, pf=pf, pff=pff)
+        "$out = all * $binary_forest")
+    gs.debug(expression)
+    mapcalc(expression,
+        out=indexfin2, binary_forest=ipl, pf=pf, pff=pff,
+        tr_limit=transitional_limit, pa_limit=patch_limit,
+        in_limit=interior_limit)
 
-    # shrink the region
-    if clip_output:
-        gs.use_temp_region()
-        reginfo = gs.parse_command("g.region", flags="gp")
-        nscor = max_index * float(reginfo['nsres'])
-        ewcor = max_index * float(reginfo['ewres'])
-        gs.run_command("g.region",
-                       n=float(reginfo['n']) - nscor,
-                       s=float(reginfo['s']) + nscor,
-                       e=float(reginfo['e']) - ewcor,
-                       w=float(reginfo['w']) + ewcor,
-                       quiet=True)
-        mapcalc("$opl = $if3", opl=opl, if3=indexfin2, quiet=True)
+    # TODO: create categories
 
-    # create categories
-    # TODO: parametrize classes (also in r.mapcalc, r.colors and desc)?
-    # TODO: translatable labels?
-    labels = LABELS
-    gs.write_command("r.category", quiet=True, map=opl,
-                     rules='-', stdin=labels, separator='space')
-
     # create color table
-    colors = COLORS_SAMBALE
-    gs.write_command("r.colors", map=opl, rules='-',
+    if options['color'] == 'riitters':
+        colors = COLORS_RIITTERS
+    if options['color'] == 'perceptual':
+        colors = COLORS_PERCEPTUAL
+    else:
+        colors = COLORS_SAMBALE
+    gs.write_command("r3.colors", map=opl, rules='-',
                      stdin=colors, quiet=True)
 
-    # write metadata for main layer
-    gs.run_command("r.support", map=opl,
-                   title="Forest fragmentation",
-                   source1="Based on %s" % ipl,
-                   description="Forest fragmentation index (6 classes)")
-    gs.raster_history(opl)
-
-    # write metadata for intermediate layers
-    if user_pf:
-        # pf layer
-        gs.run_command("r.support", map=pf,
-                       title="Proportion forested",
-                       units="Proportion",
-                       source1="Based on %s" % ipl,
-                       description="Proportion of pixels in the moving"
-                                   " window that is forested")
-        gs.raster_history(pf)
-
-    if user_pff:
-        # pff layer
-        unused, tmphist = tempfile.mkstemp()
-        text_file = open(tmphist, "w")
-        long_description = """\
-Proportion of all adjacent (cardinal directions only) pixel pairs that
-include at least one forest pixel for which both pixels are forested.
-It thus (roughly) estimates the conditional probability that, given a
-pixel of forest, its neighbor is also forest.
-"""
-        text_file.write(long_description)
-        text_file.close()
-        gs.run_command("r.support", map=pff,
-                       title="Conditional probability neighboring cell"
-                             " is forest",
-                       units="Proportion",
-                       source1="Based on %s" % ipl,
-                       description="Probability neighbor of forest cell"
-                                   " is forest",
-                       loadhistory=tmphist)
-        gs.raster_history(pff)
-        os.remove(tmphist)
-
-    # report fragmentation index and names of layers created
-
-    if flag_s:
-        gs.run_command("r.report", map=opl, units=["h", "p"],
-                       flags="n", page_width=50, quiet=True)
-
     gs.info(_("The following layers were created"))
     gs.info(_("The fragmentation index: %s") % opl)
     if user_pf:

Copied: grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r3_forestfrag_trivial.py (from rev 69853, grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_trivial.py)
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r3_forestfrag_trivial.py	                        (rev 0)
+++ grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r3_forestfrag_trivial.py	2016-11-20 22:11:34 UTC (rev 69855)
@@ -0,0 +1,197 @@
+#!/usr/bin/env python
+
+############################################################################
+#
+# MODULE:        r3_forestfrag_trivial
+# AUTHOR:        Vaclav Petras
+# PURPOSE:       Test using a small example computed by hand
+# COPYRIGHT:     (C) 2016 by Vaclav Petras and 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.
+#
+#############################################################################
+
+from grass.gunittest.case import TestCase
+from grass.gunittest.main import test
+
+
+FOREST = """\
+version: grass7
+order: nsbt
+north: 10
+south: 8
+east: 20
+west: 18
+top: 7
+bottom: 5
+rows: 3
+cols: 3
+levels: 3
+1 0 0 
+1 0 0 
+0 0 0 
+1 1 0 
+0 1 1 
+1 0 1 
+0 1 1 
+0 0 0 
+0 0 0 
+"""
+
+# output contains all cells but we check just the middle one
+PF = """\
+version: grass7
+order: nsbt
+north: 10
+south: 8
+east: 20
+west: 18
+top: 7
+bottom: 5
+rows: 3
+cols: 3
+levels: 3
+N N N
+N N N
+N N N
+N N N
+N 0.370 N
+N N N
+N N N
+N N N
+N N N
+"""
+
+PFF = """\
+version: grass7
+order: nsbt
+north: 10
+south: 8
+east: 20
+west: 18
+top: 7
+bottom: 5
+rows: 3
+cols: 3
+levels: 3
+N N N
+N N N
+N N N
+N N N
+N 0.235 N
+N N N
+N N N
+N N N
+N N N
+"""
+
+FRAG = """\
+version: grass7
+order: nsbt
+north: 10
+south: 8
+east: 20
+west: 18
+top: 7
+bottom: 5
+rows: 3
+cols: 3
+levels: 3
+N N N
+N N N
+N N N
+N N N
+N 1 N
+N N N
+N N N
+N N N
+N N N
+"""
+
+
+class TestForestFragTrivial(TestCase):
+    # TODO: replace by unified handing of maps
+    to_remove = []
+    forest = 'r3ff_test_forest_trivial'
+    forest_frag = 'r3ff_test_forest_frag_trivial'
+    forest_frag_ref = 'r3ff_test_forest_frag_trivial_ref'
+    pf_ref = 'r3ff_test_forest_frag_trivial_ref_pf'
+    pff_ref = 'r3ff_test_forest_frag_trivial_ref_pff'
+
+    def setUp(self):
+        self.use_temp_region()
+        self.runModule('r3.in.ascii', input='-', stdin_=FOREST,
+                       output=self.forest)
+        self.to_remove.append(self.forest)
+        self.runModule('g.region', raster_3d=self.forest)
+        self.runModule('r3.in.ascii', input='-', stdin_=PF,
+                       output=self.pf_ref,
+                       null_value='N')
+        self.to_remove.append(self.pf_ref)
+        self.runModule('r3.in.ascii', input='-', stdin_=PFF,
+                       output=self.pff_ref,
+                       null_value='N')
+        self.to_remove.append(self.pff_ref)
+        self.runModule('r3.in.ascii', input='-', stdin_=FRAG,
+                       output=self.forest_frag_ref,
+                       null_value='N')
+        self.to_remove.append(self.forest_frag_ref)
+        
+
+    def tearDown(self):
+        self.del_temp_region()
+        if self.to_remove:
+            self.runModule('g.remove', flags='f', type='raster',
+                           name=','.join(self.to_remove), verbose=True)
+
+    def test_simple_window(self):
+        pf = self.forest_frag + '_pf'
+        pff = self.forest_frag + '_pff'
+        self.assertModule('r3.forestfrag', input=self.forest,
+                          output=self.forest_frag, size=3,
+                          pf=pf, pff=pff, flags='t')
+        self.assertRaster3dExists(self.forest_frag)
+        self.to_remove.append(self.forest_frag)
+        self.assertRaster3dExists(pf)
+        self.to_remove.append(pf)
+        self.assertRaster3dExists(pff)
+        self.to_remove.append(pff)
+
+        # check the basic properties
+        ref_univar = dict(null_cells=0, cells=27)
+        self.assertRaster3dFitsUnivar(raster=self.forest_frag,
+                                    reference=ref_univar, precision=0)
+        # the null cells requirements for pf and pff are fuzzy
+        ref_univar = dict(cells=27)
+        self.assertRaster3dFitsUnivar(raster=pf,
+                                    reference=ref_univar, precision=0)
+        self.assertRaster3dFitsUnivar(raster=pff,
+                                    reference=ref_univar, precision=0)
+        self.assertRaster3dMinMax(pf, refmin=0, refmax=1)
+        self.assertRaster3dMinMax(pff, refmin=0, refmax=1)
+        self.assertRaster3dMinMax(self.forest_frag, refmin=0, refmax=6)
+
+        # to check just the middle cell
+        # the current implementation of assertRasters3dNoDifference
+        # ignores nulls, so it just works regardless of region
+        # there is no zoom for 3D
+        # self.runModule('g.region', zoom=self.pff_ref)
+
+        # check the values
+        self.assertRasters3dNoDifference(actual=pf,
+                                         reference=self.pf_ref,
+                                         precision=0.001)
+
+        self.assertRasters3dNoDifference(actual=pff,
+                                         reference=self.pff_ref,
+                                         precision=0.001)
+
+        self.assertRasters3dNoDifference(actual=self.forest_frag,
+                                         reference=self.forest_frag_ref,
+                                         precision=0)  # it's CELL type
+
+
+if __name__ == '__main__':
+    test()

Deleted: grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_ncspm.py
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_ncspm.py	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_ncspm.py	2016-11-20 22:11:34 UTC (rev 69855)
@@ -1,62 +0,0 @@
-#!/usr/bin/env python
-
-############################################################################
-#
-# MODULE:        r_forestfrag_ncspm
-# AUTHOR:        Vaclav Petras
-# PURPOSE:       Test using NC SPM data
-# COPYRIGHT:     (C) 2016 by Vaclav Petras and 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.
-#
-#############################################################################
-
-from grass.gunittest.case import TestCase
-from grass.gunittest.main import test
-import grass.script.raster as gr
-
-
-class TestForestFrag(TestCase):
-    # TODO: replace by unified handing of maps
-    to_remove = []
-    landclass = 'landclass96'  # in full NC SPM
-    forest = 'rff_test_forest'
-    forest_frag = 'rff_test_forest_frag'
-
-    def setUp(self):
-        self.use_temp_region()
-        self.runModule('g.region', raster=self.landclass)
-        gr.mapcalc("{f} = if({c} == 5, 1, 0)".format(c=self.landclass,
-                                                     f=self.forest))
-        self.to_remove.append(self.forest)
-
-    def tearDown(self):
-        self.del_temp_region()
-        if self.to_remove:
-            self.runModule('g.remove', flags='f', type='raster',
-                           name=','.join(self.to_remove), verbose=True)
-
-    def test_3x3(self):
-        self.assertModule('r.forestfrag', input=self.forest,
-                          output=self.forest_frag, window=3)
-        self.assertRasterExists(self.forest_frag)
-        self.to_remove.append(self.forest_frag)
-        self.assertRasterMinMax(self.forest_frag, refmin=0, refmax=6)
-        # we don't have a check sum test for raster, so we just test all
-        ref = dict(north=228527.25, south=215018.25,
-                   east=644971, west=629980,
-                   nsres=28.5, ewres=28.5, rows=474, cols=526,
-                   cells=249324, datatype='CELL', ncats=6)
-        self.assertRasterFitsInfo(raster=self.forest_frag,
-                                  reference=ref, precision=0.0006)
-        ref = dict(n=249323, null_cells=1, cells=249324, min=0, max=5,
-                   range=5, mean=2.2925, mean_of_abs=2.2925,
-                   stddev=2.3564, variance=5.5527, coeff_var=102.7879,
-                   sum=571572)
-        self.assertRasterFitsUnivar(raster=self.forest_frag,
-                                    reference=ref, precision=0.0006)
-
-if __name__ == '__main__':
-    test()

Deleted: grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_trivial.py
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_trivial.py	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_trivial.py	2016-11-20 22:11:34 UTC (rev 69855)
@@ -1,162 +0,0 @@
-#!/usr/bin/env python
-
-############################################################################
-#
-# MODULE:        r_forestfrag_trivial
-# AUTHOR:        Vaclav Petras
-# PURPOSE:       Test using example from the Riitters et al. 2000 paper
-# COPYRIGHT:     (C) 2016 by Vaclav Petras and 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.
-#
-#############################################################################
-
-# This test led to discovery of #3067 (r68717) which was an r.mapcalc
-# row indexing bug.
-
-from grass.gunittest.case import TestCase
-from grass.gunittest.main import test
-
-
-FOREST_RIITTERS = """\
-north: 10
-south: 8
-east: 20
-west: 18
-rows: 3
-cols: 3
-1 0 1
-1 1 0
-1 1 0
-"""
-
-FRAG_RIITTERS = """\
-north: 10
-south: 8
-east: 20
-west: 18
-rows: 3
-cols: 3
-null: N
-N N N
-N 4 N
-N N N
-"""
-
-PF_RIITTERS = """\
-north: 10
-south: 8
-east: 20
-west: 18
-rows: 3
-cols: 3
-null: N
-N N N
-N 0.67 N
-N N N
-"""
-
-PFF_RIITTERS = """\
-north: 10
-south: 8
-east: 20
-west: 18
-rows: 3
-cols: 3
-null: N
-N N N
-N 0.45 N
-N N N
-"""
-
-
-class TestForestFragTrivial(TestCase):
-    # TODO: replace by unified handing of maps
-    to_remove = []
-    forest = 'rff_test_forest_trivial'
-    forest_frag = 'rff_test_forest_frag_trivial'
-    forest_frag_ref = 'rff_test_forest_frag_trivial_ref'
-
-    def setUp(self):
-        self.use_temp_region()
-
-    def tearDown(self):
-        self.del_temp_region()
-        if 0 and self.to_remove:
-            self.runModule('g.remove', flags='f', type='raster',
-                           name=','.join(self.to_remove), verbose=True)
-
-    def forest_frag_general(self, forest, window, reference):
-        self.runModule('r.in.ascii', input='-', stdin_=forest,
-                       output=self.forest)
-        self.to_remove.append(self.forest)
-        self.runModule('g.region', raster=self.forest)
-        self.assertRasterMinMax(self.forest,
-                                refmin=0, refmax=1)
-        self.runModule('r.in.ascii', input='-', stdin_=reference,
-                       output=self.forest_frag_ref)
-        self.to_remove.append(self.forest_frag_ref)
-        pf_ref = self.forest_frag_ref + '_pf'
-        self.runModule('r.in.ascii', input='-', stdin_=PF_RIITTERS,
-                       output=pf_ref)
-        self.to_remove.append(pf_ref)
-        pff_ref = self.forest_frag_ref + '_pff'
-        self.runModule('r.in.ascii', input='-', stdin_=PFF_RIITTERS,
-                       output=pff_ref)
-        self.to_remove.append(pff_ref)
-        # just check if the reference is all right
-        theoretical_min = 0
-        theoretical_max = 6
-        self.assertRasterMinMax(self.forest_frag_ref,
-                                refmin=theoretical_min,
-                                refmax=theoretical_max)
-        ref_univar = dict(null_cells=8)
-        self.assertRasterFitsUnivar(raster=self.forest_frag_ref,
-                                    reference=ref_univar, precision=0)
-
-        # actually run the module
-        pf = self.forest_frag + '_pf'
-        pff = self.forest_frag + '_pff'
-        self.assertModule('r.forestfrag', input=self.forest,
-                          output=self.forest_frag, window=window,
-                          pf=pf, pff=pff)
-        self.assertRasterExists(self.forest_frag)
-        self.to_remove.append(self.forest_frag)
-        self.assertRasterExists(pf)
-        self.to_remove.append(pf)
-        self.assertRasterExists(pff)
-        self.to_remove.append(pff)
-
-        # check the basic properties
-        self.assertRasterMinMax(self.forest_frag,
-                                refmin=theoretical_min,
-                                refmax=theoretical_max)
-        ref_univar = dict(null_cells=0)
-        self.assertRasterFitsUnivar(raster=self.forest_frag,
-                                    reference=ref_univar, precision=0)
-
-        self.runModule('g.region', zoom=pff_ref)
-
-        # check cell by cell
-
-        self.assertRastersNoDifference(actual=pf,
-                                       reference=pf_ref,
-                                       precision=0.01)
-
-        self.assertRastersNoDifference(actual=pff,
-                                       reference=pff_ref,
-                                       precision=0.01)
-
-        self.assertRastersNoDifference(actual=self.forest_frag,
-                                       reference=self.forest_frag_ref,
-                                       precision=0)  # it's CELL type
-
-
-    def test_riitters(self):
-        self.forest_frag_general(FOREST_RIITTERS, 3, FRAG_RIITTERS)
-
-
-if __name__ == '__main__':
-    test()

Modified: grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_xy.py
===================================================================
--- grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_xy.py	2016-11-20 21:00:28 UTC (rev 69854)
+++ grass-addons/grass7/raster3d/r3.forestfrag/testsuite/r_forestfrag_xy.py	2016-11-20 22:11:34 UTC (rev 69855)
@@ -1,5 +1,7 @@
 #!/usr/bin/env python
 
+# TODO: create a 3D version of this to test r3.forestfrag
+
 ############################################################################
 #
 # MODULE:        r_forestfrag_xy



More information about the grass-commit mailing list