*$ CREATE MAGFLD.FOR *COPY MAGFLD * *===magfld=============================================================* * SUBROUTINE MAGFLD ( X, Y, Z, BTX, BTY, BTZ, B, NREG, IDISC ) INCLUDE '(DBLPRC)' INCLUDE '(DIMPAR)' INCLUDE '(IOUNIT)' * *----------------------------------------------------------------------* * * * Copyright (C) 1988-2010 by Alberto Fasso` & Alfredo Ferrari * * All Rights Reserved. * * * * * * Created in 1988 by Alberto Fasso` * * * * * * Last change on 06-Nov-10 by Alfredo Ferrari * * * * Input variables: * * x,y,z = current position * * nreg = current region * * Output variables: * * btx,bty,btz = cosines of the magn. field vector * * B = magnetic field intensity (Tesla) * * idisc = set to 1 if the particle has to be discarded * * * *----------------------------------------------------------------------* * Edit by Yang Yao trying to implent DIP and QUAD field to HFRS INCLUDE '(CMEMFL)' INCLUDE '(CSMCRY)' INCLUDE '(RTDFCM)' INCLUDE '(LTCLCM)' INCLUDE '(BEAMCM)' INCLUDE '(TRACKR)' * Define the string length and allocate space CHARACTER*8 REGNAM * LFIRST , to flag first call of routine LOGICAL LFIRST DATA LFIRST /.TRUE./ SAVE LFIRST * Magnetic filed data * (1)Dipole field(T) DATA DIPFLD /1.5923567D+00/ SAVE DIPFLD * (2)Quadrupole field(T/m) * Q1 DATA QUAFLD1 /4.33355D+00/ SAVE QUAFLD1 * Q2 DATA QUAFLD2 /4.665797D+00/ SAVE QUAFLD2 * Q3 DATA QUAFLD3 /2.51449D+00/ SAVE QUAFLD3 * Magnetic Regions:F:First;S:Second;T:Third SAVE QUADF, QUADS, QUADT, DIPF * Charge Correction for the magnetic field * CHFACT: Quadrupoles charge factor * Default values IDISC = 0 BTX = ZERZER BTY = ZERZER BTZ = ONEONE B = ZERZER *------------------------------------ * First Call: Initialisation *------------------------------------ IF (LFIRST)THEN * Some useful printout WRITE (LUNOUT,*) '' WRITE (LUNOUT,*) 'Dedicated magfld.f for U beam in HFRS' WRITE (LUNOUT,*) 'Initialisation started' WRITE (LUNOUT,*) '' WRITE (LUNOUT,*) '' WRITE (LUNOUT,*) 'Effective charge for U: ', Zfrttk,'+' WRITE (LUNOUT,*) 'Dipo field: ',DIPFLD,'T;' WRITE (LUNOUT,*) 'Quad field: ',QUAFLD1,'T/m;' WRITE (LUNOUT,*) 'Quad field: ',QUAFLD2,'T/m;' WRITE (LUNOUT,*) 'Quad field: ',QUAFLD3,'T/m;' WRITE (LUNOUT,*) '' * Retrieve the indices of the magnetic regions starting from their * names and store them CALL GEON2R("Q1Vac",QUADF,IERR) CALL GEON2R("Q2Vac",QUADS,IERR) CALL GEON2R("Q3Vac",QUADT,IERR) CALL GEON2R("Dip1Vac",DIPF,IERR) * Flush the .out file, in order to free the buffer: CALL FFLUSH CALL FLUSH(LUNOUT) * Avoid to initialise again the routine LFIRST =.FALSE. ENDIF * ASSIGN correct magnetic field to each REGION CALL GEOR2N (NREG,REGNAM,IERR) CHFACT = 34/Zfrttk IF(REGNAM .EQ. 'Q1Vac')THEN BTX = QUAFLD1*CHFACT*Y/1D+02 BTY = QUAFLD1*CHFACT*X/1D+02 BTZ = ZERZER ELSEIF(REGNAM .EQ. 'Q2Vac')THEN BTX = -QUAFLD2*CHFACT*Y/1D+02 BTY = -QUAFLD2*CHFACT*X/1D+02 BTZ = ZERZER ELSEIF(REGNAM .EQ. 'Q3Vac')THEN BTX = QUAFLD3*CHFACT*Y/1D+02 BTY = QUAFLD3*CHFACT*X/1D+02 BTZ = ZERZER ELSEIF(REGNAM .EQ. 'Dip1Vac')THEN BTX = ZERZER BTY = DIPFLD*CHFACT BTY = ZERZER ELSE CALL FLABRT('MAGFLD','NO MAGFLD IN THIS REGION!') ENDIF * Get field Module B = SQRT(BTX**2+BTY**2+BTZ**2) * Check against numerical precision IF (B.GT.1.0D-12)THEN * Normalised field components BTX = BTX/B BTY = BTY/B BTZ = BTZ/B ELSE * Restore default values BTX = ZERZER BTY = ZERZER BTZ = ONEONE B = ZERZER ENDIF RETURN *=== End of subroutine Magfld =========================================* END