← Back to branch summary

~ubuntu-branches/ubuntu/vivid/grass/vivid-proposed

« back to all changes in this revision

Viewing changes to temporal/t.vect.what.strds/t.vect.what.strds.py

  • Committer: Package Import Robot
  • Author(s): Bas Couwenberg
  • Date: 2015-02-20 23:12:08 UTC
  • mfrom: (8.2.6 experimental)
  • Revision ID: package-import@ubuntu.com-20150220231208-1u6qvqm84v430b10
Tags: 7.0.0-1~exp1
* New upstream release.
* Update python-ctypes-ternary.patch to use if/else instead of and/or.
* Drop check4dev patch, rely on upstream check.
* Add build dependency on libpq-dev to grass-dev for libpq-fe.h.
* Drop patches applied upstream, refresh remaining patches.
* Update symlinks for images switched from jpg to png.

Show diffs side-by-side

added added

removed removed

Lines of Context:
temporal/t.vect.what.strds/t.vect.what.strds.py
 
1
#!/usr/bin/env python
 
2
# -*- coding: utf-8 -*-
 
3
############################################################################
 
4
#
 
5
# MODULE:       t.vect.what.strds
 
6
# AUTHOR(S):    Soeren Gebbert
 
7
#
 
8
# PURPOSE:      Store raster map values at spatial and temporal positions of vector points as vector attributes.
 
9
# COPYRIGHT:    (C) 2011-2014 by the GRASS Development Team
 
10
#
 
11
#               This program is free software under the GNU General Public
 
12
#               License (version 2). Read the file COPYING that comes with GRASS
 
13
#               for details.
 
14
#
 
15
#############################################################################
 
16
 
 
17
#%module
 
18
#% description: Stores raster map values at spatial and temporal positions of vector points as vector attributes.
 
19
#% keyword: temporal
 
20
#% keyword: sampling
 
21
#% keyword: vector
 
22
#%end
 
23
 
 
24
#%option G_OPT_STVDS_INPUT
 
25
#%end
 
26
 
 
27
#%option G_OPT_STRDS_INPUT
 
28
#% key: strds
 
29
#%end
 
30
 
 
31
#%option
 
32
#% key: column
 
33
#% type: string
 
34
#% label: Name of the vector column to be created and to store sampled raster values
 
35
#% description: The use of a column name forces t.vect.what.rast to sample only values from the first map found in an interval. Otherwise the raster map names are used as column names
 
36
#% required: no
 
37
#% multiple: no
 
38
#%end
 
39
 
 
40
#%option
 
41
#% key: method
 
42
#% type: string
 
43
#% description: Aggregate operation to be performed on the raster maps
 
44
#% required: yes
 
45
#% multiple: no
 
46
#% options: disabled,average,count,median,mode,minimum,min_raster,maximum,max_raster,stddev,range,sum,variance,diversity,slope,offset,detcoeff,quart1,quart3,perc90,quantile,skewness,kurtosis
 
47
#% answer: disabled
 
48
#%end
 
49
 
 
50
#%option G_OPT_DB_WHERE
 
51
#%end
 
52
 
 
53
#%option G_OPT_T_WHERE
 
54
#% key: t_where
 
55
#%end
 
56
 
 
57
#%option G_OPT_T_SAMPLE
 
58
#%end
 
59
 
 
60
import os
 
61
import grass.script as grass
 
62
import grass.temporal as tgis
 
63
import grass.script.raster as raster
 
64
from grass.exceptions import CalledModuleError
 
65
 
 
66
############################################################################
 
67
 
 
68
 
 
69
def main():
 
70
 
 
71
    # Get the options
 
72
    input = options["input"]
 
73
    strds = options["strds"]
 
74
    where = options["where"]
 
75
    column = options["column"]
 
76
    method = options["method"]
 
77
    tempwhere = options["t_where"]
 
78
    sampling = options["sampling"]
 
79
 
 
80
    if where == "" or where == " " or where == "\n":
 
81
        where = None
 
82
 
 
83
    # Make sure the temporal database exists
 
84
    tgis.init()
 
85
    # We need a database interface
 
86
    dbif = tgis.SQLDatabaseInterfaceConnection()
 
87
    dbif.connect()
 
88
 
 
89
    sp = tgis.open_old_stds(input, "stvds", dbif)
 
90
    strds_sp = tgis.open_old_stds(strds, "strds", dbif)
 
91
 
 
92
    if strds_sp.get_temporal_type() != sp.get_temporal_type():
 
93
        dbif.close()
 
94
        grass.fatal(_("Input and aggregation dataset must "
 
95
                      "have the same temporal type"))
 
96
 
 
97
    # Check if intervals are present in the sample dataset
 
98
    if sp.get_temporal_type() == "absolute":
 
99
        map_time = sp.absolute_time.get_map_time()
 
100
    else:
 
101
        map_time = sp.relative_time.get_map_time()
 
102
 
 
103
    if map_time != "interval":
 
104
        dbif.close()
 
105
        grass.fatal(_("All registered maps of the space time vector "
 
106
                      "dataset must have time intervals"))
 
107
 
 
108
    rows = sp.get_registered_maps("name,layer,mapset,start_time,end_time",
 
109
                                  tempwhere, "start_time", dbif)
 
110
 
 
111
    if not rows:
 
112
        dbif.close()
 
113
        grass.fatal(_("Space time vector dataset <%s> is empty") % sp.get_id())
 
114
 
 
115
    # Sample the raster dataset with the vector dataset and run v.what.rast
 
116
    for row in rows:
 
117
        start = row["start_time"]
 
118
        end = row["end_time"]
 
119
        vectmap = row["name"] + "@" + row["mapset"]
 
120
        layer = row["layer"]
 
121
 
 
122
        raster_maps = tgis.collect_map_names(
 
123
            strds_sp, dbif, start, end, sampling)
 
124
 
 
125
        aggreagated_map_name = None
 
126
 
 
127
        if raster_maps:
 
128
            # Aggregation
 
129
            if method != "disabled" and len(raster_maps) > 1:
 
130
                # Generate the temporary map name
 
131
                aggreagated_map_name = "aggreagated_map_name_" + \
 
132
                    str(os.getpid())
 
133
                new_map = tgis.aggregate_raster_maps(raster_maps,
 
134
                                                     aggreagated_map_name,
 
135
                                                     start, end, 0, method,
 
136
                                                     False, dbif)
 
137
                aggreagated_map_name = aggreagated_map_name + "_0"
 
138
                if new_map is None:
 
139
                    continue
 
140
                # We overwrite the raster_maps list
 
141
                raster_maps = (new_map.get_id(), )
 
142
 
 
143
            for rastermap in raster_maps:
 
144
 
 
145
                if column:
 
146
                    col_name = column
 
147
                else:
 
148
                    # Create a new column with the SQL compliant
 
149
                    # name of the sampled raster map
 
150
                    col_name = rastermap.split("@")[0].replace(".", "_")
 
151
 
 
152
                coltype = "DOUBLE PRECISION"
 
153
                # Get raster type
 
154
                rasterinfo = raster.raster_info(rastermap)
 
155
                if rasterinfo["datatype"] == "CELL":
 
156
                    coltype = "INT"
 
157
 
 
158
                try:
 
159
                    if layer:
 
160
                        grass.run_command("v.db.addcolumn",
 
161
                                          map=vectmap, layer=layer,
 
162
                                          column="%s %s" % (col_name, coltype),
 
163
                                          overwrite=grass.overwrite())
 
164
                    else:
 
165
                        grass.run_command("v.db.addcolumn", map=vectmap,
 
166
                                          column="%s %s" % (col_name, coltype),
 
167
                                          overwrite=grass.overwrite())
 
168
                except CalledModuleError:
 
169
                    dbif.close()
 
170
                    grass.fatal(_("Unable to add column %s to vector map <%s>")
 
171
                                % (col_name, vectmap))
 
172
 
 
173
                # Call v.what.rast
 
174
                try:
 
175
                    if layer:
 
176
                        grass.run_command("v.what.rast", map=vectmap,
 
177
                                          layer=layer, raster=rastermap,
 
178
                                          column=col_name, where=where)
 
179
                    else:
 
180
                        grass.run_command("v.what.rast", map=vectmap,
 
181
                                          raster=rastermap, column=col_name,
 
182
                                          where=where)
 
183
                except CalledModuleError:
 
184
                    dbif.close()
 
185
                    grass.fatal(_("Unable to run v.what.rast for vector map "
 
186
                                  "<%s> and raster map <%s>") % (vectmap,
 
187
                                                                 rastermap))
 
188
 
 
189
                if aggreagated_map_name:
 
190
                    try:
 
191
                        grass.run_command("g.remove", flags='f', type='raster',
 
192
                                          name=aggreagated_map_name)
 
193
                    except CalledModuleError:
 
194
                        dbif.close()
 
195
                        grass.fatal(_("Unable to remove raster map <%s>")
 
196
                                    % (aggreagated_map_name))
 
197
 
 
198
                # Use the first map in case a column names was provided
 
199
                if column:
 
200
                    break
 
201
 
 
202
    dbif.close()
 
203
 
 
204
if __name__ == "__main__":
 
205
    options, flags = grass.parser()
 
206
    main()