[GRASS-user] Watershed delineation

[Hamish]

Daniel Victoria wrote:
> > I trying to delineate watersheds for several points in a very big
> > DEM
> >
> >  | Data Type: CELL |
> >  | Rows: 4803 |
> >  | Columns: 5172 |
> >  | Total Cells: 24841116
> >
> > I've successfully ran r.terraflow on the DEM but now I'm stuck.
> > How do I get the watershed area for my points of interest? I've
> > tried using the direction map produced with r.water.outlet but it
> > did not work (I used MFD flow calculation).
> >
> > Also, from what I understood, the sink watershed map is suposed to
> > be the watersheds delineated by the terraflow module, is that
> > correct?

Near-to, but not exactly.

> > Anyway, I think I'll let r.watershed run for the night

For a region that big r.watershed will take a long time to finish, but
it can be worth the wait. IIRC you get more accurate results from
r.watershed than from r.terraflow. (see historical post by ?
Helena/Laura/Charles?)

The hardest part with r.watershed is getting the threshold setting near
to the scale you wish to work with.

For watershed delimination see (new) examples here:
 http://grass.ibiblio.org/grass63/manuals/html63_user/r.watershed.html
 http://grass.gdf-hannover.de/wiki/Psmap_flooding_example


Jaime Carrera wrote:
> > that's right; you end up with a lot of sub-basins, which you can
> > aggregate according to your need. They are numbered on a logical
> > sequence, so it's really easy to aggregate them with r.mapcalc;
> > just check where a key outlet is and go upstream from there.

Daniel:
> Well, I think that something has gone wrong then because the basin
> numbering does not seam to be in a logical sequence. Or I have not
> understood what is the logic in the numbering sequence.

AFAIR Jamie is correct, there is some logical sequence there. I don't
see much about it in the help page but if you search the mailing list
archive or source code you can find it.

> Just for the record, here is my g.region -p output:
...
> nsres:      899.99917697
> ewres:      899.9166173
> rows:       4803
> cols:       5172
...
> The resolution is strange because the DEM was projected from a 30
> arcsec latlon map. I tried to set it to 900 meters but the output was
> like this...

I think there is a way to have r.proj use the target location's region
setting (region cropping flag?). The ugly way to clean the resolution
is to force a 900m res with 'g.region res=900 -a'



good luck,
Hamish



      ____________________________________________________________________________________
Never miss a thing.  Make Yahoo your home page. 
http://www.yahoo.com/r/hs