Skip to content

Commit

Permalink
Merge pull request #1 from jhamman/fix/uh_convolution
Browse files Browse the repository at this point in the history 
Fix/uh_convolution
  • Loading branch information
Joe Hamman authored Dec 23, 2016
2 parents 8ad59be + f80f9eb commit ee66188
Show file tree
Hide file tree
Showing 7 changed files with 75 additions and 111 deletions.
6 changes: 2 additions & 4 deletions src/Makefile
View file
Original file line number Diff line number Diff line change
@@ -1,14 +1,12 @@
########################################################################
### rout.f makefile ####################################################
### rout.f makefile ####################################################
########################################################################
#
# Routing algorithm written D. Lohmann
#
# This is a slightly modified code (main algotrithms unchanged -IO and
# array dimensions simplified).
# Maintained by G. O'Donnell (tempgd@hydro.washington.edu) and Andy Wood
#
# $Id: Makefile,v 1.1 2005/04/07 05:07:28 vicadmin Exp $
#

#This program uses the non-standard Fortran argument GETARG
Expand All @@ -22,7 +20,7 @@ FFLAGS = -O -C -ffixed-line-length-none
#for debugging
#FFLAGS = -C -g -lm -ffixed-line-length-none

FC=g77
FC=gfortran

HFILES= parameter.h

Expand Down
20 changes: 8 additions & 12 deletions src/init_routines.f
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,10 +10,6 @@ SUBROUTINE INIT_ARRAY( A, NROW, NCOL, VALUE )

IMPLICIT NONE

c RCS ID STRING
CHARACTER*50 RCSID
DATA RCSID/"$Id: init_routines.f,v 1.1 2005/04/07 05:07:28 vicadmin Exp $"/

INTEGER NCOL, NROW
INTEGER I, J
REAL A(NCOL,NROW)
Expand Down Expand Up @@ -65,7 +61,7 @@ SUBROUTINE CREATE_VIC_NAMES( JLOC, ILOC, EXTEN, CLEN, DPREC )
END DO

CLEN=CLEN-1

RETURN
END

Expand Down Expand Up @@ -95,27 +91,27 @@ SUBROUTINE SEARCH_CATCHMENT
II = I
JJ = J
300 CONTINUE
IF ((II .GT. ICOL) .OR. (II .LT.1) .OR.
IF ((II .GT. ICOL) .OR. (II .LT.1) .OR.
& (JJ .GT. IROW) .OR. (JJ .LT.1)) THEN
GOTO 310
END IF
IF ((II .EQ. PI) .AND. (JJ .EQ. PJ)) THEN
IF ((II .EQ. PI) .AND. (JJ .EQ. PJ)) THEN
NO_OF_BOX = NO_OF_BOX + 1
CATCHIJ(NO_OF_BOX,1) = I
CATCHIJ(NO_OF_BOX,2) = J
GOTO 310
ELSE
IF ((DIREC(II,JJ,1).NE.0) .AND. !check if the current
ELSE
IF ((DIREC(II,JJ,1).NE.0) .AND. !check if the current
& (DIREC(II,JJ,2) .NE.0)) THEN !ii,jj cell routes down
III = DIREC(II,JJ,1) !to the subbasin outlet
JJJ = DIREC(II,JJ,2) !point, following the
JJJ = DIREC(II,JJ,2) !point, following the
II = III !direction of direc(,)
JJ = JJJ !from each cell
GOTO 300
END IF !if you get there,
END IF !if you get there,
END IF !no_of_box increments
310 CONTINUE !and you try another
END DO !cell.
END DO !cell.
END DO

WRITE(*,*) 'Number of grid cells upstream of present station',
Expand Down
22 changes: 9 additions & 13 deletions src/make_convolution.f
View file
Original file line number Diff line number Diff line change
@@ -1,22 +1,18 @@

SUBROUTINE MAKE_CONVOLUTION
& (NCOL, NROW, NOB, PMAX, DAYS, CATCHIJ,
& (NCOL, NROW, NOB, PMAX, DAYS, CATCHIJ,
& BASE, RUNO, FLOW, KE, UH_DAY, UH_S, FRACTION, FACTOR_SUM,
& XC, YC, SIZE, DPREC, INPATH,ICOL,NDAY,IDAY,IMONTH,IYEAR,
& MO, YR, NYR)

IMPLICIT NONE

c RCS ID STRING
CHARACTER*50 RCSID
DATA RCSID/"$Id: make_convolution.f,v 1.1 2005/04/07 05:07:28 vicadmin Exp $"/

INTEGER N, I, J, DAYS, NDAY, II, JJ
INTEGER NCOL,NROW,ICOL,NOB,PMAX,KE,UH_DAY
INTEGER CATCHIJ(PMAX,2)
INTEGER NYR
REAL UH_S(PMAX,KE+UH_DAY-1)
REAL BASE(DAYS), RUNO(DAYS), FLOW(DAYS)
REAL BASE(DAYS), RUNO(DAYS), FLOW(DAYS)
REAL FRACTION(NCOL,NROW)

REAL PI, RERD, FACTOR, FACTOR_SUM
Expand All @@ -36,17 +32,17 @@ SUBROUTINE MAKE_CONVOLUTION

REAL STORAGE, K_CONST
REAL DUM1,DUM2

INTEGER IDAY(DAYS), IMONTH(DAYS), IYEAR(DAYS)
INTEGER MO(12*NYR),YR(12*NYR)
INTEGER MISS_NUM
C MISS_NUM is the number of grid cell output files not found

MISS_NUM=0

C *** 0 <= K_CONST = 1.0
C *** 0 <= K_CONST = 1.0
C *** K_CONST smaller 1.0 makes it a simple linear storage

K_CONST = 1.0

PI = ATAN(1.0) * 4.0
Expand All @@ -66,18 +62,18 @@ SUBROUTINE MAKE_CONVOLUTION
END DO
II = CATCHIJ(N,1)
JJ = CATCHIJ(N,2)

c the grid has been flipped left to right
c find the revised cooordinates

ILOC=XC + (ICOL-II)*SIZE + SIZE/2.0
JLOC=YC + JJ*SIZE - SIZE/2.0

AREA = RERD**2*ABS(SIZE)*PI/180* !give area of box in
AREA = RERD**2*ABS(SIZE)*PI/180* !give area of box in
& ABS(SIN((JLOC-SIZE/2.0)*PI/180)- !square meters
$ SIN((JLOC+SIZE/2.0)*PI/180))


AREA_SUM = AREA_SUM + AREA

c WRITE(*,*) N, ILOC, JLOC
Expand All @@ -87,7 +83,7 @@ SUBROUTINE MAKE_CONVOLUTION

FACTOR_SUM = FACTOR_SUM + FACTOR


call create_vic_names(jloc,iloc,loc,clen,dprec)

c print*, INPATH(1:INDEX(INPATH,' ')-1)//LOC(1:CLEN)
Expand Down
37 changes: 16 additions & 21 deletions src/read_routines.f
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
c SUBROUTINES RELATED TO READING
c SUBROUTINES RELATED TO READING
c read_diff()
c read_fraction()
c read_grid_uh()
Expand All @@ -10,10 +10,6 @@
SUBROUTINE READ_DIFF(DIFF,NCOL,NROW,FILENAME,
$ IROW, ICOL)

c RCS ID STRING
CHARACTER*50 RCSID
DATA RCSID/"$Id: read_routines.f,v 1.1 2005/04/07 05:07:28 vicadmin Exp $"/

INTEGER NCOL,NROW,IROW,ICOL,I,J
REAL DIFF(NCOL,NROW)
CHARACTER*72 FILENAME
Expand All @@ -26,8 +22,8 @@ SUBROUTINE READ_DIFF(DIFF,NCOL,NROW,FILENAME,
END DO

DO J = IROW,1,-1
READ(10,*) (DIFF(I,J), I=ICOL,1,-1)
END DO
READ(10,*) (DIFF(I,J), I=ICOL,1,-1)
END DO

CLOSE(10)

Expand All @@ -54,8 +50,8 @@ SUBROUTINE READ_FRACTION(FRACTION,NCOL,NROW,FILENAME,
END DO

DO J = IROW,1,-1
READ(22,*) (FRACTION(I,J), I=ICOL,1,-1)
END DO
READ(22,*) (FRACTION(I,J), I=ICOL,1,-1)
END DO

CLOSE(22)

Expand All @@ -69,7 +65,7 @@ SUBROUTINE READ_FRACTION(FRACTION,NCOL,NROW,FILENAME,
SUBROUTINE READ_GRID_UH
& (UH_BOX, KE, PMAX, NOB, CATCHIJ,FILENAME)


IMPLICIT NONE

INTEGER KE, PMAX, NOB
Expand Down Expand Up @@ -105,8 +101,8 @@ SUBROUTINE READ_VELO(VELO,NCOL,NROW,FILENAME,
END DO

DO J = IROW,1,-1
READ(10,*) (VELO(I,J), I=ICOL,1,-1)
END DO
READ(10,*) (VELO(I,J), I=ICOL,1,-1)
END DO

CLOSE(10)

Expand All @@ -133,8 +129,8 @@ SUBROUTINE READ_XMASK(XMASK,NCOL,NROW,FILENAME,
END DO

DO J = IROW,1,-1
READ(10,*, END=20) (XMASK(I,J), I=ICOL,1,-1)
END DO
READ(10,*, END=20) (XMASK(I,J), I=ICOL,1,-1)
END DO
CLOSE(10)

RETURN
Expand All @@ -153,11 +149,11 @@ SUBROUTINE READ_DIREC(DIREC,NCOL,NROW,H,XC,YC,SIZE
IMPLICIT NONE

INTEGER NCOL,NROW,I,J,IROW,ICOL,IMISS
INTEGER DIREC(NCOL,NROW,2)
INTEGER DIREC(NCOL,NROW,2)
INTEGER H(NCOL,NROW)
REAL XC, YC, SIZE
CHARACTER*72 FILENAME
CHARACTER*14 CDUM
CHARACTER*14 CDUM

OPEN(10, FILE = FILENAME, FORM = 'FORMATTED',
$ STATUS='OLD',ERR=9001)
Expand All @@ -174,10 +170,10 @@ SUBROUTINE READ_DIREC(DIREC,NCOL,NROW,H,XC,YC,SIZE
$ irow, icol
STOP
ENDIF

DO J = IROW,1,-1
READ(10,*) (H(I,J), I=ICOL,1,-1)
END DO
READ(10,*) (H(I,J), I=ICOL,1,-1)
END DO
CLOSE(10)

DO I = 1, ICOL
Expand All @@ -194,7 +190,7 @@ SUBROUTINE READ_DIREC(DIREC,NCOL,NROW,H,XC,YC,SIZE
ELSE IF (H(I,J) .EQ. 3) THEN
DIREC(I,J,1) = I-1
DIREC(I,J,2) = J
ELSE IF (H(I,J) .EQ. 4) THEN
ELSE IF (H(I,J) .EQ. 4) THEN
DIREC(I,J,1) = I-1
DIREC(I,J,2) = J-1
ELSE IF (H(I,J) .EQ. 5) THEN
Expand All @@ -217,4 +213,3 @@ SUBROUTINE READ_DIREC(DIREC,NCOL,NROW,H,XC,YC,SIZE
$ FILENAME
STOP
END

38 changes: 17 additions & 21 deletions src/rout.f
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,22 +9,18 @@ PROGRAM rout
c See WA Hydrology Homepage for operational details.

c Modified 5/99 to read in the uh_s array if it has already
c been generated in a previous run.
c been generated in a previous run.

c Modified 2/2001 by edm to include month and year in output
c and also check dates in VIC output files and calculate NDAYS
c
IMPLICIT NONE

c RCS ID STRING
CHARACTER*50 RCSID
DATA RCSID/"$Id: rout.f,v 1.1 2005/04/07 05:07:29 vicadmin Exp $"/

integer IARGC

integer isaleap
external isaleap

c change dimensions here
c nrow and ncol should be larger than the grid
c nyr should equal run length yrs+1
Expand All @@ -42,9 +38,9 @@ PROGRAM rout
PARAMETER (UH_DAY = 96 )
PARAMETER (TMAX = UH_DAY*24)
PARAMETER (PMAX = 10000 )

INTEGER DIREC(NCOL,NROW,2)
REAL VELO(NCOL,NROW), DIFF(NCOL,NROW)
REAL VELO(NCOL,NROW), DIFF(NCOL,NROW)
REAL XMASK(NCOL,NROW), FRACTION(NCOL,NROW)
REAL UH_BOX(PMAX,KE), UHM(NCOL,NROW,LE)
REAL UH_S(PMAX,KE+UH_DAY-1)
Expand All @@ -55,11 +51,11 @@ PROGRAM rout
INTEGER NO_OF_BOX
INTEGER CATCHIJ(PMAX,2)
INTEGER H(NCOL,NROW)

INTEGER PI, PJ
REAL UH_DAILY(PMAX,UH_DAY)
REAL FR(TMAX,2)

INTEGER NR
INTEGER IROW, ICOL
INTEGER LP,M,Y
Expand Down Expand Up @@ -174,7 +170,7 @@ PROGRAM rout
DO J=START_MO,12*(STOP_YEAR-START_YEAR)+STOP_MO
IF(M.EQ.2) THEN
LP=isaleap(Y)
ELSE
ELSE
LP=0
ENDIF
NDAY = NDAY+DAYS_IN_MONTH(M)+LP
Expand Down Expand Up @@ -205,11 +201,11 @@ PROGRAM rout
C Loop over required stations

100 CONTINUE
READ(10,*,END=110)
READ(10,*,END=110)
& NR, NAME, PI, PJ, AREA
READ(10,'(A80)',END=110) UH_STRING !new, AW: uh_string
IF (NR .EQ. 1) THEN
WRITE(*,'(I2,2X,A,I4,I4,G12.6)')
WRITE(*,'(I2,2X,A,I4,I4,G12.6)')
& NR, NAME, PI, PJ
PRINT*, 'Routing station: ', NAME
c note, the arrays are flipped left to right
Expand All @@ -220,7 +216,7 @@ PROGRAM rout
CALL SEARCH_CATCHMENT
& (PI,PJ,DIREC,NCOL,NROW,
& NO_OF_BOX,CATCHIJ,PMAX,IROW,ICOL)

print*, 'reading grid_UH...'
CALL READ_GRID_UH
& (UH_BOX, KE, PMAX, NO_OF_BOX, CATCHIJ,FILENAME)
Expand All @@ -237,25 +233,25 @@ PROGRAM rout
& CATCHIJ, BASE, RUNO, FLOW, KE, UH_DAY, UH_S, FRACTION,
& FACTOR_SUM,XC,YC,SIZE,DPREC,INPATH,ICOL,NDAY,
& IDAY,IMONTH,IYEAR, MO, YR, NYR)

print*, 'writing data...'
CALL WRITE_DATA
& (FLOW, NDAY, NAME5, FACTOR_SUM, OUTPATH,IDAY,IMONTH,IYEAR)
& (FLOW, NDAY, NAME5, FACTOR_SUM, OUTPATH,IDAY,IMONTH,IYEAR)

CALL WRITE_MONTH
& (DAYS_IN_MONTH,START_YEAR, STOP_YEAR, FIRST_YEAR,
& LAST_YEAR, START_MO, STOP_MO, FIRST_MO,
& (DAYS_IN_MONTH,START_YEAR, STOP_YEAR, FIRST_YEAR,
& LAST_YEAR, START_MO, STOP_MO, FIRST_MO,
& LAST_MO,
& NAME5, DAYS, FLOW, FACTOR_SUM, MONTHLY, MONTHLY_mm,
& NAME5, DAYS, FLOW, FACTOR_SUM, MONTHLY, MONTHLY_mm,
& YEARLY,YEARLY_mm,OUTPATH,NDAY,IMONTH,IYEAR,MO,YR,NMONTHS,NYR)


END IF
GOTO 100
110 CONTINUE

STOP
9001 WRITE(*,*) 'CANNOT OPEN: ', FILE_INPUT
9001 WRITE(*,*) 'CANNOT OPEN: ', FILE_INPUT
END
c ***********************************************
c FUNCTION ISALEAP
Expand Down
Loading

0 comments on commit ee66188

Please sign in to comment.