Added univariate statistics module for space time raster3d datasets

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@51409 15284696-431f-4ddb-bdfa-cd5b030d7da7

lib/python/temporal/Makefile

@@ -8,7 +8,7 @@ PYDIR = $(ETC)/python
 GDIR = $(PYDIR)/grass
 DSTDIR = $(GDIR)/temporal
 
-MODULES = base core abstract_dataset abstract_map_dataset abstract_space_time_dataset space_time_datasets space_time_datasets_tools metadata spatial_extent temporal_extent datetime_math temporal_granularity temporal_relationships unit_tests aggregation extract mapcalc
+MODULES = base core abstract_dataset abstract_map_dataset abstract_space_time_dataset space_time_datasets space_time_datasets_tools metadata spatial_extent temporal_extent datetime_math temporal_granularity temporal_relationships unit_tests aggregation extract mapcalc univar_statistics
 
 PYFILES := $(patsubst %,$(DSTDIR)/%.py,$(MODULES) __init__)
 PYCFILES := $(patsubst %,$(DSTDIR)/%.pyc,$(MODULES) __init__)

lib/python/temporal/__init__.py

@@ -15,3 +15,4 @@ from unit_tests import *
 from aggregation import *
 from extract import *
 from mapcalc import *
+from univar_statistics import *

lib/python/temporal/univar_statistics.py

@@ -0,0 +1,102 @@
+"""!@package grass.temporal
+
+@brief GRASS Python scripting module (temporal GIS functions)
+
+Temporal GIS related functions to be used in Python scripts.
+
+Usage:
+
+@code
+import grass.temporal as tgis
+
+tgis.print_gridded_dataset_univar_statistics(type, input, where, extended, header, fs)
+
+...
+@endcode
+
+(C) 2008-2011 by the GRASS Development Team
+This program is free software under the GNU General Public
+License (>=v2). Read the file COPYING that comes with GRASS
+for details.
+
+@author Soeren Gebbert
+"""
+
+from space_time_datasets_tools import *
+    
+def print_gridded_dataset_univar_statistics(type, input, where, extended, header, fs):
+    """!Print univariate statistics for a space time raster or raster3d dataset
+    
+       @param type Must be "strds" or "str3ds"
+       @param input The name of the space time dataset
+       @param where A temporal database where statement
+       @param extended If True compute extended statistics 
+       @param header   If True print column names as header 
+       @param fs Field separator 
+    """
+    
+    # We need a database interface
+    dbif = sql_database_interface_connection()
+    dbif.connect()
+   
+    mapset =  core.gisenv()["MAPSET"]
+
+    if input.find("@") >= 0:
+        id = input
+    else:
+        id = input + "@" + mapset
+
+    sp = dataset_factory(type, id)
+    
+    if sp.is_in_db(dbif) == False:
+        dbif.close()
+        core.fatal(_("Space time %s dataset <%s> not found") % (sp.get_new_map_instance(None).get_type(), id))
+
+    sp.select(dbif)
+
+    rows = sp.get_registered_maps("id,start_time,end_time", where, "start_time", dbif)
+
+    if not rows:
+            dbif.close()
+            core.fatal(_("Space time %s dataset <%s> is empty") % (sp.get_new_map_instance(None).get_type(), out_id))
+
+    if header == True:
+        print "id" + fs + "start" + fs + "end" + fs + "mean" + fs + "min" + fs + "max" + fs,
+        print "mean_of_abs" + fs + "stddev" + fs + "variance" + fs,
+        if extended == True:
+            print "coeff_var" + fs + "sum" + fs + "null_cells" + fs + "cells" + fs,
+            print "first_quartile" + fs + "median" + fs + "third_quartile" + fs + "percentile_90" 
+        else:
+            print "coeff_var" + fs + "sum" + fs + "null_cells" + fs + "cells" 
+
+    for row in rows:
+        id = row["id"]
+        start = row["start_time"]
+        end = row["end_time"]
+
+        flag="g"
+
+        if extended == True:
+            flag += "e"
+
+	if type == "strds":
+	    stats = core.parse_command("r.univar", map=id, flags=flag)
+	elif type == "str3ds":
+	    stats = core.parse_command("r3.univar", map=id, flags=flag)
+
+        print str(id) + fs + str(start) + fs + str(end),
+        print fs + str(stats["mean"]) + fs + str(stats["min"]) + fs + str(stats["max"]) + fs + str(stats["mean_of_abs"]),
+        print fs + str(stats["stddev"]) + fs + str(stats["variance"]) + fs + str(stats["coeff_var"]) + fs + str(stats["sum"]),
+
+        if extended == True:
+            print fs + str(stats["null_cells"]) + fs + str(stats["cells"]) + fs,
+            print str(stats["first_quartile"]) + fs + str(stats["median"]) + fs + str(stats["third_quartile"]) + fs + str(stats["percentile_90"]) 
+        else:
+            print fs + str(stats["null_cells"]) + fs + str(stats["cells"])
+        
+    dbif.close()
+
+if __name__ == "__main__":
+    options, flags = core.parser()
+    main()
+

temporal/Makefile

@@ -24,6 +24,7 @@ SUBDIRS = \
 	t.rast3d.list \
 	t.rast3d.extract \
 	t.rast3d.mapcalc \
+	t.rast3d.univar \
 	t.vect.list \
 	t.vect.extract \
 	t.vect.what.strds \

temporal/t.rast.univar/t.rast.univar.py

@@ -16,9 +16,9 @@
 
 #%module
 #% description: Calculates univariate statistics from the non-null cells for each registered raster map of a space time raster dataset
-#% keywords: spacetime raster dataset
-#% keywords: raster
+#% keywords: temporal
 #% keywords: statistics
+#% keywords: raster
 #%end
 
 #%option G_OPT_STRDS_INPUT
@@ -61,63 +61,8 @@ def main():
 
     # Make sure the temporal database exists
     tgis.create_temporal_database()
-    # We need a database interface
-    dbif = tgis.sql_database_interface_connection()
-    dbif.connect()
-   
-    mapset =  grass.gisenv()["MAPSET"]
-
-    if input.find("@") >= 0:
-        id = input
-    else:
-        id = input + "@" + mapset
-
-    sp = tgis.space_time_raster_dataset(id)
     
-    if sp.is_in_db(dbif) == False:
-        dbif.close()
-        grass.fatal(_("Space time %s dataset <%s> not found") % (sp.get_new_map_instance(None).get_type(), id))
-
-    sp.select(dbif)
-
-    rows = sp.get_registered_maps("id,start_time,end_time", where, "start_time", dbif)
-
-    if not rows:
-            dbif.close()
-            grass.fatal(_("Space time raster dataset <%s> is empty") % out_id)
-
-    if header == True:
-        print "id" + fs + "start" + fs + "end" + fs + "mean" + fs + "min" + fs + "max" + fs,
-        print "mean_of_abs" + fs + "stddev" + fs + "variance" + fs,
-        if extended == True:
-            print "coeff_var" + fs + "sum" + fs + "null_cells" + fs + "cells" + fs,
-            print "first_quartile" + fs + "median" + fs + "third_quartile" + fs + "percentile_90" 
-        else:
-            print "coeff_var" + fs + "sum" + fs + "null_cells" + fs + "cells" 
-
-    for row in rows:
-        id = row["id"]
-        start = row["start_time"]
-        end = row["end_time"]
-
-        flag="g"
-
-        if extended == True:
-            flag += "e"
-
-        stats = grass.parse_command("r.univar", map=id, flags=flag)
-
-        print str(id) + fs + str(start) + fs + str(end),
-        print fs + str(stats["mean"]) + fs + str(stats["min"]) + fs + str(stats["max"]) + fs + str(stats["mean_of_abs"]),
-        print fs + str(stats["stddev"]) + fs + str(stats["variance"]) + fs + str(stats["coeff_var"]) + fs + str(stats["sum"]),
-
-        if extended == True:
-            print fs + str(stats["null_cells"]) + fs + str(stats["cells"]) + fs,
-            print str(stats["first_quartile"]) + fs + str(stats["median"]) + fs + str(stats["third_quartile"]) + fs + str(stats["percentile_90"]) 
-        else:
-            print fs + str(stats["null_cells"]) + fs + str(stats["cells"])
-        
-    dbif.close()
+    tgis.print_gridded_dataset_univar_statistics("strds", input, where, extended, header, fs)
 
 if __name__ == "__main__":
     options, flags = grass.parser()

temporal/t.rast.univar/test.t.rast.univar.sh

@@ -20,3 +20,4 @@ t.rast.univar -he input=precip_abs1
 
 t.unregister type=rast maps=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6
 t.remove type=strds input=precip_abs1
+g.remove rast=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6

temporal/t.rast3d.univar/Makefile

@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../../
+
+PGM = t.rast3d.univar
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script $(TEST_DST)

temporal/t.rast3d.univar/t.rast3d.univar.py

@@ -0,0 +1,70 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+############################################################################
+#
+# MODULE:	t.rast3d.univar
+# AUTHOR(S):	Soeren Gebbert
+#
+# PURPOSE:	Calculates univariate statistics from the non-null cells for each registered raster3d map of a space time raster3d dataset
+# COPYRIGHT:	(C) 2011 by the GRASS Development Team
+#
+#		This program is free software under the GNU General Public
+#		License (version 2). Read the file COPYING that comes with GRASS
+#		for details.
+#
+#############################################################################
+
+#%module
+#% description: Calculates univariate statistics from the non-null cells for each registered raster3d map of a space time raster3d dataset
+#% keywords: temporal
+#% keywords: statistics
+#% keywords: raster
+#%end
+
+#%option G_OPT_STR3DS_INPUT
+#%end
+
+#%option G_OPT_T_WHERE
+#%end
+
+#%option
+#% key: fs
+#% type: string
+#% description: The field separator character between the output columns
+#% required: no
+#% answer: |
+#%end
+
+#%flag
+#% key: e
+#% description: Calculate extended statistics
+#%end
+
+#%flag
+#% key: h
+#% description: Print column names 
+#%end
+
+import grass.script as grass
+import grass.temporal as tgis
+
+############################################################################
+
+def main():
+
+    # Get the options
+    input = options["input"]
+    where = options["where"]
+    extended = flags["e"]
+    header = flags["h"]
+    fs = options["fs"]
+
+    # Make sure the temporal database exists
+    tgis.create_temporal_database()
+    
+    tgis.print_gridded_dataset_univar_statistics("str3ds", input, where, extended, header, fs)
+ 
+if __name__ == "__main__":
+    options, flags = grass.parser()
+    main()
+

temporal/t.rast3d.univar/test.t.rast3d.univar.sh

@@ -0,0 +1,24 @@
+# Univariate statitsics for space time raster3d datasets
+
+# We need to set a specific region in the
+# @preprocess step of this test. 
+# The region setting should work for UTM and LL test locations
+g.region s=0 n=80 w=0 e=120 b=0 t=50 res=10 res3=10 -p3
+# Data generation
+r3.mapcalc --o expr="prec_1 = rand(0, 550)"
+r3.mapcalc --o expr="prec_2 = rand(0, 450)"
+r3.mapcalc --o expr="prec_3 = rand(0, 320)"
+r3.mapcalc --o expr="prec_4 = rand(0, 510)"
+r3.mapcalc --o expr="prec_5 = rand(0, 300)"
+r3.mapcalc --o expr="prec_6 = rand(0, 650)"
+
+t.create --o type=str3ds temporaltype=absolute output=precip_abs1 title="A test" descr="A test"
+t.register type=rast3d --v -i input=precip_abs1 maps=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6 start="2001-01-15 12:05:45" increment="14 days"
+
+# The first @test
+t.rast3d.univar -he input=precip_abs1 
+
+t.unregister type=rast3d maps=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6
+t.remove type=str3ds input=precip_abs1
+g.remove rast3d=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6
+