Fundamental Technologies

Cassini MIMI Pages

Numerical Computation of Energy-Dependent Geometric Factors of E and F Electron Detectors of CASSINI/MIMI/LEMMS, Technical Report by Xiaodong Hong and Thomas P. Armstrong, May 10, 1997

 

*******************************************************************************
*                                                                             *
*				PROGRAM I.1				      *
*                                                                             *
*******************************************************************************
*                                                                             *
*                               TRAJ4P.FOR                                    *
* THIS IS THE MAIN TRAJECTORY PROGRAM.  ALL VARIABLES AND PARAMETERS ARE      *
* DEFINED IN THE COMMENTS FOUND ON PAGE 65 OF SHODHAN'S THESIS.               *
*                                                                             *
*******************************************************************************
 
      	PROGRAM TRAJ1PHI
      	IMPLICIT NONE
      	CHARACTER*72 FNAME,FNAME1,FNAME2,FNAME3,FNAME4
      	INTEGER I,I1,IHLF,J,MAXE,NCOUNT,NDIM,NHIT,NPHI,NTHETA,NU
      	INTEGER NXK,NYK,NZK,MAXSURF,wscatter,SCATTER,nSCATTER,startnsurf
      	PARAMETER (MAXE=180)
      	REAL*8 C,CON,ERRWT
      	PARAMETER (CON=1.0D0,NDIM=6)
      	REAL*8 STEPHI,STEPTHETA,EK,ENEK,PHIMIN,EVELO,
     & 		PHI,QMC,QMCP,
     & 	THETA,THETAMAX,THETAMID,V,VX,VY,VZ
      	REAL*8 AUX(16,NDIM),DERY(NDIM),PRMT(5),X0,Y(NDIM),Y0,Z0
      	REAL CPUTIME,TIMER,ZTIM0
        logical firstcheck,secondcheck
	Integer SEED1,SEED2,SEED3,Seed4
	COMMON /SEED1/SEED1,/SEED2/SEED2,/SEED3/SEED3,/Seed4/Seed4
      	INCLUDE 'PASS5.CMN'
      	COMMON /QMC/QMC
      	COMMON /NCOUNT/NCOUNT,/NU/NU,/NHIT/NHIT
C     	COMMON /FNAME/FNAME
      	COMMON /FNAME4/FNAME4
      	COMMON /FNAME1/FNAME1
      	COMMON /FNAME2/FNAME2,FNAME3
	common /wscatter/wscatter,/scatter/scatter,/Nscatter/Nscatter
	common /theta/theta,/phi/phi,/EK/EK,/startnsurf/startnsurf
	common /firstcheck/firstcheck,/secondcheck/secondcheck
	data scatter/0/
	DATA SEED1/1234441519/
	DATA SEED2/278161611/
	DATA SEED3/467321899/
	Data Seed4/598516711/
      	DATA C/2.998D0/,QMCP/0.175602D0/
      	EXTERNAL FCT,OUTP
C    	ZTIM0=TIMER()
	write(6,*) 'ENTER THE CHOICE OF MODEL OF THE SCATTER'
	write(6,*) 'NON SCATTER----------->0'
	write(6,*) 'SPECULAR SCATTER------>1'
	write(6,*) 'DIFFUSE SCATTER------->2'
	read(5,*) wscatter    
	if (wscatter.eq.2) then
	    write(6,*) 'ENTER NSCATTER'
	    read(5,*) Nscatter
	end if
      	CALL GEOM    
    
C     	OPEN(UNIT=7,FILE='TRAJSH.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)
      	OPEN(UNIT=8,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME1)
      	ENCODE (19,110,FNAME2)NXK,NYK,NZK,IDINT(ENEK)
      	OPEN(UNIT=1,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME2)
      	ENCODE (20,120,FNAME3)NXK,NYK,NZK,IDINT(ENEK)
      	OPEN(UNIT=2,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME3)
      	ENCODE (23,125,FNAME4)NXK,NYK,NZK,IDINT(ENEK)
      	OPEN(UNIT=3,STATUS='NEW',ACCESS='SEQUENTIAL',FILE=FNAME4)
      	V = EVELO(EK)          !COMPUTE V FOR THIS ENERGY OF e
 
      	QMC = QMCP*DSQRT(1 - (V/C)**2)  
      	WRITE (6,*) 'ENTER MID. THETA,THETAMAX,DEL THETA'
      	READ (5,*) THETAMID,THETAMAX,STEPTHETA
      	WRITE (6,*) 'ENTER MIN. PHI,NO. OF STEPS,DEL PHI'
      	READ (5,*) PHIMIN,NPHI,STEPHI
      	WRITE (8,70) EK,V
      	WRITE (8,80) THETAMID,THETAMAX,STEPTHETA
      	WRITE (8,90) PHIMIN,NPHI,STEPHI
      	WRITE (8,20) PRMT(1),PRMT(2),PRMT(3),PRMT(4)
      	WRITE (8,*)
C      	WRITE (8,95) X0,Y0,Z0
      	WRITE (1,130)
      	WRITE (2,140)
      	WRITE (3,150)
C     TO CONVERT FROM INCHES TO CMS. 
C      	X0 = X0 * CON
C      	Y0 = Y0 * CON
C      	Z0 = Z0 * CON
	write(6,*)X0,Y0,Z0
      	WRITE (8,100) X0,Y0,Z0
      	WRITE (8,*)
 
      	NPAS=0 
      	DO I=1,NPHI
            PHI=PHIMIN + DFLOAT(I-1)*STEPHI       
            DO J=-1,1,2
          	IF (J .EQ. -1) THEN
           	    THETA=THETAMID
            	ELSE
           	    THETA=THETAMID + DFLOAT(J)*STEPTHETA
          	END IF
          	NHIT=2
          	DO WHILE (NHIT .EQ. 2 .AND. THETA .LT. THETAMAX)
C            	    ENCODE (23,50,FNAME)NXK,NYK,NZK,IDINT(ENEK),
C     &                          IDINT(THETA),IDINT(PHI)
C            	    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(6,*) J
C           TO COMPUTE THE VELOCITY PROJECTIONS Vx,Vy,Vz
            	    CALL VELOPROJ(V,VX,VY,VZ,THETA,PHI)
    
C            	    WRITE(6,*) ' THETA ',THETA,' PHI ',PHI
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.
C	            write(6,*)' THETA ',THETA,' PHI ',PHI
	    	    CALL DHPCG(PRMT,Y,DERY,NDIM,IHLF,FCT,OUTP,AUX)
C            	    IF (NHIT .EQ. 2) THEN
C             	    	PAS(NPAS,1) = THETA
C             	    	PAS(NPAS,2) = PHI
C            	    END IF
C            	    WRITE(6,*) 
C            	    WRITE(6,*) ' IHLF NO. OF BISECTIONS OF STEP: ',IHLF
C            	    WRITE(6,*) 'TOTAL NO. OF POINTS IN THE TRAJECTORY: ',NU-1	    
	    	    IF ((THETA.LT.Thetamax).AND.(THETA.GT.60)) NHIT = 2
	    		THETA = THETA+DFLOAT(J)*STEPTHETA
          	END DO
            END DO 
      	END DO
	WRITE (1,*) NPAS
	WRITE (2,*) NPAS	
      	CALL PASSOUTPUT
      	WRITE(6,*) 'NPAS:',NPAS
C     	CPUTIME=TIMER()-ZTIM0
      	WRITE(6,*) 'C.P.U. TIME: ',CPUTIME
      	WRITE(1,40)
      	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,'MIDDLE THETA(in deg.)',F10.2,2X,'MAX.THETA',
     & 	F10.2,2X,'DEL THETA',F10.2)
 90   	FORMAT (1X,'INITIAL PHI  (in deg.)',F10.2,2X,'NO. OF STEPS',
     & 	  I5,2X,'DEL PHI  ',F10.2)
 95   	FORMAT (1X,'STARTING POSITION(in in.)',2X,'X(+02)',D14.7,
     & 	  2X,'Y(+02)',D14.7,2X,'Z(+02)',D14.7)
100   	FORMAT (1X,'STARTING POSITION(in 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')
125   	FORMAT (I3,I3,I4,'.E',I5,'OUTVEL')
130   	FORMAT(1X,'POLAR & AZIMUTHAL ANGLES AT THE DETECTOR')
140   	FORMAT(1X,'POLAR & AZIMUTHAL ANGLES AT THE APERTURE')
150   	FORMAT(1X,'POLAR & AZIMUTHAL ANGLES AT THE APERTURE FROM VEL.')
      	CLOSE(1)
      	CLOSE(2)
      	CLOSE(3)     
      	close(8)
      	stop
      	END           

Return to Technical Report table of contents.
Return to Cassini MIMI table of contents page.
Return to Fundamental Technologies Home Page.

Updated 6/6/02, T. Hunt-Ward
tizby@ftecs.com