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

Added new tracking analysis script.

	new file:   scripts/rec_central_electrons.cxx
parent 961c27a3
Branches
Tags
1 merge request!27Added new tracking analysis script.
#include "ROOT/RDataFrame.hxx"
#include <iostream>
#include "dd4pod/Geant4ParticleCollection.h"
#include "eicd/TrackParametersCollection.h"
#include "eicd/ClusterCollection.h"
#include "eicd/ClusterData.h"
//namespace edm4hep {
//
//std::vector<float> pt (std::vector<MCParticleData> const& in){
// std::vector<float> result;
// for (size_t i = 0; i < in.size(); ++i) {
// result.push_back(std::sqrt(in[i].momentum.x * in[i].momentum.x + in[i].momentum.y * in[i].momentum.y));
// }
// return result;
//}
//
//std::vector<float> eta(std::vector<MCParticleData> const& in){
// std::vector<float> result;
// ROOT::Math::PxPyPzMVector lv;
// for (size_t i = 0; i < in.size(); ++i) {
// lv.SetCoordinates(in[i].momentum.x, in[i].momentum.y, in[i].momentum.z, in[i].mass);
// result.push_back(lv.Eta());
// }
// return result;
//}
//
//std::vector<float> cos_theta(std::vector<MCParticleData> const& in){
// std::vector<float> result;
// ROOT::Math::PxPyPzMVector lv;
// for (size_t i = 0; i < in.size(); ++i) {
// lv.SetCoordinates(in[i].momentum.x, in[i].momentum.y, in[i].momentum.z, in[i].mass);
// result.push_back(cos(lv.Theta()));
// }
// return result;
//}
//
//}
auto p_track = [](std::vector<eic::TrackParametersData> const& in) {
std::vector<double> result;
for (size_t i = 0; i < in.size(); ++i) {
result.push_back(std::abs(1.0/(in[i].qOverP)));
}
return result;
};
std::vector<float> pt (std::vector<dd4pod::Geant4ParticleData> const& in){
std::vector<float> result;
for (size_t i = 0; i < in.size(); ++i) {
result.push_back(std::sqrt(in[i].psx * in[i].psx + in[i].psy * in[i].psy));
}
return result;
}
using ROOT::RDataFrame;
using namespace ROOT::VecOps;
auto momentum = [](std::vector<ROOT::Math::PxPyPzMVector> const& in) {
std::vector<double> result;
for (size_t i = 0; i < in.size(); ++i) {
result.push_back(in[i].E());
}
return result;
};
auto theta = [](std::vector<ROOT::Math::PxPyPzMVector> const& in) {
std::vector<double> result;
for (size_t i = 0; i < in.size(); ++i) {
result.push_back(in[i].Theta()*180/M_PI);
}
return result;
};
auto fourvec = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) {
std::vector<ROOT::Math::PxPyPzMVector> result;
ROOT::Math::PxPyPzMVector lv;
for (size_t i = 0; i < in.size(); ++i) {
lv.SetCoordinates(in[i].psx, in[i].psy, in[i].psz, in[i].mass);
result.push_back(lv);
}
return result;
};
int rec_central_electrons(const char* fname = "topside/rec_central_electrons.root") {
ROOT::EnableImplicitMT();
ROOT::RDataFrame df("events", fname);
auto df0 = df.Define("isThrown", "mcparticles2.genStatus == 1")
.Define("thrownParticles", "mcparticles2[isThrown]")
.Define("thrownP", fourvec, {"thrownParticles"})
.Define("p_thrown", momentum, {"thrownP"})
.Define("nTracks", "outputTrackParameters.size()")
.Define("p_track", p_track, {"outputTrackParameters"})
.Define("delta_p",
[](const std::vector<double> &tracks,
const std::vector<double> &thrown) {
std::vector<double> res;
for (const auto &p1 : thrown) {
for (const auto &p2 : tracks) {
res.push_back(p1 - p2);
}
}
return res;
},
{"p_track", "p_thrown"});
auto h_nTracks = df0.Histo1D({"h_nTracks", "; N tracks ", 10, 0, 10}, "nTracks");
auto h_pTracks = df0.Histo1D({"h_pTracks", "; GeV/c ", 100, 0, 10}, "p_track");
auto h_delta_p = df0.Histo1D({"h_delta_p", "; GeV/c ", 100, -10, 10}, "delta_p");
auto c = new TCanvas();
h_nTracks->DrawCopy();
c->SaveAs("results/rec_central_electrons_nTracks.png");
h_pTracks->DrawCopy();
c->SaveAs("results/rec_central_electrons_pTracks.png");
h_delta_p->DrawCopy();
c->SaveAs("results/rec_central_electrons_delta_p.png");
return 0;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment