"""
Interpolate a dem from a line shapefile
The shapefile name is passed on the command line
It must have a column named "level" with elevation values
"""

import os, sys
import grass.script as grass

if "GISBASE" not in os.environ:
    print "You must be in GRASS GIS to run this program."
    sys.exit(1)

if len(sys.argv) < 2:
    print 'No shapefile specified'
    sys.exit(1)

# Get shapefile name on command line
shp_path    = sys.argv[1]
shp         = os.path.basename(shp_path)
ctour_orig  = os.path.splitext(shp)[0].lower()
ctour       = ctour_orig.replace(' ','_')
ctour_rast  = ctour + '_rast'
dem_tile    = ctour + '_dem'

print "Processing shapefile: %s" % ctour

grass.run_command('v.in.ogr', dsn=shp_path, output=ctour, type='line', overwrite=True)
grass.run_command('g.region', flags='p', vect=ctour, overwrite=True)
grass.run_command('v.to.rast', input=ctour, output=ctour, type='line', use='attr', attrcol='level', overwrite=True)
grass.run_command('r.surf.contour', input=ctour, output=dem_tile, overwrite=True)

print "Completed creating raster tile: %s" % dem_tile

for s in *.shp; do python shp_to_dem.py $s; done