1

I use my python script to call GRASS module r.water.outlet from outside, and here's my script(I perform this script in pyscripter):

import os
import sys

#set up GRASS environment variables
sys.path.append(os.path.join(os.environ['GISBASE'], 'etc', 'python'))
import grass.script as g
import grass.script.setup as gsetup
gisbase = os.environ['GISBASE']
gisdb = 'C:\Users\Heinz\Documents\grassdata'
location = 'newLocation'
mapset = 'TC'
gsetup.init(gisbase, gisdb, location, mapset)

#csv file reading and importation
import csv
rows = list(open('tc_sta.csv'))
totalrows = len(rows) - 1

#loop through the csv(coordinates) file in r.water.outlet module
f = open('tc_sta.csv', 'r')
element = list(csv.reader(f))
i = 0
j = 0
while True:
    if i <= totalrows:
        g.run_command('r.water.outlet', drainage = 'dra', basin = 'b' + str(i + 1), east = element[i][j], north = element[i][j + 1])
        i = i + 1
    else:
        break

#print all files in grass database to check results of the module
print g.read_command('g.list', _type='rast')

and here's the output in console,

enter image description here

The output shows that the GRASS module successfully runs and generates 8 new basin area raster files(b1 to b8), but I still want to check if there are something wrong. So I launch GRASS GUI (GRASS command line can't display map under Windows platform), using d.rast command to get visual of those output maps.

But I can only see b4 in display window. Others' are blank in display window. I 've checked the attribute of every output map using r.category, and they all have the same value as b4's.

I display b4 and input coordinates of module r.water.outlet(thus, batch point) together to check if other coordinates except b4's are out of range of basin area, but it shows not,

enter image description here

I am wondering that are there better ways to check this?(Seems that I can't use python from outside to display maps.) Because somehow maybe the outputs are all right.

I work with GRASS 6.4.3 and python 2.7.6 under Windows 8.1 64-bit.

1 Answer 1

1

Going back to basics:

  1. I don't see a region setting in your script. You should first do a g.run_comand('g.region',rast='dra').
  2. Next, are you sure that the points in the csv file fall exactly on the streams? If an outlet point is even slightly off the stream channel, you will get very tiny basins.
1
  • you're right! I have to extract coordinates from splitting the output file 'drainage' then I could get all fine results of r.water.outlet, right? Thus I even don't have to check it, right?
    – Heinz
    Mar 25, 2014 at 9:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.