diff --git a/benchmarks/dvmp/analysis/dvmp.h b/benchmarks/dvmp/analysis/dvmp.h index 8c1e8a09a930af8b9e625f3ca6c154ac8e61438b..9b554e1419f81dd2e11a46bb743c688fefb848e0 100644 --- a/benchmarks/dvmp/analysis/dvmp.h +++ b/benchmarks/dvmp/analysis/dvmp.h @@ -120,16 +120,6 @@ namespace util { return momenta; } - - - inline double test_jpsi_pt(const std::vector<ROOT::Math::PxPyPzMVector>& momenta){ - double px = momenta[4].px(); - double py = momenta[4].py(); - double pt = sqrt(px*px + py*py); - return pt; - } - - //for Dummy rc inline inv_quant calc_inv_quant_rec(const std::vector<ROOT::Math::PxPyPzMVector>& parts, const double pdg_mass, const double daughter_mass){ int first = -1; @@ -196,6 +186,22 @@ namespace util { return quantities; } + inline inv_quant calc_inv_quant_both(const std::vector<ROOT::Math::PxPyPzMVector>& parts) + { + //0:e0 1:p0 2:e1 3:p1 4:recoil system (without p1) 5:l1 from 4 6:l2 from 4 + ROOT::Math::PxPyPzMVector q(parts[0] - parts[2]); + ROOT::Math::PxPyPzMVector P(parts[1]); + ROOT::Math::PxPyPzMVector Delta(parts[3] - parts[1]);//jpsi->l l' gamma, ignore gamma + + //ROOT::Math::PxPyPzMVector Delta(parts[6] - parts[3]);//exact jpsi + double nu = q.Dot(P) / P.mass(); + double Q2 = -q.Dot(q); + double t = Delta.Dot(Delta); + inv_quant quantities = {nu, Q2, Q2/2./P.mass()/nu, t}; + return quantities; + } + + inline double get_nu(inv_quant quantities) { return quantities.nu / 1000.; } inline double get_Q2(inv_quant quantities) { return quantities.Q2; } inline double get_x(inv_quant quantities) { return quantities.x; } diff --git a/benchmarks/dvmp/analysis/vm_invar.cxx b/benchmarks/dvmp/analysis/vm_invar.cxx index 647aa235537f587ac2f3860e00212a18f64c49cf..353c816dbf620bd8f6a894fcd171d65dd2ea255a 100644 --- a/benchmarks/dvmp/analysis/vm_invar.cxx +++ b/benchmarks/dvmp/analysis/vm_invar.cxx @@ -84,7 +84,8 @@ int vm_invar(const std::string& config_name) // Define analysis flow auto d_im = d.Define("p_rec", util::momenta_RC, {"DummyReconstructedParticles"}) - .Define("test_sim_ordered", momenta_ordered_sim, {"mcparticles2"}) + .Define("p_sim_ordered", momenta_ordered_sim, {"mcparticles2"}) + .Define("p_rec_ordered", momenta_ordered_rec, {"DummyReconstructedParticles"}) .Define("N", "p_rec.size()") .Define("p_sim", util::momenta_from_simulation, {"mcparticles2"}) //================================================================ @@ -98,7 +99,8 @@ int vm_invar(const std::string& config_name) .Define("Q2_sim", util::get_Q2, {"invariant_quantities_sim"}) .Define("x_sim", util::get_x, {"invariant_quantities_sim"}) .Define("t_sim", util::get_t, {"invariant_quantities_sim"}) - .Define("test_pt", util::test_jpsi_pt, {"test_sim_ordered"}); + .Define("invariant_quantities_rec", calc_inv_quant_both, {"p_rec_ordered"}) + .Define("invariant_quantities_sim", calc_inv_quant_both, {"p_sim_ordered"}); //================================================================ // Define output histograms @@ -113,8 +115,6 @@ int vm_invar(const std::string& config_name) auto h_x_rec = d_im.Histo1D({"h_x_rec", ";x;#", 100, 0., 0.1}, "x_rec"); auto h_t_rec = d_im.Histo1D({"h_t_rec", ";t;#", 100, -1., 0.}, "t_rec"); - auto h_test_pt = d_im.Histo1D({"h_test_pt", "h_test_pt", 50, 0., 10.}, "test_pt"); - // Plot our histograms. // TODO: to start I'm explicitly plotting the histograms, but want to // factorize out the plotting code moving forward. @@ -205,9 +205,9 @@ int vm_invar(const std::string& config_name) tptr3->SetTextColor(plot::kMpOrange); t3->Draw(); - // pad 3 x + // pad 4 t c.cd(4); - /*auto& ht_rec = *h_t_rec; + auto& ht_rec = *h_t_rec; auto& ht_sim = *h_t_sim; // histogram style ht_rec.SetLineColor(plot::kMpOrange); @@ -230,25 +230,7 @@ int vm_invar(const std::string& config_name) tptr4->SetTextColor(plot::kMpBlue); tptr4 = t4->AddText("rec(PlaceHolder)"); tptr4->SetTextColor(plot::kMpOrange); - t4->Draw();*/ - auto& htest = *h_test_pt; - htest.SetLineColor(plot::kMpBlue); - htest.SetLineWidth(2); - htest.GetXaxis()->CenterTitle(); - htest.GetYaxis()->CenterTitle(); - htest.DrawClone("hist"); - plot::draw_label(10, 100, detector); - TText* tptrtest; - auto ttest = new TPaveText(.6, .8417, .9, .925, "NB NDC"); - ttest->SetFillColorAlpha(kWhite, 0); - ttest->SetTextFont(43); - ttest->SetTextSize(25); - tptrtest = ttest->AddText("simulated"); - tptrtest->SetTextColor(plot::kMpBlue); - ttest->Draw(); - - - + t4->Draw(); c.Print(fmt::format("{}InvariantQuantities.png", output_prefix).c_str());