ATLAS Offline Software
SimHitTreeCreator.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 // class header include
6 #include "SimHitTreeCreator.h"
7 
8 // Athena includes
9 #include "StoreGate/ReadHandle.h"
11 
12 // HepMC includes
13 #include "AtlasHepMC/GenParticle.h"
14 
15 // ROOT includes
16 #include "TTree.h"
17 
18 ISF::SimHitTreeCreator::SimHitTreeCreator( const std::string& name, ISvcLocator* pSvcLocator )
19  : ::AthAlgorithm( name, pSvcLocator )
20  , m_thistSvc("THistSvc",name)
21  , m_validationStream("ISFSimHit")
22  , m_t_simHits(nullptr)
23  , m_pileup(-1)
24  , m_type(-1)
25  , m_id(-1)
26  , m_mother(-1)
27  , m_barcode(-1)
28  , m_time(-1.)
29  , m_drift(-1.)
30  , m_edeposit(-1.)
31  , m_momentum(-1.)
32  , m_theta(-1.)
33  , m_phi(-1.)
34  , m_eta(-1.)
35  , m_bcmHits("BCMHits")
36  , m_blmHits("BLMHits")
37  , m_pixHits("PixelHits")
38  , m_sctHits("SCT_Hits")
39  , m_trtHits("TRTUncompressedHits")
40  , m_pixPileupHits("PileupPixelHits")
41  , m_sctPileupHits("PileupSCT_Hits")
42  , m_trtPileupHits("PileupTRTUncompressedHits")
43  , m_mdtHits("MDT_Hits")
44  , m_rpcHits("RPC_Hits")
45  , m_tgcHits("TGC_Hits")
46  , m_cscHits("CSC_Hits")
47 {
48  declareProperty("ValidationStreamName",
49  m_validationStream = "ISFSimHit",
50  "Name of the output stream" );
51  declareProperty("THistService",
52  m_thistSvc,
53  "The THistSvc" );
54 
55  declareProperty("BCM_HitCollection", m_bcmHits );
56  declareProperty("BLM_HitCollection", m_blmHits );
57  declareProperty("PixelHitCollection", m_pixHits );
58  declareProperty("SCT_HitCollection", m_sctHits );
59  declareProperty("TRT_HitCollection", m_trtHits );
60  declareProperty("PileupPixelHitCollection", m_pixPileupHits );
61  declareProperty("PileupSCT_HitCollection", m_sctPileupHits );
62  declareProperty("PileupTRT_HitCollection", m_trtPileupHits );
63  declareProperty("MDT_HitCollection", m_mdtHits );
64  declareProperty("RPC_HitCollection", m_rpcHits );
65  declareProperty("TGC_HitCollection", m_tgcHits );
66  declareProperty("CSC_HitCollection", m_cscHits );
67 
68 }
69 
70 // Destructor
73 {}
74 
77  // retrieve the histogram service
78  ATH_CHECK( m_thistSvc.retrieve() );
79 
80  ATH_CHECK(this->createSimHitsTree());
81 
82  ATH_CHECK( m_bcmHits.initialize() );
83  ATH_CHECK( m_blmHits.initialize() );
84  ATH_CHECK( m_pixHits.initialize() );
85  ATH_CHECK( m_sctHits.initialize() );
86  ATH_CHECK( m_trtHits.initialize() );
87  ATH_CHECK( m_pixPileupHits.initialize() );
88  ATH_CHECK( m_sctPileupHits.initialize() );
89  ATH_CHECK( m_trtPileupHits.initialize() );
90  ATH_CHECK( m_mdtHits.initialize() );
91  ATH_CHECK( m_rpcHits.initialize() );
92  ATH_CHECK( m_tgcHits.initialize() );
93  ATH_CHECK( m_cscHits.initialize() );
94  return StatusCode::SUCCESS;
95 }
96 
98  ATH_CHECK(this->fillSimHitsTree());
99  return StatusCode::SUCCESS;
100 }
101 
104  // Create the prefix of histogram names for the THistSvc
105  //std::string prefix = "/" + m_validationStream + "/";
106  const char *treeName="simhits";
107  const std::string prefix = "/" + m_validationStream + "/"+ treeName;
108  m_t_simHits = new TTree( treeName, treeName );
109  m_t_simHits->Branch("pileup" , &m_pileup , "pileup/I" );
110  m_t_simHits->Branch("type" , &m_type , "type/I" );
111  m_t_simHits->Branch("id" , &m_id , "id/I" );
112  m_t_simHits->Branch("mother" , &m_mother , "mother/I" );
113  m_t_simHits->Branch("barcode" , &m_barcode , "barcode/I" );
114  m_t_simHits->Branch("momentum" , &m_momentum , "momentum/F" );
115  m_t_simHits->Branch("time" , &m_time , "time/F" );
116  m_t_simHits->Branch("drift" , &m_drift , "dist/F" );
117  m_t_simHits->Branch("edeposit" , &m_edeposit , "edeposit/F" );
118  m_t_simHits->Branch("theta" , &m_theta , "theta/F" );
119  m_t_simHits->Branch("phi" , &m_phi , "phi/F" );
120  m_t_simHits->Branch("eta" , &m_eta , "eta/F" );
121 
122  // register the Tree to the THistSvc and return it's StatusCodes
123  ATH_CHECK( m_thistSvc->regTree( prefix, m_t_simHits) );
124 
125  return StatusCode::SUCCESS;
126 }
127 
130 {
131  // loop over collections
132 
133  SG::ReadHandle<MDTSimHitCollection> mdtHits(m_mdtHits);
134  if (mdtHits.isValid() && mdtHits->size()) {
136  while ( ih!=mdtHits->end()) {
137  m_type = 1;
138  m_id = (*ih).MDTid();
139  m_mother = (*ih).particleEncoding();
140  m_time = (*ih).globalTime();
141  m_drift = (*ih).driftRadius();
142  m_edeposit = (*ih).energyDeposit();
143  m_barcode = (*ih).truthBarcode();
144  HepMcParticleLink HMPL = (*ih).particleLink();
145  this->addHepMcParticleLinkInfoToTree(HMPL);
146 
147  ++ih;
148  while (ih!=mdtHits->end() && m_id==(*ih).MDTid() && m_barcode==(*ih).truthBarcode() ) {
149  // merge energy deposits and move on
150  m_edeposit += (*ih).energyDeposit();
151  ++ih;
152  }
153  m_t_simHits->Fill();
154  }
155  }
156 
157  SG::ReadHandle<RPCSimHitCollection> rpcHits(m_rpcHits);
158  if (rpcHits.isValid() && rpcHits->size()) {
160  while (ih!=rpcHits->end()) {
161  m_type = 2;
162  m_id = (*ih).RPCid();
163  m_mother = (*ih).particleEncoding();
164  m_time = (*ih).globalTime();
165  m_drift = 0.;
166  m_edeposit = (*ih).energyDeposit();
167  m_barcode = (*ih).truthBarcode();
168  HepMcParticleLink HMPL = (*ih).particleLink();
169  this->addHepMcParticleLinkInfoToTree(HMPL);
170 
171  ++ih;
172  while (ih!=rpcHits->end() && m_id==(*ih).RPCid() && m_barcode==(*ih).truthBarcode() ) {
173  // merge energy deposits and move on
174  m_edeposit += (*ih).energyDeposit();
175  ++ih;
176  }
177  m_t_simHits->Fill();
178  }
179  }
180 
181  SG::ReadHandle<TGCSimHitCollection> tgcHits(m_tgcHits);
182  if (tgcHits.isValid() && tgcHits->size()) {
184  while ( ih!=tgcHits->end()) {
185  m_type = 3;
186  m_id = (*ih).TGCid();
187  m_mother = (*ih).particleEncoding();
188  m_time = (*ih).globalTime();
189  m_drift = 0.;
190  m_edeposit = (*ih).energyDeposit();
191  m_barcode = (*ih).truthBarcode();
192  HepMcParticleLink HMPL = (*ih).particleLink();
193  this->addHepMcParticleLinkInfoToTree(HMPL);
194 
195  ++ih;
196  while (ih!=tgcHits->end() && m_id==(*ih).TGCid() && m_barcode==(*ih).truthBarcode() ) {
197  // merge energy deposits and move on
198  m_edeposit += (*ih).energyDeposit();
199  ++ih;
200  }
201  m_t_simHits->Fill();
202  }
203  }
204 
205  SG::ReadHandle<CSCSimHitCollection> cscHits(m_cscHits);
206  if (cscHits.isValid() && cscHits->size()) {
208  while ( ih!=cscHits->end()) {
209  m_type = 4;
210  m_id = (*ih).CSCid();
211  m_mother = (*ih).particleID();
212  m_time = (*ih).globalTime();
213  m_drift = 0.;
214  m_edeposit = (*ih).energyDeposit();
215  m_barcode = (*ih).truthBarcode();
216  HepMcParticleLink HMPL = (*ih).particleLink();
217  this->addHepMcParticleLinkInfoToTree(HMPL);
218 
219  ++ih;
220  while (ih!=cscHits->end() && m_id==(*ih).CSCid() && m_barcode==(*ih).truthBarcode() ) {
221  // merge energy deposits and move on
222  m_edeposit += (*ih).energyDeposit();
223  ++ih;
224  }
225  m_t_simHits->Fill();
226  }
227  }
228 
229  for (int ipileup=0;ipileup<2;ipileup++) {
230  m_pileup = ipileup;
231 
232  SG::ReadHandle<SiHitCollection> pixHits((ipileup==0) ? m_pixHits : m_pixPileupHits);
233  if (pixHits.isValid() && pixHits->size()) {
235  while (ih!=pixHits->end()) {
236  m_type = 5;
237  m_id = (*ih).identify();
238  HepMcParticleLink HMPL = (*ih).particleLink();
239  if (HMPL.isValid()) m_mother = (HMPL.cptr())->pdg_id();
240  else m_mother=0;
241  m_time = (*ih).meanTime();
242  m_drift = 0.;
243  m_edeposit = (*ih).energyLoss();
244  m_barcode = (*ih).truthBarcode();
245  this->addHepMcParticleLinkInfoToTree(HMPL);
246 
247  ++ih;
248  while (ih!=pixHits->end() && ((unsigned int)m_id)==(*ih).identify() && m_barcode==(*ih).truthBarcode() ) {
249  // merge energy deposits and move on
250  m_edeposit += (*ih).energyLoss();
251  ++ih;
252  }
253  m_t_simHits->Fill();
254  }
255  }
256 
257  SG::ReadHandle<SiHitCollection> sctHits((ipileup==0) ? m_sctHits : m_sctPileupHits);
258  if (sctHits.isValid() && sctHits->size()) {
260  while (ih!=sctHits->end()) {
261  m_type = 6;
262  m_id = (*ih).identify();
263  HepMcParticleLink HMPL = (*ih).particleLink();
264  if (HMPL.isValid()) m_mother = (HMPL.cptr())->pdg_id();
265  else m_mother=0;
266  m_time = (*ih).meanTime();
267  m_drift = 0.;
268  m_edeposit = (*ih).energyLoss();
269  m_barcode = (*ih).truthBarcode();
270  this->addHepMcParticleLinkInfoToTree(HMPL);
271 
272  ++ih;
273  while (ih!=sctHits->end() && ((unsigned int)m_id)==(*ih).identify() && m_barcode==(*ih).truthBarcode() ) {
274  // merge energy deposits and move on
275  m_edeposit += (*ih).energyLoss();
276  ++ih;
277  }
278  m_t_simHits->Fill();
279  }
280  }
281 
282  SG::ReadHandle<TRTUncompressedHitCollection> trtHits((ipileup==0) ? m_trtHits : m_trtPileupHits);
283  if (trtHits.isValid() && trtHits->size()) {
285  while ( ih!=trtHits->end()) {
286  m_type = 7;
287  m_id = (*ih).GetHitID();
288  m_mother = (*ih).GetParticleEncoding();
289  HepMcParticleLink HMPL = (*ih).particleLink();
290  m_time = (*ih).GetGlobalTime();
291  m_drift = 0.;
292  m_edeposit = (*ih).GetEnergyDeposit();
293  m_barcode = (*ih).truthBarcode();
294  this->addHepMcParticleLinkInfoToTree(HMPL);
295 
296  ++ih;
297  while (ih!=trtHits->end() && m_id==(*ih).GetHitID() && m_barcode==(*ih).truthBarcode() ) {
298  // merge energy deposits and move on
299  m_edeposit += (*ih).GetEnergyDeposit();
300  ++ih;
301  }
302  m_t_simHits->Fill();
303  }
304  }
305  }
306  return StatusCode::SUCCESS;
307 }
308 
309 //** Add information from HepMcParticleLink to TTree - common for all SimHit types */
311  if (HMPL.isValid()) {
312  auto t_mom=(HMPL.cptr())->momentum();
313  m_momentum = std::sqrt( t_mom.x()*t_mom.x()+t_mom.y()*t_mom.y()+t_mom.z()*t_mom.z());
314  m_eta= (HMPL.cptr())->momentum().eta();
315  m_theta= (HMPL.cptr())->momentum().theta();
316  m_phi = (HMPL.cptr())->momentum().phi();
317  }
318  else {
319  m_momentum=-1.;
320  m_theta=-1.;
321  m_eta=-10.;
322  m_phi = -10.;
323  }
324  return;
325 }
ISF::SimHitTreeCreator::m_pixHits
SG::ReadHandleKey< SiHitCollection > m_pixHits
Definition: SimHitTreeCreator.h:71
ISF::SimHitTreeCreator::m_thistSvc
ServiceHandle< ITHistSvc > m_thistSvc
the histogram service
Definition: SimHitTreeCreator.h:50
ISF::SimHitTreeCreator::m_blmHits
SG::ReadHandleKey< SiHitCollection > m_blmHits
Definition: SimHitTreeCreator.h:70
ISF::SimHitTreeCreator::m_sctPileupHits
SG::ReadHandleKey< SiHitCollection > m_sctPileupHits
Definition: SimHitTreeCreator.h:75
ISF::SimHitTreeCreator::initialize
virtual StatusCode initialize() override final
Athena algorithm's interface method initialize()
Definition: SimHitTreeCreator.cxx:76
ISF::SimHitTreeCreator::m_sctHits
SG::ReadHandleKey< SiHitCollection > m_sctHits
Definition: SimHitTreeCreator.h:72
ISF::SimHitTreeCreator::m_rpcHits
SG::ReadHandleKey< RPCSimHitCollection > m_rpcHits
Definition: SimHitTreeCreator.h:80
SG::ReadHandle
Definition: StoreGate/StoreGate/ReadHandle.h:70
AthCommonDataStore< AthCommonMsg< Algorithm > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
ISF::SimHitTreeCreator::~SimHitTreeCreator
virtual ~SimHitTreeCreator()
Destructor.
Definition: SimHitTreeCreator.cxx:72
ISF::SimHitTreeCreator::addHepMcParticleLinkInfoToTree
void addHepMcParticleLinkInfoToTree(HepMcParticleLink &HMPL)
Definition: SimHitTreeCreator.cxx:310
ISF::SimHitTreeCreator::SimHitTreeCreator
SimHitTreeCreator(const std::string &name, ISvcLocator *pSvcLocator)
Constructor with parameters.
Definition: SimHitTreeCreator.cxx:18
ISF::SimHitTreeCreator::m_pixPileupHits
SG::ReadHandleKey< SiHitCollection > m_pixPileupHits
Definition: SimHitTreeCreator.h:74
AtlasHitsVector::begin
const_iterator begin() const
Definition: AtlasHitsVector.h:131
GenParticle.h
AtlasHitsVector::const_iterator
CONT::const_iterator const_iterator
Definition: AtlasHitsVector.h:43
ISF::SimHitTreeCreator::m_trtPileupHits
SG::ReadHandleKey< TRTUncompressedHitCollection > m_trtPileupHits
Definition: SimHitTreeCreator.h:76
ISF::SimHitTreeCreator::m_trtHits
SG::ReadHandleKey< TRTUncompressedHitCollection > m_trtHits
Definition: SimHitTreeCreator.h:73
ISF::SimHitTreeCreator::createSimHitsTree
StatusCode createSimHitsTree()
Create the simhits tree - validation mode only.
Definition: SimHitTreeCreator.cxx:103
m_type
TokenType m_type
the type
Definition: TProperty.cxx:44
ISF::SimHitTreeCreator::m_cscHits
SG::ReadHandleKey< CSCSimHitCollection > m_cscHits
Definition: SimHitTreeCreator.h:82
ParticleGun_EoverP_Config.momentum
momentum
Definition: ParticleGun_EoverP_Config.py:63
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
checkCorrelInHIST.prefix
dictionary prefix
Definition: checkCorrelInHIST.py:391
dumpFileToPlots.treeName
string treeName
Definition: dumpFileToPlots.py:20
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
AthAlgorithm
Definition: AthAlgorithm.h:47
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
ISF::SimHitTreeCreator::execute
virtual StatusCode execute() override final
Athena algorithm's interface method execute()
Definition: SimHitTreeCreator.cxx:97
ISF::SimHitTreeCreator::m_validationStream
std::string m_validationStream
validation THist stream name
Definition: SimHitTreeCreator.h:51
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:228
SimHitTreeCreator.h
ISF::SimHitTreeCreator::fillSimHitsTree
StatusCode fillSimHitsTree()
Fill the simhits tree - validation mode only.
Definition: SimHitTreeCreator.cxx:129
ISF::SimHitTreeCreator::m_mdtHits
SG::ReadHandleKey< MDTSimHitCollection > m_mdtHits
Definition: SimHitTreeCreator.h:79
AtlasHitsVector::end
const_iterator end() const
Definition: AtlasHitsVector.h:134
ISF::SimHitTreeCreator::m_bcmHits
SG::ReadHandleKey< SiHitCollection > m_bcmHits
Definition: SimHitTreeCreator.h:69
AtlasHitsVector::size
size_type size() const
Definition: AtlasHitsVector.h:143
ReadHandle.h
Handle class for reading from StoreGate.
ISF::SimHitTreeCreator::m_tgcHits
SG::ReadHandleKey< TGCSimHitCollection > m_tgcHits
Definition: SimHitTreeCreator.h:81