PROGRAM EXTRACT_MODFLOW_DATA C works on MODFLOW head files or MODFLOW-GWT conc. files CHARACTER*256 FILENAME,FILENAME2,OUTNAME CHARACTER*3 A CHARACTER*4 FTYPE character*16 text character*20 fmt integer i,j,k,kstp,kper,nc,nr,ilay,ncr,iu,kmax logical binary real pertim,totim INTEGER OUTUNIT,RANGEC,RANGER,RANGEL LOGICAL EXISTS ALLOCATABLE BUFF(:,:,:),BUFF2(:,:,:),BUFFR(:,:,:),BUFFL(:,:,:) ALLOCATABLE BUFF3(:,:,:),BUFF4(:,:,:) C open file INUNIT=9 INUNIT2=19 C read (*,*) 50 WRITE(*,*) ' Enter the name of the "EXTRACTOR" INPUT DATA FILE:' READ(*,'(A)') FILENAME INQUIRE(FILE=FILENAME,EXIST=EXISTS) IF(.NOT.EXISTS) THEN WRITE(*,*) ' File does not exist; PROGRAM STOPPING' STOP END IF OPEN(UNIT=INUNIT,FILE=FILENAME) c OPEN(UNIT=INUNIT,FILE=FILENAME,STATUS='OLD') C read input file READ(INUNIT,*) FTYPE READ(INUNIT,*) NCOL,NROW,NLAY 60 FORMAT(3I10) 62 FORMAT(E12.4) 64 FORMAT(A) IF (FTYPE.EQ.'HEAD'.OR.FTYPE.EQ.'DRAW'.OR.FTYPE.EQ.'VELO') THEN READ(INUNIT,*) KPER,KSTP 66 FORMAT (2I10) ELSE IF (FTYPE.EQ.'CONC') THEN READ(INUNIT,*) KPER,KSTP,IMOV END IF READ(INUNIT,*) ICOL1,ICOLEND READ(INUNIT,*) IROW1,IROWEND READ(INUNIT,*) ILAY1,ILAYEND C C READ DATA FILE WITH MATRIX TO BE MANIPULATED ALLOCATE (BUFF(NCOL,NROW,NLAY)) IF (FTYPE.EQ.'HEAD') THEN WRITE(*,*) ' Enter the name of the HEAD DATA FILE:' READ(*,'(A)') FILENAME2 INQUIRE(FILE=FILENAME2,EXIST=EXISTS) IF(.NOT.EXISTS) THEN WRITE(*,*) ' HEAD File does not exist; PROGRAM STOPPING' STOP END IF OPEN(UNIT=INUNIT2,FILE=FILENAME2) ELSE IF (FTYPE.EQ.'DRAW') THEN WRITE(*,*) ' Enter the name of the DRAW (drawdown) DATA FILE:' READ(*,'(A)') FILENAME2 INQUIRE(FILE=FILENAME2,EXIST=EXISTS) IF(.NOT.EXISTS) THEN WRITE(*,*) ' DRAW File does not exist; PROGRAM STOPPING' STOP END IF OPEN(UNIT=INUNIT2,FILE=FILENAME2) ELSE IF (FTYPE.EQ.'VELO') THEN ALLOCATE (BUFFR(NCOL,NROW,NLAY)) ALLOCATE (BUFFL(NCOL,NROW,NLAY)) WRITE(*,*) ' Enter the name of the VELO (velocity) DATA FILE:' READ(*,'(A)') FILENAME2 INQUIRE(FILE=FILENAME2,EXIST=EXISTS) IF(.NOT.EXISTS) THEN WRITE(*,*) ' VELO File does not exist; PROGRAM STOPPING' STOP END IF OPEN(UNIT=INUNIT2,FILE=FILENAME2) ELSE IF (FTYPE.EQ.'CONC') THEN WRITE(*,*) ' Enter the name of the CONC DATA FILE:' READ(*,'(A)') FILENAME2 INQUIRE(FILE=FILENAME2,EXIST=EXISTS) IF(.NOT.EXISTS) THEN WRITE(*,*) ' CONC File does not exist; PROGRAM STOPPING' STOP END IF OPEN(UNIT=INUNIT2,FILE=FILENAME2) END IF C C READ MATRIX TO BE MANIPULATED IF (FTYPE.EQ.'HEAD'.OR.FTYPE.EQ.'DRAW') THEN 800 format(1x,2i5,1p,2e15.6,1x,a,3i6,1x,a) C 55 k=0 kmax=nlay C 20 continue k=k+1 read(INUNIT2,800,end=80,err=90) kstp2,kper2,pertim,totim, & text,nc,nr,ilay,fmt if (nc.ne.ncol.or.nr.ne.nrow.or.ilay.ne.k) go to 90 do i=1,nr read(INUNIT2,fmt,end=80,err=90) (buff(j,I,K),j=1,nc) end do C if (k.lt.kmax) goto 20 C C CHECK IF AT CORRECT TIME LEVEL IF (KSTP.EQ.KSTP2.AND.KPER.EQ.KPER2) THEN GO TO 101 ELSE GO TO 55 END IF 101 CONTINUE ELSE IF (FTYPE.EQ.'CONC') THEN 56 k=0 kmax=nlay C 21 continue k=k+1 801 format(a13,27x,i5,7x,i5,7x,i5,7x,i5,9x,f10.0) C read(INUNIT2,801,end=80,err=90) text,ilay,imov2,kstp2,kper2,sumtch continue if (ilay.ne.k) go to 90 do i=1,nrow read(INUNIT2,777,end=80,err=90) (buff(j,I,K),j=1,ncol) end do C if (k.lt.kmax) go to 21 C C CHECK IF AT CORRECT TIME LEVEL IF (imov.eq.imov2.and.KSTP.EQ.KSTP2.AND.KPER.EQ.KPER2) THEN GO TO 102 ELSE GO TO 56 END IF 102 CONTINUE ELSE IF (FTYPE.EQ.'VELO') THEN 356 k=0 kmax=nlay C 221 continue k=k+1 802 format(52x,i5,7x,i5,7x,i5) C C READ VELOCITY IN COLUMN-DIRECTION read(INUNIT2,802,end=80,err=90) ilay,kstp2,kper2 continue if (ilay.ne.k) go to 90 do i=1,nrow read(INUNIT2,777,end=80,err=90) (buff(j,I,K),j=1,ncol) end do c if (k.lt.kmax) go to 221 C C READ VELOCITY IN ROW-DIRECTION c K=0 321 CONTINUE c K=K+1 read(INUNIT2,802,end=80,err=90) ilay,kstp2,kper2 continue if (ilay.ne.k) go to 90 do i=1,nrow read(INUNIT2,777,end=80,err=90) (buffR(j,I,K),j=1,ncol) end do c if (k.lt.kmax) go to 321 C C READ VELOCITY IN LAYER-DIRECTION c K=0 421 CONTINUE c K=K+1 read(INUNIT2,802,end=80,err=90) ilay,kstp2,kper2 continue if (ilay.ne.k) go to 90 do i=1,nrow read(INUNIT2,777,end=80,err=90) (buffL(j,I,K),j=1,ncol) end do C if (k.lt.kmax) go to 221 C C CHECK IF AT CORRECT TIME LEVEL IF (KSTP.EQ.KSTP2.AND.KPER.EQ.KPER2) THEN GO TO 402 ELSE GO TO 356 END IF 402 CONTINUE END IF C C CREATE NEW MATRIX C RANGEC=ICOLEND-ICOL1+1 RANGER=IROWEND-IROW1+1 RANGEL=ILAYEND-ILAY1+1 IF (RANGEL.EQ.1) THEN KEND=1 IEND=RANGEC JEND=RANGER ALLOCATE (BUFF2(IEND,JEND,KEND)) IF (FTYPE.EQ.'VELO') THEN ALLOCATE (BUFF3(IEND,JEND,KEND)) ALLOCATE (BUFF4(IEND,JEND,KEND)) END IF DO 200 I=1,RANGEC DO 200 J=1,RANGER I2=I+ICOL1-1 J2=J+IROW1-1 BUFF2(I,J,1)=BUFF(I2,J2,ILAY1) IF (FTYPE.EQ.'VELO') THEN BUFF3(I,J,1)=BUFFR(I2,J2,ILAY1) BUFF4(I,J,1)=BUFFL(I2,J2,ILAY1) END IF 200 CONTINUE ELSE IF (RANGEC.EQ.1) THEN KEND=1 IEND=RANGER JEND=RANGEL ALLOCATE (BUFF2(IEND,JEND,KEND)) IF (FTYPE.EQ.'VELO') THEN ALLOCATE (BUFF3(IEND,JEND,KEND)) ALLOCATE (BUFF4(IEND,JEND,KEND)) END IF DO 201 I=1,RANGER DO 201 J=1,RANGEL I2=I+IROW1-1 J2=J+ILAY1-1 BUFF2(I,J,1)=BUFF(ICOL1,I2,J2) IF (FTYPE.EQ.'VELO') THEN BUFF3(I,J,1)=BUFFR(ICOL1,I2,J2) BUFF4(I,J,1)=BUFFL(ICOL1,I2,J2) END IF 201 CONTINUE ELSE IF (RANGER.EQ.1) THEN KEND=1 IEND=RANGEC JEND=RANGEL ALLOCATE (BUFF2(IEND,JEND,KEND)) IF (FTYPE.EQ.'VELO') THEN ALLOCATE (BUFF3(IEND,JEND,KEND)) ALLOCATE (BUFF4(IEND,JEND,KEND)) END IF DO 202 I=1,RANGEC DO 202 J=1,RANGEL I2=I+ICOL1-1 J2=J+ILAY1-1 BUFF2(I,J,1)=BUFF(I2,IROW1,J2) IF (FTYPE.EQ.'VELO') THEN BUFF3(I,J,1)=BUFFR(I2,IROW1,J2) BUFF4(I,J,1)=BUFFL(I2,IROW1,J2) END IF 202 CONTINUE ELSE KEND=RANGEL IEND=RANGEC JEND=RANGER ALLOCATE (BUFF2(IEND,JEND,KEND)) IF (FTYPE.EQ.'VELO') THEN ALLOCATE (BUFF3(IEND,JEND,KEND)) ALLOCATE (BUFF4(IEND,JEND,KEND)) END IF DO 203 I=1,RANGEC DO 203 J=1,RANGER DO 203 K=1,RANGEL I2=I+ICOL1-1 J2=J+IROW1-1 K2=K+ILAY1-1 BUFF2(I,J,K)=BUFF(I2,J2,K2) IF (FTYPE.EQ.'VELO') THEN BUFF3(I,J,K)=BUFFR(I2,J2,K2) BUFF4(I,J,K)=BUFFL(I2,J2,K2) END IF 203 CONTINUE END IF C C WRITE OUTPUT FILE C C prepare file (append suffix) c WRITE(A,'(A)') FILENAME OUTNAME='CUT.'//filename OUTUNIT=10 OPEN(UNIT=OUTUNIT,FILE=OUTNAME) C IF (FTYPE.EQ.'HEAD'.OR.FTYPE.EQ.'DRAW') THEN DO 206 K=1,KEND DO 206 J=1,JEND WRITE (OUTUNIT,FMT) (BUFF2(I,J,K),I=1,IEND) 206 CONTINUE else IF (FTYPE.EQ.'CONC') THEN DO 207 K=1,KEND DO 207 J=1,JEND WRITE (OUTUNIT,777) (BUFF2(I,J,K),I=1,IEND) 777 FORMAT (1P10E12.4) 207 CONTINUE else IF (FTYPE.EQ.'VELO') THEN DO 208 K=1,KEND DO 208 J=1,JEND WRITE (OUTUNIT,777) (BUFF2(I,J,K),I=1,IEND) 208 CONTINUE DO 209 K=1,KEND DO 209 J=1,JEND WRITE (OUTUNIT,777) (BUFF3(I,J,K),I=1,IEND) 209 CONTINUE DO 211 K=1,KEND DO 211 J=1,JEND WRITE (OUTUNIT,777) (BUFF4(I,J,K),I=1,IEND) 211 CONTINUE ELSE WRITE (*,*) ' OUTPUT ERROR' END IF C go to 95 80 write (*,*) ' end of data during READ' go to 95 90 write (*,*) ' error during READ' 95 continue CLOSE(OUTUNIT) DEALLOCATE (BUFF,BUFF2) 100 CONTINUE c STOP END