SetEventListRange();
SetJetPhiRange();
SetJetEtaRange();
+ // Options determining what to keep in the stack (Heavy flavour generation)
+ fStackFillOpt = kFlavorSelection; // Keep particle with selected flavor
+ fFeedDownOpt = kTRUE; // allow feed down from higher family
+ // Fragmentation on/off
+ fFragmentation = kTRUE;
+ // Default counting mode
+ fCountMode = kCountAll;
}
AliGenPythia::AliGenPythia(const AliGenPythia & Pythia)
fPythia->SetCKIN(3,fPtHardMin);
fPythia->SetCKIN(4,fPtHardMax);
if (fNucA1 > 0 && fNucA2 > 0) fPythia->SetNuclei(fNucA1, fNucA2);
+ // Fragmentation?
+ if (fFragmentation) {
+ fPythia->SetMSTP(111,1);
+ } else {
+ fPythia->SetMSTP(111,0);
+ }
+
fPythia->ProcInit(fProcess,fEnergyCMS,fStrucFunc);
// fPythia->Pylist(0);
fParentSelect[3] = 4122;
fFlavorSelect = 4;
break;
+ case kPyD0PbMNR:
+ fParentSelect[0] = 421;
+ fFlavorSelect = 4;
+ break;
case kPyBeauty:
fParentSelect[0]= 511;
fParentSelect[1]= 521;
for (i=0; i< np; i++) {
pParent[i] = -1;
pSelected[i] = 0;
+ trackIt[i] = 0;
}
- printf("\n **************************************************%d\n",np);
- Int_t nc = 0;
+ // printf("\n **************************************************%d\n",np);
+ Int_t nc = 0; // Total n. of selected particles
+ Int_t nParents = 0; // Selected parents
+ Int_t nTkbles = 0; // Trackable particles
if (fProcess != kPyMb && fProcess != kPyJets && fProcess != kPyDirectGamma) {
for (i = 0; i<np; i++) {
kfMo = TMath::Abs(mother->GetPdgCode());
}
// printf("\n particle (all) %d %d %d", i, pSelected[i], kf);
- if (kfl >= fFlavorSelect) {
+ // What to keep in Stack?
+ Bool_t flavorOK = kFALSE;
+ Bool_t selectOK = kFALSE;
+ if (fFeedDownOpt) {
+ if (kfl >= fFlavorSelect) flavorOK = kTRUE;
+ } else {
+ if (kfl > fFlavorSelect) {
+ nc = -1;
+ break;
+ }
+ if (kfl == fFlavorSelect) flavorOK = kTRUE;
+ }
+ switch (fStackFillOpt) {
+ case kFlavorSelection:
+ selectOK = kTRUE;
+ break;
+ case kParentSelection:
+ if (ParentSelected(kf) || kf <= 10) selectOK = kTRUE;
+ break;
+ }
+ if (flavorOK && selectOK) {
//
// Heavy flavor hadron or quark
//
// Kinematic seletion on final state heavy flavor mesons
if (ParentSelected(kf) && !KinematicSelection(iparticle, 0))
{
- nc = -1;
- break;
+ continue;
}
pSelected[i] = 1;
+ if (ParentSelected(kf)) ++nParents; // Update parent count
// printf("\n particle (HF) %d %d %d", i, pSelected[i], kf);
} else {
// Kinematic seletion on decay products
if (fCutOnChild && ParentSelected(kfMo) && ChildSelected(kf)
&& !KinematicSelection(iparticle, 1))
{
- nc = -1;
- break;
+ continue;
}
//
// Decay products
}
if (pSelected[i] == -1) pSelected[i] = 0;
if (!pSelected[i]) continue;
+ // Count quarks only if you did not include fragmentation
+ if (fFragmentation && kf <= 10) continue;
nc++;
// Decision on tracking
trackIt[i] = 0;
} else {
if (ParentSelected(kf)) trackIt[i] = 0;
}
+ if (trackIt[i] == 1) ++nTkbles; // Update trackable counter
//
//
} // particle selection loop
- if (nc > -1) {
+ if (nc > 0) {
for (i = 0; i<np; i++) {
if (!pSelected[i]) continue;
TParticle * iparticle = (TParticle *) fParticles->At(i);
kf = CheckPDGCode(iparticle->GetPdgCode());
- Int_t ks = iparticle->GetStatusCode();
p[0] = iparticle->Px();
p[1] = iparticle->Py();
p[2] = iparticle->Pz();
Int_t ipa = iparticle->GetFirstMother()-1;
Int_t iparent = (ipa > -1) ? pParent[ipa] : -1;
SetTrack(fTrackIt*trackIt[i] ,
- iparent, kf, p, origin, polar, tof, kPPrimary, nt, 1, ks);
+ iparent, kf, p, origin, polar, tof, kPPrimary, nt, 1.);
pParent[i] = nt;
KeepTrack(nt);
} // SetTrack loop
if (trackIt) delete[] trackIt;
if (nc > 0) {
- jev+=nc;
+ switch (fCountMode) {
+ case kCountAll:
+ // printf(" Count all \n");
+ jev += nc;
+ break;
+ case kCountParents:
+ // printf(" Count parents \n");
+ jev += nParents;
+ break;
+ case kCountTrackables:
+ // printf(" Count trackable \n");
+ jev += nTkbles;
+ break;
+ }
if (jev >= fNpart || fNpart == -1) {
fKineBias=Float_t(fNpart)/Float_t(fTrials);
printf("\n Trials: %i %i %i\n",fTrials, fNpart, jev);
origin[2] = fOrigin[2]+iparticle->Vz()/10.;
Float_t tof=kconv*iparticle->T();
SetTrack(fTrackIt*trackIt, iparent, kf, p, origin, polar,
- tof, kPPrimary, nt, 1., ks);
+ tof, kPPrimary, nt);
KeepTrack(nt);
pParent[i] = nt;
} // select particle
class AliGenPythia : public AliGenMC
{
public:
+
+ typedef enum {kFlavorSelection, kParentSelection} StackFillOpt_t;
+ typedef enum {kCountAll, kCountParents, kCountTrackables} CountMode_t;
+
AliGenPythia();
AliGenPythia(Int_t npart);
AliGenPythia(const AliGenPythia &Pythia);
{fEtaMinGamma = etamin; fEtaMaxGamma = etamax;}
virtual void SetGammaPhiRange(Float_t phimin = -180., Float_t phimax = 180.)
{fPhiMinGamma = TMath::Pi()*phimin/180.; fPhiMaxGamma = TMath::Pi()*phimax/180.;}
+ // Set option for feed down from higher family
+ virtual void SetFeedDownHigherFamily(Bool_t opt) {
+ fFeedDownOpt = opt;
+ }
+ // Set option for selecting particles kept in stack according to flavor
+ // or to parent selection
+ virtual void SetStackFillOpt(StackFillOpt_t opt) {
+ fStackFillOpt = opt;
+ }
+ // Set fragmentation option
+ virtual void SetFragmentation(const Bool_t opt) {
+ fFragmentation = opt;
+ }
+ // Set counting mode
+ virtual void SetCountMode(const CountMode_t mode) {
+ fCountMode = mode;
+ }
// get cross section of process
virtual Float_t GetXsection() {return fXsection;}
Float_t fPhiMinGamma; // Minimum phi of triggered gamma
Float_t fPhiMaxGamma; // Maximum phi of triggered gamma
+ StackFillOpt_t fStackFillOpt; // Stack filling with all particles with
+ // that flavour or only with selected
+ // parents and their decays
+ Bool_t fFeedDownOpt; // Option to set feed down from higher
+ // quark families (e.g. b->c)
+ Bool_t fFragmentation; // Option to activate fragmentation by Pythia
+ //
+ // Options for counting when the event will be finished.
+ // fCountMode = kCountAll --> All particles that end up in the
+ // stack are counted
+ // fCountMode = kCountParents --> Only selected parents are counted
+ // fCountMode = kCountTrackabless --> Only particles flagged for tracking
+ // are counted
+ CountMode_t fCountMode;
+
private:
// adjust the weight from kinematic cuts
void AdjustWeights();