New VTK time series export modules. New option for tr.list: sample space

time datasets with granularity. Bug fixing and tests.


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

temporal/Makefile

@@ -17,6 +17,7 @@ SUBDIRS = \
 	tr3.list \
 	tr.series \
 	tr.export \
+	tr.out.vtk \
 	tr.import \
 	tr.extract \
 	tr3.register \

temporal/t.time.rel/test.t.time.rel.sh

@@ -28,9 +28,13 @@ tr.register --v input=precip_rel3 maps=prec_1,prec_2,prec_3
 # Check if the space time inputs are updated correctly
 t.time.rel --v input=prec_1,prec_2,prec_3 start=0 increment=1000
 t.info type=strds input=precip_rel1
+t.topology type=strds input=precip_rel1
+t.topology -m type=strds input=precip_rel1
 
 t.time.rel --v -i input=prec_1,prec_2,prec_3 start=0 increment=1000 -i
 t.info type=strds input=precip_rel1
+t.topology type=strds input=precip_rel1
+t.topology -m type=strds input=precip_rel1
 
 t.time.rel --v input=prec_4,prec_5 start=5000 end=6000
 t.info type=rast input=prec_4

temporal/t.topology/t.topology.py

@@ -41,28 +41,8 @@
 #%end
 
 #%flag
-#% key: c
-#% description: Check temporal topology
-#%end
-
-#%flag
-#% key: g
-#% description: Compute granularity
-#%end
-
-#%flag
-#% key: t
-#% description: Print temporal relation matrix and exit
-#%end
-
-#%flag
-#% key: r
-#% description: Print temporal relation count
-#%end
-
-#%flag
 #% key: m
-#% description: Print temporal map type count
+#% description: Print temporal relation matrix and exit
 #%end
 
 import grass.script as grass
@@ -76,11 +56,7 @@ def main():
     name = options["input"]
     type = options["type"]
     where = options["where"]
-    tmatrix = flags['t']
-    relation = flags['r']
-    map_types = flags['m']
-    check_topo = flags['c']
-    comp_gran = flags['g']
+    tmatrix = flags['m']
 
     # Make sure the temporal database exists
     tgis.create_temporal_database()
@@ -103,36 +79,79 @@ def main():
     # Get ordered map list
     maps = sp.get_registered_maps_as_objects(where=where, order="start_time", dbif=None)
 
-    if check_topo:
-        check = sp.check_temporal_topology(maps)
-        if check:
-            print "Temporal topology is valid"
-        else:
-            print "Temporal topology is invalid"
-
     if tmatrix:
         matrix = sp.print_temporal_relation_matrix(maps)
+        return
 
-    if relation:
-        dict = sp.count_temporal_relations(maps)
-
-        for key in dict.keys():
-            print key, dict[key]
+    sp.base.print_info()
 
-    if map_types:
-        dict = sp.count_temporal_types(maps)
-        
-        for key in dict.keys():
-            print key, dict[key]
-
-        if dict["invalid"] == 0:
-            print "Gaps", sp.count_gaps(maps)
-        else:
-            print "Gaps", None
+    #      0123456789012345678901234567890
+    print " +-------------------- Temporal topology -------------------------------------+"
+    if where:
+        print " | Is subset of dataset: ...... True"
+    else:
+        print " | Is subset of dataset: ...... False"
 
-    if comp_gran:
+    check = sp.check_temporal_topology(maps)
+    if check:
+        #      0123456789012345678901234567890
+        print " | Temporal topology is: ...... valid"                
+    else:
+        #      0123456789012345678901234567890
+        print " | Temporal topology is: ...... invalid"                
+
+    dict = sp.count_temporal_types(maps)
+ 
+    for key in dict.keys():
+        if key == "interval":
+            #      0123456789012345678901234567890
+            print " | Number of intervals: ....... %s" % (dict[key])
+        if key == "point":
+            print " | Number of points: .......... %s" % (dict[key])
+        if key == "invalid":
+            print " | Invalid time stamps: ....... %s" % (dict[key])
+
+    #      0123456789012345678901234567890
+    print " | Number of gaps: ............ %i" % sp.count_gaps(maps)                
+
+    if sp.is_time_absolute():
         gran = tgis.compute_absolute_time_granularity(maps)
-        print "Granularity", gran
+    else:
+        gran = tgis.compute_relative_time_granularity(maps)
+    print " | Granularity: ............... %s" % str(gran)
+
+    print " +-------------------- Topological relations ---------------------------------+"
+    dict = sp.count_temporal_relations(maps)
+
+    for key in dict.keys():
+        if key == "equivalent":
+            #      0123456789012345678901234567890
+            print " | Equivalent: ................ %s" % (dict[key])
+        if key == "during":
+            print " | During: .................... %s" % (dict[key])
+        if key == "contains":
+            print " | Contains: .................. %s" % (dict[key])
+        if key == "overlaps":
+            print " | Overlaps: .................. %s" % (dict[key])
+        if key == "overlapped":
+            print " | Overlapped: ................ %s" % (dict[key])
+        if key == "after":
+            print " | After: ..................... %s" % (dict[key])
+        if key == "before":
+            print " | Before: .................... %s" % (dict[key])
+        if key == "starts":
+            print " | Starts: .................... %s" % (dict[key])
+        if key == "finishes":
+            print " | Finishes: .................. %s" % (dict[key])
+        if key == "started":
+            print " | Started: ................... %s" % (dict[key])
+        if key == "finished":
+            print " | Finished: .................. %s" % (dict[key])
+        if key == "follows":
+            print " | Follows: ................... %s" % (dict[key])
+        if key == "precedes":
+            print " | Precedes: .................. %s" % (dict[key])
+    print " +----------------------------------------------------------------------------+"
 
 if __name__ == "__main__":
     options, flags = grass.parser()

temporal/t.topology/test.t.topology.sh

@@ -74,29 +74,29 @@ EOF
 t.create --o type=strds temporaltype=absolute output=precip_abs title="A test with input files" descr="A test with input files"
 
 tr.register -i input=precip_abs file=$n1 start="2001-01-01" increment="1 months"
-t.info type=strds input=precip_abs
 cat $n1
-t.topology -gctmr input=precip_abs
+t.topology    input=precip_abs
+t.topology -m input=precip_abs
 
 tr.register -i input=precip_abs file=$n2 start=file
-t.info type=strds input=precip_abs
 cat $n2
-t.topology -gctmr input=precip_abs
+t.topology    input=precip_abs
+t.topology -m input=precip_abs
 
 tr.register -i input=precip_abs file=$n3 start=file end=file
-t.info type=strds input=precip_abs
 cat $n3
-t.topology -gctmr input=precip_abs
+t.topology    input=precip_abs
+t.topology -m input=precip_abs
 
 tr.register -i input=precip_abs file=$n4 start=file end=file
-t.info type=strds input=precip_abs
 cat $n4
-t.topology -gctmr input=precip_abs
+t.topology    input=precip_abs
+t.topology -m input=precip_abs
 
 tr.register -i input=precip_abs file=$n5 start=file end=file
-t.info type=strds input=precip_abs
 cat $n5
-t.topology -gctmr input=precip_abs
+t.topology    input=precip_abs
+t.topology -m input=precip_abs
 
 t.remove type=strds input=precip_abs
 t.remove type=rast file=$n1

temporal/t.topology/test.t.topology.reltime.sh

@@ -0,0 +1,102 @@
+# This is a test to register and unregister raster maps in
+# space time raster input.
+# The raster maps will be registered in different space time raster
+# inputs
+
+# We need to set a specific region in the
+# @preprocess step of this test. We generate
+# raster with r.mapcalc and create two space time raster inputs
+# with relative and relolute time
+# 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
+
+r.mapcalc --o expr="prec_1 = rand(0, 550)"
+r.mapcalc --o expr="prec_2 = rand(0, 450)"
+r.mapcalc --o expr="prec_3 = rand(0, 320)"
+r.mapcalc --o expr="prec_4 = rand(0, 510)"
+r.mapcalc --o expr="prec_5 = rand(0, 300)"
+r.mapcalc --o expr="prec_6 = rand(0, 650)"
+
+n1=`g.tempfile pid=1 -d` 
+n2=`g.tempfile pid=2 -d`
+n3=`g.tempfile pid=3 -d`
+n4=`g.tempfile pid=4 -d`
+n5=`g.tempfile pid=5 -d`
+
+cat > $n1 << EOF
+prec_1
+prec_2
+prec_3
+prec_4
+prec_5
+prec_6
+EOF
+
+cat > $n2 << EOF
+prec_1|1
+prec_2|2
+prec_3|3
+prec_4|4
+prec_5|5
+prec_6|6
+EOF
+
+cat > $n3 << EOF
+prec_1|1|4
+prec_2|5|7
+prec_3|8|10
+prec_4|11|12
+prec_5|13|14
+prec_6|15|25
+EOF
+
+cat > $n4 << EOF
+prec_1|2001|2007
+prec_2|2002|2004
+prec_3|2003|2004
+prec_4|2004|2006
+prec_5|2005|2006
+prec_6|2006|2007
+EOF
+
+cat > $n5 << EOF
+prec_1|2001|2003
+prec_2|2002|2004
+prec_3|2003|2006
+prec_4|2004|2004
+prec_5|2005|2005
+prec_6|2006|2007
+EOF
+
+
+# The first @test
+# We create the space time raster inputs and register the raster maps with relolute time interval
+t.create --o type=strds temporaltype=relative output=precip_rel title="A test with input files" descr="A test with input files"
+
+tr.register -i --v input=precip_rel file=$n1 start=0 increment=1
+cat $n1
+t.topology    input=precip_rel
+t.topology -m input=precip_rel
+
+tr.register -i input=precip_rel file=$n2 start=file
+cat $n2
+t.topology    input=precip_rel
+t.topology -m input=precip_rel
+
+tr.register -i input=precip_rel file=$n3 start=file end=file
+cat $n3
+t.topology    input=precip_rel
+t.topology -m input=precip_rel
+
+tr.register -i input=precip_rel file=$n4 start=file end=file
+cat $n4
+t.topology    input=precip_rel
+t.topology -m input=precip_rel
+
+tr.register -i input=precip_rel file=$n5 start=file end=file
+cat $n5
+t.topology    input=precip_rel
+t.topology -m input=precip_rel
+
+t.remove type=strds input=precip_rel
+t.remove type=rast file=$n1

temporal/tr.list/test.tr.list.sh

@@ -64,19 +64,29 @@ prec_5|2001-05-01|2001-05-01
 prec_6|2001-06-01|2001-07-01
 EOF
 
+t.create --o type=strds temporaltype=absolute output=precip_abs0 title="A test with input files" descr="A test with input files"
 t.create --o type=strds temporaltype=absolute output=precip_abs1 title="A test with input files" descr="A test with input files"
 t.create --o type=strds temporaltype=absolute output=precip_abs2 title="A test with input files" descr="A test with input files"
 t.create --o type=strds temporaltype=absolute output=precip_abs3 title="A test with input files" descr="A test with input files"
 t.create --o type=strds temporaltype=absolute output=precip_abs4 title="A test with input files" descr="A test with input files"
 t.create --o type=strds temporaltype=absolute output=precip_abs5 title="A test with input files" descr="A test with input files"
 
-# The first @test
+# The @test
+tr.register -i input=precip_abs0 file=$n1 start="2001-01-01" increment="1 months"
+tr.list    fs=" | " method=comma     input=precip_abs0
+tr.list -h input=precip_abs1
+tr.list -h fs=" | " method=cols      input=precip_abs0
+tr.list -h fs=" | " method=delta     input=precip_abs0
+tr.list -h fs=" | " method=deltagaps input=precip_abs0
+tr.list -h fs=" | " method=gran      input=precip_abs0
+
 tr.register    input=precip_abs1 file=$n1 start="2001-01-01" increment="1 months"
 tr.list    fs=" | " method=comma     input=precip_abs1
 tr.list -h input=precip_abs1
 tr.list -h fs=" | " method=cols      input=precip_abs1
 tr.list -h fs=" | " method=delta     input=precip_abs1
 tr.list -h fs=" | " method=deltagaps input=precip_abs1
+tr.list -h fs=" | " method=gran      input=precip_abs1
 
 tr.register -i input=precip_abs2 file=$n2 start=file
 tr.list    fs=" | " method=comma     input=precip_abs2
@@ -84,16 +94,19 @@ tr.list -h input=precip_abs2
 tr.list -h fs=" | " method=cols      input=precip_abs2
 tr.list -h fs=" | " method=delta     input=precip_abs2
 tr.list -h fs=" | " method=deltagaps input=precip_abs2
+tr.list -h fs=" | " method=gran      input=precip_abs2
 
 tr.register -i input=precip_abs3 file=$n3 start=file end=file
 tr.list    fs=" | " method=comma     input=precip_abs3
 tr.list -h fs=" | " method=delta     input=precip_abs3
 tr.list -h fs=" | " method=deltagaps input=precip_abs3
+tr.list -h fs=" | " method=gran      input=precip_abs3
 
 tr.register -i input=precip_abs4 file=$n4 start=file end=file
 tr.list    fs=" | " method=comma     input=precip_abs4
 tr.list -h fs=" | " method=delta     input=precip_abs4
 tr.list -h fs=" | " method=deltagaps input=precip_abs4
+tr.list -h fs=" | " method=gran      input=precip_abs4
 
 tr.register -i input=precip_abs5 file=$n5 start=file end=file
 tr.list    fs=" | " method=comma     input=precip_abs5
@@ -101,6 +114,7 @@ tr.list -h input=precip_abs5
 tr.list -h fs=" | " method=cols      input=precip_abs5
 tr.list -h fs=" | " method=delta     input=precip_abs5
 tr.list -h fs=" | " method=deltagaps input=precip_abs5
+tr.list -h fs=" | " method=gran      input=precip_abs5
 
 t.remove type=rast input=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6
-t.remove type=strds input=precip_abs1,precip_abs2,precip_abs3,precip_abs4,precip_abs5
+t.remove type=strds input=precip_abs0,precip_abs1,precip_abs2,precip_abs3,precip_abs4,precip_abs5

temporal/tr.list/tr.list.py

@@ -61,10 +61,10 @@
 #%option
 #% key: method
 #% type: string
-#% description: Which method should be used fot data listing
+#% description: Which method should be used for data listing
 #% required: no
 #% multiple: no
-#% options: cols,comma,delta,deltagaps
+#% options: cols,comma,delta,deltagaps,gran
 #% answer: cols
 #%end
 

temporal/tr.out.vtk/Makefile

@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../../
+
+PGM = tr.out.vtk
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script $(TEST_DST)

temporal/tr.out.vtk/test.tr.out.vtk.sh

@@ -0,0 +1,53 @@
+# This is a test to export space time raster datasets as VTK time series data
+
+# We need to set a specific region in the
+# @preprocess step of this test. We generate
+# raster with r.mapcalc and create a space time raster datasets
+# The region setting should work for UTM and LL test locations
+g.region s=0 n=80 w=0 e=120 res=0.5 -p
+
+r.mapcalc --o expr="prec_1 = rand(0, 550)"
+r.mapcalc --o expr="prec_2 = rand(0, 450)"
+r.mapcalc --o expr="prec_3 = rand(0, 320)"
+r.mapcalc --o expr="prec_4 = rand(0, 510)"
+r.mapcalc --o expr="prec_5 = rand(0, 300)"
+r.mapcalc --o expr="prec_6 = rand(0, 650)"
+r.mapcalc --o expr="elevation = sin(row() + col()) * 10"
+
+n1=`g.tempfile pid=1 -d` 
+
+cat > $n1 << EOF
+prec_1|2001-01-01|2001-02-01
+prec_2|2001-02-01|2001-03-01
+prec_3|2001-03-01|2001-04-01
+prec_4|2001-05-01|2001-06-01
+prec_5|2001-07-01|2001-08-01
+prec_6|2001-08-01|2001-09-01
+EOF
+
+t.create --o type=strds temporaltype=absolute output=precip_abs1 title="A test with input files" descr="A test with input files"
+t.create --o type=strds temporaltype=absolute output=precip_abs2 title="A test with input files" descr="A test with input files"
+t.create --o type=strds temporaltype=absolute output=precip_abs3 title="A test with input files" descr="A test with input files"
+
+# The first @test
+mkdir /tmp/test1
+tr.register -i input=precip_abs1 file=$n1 start=file end=file
+tr.out.vtk input=precip_abs1 expdir=/tmp/test1 
+ls -al /tmp/test1 
+
+mkdir /tmp/test2
+tr.register -i input=precip_abs2 file=$n1 
+tr.out.vtk input=precip_abs2 expdir=/tmp/test2 elevation=elevation
+ls -al /tmp/test2 
+
+mkdir /tmp/test3
+tr.register -i input=precip_abs3 file=$n1
+tr.out.vtk -g input=precip_abs3 expdir=/tmp/test3 elevation=elevation
+ls -al /tmp/test3 
+
+rm -rf /tmp/test1
+rm -rf /tmp/test2
+rm -rf /tmp/test3
+
+t.remove type=rast input=prec_1,prec_2,prec_3,prec_4,prec_5,prec_6
+t.remove type=strds input=precip_abs1,precip_abs2,precip_abs3

temporal/tr.out.vtk/tr.out.vtk.py

@@ -0,0 +1,169 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+############################################################################
+#
+# MODULE:	tr.out.vtk
+# AUTHOR(S):	Soeren Gebbert
+#               
+# PURPOSE:	Export space time raster dataset as VTK time series
+# 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: Export space time raster dataset as VTK time series
+#% keywords: dataset
+#% keywords: spacetime
+#% keywords: raster
+#% keywords: export
+#% keywords: VTK
+#%end
+
+#%option
+#% key: input
+#% type: string
+#% description: Name of a space time raster dataset
+#% required: yes
+#% multiple: no
+#%end
+
+#%option
+#% key: expdir
+#% type: string
+#% description: Path to the export directory
+#% required: yes
+#% multiple: no
+#%end
+
+#%option G_OPT_R_INPUT
+#% key: elevation
+#% type: string
+#% description: Elevation raster map
+#% required: no
+#% multiple: no
+#%end
+
+#%option
+#% key: where
+#% type: string
+#% description: A where statement for selected listing e.g: (start_time < '2001-01-01' and end_time > '2001-01-01')
+#% required: no
+#% multiple: no
+#%end
+
+#%option
+#% key: null
+#% type: double
+#% description: Value to represent no data cell
+#% required: no
+#% multiple: no
+#% answer: -10.0
+#%end
+
+#%flag
+#% key: p
+#% description: Create VTK point data instead of VTK cell data (if no elevation map is given)
+#%end
+
+#%flag
+#% key: c
+#% description: Correct the coordinates to fit the VTK-OpenGL precision
+#%end
+
+#%flag
+#% key: g
+#% description: Export files using the space time dataset granularity for equidistant time between maps, where statement will be ignored
+#%end
+
+
+
+import os
+import grass.script as grass
+import grass.temporal as tgis
+
+############################################################################
+
+def main():
+
+    # Get the options
+    input = options["input"]
+    elevation = options["elevation"]
+    expdir = options["expdir"]
+    where = options["where"]
+    null = options["null"]
+    use_pdata = flags["p"]
+    coorcorr = flags["c"]
+    use_granularity = flags["g"]
+
+    # Make sure the temporal database exists
+    tgis.create_temporal_database()
+
+    if not os.path.exists(expdir): 
+        grass.fatal(_("Export directory <%s> not found.") % expdir)
+
+    os.chdir(expdir)
+
+    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() == False:
+        grass.fatal(_("Dataset <%s> not found in temporal database") % (id))
+
+    sp.select()
+
+    if use_granularity:
+        maps = sp.get_registered_maps_as_objects_by_granularity()
+        # Create a NULL map in case of granularity support
+        null_map = "temporary_null_map_%i" % os.getpid()
+        grass.mapcalc("%s = null()" % (null_map))
+    else:
+        maps = sp.get_registered_maps_as_objects(where, "start_time", None)
+
+    # To have scalar values with the same name, we need to copy the raster maps using a single name
+    map_name = "%s_%i" % (sp.base.get_name(), os.getpid())
+
+    count = 0
+    if maps:
+        for map in maps:
+            id = map.get_id()
+            # None ids will be replaced by NULL maps
+            if id == None:
+                id = null_map
+            
+            grass.run_command("g.copy", rast="%s,%s" % (id, map_name), overwrite=True)
+            out_name = "%6.6i_%s.vtk" % (count, sp.base.get_name()) 
+
+            mflags=""
+            if use_pdata:
+                mflags += "p"
+            if coorcorr:
+                mflags += "c"
+
+            # Export the raster map with r.out.vtk
+            if elevation:
+                ret = grass.run_command("r.out.vtk", flags=mflags, null=null, input=map_name, elevation=elevation, output=out_name, overwrite=grass.overwrite())
+            else:
+                ret = grass.run_command("r.out.vtk", flags=mflags, null=null, input=map_name, output=out_name, overwrite=grass.overwrite())
+
+            if ret != 0:
+                grass.fatal(_("Unable to export raster map <%s>" % name))
+
+            count += 1
+
+    if use_granularity:
+        grass.run_command("g.remove", rast=null_map)
+    grass.run_command("g.remove", rast=map_name)
+
+if __name__ == "__main__":
+    options, flags = grass.parser()
+    main()

temporal/tr3.list/tr3.list.py

@@ -61,10 +61,10 @@
 #%option
 #% key: method
 #% type: string
-#% description: Which method should be used fot data listing
+#% description: Which method should be used for data listing
 #% required: no
 #% multiple: no
-#% options: cols,comma,delta,deltagaps
+#% options: cols,comma,delta,deltagaps,gran
 #% answer: cols
 #%end
 

temporal/tv.list/tv.list.py

@@ -61,10 +61,10 @@
 #%option
 #% key: method
 #% type: string
-#% description: Which method should be used fot data listing
+#% description: Which method should be used for data listing
 #% required: no
 #% multiple: no
-#% options: cols,comma,delta,deltagaps
+#% options: cols,comma,delta,deltagaps,gran
 #% answer: cols
 #%end