Re-enable DIS

benchmarks/dis/analysis/dis_electrons.cxx

@@ -23,6 +23,16 @@
 #include "nlohmann/json.hpp"
 #include "eicd/ReconstructedParticleData.h"
 
+// Select only electrons
+inline auto select_electrons(const std::vector<eic::ReconstructedParticleData>& parts)
+{
+  std::vector<eic::ReconstructedParticleData> electrons;
+  for (auto &part: parts)
+    if (parts.pid == 11)
+      electrons.push_back(parts);
+  return electrons;
+}
+
 // Get a vector of 4-momenta from the reconstructed data.
 inline auto momenta_from_reconstruction(const std::vector<eic::ReconstructedParticleData>& parts)
 {
@@ -62,7 +72,7 @@ inline auto sort_momenta(const std::vector<ROOT::Math::PxPyPzEVector>& mom)
 inline auto Q2(const std::vector<ROOT::Math::PxPyPzEVector>& mom)
 {
   std::vector<double> Q2Vec(mom.size() );
-  ROOT::Math::PxPyPzEVector beamMom = {0, 0, 18, 18};
+  ROOT::Math::PxPyPzEVector beamMom = {0, 0, -18, 18};
   std::transform(mom.begin(), mom.end(), Q2Vec.begin(), [beamMom](const auto& part) {
     return -(part - beamMom).M2();
   });
@@ -203,23 +213,29 @@ int dis_electrons(const std::string& config_name)
                 .Define("mom_rec", common_bench::mom, {"p_rec"});
 */
 
-  auto d0 = d.Define("p_recon", momenta_from_reconstruction, {"ReconstructedParticles"})
-                .Define("p_recon_sort", sort_momenta, {"p_recon"})
-                .Define("Q2_recon", Q2, {"p_recon_sort"})
-                .Define("Q2_recon_rand", randomize, {"Q2_recon"})
-                .Define("elec_Q2_recon_rand", "Q2_recon_rand[0]")
-                .Define("p_sim_M", common_bench::momenta_from_simulation, {"mcparticles2"})
-                .Define("p_sim", convertMtoE, {"p_sim_M"})
+  auto d0 = d   .Define("electrons_rec", select_electrons, {"ReconstructedParticles"})
+                .Define("p_rec", momenta_from_reconstruction, {"electrons_rec"})
+                .Define("p_rec_sort", sort_momenta, {"p_rec"})
+                .Define("Q2_rec", Q2, {"p_rec_sort"})
+                .Define("Q2_rec_rand", randomize, {"Q2_rec"})
+                .Define("elec_Q2_rec_rand", "Q2_rec_rand[0]")
+
+                .Define("electrons_sim", select_electrons, {"GeneratedParticles"})
+                .Define("p_sim", momenta_from_reconstruction, {"electrons_sim"})
                 .Define("p_sim_sort", sort_momenta, {"p_sim"})
                 .Define("Q2_sim", Q2, {"p_sim_sort"})
                 .Define("elec_Q2_sim", "Q2_sim[0]")
-                .Define("Q2_diff", "(elec_Q2_recon_rand - elec_Q2_sim)/elec_Q2_sim")
-                .Define("p_diff", sub, {"p_recon","p_sim"});
+
+                .Define("Q2_diff", "(elec_Q2_rec_rand - elec_Q2_sim)/elec_Q2_sim")
+                .Define("p_diff", sub, {"p_rec","p_sim"})
+                .Define("mom_rec", common_bench::mom, {"p_rec"})
+                .Define("mom_sim", common_bench::mom, {"p_sim"})
                 /*
                 .Define("electrons_sim", elec_PID_sim, {"sorted_sim"})
                 .Define("Q2_sim_elec_pid", Q2_from_sim, {"electrons_sim"})
                 .Define("elec_Q2_sim_pid", "Q2_sim_elec_pid[0]");
                 */
+                ;
 
   //Momentum
   auto h_mom_sim = d0.Histo1D({"h_mom_sim", "; GeV; counts", 100, 0, 50}, "mom_sim");
@@ -227,7 +243,7 @@ int dis_electrons(const std::string& config_name)
   auto h_mom_res = d0.Histo1D({"h_mom_res", ";    ; counts", 100, -10, 10}, "p_diff");
   //Q2
   auto h_Q2_sim = d0.Histo1D({"h_Q2_sim", "; GeV; counts", 100, -5, 25}, "elec_Q2_sim");
-  auto h_Q2_rec = d0.Histo1D({"h_Q2_rec", "; GeV; counts", 100, -5, 25}, "elec_Q2_recon_rand");
+  auto h_Q2_rec = d0.Histo1D({"h_Q2_rec", "; GeV; counts", 100, -5, 25}, "elec_Q2_rec_rand");
   auto h_Q2_res = d0.Histo1D({"h_Q2_res", ";    ; counts", 100, -10, 10}, "Q2_diff");
 
   /*