*******************************************************************************
*                                                                             *
*			PROGRAM 2.1			   *
*                                                                             *
*******************************************************************************
*                                                                             *
*                               TRAJ3P.FOR                                    *
* THIS IS THE MODIFICATION TO THE MAIN PROGRAM NEEDED TO DO TRAJECTORY PLOTS. *
* PRIMARY DIFFERENCES BETWEEN THIS AND TRAJ6PT ARE EXPLAINED ON PAGE 108 OF   *
* SHODHAN'S THESIS.                                                           *
*                                                                             *
*******************************************************************************
 

      	PROGRAM TRAJ1PHI

      	IMPLICIT NONE
      	CHARACTER*72 FNAME,FNAME1,FNAME2,FNAME3
      	INTEGER I,I1,IHLF,J,MAXE,NCOUNT,NDIM,NHIT,NPAS,NPHI,NU,NTHETA
      	INTEGER NXK,NYK,NZK,NTRAJ,wscatter,SCATTER,nSCATTER,k,m,
     &	startnsurf
        logical firstcheck,secondcheck
      	PARAMETER (MAXE=180)
      	REAL*8 C,CON,ERRWT
      	REAL*8 STEPHI,STEPTHETA,EK,ENEK,PHIMIN,PHIMAX,EVELO,PHI(MAXE),
     &		QMC,QMCP,V,VX,VY,VZ,THETA(MAXE),THETAMIN

	INTEGER SEED1,SEED2,SEED3,seed4
	COMMON /SEED1/SEED1,/SEED2/SEED2,/SEED3/SEED3,/seed4/seed4

      	INCLUDE 'PASS.CMN'
      	COMMON /QMC/QMC
      	COMMON /NCOUNT/NCOUNT,/NU/NU,/FNAME/FNAME,/NHIT/NHIT
      	COMMON /FNAME1/FNAME1
      	COMMON /FNAME2/FNAME2,FNAME3

	common /wscatter/wscatter,/scatter/scatter,/Nscatter/Nscatter
	common /theta/theta,/phi/phi,/startnsurf/startnsurf
	common /firstcheck/firstcheck,/secondcheck/secondcheck              
	data scatter/0/,Nscatter/10/
	DATA SEED1/1234441519/
	DATA SEED2/278161611/
	DATA SEED3/467321899/
	Data Seed4/598516711/

	
      	DATA C/2.998D0/,QMCP/0.175602D0/
      	PARAMETER (CON=1.00D0,NDIM=6)
      	REAL*8 AUX(16,NDIM),DERY(NDIM),PRMT(5),X0,Y(NDIM),Y0,Z0
      	REAL CPUTIME,TIMER,ZTIM0

      	EXTERNAL FCT,OUTP
      	ZTIM0=TIMER()

c1	read(5,*) m
c1	write(6,*)m
c1      do k=1,m	

	write(6,*) 'ENTER THE CHOICE OF THE MODEL OF SCATTER'
	write(6,*) 'NON SCATTER----------->0'
	write(6,*) 'SPECULAR SCATTER------>1'
	write(6,*) 'DIFFUSE SCATTER------->2'
	read(5,*) wscatter
      

      	CALL GEOM    
    

C     	OPEN(UNIT=3,FILE='TRAJ.DAT',ACCESS='SEQUENTIAL',STATUS='NEW')
C     	OPEN(UNIT=4,FILE='TRAJCO.DAT',ACCESS='SEQUENTIAL',STATUS='NEW')
c      	OPEN(UNIT=7,FILE='TRAJSH.DAT',ACCESS='SEQUENTIAL',STATUS='NEW')
C     	OPEN(UNIT=8,FILE='PASS.DAT',ACCESS='SEQUENTIAL',STATUS='NEW')
      

      	WRITE(6,*) 'ENTER THE INITIAL COORDINATES OF THE e IN CM/100s'
      	READ(5,*) X0,Y0,Z0

 
C     TO ENTER QUANTITIES FOR PRMT
      	WRITE(6,*) 'ENTER LOWER BOUND ON TIME (t=0)'
      	READ(5,*) PRMT(1)
      	WRITE(6,*) 'ENTER THE UPPER BOUND ON TIME'
      	READ(5,*) PRMT(2)
      	WRITE(6,*) 'ENTER THE TIME STEP'
      	READ(5,*) PRMT(3)
      	WRITE(6,*) 'ENTER THE ERROR BOUND'
      	READ(5,*) PRMT(4)
   

      	WRITE(6,*) 'ENTER THE ENERGY OF THE PARTICLE IN MEV'
      	READ(5,*) EK

      	ENEK=IDINT(1000.0D0*EK)
      	NXK=IDINT(10000.0D0*dabs(X0))
      	NYK=IDINT(10000.0D0*Y0)
      	NZK=IDINT(10000.0D0*Z0)

      	ENCODE (17,60,FNAME1)NXK,NYK,NZK,IDINT(ENEK)
c      	OPEN(UNIT=8,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME1)

      	ENCODE (19,110,FNAME2)NXK,NYK,NZK,IDINT(ENEK)
c      	OPEN(UNIT=1,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME2)

      	ENCODE (20,120,FNAME3)NXK,NYK,NZK,IDINT(ENEK)
c      	OPEN(UNIT=2,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME3)

      	V = EVELO(EK)          !COMPUTE V FOR THIS ENERGY OF e
      	QMC = QMCP*DSQRT(1 - (V/C)**2)  

      	WRITE(6,*) 'ENTER MIN. THETA,NO. OF STEPS,DEL THETA'
      	READ (5,*) THETAMIN,NTHETA,STEPTHETA

      	WRITE(6,*) 'ENTER MIN. PHI,NO. OF STEPS,DEL PHI'
      	READ (5,*) PHIMIN,NPHI,STEPHI

c      	WRITE (8,70) EK,V
c      	WRITE (8,80) THETAMIN,NTHETA,STEPTHETA
c      	WRITE (8,90) PHIMIN,NPHI,STEPHI
c      	WRITE (8,20) PRMT(1),PRMT(2),PRMT(3),PRMT(4)
c      	WRITE (8,*)
c      	WRITE (8,95) X0,Y0,Z0

c      	WRITE (1,150)
c      	WRITE (2,140)

      	X0 = X0 * CON
      	Y0 = Y0 * CON
      	Z0 = Z0 * CON

c      	WRITE(8,100) X0,Y0,Z0
c      	WRITE(8,*)

      	DO I=1,NTHETA
	    THETA(I)=THETAMIN + DFLOAT(I-1)*STEPTHETA
      	END DO

      	DO I=1,NPHI
       	    PHI(I) = PHIMIN + (I-1)*STEPHI
      	END DO
 

      	NPAS=0
      	DO J=1,NTHETA 
       	    DO I=1,NPHI

        	ENCODE (23,50,FNAME)NXK,NYK,NZK,IDINT(ENEK),
     &                   IDINT(THETA(J)),IDINT(PHI(I))
        	OPEN(UNIT=4,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME)
                

        	NCOUNT = 0        
       		NU = 1
 

c        	WRITE(7,05) I
c        	WRITE(7,10) EK,V,THETA(J),PHI(I)
c        	WRITE(7,20) PRMT(1),PRMT(2),PRMT(3),PRMT(4)
c       	 	WRITE(7,30)

c        	WRITE(4,*) 
c        	WRITE(4,05) (k+54)
c		WRITE(4,*)

C       TO COMPUTE THE VELOCITY PROJECTIONS Vx,Vy,Vz
        	CALL VELOPROJ(V,VX,VY,VZ,THETA(j),PHI(I))
   

		WRITE(6,*)
       	WRITE(6,*) ' X: ',X0,' Y: ',Y0,' Z: ',Z0
       	WRITE(6,*) 'ENERGY (in mev) ',EK,' V ',V,' QMC ',QMC
c        	WRITE(6,*) ' THETA( ',J,')',THETA(J),' PHI( ',I,')',PHI(I)
C       	WRITE(6,*) ' VX ',VX,' VY ',VY,' VZ ',VZ

C       TO INITIALISE THE INITIAL VALUES 
        	Y(1) = X0
        	Y(2) = Y0
        	Y(3) = Z0
        	Y(4) = VX
        	Y(5) = VY
        	Y(6) = VZ

        	ERRWT = 1.0D0/6.0D0
        	DO I1=1,NDIM
        	    DERY(I1) = ERRWT
                END DO
		    startnsurf=1
		    firstcheck=.false.
		    secondcheck=.false.
            	CALL DHPCG(PRMT,Y,DERY,NDIM,IHLF,FCT,OUTP,AUX)

            	IF (NHIT.EQ.2) THEN
	    	    PAS(NPAS,1) = THETA(J)
            	    PAS(NPAS,2) = PHI(I)
	    	END IF

c            	WRITE(7,40)
C           	WRITE(3,40)
c            	WRITE(4,40)
 

C           	WRITE(6,*) ' Y ',(Y(I1),I1=1,NDIM)
            	WRITE(6,*) ' IHLF NO. OF BISECTIONS OF STEP: ',IHLF
            	WRITE(6,*) 'TOTAL NO. OF POINTS IN THE TRAJECTORY: ',NU-1
       	    END DO
      	END DO
c1       end do

      	CALL PASSOUTPUT

      	WRITE(6,*) 'NPAS',NPAS
      	CPUTIME=TIMER()-ZTIM0
      	WRITE(6,*) 'C.P.U. TIME: ',CPUTIME
c      	WRITE(1,40)
c      	WRITE(2,40)

 05   	FORMAT(1X,I3)
 10   	FORMAT(1X,'ENERGY(in mev)',D10.3,2X,'VELOCITY(*10+10)',F12.5,2X,
     & 	    'THETA(in deg.)',F10.3,2X,'PHI(in deg.)',F10.3)
 20   	FORMAT(1X,'INITIAL TIME(*10-08)',F12.6,2X,'FINAL TIME(*10-08)',
     & 	  F12.6,2X,'INITIAL STEP(*10-08)',F14.8,1X,'ERROR BOUND',F19.12)
 30   	FORMAT(4X,'T(-08)',10X,'X(+02)',10X,'Y(+02)',10X,'Z(+02)',10X,
     & 	  'VX(+10)',10X,'VY(+10)',10X,'VZ(+10)',10X,'(V+10)')
 40   	FORMAT(X,'---*---*---*--- END OF ENERGY ---*---*---*---')
 50   	FORMAT (I3,I3,I4,'.E',I5,I3,I3)
 60   	FORMAT (I3,I3,I4,'.E',I5)
 70   	FORMAT (1X,'ENERGY(IN MEV)',D10.3,2X,'VELOCITY(*10+10)',F12.5)
 80   	FORMAT (1X,'INITIAL THETA',F10.2,2X,'NO. OF STEPS',
     &   I5,2X,'DEL THETA',F10.2)
 90   	FORMAT (1X,'INITIAL PHI',F10.2,2X,'NO. OF STEPS',
     &   I5,2X,'DEL PHI',F10.2)
 95   	FORMAT (1X,'STARTING POSITION(INCHES)',2X,'X(+02)',D14.7,
     &   2X,'Y(+02)',D14.7,2X,'Z(+02)',D14.7)
100   	FORMAT (1X,'STARTING POSITION(CM)    ',2X,'X(+02)',D14.7,
     &   2X,'Y(+02)',D14.7,2X,'Z(+02)',D14.7)
110   	FORMAT (I3,I3,I4,'.E',I5,'IN')
120   	FORMAT (I3,I3,I4,'.E',I5,'OUT')
130   	FORMAT (1X,'POLAR & AZIM ANGLES AT DETECTOR')
140   	FORMAT (1X,'POLAR & AZIM ANGLES AT THE APERATURE')
150     format (1x,'The number of file and filename for the escaped e`s
     &		trajectory') 
      	CLOSE(1)
c      	CLOSE(2)
c      	CLOSE(7)      
      	CLOSE(4)
c      	CLOSE(8)
      	STOP
      	END

 

Return to the Table of Contents for Hong's MS Thesis