Univariate statistics for space time vector datasets

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

lib/python/temporal/space_time_datasets.py

@@ -583,9 +583,9 @@ class vector_dataset(abstract_map_dataset):
         if libvector.Vect_open_old_head2(byref(Map), name, mapset, "1") < 2:
 	    # force level 1, open fully
 	    # NOTE: number of points, lines, boundaries, centroids, faces, kernels is still available
-	    libvector.Vect_close(byref(Map))
 	    libvector.Vect_set_open_level(1) # no topology
 	    with_topo = False
+	    core.message(_("Open map without topology support"))
 	    if libvector.Vect_open_old2(byref(Map), name, mapset, "1") < 1:
 		core.fatal(_("Unable to open vector map <%s>"%(libvector.Vect_get_full_name(byref(Map)))))
 
@@ -603,26 +603,33 @@ class vector_dataset(abstract_map_dataset):
 	kvp["is_3d"] = bool(libvector.Vect_is_3d(byref(Map)))
 	
 	# Read number of features
-	kvp["points"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_POINT)
-	kvp["lines"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_LINE)
-	kvp["boundaries"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_BOUNDARY)
-	kvp["centroids"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_CENTROID)
-	kvp["faces"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_FACE)
-	kvp["kernels"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_KERNEL)
-	
-	# Summarize the primitives
-	kvp["primitives"] = kvp["points"] + kvp["lines"] + kvp["boundaries"] + kvp["centroids"]
-	if kvp["is_3d"]:
-	    kvp["primitives"] += kvp["faces"] + kvp["kernels"]
-	
-	# Read topology information
 	if with_topo:
+	    kvp["points"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_POINT)
+	    kvp["lines"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_LINE)
+	    kvp["boundaries"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_BOUNDARY)
+	    kvp["centroids"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_CENTROID)
+	    kvp["faces"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_FACE)
+	    kvp["kernels"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_KERNEL)
+	
+	    # Summarize the primitives
+	    kvp["primitives"] = kvp["points"] + kvp["lines"] + kvp["boundaries"] + kvp["centroids"]
+	    if kvp["is_3d"]:
+	        kvp["primitives"] += kvp["faces"] + kvp["kernels"]
+
+	    # Read topology information
 	    kvp["nodes"] = libvector.Vect_get_num_nodes(byref(Map))
 	    kvp["areas"] = libvector.Vect_get_num_areas(byref(Map))
 	    kvp["islands"] = libvector.Vect_get_num_islands(byref(Map))
 	    kvp["holes"] = libvector.Vect_get_num_holes(byref(Map))
 	    kvp["volumes"] = libvector.Vect_get_num_primitives(byref(Map), libvector.GV_VOLUME)
 	else:
+	    kvp["points"] = None
+	    kvp["lines"] = None
+	    kvp["boundaries"] = None
+	    kvp["centroids"] = None
+	    kvp["faces"] = None
+	    kvp["kernels"] = None
+	    kvp["primitives"] = None
 	    kvp["nodes"] = None
 	    kvp["areas"] = None
 	    kvp["islands"] = None

lib/python/temporal/univar_statistics.py

@@ -23,11 +23,13 @@ for details.
 """
 
 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 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 
@@ -57,8 +59,8 @@ def print_gridded_dataset_univar_statistics(type, input, where, extended, header
     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))
+        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,
@@ -96,7 +98,99 @@ def print_gridded_dataset_univar_statistics(type, input, where, extended, header
         
     dbif.close()
 
-if __name__ == "__main__":
-    options, flags = core.parser()
-    main()
+###############################################################################
+    
+def print_vector_dataset_univar_statistics(input, twhere, layer, type, column, where, extended, header, fs):
+    """!Print univariate statistics for a space time raster or raster3d dataset
+    
+       @param input The name of the space time dataset
+       @param twhere A temporal database where statement
+       @param layer The layer number used in case no layer is present in the temporal dataset
+       @param type options: point,line,boundary,centroid,area
+       @param column The name of the attribute column
+       @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("stvds", 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,layer", twhere, "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))
+
+    string = ""
+    if header == True:
+        string += "id" + fs + "start" + fs + "end" + fs + "n" + fs + "nmissing" + fs + "nnull" + fs
+        string += "min" + fs + "max" + fs + "range"
+	if type == "point" or type == "centroid":
+            string += fs + "mean" + fs + "mean_abs" + fs + "population_stddev" + fs + "population_variance" + fs
+	    string += "population_coeff_variation" + fs + "sample_stddev" + fs + "sample_variance" + fs
+	    string += "kurtosis" + fs + "skewness"
+            if extended == True:
+                string+= fs + "first_quartile" + fs + "median" + fs + "third_quartile" + fs + "percentile_90" 
+
+	print string
+
+    for row in rows:
+        id = row["id"]
+        start = row["start_time"]
+        end = row["end_time"]
+	mylayer = row["layer"]
+
+        flags="g"
+
+        if extended == True:
+            flags += "e"
+
+	if not mylayer:
+	    mylayer = layer
+
+        stats = core.parse_command("v.univar", map=id, where=where, column=column, layer=mylayer, type=type, flags=flags)
+
+	string = ""
+	if stats:
+            string += str(id) + fs + str(start) + fs + str(end)
+            string += fs + str(stats["n"]) + fs + str(stats["nmissing"]) + fs + str(stats["nnull"])
+	    if stats.has_key("min"):
+            	string += fs + str(stats["min"]) + fs + str(stats["max"]) + fs + str(stats["range"])
+	    else:
+            	string += fs + fs + fs
+
+	    if type == "point" or type == "centroid":
+		if stats.has_key("mean"):
+            	    string += fs + str(stats["mean"]) + fs + str(stats["mean_abs"]) + fs + str(stats["population_stddev"]) + fs + str(stats["population_variance"])
+            	    string += fs + str(stats["population_coeff_variation"]) + fs + str(stats["sample_stddev"]) + fs + str(stats["sample_variance"])
+            	    string += fs + str(stats["kurtosis"]) + fs + str(stats["skewness"])
+	        else:
+            	    string += fs + fs + fs + fs + fs + fs + fs + fs + fs
+                if extended == True:
+		    if stats.has_key("first_quartile"):
+                        string += fs + str(stats["first_quartile"]) + fs + str(stats["median"]) + fs + str(stats["third_quartile"]) + fs + str(stats["percentile_90"]) 
+		    else:
+                        string += fs + fs + fs + fs 
+        
+	    print string
+
+    dbif.close()