[GRASS-SVN] r70283 - sandbox/martinl

svn_grass at osgeo.org svn_grass at osgeo.org
Fri Jan 6 10:49:53 PST 2017


Author: martinl
Date: 2017-01-06 10:49:53 -0800 (Fri, 06 Jan 2017)
New Revision: 70283

Modified:
   sandbox/martinl/tiles-rast-stats.py
Log:
tiles-rast-stats: implement v.external where approach

Modified: sandbox/martinl/tiles-rast-stats.py
===================================================================
--- sandbox/martinl/tiles-rast-stats.py	2017-01-06 17:46:09 UTC (rev 70282)
+++ sandbox/martinl/tiles-rast-stats.py	2017-01-06 18:49:53 UTC (rev 70283)
@@ -26,6 +26,11 @@
 #% keyword: querying
 #% overwrite: yes
 #%end
+#%flag
+#% key: l
+#% label: Link tiles instead of import
+#% description: Requires datasource with multiple access capabilities (typically PostGIS). Also output_layer option is ignored since result is written directly to input layer
+#%end
 #%option
 #% key: dsn
 #% description: Data source name
@@ -40,9 +45,8 @@
 #% description: Name of key column used for input
 #%end
 #%option G_OPT_V_LAYER
-#% key: olayer
+#% key: output_layer
 #% description: Name for output layer
-#% required: yes
 #%end
 #%option G_OPT_R_INPUT
 #% key: raster
@@ -56,8 +60,8 @@
 #%option
 #% key: column_prefix
 #% type: string
-#% description: Column prefix for new attribute columns
-#% required : yes
+#% description: Column prefix for new attribute columns (default: map name)
+#% required : no
 #%end
 #%option
 #% key: method
@@ -82,6 +86,9 @@
 #% answer: 1
 #% type: integer
 #%end
+#%rules
+#% required: output_layer,-l
+#%end
 
 import os
 import sys
@@ -101,7 +108,7 @@
         grass.message("Cleaning...")
         Module('g.remove', type='vector', pattern='{}*'.format(basename), flags='f', quiet=True)
     
-def import_tiles(dsn, layer_name, keycolumn, mask):
+def filter_tiles(dsn, layer_name, keycolumn, mask):
     ds = ogr.Open(dsn, True) # write
     if ds is None:
         grass.fatal("Unable to open '{}'".format(dsn))
@@ -121,7 +128,6 @@
 
     layer.SetSpatialFilter(ogr.CreateGeometryFromWkt(wkt))
 
-    basename = 'tiles_rast_stats_{}'.format(os.getpid())
     tiles = []
     try:
         for feature in layer:
@@ -130,6 +136,11 @@
         grass.fatal("Column <{}> not found".format(keycolumn))
     grass.message("{} tiles filtered".format(len(tiles)))
 
+    return ds, layer, tiles
+
+def import_tiles(dsn, layer_name, keycolumn, mask):
+    ds, layer, tiles = filter_tiles(dsn, layer_name, keycolumn, mask)
+    
     # database for each map is required since we call v.rast.stats in paralell
     Module('db.connect', database='$GISDBASE/$LOCATION_NAME/$MAPSET/$MAP/sqlite.db')
 
@@ -141,16 +152,34 @@
         out_name = '{}_{}'.format(basename, tile)
         maps.append(out_name)
         new_import = copy.deepcopy(import_module)
-        queue.put(new_import(input=dsn, where="{}={}".format(keycolumn, tile),
+        queue.put(new_import(where="{}={}".format(keycolumn, tile),
                              output=out_name))
     queue.wait()
     grass.message("... done in {:.1f} sec".format(time.time() - start))
     
     return maps, ds, layer
 
+def link_tiles(dsn, layer_name, keycolumn, mask):
+    ds, layer, tiles = filter_tiles(dsn, layer_name, keycolumn, mask)
+    
+    grass.message("Linking tiles...")
+    start = time.time()
+    import_module = Module('v.external', input=dsn, layer=layer_name, quiet=True, overwrite=True, run_=False)
+    maps = []
+    for tile in tiles:
+        out_name = '{}_{}'.format(basename, tile)
+        maps.append(out_name)
+        new_import = copy.deepcopy(import_module)
+        queue.put(new_import(where="{}={}".format(keycolumn, tile),
+                             output=out_name))
+    queue.wait()
+    grass.message("... done in {:.1f} sec".format(time.time() - start))
+    
+    return maps, ds, layer
+
 def perform_statistics(maps, raster, column_prefix, method, percentile):
-    region_module =  Module("g.region", flags='a', raster=raster, run_=False)
-    stats_module = Module('v.rast.stats', raster=raster, column_prefix=column_prefix, method=method.split(','),
+    region_module =  Module("g.region", align=raster, run_=False)
+    stats_module = Module('v.rast.stats', flags='c', raster=raster, column_prefix=column_prefix, method=method.split(','),
                           percentile=percentile, quiet=True, run_=False)
     grass.message("Performing statistics...")
     start = time.time()
@@ -218,25 +247,38 @@
         grass.fatal("Raster map <{}> not found".format(options['raster']))
 
     start = time.time()
-    
-    # import tiles
-    maps, ds, ilayer = import_tiles(options['dsn'], options['layer'], options['keycolumn'], options['mask'])
 
-    # Remove output if it already exists
-    if ds.GetLayerByName(options['olayer']):
-        if os.getenv('GRASS_OVERWRITE', '0') == '1':
-            ds.DeleteLayer(options['olayer'])
-        else:
-            grass.fatal("option <output>: <{}> exists. To overwrite, "
-                        "use the --overwrite flag".format(options['olayer']))
+    if flags['l']:
+        if options['output_layer']:
+            grass.warning("Option <{}> will be ignored".format('output_layer'))
 
+        # link tiles
+        maps, ds, ilayer = link_tiles(options['dsn'], options['layer'], options['keycolumn'], options['mask'])
+        
+    else:
+        # import tiles
+        maps, ds, ilayer = import_tiles(options['dsn'], options['layer'], options['keycolumn'], options['mask'])
+
+        # Remove output if it already exists
+        if ds.GetLayerByName(options['output_layer']):
+            if os.getenv('GRASS_OVERWRITE', '0') == '1':
+                ds.DeleteLayer(options['output_layer'])
+            else:
+                grass.fatal("option <output>: <{}> exists. To overwrite, "
+                            "use the --overwrite flag".format(options['output_layer']))
+
+    column_prefix = options['raster'] if not options['column_prefix'] else options['column_prefix']
+    
     # perform statistics
-    perform_statistics(maps, options['raster'], options['column_prefix'], options['method'],
+    perform_statistics(maps, options['raster'],
+                       column_prefix,
+                       options['method'],
                        options['percentile'])
 
     # write output
-    write_output(maps, ds, ilayer, options['column_prefix'], options['method'], options['percentile'],
-                 options['olayer'])
+    if not flags['l']:
+        write_output(maps, ds, ilayer, column_prefix, options['method'], options['percentile'],
+                     options['output_layer'])
 
     grass.message("Done in {:.1f} sec".format(time.time() - start))
 
@@ -246,7 +288,7 @@
 if __name__ == "__main__":
     options, flags = grass.parser()
 
-    basename = None
+    basename = 'tiles_rast_stats_{}'.format(os.getpid())
 
     # queue for parallel jobs
     queue = ParallelModuleQueue(int(options['nprocs']))



More information about the grass-commit mailing list