c********************************************************************** subroutine showdigi c c c designed by J. Otwinowski c modified last on 05/12/2008 by I. Koenig c********************************************************************* implicit none #include "geant321/gcflag.inc" #include "showtups.inc" #include "user.inc" integer nhits, nhdim, nhmax, nvdim, ipl, i, nh parameter (nhdim=14, nhmax=1500, nvdim=3) character*4 iudet integer itrack(nhmax), numvs(nvdim), numbv(nvdim,nhmax) real hits(nhdim,nhmax) integer isect, cur_trk real x_in, y_in, z_in, theta_in, phi_in, pxy real p,E,beta do isect=1,nsect ! loop over sectors ntrk = 0 numvs(1) = isect numvs(2) = 0 numvs(3) = 0 do ipl=1,3 ! loop over planes iudet = 'S'//char(48+ipl)//'AI' call gfhits('SHOW',iudet,nvdim,nhdim,nhmax,0,numvs, + itrack,numbv,hits,nhits) IF (NHITS.GT.NHMAX) THEN ! overflow WRITE(6,*) 'ERROR: overflow in showdigi.F (nHits>', & NHMAX,') in Shower',ISECT,IPL RETURN ENDIF cur_trk = 0 do nh=1,nhits ! loop over hits if(hits(8,nh).eq.0.0 .or. nint(hits(11,nh)).eq.1) then ! EDEP=0 or new volume entered cur_trk = hits(7,nh) ! save current track x_in = hits(1,nh) ! save position of entrance y_in = hits(2,nh) z_in = hits(3,nh) theta_in = 57.29577*acos(hits(14,nh)) ! angle of incidence pxy = sqrt(hits(12,nh)**2+hits(13,nh)**2) ! TRAP shape OK p = hits(10,nh) ! save total momentum E = hits(4,nh) ! save total energy beta = p/E ! save beta of particle if(pxy.ne.0.0) then ! azimuthal angle phi_in = 57.29577*acos(hits(12,nh)/pxy) ! TRAP shape OK else phi_in = 0.0 endif if(hits(13,nh).lt.0.0) phi_in = 360.0 - phi_in else if (cur_trk.eq.hits(7,nh)) then ! still same track on ! exit of volume ? if(ntrk.lt.MAXTRKSHOW) then ntrk = ntrk + 1 ! nb. of tracks showtrk(ntrk) = hits(7,nh) ! track number showdet(ntrk) = hits(6,nh) ! det. plane number showe(ntrk) = hits(8,nh)*1.0e3 ! energy loss (in MeV) showx(ntrk) = lunit*x_in ! x in mm (or cm) showy(ntrk) = lunit*y_in ! y (TRAP shape !) c showy(ntrk) = lunit*z_in ! y (TRD1 shape !) showtheta(ntrk) = theta_in ! angle of incidence showphi(ntrk) = phi_in ! azimuthal angle showbeta(ntrk) = beta ! beta = v/c = p/E else WRITE(6,*) 'ERROR: too many hits (>', & MAXTRKSHOW,') in Shower sector',ISECT RETURN endif endif enddo ! end loop over hits enddo ! end loop over planes if(iroot.eq.0) then call hfnt(4000+isect) else call fillshow(isect) endif enddo ! end loop over sectors return end