ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
RootSimHitWriter.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file RootSimHitWriter.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2017-2018 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #include "ACTFW/Io/Root/RootSimHitWriter.hpp"
10 
11 #include <TFile.h>
12 #include <TTree.h>
13 #include <ios>
14 #include <stdexcept>
15 
16 #include "Acts/Utilities/Units.hpp"
17 
18 FW::RootSimHitWriter::RootSimHitWriter(const FW::RootSimHitWriter::Config& cfg,
19  Acts::Logging::Level lvl)
20  : WriterT(cfg.inputSimulatedHits, "RootSimHitWriter", lvl), m_cfg(cfg) {
21  // inputParticles is already checked by base constructor
22  if (m_cfg.filePath.empty()) {
23  throw std::invalid_argument("Missing file path");
24  }
25  if (m_cfg.treeName.empty()) {
26  throw std::invalid_argument("Missing tree name");
27  }
28 
29  // open root file and create the tree
30  m_outputFile = TFile::Open(m_cfg.filePath.c_str(), m_cfg.fileMode.c_str());
31  if (m_outputFile == nullptr) {
32  throw std::ios_base::failure("Could not open '" + m_cfg.filePath + "'");
33  }
34  m_outputFile->cd();
35  m_outputTree = new TTree(m_cfg.treeName.c_str(), m_cfg.treeName.c_str());
36  if (m_outputTree == nullptr) {
37  throw std::bad_alloc();
38  }
39 
40  // setup the branches
41  m_outputTree->Branch("event_id", &m_eventId);
42  m_outputTree->Branch("geometry_id", &m_geometryId, "geometry_id/l");
43  m_outputTree->Branch("particle_id", &m_particleId, "particle_id/l");
44  m_outputTree->Branch("tx", &m_tx);
45  m_outputTree->Branch("ty", &m_ty);
46  m_outputTree->Branch("tz", &m_tz);
47  m_outputTree->Branch("tt", &m_tt);
48  m_outputTree->Branch("tpx", &m_tpx);
49  m_outputTree->Branch("tpy", &m_tpy);
50  m_outputTree->Branch("tpz", &m_tpz);
51  m_outputTree->Branch("te", &m_te);
52  m_outputTree->Branch("deltapx", &m_deltapx);
53  m_outputTree->Branch("deltapy", &m_deltapy);
54  m_outputTree->Branch("deltapz", &m_deltapz);
55  m_outputTree->Branch("deltae", &m_deltae);
56  m_outputTree->Branch("index", &m_index);
57  m_outputTree->Branch("volume_id", &m_volumeId);
58  m_outputTree->Branch("boundary_id", &m_boundaryId);
59  m_outputTree->Branch("layer_id", &m_layerId);
60  m_outputTree->Branch("approach_id", &m_approachId);
61  m_outputTree->Branch("sensitive_id", &m_sensitiveId);
62 }
63 
64 FW::RootSimHitWriter::~RootSimHitWriter() {
65  if (m_outputFile) {
66  m_outputFile->Close();
67  }
68 }
69 
70 FW::ProcessCode FW::RootSimHitWriter::endRun() {
71  if (m_outputFile) {
72  m_outputFile->cd();
73  m_outputTree->Write();
74  ACTS_VERBOSE("Wrote hits to tree '" << m_cfg.treeName << "' in '"
75  << m_cfg.filePath << "'");
76  }
77  return ProcessCode::SUCCESS;
78 }
79 
80 FW::ProcessCode FW::RootSimHitWriter::writeT(const AlgorithmContext& ctx,
81  const FW::SimHitContainer& hits) {
82  if (not m_outputFile) {
83  ACTS_ERROR("Missing output file");
84  return ProcessCode::ABORT;
85  }
86 
87  // ensure exclusive access to tree/file while writing
88  std::lock_guard<std::mutex> lock(m_writeMutex);
89 
90  // Get the event number
91  m_eventId = ctx.eventNumber;
92  for (const auto& hit : hits) {
93  m_particleId = hit.particleId().value();
94  m_geometryId = hit.geometryId().value();
95  // write hit position
96  m_tx = hit.position4().x() / Acts::UnitConstants::mm;
97  m_ty = hit.position4().y() / Acts::UnitConstants::mm;
98  m_tz = hit.position4().z() / Acts::UnitConstants::mm;
99  m_tt = hit.position4().w() / Acts::UnitConstants::ns;
100  // write four-momentum before interaction
101  m_tpx = hit.momentum4Before().x() / Acts::UnitConstants::GeV;
102  m_tpy = hit.momentum4Before().y() / Acts::UnitConstants::GeV;
103  m_tpz = hit.momentum4Before().z() / Acts::UnitConstants::GeV;
104  m_te = hit.momentum4Before().w() / Acts::UnitConstants::GeV;
105  // write four-momentum change due to interaction
106  const auto delta4 = hit.momentum4After() - hit.momentum4Before();
107  m_deltapx = delta4.x() / Acts::UnitConstants::GeV;
108  m_deltapy = delta4.y() / Acts::UnitConstants::GeV;
109  m_deltapz = delta4.z() / Acts::UnitConstants::GeV;
110  m_deltae = delta4.w() / Acts::UnitConstants::GeV;
111  // write hit index along trajectory
112  m_index = hit.index();
113  // decoded geometry for simplicity
114  m_volumeId = hit.geometryId().volume();
115  m_boundaryId = hit.geometryId().boundary();
116  m_layerId = hit.geometryId().layer();
117  m_approachId = hit.geometryId().approach();
118  m_sensitiveId = hit.geometryId().sensitive();
119  // Fill the tree
120  m_outputTree->Fill();
121  }
122  return FW::ProcessCode::SUCCESS;
123 }