stream order

[Jochen Schmidt]

In article <{BCT-PGOTFE}@bnr.ca> you wrote:
: I am looking for a means of automating the classification of stream 
: segments by "order". An "Order 1" stream is one that has no other 
: stream flowing into it. An "Order 2" stream is one that is formed by 
: the junction of two "Order 1" streams, etc. Can anyone suggest a 
: GRASS routine that could be used to solve this problem?

What you want to do is to classify a stream net in a Strahler stream order,
I suppose. To get this classification, you can use the GRASS module
r.watershed armsed=... the output armsed file is a text file, which
determines the stream segments, that flow in each stream segments (The
stream segments are the output of the r.watershed streams=... option). 
I have written a shellscript, which calculates out of this text file the
strahler order of the segments, perhaps you can use it, but no warranties:

#!/bin/sh
if [ ! "$GISRC" ]
then 
    echo "Sorry, you are not running GRASS"
    exit 1
fi


# 
# Shellscript, was aus dem output von r.watershed
# eine Ordnung nach Strahler erzeugt
#

myname=`basename $0`

usage()
{
  echo
  echo "USAGE:  ${myname}  input=name armsed=name output=name"
  echo "        input    - input rastermap with streams or basins"
  echo "        armsed   - output-file of r.watershed"
  echo "        output   - output rastermap with strahler-ordering "
  echo 
  exit 1
}

# Bezeichnung der temporaeren Files

temp1=`g.tempfile pid=$$`
temp2=`g.tempfile pid=$$`
temp3=`g.tempfile pid=$$`

trap 'rm -f $temp1 $temp2 $temp3'

# Eingabe ueberpruefen

if [ $# -eq 0 ]
then
    g.ask type=old prompt="inputmap:" element=cell desc="raster" unixfile=/tmp/
$$
    . /tmp/$$
    rm -f /tmp/$$
    if [ ! "$name" ]
       then exit 0
    fi
    MAP=$name
 
   echo
   echo "output rastermap"
   echo "Hit RETURN to cancel request"
   echo
   OUT=`line`
   if  [ ! "$OUT" ]
     then exit 0
   fi

   echo
   echo "input armsed-file"
   echo "Hit RETURN to cancel request"
   echo
   IN=`line`
   if  [ ! "$IN" ]
     then exit 0
   fi

else

  # Ueberpruefung der maps
  if [ $# -ne 3 ]
  then 
      usage
  fi

  # test for correct form of arguments
  if [ `echo $1 | cut -c1-6` != "input=" -o `echo $2 | cut -c1-7` != "armsed=" 
-o `echo $3 | cut -c1-7` != "output=" ] 
  then
      usage
  fi

  # put args into named variables 
  MAP=`echo $1 | cut -c7-`
  IN=`echo $2 | cut -c8-`
  OUT=`echo $3 | cut -c8-`

  eval `g.findfile element=cell file="$MAP"`
     if [ ! "$file" ]
     then 
        echo $MAP": illegal mapname"
        exit 0
     fi
fi

# Erzeugung der Hilfsdatei
# Format: EZG-nr Nummern der Gebiete, die in EZG entwaessern

rm -f $temp1

for i in `awk '{print $1}' $IN`
do 
awk 'BEGIN{printf("%d ",'$i')} {if($4=='$i') printf("%d ",$1)} END{printf("\n")
}' $IN>> $temp1
done

# Extraktion der Ordnungen aus der Hilfsdatei

rm -f $temp3

while [ `wc -l $temp1|awk '{print $1}'` != 0 ]
do
nawk '
 NF==1{print $1" =  1">>"'$temp3'"} 
 NF>=4{print "Es sind uneindeutige Verzweigungen vorhanden" >"/dev/tty"; exit 6
}
 NF>1&&NF<4{s=$1
            i1=$2
            i2=$3
            while(getline <"'$temp3'") {if($1==i1) {o1=$3; break} else o1=0}
            close("'$temp3'")
            while(getline <"'$temp3'") {if($1==i2) {o2=$3; break} else o2=0}
            close("'$temp3'")
            if(o1==0 || o2==0) {print s,i1,i2; next} 
            if(o2==o1) {print s" = "o1+1>>"'$temp3'"}
            if(o2>o1) {print s" = "o2>>"'$temp3'"}
            if(o1>o2) {print s" = "o1>>"'$temp3'"}
            }' $temp1 >$temp2
if [ `echo $?` = "6" ]
   then exit 0
fi

mv $temp2 $temp1

done

# Erzeugung der output Map

r.reclass i=$MAP o=$OUT<$temp3

exit 0

Have fun with that hack, Jochen