From d86381e432459c8b95b3eb976ca2ccf8d93eb739 Mon Sep 17 00:00:00 2001 From: Whitney Armstrong <warmstrong@anl.gov> Date: Wed, 19 May 2021 22:34:09 -0500 Subject: [PATCH] modified: gem_tracker.xml new file: scripts/tutorial1_hit_position.cxx --- gem_tracker.xml | 9 ++- scripts/tutorial1_hit_position.cxx | 107 +++++++++++++++++++++++++++++ 2 files changed, 114 insertions(+), 2 deletions(-) create mode 100644 scripts/tutorial1_hit_position.cxx diff --git a/gem_tracker.xml b/gem_tracker.xml index 3339122..8351846 100644 --- a/gem_tracker.xml +++ b/gem_tracker.xml @@ -93,9 +93,14 @@ <readouts> <readout name="GEMTrackerHits"> - <segmentation type="PolarGridRPhi" grid_size_phi="1.0*degree" grid_size_r="1.0*cm"/> + <segmentation type="CartesianGridXY" grid_size_x="1*cm" grid_size_y="3*cm" /> + <id>system:8,barrel:2,layer:4,module:12,sensor:2,x:32:-16,y:-16</id> + <!-- + <segmentation type="PolarGridRPhi" grid_size_phi="5.0*degree" grid_size_r="1.0*cm"/> <id>system:5,barrel:3,layer:4,module:5,r:32:-16,phi:-16</id> - </readout> + --> + </readout> + </readouts> <!-- diff --git a/scripts/tutorial1_hit_position.cxx b/scripts/tutorial1_hit_position.cxx new file mode 100644 index 0000000..50ac7c9 --- /dev/null +++ b/scripts/tutorial1_hit_position.cxx @@ -0,0 +1,107 @@ +R__LOAD_LIBRARY(libfmt.so) +#include "fmt/core.h" + +#include "Math/Vector3D.h" +#include "Math/Vector4D.h" +#include "Math/VectorUtil.h" + +#include "ROOT/RDataFrame.hxx" +#include "TCanvas.h" +#include "TChain.h" +#include "TF1.h" +#include "TFile.h" +#include "TH1D.h" +#include "TH1F.h" +#include "TLegend.h" +#include "TMath.h" +#include "TRandom3.h" +#include "TTree.h" + +#include <algorithm> +#include <iterator> +#include <tuple> +#include <vector> +#include <cmath> + +#include "DD4hep/Detector.h" +#include "DDRec/CellIDPositionConverter.h" + +#include "dd4pod/Geant4ParticleCollection.h" +#include "dd4pod/TrackerHitCollection.h" + +/** Hit position example. + * + */ +void tutorial1_hit_position(const char* fname = "gem_tracker_sim.root") { + using namespace ROOT::Math; + + TChain* t = new TChain("events"); + t->Add(fname); + + ROOT::RDataFrame d0(*t, {"GEMTrackerHits", "MCParticles"}); + //std::cout << d0.DescribeDataset() << std::endl; + + // ------------------------- + // Get the DD4hep instance + // Load the compact XML file + // Initialize the position converter tool + dd4hep::Detector& detector = dd4hep::Detector::getInstance(); + detector.fromCompact("gem_tracker.xml"); + dd4hep::rec::CellIDPositionConverter cellid_converter(detector); + + // Simple lambda to define nhits branch + auto nhits = [](const std::vector<dd4pod::TrackerHitData>& evt) { return (int)evt.size(); }; + + auto local_position = [&](const std::vector<dd4pod::TrackerHitData>& hits) { + std::vector<std::array<double, 2>> result; + for (const auto& h : hits) { + // The actual hit position: + auto pos0 = (h.position); + + // **This is the most important part.** It provides a way of getting the + // geometry information in a flexible way. + // The segmentation hit postion: + auto pos1 = cellid_converter.position(h.cellID); + fmt::print(" Hit Position : {},{},{}\n", pos0.x / 10.0, pos0.y / 10.0, pos0.z / 10.0); + fmt::print("Segmentation-Cell Position : {},{},{}\n", pos1.x(), pos1.y(), pos1.z()); + result.push_back({pos0.x / 10.0 - pos1.x(), pos0.y / 10.0 - pos1.y()}); + } + return result; + }; + + auto x_pos = [&](const std::vector<std::array<double, 2>>& xypos) { + std::vector<double> result; + for (const auto& h : xypos) { + result.push_back(h.at(0)); + } + return result; + }; + auto y_pos = [&](const std::vector<std::array<double, 2>>& xypos) { + std::vector<double> result; + for (auto h : xypos) { + result.push_back(h.at(1)); + } + return result; + }; + + auto d1 = d0.Define("nhits", nhits, {"GEMTrackerHits"}) + .Filter([=](const std::vector<dd4pod::TrackerHitData>& hits) { + for (auto h : hits) { + auto pos = h.position; + if (std::hypot(pos.x,pos.y) > 30.0) { + return true; + } + } + return false; + }, + {"GEMTrackerHits"}) + .Define("xy_hit_pos", local_position, {"GEMTrackerHits"}) + .Define("x_pos", x_pos, {"xy_hit_pos"}) + .Define("y_pos", y_pos, {"xy_hit_pos"}); + + auto h_local_pos = d1.Histo2D(TH2D("h_local_pos", ";x [cm]; y [cm] ", 100, -5.0, 5.0, 100, -5.0, 5.0), "x_pos", "y_pos"); + + TCanvas* c = new TCanvas(); + h_local_pos->DrawClone("colz"); + c->SaveAs("results/hit_position.png"); +} -- GitLab