Skip to content

Commit

Permalink
Port modifications needed for HF simulations with Pythia8 (#1097)
Browse files Browse the repository at this point in the history 
* Correct use of decayer in AliGenPythiaPlus.
* fPythia object not persistent to avoid default constructur be called when reading galice.root
* Fix index to have the CutOnChild selection working with Pythia8
* Implement kHeavyFlavor as option for SetSTackFillOpt
* Update decay channels and BRs to force charm hadronic decays
* Change to SoftQCD the process for HF generation
* Order by flavour the partons that are first and second mother of a hadron
* debug printf removed
* Change PYTHIA8 to 8.243 (#997)
* Fix for newly imported Pythia8
* Omega_ccc added.
* Update AliGenITSULib.cxx
* Adding Xic+->pK-pi+ (including Xic+->pK0*), Xic+->Xi-pi+pi+ (including Xic+->Xi*0 pi+), Xic0->Xi-pi+ decays, OmegaC->OmegaPi to ForceHadronicD; adding Bs->Ds- pi+ to ForceBeautyUpgrade; fixes to ForceBeautyUpgrade in AliDecayerPythia8. Adding Xic0 and Bs to AliGenITSULib. Decay table modified to update Xic0, Xic+ and OmegaC decay channels.
* Add HF selection case with c and b proportions as in Pythia
  • Loading branch information
@fprino
fprino authored Feb 4, 2020
1 parent cd8461e commit 758f447
Show file tree
Hide file tree
Showing 748 changed files with 1,199,909 additions and 83 deletions.
74 changes: 47 additions & 27 deletions EVGEN/AliGenITSULib.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -70,25 +70,28 @@ GenFunc AliGenITSULib::GetPt(Int_t iPID, const char * sForm) const
case kLb : func=PtLbDist; break;
case kLc : func=PtLcDist; break;
case kXi_c : func=PtLcDist; break;
case kXi_czero : func=PtLcDist; break;
case kOmega_c : func=PtLcDist; break;
case kBplus : func=PtLbDist; break;
case kBzero : func=PtLbDist; break;
case kBs : func=PtLbDist; break;
case kDs : func=PtLcDist; break;
case kDplus : func=PtLcDist; break;
case kOmega_ccc: func=PtLbDist; break;
default : AliError(Form("Unknown particle type: %i, Pt dist is 0",iPID)); func=0;
}
}else {
AliError(Form("Unknown Pt distribution %s. Pt distribution is set to 0 ",sForm));
func=0;
}

return func;
}

GenFunc AliGenITSULib::GetY(Int_t iPID, const char *sForm) const
{
GenFunc func;

if(TMath::Abs(iPID) != kLc && TMath::Abs(iPID) != kLb && TMath::Abs(iPID) != kXi_c && TMath::Abs(iPID) != kBplus && TMath::Abs(iPID) != kBzero && TMath::Abs(iPID)!=kDplus && TMath::Abs(iPID)!=kDs) {
if(TMath::Abs(iPID) != kLc && TMath::Abs(iPID) != kLb && TMath::Abs(iPID) != kXi_c && TMath::Abs(iPID) != kBplus && TMath::Abs(iPID) != kBzero && TMath::Abs(iPID)!=kDplus && TMath::Abs(iPID)!=kDs && TMath::Abs(iPID) != kOmega_ccc && TMath::Abs(iPID) != kBs && TMath::Abs(iPID) != kXi_czero && TMath::Abs(iPID) != kOmega_c) {
AliError(Form("Unknown PID: %i, form: %s, returning 0",iPID,sForm)); //////////
func=0;
} else {
Expand All @@ -99,33 +102,50 @@ GenFunc AliGenITSULib::GetY(Int_t iPID, const char *sForm) const

GenFuncIp AliGenITSULib::GetIp(Int_t iPID, const char *sForm) const
{
AliDebug(1,Form(" %i - %s",iPID,sForm));
AliDebug(1,Form(" %i - %s",iPID,sForm));
// Return pointer to particle type parameterisation
GenFuncIp id;

if(TMath::Abs(iPID) != kLc && TMath::Abs(iPID) != kLb && TMath::Abs(iPID) != kXi_c && TMath::Abs(iPID) != kBplus && TMath::Abs(iPID) != kBzero && TMath::Abs(iPID)!=kDplus && TMath::Abs(iPID)!=kDs) {
AliError(Form("Unknown PID: %i, form: %s, return 0",iPID,sForm)); //////////
id = 0;
} else {
switch (iPID){
case kLc : return id=IpLcPlus;
case -kLc : return id=IpLcMinus;
case kLb : return id=IpLb;
case -kLb : return id=IpLbBar;
case kXi_c : return id=IpXic;
case -kXi_c : return id=IpXicBar;
case kBplus : return id=IpBPlus;
case kBzero : return id=IpB0;
case -kBzero: return id=IpB0Bar;
case -kBplus: return id=IpBMinus;
case kDs : return id=IpDsPlus;
case -kDs : return id=IpDsMinus;
case kDplus : return id=IpDPlus;
case -kDplus: return id=IpDMinus;
default : AliFatal(Form("Unknown particle type: %i",iPID)); id=0;
GenFuncIp id;

if(TMath::Abs(iPID) != kLc && TMath::Abs(iPID) != kLb && TMath::Abs(iPID) != kXi_c && TMath::Abs(iPID) != kBplus && TMath::Abs(iPID) != kBzero && TMath::Abs(iPID)!=kDplus && TMath::Abs(iPID)!=kDs && TMath::Abs(iPID)!=kOmega_ccc && TMath::Abs(iPID) != kBs && TMath::Abs(iPID) != kXi_czero && TMath::Abs(iPID) != kOmega_c) {
AliError(Form("Unknown PID: %i, form: %s, return 0",iPID,sForm)); //////////
id = 0;
} else {
switch (iPID){
case kLc : return id=IpLcPlus;
case -kLc : return id=IpLcMinus;
case kLb : return id=IpLb;
case -kLb : return id=IpLbBar;
case kXi_c : return id=IpXic;
case -kXi_c : return id=IpXicBar;
case kBplus : return id=IpBPlus;
case kBzero : return id=IpB0;
case -kBzero: return id=IpB0Bar;
case -kBplus: return id=IpBMinus;
case kDs : return id=IpDsPlus;
case -kDs : return id=IpDsMinus;
case kDplus : return id=IpDPlus;
case -kDplus: return id=IpDMinus;
case kOmega_ccc: return id=IpOmegaccc;
case kXi_czero : return id=IpXicZero;
case -kXi_czero : return id=IpXicZeroBar;
case kOmega_c : return id=IpOmegac;
case -kOmega_c : return id=IpOmegacBar;
case kBs : return id=IpBs;
case -kBs : return id=IpBsBar;

default : AliFatal(Form("Unknown particle type: %i",iPID)); id=0;
}

}

}
return id;
}

return id;
Int_t AliGenITSULib::IpOmegaccc(TRandom* ran)
{
Int_t ip = 4444;
if (ran->Rndm() < 0.5) ip = -ip;
return ip;
}


12 changes: 10 additions & 2 deletions EVGEN/AliGenITSULib.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@ class AliGenITSULib :public AliGenLib {

public:

enum EPartId {kLb=5122,kLc=4122,kXi_c = 4232,kBplus = 521,kBzero = 511,kDs=431,kDplus=411};
enum EPartId {kLb=5122,kLc=4122,kXi_c = 4232,kBplus = 521,kBzero = 511,kDs=431,kDplus=411, kOmega_ccc=4444,kXi_czero = 4132, kOmega_c=4332, kBs=531};

//Getters

Expand All @@ -31,17 +31,25 @@ class AliGenITSULib :public AliGenLib {
static Int_t IpLcMinus(TRandom * /*ran*/) {return -(int)kLc;}
static Int_t IpLb(TRandom * /*ran*/) {return (int)kLb;}
static Int_t IpLbBar(TRandom * /*ran*/) {return -(int)kLb;}
static Int_t IpXicZero(TRandom * /*ran*/) {return (int)kXi_czero;}
static Int_t IpXicZeroBar(TRandom * /*ran*/) {return -(int)kXi_czero;}
static Int_t IpXic(TRandom * /*ran*/) {return (int)kXi_c;}
static Int_t IpXicBar(TRandom * /*ran*/) {return -(int)kXi_c;}
static Int_t IpOmegac(TRandom * /*ran*/) {return (int)kOmega_c;}
static Int_t IpOmegacBar(TRandom * /*ran*/) {return -(int)kOmega_c;}
static Int_t IpBPlus(TRandom * /*ran*/) {return (int)kBplus;}
static Int_t IpBMinus(TRandom * /*ran*/) {return -(int)kBplus;}
static Int_t IpBs(TRandom * /*ran*/) {return (int)kBs;}
static Int_t IpBsBar(TRandom * /*ran*/) {return -(int)kBs;}
static Int_t IpB0(TRandom * /*ran*/) {return (int)kBzero;}
static Int_t IpB0Bar(TRandom * /*ran*/) {return -(int)kBzero;}
static Int_t IpDsPlus(TRandom * /*ran*/) {return (int)kDs;}
static Int_t IpDsMinus(TRandom * /*ran*/) {return -(int)kDs;}
static Int_t IpDPlus(TRandom * /*ran*/) {return (int)kDplus;}
static Int_t IpDMinus(TRandom * /*ran*/) {return -(int)kDplus;}
static Int_t IpDMinus(TRandom * /*ran*/) {return -(int)kDplus;}

static Int_t IpOmegaccc(TRandom * /*ran*/);

static Double_t PtFlat(const Double_t * /*px*/, const Double_t * /*dummy*/) {return 1;}
static Double_t YFlat (const Double_t * /*py*/, const Double_t * /*dummy*/) {return 1;}

Expand Down
47 changes: 31 additions & 16 deletions PYTHIA6/AliPythia6/AliDecayerPythia.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -653,15 +653,16 @@ void AliDecayerPythia::ForceBeautyUpgrade()
} else {
ForceParticleDecay( 5122, 4122, 1);
}
// Xi_c
ForceParticleDecay( 4232, 3312, 1);
// B+ -> D0pi+
const Int_t prod[2]={421,211};
Int_t mult[2]={1,1};
ForceParticleDecay(521,prod,mult,2,1);
// B0 -> D*-pi+
const Int_t prod2[2]={413,211};
ForceParticleDecay(511,prod2,mult,2,1);
// Bs -> Ds-pi+
const Int_t prod3[2]={431,211};
ForceParticleDecay(531,prod3,mult,2,1);
// force charm hadronic decays (D mesons and Lc)
ForceHadronicD(0,0,0);
}
Expand Down Expand Up @@ -691,18 +692,17 @@ Int_t AliDecayerPythia::CountProducts(Int_t channel, Int_t particle)
void AliDecayerPythia::ForceHadronicD(Int_t optUse4Bodies, Int_t optUseDtoV0, Int_t optForceLcChannel, Int_t optUsePionicD0)
{

// Xic->Xipipi
Int_t Xic = 4232;
Int_t productsX[3] = {3312, 211,211}, multX[3] = {1, 1, 1};
ForceParticleDecay(Xic, productsX, multX, 3);
//
// OmegaC -> Omega pi
Int_t iOmegaC=4332;
Int_t productsOc[2]={211, 3334},multOc[2]={1,1};
ForceParticleDecay(iOmegaC, productsOc, multOc, 2,1);


// Force golden D decay modes
//
const Int_t kNHadrons = 5;
const Int_t kNHadrons = 7;
Int_t channel;
Int_t hadron[kNHadrons] = {411, 421, 431, 4112, 4122};
Int_t hadron[kNHadrons] = {411, 421, 431, 4112, 4122, 4232, 4132};


// for D+ -> K0* (-> K- pi+) pi+
Expand All @@ -728,31 +728,39 @@ void AliDecayerPythia::ForceHadronicD(Int_t optUse4Bodies, Int_t optUseDtoV0, In
Int_t iLambda=3122;
//for Lambda_c -> antiK0 p
Int_t iK0bar=-311;

// for Xic+->Xi-pi+pi+, Xic+->Xi0*pi+,Xi0*->Xi-pi+ and Xic0->Xi-pi+ channels
Int_t iXiMinus=3312,iXiStar0=3324;


Int_t decayP1[kNHadrons][4] =
{
{kKMinus, kPiPlus, kPiPlus, 0},
{kKMinus, kPiPlus, 0 , 0},
{kKPlus , iKstarbar0, 0 , 0},
{-1 , -1 , -1 , -1},
{kProton, iKstarbar0, 0 , 0}
{kProton, iKstarbar0, 0 , 0},
{kProton, iKstarbar0, 0 , 0},
{kPiPlus, iXiMinus, 0 , 0}
};
Int_t decayP2[kNHadrons][4] =
{
{iKstarbar0, kPiPlus, 0 , 0},
{kKMinus , kPiPlus, kPiPlus, kPiMinus},
{iPhi , kPiPlus, 0 , 0},
{-1 , -1 , -1 , -1},
{iDeltaPP , kKMinus, 0 , 0}
{iDeltaPP , kKMinus, 0 , 0},
{kProton , kKMinus, kPiPlus, 0},
{-1 , -1 , -1 , -1}
};
Int_t decayP3[kNHadrons][4] =
{
{kPiPlus , iPhi , 0 , 0},
{kKMinus , kPiPlus, iRho0 , 0},
{kKPlus , iK0bar , 0 , 0},
{-1 , -1 , -1 , -1},
{kProton , kKMinus, kPiPlus , 0}
{kProton , kKMinus, kPiPlus , 0},
{kPiPlus , kPiPlus, iXiMinus, 0},
{-1 , -1 , -1 , -1}
};
// for Lambda_c -> Lambda_1520 pi+ -> p K- pi+, D0-> K*0 pi+ pi- -> K3pi
Int_t decayP4[kNHadrons][4] =
Expand All @@ -761,7 +769,9 @@ void AliDecayerPythia::ForceHadronicD(Int_t optUse4Bodies, Int_t optUseDtoV0, In
{iKstarbar0 , kPiPlus , kPiMinus, 0},
{-1 , -1 , -1 , -1},
{-1 , -1 , -1 , -1},
{iLambda1520 , kPiPlus , 0 , 0}
{iLambda1520 , kPiPlus , 0 , 0},
{kPiPlus , iXiStar0, 0 , 0},
{-1 , -1 , -1 , -1}
};
// for Lambda_c -> Lambda pi+, D0 -> pi0 pi+ pi-
Int_t decayP5[kNHadrons][4] =
Expand All @@ -770,7 +780,9 @@ void AliDecayerPythia::ForceHadronicD(Int_t optUse4Bodies, Int_t optUseDtoV0, In
{kPiPlus , kPiMinus, kPi0 , 0},
{-1 , -1 , -1 , -1},
{-1 , -1 , -1 , -1},
{iLambda , kPiPlus, 0 , 0}
{iLambda , kPiPlus, 0 , 0},
{-1 , -1 , -1 , -1},
{-1 , -1 , -1 , -1}
};

// for Lambda_c -> K0bar p, D0 -> K- pi+ pi0, D+ -> K0s pi+ pi0
Expand All @@ -780,7 +792,10 @@ void AliDecayerPythia::ForceHadronicD(Int_t optUse4Bodies, Int_t optUseDtoV0, In
{kKMinus , kPiPlus, kPi0 , 0},
{-1 , -1 , -1 , -1},
{-1 , -1 , -1 , -1},
{kProton , iK0bar, 0 , 0}
{kProton , iK0bar, 0 , 0},
{-1 , -1 , -1 , -1},
{-1 , -1 , -1 , -1}

};

if(optUse4Bodies==0){
Expand Down
20 changes: 17 additions & 3 deletions PYTHIA6/AliPythia6/AliGenPythia.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -522,6 +522,16 @@ void AliGenPythia::Init()
fParentSelect[6]= 5332;
fFlavorSelect = 5;
break;
case kPyHeavyFlavppMNRwmi:
fParentSelect[0]= 511; //settings to selct decay products
fParentSelect[1]= 521;
fParentSelect[2]= 531;
fParentSelect[3]= 5122;
fParentSelect[4]= 5132;
fParentSelect[5]= 5232;
fParentSelect[6]= 5332;
fFlavorSelect = 5;
break;
case kPyBeautyUnforced:
fParentSelect[0] = 511;
fParentSelect[1] = 521;
Expand Down Expand Up @@ -793,6 +803,7 @@ void AliGenPythia::Generate()
fProcess != kPyZgamma &&
fProcess != kPyCharmppMNRwmi &&
fProcess != kPyBeautyppMNRwmi &&
fProcess != kPyHeavyFlavppMNRwmi &&
fProcess != kPyBeautyJets &&
fProcess != kPyWPWHG &&
fProcess != kPyJetsPWHG &&
Expand Down Expand Up @@ -1125,7 +1136,7 @@ Int_t AliGenPythia::GenerateMB()
// Check if there is a ccbar or bbbar pair with at least one of the two
// in fYMin < y < fYMax

if (fProcess == kPyCharmppMNRwmi || fProcess == kPyBeautyppMNRwmi || fProcess == kPyBeautyJets) {
if (fProcess == kPyCharmppMNRwmi || fProcess == kPyBeautyppMNRwmi || fProcess == kPyHeavyFlavppMNRwmi || fProcess == kPyBeautyJets) {
TParticle *partCheck;
TParticle *mother;
Bool_t theQ=kFALSE,theQbar=kFALSE,inYcut=kFALSE;
Expand All @@ -1139,8 +1150,11 @@ Int_t AliGenPythia::GenerateMB()

for(i = 0; i < np; i++) {
partCheck = (TParticle*)fParticles.At(i);
pdg = partCheck->GetPdgCode();
if(TMath::Abs(pdg) == fFlavorSelect) { // quark
pdg = partCheck->GetPdgCode();
Bool_t flavSel=kFALSE;
if(TMath::Abs(pdg) == fFlavorSelect) flavSel=kTRUE;
if(fProcess == kPyHeavyFlavppMNRwmi && (TMath::Abs(pdg) == 4 || TMath::Abs(pdg) == 5)) flavSel=kTRUE;
if(flavSel){ // 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));
Expand Down
Loading

0 comments on commit 758f447

Please sign in to comment.