Talk:FUTURES land-change modeling for evaluating innovative conservation scenarios

From GRASS-Wiki
Jump to navigation Jump to search
The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.

Python: Loss of farmland and forest

import numpy as np
import grass.script as gscript 
 
# list maps computed by r.futures.parallelpga for one scenario
current_mapset = gscript.gisenv()['MAPSET']
# put the name you are using for the scenarios in the pattern option 
maps = gscript.list_grouped('raster', pattern="final*")[current_mapset]
developed = []
forest = []
agriculture = []
other = []
for each in maps:
    # compute landcover in 2035 by updating landcover 2011 with the simulation results:
    gscript.mapcalc("landuse_2035_{name} = if({name} >=1, 21, landuse_2011)".format(name=each), quiet=True, overwrite=True)
    # get information about how many cells there are of each land use category, use r.stats
    res = gscript.read_command('r.stats', flags='cn', input="landuse_2035_{name}".format(name=each), quiet=True).strip()
    for line in res.strip().splitlines():
        cat, cells = line.split()
        cat = int(cat)
        cells = int(cells)
        # developed
        if 21 <= cat <= 24:
            developed.append(cells)
        # agriculture
        elif cat == 81 or cat == 82:
            agriculture.append(cells)
        # forest
        elif 41 <= cat <= 43:
            forest.append(cells)
        # other land use
        else:
            other.append(cells)
 
# convert to NumPy arrays
forest = np.array(forest)
agriculture = np.array(agriculture)
other = np.array(other)
developed = np.array(developed)
 
print "forest: {mean} +- {std}".format(mean=np.mean(forest), std=np.std(forest))
print "agriculture: {mean} +- {std}".format(mean=np.mean(agriculture), std=np.std(agriculture))
print "other: {mean} +- {std}".format(mean=np.mean(other), std=np.std(other))
print "developed: {mean} +- {std}".format(mean=np.mean(developed), std=np.std(developed))

Python: Forest fragmentation

import grass.script as gscript 

def forest_fragmentation(landuse):
    gscript.mapcalc("forest_2035 = if({landuse} >= 41 && {landuse} <= 43, 1, 0)".format(landuse=landuse))
    gscript.run_command('r.forestfrag', flags='a', input='forest_2035', output='fragment', window=25)
    result = gscript.parse_command('r.stats', input='fragment', flags='cn', parse=(gscript.parse_key_val, {'sep': ' ', 'val_type': int}))
    gscript.run_command('g.remove', type='raster', name=['forest_2035', 'fragment'], flags='f')
    return result

current_mapset = gscript.gisenv()['MAPSET']
landuse_maps = gscript.list_grouped('raster', pattern="landuse_2035_*")[current_mapset]
landuse_maps = landuse_maps[:5]
for landuse in landuse_maps:
    result = forest_fragmentation(landuse)
    print result