feat: calculate PhPerp for SIDIS hadrons and dihadrons (#304)

JeffersonLab · Oct 28, 2024 · 21a6465 · 21a6465
1 parent 026c57a
commit 21a6465
Show file tree

Hide file tree

Showing 6 changed files with 43 additions and 14 deletions.
diff --git a/src/iguana/algorithms/physics/DihadronKinematics/Algorithm.cc b/src/iguana/algorithms/physics/DihadronKinematics/Algorithm.cc
@@ -25,6 +25,7 @@ namespace iguana::physics {
           "pdg_b/I",
           "Mh/D",
           "z/D",
+          "PhPerp/D",
           "MX/D",
           "xF/D",
           "phiH/D",
@@ -39,6 +40,7 @@ namespace iguana::physics {
     i_pdg_b    = result_schema.getEntryOrder("pdg_b");
     i_Mh       = result_schema.getEntryOrder("Mh");
     i_z        = result_schema.getEntryOrder("z");
+    i_PhPerp   = result_schema.getEntryOrder("PhPerp");
     i_MX       = result_schema.getEntryOrder("MX");
     i_xF       = result_schema.getEntryOrder("xF");
     i_phiH     = result_schema.getEntryOrder("phiH");
@@ -135,6 +137,10 @@ namespace iguana::physics {
       // calculate z
       double z = p_target.Dot(p_Ph) / p_target.Dot(p_q);
 
+      // calculate PhPerp
+      auto opt_PhPerp = tools::RejectVector(p_Ph.Vect(), p_q.Vect());
+      double PhPerp   = opt_PhPerp.has_value() ? opt_PhPerp.value().R() : tools::UNDEF;
+
       // calculate Mh
       double Mh = p_Ph.M();
 
@@ -190,6 +196,7 @@ namespace iguana::physics {
       result_bank.putInt(i_pdg_b,      dih_row, had_b.pdg);
       result_bank.putDouble(i_Mh,      dih_row, Mh);
       result_bank.putDouble(i_z,       dih_row, z);
+      result_bank.putDouble(i_PhPerp,  dih_row, PhPerp);
       result_bank.putDouble(i_MX,      dih_row, MX);
       result_bank.putDouble(i_xF,      dih_row, xF);
       result_bank.putDouble(i_phiH,    dih_row, phiH);

diff --git a/src/iguana/algorithms/physics/DihadronKinematics/Algorithm.h b/src/iguana/algorithms/physics/DihadronKinematics/Algorithm.h
@@ -21,6 +21,8 @@ namespace iguana::physics {
     double Mh;
     /// @brief @latex{z}: Momentum fraction of the fragmenting parton carried by the dihadron
     double z;
+    /// @brief @latex{P_h^\perp}: transverse momentum of the dihadron in the @latex{\perp}-frame (transverse to @latex{\vec{q}})
+    double PhPerp;
     /// @brief @latex{M_X(ehhX)}: Missing mass of the dihadron
     double MX;
     /// @brief @latex{x_F}: Feynman-x of the dihadron
@@ -95,6 +97,7 @@ namespace iguana::physics {
       int i_pdg_b;
       int i_Mh;
       int i_z;
+      int i_PhPerp;
       int i_MX;
       int i_xF;
       int i_phiH;

diff --git a/src/iguana/algorithms/physics/DihadronKinematics/Validator.cc b/src/iguana/algorithms/physics/DihadronKinematics/Validator.cc
@@ -41,6 +41,10 @@ namespace iguana::physics {
         new TH1D("z_dist", "z", n_bins, 0, 1),
         [](auto const& b, auto const r) { return b.getDouble("z", r); }
       },
+      {
+        new TH1D("PhPerp_dist", "P_{h}^{{}^{#perp}}", n_bins, 0, 2),
+        [](auto const& b, auto const r) { return b.getDouble("PhPerp", r); }
+      },
       {
         new TH1D("MX_dist", "missing mass M_{X} [GeV];", n_bins, 0, 4),
         [](auto const& b, auto const r) { return b.getDouble("MX", r); }

diff --git a/src/iguana/algorithms/physics/SingleHadronKinematics/Algorithm.cc b/src/iguana/algorithms/physics/SingleHadronKinematics/Algorithm.cc
@@ -23,6 +23,7 @@ namespace iguana::physics {
           "pindex/S",
           "pdg/I",
           "z/D",
+          "PhPerp/D",
           "MX/D",
           "xF/D",
           "phiH/D",
@@ -33,6 +34,7 @@ namespace iguana::physics {
     i_pindex = result_schema.getEntryOrder("pindex");
     i_pdg    = result_schema.getEntryOrder("pdg");
     i_z      = result_schema.getEntryOrder("z");
+    i_PhPerp = result_schema.getEntryOrder("PhPerp");
     i_MX     = result_schema.getEntryOrder("MX");
     i_xF     = result_schema.getEntryOrder("xF");
     i_phiH   = result_schema.getEntryOrder("phiH");
@@ -112,6 +114,10 @@ namespace iguana::physics {
         // calculate z
         double z = p_target.Dot(p_Ph) / p_target.Dot(p_q);
 
+        // calculate PhPerp
+        auto opt_PhPerp = tools::RejectVector(p_Ph.Vect(), p_q.Vect());
+        double PhPerp   = opt_PhPerp.has_value() ? opt_PhPerp.value().R() : tools::UNDEF;
+
         // calculate xi
         double xi = p_q.Dot(p_Ph) / p_target.Dot(p_q);
 
@@ -132,23 +138,25 @@ namespace iguana::physics {
         result_bank_rowlist.push_back(row);
 
         // fill the bank
-        result_bank.putShort(i_pindex, row, static_cast<int16_t>(row));
-        result_bank.putInt(i_pdg,      row, pdg);
-        result_bank.putDouble(i_z,     row, z);
-        result_bank.putDouble(i_MX,    row, MX);
-        result_bank.putDouble(i_xF,    row, xF);
-        result_bank.putDouble(i_phiH,  row, phiH);
-        result_bank.putDouble(i_xi,    row, xi);
+        result_bank.putShort(i_pindex,  row, static_cast<int16_t>(row));
+        result_bank.putInt(i_pdg,       row, pdg);
+        result_bank.putDouble(i_z,      row, z);
+        result_bank.putDouble(i_PhPerp, row, PhPerp);
+        result_bank.putDouble(i_MX,     row, MX);
+        result_bank.putDouble(i_xF,     row, xF);
+        result_bank.putDouble(i_phiH,   row, phiH);
+        result_bank.putDouble(i_xi,     row, xi);
       }
       else {
         // zero the row
-        result_bank.putShort(i_pindex, row, static_cast<int16_t>(row));
-        result_bank.putInt(i_pdg,      row, pdg);
-        result_bank.putDouble(i_z,     row, 0);
-        result_bank.putDouble(i_MX,    row, 0);
-        result_bank.putDouble(i_xF,    row, 0);
-        result_bank.putDouble(i_phiH,  row, 0);
-        result_bank.putDouble(i_xi,    row, 0);
+        result_bank.putShort(i_pindex,  row, static_cast<int16_t>(row));
+        result_bank.putInt(i_pdg,       row, pdg);
+        result_bank.putDouble(i_z,      row, 0);
+        result_bank.putDouble(i_PhPerp, row, 0);
+        result_bank.putDouble(i_MX,     row, 0);
+        result_bank.putDouble(i_xF,     row, 0);
+        result_bank.putDouble(i_phiH,   row, 0);
+        result_bank.putDouble(i_xi,     row, 0);
       }
     }
 

diff --git a/src/iguana/algorithms/physics/SingleHadronKinematics/Algorithm.h b/src/iguana/algorithms/physics/SingleHadronKinematics/Algorithm.h
@@ -14,6 +14,8 @@ namespace iguana::physics {
     int pdg;
     /// @brief @latex{z}: Momentum fraction of the fragmenting parton carried by the hadron
     double z;
+    /// @brief @latex{P_h^\perp}: transverse momentum of the hadron in the @latex{\perp}-frame (transverse to @latex{\vec{q}})
+    double PhPerp;
     /// @brief @latex{M_X(ehX)}: Missing mass of the hadron
     double MX;
     /// @brief @latex{x_F}: Feynman-x of the hadron
@@ -66,6 +68,7 @@ namespace iguana::physics {
       int i_pindex;
       int i_pdg;
       int i_z;
+      int i_PhPerp;
       int i_MX;
       int i_xF;
       int i_phiH;

diff --git a/src/iguana/algorithms/physics/SingleHadronKinematics/Validator.cc b/src/iguana/algorithms/physics/SingleHadronKinematics/Validator.cc
@@ -36,6 +36,10 @@ namespace iguana::physics {
         new TH1D("z_dist", "z", n_bins, 0, 1),
         [](auto const& b, auto const r) { return b.getDouble("z", r); }
       },
+      {
+        new TH1D("PhPerp_dist", "P_{h}^{{}^{#perp}}", n_bins, 0, 2),
+        [](auto const& b, auto const r) { return b.getDouble("PhPerp", r); }
+      },
       {
         new TH1D("MX_dist", "missing mass M_{X} [GeV];", n_bins, 0, 4),
         [](auto const& b, auto const r) { return b.getDouble("MX", r); }