Python Swig Examples

Advanced SWIG Python Interface Examples

• Simpler examples are given in the SWIG-Python README file and the GRASS and Python wiki page.

Overview

As you might expect, GRASS's C library functions preform many common GIS tasks using very well tested and mature code. Reimplementing them in Python would thus be a waste of programming effort and would introduce potentially buggy and slower code. (Well written C is hard to beat for speed) A major advantage of using the SWIG interface can be summarized as "the experts have already solved this problem for you, take advantage of that."

Passing arrays and specific data types

Some of GRASS's functions want to be passed arrays of numbers or of a specific data type (e.g. integer, double). To do this from Python the 3rd party NumPtr module is used to access a memory pointer object which can be passed to the function.

(any way to then free() the memory? Does "del SomePtr" do it?)

Setting up the NumPtr module

NumPtr setup:

# NumPtr - Numeric Pointer Module for Python  (GPL2)
# http://geosci.uchicago.edu/csc/numptr/
#   23k .tgz ; 100k installed

wget http://geosci.uchicago.edu/csc/numptr/NumPtr-1.1.tar.gz

tar xzf NumPtr-1.1.tar.gz
cd NumPtr-1.1
python setup.py build
#python setup.py install --prefix=/home/user
cp build/lib.linux-i686-2.4/*NumPtr.* /usr/src/grass63/swig/python/

Examples

Distance and area calculations

Once NumPtr is installed we can run our script.

GRASS's libgis (C API) distance and area functions automatically switch to using geodetic calculations when in a Lat/Lon location.

The following calculates the area of the default Spearfish region bounds and the distance between the region's far corners.

Spearfish uses a UTM (planimetric) projection.

m.distance

#!/usr/bin/python
############################################################################
#
# MODULE:       m.distance
#
# AUTHOR(S):    Hamish Bowman, Dunedin, New Zealand
#
# PURPOSE:      Find distance between two points
#		Demonstrates GRASS SWIG-Python interace
#
# COPYRIGHT:    (c) 2008 Hamish Bowman, and 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.
#
############################################################################

#%Module
#%  description: Find distance between two points
#%  keywords: miscellaneous, distance, measure
#%End
#%Option
#%  key: coord
#%  type: string
#%  required: yes
#%  multiple: no
#%  key_desc: x1,y1,x2,y2
#%  description: Starting and ending coordinates
#%End

import sys
import os

# run this before starting python to append module search path:
#   export PYTHONPATH=/usr/src/grass63/swig/python
#   check with "import sys; sys.path"
# or:
sys.path.append("/usr/src/grass63/swig/python")
import python_grass6 as g6lib


def main(): 

    #### add your code here ####

    #os.putenv("GIS_OPT_COORD", '609000,4913700,589980,4928010')
    coords = os.getenv("GIS_OPT_COORD").split(',')

    x1 = float(coords[0])
    y1 = float(coords[1])
    x2 = float(coords[2])
    y2 = float(coords[3])


    g6lib.G_gisinit('m.distance')
    # returns 0 on success


    ### calc distance ###
    
    g6lib.G_begin_distance_calculations()
    # returns 0 if projection has no metrix (ie. imagery)
    # returns 1 if projection is planimetric
    # returns 2 if projection is latitude-longitude
    
    distance = g6lib.G_distance(x1, y1, x2, y2)
    print "distance is", distance
    # 23802.0272246   (ok, matches above calc.)


    ### calc area ###
    
    g6lib.G_begin_polygon_area_calculations()
    # returns 0 if the projection is not measurable (ie. imagery or xy)
    # returns 1 if the projection is planimetric (ie. UTM or SP)
    # returns 2 if the projection is non-planimetric (ie. latitude-longitude)
    
    
    # we don't need this, but just to have a look
    g6lib.G_database_units_to_meters_factor()
    # 1.0
    
    # passing an array of values
    import Numeric
    import NumPtr
    
    # do not need to close polygon (but it doesn't hurt if you do)
    x = [ x1, x2, x2, x1 ]
    y = [ y1, y1, y2, y2 ]
    npoints = len(x)
    
    # unset variables:
    #del [Xs, Xptr, Ys, Yptr]
    #   or
    #Xs = Xptr = Ys = Yptr = None
    
    Xs = Numeric.array(x, Numeric.Float64)
    Xptr = NumPtr.getpointer(Xs)
    Ys = Numeric.array(y, Numeric.Float64)
    Yptr = NumPtr.getpointer(Ys)
    
    area = g6lib.G_area_of_polygon(Xptr, Yptr, npoints)
    print "area is", area
    # 272176200.0   (ok, matches above calc)

    #### end of your code ####
    return 

if __name__ == "__main__":
    if ( len(sys.argv) <= 1 or sys.argv[1] != "@ARGS_PARSED@" ):
        os.execvp("g.parser", [sys.argv[0]] + sys.argv)
    else:
	main();