*********************************************************************** SUBROUTINE MDCSTEP c c store hit info for MDC set DRIF as: c c HITS(1) = x of current particle in cave system c HITS(2) = y " c HITS(3) = z " c HITS(4) = px of current particle in cave system c HITS(5) = py " c HITS(6) = pz " c HITS(7) = x in MDC system (TRD1 shape used for drift volumes, ie x changes c HITS(8) = y in MDC system along z --> y and z are permuted) c HITS(9) = kinetic energy of current particle (in GeV) c HITS(10) = current particle ID c HITS(11) = TOF of particle (in sec) c HITS(12) = distance along wire (= xprime) c HITS(13) = perpendicular distance to central wire in plane (= 0-wire) c HITS(14) = volume copy number (of level=3) c HITS(15) = angle of incidence on plane (0-180 deg) c HITS(16) = azimuthal angle (0-360 deg) c ->needed to compute the shortest distance to nearest sense wire c c c last modified on: 29/11/06 by R.Holzmann *********************************************************************** IMPLICIT NONE #include "drtups.inc" #include "mdcdconst.inc" #include "geant321/gcflag.inc" #include "geant321/gckine.inc" #include "geant321/gcking.inc" #include "geant321/gclist.inc" #include "geant321/gcsets.inc" #include "geant321/gctmed.inc" #include "geant321/gctrak.inc" #include "geant321/gcvolu.inc" #ifdef CALC_REF_MATRIX #include "guref.inc" #endif #include "kinetups.inc" INTEGER NH PARAMETER (NH=16) REAL HITS(NH) CHARACTER*4 CHVOL(10), CHSET EQUIVALENCE (IHSET,CHSET) EQUIVALENCE (NAMES,CHVOL) REAL XM(3), XD(3), PM(3), PD(3), PDD(3), XDD(2), PDY REAL DEG, ANG, FP3, PXY INTEGER I, J, IHIT, ID1, ID2, IP3 LOGICAL LVOL INTEGER ISTEP SAVE ISTEP IF((INWVOL.EQ.1 .OR. INWVOL.EQ.2) + .AND. CHSET.EQ.'DRIF' .AND. CHARGE.NE.0.) THEN ! charged track has entered or left MDC drift volume c IF(VECT(7)<0.001) RETURN ! momentum is below 1 MeV, abort IF(VECT(7)<0.0002) RETURN ! momentum is below 0.2 MeV, abort ID1 = IDTYPE/10 ! identify wire plane via IDTYPE... (gcsets) ID2 = IDTYPE - 10*ID1 c IF(ID2.LE.6 .AND. INWVOL.EQ.2) RETURN ! a sense plane is left ANG = PHI(ID1,ID2) ! ...and select appropriate wire angle IF(MDCSHAPE(1:4).NE.'TRD1') + WRITE(6,*) ' MDCSTEP: Wrong shape: ', mdcshape DO I = 1,3 HITS(I) = VECT(I) ! position of hit and particle momentum HITS(I+3) = VECT(I+3)*VECT(7) XM(I) = VECT(I) PM(I) = VECT(I+3) ENDDO CALL GMTOD(XM,XD,1) ! transform coordinates into MDC system CALL GMTOD(PM,PD,2) ! transform momenta into MDC system c rotate momentum into wire plane, i.e. xprime-yprime where xprime is c parallel to wires PDD(1) = -PD(1)*COS(ANG) - PD(3)*SIN(ANG) PDD(2) = PD(1)*SIN(ANG) - PD(3)*COS(ANG) PDD(3) = PD(2) ! pz is along GEANT y axis of drift volumes FP3 = PDD(3)/ABS(PDD(3)) ! particle goes forwards or backwards ? c IF(ID2.EQ.7) THEN ! a cathode plane is IF((INWVOL.EQ.1 .AND. FP3.LT.0) .OR. ! entered backwards + (INWVOL.EQ.2 .AND. FP3.GT.0)) RETURN ! or left forwards c ENDIF ! these gymnastics are necessary to guarantee that the ! reference hit is always registered on the front of the cathode IF (NSTEP-ISTEP.EQ.1) THEN ! next hit follows to close to previous ISTEP = NSTEP RETURN ENDIF ISTEP = NSTEP PDY = PDD(2)*FP3/SQRT(PDD(2)**2+PDD(3)**2) ! cos of angle in yprime-z c rotate positions into wire plane (and permute y and z axes) XDD(1) = (-XD(1)+MDCDX(ID2,ID1))*COS(ANG) - + (XD(3)+MDCDY(ID2,ID1))*SIN(ANG) XDD(2) = (-XD(1)+MDCDX(ID2,ID1))*SIN(ANG) - + (XD(3)+MDCDY(ID2,ID1))*COS(ANG) HITS(7) = -XD(1)+MDCDX(ID2,ID1) HITS(8) = XD(3)+MDCDY(ID2,ID1) HITS(9) = GEKIN HITS(10) = FLOAT(IPART) HITS(11) = TOFG HITS(12) = XDD(1) ! xprime HITS(13) = XDD(2) ! yprime HITS(14) = FLOAT(NUMBER(3)) c HITS(15) = PDY ! cos yprime-z HITS(15) = 57.29577*ACOS(PD(2)) ! angle of incidence (with y <-> z) PXY = SQRT(PD(1)**2+PD(3)**2) ! with y <-> z IF(PXY.NE.0.0) THEN HITS(16) = 57.29577*ACOS(-PD(1)/PXY) ! x axis flipped ELSE HITS(16) = 0.0 ENDIF IF(PD(3).LT.0.0) HITS(16) = 360.0 - HITS(16) CALL GSAHIT(ISET,IDET,ITRA,NUMBV,HITS,IHIT) ENDIF RETURN END