radgen_event.f

Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file radgen_event.f

      subroutine radgen_event

      implicit none

       include "mc_set.inc"
       include "mconsp.inc"
       include "phiout.inc"
       include "tailcom.inc"
       include "cmpcom.inc"
       include "radgen.inc"
       include "mcRadCor.inc"

      integer iphot
      real phrad(4),q2true,nutrue,radweight
      real nu,q2,phi,yys,xxs

! calculate radiative corrections
      nu=sngl(gennu)
      q2=sngl(genq2)
      phi=sngl(genphi)
      yys=sngl(geny)
      xxs=sngl(genx)
      call radgen(mcset_enebeam,q2,nu,yys,xxs,phi,phrad,q2true,nutrue,
     +     radweight)

! fill mcRadCor WCB with ...

      mcradcor_id=1

! ... true kinematics

* by definition we calculate xbj using the proton mass
* ---> xbj for elastic events is A
*     mcRadCor_XTrue=q2true/(2.0d0*amp*nutrue)
      mcradcor_xtrue=q2true/(2.0d0*0.938272d0*nutrue)
      mcradcor_ytrue=nutrue/mcset_enebeam
      mcradcor_nutrue=nutrue
      mcradcor_q2true=q2true
      mcradcor_w2true=amp2 - q2true + 2.*amp*nutrue
*     print *,amp,amp2,mcRadCor_XTrue,nutrue,q2true,mcRadCor_W2True

! ... kinematics of real photon

      mcradcor_ebrems=phrad(4)
      mcradcor_thetabrems=0.
      if(phrad(4).gt.0.)
     +    mcradcor_thetabrems = acos(phrad(3)/phrad(4))
      mcradcor_phibrems=0.
      if (.not.(phrad(1).eq.0..and.phrad(2).eq.0.)) then
        mcradcor_phibrems = atan2(phrad(2),phrad(1))
        if (mcradcor_phibrems.lt.0.)
     +       mcradcor_phibrems = mcradcor_phibrems + twopi
      endif

c...if we would like to have the TSAI system angles for the real gamma than
c     mcRadCor_ThetaBrems = dthg
c     mcRadCor_PhiBrems = dphig

! ... radiative contributions

      mcradcor_sigrad=sigrad
      mcradcor_sigcor=sigcor
      mcradcor_sigcorerr=0.
      mcradcor_tailine=tine
      mcradcor_tailela=tpro
      mcradcor_tailcoh=tnuc
      mcradcor_vacuum=vac
      mcradcor_vertex=vertex
      mcradcor_small=small
      mcradcor_redfac=redfac

! ... radiative correction type

      if (ita.eq.2) then
        mcradcor_ctype='elas'
      else if (ita.eq.3) then
        mcradcor_ctype='qela'
      else
        mcradcor_ctype='inel'
      endif

      end