[GRASS-SVN] r67604 - grass-addons/grass7/raster/r.sample.category
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon Jan 18 18:50:20 PST 2016
Author: annakrat
Date: 2016-01-18 18:50:20 -0800 (Mon, 18 Jan 2016)
New Revision: 67604
Modified:
grass-addons/grass7/raster/r.sample.category/r.sample.category.py
Log:
r.sample.category: add labels from sample map; type of column is now based on sampled raster; check minimum number of points in each category (author: Paulo van Breugel)
Modified: grass-addons/grass7/raster/r.sample.category/r.sample.category.py
===================================================================
--- grass-addons/grass7/raster/r.sample.category/r.sample.category.py 2016-01-18 21:39:00 UTC (rev 67603)
+++ grass-addons/grass7/raster/r.sample.category/r.sample.category.py 2016-01-19 02:50:20 UTC (rev 67604)
@@ -44,11 +44,14 @@
#% multiple: yes
#% type: integer
#%end
+#%flag
+#% key: s
+#% description: If number of cells in category < npoints, skip category
+#%end
# TODO: Python tests for more advanced things such as overwrite or attributes
# TODO: only optional sampling of the category raster
# TODO: preserver original raster categories as vector point categories
-# TODO: ensure/check minimum number of points in each category
# TODO: specify number of points and distribute them uniformly
# TODO: specify number of points and distribute them according to histogram
# TODO: ensure/check minimum and maximum number of of points when doing histogram
@@ -66,9 +69,9 @@
def cleanup():
if gscript.find_file(name='MASK', element='cell', mapset=gscript.gisenv()['MAPSET'])['name']:
- gscript.run_command('r.mask', flags='r')
+ gscript.run_command('r.mask', flags='r', quiet=True)
if TMP:
- gscript.run_command('g.remove', flags='f', type=['raster', 'vector'], name=TMP)
+ gscript.run_command('g.remove', flags='f', type=['raster', 'vector'], name=TMP, quiet=True)
def escape_sql_column(name):
@@ -102,6 +105,7 @@
else:
sampled_rasters = []
npoints = [int(num) for num in options['npoints'].split(',')]
+ flag_s = flags['s']
if gscript.find_file(name='MASK', element='cell', mapset=gscript.gisenv()['MAPSET'])['name']:
gscript.fatal(_("MASK is active. Please remove it before proceeding."))
@@ -114,19 +118,17 @@
TMP.append(points_nocats)
# input must be CELL
- rdescribe = gscript.read_command('r.describe', map=input_raster, flags='d1', quiet=True)
- categories = []
- for line in rdescribe.splitlines():
- try:
- categories.append(int(line))
- except ValueError:
- pass
+ rdescribe = gscript.read_command('r.stats', flags='ln', input=input_raster, separator='pipe')
+ catlab = rdescribe.splitlines()
+ categories = map(int, [z.split('|')[0] for z in catlab])
+ catlab = dict([z.split('|') for z in catlab])
if len(npoints) == 1:
npoints = npoints * len(categories)
else:
if len(categories) != len(npoints):
gscript.fatal(_("Number of categories in raster does not match the number of provided sampling points numbers."))
+ # Create sample points per category
vectors = []
for i, cat in enumerate(categories):
# skip generating points if none are required
@@ -136,32 +138,57 @@
# change mask to sample zeroes and then change again to sample ones
# overwrite mask for an easy loop
gscript.run_command('r.mask', raster=input_raster, maskcats=cat, overwrite=True, quiet=True)
+
+ # Check number of cells in category
+ nrc = int(gscript.parse_command('r.univar', map=input_raster, flags='g')['n'])
+ if nrc < npoints[i]:
+ if flag_s:
+ gscript.info(_("Not enough points in category {cat}. Skipping").format(cat=categories[i]))
+ continue
+ gscript.warning(_("Number of raster cells in category {cat} < {np}. Sampling {n} points").format(cat=categories[i], np=npoints[i], n=nrc))
+ npoints[i] = nrc
+
+ # Create the points
vector = temp_name + str(cat)
vectors.append(vector)
gscript.run_command('r.random', input=input_raster, npoints=npoints[i], vector=vector, quiet=True)
TMP.append(vector)
+ gscript.run_command('r.mask', flags='r', quiet=True)
- gscript.run_command('r.mask', flags='r', quiet=True)
-
gscript.run_command('v.patch', input=vectors, output=points, quiet=True)
# remove and add gain cats so that they are unique
gscript.run_command('v.category', input=points, option='del', cat=-1, output=points_nocats, quiet=True)
# overwrite to reuse the map
gscript.run_command('v.category', input=points_nocats, option='add', output=points, overwrite=True, quiet=True)
+ # Sample layers
columns = []
column_names = []
sampled_rasters.insert(0, input_raster)
for raster in sampled_rasters:
column = escape_sql_column(strip_mapset(raster).lower())
column_names.append(column)
- # TODO: column type according to map type
- columns.append("{column} double precision".format(column=column))
+ datatype = gscript.parse_command('r.info', flags='g', map=raster)['datatype']
+ if datatype == 'CELL':
+ datatype = 'integer'
+ else:
+ datatype = 'double precision'
+ columns.append("{column} {datatype}".format(column=column, datatype=datatype))
gscript.run_command('v.db.addtable', map=points, columns=','.join(columns), quiet=True)
for raster, column in zip(sampled_rasters, column_names):
gscript.info(_("Sampling raster map %s...") % raster)
gscript.run_command('v.what.rast', map=points, type='point', raster=raster, column=column, quiet=True)
+ # Add category labels
+ if not list(set(catlab.values()))[0] and len(set(catlab.values())) == 1:
+ gscript.verbose(_("There are no category labels in the raster to add"))
+ else:
+ gscript.run_command("v.db.addcolumn", map=points, columns="label varchar(250)")
+ table_name = escape_sql_column(strip_mapset(points).lower())
+ for i in categories:
+ sqlstat = "UPDATE " + table_name + " SET label='" + catlab[str(i)] + "' WHERE " + column_names[0] + " == " + str(i)
+ gscript.run_command("db.execute", sql=sqlstat)
+
if __name__ == '__main__':
main()
More information about the grass-commit
mailing list