App-SeismicUnixGui
view release on metacpan or search on metacpan
lib/App/SeismicUnixGui/fortran/src/mmodpg4L_SU_Aug27_20_emilio/mmodpg2.for view on Meta::CPAN
c **** PARAMETERS FOR THE X-T PLOT ****
c
write(*,*) 'TO DEFINE PLOTTING AREA ENTER : '
write(*,*) ' '
rvinv=0.
call read_par_r4('Reducing Velocity (km/s),(0-For none)',rv)
if (rv.eq.0.) go to 40
rvinv=1./rv
40 continue
if(idrdtr.eq.1) then
call read_par_i4('1-Traces are red. by this vel., 0-No',idred)
endif
c
call read_par_r4('MINIMUN DISTANCE (KM)?',xmin)
call read_par_r4('MAXIMUN DISTANCE (KM)?',xmax)
call read_par_r4('MINIMUN TIME (SEC)?',tmin)
call read_par_r4('MAXIMUN TIME (SEC)?',tmax)
c **** RE-CHECK PARAMETERS ****
ID=0
call read_par_i4('1- RECHECK PARAMETERS,0- NO',ID)
if(ID.eq.1) go to 15
c
c
rvinvd = 0.0
if(idred.eq.0) rvinvd = rvinv
c *******************
c ***** READ VELOCITY DEPTH MODEL ******
c
call READMMOD(VT,VB,DZ,VST,VSB,RHOT,RHOB,nl)
c
call read_par_i4('Working layer number ?',lch)
if(lch.lt.1) lch = 1
if(lch.gt.nl) lch = nl
c
c ***************** Open X-T Window ****************
c
call pgbegin(0,'?',1,1)
c call pgpaper(12.0,8.5/11.0)
c call pgpaper(17.0,0.6)
call pgpaper(14.0,0.6)
call pgask(flag)
************************************
c
10 continue
c
c Check for vt-vb too small.
c
A3 = 0.001
c
do 20 i = 1, nl
20 if(ABS(VT(I)-VB(I)).le.A3) vb(i) = vt(i)
c
c Reflections at the bottom of layers decided automatically
c based on velocity discontinuities
c
do 34 I=1,nl-1
IR(I)=0
34 if(ABS(VT(I+1)-VB(I)).GT.A3) IR(I)=1
IR(nl) = 0 !** No reflection at the bottom of model
c
c *** COMPUTATIONS ***
c
DZ1TEM=DZ(1)
DZ(1)=DZ(1)-(SDEPTH+RDEPTH)/2. ! Correct only if 1rst layer
c is a constant vel. layer.
call txpr(nl,VT,VB,DZ,PMIN,PMAX,DP,IR,multin,ntp,ILA,P,X,T)
DZ(1)=DZ1TEM
c write(*,*) 'total # of computed points ',ntp
c
c *** PLOTTING ***
c
55 continue
c call pgvport(0.1,0.9,0.1,0.9)
call pgvport(0.075,1.0,0.08,0.925)
call pgwindow(xmin,xmax,tmax,tmin)
call pgbox('BCTN',0.0,0,'BCTN',0.0,0)
call pglabel('X(km)','T - X/Vred (s)',
c call pglabel('X(km)','T - X/8.0 (seg)',
+ 'Curvas camino-tiempo (X-T), Ondas P, Modelo 1-D')
c + 'Primary P-wave T-X Curves, 1-D Model')
c
if(idrdtr.eq.1) then
call pgsci(1)
do j=1,ntr
xline = datax1 + (j-1) * datadx
c datax = sqrt(psld2 + xline*xline)
datax = xline
a4 = datax + clip * datadx
a5 = datax - clip * datadx
do i=1,ns
a3 = f(j,i)
xa1(i) = datax + gain * a3
if(xa1(i).gt.a4) xa1(i) = a4
if(xa1(i).lt.a5) xa1(i) = a5
xa3(i) = xa1(i)
if(a3.lt.0.0) xa3(i) = datax
xa2(i) = datat1 + (i-1)*datadt - abs(datax)*rvinvd
enddo
c
call pgline(ns,xa1,xa2)
xa3(1) = datax
xa3(ns) = datax
c pgsci(1) selecciona blanco (negro en papel) para rellenar trazas
c pgsci(3) relleno en verde
call pgsci(3)
call pgpolyev(ns,xa3,xa2)
c call pgpoly(ns,xa3,xa2)
call pgsci(1)
enddo
endif
c
do 140 j=1,lch
nplj = 0
do 120 i=1,ntp
c ** SELECT LAYER J **
if(ABS(ILA(I)).NE.J) go to 120
nplj = nplj+1
xa1(nplj) = X(I)
xa5(nplj) = -X(I)
xa2(nplj) = T(I) - X(I) * rvinv !** RED. TIME
xa3(nplj) = P(I)
lib/App/SeismicUnixGui/fortran/src/mmodpg4L_SU_Aug27_20_emilio/mmodpg2.for view on Meta::CPAN
call read_par_r4('Multiplicative constant ??',vmf)
do i = 1,nl
VT(i) = vmf * VT(i)
VB(i) = vmf * VB(i)
enddo
endif
c
c **************************
c
icolor = 0
call pgpage
call pgsci(1)
go to 10
c
255 continue
c
c write modified model to terminal
c
write(*,*) ' '
write(*,*) 'MODIFIED MODEL:'
call WRIMOD2(nl,VT,VB,DZ,VST,VSB,RHOT,RHOB)
c
c write modified model to file mmodpg.out
c
OPEN(UNIT=IOUT,FILE='mmodpg.out',STATUS='UNKNOWN',
+ 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)
c
write(*,*) 'This model has been written to file:'
write(*,*) '*** mmodpg.out ***'
c
call pgend
500 format(' 1- VTOP = ',f6.3,', 2- VBOT = ',f6.3,' (km
+/s)')
502 format(' 3- DZ = ',f7.4,' (km)',', 4- VTOP and VBOT')
503 format(' 7- Increment = ',f7.4,' (km or km/s)')
504 format(' 9- Gain for data, 10- DP = ',f9.6,
+' (s/km)')
end
c
c ****************************************************
c
SUBROUTINE rdata(f, m, n, ntr, ns, datadt, fmin, fmax)
INTEGER m,n,ntr,ns,idtusec
REAL f(m,m), fmin, fmax
CHARACTER*40 FIN
LOGICAL EX
c
c Read data parameters
c
iin=26
OPEN(UNIT=iin,FILE='parmmod',STATUS='OLD')
read(iin,*) ntr,ns,idtusec
close(iin)
datadt = float(idtusec) * 1e-6
c
c115 write(*,*) 'INPUT DATA FILE NAME ?? '
c READ(5,'(A)') FIN
c INQUIRE(FILE=FIN,EXIST=EX)
c if(.NOT.EX) then
c write(*,*)'FILE DOES NOT EXIST, TRY AGAIN WITH A NEW NAME'
c go to 115
c endif
c
c OPEN(UNIT=iin,FILE=FIN,STATUS='OLD',FORM='UNFORMATTED')
c
c Read data File
c
OPEN(UNIT=iin,FILE='datammod',STATUS='OLD',FORM='UNFORMATTED')
k=1
120 READ(iin) (f(k,i), i=1,ns)
if(k.GE.ntr) go to 125
k=k+1
go to 120
125 CLOSE(UNIT=iin)
c
fmin = 1e30
fmax = -1e30
do 20 i=1,ntr
do 10 j=1,ns
fmin = min(f(i,j),fmin)
fmax = max(f(i,j),fmax)
10 continue
20 continue
c
write(*,*) 'Data min, max = ',fmin,fmax
write(*,*)
c
END
c
c ************************************************
c
subroutine read_datxy(x,y,n,fin,iin,iwrit)
c
c Reads file containing (x,y) pairs and an arbitrary number of comment lines
c in between. The (x,y) pairs are returned in arrays x and y. The
c comment lines are written on the terminal.
c
c n = number of (x,y) pairs in the file (output).
c
c fin = Default input file name (input)
c iin = reading input unit (input).
c iwrit = if iwrit.eq.1, (x,y) pairs are written on screen (input).
c
dimension x(*),y(*)
character*40 fin
character*40 comment
LOGICAL EX
c
c ********* read input file *********
c
c Check for default file "fin". If file does not exist, then ask for
c a file name
c
c write(*,*) fin
go to 117
115 write(*,*) 'Input File Name ?? '
( run in 1.366 second using v1.01-cache-2.11-cpan-99c4e6809bf )