#include "AliMC.h"
#include "AliLog.h"
#include "PyquenCommon.h"
-#include "AliLog.h"
ClassImp(AliGenPythia)
fFlavorSelect(0),
fXsection(0.),
fPythia(0),
+ fWeightPower(0.),
fPtHardMin(0.),
fPtHardMax(1.e4),
fYHardMin(-1.e10),
fYHardMax(1.e10),
fGinit(1),
fGfinal(1),
+ fCRoff(0),
fHadronisation(1),
fPatchOmegaDalitz(0),
+ fDecayerExodus(0),
fNpartons(0),
fReadFromFile(0),
fReadLHEF(0),
fEtaMaxGamma(20.),
fPhiMinGamma(0.),
fPhiMaxGamma(2. * TMath::Pi()),
+ fUseYCutHQ(kFALSE),
+ fYMinHQ(-20.),
+ fYMaxHQ(20.),
fPycellEtaMax(2.),
fPycellNEta(274),
fPycellNPhi(432),
fFeedDownOpt(kTRUE),
fFragmentation(kTRUE),
fSetNuclei(kFALSE),
+ fUseNuclearPDF(kFALSE),
+ fUseLorentzBoost(kTRUE),
fNewMIS(kFALSE),
fHFoff(kFALSE),
fNucPdf(0),
fTriggerParticle(0),
fTriggerEta(0.9),
+ fTriggerY(999.),
+ fTriggerEtaMin(0.9),
fTriggerMinPt(-1),
fTriggerMaxPt(1000),
fTriggerMultiplicity(0),
fTriggerMultiplicityEta(0),
+ fTriggerMultiplicityEtaMin(0),
+ fTriggerMultiplicityEtaMax(0),
fTriggerMultiplicityPtMin(0),
fCountMode(kCountAll),
fHeader(0),
fFlavorSelect(0),
fXsection(0.),
fPythia(0),
+ fWeightPower(0.),
fPtHardMin(0.),
fPtHardMax(1.e4),
fYHardMin(-1.e10),
fYHardMax(1.e10),
fGinit(kTRUE),
fGfinal(kTRUE),
+ fCRoff(kFALSE),
fHadronisation(kTRUE),
- fPatchOmegaDalitz(0),
+ fPatchOmegaDalitz(0),
+ fDecayerExodus(0),
fNpartons(0),
fReadFromFile(kFALSE),
fReadLHEF(0),
fEtaMaxGamma(20.),
fPhiMinGamma(0.),
fPhiMaxGamma(2. * TMath::Pi()),
+ fUseYCutHQ(kFALSE),
+ fYMinHQ(-20.),
+ fYMaxHQ(20.),
fPycellEtaMax(2.),
fPycellNEta(274),
fPycellNPhi(432),
fFeedDownOpt(kTRUE),
fFragmentation(kTRUE),
fSetNuclei(kFALSE),
+ fUseNuclearPDF(kFALSE),
+ fUseLorentzBoost(kTRUE),
fNewMIS(kFALSE),
fHFoff(kFALSE),
fNucPdf(0),
fTriggerParticle(0),
fTriggerEta(0.9),
+ fTriggerY(999.),
+ fTriggerEtaMin(0.9),
fTriggerMinPt(-1),
fTriggerMaxPt(1000),
fTriggerMultiplicity(0),
fTriggerMultiplicityEta(0),
+ fTriggerMultiplicityEtaMin(0),
+ fTriggerMultiplicityEtaMax(0),
fTriggerMultiplicityPtMin(0),
fCountMode(kCountAll),
fHeader(0),
//
fParentWeight=1./Float_t(fNpart);
//
-
-
+ if (fWeightPower != 0)
+ fPythia->SetWeightPower(fWeightPower);
fPythia->SetCKIN(3,fPtHardMin);
fPythia->SetCKIN(4,fPtHardMax);
fPythia->SetCKIN(7,fYHardMin);
fPythia->SetCKIN(8,fYHardMax);
- if (fAProjectile > 0 && fATarget > 0) fPythia->SetNuclei(fAProjectile, fATarget, fNucPdf);
if (fProjectile != "p" || fTarget != "p") fPythia->SetCollisionSystem(fProjectile,fTarget);
+ if(fUseNuclearPDF)
+ fPythia->SetNuclei(fAProjectile, fATarget, fNucPdf);
// Fragmentation?
if (fFragmentation) {
fPythia->SetMSTP(111,1);
fPythia->SetMSTP(61,fGinit);
// final state radiation
fPythia->SetMSTP(71,fGfinal);
+ //color reconnection strength
+ if(fCRoff==1)fPythia->SetMSTP(95,0);
// pt - kick
if (fPtKick > 0.) {
fPythia->SetMSTP(91,1);
case kPyDirectGamma:
case kPyLhwgMb:
break;
+ case kPyWPWHG:
case kPyW:
case kPyZ:
+ case kPyZgamma:
case kPyMBRSingleDiffraction:
case kPyMBRDoubleDiffraction:
case kPyMBRCentralDiffraction:
//
//
AliGenMC::Init();
+
+ // Reset Lorentz boost if demanded
+ if(!fUseLorentzBoost) {
+ fDyBoost = 0;
+ Warning("Init","Demand to discard Lorentz boost.\n");
+ }
//
//
//
if (fSetNuclei) {
- fDyBoost = 0;
- Warning("Init","SetNuclei used. Use SetProjectile + SetTarget instead. fDyBoost has been reset to 0\n");
+ fDyBoost = 0;
+ Warning("Init","Deprecated function SetNuclei() used (nPDFs + no boost). Use SetProjectile + SetTarget + SetUseNuclearPDF + SetUseLorentzBoost instead.\n");
}
-
fPythia->SetPARJ(200, 0.0);
fPythia->SetPARJ(199, 0.0);
fPythia->SetPARJ(198, 0.0);
fPythia->SetMSTJ(44, 2); // option to run alpha_s
fPythia->SetPARJ(82, 1.); // Cut off for parton showers
}
+
+ if ( AliLog::GetDebugLevel("","AliGenPythia") >= 1 ) {
+ fPythia->Pystat(4);
+ fPythia->Pystat(5);
+ }
}
void AliGenPythia::Generate()
fPythia->DalitzDecays();
fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 1);
}
- fPythia->Pyexec();
- }
+
+ else if (fDecayerExodus) {
+
+ fPythia->SetMDCY(fPythia->Pycomp(22) ,1, 0);
+ fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 0);
+ fPythia->SetMDCY(fPythia->Pycomp(221) ,1, 0);
+ fPythia->Pyexec();
+ fPythia->OmegaDalitz();
+ fPythia->SetMDCY(fPythia->Pycomp(111) ,1, 1);
+ fPythia->PizeroDalitz();
+ fPythia->PhiDalitz();
+ fPythia->SetMDCY(fPythia->Pycomp(221) ,1, 1);
+ fPythia->EtaDalitz();
+ fPythia->EtaprimeDalitz();
+ fPythia->SetMDCY(fPythia->Pycomp(22) ,1, 1);
+ fPythia->RhoDirect();
+ fPythia->OmegaDirect();
+ fPythia->PhiDirect();
+ fPythia->JPsiDirect();
+ }
+
+ fPythia->Pyexec();
+ }
fTrials++;
fPythia->ImportParticles(&fParticles,"All");
fProcess != kPyMbMSEL1 &&
fProcess != kPyW &&
fProcess != kPyZ &&
+ fProcess != kPyZgamma &&
fProcess != kPyCharmppMNRwmi &&
fProcess != kPyBeautyppMNRwmi &&
fProcess != kPyBeautyJets &&
- fProcess != kPyJetsPWHG &&
- fProcess != kPyCharmPWHG &&
- fProcess != kPyBeautyPWHG) {
+ fProcess != kPyWPWHG &&
+ fProcess != kPyJetsPWHG &&
+ fProcess != kPyCharmPWHG &&
+ fProcess != kPyBeautyPWHG) {
for (i = 0; i < np; i++) {
TParticle* iparticle = (TParticle *) fParticles.At(i);
// eta cut
if (fTriggerMultiplicityEta > 0 && TMath::Abs(iparticle->Eta()) > fTriggerMultiplicityEta)
continue;
+ //multiplicity check for a given eta range
+ if ((fTriggerMultiplicityEtaMin != fTriggerMultiplicityEtaMax) &&
+ (iparticle->Eta() < fTriggerMultiplicityEtaMin || iparticle->Eta() > fTriggerMultiplicityEtaMax))
+ continue;
// pt cut
if (iparticle->Pt() < fTriggerMultiplicityPtMin)
continue;
kf = CheckPDGCode(iparticle->GetPdgCode());
if (kf != fTriggerParticle) continue;
if (iparticle->Pt() == 0.) continue;
- if (TMath::Abs(iparticle->Eta()) > fTriggerEta) continue;
+ if (TMath::Abs(iparticle->Y()) > fTriggerY) continue;
+ if (fTriggerEtaMin == fTriggerEta) {
+ if (TMath::Abs(iparticle->Eta()) > fTriggerEta) continue;
+ } else {
+ if (iparticle->Eta() > fTriggerEta || iparticle->Eta() < fTriggerEtaMin) continue;
+ }
if ( iparticle->Pt() > fTriggerMaxPt || iparticle->Pt() < fTriggerMinPt ) continue;
triggered = kTRUE;
break;
Bool_t theChild=kFALSE;
Bool_t triggered=kFALSE;
Float_t y;
- Int_t pdg,mpdg,mpdgUpperFamily;
- for(i=0; i<np; i++) {
+ Int_t pdg, mpdg, mpdgUpperFamily;
+ for(i = 0; i < np; i++) {
partCheck = (TParticle*)fParticles.At(i);
pdg = partCheck->GetPdgCode();
if(TMath::Abs(pdg) == fFlavorSelect) { // quark
if(pdg>0) { theQ=kTRUE; } else { theQbar=kTRUE; }
y = 0.5*TMath::Log((partCheck->Energy()+partCheck->Pz()+1.e-13)/
(partCheck->Energy()-partCheck->Pz()+1.e-13));
- if(y>fYMin && y<fYMax) inYcut=kTRUE;
+ if(fUseYCutHQ && y>fYMinHQ && y<fYMaxHQ) inYcut=kTRUE;
+ if(!fUseYCutHQ && y>fYMin && y<fYMax) inYcut=kTRUE;
}
if(fTriggerParticle) {
if(TMath::Abs(pdg)==fTriggerParticle) triggered=kTRUE;
}
//Introducing child cuts in case kPyW, kPyZ, kPyMb, and kPyMbNonDiff
- if ( (fProcess == kPyW ||
+ if ( (
+ fProcess == kPyWPWHG ||
+ fProcess == kPyW ||
fProcess == kPyZ ||
+ fProcess == kPyZgamma ||
fProcess == kPyMbDefault ||
fProcess == kPyMb ||
fProcess == kPyMbAtlasTuneMC09 ||
return 0;
}
}
-
+
for (i = 0; i < np; i++) {
Int_t trackIt = 0;
{
// Treat protons as inside nuclei with mass numbers a1 and a2
- fAProjectile = a1;
- fATarget = a2;
- fNucPdf = pdfset; // 0 EKS98 9 EPS09LO 19 EPS09NLO
- fSetNuclei = kTRUE;
+ fAProjectile = a1;
+ fATarget = a2;
+ fNucPdf = pdfset; // 0 EKS98 9 EPS09LO 19 EPS09NLO
+ fUseNuclearPDF = kTRUE;
+ fSetNuclei = kTRUE;
}
// Builds the event header, to be called after each event
if (fHeader) delete fHeader;
fHeader = new AliGenPythiaEventHeader("Pythia");
+ fHeader->SetTitle(GetTitle());
//
// Event type
if(fProcDiff>0){
// Number of primaries
fHeader->SetNProduced(fNprimaries);
//
+// Event weight
+ fHeader->SetEventWeight(fPythia->GetVINT(97));
+//
// Jets that have triggered
//Need to store jets for b-jet studies too!
((AliGenPythiaEventHeader*) fHeader)->SetZQuench(z);
}
//
-// Store pt^hard
+// Store pt^hard and cross-section
((AliGenPythiaEventHeader*) fHeader)->SetPtHard(fPythia->GetVINT(47));
+ ((AliGenPythiaEventHeader*) fHeader)->SetXsection(fPythia->GetPARI(1));
//
// Store Event Weight
- ((AliGenPythiaEventHeader*) fHeader)->SetEventWeight(fPythia->GetPARI(7));
-
+ ((AliGenPythiaEventHeader*) fHeader)->SetEventWeight(fPythia->GetPARI(7)*fPythia->GetPARI(10));
+ // PARI(7) is 1 or -1, for weighted generation with accept/reject, e.g. POWHEG
+ // PARI(10) is a weight associated with reweighted generation, using Pyevwt
//
// Pass header
//
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff