r.hydro.CASC2D Reference Manual

Copyright 1994, 1995, 1996, F.L. Ogden and B. Saghafian

All Rights Reserved

This reference manual may be copied without the permission of the authors.
This document may be included in part or wholly in any other document provided
that this document, F.L. Ogden, and B. Saghafian are referenced as the original
source.

F.L. OGDEN AND B. SAGHAFIAN MAKE NO WARRANTIES EITHER EXPRESS OR IMPLIED
REGARDING THE SUITABILITY OF THE PROGRAM r.hydro.CASC2D AND ITS FITNESS FOR ANY
PARTICULAR PURPOSE OR THE VALIDITY OF THE INFORMATION CONTAINED IN THIS
REFERENCE MANUAL.



NAME r.hydro.CASC2D - GRASS raster command to execute fully integrated
distributed hydrologic modeling.

SYNOPSIS: r.hydro.CASC2D, r.hydro.CASC2D help, r.hydro.CASC2D [-toepidbuq]
elevation=mapname time_step=value tot_time=value discharge=mapname
outlet_east&north&slope=east,north,bedslope rain_duration=value
[watershed_mask=mapname] [initial_depth=mapname] [storage_capacity=mapname]
[interception_coefficient=mapname] [roughness_map=mapname] [Manning_n=value]
[conductivity=mapname] [capillary=mapname] [porosity=mapname]
[moisture=mapname] [pore_index=mapname] [residual_sat=mapname]
[lake_map=mapname] [lake_elev=mapname] [radar_intensity_map=mapname]
[links_map=mapname] [nodes_map=mapname] [channel_input=mapname]
[table_input=mapname] [dis_profile=mapname] [wat_surf_profile=mapname]
[hyd_location=mapname] [r_gage_file=mapname] [unif_rain_int=value]
[num_of_raingages=value] [gage_time_step=value] [radar_time_step=value]
[write_time_step=value] [unit_el_conv=value] [unit_lake=value]
[unit_space=value] [d_thresh=value] [dis_hyd_location=mapname]
[depth_map=mapname] [inf_depth_map=mapname] [surf_moist_map=mapname]
[rate_of_infil_map=mapname] [dis_rain_map=mapname]

NOTE: In above, the command line arguments have been rearranged so that the
required parameters are placed first. When doing "r.hydro.CASC2D help", the
user will see a different sequence of parameters.

DISCLAIMER: The r.hydro.CASC2D model is Copyleft 1994 by B. Saghafian and F. L.
Ogden. There are absolutely no warranties, expressed or implied, made regarding
the accuracy of the results produced by this program. Use of this program is
allowed as long as the original developers are properly acknowledged.  No part
of this document describing the use of the model may be used or duplicated
without a complete citation.  Users with no background in or understanding of
distributed hydrology are strongly advised against using this code in any mode,
particularly in operational mode.  Besides knowledge of basic hydrology,
experience with typical numerical techniques used in physically-based
hydrodynamic models is recommended as it will help the user grasp capabilities
and limitations of this model.  This manual is significantly condensed for
electronic distribution and is  in no way comprehensive.  Users are encouraged
to experiment with the model and venture in hydrology textbooks and journal
papers to learn more about the topics touched upon in this manual.  The
r.hydro.CASC2D code is continuously being  improved.  Changes in the source
code of r.hydro.CASC2D may be made at any time, without notification.  No
claims are made regarding the suitability of r.hydro.CASC2D for any particular
purpose.  The model r.hydro.CASC2D was written for research and educational
purposes. Use r.hydro.CASC2D at your own risk.


HISTORY: CASC2D originally began with the two-dimensional overland flow routing
algorithm developed in APL by Prof. P.Y. Julien at Colorado State University.
The overland flow routing module was converted from APL to FORTRAN by Bahram
Saghafian, then at Colorado State University, with the addition of Green & Ampt
infiltration and explicit channel routing (Julien and Saghafian, 1991,
Saghafian, 1992, and Julien et al., 1995). The Fortran version was
reformulated, significantly enhanced, and re-written in the C programming
language by Bahram Saghafian at the U.S. Army Construction Engineering Research
Laboratories.  Implicit channel routing code was developed and added to
r.hydro.CASC2D by Fred L. Ogden (Ogden, 1994), formerly at Colorado State
University, now Asst. Professor, Department of Civil and Environmental
Engineering, University of Connecticut, Storrs, Connecticut.  This version
became known as r.hydro.CASC2D, part of the GRASS GIS for hydrologic
simulations (Saghafian, 1993).

DESCRIPTION: r.hydro.CASC2D is a physically-based, distributed, raster
hydrologic model which simulates the hydrologic response of a watershed subject
to a given rainfall field. Input rainfall is allowed to vary in space and time.
Major components of the model include interception, infiltration, and surface
runoff routing.  Interception is a process whereby rainfall is retained by
vegetation.  Interception is estimated using an empirical three parameter
model.  Infiltration is the process whereby rainfall or surface water is pulled
into the soil by capillary and gravity forces. The Green and Ampt equation with
four parameters is applied to model the event-based infiltration. For
continuous soil moisture accounting, redistribution of soil moisture can also
be simulated whenever the non-intercepted  rainfall intensity falls below the
saturated hydraulic conductivity of the soil.  The redistribution option
requires two more soil hydraulic parameters.  Excess rainfall becomes surface
runoff and is routed as overland flow and subsequently as channel flow. The
overland flow routing formulation is based on a two-dimensional explicit finite
difference (FD)  technique, while two different FD techniques, one explicit and
one implicit, provide options for routing one-dimensional channel flow. Through
a step function, a depression depth may be specified, below which no overland
flow will be routed.

The following sections describe various aspects of the model. In executing the
model, the likelihood of making a serious mistake by an inexperienced user is
unfortunately very high due to complexity of input and output options. Thus the
user must thoroughly read all the details of this manual and then select
appropriate options for his or her needs.  Since at this time GRASS is an
integer GIS, all input maps are stored in integer form.  This requires
conversion from integer to floating point representation via a scaling factor.
The scaling factor for all

floating point values is fixed.  Therefore, the integer maps of parameter
values are created by multiplying floating point values (e.g. volumetric water
content) by a multiplier.

PARAMETERS/OPTIONS: The following input/output parameters/options control
complexity of the simulation.  Map and file names in square brackets [] are
optional.  Some maps are mutually exclusive ( logical -or- ), while some maps
require other maps to enable proper function ( logical -and- ).   Carefully
read the NOTES section.

INPUT


           TOPOGRAPHY


elevation=mapname                map of elevation (DEM).

outlet_east&north&slope=value,   easting, northing, and bed slope at the outlet
value,value (comma delimited)

              TIME


time_step=value                  computational time step duration in seconds.
                                 (typ. 1 to 30 seconds)

tot_time=value                   total simulation time in sec.

         OVERLAND FLOW


[Manning_n=value]                spatially uniform Manning's n roughness value
                                 for overland flow.

    -or-

[roughness_map=mapname]          spatially varied map of Manning's n roughness
                                 coefficient (values in 1000*Manning's n).

[watershed_mask=mapname]         map of watershed boundary (or mask).  This
                                 option is recommended, as it speeds execution
                                 greatly.

[d_thresh=value]                 threshold overland depth, in meters, below
                                 which overland routing will not be performed
                                 (i.e. average depression storage).

[initial_depth=mapname]          map of initial overland (not lakes) depth in
                                 mm.

            RAINFALL


rain_duration=value              total rainfall duration in sec.

[unif_rain_int=value]            spatially uniform rainfall intensity in mm/hr.


    -or-

[r_gage_file=filename]           raingage rainfall input file name (ASCII).

    -and-

[num_of_raingages=value]         number of recording raingages.

    -and-

[gage_time_step=value]           time step (temporal resolution) of recorded
                                 raingage data in sec.

    -or-

[radar_intensity_map=mapname]    prefix of time series of maps of radar- (or
                                 otherwise-) generated rainfall intensities in
                                 mm/hr.

    -and-

[radar_time_step=value]          time increment between radar- (or otherwise-)
                                 generated rainfall maps in sec.

          INTERCEPTION


[storage_capacity=mapname]       map of vegetation storage capacity in tenths
                                 of mm.

    -and-

[interception_coefficient=mapnam map of interception coefficient (values in
e]                               1000*actual coefficient).

          INFILTRATION


[conductivity=mapname]           map of soil saturated hydraulic conductivity
                                 in tenths of mm/hr (Req'd for G&A and Redist).

[capillary=mapname]              map of soil capillary pressure head at the
                                 wetting front in tenths of mm (Req'd for G&A
                                 and Redist).

[porosity=mapname]               map of soil effective porosity (values in
                                 1000*porosity) (Req'd for G&A and Redist).

[moisture=mapname]               map of initial soil moisture (values in
                                 1000*moisture) (Req'd for G&A and Redist).

[pore_index=mapname]             map of soil pore-size distribution index
                                 (Brooks & Corey lambda) in 1000*index (Req'd
                                 for Redist).


[residual_sat=mapname]           map of soil residual saturation (values in
                                 1000*residual saturation) (Req'd for Redist).

             LAKES


[lake_map=mapname]               map of lakes categories.

[lake_elev=mapname]              map of lakes initial water surface elevation
                                 (also see unit_lake).

        CHANNEL ROUTING


[channel_input=filename]         channel input data file name (ASCII), required
                                 for explicit (EX) and implicit (IM) channel
                                 routing methods.

[links_map=mapname]              map of channel network link numbers. (EX & IM)

[nodes_map=mapname]              map of channel network node numbers. (EX & IM)

[table_input=filename]           look-up table file for links with breakpoint
                                 cross section, link type 8, (ASCII) (IM)

[dis_profile=filename]           channel initial discharge profile file name
                                 (ASCII). (IM)

[wat_surf_profile=filename]      channel initial water surface profile file
                                 name (ASCII). (IM)

[hyd_location=filename]          file name containing link and node numbers of
                                 internal locations where discharge hydrographs
                                 are to be saved (ASCII).

             UNITS


[unit_el_conv=value]             unit conversion factor by which the values in
                                 elevation must be DIVIDED to convert them into
                                 meters.

[unit_lake=value]                unit conversion factor by which the values in
                                 lake surface  elevation map must be DIVIDED to
                                 convert them into meters.

[unit_space=value]               unit conversion factor by which all region
                                 easting and northing values must be DIVIDED to
                                 convert them into meters.

OUTPUT


discharge=filename               outlet hydrograph file name (ASCII).


[dis_hyd_location=filename]      output file name for discharge hydrograph at
                                 internal locations (ASCII)

[write_time_step=value]          time increment for writing output raster maps
                                 in seconds.

    -and-

  [depth_map=mapname]            output maps of surface depth in mm.

  [inf_depth_map=mapname]        output maps of cumulative infiltration depth
                                 in tenths of mm.

  [surf_moist_map=mapname]       output maps of surface soil moisture in number
                                 of fractions of a thousand.

  [rate_of_infil_map=mapname]    output maps of infiltration rate in mm/hr.

  [dis_rain_map=mapname]         output maps of distributed rainfall intensity
                                 in mm.



FLAGS: There are several flags whose utility is driven by data availability
and/or user's choice.


-t     spatially interpolates raingage rainfall intensities using Thiessen
       polygon technique.  The default technique uses inverse-distance-squared
       proportionality for interpolation of rainfall intensity over space.

-o     routes edge-accumulated overland flow out of active region (ONLY when no
       mask is specified). Often the surface water accumulated at the edges of
       the current region creates severe backwater effects and may limit the
       use of longer computational time steps.

-e     performs one-dimensional explicit finite difference channel routing. May
       be suitable for low- to medium-intensity rainstorms over small arid and
       semi-arid watersheds with no base flow discharge.  This option often
       limits the computational time step to small values (<10 seconds).
       Channel bed smoothing is recommended to eliminate adverse slopes. No
       hydraulic structures, except reservoir spillways, can be simulated.

-p     assumes uniform channel geometry in each link (requires -e option). If
       this option is chosen, the channel input file (channel_input) must have
       only one line per fluvial link rather than the default of one line per
       node.

-i     performs Preissmann double sweep implicit channel routing. Particularly
       suitable for watersheds with some base flow to avoid dry-bed condition
       with channel slopes less than 1%.  Supercritical slopes cannot be
       handled; a warning message will be printed if supercritical flow is
       encountered.


-b     initializes the channel depth and discharge files (similar to -d,
       requires -i option) using the standard step backwater method. This
       option must be used with -i flag and replaces -d option to write
       "dis_profile" and "wat_surf_profile" files.  At present, only link types
       1, 2, and 8 may exist in the channel network.

-d     performs channel initialization for implicit routing (similar to -b,
       requires -i option) by flooding the entire channel network with a
       horizontal water surface and draining down to normal depth using a y(t)
       outlet boundary condition (similar to -b).  It is essential for implicit
       channel routing technique that a minimum initial base flow discharge
       exist in the channels and also the corresponding initial water surface
       profile at each node in the channel network have a realistic value.
       When the depth at the outlet reaches normal depth, the values of depth
       and discharge at each node is written to "dis_profile" and
       "wat_surf_profile" files for use in start up of actual simulations.
       With this option no other component of the model is executed; only the
       implicit channel routing is performed to create initial depth and
       discharge files necessary for start up of actual simulations.

-u     writes discharges in cubic feet per second (cfs) and volumes in cubic
       feet in "discharge" file.  The internal calculations are all in SI units
       regardless of this flag. The default SI unit prints the discharges in
       cubic meters per second (cms) and volumes in cubic meters (m3).

-q     quietly skips printing iteration, time, and discharge values to the
       screen.  No status report is printed.

IMPORTANT NOTES

The user must pay close attention to the following notes prior to any
simulation:

1) watershed_mask: Although optional, preparation of this map is highly
recommended as it cuts down on the memory requirements by the amount directly
proportional to the ratio of mask area over the region area  (also see -o
flag). If the basin (mask) delineation has not been performed correctly,
surface water may accumulate near the edges of the watershed. This excess water
has no way out of the watershed and will accumulate, creating undesirable
backwater effects within the watershed which eventually dictate use of shorter
time step in order to accommodate such effects. It is recommended, in such
instances, to re-delineate the watershed along the edges. An alternative would
be to lower the elevation of the cells being flooded by excessive surface
runoff near the edges such that no artificial backwater is created.


2) elevation: The elevation map in the form of Digital Elevation Model (DEM) is
undoubtedly one of the most important inputs for distributed modeling.  Thus
the quality of the DEM plays a major role in success of distributed hydrologic
simulations.  DEMs almost always contain errors.  Large flat areas in the DEM
may be due to the limited vertical resolution of elevation data. Routing over
such flat areas usually creates problems for the numerical techniques used in
distributed physically-based models.  Unreal pits in the DEM may be artifacts
of interpolation scheme used to rasterize digitized contours, or due to coarse
resolution in areas of concave topography.  As a rule of thumb, the user must
cross check the DEM with topographic maps of the area.  One way to discover
potential errors in the DEM is to run r.hydro.CASC2D iteratively with no other
option except uniform rainfall and a relatively short time step.  Writing
surface depth maps should also be selected (see depth_map and write_time_step
options).  The simulation may terminate normally, at which case the surface
depth maps must be examined in order to determine where most water accumulation
occurred and whether such accumulation areas is justified by the topographic
map of the area.  Also  stream network may be checked for proper connectivity
via depth maps.  Alternatively, if the model crashed at a certain location
(whose address is printed on the screen as well as at the bottom of discharge
file) the user must zoom in and double check the DEM  with the topographic map.
Often times some manual editing of the DEM is necessary in order to impose the
actual drainage trend of the topo map. Prior to editing it is recommended that
the DEM be multiplied by 10 or 100, particularly if the original vertical
resolution of the DEM is also a concern. To account for this multiplication
factor use unit_el_conv=10 or unit_el_conv=100, whichever appropriate. Note
that only the unreal pits and depressions must be removed since they most
likely trap surface runoff which would have otherwise contributed to the
outlet.  The real lakes and reservoirs may be simulated if delineated (see
lake_map).  In any event, non-smoothed DEMs require short computational time
steps, while properly smoothed DEMs, particularly those with coarser grid space
resolution, allow use of longer time steps.  Another source of concern, which
was briefly mentioned in above paragraph, about the quality of the DEMs is that
a nicely connected stream network cannot be derived from non-smoothed DEM.  If
no delineated network (in the form of a vector file, for example) is available,
an approach similar to what was described before may be taken in order to edit
the DEM. However where a network has been delineated independent of the DEM
then the elevation of stream cells should be checked so that they are not
higher than those of the surrounding cells. Otherwise the stream cells will not
properly collect the surface runoff.


3) initial_depth: This is a map which contains initial overland depth values,
if any.  Rarely used since prior to the storm overland planes are often dry.
For channel initial depth see dis_profile and wat_surf_profile options.

4) storage_capacity & interception_coefficient: In r.hydro.CASC2D, the
interception rate (i) is expressed as:

        i(t)=r(t)    while I < a
        i(t)=b*r(t)  while I > a

where r(t) denotes rainfall intensity at time t; a  is the
storage capacity; b  is the interception coefficient; and I  is
the cumulative interception depth. Storage capacity map, as well as
interception coefficient map, are usually a reclass of vegetation map.  For
table of storage capacity and interception coefficient values see Gray (1970,
section 4.6) or Bras (1990) p. 233.

5) roughness_map: This map represents the spatial distribution of overland
Manning roughness coefficient n.  This map could be a reclass of vegetation
cover map and/or the land use map.  Tables of Manning's n are available in
most hydrology textbooks.  By using Manning resistance equation it is assumed
that the overland flow over watersheds satisfies the conditions of turbulent
flow over rough surfaces.

6) conductivity, capillary, porosity, & moisture: Four soil property maps are
needed for modeling infiltration process using the Green-Ampt technique.  These
maps respectively contain saturated hydraulic conductivity, capillary suction
head, effective porosity, and initial soil moisture content.  The first three
maps may be reclasses of soil texture map and the last one must be prepared
considering antecedent soil moisture conditions.  Based on soil textural
classifications, tables of estimates of the first three parameters can be found
in Rawls et al. (1983).  Wherever reliable field measurements of such
parameters are available they may be substituted for table values.  Note that
the saturated hydraulic conductivity map must be adjusted for the percent of
rock fraction within the soil, if known.  If these four maps are specified,
Green & Ampt infiltration calculations are performed.  If one or more of the
maps is not specified, no infiltration calculations are made.

7) pore_index & residual_sat: These two maps are required when continuous soil
moisture accounting is of primary interest. The model is capable of
redistributing the soil moisture during periods of no- or low- intensity
rainfall, over which infiltration capacity may recover for the next burst of
storm intensity. The technique used herein for hiatus and post-hiatus stages is
primarily based on the method by Smith el al. (1993).  In this model Green-Ampt
equation is used for post-hiatus stage. Tables of pore-size distribution index
(Brooks & Corey lambda) and residual saturation are given by Rawls et al.
(1982).  As in

(6), measured values should be substituted where available.  When the four maps
specified in (6) as well as the two maps specified here are given as command
line arguments, the redistribution infiltration routine is used.  All six maps
must be specified to enable redistribution infiltration calculations.

8) lake_map & lake_elev: The first map is the delineated lake map with
different categories corresponding to individual lakes and the second map holds
the initial water surface elevation in the lakes.  It is best to number the
lakes  categories sequentially. If a lake or reservoir connects to the implicit
channel network at its outlet, then the lake must also be numbered as a link as
reflected in the channel_input file. See channel_input and links_map options
for more details.  However isolated pond areas may be simulated as well; in
fact such ponds are recommended to be delineated as part of the lake_map.
Such  isolated lakes should not be present in the links map since their outlet
is not  connected to the channel network.  The routing of the lakes is
performed by linear  reservoir technique.  Rainfall is added directly to the
lake inflow.  No  interception or infiltration is abstracted from the lake
cells.

9) radar_intensity_map: This is the common prefix of radar- or otherwise-
generated rainfall maps. One possible source of this type of data could be  the
NEXRAD system (Crum and Albery, 1993).  Alternatively, such time series  of
raster maps may come from the output of an interpolation software, which takes
in the rainfall site data corresponding to different times and generates the
time series of raster maps.  For instance, where raingage rainfall data is
available, one could run s.surf.tps of GRASS, one time for every recording
time, and save the interpolated rainfall maps to be used as  input for
r.hydro.CASC2D. In any case one must pay attention to the unit of  intensity
which is in integer millimeters per hour for map values. Also see
radar_time_step, and rain_duration options.  Example: If one sets
radar_intensity_map=rain.map and there are a total of eight maps in regular
intervals of 10 minutes (radar_time_step=600 and rain_duration= 4800, both in
seconds), then the actual name of the maps in the current mapset would have  to
be: "rain.map.1", "rain.map.2" through "rain.map.8", respectively
corresponding to time periods of 0-10, 10-20 through 70-80 minutes.

10) links_map & nodes_map: A link is a channel segment, of finite length, which
is comprised of two or more computational nodes placed at the center of each
grid cell. Any internal boundary condition (weirs and lakes/reservoirs, for
example) is also considered a link with only two nodes, one node upstream of
the  internal link and one downstream.  Although all links and nodes must be
numbered as discussed below, the link and node maps only contain trapezoidal or
look-up table cross sections (link types 1 and 8).  Internal boundary condition
link types (e.g. weirs) do not appear in the link or node maps.  The internal
boundary condition link types do appear in the channel_input file (see

attachment).  The purpose of the link and node maps is to provide connectivity
information between overland flow grid cells and channel (link,node) pairs.
This information is used to pass overland flow to the 1-D channel model as
lateral inflow.  The nodes map must contain the node numbers corresponding to
the links map.  Thus nodes corresponding to internal boundary condition links
are not present in the nodes map.  At a cell where a link terminates and
another link begins (junctions for example) the node number must comply with
the downstream link; i.e. this cell is  assigned a value of 1 in the nodes map.
However one must note that the last node of the upstream link also lies at such
junction cells.  This last node number must be accounted for in the main ascii
channel data file (channel_input), in which the total number of nodes per link
must be given.  An exception to this rule is the most downstream link leading
to the outlet. The links map essentially describes topology of the stream
network and its strict conformation with the stream numbering conventions is
vital.  Output from the GRASS command r.watershed cannot be directly used in
conjunction with r.hydro.CASC2D.  The link numbers assigned by r.watershed must
be re-numbered in accordance with the following rules.  Additionally, the node
map may be constructed by re-numbering the link map.  The general rules for
link and node numbering are:

A) The first link is numbered 1. B) Upstream links must have smaller link
numbers than downstream links. C) Link numbers must change at junctions. D)
Link numbers may not be skipped. E) Node numbers increase in the downstream
direction. F) Looped reaches cannot be simulated. G) Streams may run only in
x-y directions and not diagonally. H) Link types cannot be mixed within a link.

A small example network with three links is given below:

             LINK MAP                        NODE MAP

     0 0 0 0 0 0 0 0 0 0 0 0 0       0 0 0 0 0 0 0 0 0 0 0 0 0
     0 0 0 0 1 0 0 0 0 0 0 0 0       0 0 0 0 1 0 0 0 0 0 0 0 0
     0 0 0 1 1 0 0 0 0 0 2 0 0       0 0 0 3 2 0 0 0 0 0 1 0 0
     0 0 0 1 0 0 0 0 0 2 2 0 0       0 0 0 4 0 0 0 0 0 3 2 0 0
     0 0 0 1 1 0 0 0 2 2 0 0 0       0 0 0 5 6 0 0 0 5 4 0 0 0
     0 0 0 0 1 0 0 2 2 0 0 0 0       0 0 0 0 7 0 0 7 6 0 0 0 0
     0 0 0 0 1 1 3 2 0 0 0 0 0       0 0 0 0 8 9 1 8 0 0 0 0 0
     0 0 0 0 0 0 3 0 0 0 0 0 0       0 0 0 0 0 0 2 0 0 0 0 0 0
     0 0 0 0 0 0 3 3 0 0 0 0 0       0 0 0 0 0 0 3 4 0 0 0 0 0
     0 0 0 0 0 0 0 3 0 0 0 0 0       0 0 0 0 0 0 0 5 0 0 0 0 0

This network has three fluvial links, link 1 has 9 nodes in the node map, while
link 2 has 8 and link 3 has 5.  Note that the junction of links 1 and 2 is
labeled link 3.  However, links 1 and 2 have hidden nodes at this junction,
which must have the same bed elevation as the first node in link 3, but
different cross-sections.  The junction is really not the same point in space
for all three links, but really represents a short distance from the
confluence.  In the channel_input file, link 1 really

has 10 nodes and link 2 has 9.  Link 3 would only have 5 nodes, because its
downstream end is not a junction, rather it is the watershed outlet.  The link
and node maps tell r.hydro.CASC2D that the overland flow in the grid cell at
row 5, column 4 is passed to link 1, node 5.

For more discussion of links and nodes maps, see channel_input description and
the last section of this document titled "ATTACHMENT for channel_input  and
table_file File Options".

11) channel_input, table_file: For a detailed description of channel data  file
set-up and requirements: see the attachment titled "ATTACHMENT for
channel_input and table_file File Options", which is a modified version of the
document presented at the CASC2D Workshop by F.L. Ogden, June 9-10, 1994, at
the University of Memphis, Tennessee (Revision 2, 10 January 1995).  Also it is
recommended for the user to refer to Ogden (1994) for an overview of the
implicit channel routing formulation and numerical scheme.  Note that for
explicit channel routing, the same channel input file may be used provided that
only link types 1 and 4 are present.  If explicit channel routing is used, the
channel talweg elevations are taken from the channel_input file, not from the
DEM.

12) dis_profile & wat_surf_profile: These two files hold, respectively, the
discharge profile and the water surface profile of all nodes within the stream
network, including internal boundary condition nodes.  These files are created
by running r.hydro.CASC2D with the -i option AND either -d or -b.  After
creation, these files are used to initialize the implicit channel routing
scheme for actual watershed simulations.  They must be present for actual
simulations.  The -d flag creates initial discharge and depth files by flooding
the entire channel network with a horizontal water surface, and draining the
network using a depth vs. time relation at the watershed outlet.  This draining
proceeds until normal depth is reached at the outlet.  At this time, the
dis_profile and wat_surf_profile files are written for use in channel
initialization for later simulations.  The -b flag creates the initialization
files using the standard step backwater method.  Note that in dis_profile the
unit of discharge is in cubic meters per second and the unit of depth in
wat_surf_profile is in meters. Except for the reservoirs, the water surface
profile file actually holds the flow depth values and thus is measured relative
to the bed.  CAUTION:  The implicit channel routing routine must have good
first approximation of depth and discharge at each node in the network.
Inaccurate values (guesses) will surely lead to a crash.  Additionally, the
implicit channel routine is for SUBCRITICAL (Froude number < 1) flow only.  A
newly created channel input file will almost always crash when first run.  It
is then the task of the user to identify the location of the problem.  Usually,
there are regions of hydraulically steep slopes which cause supercritical flow.
The channel code will warn the user of supercritical flow, including

node and link number when it occurs, and then exits.  The user should look at
that node/link combination, and alter the data file to eliminate the reach of
steep slopes.  One workable solution to prevent stability problems at low flows
is to use the so-called "Preissmann Slot" (see Cunge et al., 1980, for
example).  Be patient, as getting the channel network running is a time
consuming process.  There are a host of possible errors including: abrupt
changes in cross-section, DEM-error, etc., which can cause problems.  The
implicit routing method is not applicable to steep (mountainous or upland)
streams.  If supercritical flow is encountered in upland (1st order) streams,
they can be eliminated from the network.  The user should verify the
suitability of this approach on each watershed.

13) hyd_location: If a filename is specified, this input ascii file contains
the link and node numbers of the locations at which discharge hydrographs are
to be written to the file named in the dis_hyd_location option. The first
column in this file should contain link numbers and the second column must be
filled with the corresponding node numbers.  For a description of the output,
see the dis_hyd_location option.

14) r_gage_file: This is an ascii file which must be provided when raingage
rainfall data is being simulated.  At the top of the file, two-column lines
hold the easting and the northing for each raingage.  The number of such lines
is determined by the total number of recording raingages (num_of_raingages).
The location of any of the gages does not have to be within the current region
nor within the current watershed mask as long as the easting and the northing
are not specified relative to the current region, but are based on absolute
values in the UTM or SP coordinates.  If a gage falls outside in a different
zone to the  left of the active region's zone, then negative values are also
acceptable.  Note however that raingages well outside the watershed under
analysis generally  provide poor rainfall estimates.  The subsequent lines at
the bottom of the  raingage file must reflect temporal variation of rainfall
intensity.  The number of columns per line is equal to the number of raingages
(specified via num_of_raingages).  The columns are separated by space.  The
number of lines in this lower portion is equal to number of instances,
separated by a constant time interval, the raingages have made a recording.  As
usually is the case,  the unit of rainfall intensity for raingage data must be
in inches per hour.  Example: For three raingages, each recording rainfall
intensity every 2 minutes for a total duration of, say, 10 minutes, a file
called "rain.inp" may look  like this:


205150.0    750212.0
20545.0     750104.0
205320.0    750173.0

0.0     0.0     0.55
1.75    2.25    0.80
1.00    1.80    1.50
0.65    0.90    0.70
0.0     0.50    0.30

In above example, the eastings and northing of the first, second, and third
raingages are (205150.0,750212.0), (20545.0,750104.0), and (205320.0,750173.0)
respectively.  The intensities recorded by the first gage, for example, are
0.0,  1.75, 1.00, 0.65, and 0.0 inches per hour, respectively, over 0-2, 2-4,
4-6, 6-8,  and 8-10 minutes, etc. For this example r_gage_file=rain.inp,
num_of_raingages=3,  gage_time_step=120, and rain_duration=600. For state plane
(SP) coordinate system  the eastings and northing will have to be in feet
rather than meters.

15) outlet_east&north&slope: These three values determine the location of the
outlet, in terms of its easting and northing, and the outlet bed slope.  One
needs to make sure that the outlet described by its easting and northing is not
only within the active region but also inside watershed mask.  Often times the
region is not set to its original settings after zoom-in operations (d.zoom)
are performed and this may put the outlet outside the active region thus
causing the model to eventually crash.  The bed slope is equal to tangent of
the angle which is made between the bed profile at the outlet and horizontal
plane.  This slope is primarily used to calculate the outflow overland
discharge at the outlet, if any, based on normal depth boundary condition, when
no channel routing is performed (all surface flow treated as overland flow and
channels essentially assumed wide).  Normal depth is also assumed to prevail at
the outlet when explicit channel routing (-e) has been selected.

16) Manning_n: The alternative to simulating spatially varied roughness
coefficient is to provide a single value of Manning roughness coefficient n via
Manning_n parameter.  The user is warned if both roughness_map and Manning_n
have been specified.

17) unif_rain_int: This option represents the intensity of the spatially-
uniform temporally-constant (up to the rainfall duration) rainfall.  Therefore
this option may replace r_gage_file and radar_intensity_map options.  The unit
is in millimeters per hour.

18) num_of_raingages: The total number of recording raingage.  If this variable
is set to one, no spatial interpolation is performed and the rainfall is
treated as uniform is space but could vary in time.  The temporal variation
will be then provided by r_gage_file.  Note that when num_of_raingages is set
to one, t

he easting and northing of the gage is irrelevant although two (arbitrary)
values must still be provided at the top of raingage file (r_gage_file).

19) time_step: This represents the duration of computational time step in
seconds and is a critical variable determining the total execution time for  a
particular simulation.  These is no firm guide for the selection of the time
step; it comes with experience and strongly depends on watershed and rainfall
characteristics.  The general rule for overland routing and explicit channel
routing is that shorter time steps must be used for higher intensity storms,
finer horizontal grid resolution (grid spacing), steeper watershed slopes,
larger watershed areas, and smoother surfaces.  Stability of explicit routing
depends upon Courant number.  Unfortunately, the critical condition for Courant
number limits the length of the computational time step which must be used for
the entire simulation unless variable time step algorithm is implemented in
the future.  Shorter time steps must be used when backwater effects are
generated, mainly in flat areas which are not part of the lake map.  If the
time step is too long for any particular simulation the surface water depth  in
completely flat areas may show a checker-board pattern, i.e. oscillations   are
observed in the water surface level.  This eventually results in a crash. As
such, the time step should be decreased and the simulation repeated; or the
flat areas be delineated within the lake map.

20) gage_time_step: This variable represents the time interval, in seconds,
between recording instances of rainfall intensities by raingage(s).  It is
implied that the rainfall data has been recorded in regular intervals.  See
the  notes under r_gage_file for an example.

21) radar_time_step: A value expressing the time interval between consecutive
rainfall maps.  A uniform time step is implied.

22) rain_duration: This is the total duration of rainfall in seconds.  If
multistorm events are simulated, the time from the beginning of the first storm
to the end of the last storm constitutes the total rainfall duration.  As such,
selection of soil moisture redistribution capability is recommended via
specifying pore index (pore_index) and residual saturation (residual_sat) maps.

23) tot_time: The total simulation time in seconds.  If the falling limb of the
discharge hydrograph is of particular interest the total simulation time must
be set to a value greater than total rainfall duration plus the expected
recession time.

24) write_time_step: This time parameter determines the frequency of writing
output raster maps and is equal to the time interval, in seconds, at which
output raster maps are saved.  Also see depth_map, inf_depth_map,
surf_moist_map, rate_of_infil_map, and dis_rain_map output options.


25) unit_el_conv: This conversion factor is used to convert elevation values in
DEM (elevation) map to meters.  This parameter doesn't need to be specified if
the DEM is already in meters since the default is one.  For DEM units in cm or
ft, unit_el_conv must be respectively set to 100 or 3.281.

26) unit_lake: This conversion factor is used to convert the water surface
elevation values in the lake_elev map to meters.  The default is 1.0.

27) unit_space: This conversion factor is used to convert the horizontal grid
spacing resolution to meters.  The default is 1.0.  For state plane coordinate
system in ft, set unit_space equal to 3.281.

28) d_thresh: This parameter represents an average step retention storage
below which the overland depth will not be routed.  Another words, all
depressions less than or equal to d_thresh are filled before any overland flow
could begin.  The unit must be in meters so for 2 mm of depression storage set
d_thresh=0.002.  Higher values of depression storage would reduce the total
execution time of the model as the overland routing consumes most of the CPU
effort.

29) discharge: The discharge hydrograph computed at the outlet will be saved in
this ascii file under the current directory.  Some other information, such as
peak discharge, is also printed in this file.

30) dis_hyd_location: Whenever the hyd_location option is selected, the
discharge at individual node/link pairs will be saved in this file.  The
discharge hydrographs at the (link,node) locations specified in hyd_location
file are grouped in columns.  The number of lines in this file is determined by
the total number of iterations equal to tot_time divided by time_step.

31) depth_map, inf_depth_map, surf_moist_map, rate_of_infil_map, &
dis_rain_map: depth_map is the common filename prefix given to the time series
of output raster maps containing surface depth values in millimeters. At the
channel cells the channel flow depth will be recorded in the maps. The surface
depth maps, and all other output maps, are saved at regular intervals
determined by write_time_step option. If write_time_step is not set, no output
raster map will be saved. The first map always corresponds to the initial
condition and  naturally shows the water surface profile corresponding to the
base flow discharge within the channel network when implicit channel routing is
performed.  Similarly the last depth map corresponds to the end-of-simulation
time, or to the time at which the program finished abnormally, for example due
to selection of a long step which generated oscillations leading to a crash.
Abnormal program termination caused by oscillating depths may show negative
depths in the overland plane.  There is a hard-coded limit of 2000 output
raster maps for each simulation.  inf_depth_map and rate_of_infil_map are two
options to save the output raster maps

of cumulative infiltration depth, in tenth of mm, and rate of infiltration in
mm/hr, respectively.  These two options can only be selected when infiltration
is being computed via either the Green & Ampt or redistribution methods.
surf_moist_map output contains the soil moisture values at the soil surface and
it may solely be selected with continuous infiltration option (pore_index and
residual_sat maps are specified).  dis_rain_map option saves the instantaneous
rainfall intensities, in mm/hr,  at write_time_step intervals and may be set
for the simulations involving spatially distributed rainstorms; i.e. raingage
rainfall data or radar (or other sources of rainfall raster maps) data.
Example: If one sets depth_map=depth.out, tot_time=1000, and write_time_step=
100, then a total of 11 raster depth maps will be saved in a
normally-terminated simulation: depth.out.00, depth.out.01, depth.out.02
through  depth.out.10, respectively corresponding to times equal to 0, 100,
200,  through 1000 seconds.  These maps may be viewed sequentially in animation
form using the GRASS animator program xganim.

REFERENCES:

Bras, R. L., 1990, Hydrology: An introduction to hydrologic science,
     Addison-Wesley, Reading, Mass., 643 p.

Crum, T. D., and R. L. Alberty, 1993, The WSR-88D and the WSR-88D operational
     support facility, Bulletin of the American Meteorological Soc., 74(9),
     pp.  1669-1687.

Cunge, J.A., F.M. Holly, and A. Verwey, 1980, Practical Aspects of
     Computational River Hydraulics, Iowa Institute of Hydraulic Research,
     404HL The University of Iowa, Iowa City, IA 52242. 420 p.

Gray, D.M., 1970, Handbook on the Principles of Hydrology, National Research
     Council of Canada, Water Information Center Inc., Water Research Building,
     Manhasset Isle, Port Washington, N.Y., 11050.

Julien, P. Y., and B. Saghafian, 1991, CASC2D users manual - A two dimensional
     watershed rainfall-runoff model, Civil Engr. Report, CER90-91PYJ-BS-12,
     Colorado State University, Fort Collins, CO.

Julien, P. Y., Saghafian, B., and F. L. Ogden, 1995, "Raster-Based Hydrologic
     Modeling of Spatially-Varied Surface Runoff", Water Resources Bulletin,
     AWRA, 31(3), pp. 523-536.

Ogden, F.L., 1994, de-St Venant channel routing in distributed hydrologic
     modeling., Proc. Hydraulic Engineering `94, ASCE Hydraulics Specialty
     Conference, G.V. Cotroneo and R.R. Rumer, eds., Vol. 1, pp. 492-496.

Rawls, W. J., Brakensiek, D. L., and N. Miller, 1983, Green-Ampt infiltration
     parameters from soils data, J. of Hydraulic Engineering, ASCE, 109(1), pp.
     62-70.


     Rawls, W. J., Brakensiek, D. L., and K. E. Saxton, 1982, Estimation of
     soil water properties, Trans. of ASAE, pp. 1316-1320.

Saghafian, B., 1992, Hydrologic analysis of watershed response to spatially
     varied infiltration, Ph.D. Dissertation, Civil Engr. Dept., Colorado State
     University, Fort Collins, CO.

Saghafian, B., 1993, Implementation of a distributed hydrologic model within
     Geographic Resources Analysis Support System (GRASS), Proceedings of the
     Second International Conference on Integrating Environmental Models and
     GIS,  Breckenridge, CO.

Smith, R. E., Corradini, C., and F. Melone, 1993, Modeling infiltration for
     multistorm runoff events, Water Resources Research, 29(1), pp. 133-144.



           ATTACHMENT for channel_input and table_file File Options

GENERAL NOTES

     The naming convention associated with the Preissmann double-sweep 1-D
implicit channel routing method is based on the concept to links and nodes. A
link is a channel segment, or an internal boundary condition, which is
comprised of two or more computational nodes. All internal boundary conditions
contain two nodes, while fluvial reaches may be of any size greater than or
equal to two nodes.  The following discussion of input file format must first
distinguish between different link types.  A few of the possible link types are
presented below in Table 1.  As of August 1995, only link types 1,  2, 4, and 8
are supported.  Development is continuing on link types 3 and 7.


                              Table 1.  Link Types



  Link                                                Number of      Number of
  Type            Description of Link Type            Parameters       Nodes

    1    Fluvial Link, Trapezoidal Cross-Section          5             >=2

    2    Overflow Weir                                    8              2

    3    Culvert/Weir (not yet supported)                 8              2

    4    Reservoir                                        5              2

    7    Bridge Crossing (not yet supported)              8              2

    8    Fluvial Link, Look-up Table (Breakpoint)         4             >=2
         Cross-Section


     At present the channel model formulation accepts cross-sectional input for
only two different channel geometries, namely trapezoidal, and breakpoint via
look-up table.  Trapezoidal channel parameters which must be provided as input
at each computation node include: Manning's n, bottom width, channel  depth,
side slope, and bed elevation.  Look-up tables of cross sectional properties
must include cross-sectional area, top width, and conveyance at  equal depth
intervals.  Smooth transition in channel cross sectional properties within
links of type 8 and between all connecting fluvial links often plays a vital
role in the success of simulations.  Abrupt changes in cross-sections can lead
to numerical errors in mass conservation.


     As far as handling the reservoirs (link type 4), flow is not routed
through reservoirs in this version of the CASC2D channel routing code.  Instead
the linear reservoir approximation is used.  Among other internal boundary
conditions link types, only weirs can be simulated at present.

     The channel_input file contains the channel bed elevation at each node,
which constructs the longitudinal profile of channels in the drainage network.
Ideally, the modeler should use surveyed cross sections and talweg profiles of
the channel network.  However, extensive surveys are often impossible for the
entire drainage networks on large watersheds.  In lieu of an extensive survey,
there are existing tools for extracting channel topology from digital elevation
models (DEM).  However, if the channel network is extracted from a DEM, a
smoothing algorithm must be applied (e.g. Ogden et al., 1994) to produce
physically realistic longitudinal profiles because of errors inherent in any
DEM.


EXAMPLE OF INPUT FILE FORMAT AND CONSTRUCTION

NOTE: The channel_input file and table_file both follow SI unit convention.
All units of length are in meters, including elevations.

     In this example, a channel input data file is constructed for a fictitious
watershed.  Note that this particular example is merely for demonstration and
has not been tried in any simulation run.  The stream network will consist of
trapezoidal and look-up table fluvial links, an internal reservoir, and an
overflow weir.  Regarding the MASK map (watershed_mask option), note that all
watershed cells must be marked with 1, while all raster cells outside the
watershed boundary are marked with 0. Also all the cells in the lake should
carry the value 1 on the LAKE map (lake_map  option).  If there was another
lake in this example, it would be denoted with  the number 2, etc.

     The channel link numbering scheme typically employed in  double-sweep
routing was explained in the main text.  With reference to the LINK map shown
below, (links_map option), assume that links 1  and 2 drain into link 3, which
in-turn drains in to the lake, which is assigned link number 4.  Also links 7
and 8 drain into link 9.  Link 5 originates at the outlet of the reservoir.
Link 5 is split into link 6 because of a change in some  cross-section
property.  Links 6 and 9 flow into link 10, which is immediately  upstream from
a weir which must be considered link 11.  Finally, link 12, the  most
downstream link, lies below the weir (link 11).  Note that the LINK map
contains continuous sequence  of link numbers, except for the internal
boundary conditions such as the weir (link 11) and lakes (link 4) which does
not appear in the LINK map.



                      EXAMPLE LINK MAP

0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  1  0  0  0  0  0  0  0  0  0  0  2  0  0  0  0  0
0  0  0  1  0  0  0  0  0  0  0  0  0  2  2  0  0  0  0  0
0  0  0  1  1  0  0  0  0  0  0  0  2  2  0  0  0  0  0  0
0  0  0  0  1  1  1  0  0  0  0  0  2  0  0  0  0  0  0  0
0  0  0  0  0  0  1  1  1  0  0  0  2  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  1  0  2  2  2  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  1  1  3  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  3  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  3  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  7  7  7  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  7  7  0  0  0  0  5  0  0  0  0  0  0  0  0  0
0  0  0  0  0  7  7  0  0  0  5  5  0  0  0  0  0  0  0  0
0  8  8  8  0  0  7  0  0  0  0  5  6  0  0  0  0  0  0  0
0  0  0  8  8  8  9  9  0  0  0  0  6  6  0  0  0  0  0  0
0  0  0  0  0  0  0  9  9  9  0  0  0  6  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  9  9  9  9 10 10 10  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 12  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 12 12 12 12  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 12  0


     Now, refer to the NODE map (nodes_map option) shown below.  In the NODE
map, each link is numbered from 1 to the number of grid cells spanned by that
link, with the exception of the internal boundary conditions.  Conceptually,
internal boundary conditions (including reservoirs) have two nodes but they are
not given any node numbers in the NODE map.  This is how the program recognizes
internal boundary conditions.  Furthermore, all links except the one leading to
the watershed outlet must have an extra node to provide connectivity to the
downstream link.  Therefore, even though the NODE map may show link 1 to have
13 nodes, it actually has 14.  This implied extra node for link 1, shared
between the most downstream node of link 1 and the most upstream node of link
3, provides the connection between link 1 and link 3. The number of nodes in
each link in this example are shown below in Table 2.  Note that the node map
entries for the lake (link 4) are 0.

                      EXAMPLE NODE MAP

0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  1  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0
0  0  0  2  0  0  0  0  0  0  0  0  0  3  2  0  0  0  0  0
0  0  0  3  4  0  0  0  0  0  0  0  5  4  0  0  0  0  0  0
0  0  0  0  5  6  7  0  0  0  0  0  6  0  0  0  0  0  0  0
0  0  0  0  0  0  8  9 10  0  0  0  7  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0 11  0 10  9  8  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0 12 13  1  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  2  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  3  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  1  2  3  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0
0  0  0  0  4  5  0  0  0  0  1  0  0  0  0  0  0  0  0  0
0  0  0  0  0  6  7  0  0  0  2  3  0  0  0  0  0  0  0  0
0  1  2  3  0  0  8  0  0  0  0  4  1  0  0  0  0  0  0  0
0  0  0  4  5  6  1  2  0  0  0  0  2  3  0  0  0  0  0  0
0  0  0  0  0  0  0  3  4  5  0  0  0  4  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  6  7  8  9  1  2  3  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  1  0  0  0  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  2  3  4  5  0
0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  6  0



   Table 2. Number of Nodes in each Link for Example Watershed Stream Network



         Link Number          Number of Nodes             Number of Nodes
                              as in nodes_map          in channel-input file

              1                      13                          14

              2                      10                          11

              3                      3                           4

              4                      0                           2
         (reservoir)

              5                      4                           5

              6                      4                           5

              7                      8                           9

              8                      6                           7

              9                      9                           10

              10                     3                           4

              11                     0                           2
            (weir)

              12                     6                           6
        (outlet link)


     The first portion of the channel input file is used to pass physical
constants and simulation  parameters  to  the  model. These  include  the
gravitational acceleration "g", kinetic energy  correction factor "alpha", the
friction slope weighting factor "beta", the spatial derivative weighting
coefficient "theta", the length  of each node "dx", the computational time step
"dt" (seconds), the total simulation time "tt" (seconds), and the discharge in
all first order  streams "qmin" (m3/s).  At present, the time step "dt" must be
identical to the  computational time step used in the overland flow routing
portion of  r.hydro.CASC2D.  The program must be told the number of links, and
the largest  number of nodes in any link in the network for dynamic memory
allocation.  In our example problem, the number  of links is 12, and the
maximum number of nodes is 14 (in link 1).  Remember that all links which are
not at the outlet or not immediately upstream of a reservoir must have an extra
node for connectivity purposes.  The total number of links is called  "nlinks",
and the largest number of nodes in any link in the network is called
"maxnodes".  In this example, we will use the constants and parameters given
in Table 3:



             Table 3. physical constants and simulation  parameters


                     "g"               9.81 m/s2

                     "alpha"           1.0

                     "beta"            0.5

                     "theta"           0.55

                     "dx"              100.0 m

                     "dt"              30.0 s

                     "tt"              3600.0 s

                     "qmin"            0.07 cms

                     "nlinks"          12

                     "maxnodes"        14


     This data constitutes the first portion of channel_input file and is
arranged into a header which must have the form (note floating point and
integers):



9.81
1.0  0.5  0.55
100.0
30.0  3600.0
0.07
12
14


     The second portion of the input file describes link types as well as  the
network topology and connectivity.  This is accomplished using a  line-input
format, one line for each link in the network.  Each line contains 6 values
arranged in columns, as shown in Table 4.



                   Table 4. Connectivity information format


    Column     Data Description

       1       Link number (INT)

       2       Link type (INT)

       3       Number of links upstream from this link (max. 2, min. 0) (INT)

       4       Upstream link #1 (INT) (0 if no upstream links)

       5       Upstream link #2 (INT) (0 if no or only one upstream link)

       6       Downstream link (0 for outlet ) (INT)


     As far as link types in the current example, assume that links 1 through
9, except reservoir link 4 (link type 4), are well-described as trapezoidal
(type 1), and we are  using look-up table data to describe the cross-sections
of links 10 and 12  (type 8).  The weir (link 11) is of link type 2.  Thus the
second portion of input file for describing link types and  connectivity looks
like:

1   1   0   0   0   3
2   1   0   0   0   3
3   1   2   1   2   4
4   4   1   3   0   5
5   1   1   4   0   6
6   1   1   5   0  10
7   1   0   0   0   9
8   1   0   0   0   9
9   1   2   7   8  10
10  8   2   6   9  11
11  2   1  10   0  12
12  8   1  11   0   0

Note that only the outlet link has 0 as a downstream dependency.  It is
important that this be the only link with 0 as a downstream dependency.  Also,
the location and topology of internal boundary condition links (4 and 11) with
respect to other links are known via above connectivity information in the
channel_input file.  As an interpretation, examine link number 10 above.   It
is of link type 8 (look-up table cross-section data), has 2 upstream
dependencies (links 6 and 9), and flows into link 11.

     The third portion of the  input  file contains the  individual  hydraulic
property information for each node in the network.  Assume for now that all the
trapezoidal channels in the network have the properties shown in Table 5.



      Table 5. Cross-Sectional Properties of Example Trapezoidal Channels


     Manning roughness coefficient   varies

     Bottom width (m)                varies

     Channel depth (m)               1.75

     Side slope H:V                  varies

     Talweg elevation                dependent upon location


Now, build the input file section for link #1 which has 13 nodes on the NODE
map, plus an extra for connectivity at the upstream end of link 3  (total of 14
nodes).  This input file section will look like:

1  14
0.035  2.10  1.75  2.00  264.40
0.035  2.10  1.75  2.00  263.85
0.035  2.10  1.75  2.00  263.41
0.035  2.10  1.75  2.00  262.98
0.035  2.10  1.75  2.00  262.75
0.035  2.10  1.75  2.00  262.54
0.035  2.10  1.75  2.00  262.37
0.035  2.10  1.75  2.05  262.02
0.035  2.20  1.75  2.05  261.87
0.035  2.20  1.75  2.10  261.75
0.035  2.20  1.75  2.10  261.67
0.035  2.20  1.75  2.15  261.52
0.035  2.30  1.75  2.15  261.25
0.035  2.40  1.75  2.15  261.00

The 1 and 14 on the first line indicate that it is link 1, with 14 nodes.  The
columns represent Manning's n, bottom width, channel depth, trapezoidal side
slope H/V, and talweg elevation, respectively.

Link 2 might have an entry which looks like:

2  11
0.035  1.80  1.75  1.80  264.21
0.035  1.80  1.75  1.80  264.10
0.035  1.90  1.75  1.90  263.81
0.035  1.90  1.75  1.90  263.56
0.035  2.00  1.75  1.90  263.34
0.035  2.10  1.75  1.90  262.86
0.035  2.10  1.75  2.00  262.24
0.035  2.10  1.75  2.00  261.76
0.035  2.20  1.75  2.10  261.48
0.035  2.30  1.75  2.10  261.26
0.035  2.30  1.75  2.10  261.00


The cross-section entry for link 3 would then look like:

3  4
0.035  2.60  1.75  2.20  261.00
0.035  2.60  1.75  2.40  260.20
0.035  2.80  1.75  2.90  259.67
0.035  2.80  1.75  2.90  259.17

     The elevation of the downstream end of links 1 and 2 must be equal to the
bed elevation of the upstream end of link 3. It is required that the bed
elevation of all channel inverts at each junction be equal.  Even though they
are not exactly the same points in space, they are assumed close enough to have
the same elevation.


                    Flow                 Flow
                       \                    /
                    \   \    \         /   /     /
                     \   v    \       /   v     /
                      \        \     /         /
                       \ link 1 \   / link 2  /
                        \ node 14\ / node 11 /
                         \                  /
                          \                /
                           \              /
                            \            /
                             |          |
                             |  link 3  |
                             |  node 1  |
                             |          |


     The input for other links with trapezoidal cross section is similar to
the input for link 1 (see the next section for the complete input file).

     Additionally, assume that the reservoir (link 4) has a surface area of
0.498 square km, a rectangular spillway width of 12m, a spillway discharge
coefficient of 0.97, an initial water surface elevation of  260.10m, and a
spillway crest elevation of 260.50m.  The resulting input file entry is shown
below.

4  2
0.498  12.0  0.97  260.10  260.50
0.000   0.0  0.00    0.00    0.00

The number of nodes (line entries) in the channel_input file per reservoir link
is two, but the numbers entered in the second row entry are for further
improvements and for now their value is irrelevant.  Also note that the bed
elevation at the downstream limit of link 3 MUST be lower than the initial
water surface elevation in the reservoir (link 4). This is mandatory at all
downstream boundaries between channels flowing into the reservoirs.  Reservoirs
serve as downstream boundary conditions a

t upstream links.  The elevation of reservoirs should therefore never be
allowed to fall below critical depth in any channels which flow into the
reservoir.  If this happens, the code will crash.

     The input for links 5, 6, 7, 8, and 9 will be similar to links 1, 2, and
3.  See the constructed input file at the end of this text.  Trapezoidal
channel links are simple in terms of input.  The channel depth field is
intended to represent the average bank-full depth of the channel.  This number
is not used in the present version of the code, but will be in the future when
the overland flood plane and channel flows are coupled iteratively.  At present
the flows are not fully coupled.  Water can flow only from the overland flow
plane into the channels, not from the channels back to the flood plane.

     If look-up table data are used (link type 8), the look-up tables are
stored in a separate file.  The GRASS command line option for this look-up
table file is table_file.  If there are any link type 8 in the channel_input
data file, the program will attempt to open the file table_file (typically
called "table.dat").  If this file is not present, the program will terminate.
In channel_input file, all that is needed is the table number for a particular
cross section, and the talweg elevation.  Assume that the 4 nodes in link 10
are approximated by one cross-section, which is given as table entry 1, then
the channel_input file entry for link 10 would look like:

10  4
1  244.0 1  243.5 1  243.0 1  242.5

In this format, the first column represents the look-up table number for the
given node, and the second column represents the bed elevation for that node.
Table numbers can be mixed within a link.  However, abrupt changes in
cross-section should be avoided because they cause significant continuity
errors in the formulation.  The format of the table_file is discussed later in
this document.

     Link 11 is an overflow weir link.  While this link has no nodes which
appear in the node map for overland flow connectivity, it does have two nodes
in terms of channel topology.   The meaning of the columns in the weir node
input section is presented in Table 6.



                   Table 6.  Parameter Description for Weirs


   Line Entry      Column    Parameter Description

        1             1      forward direction weir discharge coefficient
        1             2      reverse direction weir discharge coefficient
        1             3      0.0 (reserved for future use)
        1             4      crest length, meters
        1             5      elevation of weir crest, meters
        2            1-4     0.0 (reserved for future use)
        2             5      downstream bed elevation, meters



The first (upstream) node represents the crest of the rectangular weir.  The
second (downstream) node represents the bed of the channel just downstream from
the weir.  Assume that this weir has a crest elevation of 244.4 m, a crest
width of 8 m, a discharge coefficient in the forward direction of 0.92, a
discharge coefficient for flow in the reverse direction equal to 0.85, and the
channel bed elevation immediately downstream from the weir is 239.8m.  The
entry for this weir in the channel_input file would thus appear as:

11  2
0.92   0.85  0.00  8.00  244.4
0.00   0.00  0.00  0.00  239.8

     Link 12 uses table entry 2 for the its first four nodes and entry 4 for
the remaining two nodes.  Note that link 12 has 6 nodes in the node map, and 6
nodes in the channel_input file.  Because link 12 is the outlet link, we do not
add an extra node at the downstream end for connectivity purposes.  The portion
of the channel_input file for link 12 looks like:

12  6
2  239.80
2  239.47
2  239.33
2  239.05
3  238.54
3  238.44

Hence, 238.44 is the talweg elevation at the outlet of the catchment.  Also
notice that the upstream end of link 12 and the downstream node of link 11 have
the same bed elevation.

      Weirs cannot be used as the outlet boundary condition for this channel
model.  If your watershed has a weir at the outlet, just place a short link
downstream from the weir with appropriate characteristics.



FINAL CHANNEL FILE FORMAT

     The example channel data file (typically called chn.dat) developed for
this example looks like:



----------------------BEGIN channel_input FILE HERE----------------------------
9.81
1.0  0.5  0.55
100.0
30.0  3600.0
0.07
12
14

1   1   0   0   0   3
2   1   0   0   0   3
3   1   2   1   2   4
4   4   1   3   0   5
5   1   1   4   0   6
6   1   1   5   0  10
7   1   0   0   0   9
8   1   0   0   0   9
9   1   2   7   8  10
10  8   2   6   9  11
11  2   1  10   0  12
12  8   1  11   0   0

1  14
0.035  2.10  1.75  2.00  264.40
0.035  2.10  1.75  2.00  263.85
0.035  2.10  1.75  2.00  263.41
0.035  2.10  1.75  2.00  262.98
0.035  2.10  1.75  2.00  262.75
0.035  2.10  1.75  2.00  262.54
0.035  2.10  1.75  2.00  262.37
0.035  2.10  1.75  2.05  262.02
0.035  2.20  1.75  2.05  261.87
0.035  2.20  1.75  2.10  261.75
0.035  2.20  1.75  2.10  261.67
0.035  2.20  1.75  2.15  261.52
0.035  2.30  1.75  2.15  261.25
0.035  2.40  1.75  2.15  261.00

2  11
0.035  1.80  1.75  1.80  264.21
0.035  1.80  1.75  1.80  264.10
0.035  1.90  1.75  1.90  263.81
0.035  1.90  1.75  1.90  263.56
0.035  2.00  1.75  1.90  263.34
0.035  2.10  1.75  1.90  262.86
0.035  2.10  1.75  2.00  262.24
0.035  2.10  1.75  2.00  261.76

0.035  2.20  1.75  2.10  261.48
0.035  2.30  1.75  2.10  261.26
0.035  2.30  1.75  2.10  261.00

3  4
0.035  2.60  1.75  2.20  261.00
0.035  2.60  1.75  2.40  260.20
0.035  2.80  1.75  2.90  259.67
0.035  2.80  1.75  2.90  259.17

4  2
0.498  12.0  0.97  260.10  260.50
0.000   0.0  0.00    0.00    0.00

5  5
0.035  2.50  1.75  2.10  249.20
0.035  2.50  1.75  2.10  248.65
0.035  2.50  1.75  2.15  247.94
0.035  2.50  1.75  2.20  247.31
0.035  2.50  1.75  2.20  246.97
0.035  2.50  1.75  2.25  246.42

6  5
0.035  2.50  1.75  2.00  246.42
0.035  2.50  1.75  2.00  245.97
0.035  2.50  1.75  2.00  245.24
0.035  2.50  1.75  2.00  244.73
0.035  2.50  1.75  2.00  244.00

7  9
0.035  1.75  1.75  2.00  254.20
0.035  1.75  1.75  2.00  253.55
0.035  1.75  1.75  2.00  252.92
0.035  1.90  1.75  2.00  252.29
0.035  1.90  1.75  2.00  251.71
0.035  1.90  1.75  2.00  251.24
0.035  2.00  1.75  2.00  250.82
0.035  2.00  1.75  2.00  250.34
0.035  2.00  1.75  2.00  250.08

8  7
0.035  1.60  1.75  2.00  253.71
0.035  1.90  1.75  2.00  253.21
0.035  1.90  1.75  2.00  252.94
0.035  2.00  1.75  2.00  252.55
0.035  2.10  1.75  2.00  251.80
0.035  2.10  1.75  2.00  251.11
0.035  2.10  1.75  2.00  250.08

9  10
0.035  3.10  1.75  2.90  250.08
0.035  3.20  1.75  2.90  249.16
0.035  3.30  1.75  3.00  248.63
0.035  3.40  1.75  3.00  248.02
0.035  3.40  1.75  3.10  247.33

0.035  3.20  1.75  3.10  246.65
0.035  3.10  1.75  3.10  246.05
0.035  3.20  1.75  3.10  245.88
0.035  3.20  1.75  3.10  245.52
0.035  3.50  1.75  3.10  245.24
0.035  3.70  1.75  3.40  244.74
0.035  3.50  1.75  3.40  244.51
0.035  3.10  1.75  3.40  244.00

10  4
1  244.00
1  243.50
1  243.00
1  242.50

11  2
0.92   0.85  0.00  8.00  244.40
0.00   0.00  0.00  0.00  239.80

12  6
2  239.80
2  239.47
2  239.33
2  239.05
3  238.54
3  238.44
--------------------- END OF channel_input FILE ---------------------------



TABLE_FILE INPUT FORMAT

     If link type 8 (look-up table x-sections) are present in the channel_input
file, we require a file typically called "table.dat", to store all look-up
table values.  This file must begin with an integer equal to the number of
tables contained within the file.  In our example this number is 3.  Then we
require 3 tables.  The first line of each table is an integer equal to the
table number. The second line in the table entry contains two numbers, an
integer, and a floating point (real) value.  The first number is equal to the
number of entries (rows) in each particular table, designated by "numhts". The
second number is the vertical distance between hydraulic property points, which
must be a constant for each table.  For instance, if you describe the variation
of cross-sectional area, topwidth, and conveyance at 0.5 m vertical intervals,
this number will be 0.5.

     Each table entry must then contain "numhts" lines, each with four
entries.  The first line must be:

0.0   0.0   0.0   0.0

subsequent lines must contain the following:


height   topwidth   x-section_area  x-section_conveyance.

     For instance, assume we know geometric variables for table entries 1, 2,
and 3.  The file "table.dat" (table_file option) may look like:

3
1 5  1.0
0.0   0.0    0.0    0.0
1.0   2.1    5.9   13.5
2.0   6.9   22.3  122.2
3.0  14.7   72.1  312.8
4.0  19.6  145.6  789.4

2 8  0.5
0.0   0.0    0.0    0.0
0.5   1.1    1.5    7.5
1.0   1.7    3.4   17.1
1.5   3.1   11.2   43.9
2.0   5.2   32.4   98.5
2.5   7.2   49.2  187.4
3.0   9.2   86.4  312.5
3.5  12.1  143.1  624.9

3 5  1.0
0.0   0.0    0.0    0.0
1.0   2.2    5.4   15.5
2.0   6.6   15.3  132.2
3.0  15.7   81.8  352.8
4.0  21.4  155.2  839.6


For further clarification, consider table 2 above, corresponding to cross
sectional properties of the four most upstream nodes of link 12.  Table 2 has 8
entries, each separated by 0.5 meter of depth.  At a depth of 3.0 m, for
example, the channel has a top-width of 9.2m, a cross-sectional area of 86.4
m2, and a conveyance of 312.5 m3/s.  The conveyance K is used to calculate the
discharge using:

                    Q=K*sqrt(Slope)

Where using Manning's equation, K is defined (in SI units) as:

       K = 1/n * (Area) * (Hydraulic_radius)^(2/3)

These tables can be produced easily using a spreadsheet using measured channel
cross-section data.



REFERENCES

Cunge, J.A., F.M. Holly, and A. Verwey, 1980, Practical Aspects of
     Computational River Hydraulics, Iowa Institute of Hydraulic Research,
     404HL The University of Iowa, Iowa City, IA 52242. 420 p.

Ogden, F.L., 1994, de-St Venant channel routing in distributed hydrologic
     modeling., Proc. Hydraulic Engineering `94, ASCE Hydraulics Specialty
     Conference, G.V. Cotroneo and R.R. Rumer eds., August  1-5,  1994,
     Buffalo, N.Y., pp. 492-496.

Ogden, F.L., Saghafian, B., and W.F.  Krajewski, 1994, GIS-based  channel
     extraction and smoothing  algorithm for  distributed hydrologic  modeling,
     Proc. Hydraulic Engineering `94,  ASCE Hydraulics Specialty  Conference,
     G.V. Cotroneo and R.R. Rumer eds., August  1-5,  1994, Buffalo,  N.Y.,
     pp. 237-241.
