Skip to content

Commit 0027a9c

Browse files
Adding code to run 5 TeV simulation
1 parent a908794 commit 0027a9c

File tree

3 files changed

+101
-0
lines changed

3 files changed

+101
-0
lines changed
Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,7 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8Onia_PromptSignals_gaptriggered.C
4+
funcName=GeneratorPromptJpsi_EvtGenFwdy(5,-4.3,-2.3)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_onia_triggerGap_pp5TeV.cfg
Lines changed: 82 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,82 @@
1+
int External()
2+
{
3+
int checkPdgSignal[] = {443};
4+
int checkPdgDecay = 13;
5+
double rapiditymin = -4.3; double rapiditymax = -2.3;
6+
std::string path{"o2sim_Kine.root"};
7+
std::cout << "Check for\nsignal PDG " << checkPdgSignal[0] << "\n decay PDG " << checkPdgDecay << "\n";
8+
TFile file(path.c_str(), "READ");
9+
if (file.IsZombie()) {
10+
std::cerr << "Cannot open ROOT file " << path << "\n";
11+
return 1;
12+
}
13+
14+
auto tree = (TTree*)file.Get("o2sim");
15+
std::vector<o2::MCTrack>* tracks{};
16+
tree->SetBranchAddress("MCTrack", &tracks);
17+
18+
int nLeptons{};
19+
int nAntileptons{};
20+
int nLeptonPairs{};
21+
int nLeptonPairsToBeDone{};
22+
int nSignalJpsi{};
23+
int nSignalJpsiWithinAcc{};
24+
auto nEvents = tree->GetEntries();
25+
o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine);
26+
Bool_t isInjected = kFALSE;
27+
28+
for (int i = 0; i < nEvents; i++) {
29+
tree->GetEntry(i);
30+
for (auto& track : *tracks) {
31+
auto pdg = track.GetPdgCode();
32+
auto rapidity = track.GetRapidity();
33+
auto idMoth = track.getMotherTrackId();
34+
if (pdg == checkPdgDecay) {
35+
// count leptons
36+
nLeptons++;
37+
} else if(pdg == -checkPdgDecay) {
38+
// count anti-leptons
39+
nAntileptons++;
40+
} else if (pdg == checkPdgSignal[0]) {
41+
if(idMoth < 0){
42+
// count signal PDG
43+
nSignalJpsi++;
44+
// count signal PDG within acceptance
45+
if(rapidity > rapiditymin && rapidity < rapiditymax) nSignalJpsiWithinAcc++;
46+
}
47+
auto child0 = o2::mcutils::MCTrackNavigator::getDaughter0(track, *tracks);
48+
auto child1 = o2::mcutils::MCTrackNavigator::getDaughter1(track, *tracks);
49+
if (child0 != nullptr && child1 != nullptr) {
50+
// check for parent-child relations
51+
auto pdg0 = child0->GetPdgCode();
52+
auto pdg1 = child1->GetPdgCode();
53+
std::cout << "First and last children of parent " << checkPdgSignal << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n";
54+
if (std::abs(pdg0) == checkPdgDecay && std::abs(pdg1) == checkPdgDecay && pdg0 == -pdg1) {
55+
nLeptonPairs++;
56+
if (child0->getToBeDone() && child1->getToBeDone()) {
57+
nLeptonPairsToBeDone++;
58+
}
59+
}
60+
}
61+
}
62+
}
63+
}
64+
std::cout << "#events: " << nEvents << "\n"
65+
<< "#leptons: " << nLeptons << "\n"
66+
<< "#antileptons: " << nAntileptons << "\n"
67+
<< "#signal (prompt Jpsi): " << nSignalJpsi << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalJpsiWithinAcc << "\n"
68+
<< "#lepton pairs: " << nLeptonPairs << "\n"
69+
<< "#lepton pairs to be done: " << nLeptonPairs << "\n";
70+
71+
72+
if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) {
73+
std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n";
74+
return 1;
75+
}
76+
if (nLeptonPairs != nLeptonPairsToBeDone) {
77+
std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n";
78+
return 1;
79+
}
80+
81+
return 0;
82+
}
Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,12 @@
1+
### beams
2+
Beams:idA 2212 # proton
3+
Beams:idB 2212 # proton
4+
Beams:eCM 5360. # GeV
5+
6+
### processes
7+
SoftQCD:inelastic on # all inelastic processes
8+
CharmoniumShower:all = on
9+
10+
### decays
11+
ParticleDecays:limitTau0 on
12+
ParticleDecays:tau0Max 10.

0 commit comments

Comments
 (0)