Dear list,<br><br>I am interested in finding basin area limited to say 300sq km.<br><br>1. Easier approach: is As suggested by Mr. M Metz earlier. It's quite close.<br><br>Use SFD in r.watershed. Keep threshold in terms of number of cells equal to the area required. Use a raster dem 1 cells (r.grow) bigger than the basin boundary to get positive accumulation values.<br>
<br>Then run r.stream.extract. Keep threshold equal to the number used above.<br>Find starting points . Use layer 2 option and type_code=0.<br><br>These points are having the catchment equal to the threshold. (at least, as far as I have checked. I am still on job.)<br>
<br>The river starts forming when the given threshold is reached.<br><br><br>2. Not so easy approach:<br><br>--------------------------------------------------------------------------------------------------<br>echo "delete from cat_temp_done;" > file.sql;<br>
db.execute input=file.sql;<br><br>r.watershed -s --overwrite -m memory=12000 elevation=ybdem@PERMANENT threshold=12345 accumulation=step1_acc drainage=step1_fdir basin=step1_basin stream=step1_stream<br><br>r.stream.extract --overwrite memory=12000 elevation=ybdem@PERMANENT accumulation=step1_acc threshold=12345 stream_length=2 stream_rast=step1_stream_rse stream_vect=step1_stream_rse direction=step1_stream_fdir<br>
<br>v.db.addcolumn map=step1_stream_rse@PERMANENT layer=1 columns="basin_area double precision,basin_area_point double precision"<br><br>v.db.addtable map=step1_stream_rse@PERMANENT layer=2<br><br><br>#find start points<br>
<br>v.extract -t --verbose input=step1_stream_rse@PERMANENT layer=2 type=point cats=0 output=zero_points_delete<br>v.category --overwrite --verbose input=zero_points_delete@PERMANENT output=zero_points_no_cat_delete option=del<br>
v.category --overwrite --verbose input=zero_points_no_cat_delete@PERMANENT output=zero_points_with_cat_delete option=add<br>v.select --overwrite --verbose ainput=step1_stream_rse@PERMANENT alayer=1 atype=point binput=zero_points_with_cat_delete@PERMANENT blayer=1 btype=point output=zero_points@PERMANENT operator=equals<br>
<br>#find outlets<br><br>v.extract -t --verbose input=step1_stream_rse@PERMANENT layer=2 type=point cats=2 output=outlet_points_delete<br>v.category --overwrite --verbose input=outlet_points_delete@PERMANENT output=outlet_points_no_cat_delete option=del<br>
v.category --overwrite --verbose input=outlet_points_no_cat_delete@PERMANENT output=outlet_points_with_cat_delete option=add<br>v.select --overwrite --verbose ainput=step1_stream_rse@PERMANENT alayer=1 atype=point binput=outlet_points_with_cat_delete@PERMANENT blayer=1 btype=point output=outlet_points@PERMANENT operator=equals<br>
<br><br><br>start_cats=`echo "select cat from zero_points_step1_stream_rse;" | db.select -c;`<br>end_cat=`echo "select cat from outlet_points_step1_stream_rse;" | db.select -c;`<br>for i in $start_cats<br>
do<br>for j in $end_cat<br>do<br>echo "i:$i, j:$j"<br>v.net.salesman --overwrite input=step1_stream_rse@PERMANENT output=river_path_$i type=line alayer=1 nlayer=1 ccats=$i,$j<br>v.db.addtable map=river_path_$i@PERMANENT columns="cat integer, basin_area double precision, basin_area_point double precision"<br>
path_line_cats=`echo "select cat from river_path_$i;" | db.select -c;`<br>prev_area_value=0;<br>for x in $path_line_cats<br>do<br>echo "i:$i, j:$j, x: $x"<br>r.stream.basins dirs=step1_stream_fdir@PERMANENT streams=step1_stream_rse@PERMANENT cats=$x basins=river_path_basins_$i_$x<br>
area_value=`r.stats -a -n --verbose input=river_path_basins_$i_$x@PERMANENT | awk '{print $2}'`<br>echo "Update river_path_$i set basin_area=$area_value where cat=$x;" > file.sql;<br>db.execute input=file.sql;<br>
area_value_point=$area_value-$prev_area_value<br>echo "Update river_path_$i set basin_area_point=$area_value_point where cat=$x;" > file.sql;<br>db.execute input=file.sql;<br>$prev_area_value=$area_value<br>
done<br>line_cats=`echo "select cat from river_path_$i and cat not in (select cat from cat_temp_done);" | db.select -c;`<br>for k in $line_cats<br>do<br>echo "i:$i, j:$j, x: $k"<br>r.stream.basins dirs=step1_stream_fdir@PERMANENT streams=step1_stream_rse@PERMANENT cats=$k basins=river_path_basins_$i_$k<br>
area_value=`r.stats -a -n --verbose input=river_path_basins_$i_$k@PERMANENT | awk '{print $2}'`<br>echo "Update step1_stream_rse set basin_area=$area_value where cat=$k;" > file.sql;<br>db.execute input=file.sql;<br>
area_value_point=$area_value-$prev_area_value<br>echo "Update step1_stream_rse set basin_area_point=$area_value_point where cat=$k;" > file.sql;<br>db.execute input=file.sql;<br>$prev_area_value=$area_value<br>
echo "insert into gis_schema.cat_temp_done (cat) values($k);" > file.sql;<br>db.execute input=file.sql;<br>done<br>done<br>done<br><br>___________________________________________________________________________________________________--<br>
<br>The above codes will neatly arrange the required values in the columns. Which I shall later use on as label on points on river path.<br><br>3. Not so easy -- Not so easy approach<br><br>Replace the one line code of v.net.salesman used above with this:<br>
<br>--------------------------------------------------------------------------------------------------------------------<br>elcat=1972;<br>echo "delete from work.cat_temp;" > file.sql;<br>db.execute input=file.sql;<br>
echo "delete from work.cat_temp_done;" > file.sql;<br>db.execute input=file.sql;<br>echo "delete from work.yrse100sqkm_stream_path1;" > file.sql;<br>db.execute input=file.sql;<br>while [ "$elcat" ]<br>
do<br>cats_to_write=`echo "select cat from work.yrse100sqkm_stream where end_long=(select start_long from work.yrse100sqkm_stream where cat=$elcat) and end_lat=(select start_lat from work.yrse100sqkm_stream where cat=$elcat);" | db.select -c;`<br>
echo "delete from work.cat_temp where cat= $elcat;" > file.sql;<br>db.execute input=file.sql; <br>echo "cats_to_write: $cats_to_write";<br>if [ ! -n "$cats_to_write" ]<br>then<br>echo "Update work.yrse100sqkm_stream set typecode=0 where cat=$elcat;" > file.sql;<br>
db.execute input=file.sql; <br>else<br>echo "Update work.yrse100sqkm_stream set typecode=1 where cat=$elcat;" > file.sql;<br>db.execute input=file.sql; <br>for k in $cats_to_write<br>do<br>echo "insert into work.cat_temp (cat) values($k);" > file.sql;<br>
db.execute input=file.sql; <br>prev_path=9999999;<br>prev_path=`echo "select cat from work.yrse100sqkm_stream_path1 where path_cat=$elcat; " | db.select -c;`<br>for p in $prev_path<br>do<br>echo "p:k: $p,$k";<br>
echo "insert into work.yrse100sqkm_stream_path1 (cat,path_cat) values($p,$k);" > file.sql;<br>db.execute input=file.sql; <br>done<br>echo "insert into work.yrse100sqkm_stream_path1 (cat,path_cat) values($elcat,$k);" > file.sql;<br>
db.execute input=file.sql; <br>done<br>fi<br>elcat=`echo "select cat from work.cat_temp limit 1;" | db.select -c;`<br>done<br>------------------------------------------------------------------------------------------------------------------------------------------------------------------------<br>
<br>or with this<br><br><br>--------------------------------------------------------------------------------------------------------------------------------------------------------------------------<br><br>elcat=1972;<br>
echo "delete from work.cat_temp;" > file.sql;<br>db.execute input=file.sql;<br>echo "delete from work.cat_temp_done;" > file.sql;<br>db.execute input=file.sql;<br>while [ $elcat ]<br>do<br>cats=`echo "select cat from work.yrse100sqkm_stream where end_long=(select start_long from work.yrse100sqkm_stream where cat=$elcat) and end_lat=(select start_lat from work.yrse100sqkm_stream where cat=$elcat) and cat <> $elcat; " | db.select -c;`<br>
if [ ! -n "$cats" ]<br> then<br> echo "insert into work.river_startline(cat) values($elcat);" > file.sql;<br> db.execute input=file.sql;<br>fi<br>echo "insert into work.cat_temp_done (cat) values($elcat);" > file.sql;<br>
db.execute input=file.sql; <br>cats_to_write=`echo "select cat from work.yrse100sqkm_stream where end_long=(select start_long from work.yrse100sqkm_stream where cat=$elcat) and end_lat=(select start_lat from work.yrse100sqkm_stream where cat=$elcat) and cat not in (select cat from work.cat_temp_done); " | db.select -c;`<br>
echo "delete from work.cat_temp where cat= $elcat;" > file.sql;<br>db.execute input=file.sql; <br>for k in $cats_to_write<br>do<br>echo "insert into work.cat_temp (cat) values($k);" > file.sql;<br>
db.execute input=file.sql; <br>done<br>elcat=`echo "select cat from work.cat_temp limit 1;" | db.select -c;`<br>done<br><br>and then<br><br>echo "delete from yrse100sqkm_stream_path;" > file.sql;<br>
db.execute input=file.sql;<br>elcats=`echo "select cat from work.river_startline;" | db.select -c;`<br>for i in $elcats<br>do<br>cats=$i<br>declare -i order;<br>order=1;<br>while [ $cats ]<br>do<br>echo "insert into work.yrse100sqkm_stream_path (cat,path_cat,orderno) values($cats,$i,$order);" > file.sql;<br>
db.execute input=file.sql;<br>echo "cat: $cats,path: $i,order: $order";<br>cats=`echo "select cat from work.yrse100sqkm_stream where start_long=(select end_long from work.yrse100sqkm_stream where cat=$cats) and start_lat=(select end_lat from work.yrse100sqkm_stream where cat=$cats) limit 1; " | db.select -c;`<br>
order=order+1;<br>done<br>done<br><br><br>_________________________________________________________________________________<br><br>Definitely, knowledge of python would have made things easier.<br><br>Thanks to all for the help.<br>
<br>