*$ CREATE USREOU.FOR *COPY USREOU * *=== Usreou ===========================================================* * SUBROUTINE USREOU INCLUDE '(DBLPRC)' INCLUDE '(DIMPAR)' INCLUDE '(IOUNIT)' * *----------------------------------------------------------------------* * * * Copyright (C) 1991-2005 by Alfredo Ferrari & Paola Sala * * All Rights Reserved. * * * * * * USeR Event OUtput: this routine is called at the end of each * * event * * * * Created on 01 january 1991 by Alfredo Ferrari & Paola Sala * * Infn - Milan * * * * Last change on 09-apr-99 by Alfredo Ferrari * * * * * *----------------------------------------------------------------------* * INCLUDE '(BLNKCM)' INCLUDE '(CASLIM)' INCLUDE '(FLKMAT)' INCLUDE '(SCEXFL)' INCLUDE '(PAPROP)' * logical first integer edepos(1000), ephot(200), nphot(10) real*4 xtup(10,500) real*4 regengy real*4 energy, esmall, elarge, etot real*4 egamma(10), partmom dimension iscmem(4) dimension idepos(600),endepos(600) common /pextra/ xtup, nevt, ntpevent, iicas common /deposit/ store(5000), regengy(5000), edepos, ephot, nphot common /csihits/ xgamhit(40),ygamhit(40),egamhit(40),agamhit(40) & ,pxgamhit(40), pygamhit(40), pzgamhit(40) & ,xparhit(40),yparhit(40),pparhit(40),eparhit(40) & ,aparhit(40), idhit(40), ngam, nothr data first/.true./ save * if (first) then store = 0.0d0 edepos = 0 ephot = 0 nphot = 0 first = .false. endif * * Main entry point * * Extract the energy deposits in the detector regions do ir=1,Nregs kk = icmbgn + 2*nregs + ir regengy(ir) = 1000.0d0*(comsco(kk)-store(ir)) enddo * if (iicas.eq.0) go to 190 * c if (iicas.gt.1) write(49,*) ncase,iicas c if (ncase.eq.88319) then c do ii=1,iicas c idpx = xtup(2,ii) + 0.01 c write(49,10) ncase,idpx,(xtup(ik,ii),ik=3,10) c 10 format(i8,i4,1pe13.5,3(0pf12.6),3(1pe14.5),0pf8.0) c enddo c endif * * Test for hits in the CsI and store information ngam = 0 nothr = 0 etphot = 0.0 etpart = 0.0 do ii=1,iicas idpx = xtup(2,ii) + 0.01 newreg = xtup(10,ii) + 0.01 if (newreg.eq.3) then if (ngam.lt.40) then partmom = 1000.0*sqrt( & xtup(7,ii)**2 + xtup(8,ii)**2 + xtup(9,ii)**2) if (idpx.eq.7) then if (partmom.gt.0.50) then ngam = ngam + 1 xgamhit(ngam) = xtup(4,ii) ygamhit(ngam) = xtup(5,ii) pxgamhit(ngam) = 1000.0*xtup(7,ii) pygamhit(ngam) = 1000.0*xtup(8,ii) pzgamhit(ngam) = 1000.0*xtup(9,ii) egamhit(ngam) = partmom agamhit(ngam) = xtup(3,ii) etphot = etphot + partmom endif else nothr = nothr + 1 eparhit(nothr) = sqrt(partmom*partmom & + am(idpx)*am(idpx)*1.0d6) - am(idpx)*1.0d3 etpart = etpart + eparhit(nothr) xparhit(nothr) = xtup(4,ii) yparhit(nothr) = xtup(5,ii) aparhit(nothr) = xtup(3,ii) idhit(nothr) = idpx endif endif endif enddo etcsi = etphot + etpart c if (ngam.gt.0) write(49,*) ngam,etphot,etcsi * if ((ngam.gt.0).and.(ngam.le.10)) then nphot(ngam) = nphot(ngam) + 1 do i=1,ngam ie = egamhit(i)/10.0 + 1.0 if (ie.le.200) then ephot(ie) = ephot(ie) + 1 endif enddo endif * if (ngam.eq.2) then gammom2 = (pxgamhit(1)+pxgamhit(2))**2 & + (pygamhit(1)+pygamhit(2))**2 & + (pzgamhit(1)+pzgamhit(2))**2 gamengy2 = (egamhit(1)+egamhit(2) )**2 gammass2 = gamengy2 - gammom2 if (gammass2.gt.0.0) then gammass = sqrt(gammass2) else gammass = 0.0 endif im = gammass/5.0 + 1.0 if (im.le.1000) edepos(im) = edepos(im) + 1 * write(49,*) ncase,sngl(gammom2),sngl(gamengy2),sngl(gammass) if (gammass.gt.400.0) write(49,*) ncase,sngl(gammass) endif * 190 continue * Save the previous total energy deposits iicas = 0 do ir=1,Nregs kk = icmbgn + 2*nregs + ir store(ir) = comsco(kk) enddo RETURN *=== End of subroutine Usreou =========================================* END