Skip to content
Snippets Groups Projects

EDM4hep: mcparticles -> MCParticles

Merged Wouter Deconinck requested to merge edm4hep-mcparticles into master
1 file
+ 4
5
Compare changes
  • Side-by-side
  • Inline
@@ -10,9 +10,9 @@
R__LOAD_LIBRARY(libeicd.so)
R__LOAD_LIBRARY(libDD4pod.so)
#include "dd4pod/TrackerHitCollection.h"
#include "dd4pod/TrackerHitData.h"
#include "dd4pod/Geant4ParticleCollection.h"
#include "edm4hep/SimTrackerHitCollection.h"
#include "edm4hep/SimTrackerHitData.h"
#include "edm4hep/MCParticleCollection.h"
#include "eicd/TrackParametersCollection.h"
#include "eicd/ClusterCollection.h"
#include "eicd/ClusterData.h"
@@ -30,10 +30,10 @@ auto p_track = [](std::vector<eic::TrackParametersData> const& in) {
};
std::vector<float> pt (std::vector<dd4pod::Geant4ParticleData> const& in){
std::vector<float> pt (std::vector<edm4hep::MCParticleData> const& in){
std::vector<float> result;
for (size_t i = 0; i < in.size(); ++i) {
result.push_back(std::sqrt(in[i].ps.x * in[i].ps.x + in[i].ps.y * in[i].ps.y));
result.push_back(std::sqrt(in[i].momentum.x * in[i].momentum.x + in[i].momentum.y * in[i].momentum.y));
}
return result;
}
@@ -52,11 +52,11 @@ auto theta = [](std::vector<ROOT::Math::PxPyPzMVector> const& in) {
}
return result;
};
auto fourvec = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) {
auto fourvec = [](ROOT::VecOps::RVec<edm4hep::MCParticleData> 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].ps.x, in[i].ps.y, in[i].ps.z, in[i].mass);
lv.SetCoordinates(in[i].momentum.x, in[i].momentum.y, in[i].momentum.z, in[i].mass);
result.push_back(lv);
}
return result;
@@ -82,14 +82,14 @@ auto delta_p_over_p = [](const std::vector<double>& tracks, const std::vector<do
return res;
};
int hits_far_forward_protons(const char* fname = "sim_far_forward_protons.root")
int hits_far_forward_protons(const char* fname = "sim_far_forward_protons.edm4hep.root")
{
ROOT::EnableImplicitMT();
ROOT::RDataFrame df("events", fname);
auto df0 = df.Define("isThrown", "mcparticles.genStatus == 1")
.Define("thrownParticles", "mcparticles[isThrown]")
auto df0 = df.Define("isThrown", "MCParticles.generatorStatus == 1")
.Define("thrownParticles", "MCParticles[isThrown]")
.Define("thrownP", fourvec, {"thrownParticles"})
.Define("p_thrown", momentum, {"thrownP"})
.Define("theta_thrown", theta, {"thrownP"})
@@ -105,8 +105,8 @@ int hits_far_forward_protons(const char* fname = "sim_far_forward_protons.root")
//.Define("delta_p_over_p1",delta_p_over_p, {"p_track1", "p_thrown"})
//.Define("delta_p_over_p2",delta_p_over_p, {"p_track2", "p_thrown"})
//.Define("N_VtxBarrelHits",[](std::vector<eic::TrackerHitData> hits) { return hits.size();},{"VertexBarrelRecHits"})
.Define("N_BarrelHits", [](std::vector<dd4pod::TrackerHitData> hits) { return hits.size();}, {"TrackerBarrelHits"})
.Define("N_EndcapHits", [](std::vector<dd4pod::TrackerHitData> hits) { return hits.size();}, {"TrackerEndcapHits"})
.Define("N_BarrelHits", [](std::vector<edm4hep::SimTrackerHitData> hits) { return hits.size();}, {"TrackerBarrelHits"})
.Define("N_EndcapHits", [](std::vector<edm4hep::SimTrackerHitData> hits) { return hits.size();}, {"TrackerEndcapHits1"})
;
auto hBarrel_x_vs_y = df0.Histo2D({"hBarrel_x_vs_y", "; x ; y ", 100, -900, 900,100, -900, 900 }, "TrackerBarrelHits.position.x", "TrackerBarrelHits.position.y");
Loading