[PWGCF] FemtoUniverse: Add PRF boost to identical pairs, remove qTmin selection (#14615)

Co-authored-by: ALICE Action Bot <[email protected]>

Author: prchakra

PWGCF/FemtoUniverse/Core/FemtoUniverse3DContainer.h

@@ -216,22 +216,36 @@ class FemtoUniverse3DContainer
   /// \param isiden Choosing identical or non-identical pairs
   /// \param islcm Choosing LCMS or PRF
   template <bool isMC, typename T>
-  void setPair(T const& part1, T const& part2, const int mult, bool use3dplots, bool isiden)
+  void setPair(T const& part1, T const& part2, const int mult, bool use3dplots, bool isIdenLCMS, bool isWeight, bool isIdenPRF)
   {
     std::vector<double> f3d;
     const float kT = FemtoUniverseMath::getkT(part1, mMassOne, part2, mMassTwo);
     const float mT = FemtoUniverseMath::getmT(part1, mMassOne, part2, mMassTwo);
 
-    f3d = FemtoUniverseMath::newpairfunc(part1, mMassOne, part2, mMassTwo, isiden);
+    f3d = FemtoUniverseMath::newpairfunc(part1, mMassOne, part2, mMassTwo, isIdenLCMS, isWeight, isIdenPRF);
+
+    float varout = 0.0;
+    float varside = 0.0;
+    float varlong = 0.0;
+
+    if (isIdenLCMS) {
+      varout = f3d[1];
+      varside = f3d[2];
+      varlong = f3d[3];
+    } else if (isIdenPRF) {
+      varout = f3d[6];
+      varside = f3d[7];
+      varlong = f3d[8];
+    }
 
     const float femtoObs1D = f3d[0];
-    const float femtoObsKout = f3d[1];
-    const float femtoObsKside = f3d[2];
-    const float femtoObsKlong = f3d[3];
+    const float femtoObsKout = varout;
+    const float femtoObsKside = varside;
+    const float femtoObsKlong = varlong;
 
     if (mHistogramRegistry) {
-      setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObsKout, femtoObsKside, femtoObsKlong, femtoObs1D, kT, mT, part1, part2, mult, use3dplots, isiden);
-      if (!isiden) {
+      setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObsKout, femtoObsKside, femtoObsKlong, femtoObs1D, kT, mT, part1, part2, mult, use3dplots, isIdenLCMS);
+      if (!isIdenLCMS) {
         if (femtoObsKout > 0.0) {
           mHistogramRegistry->fill(HIST(FolderSuffix[EventType]) + HIST(o2::aod::femtouniverse_mc_particle::MCTypeName[o2::aod::femtouniverse_mc_particle::MCType::kRecon]) + HIST("/KStarOutP"), femtoObs1D);
         } else {

PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h

@@ -177,7 +177,7 @@ class FemtoUniverseMath
   /// \param mass2 Mass of particle 2
   /// \param isiden Identical or non-identical particle pair
   template <typename T>
-  static std::vector<double> newpairfunc(const T& part1, const float mass1, const T& part2, const float mass2, bool isiden, bool isWeight = 0)
+  static std::vector<double> newpairfunc(const T& part1, const float mass1, const T& part2, const float mass2, bool isIdenLCMS, bool isWeight, bool isIdenPRF)
   {
     const double e1 = std::sqrt(std::pow(part1.px(), 2) + std::pow(part1.py(), 2) + std::pow(part1.pz(), 2) + std::pow(mass1, 2));
     const double e2 = std::sqrt(std::pow(part2.px(), 2) + std::pow(part2.py(), 2) + std::pow(part2.pz(), 2) + std::pow(mass2, 2));
@@ -240,7 +240,7 @@ class FemtoUniverseMath
     const double qinv = std::sqrt(fDKOutPRF * fDKOutPRF + fDKSidePRF * fDKSidePRF + fDKLongPRF * fDKLongPRF);
     const double kstar = std::sqrt(fKOut * fKOut + fDKSide * fDKSide + fDKLong * fDKLong);
 
-    if (isiden) {
+    if (isIdenLCMS) {
       vect.push_back(qinv);
       vect.push_back(fDKOutLCMS);
       vect.push_back(fDKSideLCMS);
@@ -251,9 +251,10 @@ class FemtoUniverseMath
       vect.push_back(fDKOut);
       vect.push_back(fDKSide);
       vect.push_back(fDKLong);
+      vect.push_back(1.0);
     }
 
-    if (isiden && isWeight) {
+    if (isIdenLCMS && isWeight) {
       const double x1_lcms = (vecspace_part1.x() * tPx + vecspace_part1.y() * tPy) / tPt;
       const double y1_lcms = (-vecspace_part1.x() * tPy + vecspace_part1.y() * tPx) / tPt;
       const double x2_lcms = (vecspace_part2.x() * tPx + vecspace_part2.y() * tPy) / tPt;
@@ -277,6 +278,13 @@ class FemtoUniverseMath
     } else {
       vect.push_back(1.0);
     }
+
+    if (isIdenPRF) {
+      vect.push_back(fDKOutPRF);
+      vect.push_back(fDKSidePRF);
+      vect.push_back(fDKLongPRF);
+    }
+
     return vect;
   }
 };

PWGCF/FemtoUniverse/Core/FemtoUniversePairSHCentMultKt.h

@@ -180,7 +180,7 @@ class PairSHCentMultKt
   /// \param ktval kT value
   template <typename T>
   void fillMultNumDen(T const& part1, T const& part2, uint8_t ChosenEventType,
-                      int maxl, int multval, float ktval, bool isiden, bool isqinvfill)
+                      int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
   {
     int multbinval;
     int absmultval = multval;
@@ -197,7 +197,7 @@ class PairSHCentMultKt
       return;
     }
     // std::cout<<"multbinval "<<multbinval<<std::endl;
-    fillkTNumDen(part1, part2, ChosenEventType, maxl, multbinval, ktval, isiden, isqinvfill);
+    fillkTNumDen(part1, part2, ChosenEventType, maxl, multbinval, ktval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF);
   }
 
   /// Templated function to access different kT directory and call addEventPair
@@ -209,7 +209,7 @@ class PairSHCentMultKt
   /// \param ktval kT value
   template <typename T>
   void fillkTNumDen(T const& part1, T const& part2, uint8_t ChosenEventType,
-                    int maxl, int multval, float ktval, bool isiden, bool isqinvfill)
+                    int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
   {
     int ktbinval = -1;
     if (ktval >= ktBins[0] && ktval < ktBins[1]) {
@@ -229,7 +229,7 @@ class PairSHCentMultKt
     } else {
       return;
     }
-    addEventPair(part1, part2, ChosenEventType, maxl, multval, ktbinval, isiden, isqinvfill);
+    addEventPair(part1, part2, ChosenEventType, maxl, multval, ktbinval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF);
   }
 
   /// Set the PDG codes of the two particles involved
@@ -262,20 +262,33 @@ class PairSHCentMultKt
   /// \param ktval kT value
   template <typename T>
   void addEventPair(T const& part1, T const& part2, uint8_t ChosenEventType,
-                    int /*maxl*/, int multval, int ktval, bool isiden, bool isqinvfill)
+                    int /*maxl*/, int multval, int ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
   {
     int fMultBin = multval;
     int fKtBin = ktval;
     std::vector<std::complex<double>> fYlmBuffer(kMaxJM);
     std::vector<double> f3d;
     setPionPairMass();
-    f3d = FemtoUniverseMath::newpairfunc(part1, mMassOne, part2, mMassTwo,
-                                         isiden);
 
-    const float qout = f3d[1];
-    const float qside = f3d[2];
-    const float qlong = f3d[3];
+    f3d = FemtoUniverseMath::newpairfunc(part1, mMassOne, part2, mMassTwo,
+                                         isIdenLCMS, isWeight, isIdenPRF);
+    double varout = 0.0;
+    double varside = 0.0;
+    double varlong = 0.0;
+
+    if (isIdenLCMS) {
+      varout = f3d[1];
+      varside = f3d[2];
+      varlong = f3d[3];
+    } else if (isIdenPRF) {
+      varout = f3d[6];
+      varside = f3d[7];
+      varlong = f3d[8];
+    }
 
+    const double qout = varout;
+    const double qside = varside;
+    const double qlong = varlong;
     double kv = std::sqrt(qout * qout + qside * qside + qlong * qlong);
 
     // int nqbin = fbinctn[0][0]->GetXaxis()->FindFixBin(kv);