[GRASS-SVN] r72655 - grass-addons/grass7/raster/r.valley.bottom
svn_grass at osgeo.org
svn_grass at osgeo.org
Sun Apr 29 22:26:28 PDT 2018
Author: spawley
Date: 2018-04-29 22:26:27 -0700 (Sun, 29 Apr 2018)
New Revision: 72655
Modified:
grass-addons/grass7/raster/r.valley.bottom/r.valley.bottom.py
Log:
r.valley.bottom deal with nodata values at edges as part of mapcalc expression
Modified: grass-addons/grass7/raster/r.valley.bottom/r.valley.bottom.py
===================================================================
--- grass-addons/grass7/raster/r.valley.bottom/r.valley.bottom.py 2018-04-29 22:17:12 UTC (rev 72654)
+++ grass-addons/grass7/raster/r.valley.bottom/r.valley.bottom.py 2018-04-30 05:26:27 UTC (rev 72655)
@@ -231,11 +231,7 @@
PCTL = "PCTL{L}".format(L=L+1)
TMP_RAST.append(PCTL)
-
- # cell padding of input raster
- input_grown = '_'.join(['tmp', input.split('@')[0], 'cell_padded']) + \
- ''.join([random.choice(string.ascii_letters + string.digits) for n in range(4)])
- region = cell_padding(input=input, output=input_grown, radius=radius)
+ input_grown=input
# get offsets for given neighborhood radius
offsets = focal_expr(radius=radius, window_square=window_square)
@@ -243,15 +239,20 @@
# generate grass mapcalc terms and execute
n_pixels = float(len(offsets))
terms = ["({input}[%d,%d] <= {input})".format(input=input_grown) % d for d in offsets]
+
+ # create mapcalc expr
+ # if pixel in neighborhood contains nodata, attempt to use opposite neighbor
+ # if opposite neighbor is also nodata, then use center pixel
+ terms = []
+ for d in offsets:
+ valid = ','.join(map(str, d))
+ invalid = ','.join([str(-d[0]), str(-d[1])])
+ terms.append("if(isnull({input}[{d}]), if(isnull({input}[{e}] <= {input}), {input}, {input}[{e}] <= {input}), {input}[{d}] <= {input})".format(
+ input=input_grown, d=valid, e=invalid))
+
expr = "{x} = ({s}) / {n}".format(x=PCTL, s=" + ".join(terms), n=n_pixels)
r.mapcalc(expr)
- # remove padded map
- g.remove(type='raster', name=input_grown, flags='f', quiet=True)
-
- # return to original region setting
- Region.write(region)
-
return(PCTL)
More information about the grass-commit
mailing list