Add injected nuclei generator for HMPID; extend long-lived generator with phi range selection (#1995)

Author: alcaliva

MC/config/PWGLF/ini/GeneratorLFNucleiPbPbInHMPIDAcceptance.ini

@@ -0,0 +1,6 @@
+[GeneratorExternal]
+fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_longlived.C
+funcName=generateLongLived(1000010020,10,0.1,3.0,-0.6,0.6,-0.7,1.7,-1)
+
+[GeneratorPythia8]
+config=${O2_ROOT}/share/Generators/egconfig/pythia8_hi.cfg

MC/config/PWGLF/ini/tests/GeneratorLFNucleiPbPbInHMPIDAcceptance.C

@@ -0,0 +1,53 @@
+int External()
+{
+  std::string path{"o2sim_Kine.root"};
+  std::vector<int> possiblePDGs = {1000010020, -1000010020};
+
+  int nPossiblePDGs = possiblePDGs.size();
+
+  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");
+  if (!tree)
+  {
+    std::cerr << "Cannot find tree o2sim in file " << path << "\n";
+    return 1;
+  }
+  std::vector<o2::MCTrack> *tracks{};
+  tree->SetBranchAddress("MCTrack", &tracks);
+
+  std::vector<int> injectedPDGs;
+
+  auto nEvents = tree->GetEntries();
+  for (int i = 0; i < nEvents; i++)
+  {
+    auto check = tree->GetEntry(i);
+    for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack)
+    {
+      auto track = tracks->at(idxMCTrack);
+      auto pdg = track.GetPdgCode();
+      auto it = std::find(possiblePDGs.begin(), possiblePDGs.end(), pdg);
+      if (it != possiblePDGs.end() && track.isPrimary())  // found
+      {
+        injectedPDGs.push_back(pdg);
+      }
+    }
+  }
+  std::cout << "--------------------------------\n";
+  std::cout << "# Events: " << nEvents << "\n";
+  if(injectedPDGs.empty()){
+    std::cerr << "No injected particles\n";
+    return 1; // At least one of the injected particles should be generated
+  }
+  for (int i = 0; i < nPossiblePDGs; i++)
+  {
+    std::cout << "# Injected nuclei \n";
+    std::cout << possiblePDGs[i] << ": " << std::count(injectedPDGs.begin(), injectedPDGs.end(), possiblePDGs[i]) << "\n";
+  }
+  return 0;
+}

MC/config/PWGLF/pythia8/generator_pythia8_longlived.C

@@ -16,7 +16,7 @@ class GeneratorPythia8LongLivedGun : public o2::eventgen::GeneratorPythia8
 {
 public:
   /// constructor
-  GeneratorPythia8LongLivedGun(int input_pdg, int nInject = 1, float ptMin = 1, float ptMax = 10, int input_pdg2 = -1) : pdg{input_pdg}, nParticles{nInject}, genMinPt{ptMin}, genMaxPt{ptMax}, m{getMass(input_pdg)}, pdg2{input_pdg2}
+  GeneratorPythia8LongLivedGun(int input_pdg, int nInject = 1, float ptMin = 1, float ptMax = 10, float etaMin = -1.0, float etaMax = 1.0, float phiMin = 0.0, float phiMax = TMath::Pi(), int input_pdg2 = -1) : pdg{input_pdg}, nParticles{nInject}, genMinPt{ptMin}, genMaxPt{ptMax}, genMinEta{etaMin}, genMaxEta{etaMax}, genMinPhi{phiMin}, genMaxPhi{phiMax},  m{getMass(input_pdg)}, pdg2{input_pdg2}
   {
   }
 
@@ -51,7 +51,7 @@ public:
     {
       const double pt = gRandom->Uniform(genMinPt, genMaxPt);
       const double eta = gRandom->Uniform(genMinEta, genMaxEta);
-      const double phi = gRandom->Uniform(0, TMath::TwoPi());
+      const double phi = gRandom->Uniform(genMinPhi, genMaxPhi);
       const double px{pt * std::cos(phi)};
       const double py{pt * std::sin(phi)};
       const double pz{pt * std::sinh(eta)};
@@ -69,7 +69,7 @@ public:
       {
         const double pt = gRandom->Uniform(genMinPt, genMaxPt);
         const double eta = gRandom->Uniform(genMinEta, genMaxEta);
-        const double phi = gRandom->Uniform(0, TMath::TwoPi());
+        const double phi = gRandom->Uniform(genMinPhi, genMaxPhi);
         const double px{pt * std::cos(phi)};
         const double py{pt * std::sin(phi)};
         const double pz{pt * std::sinh(eta)};
@@ -88,8 +88,10 @@ public:
 private:
   double genMinPt = 0.5;  /// minimum 3-momentum for generated particles
   double genMaxPt = 12.;  /// maximum 3-momentum for generated particles
-  double genMinEta = -1.; /// minimum pseudorapidity for generated particles
-  double genMaxEta = +1.; /// maximum pseudorapidity for generated particles
+  double genMinEta = -1.0; /// minimum pseudorapidity for generated particles
+  double genMaxEta = +1.0; /// maximum pseudorapidity for generated particles
+  double genMinPhi = 0.0; /// minimum pseudorapidity for generated particles
+  double genMaxPhi = TMath::Pi(); /// maximum pseudorapidity for generated particles
 
   double m = 0;       /// particle mass [GeV/c^2]
   int pdg = 0;        /// particle pdg code
@@ -101,7 +103,8 @@ private:
 };
 
 ///___________________________________________________________
-FairGenerator *generateLongLived(int pdg, int nInject, float ptMin = 1, float ptMax = 10, int pdg2 = -1)
+FairGenerator *generateLongLived(int pdg, int nInject, float ptMin = 1, float ptMax = 10, float etaMin = -1.0, float etaMax = 1.0, float phiMin = 0.0, float phiMax = TMath::Pi(), int pdg2 = -1)
 {
-  return new GeneratorPythia8LongLivedGun(pdg, nInject, ptMin, ptMax, pdg2);
+  return new GeneratorPythia8LongLivedGun(pdg, nInject, ptMin, ptMax, etaMin, etaMax, phiMin, phiMax, pdg2);
 }
+