Appendix A C-codes
Appendix A1 Program newnx.c--reconstructing the flow
connectivity between modeling units after some of units have been
removed
/* newnx.c
* from new1.c, 9/21/95
* Pawel Mizgalewicz, CRWR UT at Austin */
/****************************************************************
* newnx.c -- reconstructs the flow connectivity between units
* after some of them have been removed. Requires two ASCII,comma
* delimited files:
* file 1) full set of units. Each line should contain:
* unit_id, next unit_id
* file 2) reduced set of units. Each line should contain the
* unit_id number
* Output file: unit_id, id of downstream unit, unit order
*
* USAGE: newnx infile1 infile2 outfile
*
****************************************************************/
#include <stdio.h>
main (argc, argv)
int argc;
char *argv[];
/* arg 1) input_file (oid, onx => old ID, old NEXT)
2) input_file (nid => new ID)
3) output_file (nid, nnx, nor => ID, new NEXT, order)
*/
{
char loop, found;
int i, j, k, onn, nnn, nni, nor[1200], nix[1200];
long id, nx, ix,
oid[1200], onx[1200], nid[1200], nnx[1200];
FILE *fgsin, *funin, *ftable2;
if (argc < 3 )
{
printf ("wrong number of arguments = %d \n",argc);
exit(1);
}
fgsin = fopen (argv[1], "r");
i = 0;
while (fscanf(fgsin, "%li,%li", &id,&nx) != EOF)
{
oid[i] = id;
onx[i] = nx;
++i;
}
onn = i-1;
fclose(fgsin);
funin = fopen (argv[2], "r");
i = 0;
while (fscanf(funin, "%li", &id) != EOF)
{
nid[i] = id;
nnx[i] = 0;
++i;
}
nnn = i-1;
nni = i;
fclose(funin);
/* To each new unit assign old_next */
for (i = 0; i <= nnn; ++i)
for (k = 0; k <= onn; ++k )
if ( nid[i] == oid[k] )
{
nnx[i] = onx[k];
break;
}
/* check if all next are valid (have new_id specified)*/
for (i = 0; i <= nnn; ++i)
{
found = 0;
ix = 0;
while (!found)
{
++ix;
for (j = 0; j <= nnn; ++j)
if ( ( nnx[i] == nid[j] ) || ( nnx[i] == 0 ) )
{
found = 1;
break;
}
if (found) continue;
/* if not found, go to old table and read next_id */
for (j = 0; j <= onn; ++j)
if ( nnx[i] == oid[j] )
{
nnx[i] = onx[j];
break;
}
if(ix > 1000000) printf("error %i\n", ix);
}
}
/* =============== WATERSHED ORDER ================= */
/* Set initial value of order array */
for (i = 0; i <= nnn; ++i)
{
nor[i] = 1;
nix[i] = nni;
}
/* Select first order items */
for (i = 0; i <= nnn; ++i)
{
j = 0; loop = 1;
while ( loop )
{
if (nnx[i] == nid[j])
{
nor[j] = 0;
nix[i] = j;
loop = 0;
}
++j;
if (j == nni) loop = 0;
}
}
/* order of the remaining streams */
for (i=0; i<=nnn; ++i)
if(nor[i] == 1)
{
j = nix[i];
k = 2;
while ((nor[j] < k) && (j != nni))
{
nor[j] = k;
j = nix[j];
k = k + 1;
}
}
/* Write results to output file */
ftable2 = fopen (argv[3], "w");
for (j = 0; j <= nnn; ++j )
fprintf(ftable2, "%li,%li,%i\n", nid[j], nnx[j], nor[j] );
fclose(ftable2);
return 0;
}
Appendix A2 Program fdy4.c--calculating the flow
rate in ungauged streams from the available record
/************************************************************/
/* fdy4.c -- calculates flow rate in all modeling units from
* the available record. The average precipitation
* depth is used as a weight.
* arguments:
* 1) file name that contains USGS flow record
* 2) file name that contains modeling unit data
* 3) name of the output file.
* input: two ASCII files, each record contains the following
* values, coma delimited:
* file one specified by the first argument:
* gsid = USGS station identification number,
* gsnx = ID number of the downstream USGS station,
* gsqo = flow rate;
* file two specified by the second argument:
* unid = modeling unit identification number,
* unnx = ID number of the downstream modeling unit,
* ungsid = ID number of the USGS station that
* determines the zone in which the unit
* is located (ID of the gauging station that
* is downstream of the modeling unit),
* unar = area of the modeling unit,
* unor = modeling unit order in the flow system,
* unpr = average precipitation depth;
* output: ASCII file specified by the third argument:
* unid = modeling unit identification number,
* unqc = flow rate estimated in modeling unit
* (this flow is cumulative, i.e., it
* is a runoff from the drainage area
* determined by the modeling unit outlet. */
/************************************************************/
#include <stdio.h>
#define size 512
main (argc, argv)
int argc;
char *argv[];
/* arg 1) input_file (gs) must have: id, next, and flow.
2) input (modeling units--un)
must have: id, next, gsid, area, order, precip.
3) output_file: unid, qc (cumulative)
used in ArcView version:
4) name_id (copied from first column)
5) name_out */
{
char loop, notfound;
int i, j, k, l, n, unormx, d;
int unor[2000], unixx[2000], ungsix[2000], gsused[100], or,
gsni, gsnn, unni, unnn, unormxi;
long gsid[100], gsnx[100], id, nx, idgs,
unid[1200], unnx[1200], ungsid[1200];
float gsqo[100], gsqi[100], gsvl[100], gscf[100],
unar[1200], unqo[1200], unpr[1200], unqc[1200], unqct[1200],
ar, qo, pr, rncf, sum1, sum2, xxx, x2;
/* gsor[100], gsar[100], gspr[100] deleted, gsvl[100] created */
char buf[size];
FILE *fgsin, *funin, *ftable2;
if (argc < 3 )
{
printf ("wrong number of arguments = %d \n", argc);
exit(1);
}
fgsin = fopen (argv[1], "r");
/* fgets(buf, size, fgsin); */
i = 0;
while (fscanf(fgsin, "%li,%li,%f", &id,&nx,&qo) != EOF)
{
gsid[i] = id;
gsnx[i] = nx;
gsqo[i] = qo;
++i;
}
gsnn = i-1;
gsni = i;
fclose(fgsin);
/* look for records that has Q = 0. All records that has Q=0 will
have instead of next-id, the id they should have as a new unit.
*/
for ( i = 0; i <= gsnn; ++i )
{
if ( gsqo[i] > 0.0 ) continue;
for ( j = 0; j <= gsnn; ++j )
if ( gsid[i] == gsnx[j]) gsnx[j] = gsnx[i];
}
/* =========================================================*/
funin = fopen (argv[2], "r");
fgets(buf, size, funin); /* tables unload first, dummy record */
i = 0; /* that have negative area !!!! */
while (fscanf(funin, "%li,%li,%li,%f,%i,%f",
&id,&nx,&idgs,&ar,&or,&pr) != EOF)
{
unid[i] = id;
unnx[i] = nx;
ungsid[i] = idgs;
unar[i] = ar;
unor[i] = or;
unpr[i] = pr;
++i;
}
unnn = i-1;
unni = i;
fclose(funin);
/* update ungsid[i] */
for ( i = 0; i <= gsnn; ++i )
{
if ( gsqo[i] > 0.0 ) continue;
for ( j = 0; j <= unnn; ++j )
if ( gsid[i] == ungsid[j]) ungsid[j] = gsnx[i];
}
/* rebuild GS arrays, i.e. remove records GSQO <= 0 */
i = 0;
for ( j = 0; j <= gsnn; ++j )
{
if ( gsqo[j] <= 0.0 ) continue;
gsid[i] = gsid[j];
gsnx[i] = gsnx[j];
gsqo[i] = gsqo[j];
gsvl[i] = 0.0;
gsqi[i] = 0.0; /* initialization of the inflow table */
gsused[i] = 0;
++i;
}
gsnn = i-1;
gsni = i;
/* Set UNGSIX index, i.e. index that relates UNID with the GSID
= record in UN table with the record in GS table,
if ungsid[] = 0, i.e.,. no gauging station assigned to this record
index = gsnn + 1 = gsnx */
for (i = 0; i <= unnn; ++i)
{
if ( ungsid[i] == 0 )
{
ungsix[i] = gsni;
continue;
}
for (k = 0; k <= gsnn; ++k )
{
if ( ungsid[i] == gsid[k] )
ungsix[i] = k;
}
}
for (i = 0; i <= gsnn; ++i )
{
for (j = 0; j <= unnn; ++j )
if ( ungsid[j] == gsid[i] )
{
gsvl[i] = gsvl[i] + ( unpr[j] * unar[j] );
}
}
/* Calculate runoff coefficient, use first order GS watersheds only */
sum1 = 0.0;
sum2 = 0.0;
xxx = 1.0 / 8640.0 ; /* area in km2, prec in 0.001 cm/d, -> m3/s */
notfound = 1;
for (i = 0; i <= gsnn; ++i)
{
notfound = 1;
for (j = 0; j <= gsnn; ++j )
{
if (gsnx[j] == gsid[i] )
{
notfound = 0;
break;
}
}
if ( notfound == 1 )
{
++n;
sum1 = sum1 + gsqo[i];
sum2 = sum2 + ( gsvl[i] * xxx);
}
}
/* printf("s1 s2 %f,%f\n",sum1,sum2); */
rncf = sum1 / sum2;
/* 1) Calculate runoff from all UN watersheds */
/* 2) assign initial values to both, unqc[] and unqct[] */
/* 3) fill array unixx[], index of the next/downstream unit */
/* is set in "find maximum order" module */
for (i = 0; i <= unnn; ++i)
{
unqo[i] = rncf * unpr[i] * unar[i] * xxx ;
unqc[i] = unqo[i];
unqct[i] = unqo[i];
for ( j = 0; j <= unnn; ++j)
{
if (unid[j] == unnx[i])
{
unixx[i] = j;
break;
}
}
}
/* Find maximum UN order */
unormx = 0;
for (i = 0; i <= unnn; ++i)
{
if ( unormx < unor[i] )
{
unormx = unor[i];
unormxi = i;
}
}
unixx[unormxi] = unni;
/* Calculate zonal cumulative UN flow (zonal = within GS zone) */
for (k = 1; k < unormx; ++k)
{
for (i = 0; i <= unnn; ++i)
{
if ( unor[i] != k )
continue;
j = unixx[i];
if ( ungsid[j] != ungsid[i] )
continue;
unqc[j] = unqc[j] + unqc[i];
}
}
for (i = 0; i <= unnn; ++i)
{
unqct[i] = unqc[i];
}
/* Calculate sum of GS inflow (inflows ? */
for (i = 0; i <= gsnn; ++i)
{
for (j = 0; j <= gsnn; ++j)
{
if ( gsnx[i] == gsid[j] )
{
gsqi[j] = gsqi[j] + gsqo[i];
}
}
}
/* 1) Estimate correction factors GSCF */
/* 2) Assign correction factor for the UN watersheds which are
/* outside GS zones, i.e., most downstream UN watersheds
/* (next wsh for the last unit has GSID = 0) */
for (i = 0; i <= gsnn; ++i)
{
x2 = rncf * gsvl[i] * xxx;
gscf[i] = ( gsqo[i] / ( gsqi[i] + x2 ) - 1.0 ) /x2 ;
if (gsnx[i] == 0 )
{
gscf[gsni] = gscf[i];
gsid[gsni] = 0;
}
}
/* Calculate cumulative flow (total= observed inflow + calculated
cumulative flow */
for (i = 0; i <= unnn; ++i )
{
k = i;
j = unixx[k];
l = ungsix[k];
if ( ungsid[k] == ungsid[j] ) continue;
if ( gsused[l] == 1 ) continue;
xxx = gsqo[l];
gsused[l] = 1;
while (j < unni )
{
unqct[j] = unqct[j] + xxx;
k = j;
j = unixx[k];
if ( ungsid[j] != ungsid[k] ) break;
}
}
/* Final distribution of flow */
for ( i = 0; i <= unnn; ++i)
{
j = ungsix[i];
unqc[i] = unqct[i] * (1. + ( gscf[j] *unqc[i] ) );
}
/* Write results to output file */
ftable2 = fopen (argv[3], "w");
/* arc view version:
fprintf ( ftable2, "\"%s\",\"%s\"\n", argv[4], argv[5]);
*/
for (j = 0; j <= unnn; ++j )
fprintf(ftable2, "%li,%f\n", unid[j], unqc[j] );
fclose(ftable2);
return 0;
}