pysave.f

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

 
C*********************************************************************
 
C...PYSAVE
C...Saves and restores parameter and cross section values for the
C...3 gamma-p and 6 (or 4, or 9, or 13) gamma-gamma alternatives.
C...Also makes random choice between alternatives.
 
      SUBROUTINE pysave(ISAVE,IGA)
 
C...Double precision and integer declarations.
      IMPLICIT DOUBLE PRECISION(a-h, o-z)
      IMPLICIT INTEGER(i-n)
      INTEGER pyk,pychge,pycomp
C...Commonblocks.
      common/pysubs/msel,mselpd,msub(500),kfin(2,-40:40),ckin(200)
      common/pypars/mstp(200),parp(200),msti(200),pari(200)
      common/pyint1/mint(400),vint(400)
      common/pyint2/iset(500),kfpr(500,2),coef(500,20),icol(40,4,2)
      common/pyint5/ngenpd,ngen(0:500,3),xsec(0:500,3)
      common/pyint7/sigt(0:6,0:6,0:5)
      SAVE /pysubs/,/pypars/,/pyint1/,/pyint2/,/pyint5/,/pyint7/
C...Local arrays and saved variables.
      dimension ncp(15),nsubcp(15,20),msubcp(15,20),coefcp(15,20,20),
     &ngencp(15,0:20,3),xseccp(15,0:20,3),sigtcp(15,0:6,0:6,0:5),
     &intcp(15,20),recp(15,20)
      SAVE ncp,nsubcp,msubcp,coefcp,ngencp,xseccp,sigtcp,intcp,recp
 
C...Save list of subprocesses and cross-section information.
      IF(isave.EQ.1) THEN
        icp=0
        DO 120 i=1,500
          IF(msub(i).EQ.0.AND.i.NE.96.AND.i.NE.97) goto 120
          icp=icp+1
          nsubcp(iga,icp)=i
          msubcp(iga,icp)=msub(i)
          DO 100 j=1,20
            coefcp(iga,icp,j)=coef(i,j)
  100     CONTINUE
          DO 110 j=1,3
            ngencp(iga,icp,j)=ngen(i,j)
            xseccp(iga,icp,j)=xsec(i,j)
  110     CONTINUE
  120   CONTINUE
        ncp(iga)=icp
        DO 130 j=1,3
          ngencp(iga,0,j)=ngen(0,j)
          xseccp(iga,0,j)=xsec(0,j)
  130   CONTINUE
        DO 160 i1=0,6
          DO 150 i2=0,6
            DO 140 j=0,5
              sigtcp(iga,i1,i2,j)=sigt(i1,i2,j)
  140       CONTINUE
  150     CONTINUE
  160   CONTINUE
 
C...Save various common process variables.
        DO 170 j=1,10
          intcp(iga,j)=mint(40+j)
  170   CONTINUE
        intcp(iga,11)=mint(101)
        intcp(iga,12)=mint(102)
        intcp(iga,13)=mint(107)
        intcp(iga,14)=mint(108)
        intcp(iga,15)=mint(123)
        recp(iga,1)=ckin(3)
        recp(iga,2)=vint(318)
 
C...Save cross-section information only.
      ELSEIF(isave.EQ.2) THEN
        DO 190 icp=1,ncp(iga)
          i=nsubcp(iga,icp)
          DO 180 j=1,3
            ngencp(iga,icp,j)=ngen(i,j)
            xseccp(iga,icp,j)=xsec(i,j)
  180     CONTINUE
  190   CONTINUE
        DO 200 j=1,3
          ngencp(iga,0,j)=ngen(0,j)
          xseccp(iga,0,j)=xsec(0,j)
  200   CONTINUE
 
C...Choose between allowed alternatives.
      ELSEIF(isave.EQ.3.OR.isave.EQ.4) THEN
        IF(isave.EQ.4) THEN
          xsumcp=0d0
          DO 210 ig=1,mint(121)
            xsumcp=xsumcp+xseccp(ig,0,1)
  210     CONTINUE
          xsumcp=xsumcp*pyr(0)
          DO 220 ig=1,mint(121)
            iga=ig
            xsumcp=xsumcp-xseccp(ig,0,1)
            IF(xsumcp.LE.0d0) goto 230
  220     CONTINUE
  230     CONTINUE
        ENDIF
 
C...Restore cross-section information.
        DO 240 i=1,500
          msub(i)=0
  240   CONTINUE
        DO 270 icp=1,ncp(iga)
          i=nsubcp(iga,icp)
          msub(i)=msubcp(iga,icp)
          DO 250 j=1,20
            coef(i,j)=coefcp(iga,icp,j)
  250     CONTINUE
          DO 260 j=1,3
            ngen(i,j)=ngencp(iga,icp,j)
            xsec(i,j)=xseccp(iga,icp,j)
  260     CONTINUE
  270   CONTINUE
        DO 280 j=1,3
          ngen(0,j)=ngencp(iga,0,j)
          xsec(0,j)=xseccp(iga,0,j)
  280   CONTINUE
        DO 310 i1=0,6
          DO 300 i2=0,6
            DO 290 j=0,5
              sigt(i1,i2,j)=sigtcp(iga,i1,i2,j)
  290       CONTINUE
  300     CONTINUE
  310   CONTINUE
 
C...Restore various common process variables.
        DO 320 j=1,10
          mint(40+j)=intcp(iga,j)
  320   CONTINUE
        mint(101)=intcp(iga,11)
        mint(102)=intcp(iga,12)
        mint(107)=intcp(iga,13)
        mint(108)=intcp(iga,14)
        mint(123)=intcp(iga,15)
        ckin(3)=recp(iga,1)
        ckin(1)=2d0*ckin(3)
        vint(318)=recp(iga,2)
 
C...Sum up cross-section info (for PYSTAT).
      ELSEIF(isave.EQ.5) THEN
        DO 330 i=1,500
          msub(i)=0
          ngen(i,1)=0
          ngen(i,3)=0
          xsec(i,3)=0d0
  330   CONTINUE
        ngen(0,1)=0
        ngen(0,2)=0
        ngen(0,3)=0
        xsec(0,3)=0
        DO 350 ig=1,mint(121)
          DO 340 icp=1,ncp(ig)
            i=nsubcp(ig,icp)
            IF(msubcp(ig,icp).EQ.1) msub(i)=1
            ngen(i,1)=ngen(i,1)+ngencp(ig,icp,1)
            ngen(i,3)=ngen(i,3)+ngencp(ig,icp,3)
            xsec(i,3)=xsec(i,3)+xseccp(ig,icp,3)
  340     CONTINUE
          ngen(0,1)=ngen(0,1)+ngencp(ig,0,1)
          ngen(0,2)=ngen(0,2)+ngencp(ig,0,2)
          ngen(0,3)=ngen(0,3)+ngencp(ig,0,3)
          xsec(0,3)=xsec(0,3)+xseccp(ig,0,3)
  350   CONTINUE
      ENDIF
 
      RETURN
      END