[GRASS-user] latest r.width.funct.py appears to be broken
Bartolomei.Chris
Bartolomei.Chris at ensco.com
Thu Nov 1 09:58:50 PDT 2018
Good afternoon Markus & Co. 😊
It appears that the latest version of r.width.funct.py is not functioning properly ... We have a few Windows10 systems here running GRASS; I'm on v7.4.0 and another person has (had) 7.2.2. She was trying to run r.basin on a 30m resolution DEM and kept getting the following error:
Traceback (most recent call last):
File "C:\Users\XXX\AppData\Roaming\GRASS7\addons/scripts/r.width.funct.py", line 135, in <module>
sys.exit(main())
File "C:\Users\XXX\AppData\Roaming\GRASS7\addons/scripts/r.width.funct.py", line 84, in main
prc[0,0] , prc[0,1] = findint(kl,0.05) , 0.05
File "C:\Users\XXX\AppData\Roaming\GRASS7\addons/scripts/r.width.funct.py", line 128, in findint
z1, z2, f1, f2 = kl[int(Xf[0])][0], kl[int(Xf[0]-1)][0], kl[int(Xf[0])][1], kl[int(Xf[0]-1)][1]
IndexError: invalid index to scalar variable.
I tried running the same tool on my 7.4.0 version and also had the problem.
Trying to troubleshoot we initially discovered that the resolution of our region wasn't uniform NS and EW ... so resetting the extents to something divisible exactly by 30 fixed the problem on the 7.4.0 system however the 7.2.2 system was still getting the errors.
We upgraded that system to the latest version of GRASS (7.4.2), rebuilt the UTM location and mapset from scratch making sure that the difference of the extents were divisible by 30 and reinstalled the addons (r.hypso, r.width.funct, and all of the r.basin.* addons) for the whole system. We found that for some reason r.basin wants r.hypso and r.width.funct installed locally so we reinstalled those two locally as well.
We imported the DEM (a geoTIFF) same as we did on the 7.4.0 system, resampled it to a smaller region which was also adjusted to something with extents divisible by 30 same as on the 7.4.0 system.
We then ran r.basin with the outlet coordinates and while the 7.4.0 system works fine, he now 7.4.2 system still has the same errors as listed above.
I looked into r.width.funct.py as it was the module with the error and noticed the copy I have installed on the 7.4.0 system is different than the one installed on the 7.4.2 system which should be newer as it's a fresh installation of GRASS. I copied over r.width.funct.py from the 7.4.0 machine to the .4.2 machine and guess what ... it works now.
Based on cause and effect, it appears r.width.funct.py has an issue. I hadn't delved into the code itself but I'm including it here (it's not too long) ... this is the one that works:
I'm not sure what you would like to do but if other users are having a problem, then replace your r.width.funct.py with this and you should be okay.
😊
Chris
#!/usr/bin/env python
################################################################################
#
# MODULE: r.width.funct
#
# AUTHOR(S): Massimo Di Stefano, Francesco Di Stefano, Margherita Di Leo
#
# PURPOSE: The module produces the Width Function of a basin. The Width
# Function W(x) gives the number of the cells in a basin at a
# flow distance x from the outlet (distance-area function)
#
# COPYRIGHT: (c) 2010 Massimo Di Stefano, Francesco Di Stefano, Margherita Di Leo
#
# This program is free software under the GNU General Public
# License (>=v2). Read the file COPYING that comes with GRASS
# for details.
#
# REQUIRES: Matplotlib
# http://matplotlib.sourceforge.net/
#
#
################################################################################
#%module
#% description: Calculates the Width Function of a watershed basin.
#% keyword: raster
#% keyword: hydrology
#%end
#%option G_OPT_R_INPUT
#% key: map
#% description: Distance to outlet map (from r.stream.distance)
#% required: yes
#%end
#%option G_OPT_F_OUTPUT
#% key: image
#% key_desc: image
#% description: Name for output graph file (png)
#% required: yes
#%END
import sys
import os
import matplotlib #required by windows
matplotlib.use('wx') #required by windows
import matplotlib.pyplot as plt
import grass.script as grass
import numpy as np
def main():
stats = grass.read_command('r.stats', input = options['map'], sep = 'space', nv = '*', nsteps = '255', flags = 'Anc').split('\n')[:-1]
# res = cellsize
res = grass.region()['nsres']
zn = np.zeros((len(stats),4),float)
kl = np.zeros((len(stats),2),float)
prc = np.zeros((9,2),float)
for i in range(len(stats)):
if i == 0:
zn[i,0], zn[i,1] = map(float, stats[i].split(' '))
zn[i,1] = zn[i,1]
zn[i,2] = zn[i,1] * res
if i != 0:
zn[i,0], zn[i,1] = map(float, stats[i].split(' '))
zn[i,2] = zn[i,1] + zn[i-1,2]
zn[i,3] = zn[i,1] * (res**2)
totcell = sum(zn[:,1])
print "Tot. cells", totcell
totarea = totcell * (res**2)
print "Tot. area", totarea
maxdist = max(zn[:,0])
print "Max distance", maxdist
for i in range(len(stats)):
kl[i,0] = zn[i,0]
kl[i,1] = zn[i,2] / totcell
# quantiles
prc[0,0] , prc[0,1] = findint(kl,0.05) , 0.05
prc[1,0] , prc[1,1] = findint(kl,0.15) , 0.15
prc[2,0] , prc[2,1] = findint(kl,0.3) , 0.3
prc[3,0] , prc[3,1] = findint(kl,0.4) , 0.4
prc[4,0] , prc[4,1] = findint(kl,0.5) , 0.5
prc[5,0] , prc[5,1] = findint(kl,0.6) , 0.6
prc[6,0] , prc[6,1] = findint(kl,0.7) , 0.7
prc[7,0] , prc[7,1] = findint(kl,0.85) , 0.85
prc[8,0] , prc[8,1] = findint(kl,0.95) , 0.95
# plot
plotImage(zn[:,0], zn[:,3], options['image']+'_width_function.png','-','x','W(x)','Width Function')
print "==========================="
print "Width Function | quantiles"
print "==========================="
print '%.0f' %findint(kl,0.05) , "|", 0.05
print '%.0f' %findint(kl,0.15) , "|", 0.15
print '%.0f' %findint(kl,0.3) , "|", 0.3
print '%.0f' %findint(kl,0.4) , "|", 0.4
print '%.0f' %findint(kl,0.5) , "|", 0.5
print '%.0f' %findint(kl,0.6) , "|", 0.6
print '%.0f' %findint(kl,0.7) , "|", 0.7
print '%.0f' %findint(kl,0.85) , "|", 0.85
print '%.0f' %findint(kl,0.95) , "|", 0.95
print '\n'
print 'Done!'
def plotImage(x,y,image,type,xlabel,ylabel,title):
plt.plot(x, y, type)
plt.ylabel(ylabel)
plt.xlabel(xlabel)
plt.xlim( min(x), max(x) )
plt.ylim( min(y), max(y) )
plt.title(title)
plt.grid(True)
plt.savefig(image)
plt.close('all')
def findint(kl,f):
Xf = np.abs(kl-f); Xf = np.where(Xf==Xf.min())
#z1, z2, f1, f2 = kl[float(Xf[0])][0], kl[float(Xf[0]-1)][0], kl[float(Xf[0])][1], kl[float(Xf[0]-1)][1]
z1, z2, f1, f2 = kl[int(Xf[0])][0], kl[int(Xf[0]-1)][0], kl[int(Xf[0])][1], kl[int(Xf[0]-1)][1]
z = z1 + ((z2 - z1) / (f2 - f1)) * (f - f1)
return z
if __name__ == "__main__":
options, flags = grass.parser()
sys.exit(main())
More information about the grass-user
mailing list