<http://savannah.cern.ch/bugs/?102190>
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 7 Aug 2013 08:32:54 +0000 (08:32 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 7 Aug 2013 08:32:54 +0000 (08:32 +0000)
[bug #102190] POWHEG interface to AliRoot]
Rachid Guernane

12 files changed:
PYTHIA6/AliGenPythia.cxx
PYTHIA6/AliGenPythia.h
PYTHIA6/AliGenPythiaPlus.cxx
PYTHIA6/AliPythia.cxx
PYTHIA6/AliPythia6.cxx
PYTHIA6/CMakelibpythia6.4.25.pkg
PYTHIA6/POWHEG-BOX/LesHouches.h [new file with mode: 0644]
PYTHIA6/POWHEG-BOX/hepevt.h [new file with mode: 0644]
PYTHIA6/POWHEG-BOX/lhefread.f [new file with mode: 0644]
PYTHIA6/POWHEG-BOX/setup-PYTHIA-lhef.f [new file with mode: 0644]
PYTHIA6/PythiaProcesses.h
PYTHIA6/pythia6.4.25/pythia-6.4.25.f

index 5872c40..e490e6f 100644 (file)
@@ -81,6 +81,8 @@ AliGenPythia::AliGenPythia():
     fPatchOmegaDalitz(0), 
     fNpartons(0),
     fReadFromFile(0),
+    fReadLHEF(0),
+    fPartonShower(1),
     fQuench(0),
     fQhat(0.),
     fLength(0.),
@@ -131,6 +133,7 @@ AliGenPythia::AliGenPythia():
     fHeader(0),  
     fRL(0),      
     fkFileName(0),
+    fkNameLHEF(0),
     fFragPhotonInCalo(kFALSE),
     fHadronInCalo(kFALSE) ,
     fPi0InCalo(kFALSE) ,
@@ -193,6 +196,8 @@ AliGenPythia::AliGenPythia(Int_t npart)
      fPatchOmegaDalitz(0), 
      fNpartons(0),
      fReadFromFile(kFALSE),
+     fReadLHEF(0),
+     fPartonShower(1),
      fQuench(kFALSE),
      fQhat(0.),
      fLength(0.),
@@ -243,6 +248,7 @@ AliGenPythia::AliGenPythia(Int_t npart)
      fHeader(0),  
      fRL(0),      
      fkFileName(0),
+     fkNameLHEF(0),
      fFragPhotonInCalo(kFALSE),
      fHadronInCalo(kFALSE) ,
      fPi0InCalo(kFALSE) ,
@@ -403,7 +409,8 @@ void AliGenPythia::Init()
        fPythia->SetMSTP(91,0);
     }
 
-
+          if (fReadLHEF) fPythia->OpenFortranFile(97, const_cast<char*>(fkNameLHEF));
+       
     if (fReadFromFile) {
        fRL  =  AliRunLoader::Open(fkFileName, "Partons");
        fRL->LoadKinematics();
@@ -441,6 +448,7 @@ void AliGenPythia::Init()
     case kPyCharmpPbMNR:
     case kPyCharmppMNR:
     case kPyCharmppMNRwmi:
+    case kPyCharmPWHG:
        fParentSelect[0] =   411;
        fParentSelect[1] =   421;
        fParentSelect[2] =   431;
@@ -478,6 +486,7 @@ void AliGenPythia::Init()
     case kPyBeautypPbMNR:
     case kPyBeautyppMNR:
     case kPyBeautyppMNRwmi:
+    case kPyBeautyPWHG:
        fParentSelect[0]=  511;
        fParentSelect[1]=  521;
        fParentSelect[2]=  531;
@@ -508,6 +517,7 @@ void AliGenPythia::Init()
     case kPyMbNonDiffr:
     case kPyMbMSEL1:
     case kPyJets:
+    case kPyJetsPWHG:
     case kPyDirectGamma:
     case kPyLhwgMb:    
        break;
@@ -720,7 +730,10 @@ void AliGenPythia::Generate()
            fProcess != kPyZ &&
            fProcess != kPyCharmppMNRwmi && 
            fProcess != kPyBeautyppMNRwmi &&
-           fProcess != kPyBeautyJets) {
+           fProcess != kPyBeautyJets &&
+     fProcess != kPyJetsPWHG &&
+     fProcess != kPyCharmPWHG &&
+     fProcess != kPyBeautyPWHG) {
            
            for (i = 0; i < np; i++) {
                TParticle* iparticle = (TParticle *) fParticles.At(i);
@@ -949,7 +962,7 @@ Int_t  AliGenPythia::GenerateMB()
 
     Int_t* pParent = new Int_t[np];
     for (i=0; i< np; i++) pParent[i] = -1;
-    if ((fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi) 
+    if ((fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi || fProcess == kPyJetsPWHG || fProcess == kPyCharmPWHG || fProcess == kPyBeautyPWHG)
        && fEtMinJet > 0.) {
        TParticle* jet1 = (TParticle *) fParticles.At(6);
        TParticle* jet2 = (TParticle *) fParticles.At(7);
@@ -1111,7 +1124,7 @@ Int_t  AliGenPythia::GenerateMB()
        Int_t km = iparticle->GetFirstMother();
        if (
            (((ks == 1  && kf!=0 && KinematicSelection(iparticle, 0)) || (ks !=1)) && ((fStackFillOpt != kHeavyFlavor) || IsFromHeavyFlavor(i))) ||
-           ((fProcess == kPyJets || fProcess == kPyBeautyJets) && ks == 21 && km == 0 && i>1)
+           ((fProcess == kPyJets || fProcess == kPyBeautyJets || fProcess == kPyJetsPWHG || fProcess == kPyCharmPWHG || fProcess == kPyBeautyPWHG) && ks == 21 && km == 0 && i>1)
            ) 
          {
             nc++;
@@ -1229,7 +1242,7 @@ void AliGenPythia::MakeHeader()
 // Jets that have triggered
 
     //Need to store jets for b-jet studies too!
-    if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi)
+    if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyBeautyJets || fProcess == kPyBeautyppMNRwmi || fProcess == kPyJetsPWHG || fProcess == kPyCharmPWHG || fProcess == kPyBeautyPWHG)
     {
        Int_t ntrig, njet;
        Float_t jets[4][10];
@@ -1295,6 +1308,11 @@ void AliGenPythia::MakeHeader()
 //
 // Store pt^hard 
     ((AliGenPythiaEventHeader*) fHeader)->SetPtHard(fPythia->GetVINT(47));
+               
+//
+// Store Event Weight
+               ((AliGenPythiaEventHeader*) fHeader)->SetEventWeight(fPythia->GetPARI(7));
+                               
 //
 //  Pass header
 //
@@ -1317,7 +1335,7 @@ Bool_t AliGenPythia::CheckTrigger(const TParticle* jet1, const TParticle* jet2)
     pdg[1] = jet2->GetPdgCode();    
     Bool_t   triggered = kFALSE;
 
-    if (fProcess == kPyJets || fProcess == kPyBeautyJets || fProcess ==  kPyBeautyppMNRwmi) {
+    if (fProcess == kPyJets || fProcess == kPyBeautyJets || fProcess ==  kPyBeautyppMNRwmi || fProcess == kPyJetsPWHG || fProcess == kPyCharmPWHG || fProcess == kPyBeautyPWHG) {
        Int_t njets = 0;
        Int_t ntrig = 0;
        Float_t jets[4][10];
index 7fc76ad..1f326fb 100644 (file)
@@ -198,6 +198,7 @@ class AliGenPythia : public AliGenMC
     virtual void SetHadronisation(Int_t flag = 1) {fHadronisation = flag;}
     virtual void SetPatchOmegaDalitz(Int_t flag = 1) {fPatchOmegaDalitz = flag;}
     virtual void SetReadFromFile(const Text_t *filname) {fkFileName = filname;  fReadFromFile = 1;}    
+    virtual void SetReadLHEF(const Text_t *filename) {fkNameLHEF = filename; fReadLHEF = 1;}    
 
     //
     // Pile-up
@@ -281,6 +282,7 @@ class AliGenPythia : public AliGenMC
     Bool_t      fPatchOmegaDalitz;  //flag for omega dalitz decay patch
     Int_t       fNpartons;          //Number of partons before hadronisation
     Int_t       fReadFromFile;      //read partons from file
+    Int_t       fReadLHEF;          //read lhef file
     Int_t       fQuench;            //Flag for quenching
     Float_t     fQhat;              //Transport coefficient (GeV^2/fm)
     Float_t     fLength;            //Medium length (fm)
@@ -342,7 +344,8 @@ class AliGenPythia : public AliGenMC
     AliGenPythiaEventHeader* fHeader;  //! Event header
     AliRunLoader*            fRL;      //! Run Loader
     const Text_t* fkFileName;          //! Name of file to read from
-
+    const Text_t* fkNameLHEF;          //! Name of lhef file to read from
+    Int_t fPartonShower;
 
     Bool_t fFragPhotonInCalo; // Option to ask for Fragmentation Photon in calorimeters acceptance
     Bool_t fHadronInCalo;     // Option to ask for hadron (not pi0) in calorimeters acceptance
index 6131d5f..fc69250 100644 (file)
@@ -374,6 +374,7 @@ void AliGenPythiaPlus::Init()
     case kPyCharmpPbMNR:
     case kPyCharmppMNR:
     case kPyCharmppMNRwmi:
+    case kPyCharmPWHG:
        fParentSelect[0] =   411;
        fParentSelect[1] =   421;
        fParentSelect[2] =   431;
@@ -411,6 +412,7 @@ void AliGenPythiaPlus::Init()
     case kPyBeautypPbMNR:
     case kPyBeautyppMNR:
     case kPyBeautyppMNRwmi:
+    case kPyBeautyPWHG:
        fParentSelect[0]=  511;
        fParentSelect[1]=  521;
        fParentSelect[2]=  531;
@@ -441,6 +443,7 @@ void AliGenPythiaPlus::Init()
     case kPyMbNonDiffr:
     case kPyMbMSEL1:
     case kPyJets:
+    case kPyJetsPWHG:
     case kPyDirectGamma:
     case kPyLhwgMb:    
        break;
@@ -604,7 +607,10 @@ void AliGenPythiaPlus::Generate()
            fProcess != kPyW && 
            fProcess != kPyZ &&
            fProcess != kPyCharmppMNRwmi && 
-           fProcess != kPyBeautyppMNRwmi) {
+           fProcess != kPyBeautyppMNRwmi &&
+           fProcess != kPyJetsPWHG &&
+           fProcess != kPyCharmPWHG &&
+     fProcess != kPyBeautyPWHG) {
            
            for (i = 0; i < np; i++) {
                TParticle* iparticle = (TParticle *) fParticles.At(i);
@@ -827,7 +833,7 @@ Int_t  AliGenPythiaPlus::GenerateMB()
     
     Int_t* pParent = new Int_t[np];
     for (i=0; i< np; i++) pParent[i] = -1;
-    if (fProcess == kPyJets || fProcess == kPyDirectGamma) {
+    if (fProcess == kPyJets || fProcess == kPyDirectGamma || fProcess == kPyJetsPWHG || fProcess == kPyCharmPWHG || fProcess == kPyBeautyPWHG ) {
        TParticle* jet1 = (TParticle *) fParticles.At(6);
        TParticle* jet2 = (TParticle *) fParticles.At(7);
        if (!CheckTrigger(jet1, jet2)) {
@@ -837,7 +843,7 @@ Int_t  AliGenPythiaPlus::GenerateMB()
     }
 
     // Select jets with fragmentation photon or pi0 going to PHOS or EMCAL
-    if (fProcess == kPyJets && (fFragPhotonInCalo || fPi0InCalo) ) {
+    if ((fProcess == kPyJets || fProcess == kPyJetsPWHG) && (fFragPhotonInCalo || fPi0InCalo) ) {
 
       Bool_t ok = kFALSE;
 
@@ -870,7 +876,7 @@ Int_t  AliGenPythiaPlus::GenerateMB()
     
     
      // Select events with a photon  pt > min pt going to PHOS eta acceptance or exactly PHOS eta phi
-    if ((fProcess == kPyJets || fProcess == kPyDirectGamma) && fPhotonInCalo && (fCheckPHOSeta || fCheckPHOS)){
+    if ((fProcess == kPyJets || fProcess == kPyJetsPWHG || fProcess == kPyDirectGamma) && fPhotonInCalo && (fCheckPHOSeta || fCheckPHOS)){
 
       Bool_t okd = kFALSE;
 
@@ -990,7 +996,7 @@ Int_t  AliGenPythiaPlus::GenerateMB()
        Int_t km = iparticle->GetFirstMother();
        if ((ks == 1  && kf!=0 && KinematicSelection(iparticle, 0)) ||
            (ks != 1) ||
-           (fProcess == kPyJets && ks == 21 && km == 0 && i>1)) {
+           ((fProcess == kPyJets || fProcess == kPyJetsPWHG) && ks == 21 && km == 0 && i>1)) {
            nc++;
            if (ks == 1) trackIt = 1;
 
@@ -1116,7 +1122,7 @@ void AliGenPythiaPlus::MakeHeader()
 //
 // Jets that have triggered
 
-    if (fProcess == kPyJets)
+    if (fProcess == kPyJets || fProcess == kPyJetsPWHG)
     {
        Int_t ntrig, njet;
        Float_t jets[4][10];
@@ -1193,7 +1199,7 @@ Bool_t AliGenPythiaPlus::CheckTrigger(const TParticle* jet1, const TParticle* je
     pdg[1] = jet2->GetPdgCode();    
     Bool_t   triggered = kFALSE;
 
-    if (fProcess == kPyJets) {
+    if (fProcess == kPyJets || fProcess == kPyJetsPWHG) {
        Int_t njets = 0;
        Int_t ntrig = 0;
        Float_t jets[4][10];
index f8db4a6..918e67c 100644 (file)
@@ -616,6 +616,28 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
       case kPyMBRDoubleDiffraction:
       case kPyMBRCentralDiffraction:
       break;
+      case kPyJetsPWHG:
+      //    N.B.
+      //    ====
+      //    For the case of jet production the following parameter setting
+      //    limits the transverse momentum of secondary scatterings, due
+      //    to multiple parton interactions, to be less than that of the
+      //    primary interaction (see POWHEG Dijet paper arXiv:1012.3380
+      //    [hep-ph] sec. 4.1 and also the PYTHIA Manual).
+      SetMSTP(86,1);
+      
+      //    maximum number of errors before pythia aborts (def=10)
+      SetMSTU(22,10);
+      //    number of warnings printed on the shell
+      SetMSTU(26,20);
+      break;
+    
+      case kPyCharmPWHG:
+      case kPyBeautyPWHG:
+      //    number of warnings printed on the shell
+      SetMSTU(26,20);
+            
+      break;
     }
 //
 //  Initialize PYTHIA
@@ -625,7 +647,11 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
     
 //  
     SetMSTP(41,1);   // all resonance decays switched on
-    Initialize("CMS",fProjectile,fTarget,fEcms);
+    if (process == kPyJetsPWHG || process == kPyCharmPWHG || process == kPyBeautyPWHG) {
+      Initialize("USER","","",0.);
+    } else {   
+      Initialize("CMS",fProjectile,fTarget,fEcms);
+    }
     fOmegaDalitz.Init();
 }
 
index 2f631f9..4196feb 100644 (file)
@@ -613,13 +613,29 @@ void AliPythia6::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfu
       case kPyMBRDoubleDiffraction:
       case kPyMBRCentralDiffraction:
       break;  
-
+      case kPyJetsPWHG:
+      //    N.B.
+      //    ====
+      //    For the case of jet production the following parameter setting
+      //    limits the transverse momentum of secondary scatterings, due
+      //    to multiple parton interactions, to be less than that of the
+      //    primary interaction (see POWHEG Dijet paper arXiv:1012.3380
+      //    [hep-ph] sec. 4.1 and also the PYTHIA Manual).
+      SetMSTP(86,1);
+      //    maximum number of errors before pythia aborts (def=10)
+      SetMSTU(22,10);
+      //    number of warnings printed on the shell
+      SetMSTU(26,20);
+      break;
     }
 //
 //  Initialize PYTHIA
     SetMSTP(41,1);   // all resonance decays switched on
-    Initialize("CMS",fProjectile,fTarget,fEcms);
-    
+    if (process == kPyJetsPWHG || process == kPyCharmPWHG || process == kPyBeautyPWHG) {
+      Initialize("USER","","",0.);
+    } else {   
+      Initialize("CMS",fProjectile,fTarget,fEcms);
+    }
 }
 
 Int_t AliPythia6::CheckedLuComp(Int_t kf)
index e8103be..ed0ddc4 100644 (file)
@@ -35,7 +35,7 @@ set ( EXPORT )
 
 set ( CSRCS  pythia6.4.25/main.c pythia6.4.25/pythia6_common_address.c)
 
-set ( FSRCS  pythia6.4.25/pythia6_common_block_address.F pythia6.4.25/tpythia6_called_from_cc.F pythia6.4.25/pythia-6.4.25.f pythia6.4.25/pydummy_6.4.25.f)
+set ( FSRCS  pythia6.4.25/pythia6_common_block_address.F pythia6.4.25/tpythia6_called_from_cc.F pythia6.4.25/pythia-6.4.25.f pythia6.4.25/pydummy_6.4.25.f POWHEG-BOX/setup-PYTHIA-lhef.f POWHEG-BOX/lhefread.f)
 
 if( ALICE_TARGET STREQUAL "win32gcc")
        
diff --git a/PYTHIA6/POWHEG-BOX/LesHouches.h b/PYTHIA6/POWHEG-BOX/LesHouches.h
new file mode 100644 (file)
index 0000000..d3bbffd
--- /dev/null
@@ -0,0 +1,18 @@
+c -*-Fortran-*-
+
+      integer maxpup
+      parameter(maxpup=100)
+      integer idbmup,pdfgup,pdfsup,idwtup,nprup,lprup
+      double precision ebmup,xsecup,xerrup,xmaxup
+      common /heprup/ idbmup(2),ebmup(2),pdfgup(2),pdfsup(2),
+     &                idwtup,nprup,xsecup(maxpup),xerrup(maxpup),
+     &                xmaxup(maxpup),lprup(maxpup)
+      integer maxnup
+      parameter (maxnup=500)
+      integer nup,idprup,idup,istup,mothup,icolup
+      double precision xwgtup,scalup,aqedup,aqcdup,pup,vtimup,spinup
+      common/hepeup/nup,idprup,xwgtup,scalup,aqedup,aqcdup,
+     &              idup(maxnup),istup(maxnup),mothup(2,maxnup),
+     &              icolup(2,maxnup),pup(5,maxnup),vtimup(maxnup),
+     &              spinup(maxnup)
+      save /hepeup/
diff --git a/PYTHIA6/POWHEG-BOX/hepevt.h b/PYTHIA6/POWHEG-BOX/hepevt.h
new file mode 100644 (file)
index 0000000..5b49922
--- /dev/null
@@ -0,0 +1,9 @@
+c -*-Fortran-*-
+
+      integer NMXHEP,NEVHEP,NHEP,ISTHEP,IDHEP,
+     &     JMOHEP,JDAHEP
+      double precision phep,vhep
+      PARAMETER (NMXHEP=4000)
+      COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
+     &   JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
+      save /HEPEVT/
diff --git a/PYTHIA6/POWHEG-BOX/lhefread.f b/PYTHIA6/POWHEG-BOX/lhefread.f
new file mode 100644 (file)
index 0000000..5b5cd2f
--- /dev/null
@@ -0,0 +1,64 @@
+c...lhefheader(nlf)
+c...reads initialization information from a les houches events file on unit nlf. 
+      subroutine lhefreadhdr(nlf)
+      implicit none
+      integer nlf
+      character * 100 string
+      integer ipr
+      include 'LesHouches.h'
+ 1    read(nlf,fmt='(a)',err=998,end=998) string
+      if(string(1:5).eq.'<init') then
+         read(nlf,*) idbmup(1),idbmup(2),ebmup(1),ebmup(2),
+     &        pdfgup(1),pdfgup(2),pdfsup(1),pdfsup(2),idwtup,nprup
+         do ipr=1,nprup
+            read(nlf,*) xsecup(ipr),xerrup(ipr),xmaxup(ipr),
+     &           lprup(ipr)
+         enddo
+         goto 999
+      else
+         goto 1
+      endif
+ 998  write(*,*) 'lhefreadhdr: could not find <init> data'
+      call exit(1)
+ 999  end
+
+
+c...reads event information from a les houches events file on unit nlf. 
+      subroutine lhefreadev(nlf)
+      implicit none
+      integer nlf
+      character * 100 string
+      include 'LesHouches.h'
+      integer i,j
+ 1    continue
+      string=' '
+      read(nlf,fmt='(a)',err=777,end=666) string
+      if(string.eq.'</LesHouchesEvents>') then
+         goto 998
+      endif
+      if(string(1:6).eq.'<event') then
+c on error try next event. The error may be cause by merging
+c truncated event files. On EOF return with no event found
+         read(nlf,*,end=998,err=1)nup,idprup,xwgtup,scalup,aqedup,aqcdup
+         do i=1,nup
+            read(nlf,*,end=998,err=1) idup(i),istup(i),mothup(1,i),
+     &           mothup(2,i),icolup(1,i),icolup(2,i),(pup(j,i),j=1,5),
+     &           vtimup(i),spinup(i)
+         enddo
+         goto 999
+      else
+         goto 1
+      endif
+c no event found:
+ 777  continue
+      print *,"Error in reading"
+      print *,string
+      stop
+ 666  continue
+      print *,"reached EOF"
+      print *,string
+      stop
+ 998  continue
+      print *,"read </LesHouchesEvents>"
+      nup=0      
+ 999  end
diff --git a/PYTHIA6/POWHEG-BOX/setup-PYTHIA-lhef.f b/PYTHIA6/POWHEG-BOX/setup-PYTHIA-lhef.f
new file mode 100644 (file)
index 0000000..5b5f688
--- /dev/null
@@ -0,0 +1,30 @@
+      subroutine UPINIT
+      implicit none
+      include 'hepevt.h'
+      include 'LesHouches.h'
+      double precision parp,pari
+      integer mstp,msti
+      common/pypars/mstp(200),parp(200),msti(200),pari(200)
+      integer MSTU,MSTJ
+      double precision PARU,PARJ
+      COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
+      integer MDCY,MDME,KFDP
+      double precision brat
+      COMMON/PYDAT3/MDCY(500,3),MDME(8000,2),BRAT(8000),KFDP(8000,5)
+      integer pycomp
+      external pycomp
+      integer maxev
+      common/mcmaxev/maxev
+      nevhep=0
+c read the header first, so lprup is set
+      call lhefreadhdr(97)
+      end
+
+      subroutine UPEVNT
+      implicit none
+      call lhefreadev(97)
+      end
+
+      subroutine UPVETO
+c pythia routine to abort event
+      end
index 920fb03..519e7f9 100644 (file)
@@ -10,7 +10,7 @@ typedef enum
     kPyCharmppMNR, kPyCharmppMNRwmi, kPyD0ppMNR, kPyDPlusppMNR, kPyDPlusStrangeppMNR, 
     kPyBeautyppMNR, kPyBeautyppMNRwmi, kPyBeautyJets, kPyW, kPyZ, kPyLambdacppMNR, kPyMbMSEL1,
     kPyOldUEQ2ordered, kPyOldUEQ2ordered2, kPyOldPopcorn,
-    kPyLhwgMb, kPyMbDefault, kPyMbAtlasTuneMC09, kPyMBRSingleDiffraction, kPyMBRDoubleDiffraction, kPyMBRCentralDiffraction
+    kPyLhwgMb, kPyMbDefault, kPyMbAtlasTuneMC09, kPyMBRSingleDiffraction, kPyMBRDoubleDiffraction, kPyMBRCentralDiffraction, kPyJetsPWHG, kPyCharmPWHG, kPyBeautyPWHG
 }
 Process_t;
 
index cbaef1e..c9fee3e 100644 (file)
@@ -80120,7 +80120,7 @@ C...on incoming beams and allowed processes.
 
 C...New example: handles a standard Les Houches Events File.
 
-      SUBROUTINE UPINIT
+      SUBROUTINE UPINITD
  
 C...Double precision and integer declarations.
       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
@@ -80240,7 +80240,7 @@ C...HEPEUP commonblock, including (often) an event weight.
 
 C...New example: handles a standard Les Houches Events File.
 
-      SUBROUTINE UPEVNT
+      SUBROUTINE UPEVNTD
  
 C...Double precision and integer declarations.
       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
@@ -80371,7 +80371,7 @@ C...The user decision is to be conveyed by the IVETO value.
 C...IVETO = 0 : retain current event and generate in full;
 C...      = 1 : abort generation of current event and move to next.
  
-      SUBROUTINE UPVETO(IVETO)
+      SUBROUTINE UPVETOD(IVETO)
  
 C...HEPEVT commonblock.
       PARAMETER (NMXHEP=4000)