pyrvgw.f

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

 
C*********************************************************************
 
C...PYRVGW
C...Generalized Matrix Element for R-Violating 3-body widths.
C...P. Z. Skands
      SUBROUTINE pyrvgw(KFIN,ID1,ID2,ID3,XLAM)
 
      IMPLICIT DOUBLE PRECISION (a-h,o-z)
      IMPLICIT integer(i-n)
      parameter(ksusy1=1000000,ksusy2=2000000,ktechn=3000000,
     &kexcit=4000000,kdimen=5000000)
      parameter(eps=1d-4)
      common/pydat2/kchg(500,4),pmas(500,4),parf(2000),vckm(4,4)
      common/pyrvnv/ab(2,16,2),rms(0:3),res(6,2),intres(6,3),idr,idr2
     &     ,dcmass,kfr(3)
      common/pyssmt/zmix(4,4),umix(2,2),vmix(2,2),smz(4),smw(2),
     & sfmix(16,4),zmixi(4,4),umixi(2,2),vmixi(2,2)
      DOUBLE PRECISION xlim(3,3)
      INTEGER kc(0:3), pycomp
      LOGICAL dcmass, dcheck(6)
      SAVE /pydat2/,/pyrvnv/,/pyssmt/
 
      xlam   = 0d0
 
      kc(0)  = pycomp(kfin)
      kc(1)  = pycomp(id1)
      kc(2)  = pycomp(id2)
      kc(3)  = pycomp(id3)
      rms(0) = pmas(kc(0),1)
      rms(1) = pymrun(id1,pmas(kc(1),1)**2)
      rms(2) = pymrun(id2,pmas(kc(2),1)**2)
      rms(3) = pymrun(id3,pmas(kc(3),1)**2)
C...INITIALIZE OUTER INTEGRATION LIMITS AND KINEMATICS CHECK
      xlim(1,1)=(rms(1)+rms(2))**2
      xlim(1,2)=(rms(0)-rms(3))**2
      xlim(1,3)=xlim(1,2)-xlim(1,1)
      xlim(2,1)=(rms(2)+rms(3))**2
      xlim(2,2)=(rms(0)-rms(1))**2
      xlim(2,3)=xlim(2,2)-xlim(2,1)
      xlim(3,1)=(rms(1)+rms(3))**2
      xlim(3,2)=(rms(0)-rms(2))**2
      xlim(3,3)=xlim(3,2)-xlim(3,1)
C...Check Phase Space
      IF (xlim(1,3).LT.0d0.OR.xlim(2,3).LT.0d0.OR.xlim(3,3).LT.0d0) THEN
        RETURN
      ENDIF
 
C...INITIALIZE RESONANCE INFORMATION
      DO 110 jres = 1,3
        DO 100 imass = 1,2
          ires = 2*(jres-1)+imass
          intres(ires,1) = 0
          dcheck(ires)   =.false.
C...NO RIGHT-HANDED NEUTRINOS
          IF (((imass.EQ.2).AND.((iabs(kfr(jres)).EQ.12).or
     &         .(iabs(kfr(jres)).EQ.14).OR.(iabs(kfr(jres)).EQ.16))).or
     &         .kfr(jres).EQ.0) goto 100
          res(ires,1) = pmas(pycomp(imass*ksusy1+iabs(kfr(jres))),1)
          res(ires,2) = pmas(pycomp(imass*ksusy1+iabs(kfr(jres))),2)
          intres(ires,1) = iabs(kfr(jres))
          intres(ires,2) = imass
          IF (kfr(jres).LT.0) intres(ires,3) = 1
          IF (kfr(jres).GT.0) intres(ires,3) = 0
  100   CONTINUE
  110 CONTINUE
 
C...SUM OVER DIAGRAMS AND INTEGRATE OVER PHASE SPACE
 
C...RESONANCE CONTRIBUTIONS
C...(Only sum contributions where the resonance is off shell).
C...Store whether diagram on/off in DCHECK.
C...LOOP OVER MASS STATES
      DO 120 j=1,2
        idr=j
        tmix = sfmix(intres(idr,1),2*j+intres(idr,3)-1)**2
        IF ((rms(0).LT.(rms(1)+res(idr,1)).OR.(res(idr,1).LT.(rms(2)
     &       +rms(3)))).AND.tmix.GT.eps.AND.intres(idr,1).NE.0) THEN
          dcheck(idr) =.true.
          xlam = xlam + tmix * pyrvi1(2,3,1)
        ENDIF
 
        idr=j+2
        tmix = sfmix(intres(idr,1),2*j+intres(idr,3)-1)**2
        IF ((rms(0).LT.(rms(2)+res(idr,1)).OR.(res(idr,1).LT.(rms(1)
     &       +rms(3)))).AND.tmix.GT.eps.AND.intres(idr,1).NE.0) THEN
          dcheck(idr) =.true.
          xlam = xlam + tmix * pyrvi1(1,3,2)
        ENDIF
 
        idr=j+4
        tmix = sfmix(intres(idr,1),2*j+intres(idr,3)-1)**2
        IF ((rms(0).LT.(rms(3)+res(idr,1)).OR.(res(idr,1).LT.(rms(1)
     &       +rms(2)))).AND.tmix.GT.eps.AND.intres(idr,1).NE.0) THEN
          dcheck(idr) =.true.
          xlam = xlam + tmix * pyrvi1(1,2,3)
        ENDIF
  120 CONTINUE
C... L-R INTERFERENCES
C... (Only add contributions where both contributing diagrams
C... are non-resonant).
      idr=1
      IF (dcheck(1).AND.dcheck(2)) THEN
C...Bug corrected 11/12 2001. Skands.
        xlam  = xlam + 2d0 * pyrvi2(2,3,1)
     &     * sfmix(intres(1,1),2+intres(1,3)-1)
     &     * sfmix(intres(2,1),4+intres(2,3)-1)
      ENDIF
 
      idr=3
      IF (dcheck(3).AND.dcheck(4)) THEN
        xlam  = xlam + 2d0 * pyrvi2(1,3,2)
     &     * sfmix(intres(3,1),2+intres(3,3)-1)
     &     * sfmix(intres(4,1),4+intres(4,3)-1)
      ENDIF
 
      idr=5
      IF (dcheck(5).AND.dcheck(6)) THEN
        xlam  = xlam + 2d0 * pyrvi2(1,2,3)
     &     * sfmix(intres(5,1),2+intres(5,3)-1)
     &     * sfmix(intres(6,1),4+intres(6,3)-1)
      ENDIF
C... TRUE INTERFERENCES
C... (Only add contributions where both contributing diagrams
C... are non-resonant).
      pref=-2d0
      IF ((kfin-ksusy1).EQ.24.OR.(kfin-ksusy1).EQ.37) pref=2d0
      DO 140 ikr1 = 1,2
        DO 130 ikr2 = 1,2
          idr  = ikr1+2
          idr2 = ikr2
          IF (dcheck(idr).AND.dcheck(idr2)) THEN
            xlam = xlam + pref*pyrvi3(1,3,2) *
     &           sfmix(intres(idr,1),2*ikr1+intres(idr,3)-1)
     &           *sfmix(intres(idr2,1),2*ikr2+intres(idr2,3)-1)
          ENDIF
 
          idr  = ikr1+4
          idr2 = ikr2
          IF (dcheck(idr).AND.dcheck(idr2)) THEN
            xlam = xlam + pref*pyrvi3(1,2,3) *
     &           sfmix(intres(idr,1),2*ikr1+intres(idr,3)-1)
     &           *sfmix(intres(idr2,1),2*ikr2+intres(idr2,3)-1)
          ENDIF
 
          idr  = ikr1+4
          idr2 = ikr2+2
          IF (dcheck(idr).AND.dcheck(idr2)) THEN
            xlam = xlam + pref*pyrvi3(2,1,3) *
     &           sfmix(intres(idr,1),2*ikr1+intres(idr,3)-1)
     &           *sfmix(intres(idr2,1),2*ikr2+intres(idr2,3)-1)
          ENDIF
  130   CONTINUE
  140 CONTINUE
 
      RETURN
      END