From 69a7be56252e8875f7012831c54a73cbcec3fa0b Mon Sep 17 00:00:00 2001 From: Wouter Deconinck <wouter.deconinck@umanitoba.ca> Date: Sat, 23 Oct 2021 02:40:00 +0000 Subject: [PATCH] Inclusive Kinematics JB DA --- benchmarks/dis/analysis/dis_electrons.cxx | 86 ++++++++++++++++++----- options/reconstruction.py | 20 +++++- 2 files changed, 86 insertions(+), 20 deletions(-) diff --git a/benchmarks/dis/analysis/dis_electrons.cxx b/benchmarks/dis/analysis/dis_electrons.cxx index 8f89147f..e59c69a5 100644 --- a/benchmarks/dis/analysis/dis_electrons.cxx +++ b/benchmarks/dis/analysis/dis_electrons.cxx @@ -73,26 +73,42 @@ int dis_electrons(const std::string& config_name) }; auto d0 = d.Define("Q2_sim", "InclusiveKinematicsTruth.Q2") - .Define("Q2_rec", "InclusiveKinematicsElectron.Q2") - .Define("Q2_res", combinatorial_diff_ratio, {"Q2_sim", "Q2_rec"}) + .Define("Q2_el", "InclusiveKinematicsElectron.Q2") + .Define("Q2_jb", "InclusiveKinematicsJB.Q2") + .Define("Q2_da", "InclusiveKinematicsDA.Q2") + .Define("Q2_el_res", combinatorial_diff_ratio, {"Q2_sim", "Q2_el"}) + .Define("Q2_jb_res", combinatorial_diff_ratio, {"Q2_sim", "Q2_jb"}) + .Define("Q2_da_res", combinatorial_diff_ratio, {"Q2_sim", "Q2_da"}) .Define("x_sim", "InclusiveKinematicsTruth.x") - .Define("x_rec", "InclusiveKinematicsElectron.x") - .Define("x_res", combinatorial_diff_ratio, {"x_sim", "x_rec"}) + .Define("x_el", "InclusiveKinematicsElectron.x") + .Define("x_jb", "InclusiveKinematicsJB.x") + .Define("x_da", "InclusiveKinematicsDA.x") + .Define("x_el_res", combinatorial_diff_ratio, {"x_sim", "x_el"}) + .Define("x_jb_res", combinatorial_diff_ratio, {"x_sim", "x_jb"}) + .Define("x_da_res", combinatorial_diff_ratio, {"x_sim", "x_da"}) ; //Q2 auto h_Q2_sim = d0.Histo1D({"h_Q2_sim", "; GeV^2; counts", 100, -5, 25}, "Q2_sim"); - auto h_Q2_rec = d0.Histo1D({"h_Q2_rec", "; GeV^2; counts", 100, -5, 25}, "Q2_rec"); - auto h_Q2_res = d0.Histo1D({"h_Q2_res", "; ; counts", 100, -1, 1}, "Q2_res"); + auto h_Q2_el = d0.Histo1D({"h_Q2_el", "; GeV^2; counts", 100, -5, 25}, "Q2_el"); + auto h_Q2_jb = d0.Histo1D({"h_Q2_jb", "; GeV^2; counts", 100, -5, 25}, "Q2_jb"); + auto h_Q2_da = d0.Histo1D({"h_Q2_da", "; GeV^2; counts", 100, -5, 25}, "Q2_da"); + auto h_Q2_el_res = d0.Histo1D({"h_Q2_el_res", "; ; counts", 100, -1, 1}, "Q2_el_res"); + auto h_Q2_jb_res = d0.Histo1D({"h_Q2_jb_res", "; ; counts", 100, -1, 1}, "Q2_jb_res"); + auto h_Q2_da_res = d0.Histo1D({"h_Q2_da_res", "; ; counts", 100, -1, 1}, "Q2_da_res"); //x auto h_x_sim = d0.Histo1D({"h_x_sim", "; ; counts", 100, 0, +1}, "x_sim"); - auto h_x_rec = d0.Histo1D({"h_x_rec", "; ; counts", 100, 0, +1}, "x_rec"); - auto h_x_res = d0.Histo1D({"h_x_res", "; ; counts", 100, -1, 1}, "x_res"); + auto h_x_el = d0.Histo1D({"h_x_el", "; ; counts", 100, 0, +1}, "x_el"); + auto h_x_jb = d0.Histo1D({"h_x_jb", "; ; counts", 100, 0, +1}, "x_jb"); + auto h_x_da = d0.Histo1D({"h_x_da", "; ; counts", 100, 0, +1}, "x_da"); + auto h_x_el_res = d0.Histo1D({"h_x_el_res", "; ; counts", 100, -1, 1}, "x_el_res"); + auto h_x_jb_res = d0.Histo1D({"h_x_jb_res", "; ; counts", 100, -1, 1}, "x_jb_res"); + auto h_x_da_res = d0.Histo1D({"h_x_da_res", "; ; counts", 100, -1, 1}, "x_da_res"); - TFitResultPtr f_Q2_res = h_Q2_res->Fit("gaus", "S"); - if (f_Q2_res == 0) f_Q2_res->Print("V"); - TFitResultPtr f_x_res = h_x_res->Fit("gaus", "S"); - if (f_x_res == 0) f_x_res->Print("V"); + TFitResultPtr f_Q2_el_res = h_Q2_el_res->Fit("gaus", "S"); + if (f_Q2_el_res == 0) f_Q2_el_res->Print("V"); + TFitResultPtr f_x_el_res = h_x_el_res->Fit("gaus", "S"); + if (f_x_el_res == 0) f_x_el_res->Print("V"); // Plot our histograms. // TODO: to start I'm explicitly plotting the histograms, but want to @@ -105,18 +121,26 @@ int dis_electrons(const std::string& config_name) gPad->SetLogx(false); gPad->SetLogy(true); auto& h1 = *h_Q2_sim; - auto& h2 = *h_Q2_rec; + auto& h2 = *h_Q2_el; + auto& h3 = *h_Q2_jb; + auto& h4 = *h_Q2_da; // histogram style h1.SetLineColor(common_bench::plot::kMpBlue); h1.SetLineWidth(2); h2.SetLineColor(common_bench::plot::kMpOrange); h2.SetLineWidth(2); + h3.SetLineColor(common_bench::plot::kMpRed); + h3.SetLineWidth(2); + h4.SetLineColor(common_bench::plot::kMpGreen); + h4.SetLineWidth(2); // axes h1.GetXaxis()->CenterTitle(); h1.GetYaxis()->CenterTitle(); // draw everything h1.DrawClone("hist"); h2.DrawClone("hist same"); + h3.DrawClone("hist same"); + h4.DrawClone("hist same"); // FIXME hardcoded beam configuration common_bench::plot::draw_label(18, 275, detector); TText* tptr1; @@ -126,8 +150,12 @@ int dis_electrons(const std::string& config_name) t1.SetTextSize(25); tptr1 = t1.AddText("simulated"); tptr1->SetTextColor(common_bench::plot::kMpBlue); - tptr1 = t1.AddText("reconstructed"); + tptr1 = t1.AddText("electron"); tptr1->SetTextColor(common_bench::plot::kMpOrange); + tptr1 = t1.AddText("JB"); + tptr1->SetTextColor(common_bench::plot::kMpRed); + tptr1 = t1.AddText("DA"); + tptr1->SetTextColor(common_bench::plot::kMpGreen); t1.Draw(); c.Print(fmt::format("{}Q2.png", output_prefix).c_str()); } @@ -138,15 +166,23 @@ int dis_electrons(const std::string& config_name) c.cd(); gPad->SetLogx(false); gPad->SetLogy(true); - auto& h1 = *h_Q2_res; + auto& h1 = *h_Q2_el_res; + auto& h2 = *h_Q2_jb_res; + auto& h3 = *h_Q2_da_res; // histogram style h1.SetLineColor(common_bench::plot::kMpBlue); h1.SetLineWidth(2); + h2.SetLineColor(common_bench::plot::kMpOrange); + h2.SetLineWidth(2); + h3.SetLineColor(common_bench::plot::kMpRed); + h3.SetLineWidth(2); // axes h1.GetXaxis()->CenterTitle(); h1.GetYaxis()->CenterTitle(); // draw everything h1.DrawClone("hist"); + h2.DrawClone("hist same"); + h3.DrawClone("hist same"); // FIXME hardcoded beam configuration common_bench::plot::draw_label(18, 275, detector); c.Print(fmt::format("{}Q2resolution.png", output_prefix).c_str()); @@ -159,18 +195,26 @@ int dis_electrons(const std::string& config_name) gPad->SetLogx(true); gPad->SetLogy(true); auto& h1 = *h_x_sim; - auto& h2 = *h_x_rec; + auto& h2 = *h_x_el; + auto& h3 = *h_x_jb; + auto& h4 = *h_x_da; // histogram style h1.SetLineColor(common_bench::plot::kMpBlue); h1.SetLineWidth(2); h2.SetLineColor(common_bench::plot::kMpOrange); h2.SetLineWidth(2); + h3.SetLineColor(common_bench::plot::kMpRed); + h3.SetLineWidth(2); + h4.SetLineColor(common_bench::plot::kMpGreen); + h4.SetLineWidth(2); // axes h1.GetXaxis()->CenterTitle(); h1.GetYaxis()->CenterTitle(); // draw everything h1.DrawClone("hist"); h2.DrawClone("hist same"); + h3.DrawClone("hist same"); + h4.DrawClone("hist same"); // FIXME hardcoded beam configuration common_bench::plot::draw_label(18, 275, detector); TText* tptr1; @@ -192,15 +236,23 @@ int dis_electrons(const std::string& config_name) c.cd(); gPad->SetLogx(false); gPad->SetLogy(true); - auto& h1 = *h_x_res; + auto& h1 = *h_x_el_res; + auto& h2 = *h_x_jb_res; + auto& h3 = *h_x_da_res; // histogram style h1.SetLineColor(common_bench::plot::kMpBlue); h1.SetLineWidth(2); + h2.SetLineColor(common_bench::plot::kMpOrange); + h2.SetLineWidth(2); + h3.SetLineColor(common_bench::plot::kMpRed); + h3.SetLineWidth(2); // axes h1.GetXaxis()->CenterTitle(); h1.GetYaxis()->CenterTitle(); // draw everything h1.DrawClone("hist"); + h2.DrawClone("hist same"); + h3.DrawClone("hist same"); // FIXME hardcoded beam configuration common_bench::plot::draw_label(18, 275, detector); c.Print(fmt::format("{}xresolution.png", output_prefix).c_str()); diff --git a/options/reconstruction.py b/options/reconstruction.py index 758c8649..12e59c26 100644 --- a/options/reconstruction.py +++ b/options/reconstruction.py @@ -93,6 +93,8 @@ from Configurables import Jug__Reco__CKFTracking as CKFTracking from Configurables import Jug__Reco__ParticlesFromTrackFit as ParticlesFromTrackFit from Configurables import Jug__Reco__TrajectoryFromTrackFit as TrajectoryFromTrackFit from Configurables import Jug__Reco__InclusiveKinematicsElectron as InclusiveKinematicsElectron +from Configurables import Jug__Reco__InclusiveKinematicsDA as InclusiveKinematicsDA +from Configurables import Jug__Reco__InclusiveKinematicsJB as InclusiveKinematicsJB from Configurables import Jug__Reco__FarForwardParticles as FFRecoRP from Configurables import Jug__Reco__FarForwardParticlesOMD as FFRecoOMD @@ -679,13 +681,25 @@ if 'acadia' in detector_version: outputHitCollection="MRICHRecHits") algorithms.append(mrich_reco) -# Electron kinematics -electron_incl_kin = InclusiveKinematicsElectron("electron_incl_kin", +# Inclusive kinematics +incl_kin_electron = InclusiveKinematicsElectron("incl_kin_electron", inputMCParticles="mcparticles", inputParticles="ReconstructedParticles", outputData="InclusiveKinematicsElectron" ) -algorithms.append(electron_incl_kin) +algorithms.append(incl_kin_electron) +incl_kin_jb = InclusiveKinematicsJB("incl_kin_jb", + inputMCParticles="mcparticles", + inputParticles="ReconstructedParticles", + outputData="InclusiveKinematicsJB" +) +algorithms.append(incl_kin_jb) +incl_kin_da = InclusiveKinematicsDA("incl_kin_da", + inputMCParticles="mcparticles", + inputParticles="ReconstructedParticles", + outputData="InclusiveKinematicsDA" +) +algorithms.append(incl_kin_da) # Output podout = PodioOutput("out", filename=output_rec) -- GitLab