view release on metacpan or  search on metacpan

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

	SUBROUTINE TXGRD(P,VT,VB,DZ,T,X,IN)
C THIS SUBROUTINE COMPUTES 2-WAY TRAVELTIME (T) AND HORIZONTAL DISTANCE (X)
C FOR A RAY OF RAY PARAMETER P TRAVELLING IN A CONSTANT VELOCITY GRADIENT
C LAYER OF THICKNESS DZ AND TOP AND BOTTOM VELOCITIES VT AND VB
C RESPECTIVELY.THE SORT OF COMPUTATION DONE IN THE SUBROUTINE IS TRANFERRED
C TO THE MAIN PROGRAM BY THE INTEGER  IN . WE HAVE THE FOLLOWING CODE:
C
C IN=0  IF 1/P L.E. VT  NO CONTUTATION IS CARRIED OUT IN THIS CASE
C
C IN=1 IF VB.L.E. 1/P.GT.VT. RAY IS TURNED WITHIN THE LAYER
C
C IN=2 IF 1/P.G.T.VB  RAY NOT TURNED WITHIN THE LAYER
C
C IN=3  IF  VT=VB (HOMOGENEOUS LAYER)
C
C TEST FOR P=0 *****
	IF(P.EQ.0.) GO TO 10
C
	PM1=1./P
C *** TEST FOR 1/P ****
	IF(PM1.GT.VT) GO TO 10
	IN=0
	RETURN
C
C *** TEST FOR HOMOGENEUOS LAYER  ****
C
 10	IF(VT.NE.VB) GO TO 30
C ***COMPUTE X AND T FOR HOMOGENEOUS LAYER  ****
	A1=P*VT
	A2=2.*DZ/SQRT(1.-A1**2)
	T=A2/VT
	X=A2*A1
	IN=3
	RETURN
C
C **COMPUTE 2./GRADIENT :::
 30	GM1T2=2.*DZ/(VB-VT)
C
	A1=1./VB
	IF(P.GE.A1) GO TO 50
C ***COMPUTATION FOR RAY THAT PASSES THROUGH THE LAYER ***
C
	A1=VB/VT
	A2=SQRT(1.-(P*VT)**2)
	A3=SQRT(1.-(P*VB)**2)
	T=GM1T2*ALOG(A1*(1.+A2)/(1.+A3))
	X=0.
	IF(P.EQ.0.) GO TO 45
	X=GM1T2*(A2-A3)/P
 45	IN=2
	RETURN
C *** COMPUTATION FOR TURNING RAY ***
C
 50	A1=P*VT
	A2=SQRT(1.-A1**2)
	T=GM1T2*ALOG((1.+A2)/A1)
	X=GM1T2*A2/P
	IN=1
	RETURN
	END