Skip to content
Snippets Groups Projects
Commit 0773c220 authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Fix dvmp

parent ade23046
No related branches found
No related tags found
1 merge request!64Fix dvmp
Hide whitespace changes
Inline Side-by-side
benchmarks/dvmp/analysis/dvmp.h
+ 7
− 7
View file @ 0773c220
#ifndef DVMP_H
#define DVMP_H
#include <util.h>
#include "common_bench/util.h"
#include <algorithm>
#include <cmath>
#include <exception>
#include <fmt/core.h>
#include <limits>
#include <string>
#include <vector>
#include <Math/Vector4D.h>
#include "fmt/core.h"
#include "Math/Vector4D.h"
// Additional utility functions for DVMP benchmarks. Where useful, these can be
// promoted to the top-level util library
......@@ -26,7 +26,7 @@ namespace util {
//0:e0 1:p0 2:e1 3:p1 4:recoil system (without p1) 5:l1 from 4 6:l2 from 4
inline auto momenta_sort_sim(const std::vector<dd4pod::Geant4ParticleData>& parts, std::string_view mother, std::string_view daughter){//mother and daughter are not used yet; will be useful when generater is different and/or when the mcparticles doesn't follow the same order in all events
std::vector<ROOT::Math::PxPyPzMVector> momenta{7};
int order_map[7] = {0, 3, 2, 6, 5, 7, 8};
int order_map[7] = {0, 3, 2, 6, 5, 7, 8};
for(int i = 0 ; i < 7 ; i++){
double px = parts[order_map[i]].psx;
double py = parts[order_map[i]].psy;
......@@ -49,10 +49,10 @@ namespace util {
for(int i = 0 ; i < 7 ; i++) momenta[i].SetPxPyPzE(0., 0., 0., 0.); //initialize as all 0
//manually set incoming electron and proton;
double e0_mass = get_pdg_mass("electron");
double e0_mass = common_bench::get_pdg_mass("electron");
double e0_pz = 1.305e-8 - 10.;
momenta[0].SetPxPyPzE(0., 0., e0_pz, sqrt(e0_mass*e0_mass + e0_pz*e0_pz));
double p0_mass = get_pdg_mass("proton");
double p0_mass = common_bench::get_pdg_mass("proton");
double p0_pz = 99.995598 + 1.313e-7 + 8.783e-11;
momenta[1].SetPxPyPzE(0., 0., p0_pz, sqrt(p0_mass*p0_mass + p0_pz*p0_pz));
......@@ -88,7 +88,7 @@ namespace util {
double energy_tmp2 = parts[j].energy;
lpt_2.SetPxPyPzE(parts[j].p.x, parts[j].p.y, parts[j].p.z, energy_tmp2);
const double new_mass{(lpt_1 + lpt_2).mass()};
if (fabs(new_mass - get_pdg_mass(mother)) < fabs(best_mass - get_pdg_mass(mother))) {
if (fabs(new_mass - common_bench::get_pdg_mass(mother)) < fabs(best_mass - common_bench::get_pdg_mass(mother))) {
first = i;
second = j;
best_mass = new_mass;
......
benchmarks/dvmp/analysis/vm_invar.cxx
+ 32
− 30
View file @ 0773c220
#include "dvmp.h"
#include "plot.h"
#include <common_bench/benchmark.h>
#include <common_bench/mt.h>
#include <common_bench/util.h>
#include "common_bench/plot.h"
#include "common_bench/benchmark.h"
#include "common_bench/mt.h"
#include "common_bench/util.h"
#include "ROOT/RDataFrame.hxx"
#include <ROOT/RDataFrame.hxx>
#include <cmath>
#include <fmt/color.h>
#include <fmt/core.h>
#include <fstream>
#include <iostream>
#include <nlohmann/json.hpp>
#include <string>
#include <vector>
#include "nlohmann/json.hpp"
#include "fmt/color.h"
#include "fmt/core.h"
// Run VM invariant-mass-based benchmarks on an input reconstruction file for
// a desired vector meson (e.g. jpsi) and a desired decay particle (e.g. muon)
// Output figures are written to our output prefix (which includes the output
......@@ -57,8 +59,8 @@ int vm_invar(const std::string& config_name)
ROOT::EnableImplicitMT(kNumThreads);
// The particles we are looking for. E.g. J/psi decaying into e+e-
const double vm_mass = util::get_pdg_mass(vm_name);
const double decay_mass = util::get_pdg_mass(decay_name);
const double vm_mass = common_bench::get_pdg_mass(vm_name);
const double decay_mass = common_bench::get_pdg_mass(decay_name);
// Ensure our output prefix always ends on a dot, a slash or a dash
if (output_prefix.back() != '.' && output_prefix.back() != '/' && output_prefix.back() != '-') {
......@@ -125,9 +127,9 @@ int vm_invar(const std::string& config_name)
auto& hy_rec = *h_y_rec;
auto& hy_sim = *h_y_sim;
// histogram style
hy_rec.SetLineColor(plot::kMpOrange);
hy_rec.SetLineColor(common_bench::plot::kMpOrange);
hy_rec.SetLineWidth(1);
hy_sim.SetLineColor(plot::kMpBlue);
hy_sim.SetLineColor(common_bench::plot::kMpBlue);
hy_sim.SetLineWidth(2);
// axes
hy_sim.GetXaxis()->CenterTitle();
......@@ -135,16 +137,16 @@ int vm_invar(const std::string& config_name)
hy_sim.DrawClone("hist");
hy_rec.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr1;
auto t1 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t1->SetFillColorAlpha(kWhite, 0);
t1->SetTextFont(43);
t1->SetTextSize(25);
tptr1 = t1->AddText("simulated");
tptr1->SetTextColor(plot::kMpBlue);
tptr1->SetTextColor(common_bench::plot::kMpBlue);
tptr1 = t1->AddText("rec(PlaceHolder)");
tptr1->SetTextColor(plot::kMpOrange);
tptr1->SetTextColor(common_bench::plot::kMpOrange);
t1->Draw();
// pad 2 Q2
......@@ -152,9 +154,9 @@ int vm_invar(const std::string& config_name)
auto& hQ2_rec = *h_Q2_rec;
auto& hQ2_sim = *h_Q2_sim;
// histogram style
hQ2_rec.SetLineColor(plot::kMpOrange);
hQ2_rec.SetLineColor(common_bench::plot::kMpOrange);
hQ2_rec.SetLineWidth(1);
hQ2_sim.SetLineColor(plot::kMpBlue);
hQ2_sim.SetLineColor(common_bench::plot::kMpBlue);
hQ2_sim.SetLineWidth(2);
// axes
hQ2_sim.GetXaxis()->CenterTitle();
......@@ -162,16 +164,16 @@ int vm_invar(const std::string& config_name)
hQ2_sim.DrawClone("hist");
hQ2_rec.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr2;
auto t2 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t2->SetFillColorAlpha(kWhite, 0);
t2->SetTextFont(43);
t2->SetTextSize(25);
tptr2 = t2->AddText("simulated");
tptr2->SetTextColor(plot::kMpBlue);
tptr2->SetTextColor(common_bench::plot::kMpBlue);
tptr2 = t2->AddText("rec(PlaceHolder)");
tptr2->SetTextColor(plot::kMpOrange);
tptr2->SetTextColor(common_bench::plot::kMpOrange);
t2->Draw();
// pad 3 x
......@@ -179,9 +181,9 @@ int vm_invar(const std::string& config_name)
auto& hx_rec = *h_x_rec;
auto& hx_sim = *h_x_sim;
// histogram style
hx_rec.SetLineColor(plot::kMpOrange);
hx_rec.SetLineColor(common_bench::plot::kMpOrange);
hx_rec.SetLineWidth(1);
hx_sim.SetLineColor(plot::kMpBlue);
hx_sim.SetLineColor(common_bench::plot::kMpBlue);
hx_sim.SetLineWidth(2);
// axes
hx_sim.GetXaxis()->CenterTitle();
......@@ -189,16 +191,16 @@ int vm_invar(const std::string& config_name)
hx_sim.DrawClone("hist");
hx_rec.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr3;
auto t3 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t3->SetFillColorAlpha(kWhite, 0);
t3->SetTextFont(43);
t3->SetTextSize(25);
tptr3 = t3->AddText("simulated");
tptr3->SetTextColor(plot::kMpBlue);
tptr3->SetTextColor(common_bench::plot::kMpBlue);
tptr3 = t3->AddText("rec(PlaceHolder)");
tptr3->SetTextColor(plot::kMpOrange);
tptr3->SetTextColor(common_bench::plot::kMpOrange);
t3->Draw();
// pad 4 t
......@@ -206,9 +208,9 @@ int vm_invar(const std::string& config_name)
auto& ht_rec = *h_t_rec;
auto& ht_sim = *h_t_sim;
// histogram style
ht_rec.SetLineColor(plot::kMpOrange);
ht_rec.SetLineColor(common_bench::plot::kMpOrange);
ht_rec.SetLineWidth(1);
ht_sim.SetLineColor(plot::kMpBlue);
ht_sim.SetLineColor(common_bench::plot::kMpBlue);
ht_sim.SetLineWidth(2);
// axes
ht_sim.GetXaxis()->CenterTitle();
......@@ -216,16 +218,16 @@ int vm_invar(const std::string& config_name)
ht_sim.DrawClone("hist");
ht_rec.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr4;
auto t4 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t4->SetFillColorAlpha(kWhite, 0);
t4->SetTextFont(43);
t4->SetTextSize(25);
tptr4 = t4->AddText("simulated");
tptr4->SetTextColor(plot::kMpBlue);
tptr4->SetTextColor(common_bench::plot::kMpBlue);
tptr4 = t4->AddText("rec(PlaceHolder)");
tptr4->SetTextColor(plot::kMpOrange);
tptr4->SetTextColor(common_bench::plot::kMpOrange);
t4->Draw();
c.Print(fmt::format("{}InvariantQuantities.png", output_prefix).c_str());
......
benchmarks/dvmp/analysis/vm_mass.cxx
+ 36
− 34
View file @ 0773c220
#include "dvmp.h"
#include "plot.h"
#include "common_bench/plot.h"
#include <common_bench/benchmark.h>
#include <common_bench/mt.h>
#include <common_bench/util.h>
#include "common_bench/benchmark.h"
#include "common_bench/mt.h"
#include "common_bench/util.h"
#include <ROOT/RDataFrame.hxx>
#include "ROOT/RDataFrame.hxx"
#include <cmath>
#include <fmt/color.h>
#include <fmt/core.h>
#include <fstream>
#include <iostream>
#include <nlohmann/json.hpp>
#include <string>
#include <vector>
#include "fmt/color.h"
#include "fmt/core.h"
#include "nlohmann/json.hpp"
#include "eicd/ReconstructedParticleCollection.h"
#include "eicd/ReconstructedParticleData.h"
......@@ -74,8 +76,8 @@ int vm_mass(const std::string& config_name)
ROOT::EnableImplicitMT(kNumThreads);
// The particles we are looking for. E.g. J/psi decaying into e+e-
const double vm_mass = util::get_pdg_mass(vm_name);
const double decay_mass = util::get_pdg_mass(decay_name);
const double vm_mass = common_bench::get_pdg_mass(vm_name);
const double decay_mass = common_bench::get_pdg_mass(decay_name);
// Ensure our output prefix always ends on a dot, a slash or a dash
if (output_prefix.back() != '.' && output_prefix.back() != '/' && output_prefix.back() != '-') {
......@@ -89,14 +91,14 @@ int vm_mass(const std::string& config_name)
// types
auto find_decay_pair = [vm_mass, decay_mass](const std::vector<ROOT::Math::PxPyPzMVector>& parts) {
return util::find_decay_pair(parts, vm_mass, decay_mass);
return common_bench::find_decay_pair(parts, vm_mass, decay_mass);
};
// util::PrintGeant4(mcparticles2);
// common_bench::PrintGeant4(mcparticles2);
// Define analysis flow
auto d_im = d.Define("p_rec", util::momenta_RC, {"DummyReconstructedParticles"}) //using dummy rc
auto d_im = d.Define("p_rec", common_bench::momenta_RC, {"DummyReconstructedParticles"}) //using dummy rc
.Define("N", "p_rec.size()")
.Define("p_sim", util::momenta_from_simulation, {"mcparticles2"})
.Define("p_sim", common_bench::momenta_from_simulation, {"mcparticles2"})
.Define("decay_pair_rec", find_decay_pair, {"p_rec"})
.Define("decay_pair_sim", find_decay_pair, {"p_sim"})
.Define("p_vm_rec", "decay_pair_rec.first + decay_pair_rec.second")
......@@ -137,9 +139,9 @@ int vm_mass(const std::string& config_name)
auto& h11 = *h_im_sim;
auto& h12 = *h_im_rec;
// histogram style
h11.SetLineColor(plot::kMpBlue);
h11.SetLineColor(common_bench::plot::kMpBlue);
h11.SetLineWidth(2);
h12.SetLineColor(plot::kMpOrange);
h12.SetLineColor(common_bench::plot::kMpOrange);
h12.SetLineWidth(1);
// axes
h11.GetXaxis()->CenterTitle();
......@@ -162,16 +164,16 @@ int vm_mass(const std::string& config_name)
mfMass->Draw("same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr1;
auto t1 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t1->SetFillColorAlpha(kWhite, 0);
t1->SetTextFont(43);
t1->SetTextSize(25);
tptr1 = t1->AddText("simulated");
tptr1->SetTextColor(plot::kMpBlue);
tptr1->SetTextColor(common_bench::plot::kMpBlue);
tptr1 = t1->AddText("reconstructed");
tptr1->SetTextColor(plot::kMpOrange);
tptr1->SetTextColor(common_bench::plot::kMpOrange);
t1->Draw();
// pad 2 pt
......@@ -181,9 +183,9 @@ int vm_mass(const std::string& config_name)
auto& h21 = *h_pt_sim;
auto& h22 = *h_pt_rec;
// histogram style
h21.SetLineColor(plot::kMpBlue);
h21.SetLineColor(common_bench::plot::kMpBlue);
h21.SetLineWidth(2);
h22.SetLineColor(plot::kMpOrange);
h22.SetLineColor(common_bench::plot::kMpOrange);
h22.SetLineWidth(1);
// axes
h21.GetXaxis()->CenterTitle();
......@@ -193,16 +195,16 @@ int vm_mass(const std::string& config_name)
h22.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr2;
auto t2 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t2->SetFillColorAlpha(kWhite, 0);
t2->SetTextFont(43);
t2->SetTextSize(25);
tptr2 = t2->AddText("simulated");
tptr2->SetTextColor(plot::kMpBlue);
tptr2->SetTextColor(common_bench::plot::kMpBlue);
tptr2 = t2->AddText("reconstructed");
tptr2->SetTextColor(plot::kMpOrange);
tptr2->SetTextColor(common_bench::plot::kMpOrange);
t2->Draw();
// pad 3 phi
......@@ -212,9 +214,9 @@ int vm_mass(const std::string& config_name)
auto& h31 = *h_phi_sim;
auto& h32 = *h_phi_rec;
// histogram style
h31.SetLineColor(plot::kMpBlue);
h31.SetLineColor(common_bench::plot::kMpBlue);
h31.SetLineWidth(2);
h32.SetLineColor(plot::kMpOrange);
h32.SetLineColor(common_bench::plot::kMpOrange);
h32.SetLineWidth(1);
// axes
h31.GetXaxis()->CenterTitle();
......@@ -223,16 +225,16 @@ int vm_mass(const std::string& config_name)
h31.DrawClone("hist");
h32.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr3;
auto t3 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t3->SetFillColorAlpha(kWhite, 0);
t3->SetTextFont(43);
t3->SetTextSize(25);
tptr3 = t3->AddText("simulated");
tptr3->SetTextColor(plot::kMpBlue);
tptr3->SetTextColor(common_bench::plot::kMpBlue);
tptr3 = t3->AddText("reconstructed");
tptr3->SetTextColor(plot::kMpOrange);
tptr3->SetTextColor(common_bench::plot::kMpOrange);
t3->Draw();
// pad 4 rapidity
......@@ -242,9 +244,9 @@ int vm_mass(const std::string& config_name)
auto& h41 = *h_eta_sim;
auto& h42 = *h_eta_rec;
// histogram style
h41.SetLineColor(plot::kMpBlue);
h41.SetLineColor(common_bench::plot::kMpBlue);
h41.SetLineWidth(2);
h42.SetLineColor(plot::kMpOrange);
h42.SetLineColor(common_bench::plot::kMpOrange);
h42.SetLineWidth(1);
// axes
h41.GetXaxis()->CenterTitle();
......@@ -253,16 +255,16 @@ int vm_mass(const std::string& config_name)
h41.DrawClone("hist");
h42.DrawClone("hist same");
// FIXME hardcoded beam configuration
plot::draw_label(10, 100, detector);
common_bench::plot::draw_label(10, 100, detector);
TText* tptr4;
auto t4 = new TPaveText(.6, .8417, .9, .925, "NB NDC");
t4->SetFillColorAlpha(kWhite, 0);
t4->SetTextFont(43);
t4->SetTextSize(25);
tptr4 = t4->AddText("simulated");
tptr4->SetTextColor(plot::kMpBlue);
tptr4->SetTextColor(common_bench::plot::kMpBlue);
tptr4 = t4->AddText("reconstructed");
tptr4->SetTextColor(plot::kMpOrange);
tptr4->SetTextColor(common_bench::plot::kMpOrange);
t4->Draw();
c.Print(fmt::format("{}vm_mass_pt_phi_rapidity.png", output_prefix).c_str());
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment