Add extra files for HF -> mumu simulations in DQ

MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C

@@ -162,7 +162,7 @@ FairGenerator*
   auto gen = new o2::eventgen::GeneratorEvtGen<o2::eventgen::GeneratorHFToMu_EvtGenFwdY_gaptriggered>();
   if (isbb == false) {
     gen->setPDG(4);
-    pdgs = "411;421;431;4122;4232;4332";
+    pdgs = "411;421;431;4122;4132;4232;4332";
   } else {
     gen->setPDG(5);
     pdgs = "511;521;531;541;5112;5122;5232;5132;5332";
@@ -171,7 +171,7 @@ FairGenerator*
 
   gen->setRapidityHadron(rapidityMin,rapidityMax);
   gen->setHadronMultiplicity(1);
-  TString pathO2table = gSystem->ExpandPathName("$O2DPG_MC_CONFIG_ROOT/MC/config/PWGDQ/pythia8/decayer/switchOffChadrons.cfg");
+  TString pathO2table = gSystem->ExpandPathName("$O2DPG_MC_CONFIG_ROOT/MC/config/PWGDQ/pythia8/decayer/switchOffCBhadrons.cfg");
   gen->readFile(pathO2table.Data());
   gen->setConfigMBdecays(pathO2table);
   gen->setVerbose(verbose);

MC/config/PWGDQ/ini/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap.ini

@@ -3,12 +3,12 @@
 
 [GeneratorExternal]
 fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C 
-funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, true, false)
+funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, true, true)
 
 ### The external generator derives from GeneratorPythia8.
 ### This part configures the bits of the interface: configuration and user hooks
 
 [GeneratorPythia8]
-config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_hf.cfg
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_hard_cr2.cfg
 hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
 hooksFuncName = pythia8_userhooks_bbbar(-4.3, -2.3)

MC/config/PWGDQ/ini/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_inel.ini

@@ -0,0 +1,14 @@
+### The setup uses an external event generator
+### This part sets the path of the file and the function call to retrieve it
+
+[GeneratorExternal]
+fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C 
+funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, true, true)
+
+### The external generator derives from GeneratorPythia8.
+### This part configures the bits of the interface: configuration and user hooks
+
+[GeneratorPythia8]
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_cr2.cfg
+hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
+hooksFuncName = pythia8_userhooks_bbbar(-4.3, -2.3)

MC/config/PWGDQ/ini/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_natural.ini

@@ -0,0 +1,14 @@
+### The setup uses an external event generator
+### This part sets the path of the file and the function call to retrieve it
+
+[GeneratorExternal]
+fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C 
+funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, true, false)
+
+### The external generator derives from GeneratorPythia8.
+### This part configures the bits of the interface: configuration and user hooks
+
+[GeneratorPythia8]
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_hard_cr2.cfg
+hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
+hooksFuncName = pythia8_userhooks_bbbar(-4.3, -2.3)

MC/config/PWGDQ/ini/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_natural_inel.ini

@@ -0,0 +1,14 @@
+### The setup uses an external event generator
+### This part sets the path of the file and the function call to retrieve it
+
+[GeneratorExternal]
+fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C 
+funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, true, false)
+
+### The external generator derives from GeneratorPythia8.
+### This part configures the bits of the interface: configuration and user hooks
+
+[GeneratorPythia8]
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_cr2.cfg
+hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
+hooksFuncName = pythia8_userhooks_bbbar(-4.3, -2.3)

MC/config/PWGDQ/ini/GeneratorHF_ccbarToMuonsSemileptonic_fwdy_inel_TriggerGap.ini

@@ -9,6 +9,6 @@ funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, false,
 ### This part configures the bits of the interface: configuration and user hooks
 
 [GeneratorPythia8]
-config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_inel.cfg
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_cr2.cfg
 hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
 hooksFuncName = pythia8_userhooks_ccbar(-4.3, -2.3)

MC/config/PWGDQ/ini/GeneratorHF_ccbarToMuonsSemileptonic_fwdy_inel_TriggerGap_natural.ini

@@ -0,0 +1,14 @@
+### The setup uses an external event generator
+### This part sets the path of the file and the function call to retrieve it
+
+[GeneratorExternal]
+fileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq.C 
+funcName = GeneratorHFToMu_EvtGenFwdY_gaptriggered_dq(-4.3, -2.3, false, false, false)
+
+### The external generator derives from GeneratorPythia8.
+### This part configures the bits of the interface: configuration and user hooks
+
+[GeneratorPythia8]
+config = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_cr2.cfg
+hooksFileName = ${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
+hooksFuncName = pythia8_userhooks_ccbar(-4.3, -2.3)

MC/config/PWGDQ/ini/tests/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_inel.C

@@ -0,0 +1,77 @@
+int External()
+{
+  int checkPdgDecay = 13;
+  std::string path{"o2sim_Kine.root"};
+  TFile file(path.c_str(), "READ");
+  if (file.IsZombie()) {
+    std::cerr << "Cannot open ROOT file " << path << "\n";
+    return 1;
+  }
+  auto tree = (TTree*)file.Get("o2sim");
+  std::vector<o2::MCTrack>* tracks{};
+  tree->SetBranchAddress("MCTrack", &tracks);
+
+  int nLeptons{};
+  int nLeptonsInAcceptance{};
+  int nLeptonsToBeDone{};
+  int nSignalPairs{};
+  int nLeptonPairs{};
+  int nLeptonPairsInAcceptance{};
+  int nLeptonPairsToBeDone{};
+  auto nEvents = tree->GetEntries();
+
+  for (int i = 0; i < nEvents; i++) {
+    tree->GetEntry(i);
+    int nleptonseinacc = 0;
+    int nleptonse = 0;
+    int nleptonseToBeDone = 0;
+    int nopenHeavy = 0;
+    for (auto& track : *tracks) {
+      auto pdg = track.GetPdgCode();
+      auto y = track.GetRapidity();
+      if (std::abs(pdg) == checkPdgDecay) {
+        int igmother = track.getMotherTrackId();
+        if (igmother > 0) {
+          auto gmTrack = (*tracks)[igmother];
+          int gmpdg = gmTrack.GetPdgCode();
+          if ( int(std::abs(gmpdg)/100.) == 4 || int(std::abs(gmpdg)/1000.) == 4 || int(std::abs(gmpdg)/100.) == 5 || int(std::abs(gmpdg)/1000.) == 5 ) {
+            nLeptons++;
+            nleptonse++;
+            if (-4.3 < y && y < -2.2) {
+              nleptonseinacc++;
+              nLeptonsInAcceptance++;
+            }
+            if (track.getToBeDone()) {
+              nLeptonsToBeDone++;
+              nleptonseToBeDone++;
+            }
+          }
+        }    
+      } else if (std::abs(pdg) == 411 || std::abs(pdg) == 421 || std::abs(pdg) == 431 || std::abs(pdg) == 4122 || std::abs(pdg) == 4132 || std::abs(pdg) == 4232 || std::abs(pdg) == 4332 || std::abs(pdg) == 511 || std::abs(pdg) == 521 || std::abs(pdg) == 531 || std::abs(pdg) == 541 || std::abs(pdg) == 5112 || std::abs(pdg) == 5122 || std::abs(pdg) == 5232 || std::abs(pdg) == 5132 || std::abs(pdg) == 5332) {
+        nopenHeavy++;
+      }
+    }
+    if (nopenHeavy > 1) nSignalPairs++;
+    if (nleptonse > 1) nLeptonPairs++;
+    if (nleptonseToBeDone > 1) nLeptonPairsToBeDone++;
+    if (nleptonseinacc > 1) nLeptonPairsInAcceptance++;
+  }
+  std::cout << "#events: " << nEvents << "\n"
+  << "#muons in acceptance: " << nLeptonsInAcceptance << "\n"
+  << "#muon pairs in acceptance: " << nLeptonPairsInAcceptance << "\n"
+  << "#muons: " << nLeptons << "\n"
+  << "#muons to be done: " << nLeptonsToBeDone << "\n"
+  << "#signal pairs: " << nSignalPairs << "\n"
+  << "#muon pairs: " << nLeptonPairs << "\n"
+  << "#muon pairs to be done: " << nLeptonPairsToBeDone << "\n";
+  if (nLeptonPairs != nLeptonPairsToBeDone) {
+    std::cerr << "The number of muon pairs should be the same as the number of muon pairs which should be transported.\n";
+    return 1;
+  }
+  if (nLeptons != nLeptonsToBeDone) {
+    std::cerr << "The number of muons should be the same as the number of muons which should be transported.\n";
+    return 1;
+  }
+
+  return 0;
+}

MC/config/PWGDQ/ini/tests/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_natural.C

@@ -0,0 +1,77 @@
+int External()
+{
+  int checkPdgDecay = 13;
+  std::string path{"o2sim_Kine.root"};
+  TFile file(path.c_str(), "READ");
+  if (file.IsZombie()) {
+    std::cerr << "Cannot open ROOT file " << path << "\n";
+    return 1;
+  }
+  auto tree = (TTree*)file.Get("o2sim");
+  std::vector<o2::MCTrack>* tracks{};
+  tree->SetBranchAddress("MCTrack", &tracks);
+
+  int nLeptons{};
+  int nLeptonsInAcceptance{};
+  int nLeptonsToBeDone{};
+  int nSignalPairs{};
+  int nLeptonPairs{};
+  int nLeptonPairsInAcceptance{};
+  int nLeptonPairsToBeDone{};
+  auto nEvents = tree->GetEntries();
+
+  for (int i = 0; i < nEvents; i++) {
+    tree->GetEntry(i);
+    int nleptonseinacc = 0;
+    int nleptonse = 0;
+    int nleptonseToBeDone = 0;
+    int nopenHeavy = 0;
+    for (auto& track : *tracks) {
+      auto pdg = track.GetPdgCode();
+      auto y = track.GetRapidity();
+      if (std::abs(pdg) == checkPdgDecay) {
+        int igmother = track.getMotherTrackId();
+        if (igmother > 0) {
+          auto gmTrack = (*tracks)[igmother];
+          int gmpdg = gmTrack.GetPdgCode();
+          if ( int(std::abs(gmpdg)/100.) == 4 || int(std::abs(gmpdg)/1000.) == 4 || int(std::abs(gmpdg)/100.) == 5 || int(std::abs(gmpdg)/1000.) == 5 ) {
+            nLeptons++;
+            nleptonse++;
+            if (-4.3 < y && y < -2.2) {
+              nleptonseinacc++;
+              nLeptonsInAcceptance++;
+            }
+            if (track.getToBeDone()) {
+              nLeptonsToBeDone++;
+              nleptonseToBeDone++;
+            }
+          }
+        }    
+      } else if (std::abs(pdg) == 411 || std::abs(pdg) == 421 || std::abs(pdg) == 431 || std::abs(pdg) == 4122 || std::abs(pdg) == 4132 || std::abs(pdg) == 4232 || std::abs(pdg) == 4332 || std::abs(pdg) == 511 || std::abs(pdg) == 521 || std::abs(pdg) == 531 || std::abs(pdg) == 541 || std::abs(pdg) == 5112 || std::abs(pdg) == 5122 || std::abs(pdg) == 5232 || std::abs(pdg) == 5132 || std::abs(pdg) == 5332) {
+        nopenHeavy++;
+      }
+    }
+    if (nopenHeavy > 1) nSignalPairs++;
+    if (nleptonse > 1) nLeptonPairs++;
+    if (nleptonseToBeDone > 1) nLeptonPairsToBeDone++;
+    if (nleptonseinacc > 1) nLeptonPairsInAcceptance++;
+  }
+  std::cout << "#events: " << nEvents << "\n"
+  << "#muons in acceptance: " << nLeptonsInAcceptance << "\n"
+  << "#muon pairs in acceptance: " << nLeptonPairsInAcceptance << "\n"
+  << "#muons: " << nLeptons << "\n"
+  << "#muons to be done: " << nLeptonsToBeDone << "\n"
+  << "#signal pairs: " << nSignalPairs << "\n"
+  << "#muon pairs: " << nLeptonPairs << "\n"
+  << "#muon pairs to be done: " << nLeptonPairsToBeDone << "\n";
+  if (nLeptonPairs != nLeptonPairsToBeDone) {
+    std::cerr << "The number of muon pairs should be the same as the number of muon pairs which should be transported.\n";
+    return 1;
+  }
+  if (nLeptons != nLeptonsToBeDone) {
+    std::cerr << "The number of muons should be the same as the number of muons which should be transported.\n";
+    return 1;
+  }
+
+  return 0;
+}

MC/config/PWGDQ/ini/tests/GeneratorHF_bbbarToDDbarToMuons_fwdy_TriggerGap_natural_inel.C

@@ -0,0 +1,77 @@
+int External()
+{
+  int checkPdgDecay = 13;
+  std::string path{"o2sim_Kine.root"};
+  TFile file(path.c_str(), "READ");
+  if (file.IsZombie()) {
+    std::cerr << "Cannot open ROOT file " << path << "\n";
+    return 1;
+  }
+  auto tree = (TTree*)file.Get("o2sim");
+  std::vector<o2::MCTrack>* tracks{};
+  tree->SetBranchAddress("MCTrack", &tracks);
+
+  int nLeptons{};
+  int nLeptonsInAcceptance{};
+  int nLeptonsToBeDone{};
+  int nSignalPairs{};
+  int nLeptonPairs{};
+  int nLeptonPairsInAcceptance{};
+  int nLeptonPairsToBeDone{};
+  auto nEvents = tree->GetEntries();
+
+  for (int i = 0; i < nEvents; i++) {
+    tree->GetEntry(i);
+    int nleptonseinacc = 0;
+    int nleptonse = 0;
+    int nleptonseToBeDone = 0;
+    int nopenHeavy = 0;
+    for (auto& track : *tracks) {
+      auto pdg = track.GetPdgCode();
+      auto y = track.GetRapidity();
+      if (std::abs(pdg) == checkPdgDecay) {
+        int igmother = track.getMotherTrackId();
+        if (igmother > 0) {
+          auto gmTrack = (*tracks)[igmother];
+          int gmpdg = gmTrack.GetPdgCode();
+          if ( int(std::abs(gmpdg)/100.) == 4 || int(std::abs(gmpdg)/1000.) == 4 || int(std::abs(gmpdg)/100.) == 5 || int(std::abs(gmpdg)/1000.) == 5 ) {
+            nLeptons++;
+            nleptonse++;
+            if (-4.3 < y && y < -2.2) {
+              nleptonseinacc++;
+              nLeptonsInAcceptance++;
+            }
+            if (track.getToBeDone()) {
+              nLeptonsToBeDone++;
+              nleptonseToBeDone++;
+            }
+          }
+        }    
+      } else if (std::abs(pdg) == 411 || std::abs(pdg) == 421 || std::abs(pdg) == 431 || std::abs(pdg) == 4122 || std::abs(pdg) == 4132 || std::abs(pdg) == 4232 || std::abs(pdg) == 4332 || std::abs(pdg) == 511 || std::abs(pdg) == 521 || std::abs(pdg) == 531 || std::abs(pdg) == 541 || std::abs(pdg) == 5112 || std::abs(pdg) == 5122 || std::abs(pdg) == 5232 || std::abs(pdg) == 5132 || std::abs(pdg) == 5332) {
+        nopenHeavy++;
+      }
+    }
+    if (nopenHeavy > 1) nSignalPairs++;
+    if (nleptonse > 1) nLeptonPairs++;
+    if (nleptonseToBeDone > 1) nLeptonPairsToBeDone++;
+    if (nleptonseinacc > 1) nLeptonPairsInAcceptance++;
+  }
+  std::cout << "#events: " << nEvents << "\n"
+  << "#muons in acceptance: " << nLeptonsInAcceptance << "\n"
+  << "#muon pairs in acceptance: " << nLeptonPairsInAcceptance << "\n"
+  << "#muons: " << nLeptons << "\n"
+  << "#muons to be done: " << nLeptonsToBeDone << "\n"
+  << "#signal pairs: " << nSignalPairs << "\n"
+  << "#muon pairs: " << nLeptonPairs << "\n"
+  << "#muon pairs to be done: " << nLeptonPairsToBeDone << "\n";
+  if (nLeptonPairs != nLeptonPairsToBeDone) {
+    std::cerr << "The number of muon pairs should be the same as the number of muon pairs which should be transported.\n";
+    return 1;
+  }
+  if (nLeptons != nLeptonsToBeDone) {
+    std::cerr << "The number of muons should be the same as the number of muons which should be transported.\n";
+    return 1;
+  }
+
+  return 0;
+}