view release on metacpan or  search on metacpan

lib/App/SeismicUnixGui/fortran/src/readmmod.for  view on Meta::CPAN

	VSVPR=0.0
	CALL READ_PAR_R4('VS/VP RATIO TO DEFAULT VS?',VSVPR)
	CALL READ_PAR_I4(' IOP TO DEFINE DENSITIES ?',IOP)
c
	I=1
	DE1=1.5
	AU1=0.
 10	write(*,*) 'LAYER ',I
	write(*,*)
	VT(I)=DE1
	CALL READ_PAR_R4('VPTOP (KM/SEC) ?',VT(I))
	IF(VT(I).LT.0) GO TO 100
	CALL READ_PAR_R4('VPBOTTOM (KM/SEC) ?',VB(I))
	CALL READ_PAR_R4('DEPTH OR TWTT OR THICKNESS (KM OR SEC)?',AU2)
	VST(I)=VSVPR*VT(I)
	VSB(I)=VSVPR*VB(I)
	CALL DENFVP(VT(I),RHOT(I),IOP)
	CALL DENFVP(VB(I),RHOB(I),IOP)
	IF(I.EQ.1) THEN
	VST(I)=0.0
	VSB(I)=0.0
	RHOT(I)=1.0
	RHOB(I)=1.0
	ENDIF
	CALL READ_PAR_R4(
     +  'VSTOP (KM/SEC)? (ENTER 0.0 TO DEFINE A FLUID)',VST(I))
	CALL READ_PAR_R4(
     +  'VSBOTTOM (KM/SEC)? (ENTER 0.0 TO DEFINE A FLUID)',VSB(I))
	CALL READ_PAR_R4('TOP DENSITY (G/CC) ?',RHOT(I))
	CALL READ_PAR_R4('BOTTOM DENSITY (G/CC) ?',RHOB(I))
c
	IF(ID.EQ.0) THEN
	DZ(I)=AU2-AU1
	AU1=AU2
	ELSE
		IF(ID.EQ.1) THEN
		AU3=AU2-AU1 ! LAYER I TWO WAY TRAVEL TIME
		AU1=AU2
		CALL THI(AU3,VT(I),VB(I),DZ(I))
		ELSE
		DZ(I)=AU2
		ENDIF
	ENDIF
	I=I+1
	DE1=VB(I-1)
	GO TO 10
100	NL=I-1
   	VB(I)=0.0
	DZ(I)=0.0
	VST(I)=0.0
	VSB(I)=0.0
	RHOT(I)=0.0
	RHOB(I)=0.0
c
c *** WRITE MODEL IN FILE ***
c
97      ID=1
  	CALL READ_PAR_I4('1- WRITE MODEL IN FILE, 0- NO',ID)
	IF(ID.EQ.0) GO TO 135
	STAT='NEW'
103	write(*,*)'OUTPUT FILE NAME ?? '
	READ(*,'(A)') FOUT
	INQUIRE(FILE=FOUT,EXIST=EX)
	IF(EX) THEN
	write(*,*) 'FILE ALREADY EXISTS: 1- OVERWRITE IT,0- TRY ',
     +  'AGAIN WITH A NEW NAME'
	read(5,*) ID
		IF(ID.EQ.0) GO TO 103
		IF(ID.EQ.1) STAT='UNKNOWN'
	ENDIF
	OPEN(UNIT=IOUT,FILE=FOUT,STATUS=STAT,FORM='UNFORMATTED')
	DO K=1,NL+1
	write(IOUT) VT(K),VB(K),DZ(K),VST(K),VSB(K),RHOT(K),RHOB(K)
	ENDDO
	CLOSE(UNIT=IOUT)
	GO TO 135

! *** READ OLD MODEL FROM DISK ***

110	INQUIRE(FILE=inbound_model,EXIST=EX)
!      print *, 'readmmod.for,inbound_model=',inbound_model
        IF(EX) then
!		write(*,*) ''
!		write(*,*) 'readmmod,Default file immodpg.out exists.'
!		write(*,*) ''
!      call read_par_i4('0- use it , 1- No',idef)
!              if(idef.eq.1) 117
              go to 117
	 else
	       print *, 'Default file immodpg.out is missing!'
 !             if(idef.eq.1) go to 115
              go to 115

	end if

c
115     write(*,*) 'INPUT FILE NAME ?? '
	READ(5,'(A)') inbound_model
	INQUIRE(FILE=inbound_model,EXIST=EX)
	IF(.NOT.EX) THEN
	write(*,*)'FILE DOES NOT EXIST, TRY AGAIN WITH A NEW NAME'
	GO TO 115
	ENDIF
c
117	continue
	OPEN(UNIT=IIN,FILE=inbound_model,STATUS='OLD',FORM='UNFORMATTED')
	K=1
120     READ(IIN) VT(K),VB(K),DZ(K),VST(K),VSB(K),RHOT(K),RHOB(K)
	IF(VT(K).LT.0.) GO TO 125
	K=K+1
	GO TO 120
125	NL=K-1
c
	CLOSE(UNIT=IIN)
!   	write(*,*)'readmmod,FILE = ',inbound_model
	CALL WRIMOD2(NL,VT,VB,DZ,VST,VSB,RHOT,RHOB)
c
c *** MODIFICATIONS ***
c
	ID=0
!	Juan's modification for immodpg.for July 25 2020
!       write(*,*) 'readmmod.for L207'
!	CALL READ_PAR_I4('CHANGE THIS FILE? 1-YES 0-NO',ID)
!        write(*,*) 'readmmod.for L209 ID=', ID
!       FORCE NO READING OF MODEL
	ID=0
	IF(ID.NE.1) GO TO 135
127	write(*,*)'1-DELETE OR 2-INSERT AFTER LAYER,3- CHANGE LAYER '
	read(5,*) IMOD
	write(*,*)'LAYER NUMBER ?? '
	read(5,*) LNU
c
	IF(IMOD.EQ.1) THEN
	DO K=LNU,NL
   	VT(K)=VT(K+1)
	VB(K)=VB(K+1)
	DZ(K)=DZ(K+1)
	VST(K)=VST(K+1)
	VSB(K)=VSB(K+1)
	RHOT(K)=RHOT(K+1)
	RHOB(K)=RHOB(K+1)
	ENDDO
	NL=NL-1
	GO TO 130
	ENDIF
c
	I=LNU
	IF(IMOD.EQ.2) THEN
	K=NL+1
	I=LNU+1
57	CONTINUE
   	VT(K+1)=VT(K)
	VB(K+1)=VB(K)
	DZ(K+1)=DZ(K)
	VST(K+1)=VST(K)
	VSB(K+1)=VSB(K)
	RHOT(K+1)=RHOT(K)