Python Ctypes Examples: Difference between revisions
Jump to navigation
Jump to search
(Created page with "== Accessing geometry == This script switches X, Y coordinates and multiple them by -1. <source lang="python"> import sys import grass.script as grass from grass.lib.vector i...") |
(note what's grass7) |
||
(10 intermediate revisions by 4 users not shown) | |||
Line 1: | Line 1: | ||
== Accessing geometry == | * See [[GRASS and Python]] for more information. | ||
* See also [http://trac.osgeo.org/grass/browser/grass/trunk/doc/python/raster_example_ctypes.py raster], [http://trac.osgeo.org/grass/browser/grass/trunk/doc/python/vector_example_ctypes.py vector] sample scripts. (GRASS 7) | |||
__TOC__ | |||
== Raster == | |||
=== Working with (geographic) coordinates for each raster cell and the computational region === | |||
The ctypes interfaces are thin wrappers around the corresponding C | |||
functions. Any script which uses them should have essentially the same | |||
structure as a module written in C. Structures can be allocated by | |||
using the structure name as a function, and passed by reference using | |||
[http://docs.python.org/library/ctypes.html ctypes]' byref() function. | |||
The second argument to G_col_to_easting() and G_row_to_northing() is a | |||
"const struct Cell_head *". The information can be obtained from | |||
G_get_window(), e.g.: | |||
<source lang="python"> | |||
import sys | |||
import grass.lib.gis as gis | |||
from ctypes import * | |||
gis.G_gisinit(sys.argv[0]) | |||
region = gis.Cell_head() | |||
gis.G_get_window(byref(region)) | |||
coor_col = gis.G_col_to_easting((thecolumn + 0.5), byref(region)) | |||
coor_row = gis.G_row_to_northing((therow + 0.5), byref(region)) | |||
</source> | |||
=== Basic raster metadata access methods === | |||
This is a simple demonstration script to show the ctypes style access to some of the raster metadata. | |||
(''GRASS 7'') | |||
<source lang="python"> | |||
from ctypes import * | |||
import grass.lib.gis as libgis | |||
import grass.lib.raster as libraster | |||
import grass.script as grass | |||
import grass.temporal as tgis | |||
def print_raster_info(name): | |||
grass.run_command("r.mapcalc", expr=name + " = 5.0", overwrite=True) | |||
grass.run_command("r.timestamp", map=name, date="13 Jan 2003 / 15 Jan 2003") | |||
# Create a region struct | |||
region = libgis.Cell_head() | |||
# Read the raster region | |||
libraster.Rast_get_cellhd("test", "", byref(region)) | |||
print region.cols | |||
print region.rows | |||
print region.north | |||
print region.south | |||
print region.east | |||
print region.west | |||
print region.ns_res | |||
print region.ew_res | |||
print region.top | |||
print region.bottom | |||
print region.tb_res | |||
print region.depths | |||
print region.proj | |||
# Get Map type | |||
maptype = libraster.Rast_map_type(name, "") | |||
if maptype == libraster.DCELL_TYPE: | |||
print "DCELL" | |||
elif maptype == libraster.FCELL_TYPE: | |||
print "FCELL" | |||
elif maptype == libraster.CELL_TYPE: | |||
print "CELL" | |||
else: | |||
print "Unknown type" | |||
# Read range | |||
if libraster.Rast_map_is_fp(name, ""): | |||
print "Floating point map" | |||
range = libraster.FPRange() | |||
libraster.Rast_read_fp_range(name, "", byref(range)) | |||
min = libgis.DCELL() | |||
max = libgis.DCELL() | |||
libraster.Rast_get_fp_range_min_max(byref(range), byref(min), byref(max)) | |||
else: | |||
print "Cell map" | |||
range = libraster.Range() | |||
libraster.Rast_read_range(name, "", byref(range)) | |||
min = libgis.CELL() | |||
max = libgis.CELL() | |||
libraster.Rast_get_fp_range_min_max(byref(range), byref(min), byref(max)) | |||
minimum = min.value | |||
maximum = max.value | |||
print minimum | |||
print maximum | |||
# Read time stamp | |||
tstmp = libgis.TimeStamp() | |||
libgis.G_read_raster_timestamp(name, "", byref(tstmp)) | |||
print tstmp.count | |||
dt1 = tstmp.dt[0] | |||
dt2 = tstmp.dt[1] | |||
if tstmp.count >= 1: | |||
print dt1.mode | |||
print dt1.year | |||
print dt1.month | |||
print dt1.day | |||
print dt1.hour | |||
print dt1.minute | |||
print dt1.second | |||
if tstmp.count > 1: | |||
print dt2.mode | |||
print dt2.year | |||
print dt2.month | |||
print dt2.day | |||
print dt2.hour | |||
print dt2.minute | |||
print dt2.second | |||
if __name__ == "__main__": | |||
print_raster_info("test") | |||
</source> | |||
== Vector == | |||
=== Accessing feature geometry === | |||
This script switches X, Y coordinates and multiple them by -1. | This script switches X, Y coordinates and multiple them by -1. | ||
<source lang="python"> | <source lang="python"> | ||
#!/usr/bin/python | |||
import sys | import sys | ||
from grass.lib.vector import * | from grass.lib.vector import * | ||
Line 26: | Line 165: | ||
nlines = Vect_get_num_lines(map_info) | nlines = Vect_get_num_lines(map_info) | ||
for line in range(1, nlines): | for line in range(1, nlines + 1): | ||
ltype = Vect_read_line(map_info, points, cats, line) | ltype = Vect_read_line(map_info, points, cats, line) | ||
if line % 100 == 0: | if line % 100 == 0: | ||
Line 39: | Line 178: | ||
points, cats) | points, cats) | ||
line += 1 | line += 1 | ||
sys.stderr.write("\r") | |||
Vect_destroy_line_struct(points) | Vect_destroy_line_struct(points) | ||
Vect_destroy_cats_struct(cats) | Vect_destroy_cats_struct(cats) | ||
Vect_build_partial(map_info, GV_BUILD_NONE) | |||
Vect_build(map_info) | |||
Vect_close(map_info) | Vect_close(map_info) | ||
</source> | </source> | ||
== See also == | |||
* http://docs.python.org/library/ctypes.html | |||
{{Python}} |
Latest revision as of 10:02, 20 May 2012
- See GRASS and Python for more information.
Raster
Working with (geographic) coordinates for each raster cell and the computational region
The ctypes interfaces are thin wrappers around the corresponding C functions. Any script which uses them should have essentially the same structure as a module written in C. Structures can be allocated by using the structure name as a function, and passed by reference using ctypes' byref() function.
The second argument to G_col_to_easting() and G_row_to_northing() is a "const struct Cell_head *". The information can be obtained from G_get_window(), e.g.:
import sys
import grass.lib.gis as gis
from ctypes import *
gis.G_gisinit(sys.argv[0])
region = gis.Cell_head()
gis.G_get_window(byref(region))
coor_col = gis.G_col_to_easting((thecolumn + 0.5), byref(region))
coor_row = gis.G_row_to_northing((therow + 0.5), byref(region))
Basic raster metadata access methods
This is a simple demonstration script to show the ctypes style access to some of the raster metadata.
(GRASS 7)
from ctypes import *
import grass.lib.gis as libgis
import grass.lib.raster as libraster
import grass.script as grass
import grass.temporal as tgis
def print_raster_info(name):
grass.run_command("r.mapcalc", expr=name + " = 5.0", overwrite=True)
grass.run_command("r.timestamp", map=name, date="13 Jan 2003 / 15 Jan 2003")
# Create a region struct
region = libgis.Cell_head()
# Read the raster region
libraster.Rast_get_cellhd("test", "", byref(region))
print region.cols
print region.rows
print region.north
print region.south
print region.east
print region.west
print region.ns_res
print region.ew_res
print region.top
print region.bottom
print region.tb_res
print region.depths
print region.proj
# Get Map type
maptype = libraster.Rast_map_type(name, "")
if maptype == libraster.DCELL_TYPE:
print "DCELL"
elif maptype == libraster.FCELL_TYPE:
print "FCELL"
elif maptype == libraster.CELL_TYPE:
print "CELL"
else:
print "Unknown type"
# Read range
if libraster.Rast_map_is_fp(name, ""):
print "Floating point map"
range = libraster.FPRange()
libraster.Rast_read_fp_range(name, "", byref(range))
min = libgis.DCELL()
max = libgis.DCELL()
libraster.Rast_get_fp_range_min_max(byref(range), byref(min), byref(max))
else:
print "Cell map"
range = libraster.Range()
libraster.Rast_read_range(name, "", byref(range))
min = libgis.CELL()
max = libgis.CELL()
libraster.Rast_get_fp_range_min_max(byref(range), byref(min), byref(max))
minimum = min.value
maximum = max.value
print minimum
print maximum
# Read time stamp
tstmp = libgis.TimeStamp()
libgis.G_read_raster_timestamp(name, "", byref(tstmp))
print tstmp.count
dt1 = tstmp.dt[0]
dt2 = tstmp.dt[1]
if tstmp.count >= 1:
print dt1.mode
print dt1.year
print dt1.month
print dt1.day
print dt1.hour
print dt1.minute
print dt1.second
if tstmp.count > 1:
print dt2.mode
print dt2.year
print dt2.month
print dt2.day
print dt2.hour
print dt2.minute
print dt2.second
if __name__ == "__main__":
print_raster_info("test")
Vector
Accessing feature geometry
This script switches X, Y coordinates and multiple them by -1.
#!/usr/bin/python
import sys
from grass.lib.vector import *
if len(sys.argv) < 2:
sys.exit("Usage: %s: vector" % sys.argv[0])
name = sys.argv[1]
map_info = pointer(Map_info())
points = Vect_new_line_struct()
cats = Vect_new_cats_struct()
c_points = points.contents
level = Vect_open_update(map_info, name, "")
if level < 2:
sys.exit("Topology not available")
nlines = Vect_get_num_lines(map_info)
for line in range(1, nlines + 1):
ltype = Vect_read_line(map_info, points, cats, line)
if line % 100 == 0:
sys.stderr.write("\r%d" % line)
for i in range(c_points.n_points):
x = c_points.x[i]
c_points.x[i] = -1 * c_points.y[i]
c_points.y[i] = -1 * x
Vect_rewrite_line(map_info, line, ltype,
points, cats)
line += 1
sys.stderr.write("\r")
Vect_destroy_line_struct(points)
Vect_destroy_cats_struct(cats)
Vect_build_partial(map_info, GV_BUILD_NONE)
Vect_build(map_info)
Vect_close(map_info)