|
|
@@ -40,7 +40,7 @@
|
|
|
#%end
|
|
|
#%option G_OPT_V_FIELD
|
|
|
#%end
|
|
|
-#%option G_OPT_R_INPUT
|
|
|
+#%option G_OPT_R_INPUTS
|
|
|
#% key: raster
|
|
|
#% description: Name of input raster map to calculate statistics from
|
|
|
#%end
|
|
|
@@ -49,6 +49,7 @@
|
|
|
#% type: string
|
|
|
#% description: Column prefix for new attribute columns
|
|
|
#% required : yes
|
|
|
+#% multiple: yes
|
|
|
#%end
|
|
|
#%option
|
|
|
#% key: method
|
|
|
@@ -99,8 +100,8 @@ def main():
|
|
|
|
|
|
nuldev = open(os.devnull, 'w')
|
|
|
|
|
|
- raster = options['raster']
|
|
|
- colprefix = options['column_prefix']
|
|
|
+ rasters = options['raster'].split(',')
|
|
|
+ colprefixes = options['column_prefix'].split(',')
|
|
|
vector = options['map']
|
|
|
layer = options['layer']
|
|
|
percentile = options['percentile']
|
|
|
@@ -129,20 +130,27 @@ def main():
|
|
|
# we need this for non-DBF driver:
|
|
|
dbfdriver = fi['driver'] == 'dbf'
|
|
|
|
|
|
+ # colprefix for every raster map?
|
|
|
+ if len(colprefixes) != len(rasters):
|
|
|
+ grass.fatal(_("Number of raster maps ({0}) different from \
|
|
|
+ number of column prefixes ({1})". format(len(rasters),
|
|
|
+ len(colprefixes))))
|
|
|
+
|
|
|
vector = vs[0]
|
|
|
|
|
|
rastertmp = "%s_%s" % (vector, tmpname)
|
|
|
|
|
|
- # check the input raster map
|
|
|
- if not grass.find_file(raster, 'cell')['file']:
|
|
|
- grass.fatal(_("Raster map <%s> not found") % raster)
|
|
|
+ for raster in rasters:
|
|
|
+ # check the input raster map
|
|
|
+ if not grass.find_file(raster, 'cell')['file']:
|
|
|
+ grass.fatal(_("Raster map <%s> not found") % raster)
|
|
|
|
|
|
# save current settings:
|
|
|
grass.use_temp_region()
|
|
|
|
|
|
# Temporarily aligning region resolution to $RASTER resolution
|
|
|
# keep boundary settings
|
|
|
- grass.run_command('g.region', align=raster)
|
|
|
+ grass.run_command('g.region', align=rasters[0])
|
|
|
|
|
|
# prepare base raster for zonal statistics
|
|
|
try:
|
|
|
@@ -170,6 +178,29 @@ def main():
|
|
|
if number < 1:
|
|
|
grass.fatal(_("No categories found in raster map"))
|
|
|
|
|
|
+ # Check if all categories got converted
|
|
|
+ # Report categories from vector map
|
|
|
+ vect_cats = grass.read_command('v.category', input=vector, option='report',
|
|
|
+ flags='g').rstrip('\n').split('\n')
|
|
|
+
|
|
|
+ # get number of all categories in selected layer
|
|
|
+ for vcl in vect_cats:
|
|
|
+ if vcl.split(' ')[0] == layer and vcl.split(' ')[1] == 'all':
|
|
|
+ vect_cats_n = int(vcl.split(' ')[2])
|
|
|
+
|
|
|
+ if vect_cats_n != number:
|
|
|
+ grass.warning(_("Not all vector categories converted to raster. \
|
|
|
+ Converted {0} of {1}.".format(number, vect_cats_n)))
|
|
|
+
|
|
|
+ # check if DBF driver used, in this case cut to 10 chars col names:
|
|
|
+ try:
|
|
|
+ fi = grass.vector_db(map=vector)[int(layer)]
|
|
|
+ except KeyError:
|
|
|
+ grass.fatal(
|
|
|
+ _('There is no table connected to this map. Run v.db.connect or v.db.addtable first.'))
|
|
|
+ # we need this for non-DBF driver:
|
|
|
+ dbfdriver = fi['driver'] == 'dbf'
|
|
|
+
|
|
|
# Find out which table is linked to the vector map on the given layer
|
|
|
if not fi['table']:
|
|
|
grass.fatal(
|
|
|
@@ -178,132 +209,135 @@ def main():
|
|
|
# replaced by user choiche
|
|
|
#basecols = ['n', 'min', 'max', 'range', 'mean', 'stddev', 'variance', 'cf_var', 'sum']
|
|
|
|
|
|
- # we need at least three chars to distinguish [mea]n from [med]ian
|
|
|
- # so colprefix can't be longer than 6 chars with DBF driver
|
|
|
- if dbfdriver:
|
|
|
- colprefix = colprefix[:6]
|
|
|
- variables_dbf = {}
|
|
|
-
|
|
|
- # by default perccol variable is used only for "variables" variable
|
|
|
- perccol = "percentile"
|
|
|
- perc = None
|
|
|
- for b in basecols:
|
|
|
- if b.startswith('p'):
|
|
|
- perc = b
|
|
|
- if perc:
|
|
|
- # namespace is limited in DBF but the % value is important
|
|
|
+ for i in xrange(len(rasters)):
|
|
|
+ raster = rasters[i]
|
|
|
+ colprefix = colprefixes[i]
|
|
|
+ # we need at least three chars to distinguish [mea]n from [med]ian
|
|
|
+ # so colprefix can't be longer than 6 chars with DBF driver
|
|
|
if dbfdriver:
|
|
|
- perccol = "per" + percentile
|
|
|
- else:
|
|
|
- perccol = "percentile_" + percentile
|
|
|
- percindex = basecols.index(perc)
|
|
|
- basecols[percindex] = perccol
|
|
|
-
|
|
|
- # dictionary with name of methods and position in "r.univar -gt" output
|
|
|
- variables = {'number': 2, 'null_cells': 3, 'minimum': 4, 'maximum': 5, 'range': 6,
|
|
|
- 'average': 7, 'stddev': 9, 'variance': 10, 'coeff_var': 11,
|
|
|
- 'sum': 12, 'first_quartile': 14, 'median': 15,
|
|
|
- 'third_quartile': 16, perccol: 17}
|
|
|
- # this list is used to set the 'e' flag for r.univar
|
|
|
- extracols = ['first_quartile', 'median', 'third_quartile', perccol]
|
|
|
- addcols = []
|
|
|
- colnames = []
|
|
|
- extstat = ""
|
|
|
- for i in basecols:
|
|
|
- # this check the complete name of out input that should be truncated
|
|
|
- for k in variables.keys():
|
|
|
- if i in k:
|
|
|
- i = k
|
|
|
- break
|
|
|
- if i in extracols:
|
|
|
- extstat = 'e'
|
|
|
- # check if column already present
|
|
|
- currcolumn = ("%s_%s" % (colprefix, i))
|
|
|
- if dbfdriver:
|
|
|
- currcolumn = currcolumn[:10]
|
|
|
- variables_dbf[currcolumn.replace("%s_" % colprefix, '')] = i
|
|
|
-
|
|
|
- colnames.append(currcolumn)
|
|
|
- if currcolumn in grass.vector_columns(vector, layer).keys():
|
|
|
- if not flags['c']:
|
|
|
- grass.fatal((_("Cannot create column <%s> (already present). ") % currcolumn) +
|
|
|
- _("Use -c flag to update values in this column."))
|
|
|
- else:
|
|
|
- if i == "n":
|
|
|
- coltype = "INTEGER"
|
|
|
+ colprefix = colprefix[:6]
|
|
|
+ variables_dbf = {}
|
|
|
+
|
|
|
+ # by default perccol variable is used only for "variables" variable
|
|
|
+ perccol = "percentile"
|
|
|
+ perc = None
|
|
|
+ for b in basecols:
|
|
|
+ if b.startswith('p'):
|
|
|
+ perc = b
|
|
|
+ if perc:
|
|
|
+ # namespace is limited in DBF but the % value is important
|
|
|
+ if dbfdriver:
|
|
|
+ perccol = "per" + percentile
|
|
|
else:
|
|
|
- coltype = "DOUBLE PRECISION"
|
|
|
- addcols.append(currcolumn + ' ' + coltype)
|
|
|
-
|
|
|
- if addcols:
|
|
|
- grass.verbose(_("Adding columns '%s'") % addcols)
|
|
|
- try:
|
|
|
- grass.run_command('v.db.addcolumn', map=vector, columns=addcols,
|
|
|
- layer=layer)
|
|
|
- except CalledModuleError:
|
|
|
- grass.fatal(_("Adding columns failed. Exiting."))
|
|
|
-
|
|
|
- # calculate statistics:
|
|
|
- grass.message(_("Processing input data (%d categories)...") % number)
|
|
|
-
|
|
|
- # get rid of any earlier attempts
|
|
|
- grass.try_remove(sqltmp)
|
|
|
-
|
|
|
- f = open(sqltmp, 'w')
|
|
|
-
|
|
|
- # do the stats
|
|
|
- p = grass.pipe_command('r.univar', flags='t' + extstat, map=raster,
|
|
|
- zones=rastertmp, percentile=percentile, sep=';')
|
|
|
-
|
|
|
- first_line = 1
|
|
|
-
|
|
|
- f.write("{0}\n".format(grass.db_begin_transaction(fi['driver'])))
|
|
|
- for line in p.stdout:
|
|
|
- if first_line:
|
|
|
- first_line = 0
|
|
|
- continue
|
|
|
-
|
|
|
- vars = decode(line).rstrip('\r\n').split(';')
|
|
|
-
|
|
|
- f.write("UPDATE %s SET" % fi['table'])
|
|
|
- first_var = 1
|
|
|
- for colname in colnames:
|
|
|
- variable = colname.replace("%s_" % colprefix, '', 1)
|
|
|
+ perccol = "percentile_" + percentile
|
|
|
+ percindex = basecols.index(perc)
|
|
|
+ basecols[percindex] = perccol
|
|
|
+
|
|
|
+ # dictionary with name of methods and position in "r.univar -gt" output
|
|
|
+ variables = {'number': 2, 'null_cells': 2, 'minimum': 4, 'maximum': 5, 'range': 6,
|
|
|
+ 'average': 7, 'stddev': 9, 'variance': 10, 'coeff_var': 11,
|
|
|
+ 'sum': 12, 'first_quartile': 14, 'median': 15,
|
|
|
+ 'third_quartile': 16, perccol: 17}
|
|
|
+ # this list is used to set the 'e' flag for r.univar
|
|
|
+ extracols = ['first_quartile', 'median', 'third_quartile', perccol]
|
|
|
+ addcols = []
|
|
|
+ colnames = []
|
|
|
+ extstat = ""
|
|
|
+ for i in basecols:
|
|
|
+ # this check the complete name of out input that should be truncated
|
|
|
+ for k in variables.keys():
|
|
|
+ if i in k:
|
|
|
+ i = k
|
|
|
+ break
|
|
|
+ if i in extracols:
|
|
|
+ extstat = 'e'
|
|
|
+ # check if column already present
|
|
|
+ currcolumn = ("%s_%s" % (colprefix, i))
|
|
|
if dbfdriver:
|
|
|
- variable = variables_dbf[variable]
|
|
|
- i = variables[variable]
|
|
|
- value = vars[i]
|
|
|
- # convert nan, +nan, -nan, inf, +inf, -inf, Infinity, +Infinity,
|
|
|
- # -Infinity to NULL
|
|
|
- if value.lower().endswith('nan') or 'inf' in value.lower():
|
|
|
- value = 'NULL'
|
|
|
- if not first_var:
|
|
|
- f.write(" , ")
|
|
|
+ currcolumn = currcolumn[:10]
|
|
|
+ variables_dbf[currcolumn.replace("%s_" % colprefix, '')] = i
|
|
|
+
|
|
|
+ colnames.append(currcolumn)
|
|
|
+ if currcolumn in grass.vector_columns(vector, layer).keys():
|
|
|
+ if not flags['c']:
|
|
|
+ grass.fatal((_("Cannot create column <%s> (already present). ") % currcolumn) +
|
|
|
+ _("Use -c flag to update values in this column."))
|
|
|
else:
|
|
|
- first_var = 0
|
|
|
- f.write(" %s=%s" % (colname, value))
|
|
|
-
|
|
|
- f.write(" WHERE %s=%s;\n" % (fi['key'], vars[0]))
|
|
|
- f.write("{0}\n".format(grass.db_commit_transaction(fi['driver'])))
|
|
|
- p.wait()
|
|
|
- f.close()
|
|
|
-
|
|
|
- grass.message(_("Updating the database ..."))
|
|
|
- exitcode = 0
|
|
|
- try:
|
|
|
- grass.run_command('db.execute', input=sqltmp,
|
|
|
- database=fi['database'], driver=fi['driver'])
|
|
|
- grass.verbose((_("Statistics calculated from raster map <{raster}>"
|
|
|
- " and uploaded to attribute table"
|
|
|
- " of vector map <{vector}>."
|
|
|
- ).format(raster=raster, vector=vector)))
|
|
|
- except CalledModuleError:
|
|
|
- grass.warning(
|
|
|
- _("Failed to upload statistics to attribute table of vector map <%s>.") %
|
|
|
- vector)
|
|
|
- exitcode = 1
|
|
|
+ if i == "n":
|
|
|
+ coltype = "INTEGER"
|
|
|
+ else:
|
|
|
+ coltype = "DOUBLE PRECISION"
|
|
|
+ addcols.append(currcolumn + ' ' + coltype)
|
|
|
+
|
|
|
+ if addcols:
|
|
|
+ grass.verbose(_("Adding columns '%s'") % addcols)
|
|
|
+ try:
|
|
|
+ grass.run_command('v.db.addcolumn', map=vector, columns=addcols,
|
|
|
+ layer=layer)
|
|
|
+ except CalledModuleError:
|
|
|
+ grass.fatal(_("Adding columns failed. Exiting."))
|
|
|
+
|
|
|
+ # calculate statistics:
|
|
|
+ grass.message(_("Processing input data (%d categories)...") % number)
|
|
|
+
|
|
|
+ # get rid of any earlier attempts
|
|
|
+ grass.try_remove(sqltmp)
|
|
|
+
|
|
|
+ f = open(sqltmp, 'w')
|
|
|
+
|
|
|
+ # do the stats
|
|
|
+ p = grass.pipe_command('r.univar', flags='t' + extstat, map=raster,
|
|
|
+ zones=rastertmp, percentile=percentile, sep=';')
|
|
|
+
|
|
|
+ first_line = 1
|
|
|
+
|
|
|
+ f.write("{0}\n".format(grass.db_begin_transaction(fi['driver'])))
|
|
|
+ for line in p.stdout:
|
|
|
+ if first_line:
|
|
|
+ first_line = 0
|
|
|
+ continue
|
|
|
+
|
|
|
+ vars = decode(line).rstrip('\r\n').split(';')
|
|
|
+
|
|
|
+ f.write("UPDATE %s SET" % fi['table'])
|
|
|
+ first_var = 1
|
|
|
+ for colname in colnames:
|
|
|
+ variable = colname.replace("%s_" % colprefix, '', 1)
|
|
|
+ if dbfdriver:
|
|
|
+ variable = variables_dbf[variable]
|
|
|
+ i = variables[variable]
|
|
|
+ value = vars[i]
|
|
|
+ # convert nan, +nan, -nan, inf, +inf, -inf, Infinity, +Infinity,
|
|
|
+ # -Infinity to NULL
|
|
|
+ if value.lower().endswith('nan') or 'inf' in value.lower():
|
|
|
+ value = 'NULL'
|
|
|
+ if not first_var:
|
|
|
+ f.write(" , ")
|
|
|
+ else:
|
|
|
+ first_var = 0
|
|
|
+ f.write(" %s=%s" % (colname, value))
|
|
|
+
|
|
|
+ f.write(" WHERE %s=%s;\n" % (fi['key'], vars[0]))
|
|
|
+ f.write("{0}\n".format(grass.db_commit_transaction(fi['driver'])))
|
|
|
+ p.wait()
|
|
|
+ f.close()
|
|
|
+
|
|
|
+ grass.message(_("Updating the database ..."))
|
|
|
+ exitcode = 0
|
|
|
+ try:
|
|
|
+ grass.run_command('db.execute', input=sqltmp,
|
|
|
+ database=fi['database'], driver=fi['driver'])
|
|
|
+ grass.verbose((_("Statistics calculated from raster map <{raster}>"
|
|
|
+ " and uploaded to attribute table"
|
|
|
+ " of vector map <{vector}>."
|
|
|
+ ).format(raster=raster, vector=vector)))
|
|
|
+ except CalledModuleError:
|
|
|
+ grass.warning(
|
|
|
+ _("Failed to upload statistics to attribute table of vector map <%s>.") %
|
|
|
+ vector)
|
|
|
+ exitcode = 1
|
|
|
|
|
|
- sys.exit(exitcode)
|
|
|
+ sys.exit(exitcode)
|
|
|
|
|
|
if __name__ == "__main__":
|
|
|
options, flags = grass.parser()
|
Stefan Blumentrath