Introduction to GRASS GIS and Scripting in GRASS using Python

FOSS4G-PH 2016 workshop (April 22, 2016)

MH 210, College of Engineering, University of the Philippines Diliman, Quezon City 1191

Engr. Ben Hur S. Pintor

     https://github.com/bnhr07b

     [email protected]

Python Scripting in GRASS

Let's say I want to compute for the daily solar radiation values for January then average them to get the monthly average value. Instead of running r.sun manually 31 times, we can create a script that only needs to be run once.

There are two (2) libraries that will allow us to script using Python in GRASS. These are the Python Scripting Library and the PyGRASS module. Let's try to create the same script using both methods.

The Python Scripting Library

#!/usr/bin/env python

import grass.script as gscript      # import the Python Scripting Library

def main():
    startday = 1    # the Julian date when to start the computations
    endday = 31     # the Julian date when to end the computations
    basename = "psl_clear_GHI"      # basename of the created solar radiation maps
    output = "psl_clear_GHI_ave"    # name of the average map
    calc = ""   # variable to hold the mapcalc calculation

    for day in range(startday, endday+1):   # run r.sun from startday to endday

        ghi = "{}_{}".format(basename, day)    # create a new name for the output ghi map

        gscript.run_command('r.sun', elevation='elev', aspect='aspect', 
            slope='slope', linke='linke_Jan',horizon_basename='horizon', 
            horizon_step='15', glob_rad=ghi, day=day)   # run r.sun

    '''Create the mapcalc calculation string'''
    for day in range(startday, endday):
        calc += "{}_{} + ".format(basename, day)

    calc += "{}_{}".format(basename, endday)

    gscript.mapcalc('{} = {}'.format(output, calc))     # Perform the mapcalc


if __name__ == "__main__":
    main()

PyGRASS

#!/usr/bin/env python

from grass.pygrass.modules import Module    # import Module from PyGRASS

r_sun = Module('r.sun')     # r_sun is the r.sun module
r_mapcalc = Module('r.mapcalc')     #r_mapcalc is the r.mapcalc module

def main():
    startday = 1    # the Julian date when to start the computations
    endday = 31     # the Julian date when to end the computations
    basename = "pygrass_clear_GHI"      # basename of the created solar radiation maps
    output = "pygrass_clear_GHI_ave"    # name of the average map
    calc = ""   # variable to hold the mapcalc calculation

    for day in range(startday, endday+1):   # run r.sun from startday to endday

        ghi = "{}_{}".format(basename, day)    # create a new name for the output ghi map

        r_sun(elevation='elev', aspect='aspect', slope='slope', 
            linke='linke_Jan', horizon_basename='horizon', horizon_step='15', 
            glob_rad=ghi, day=day, verbose=True, overwrite=True)    # run r.sun

    '''Create the mapcalc calculation string'''
    for day in range(startday, endday):
        calc += "{}_{} + ".format(basename, day)

    calc += "{}_{}".format(basename, endday)

    r_mapcalc(expression="{} = {}".format(output, calc))


if __name__ == "__main__":
    main()

Adding a GRASS module interface to the script

#!/usr/bin/env python

'''This script doesn't work when run from the terminal. 
It must be run from File -> Launch script menu'''

'''The parts below define the module interface.'''
#%module
#% description: Computes daily solar radiation from start day to end day and averages them
#% keyword: raster
#% keyword: r.sun
#% keyword: average
#%end
#%option
#% key: startday
#% type: integer
#% description: The Julian day to start computations
#% required : yes
#% answer: 1
#%end
#%option
#% key: endday
#% type: integer
#% description: The Julian day to stop computations
#% required : yes
#% answer: 31
#%end
#%option
#% key: basename
#% type: string
#% description: basename for daily GHI maps
#% required : no
#% answer: clear_GHI
#%end
#%option
#% key: output
#% type: string
#% description: The name of the output average raster
#% required : no
#% answer: clear_GHI_ave
#%end


import sys
import grass.script as gscript      # import the Python Scripting Library
from grass.pygrass.modules import Module    # import Module from PyGRASS

r_sun = Module('r.sun')     # r_sun is the r.sun module
r_mapcalc = Module('r.mapcalc')     #r_mapcalc is the r.mapcalc module

def main():
    '''Set options (corresponds to the options stated above)'''
    startday = int(options['startday'])     # the Julian date when to start the computations
    endday = int(options['endday'])         # the Julian date when to end the computations 
    basename = options['basename']          # basename of the created solar radiation maps
    output = options['output']              # name of the average map
    calc = ""   # variable to hold the mapcalc calculation

    for day in range(startday, endday+1):   # run r.sun from startday to endday

        ghi = "{}_{}".format(basename, day)    # create a new name for the output ghi map

        r_sun(elevation='elev', aspect='aspect', slope='slope', 
            linke='linke_Jan', horizon_basename='horizon', horizon_step='15', 
            glob_rad=ghi, day=day, verbose=True, overwrite=True)    #run r.sun

    '''Create the mapcalc calculation string'''
    for day in range(startday, endday):
        calc += "{}_{} + ".format(basename, day)

    calc += "{}_{}".format(basename, endday)

    r_mapcalc(expression="{} = {}".format(output, calc))


if __name__ == "__main__":
    options, flags = gscript.parser()
    main()

For more information, visit:

• GRASS GIS Python Library Documentation

---

Sitemap

Introduction | Raster Processing and Analysis | Vector Processing and Analysis | WORKSHOP 1 | WORKSHOP 2 (Scripting)