BOSS 7.0.1
BESIII Offline Software System
Loading...
Searching...
No Matches
NavigationTestAlg Class Reference

#include <NavigationTestAlg.h>

+ Inheritance diagram for NavigationTestAlg:

Public Member Functions

 NavigationTestAlg (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 
 NavigationTestAlg (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 

Detailed Description

Constructor & Destructor Documentation

◆ NavigationTestAlg() [1/2]

NavigationTestAlg::NavigationTestAlg ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 28 of file NavigationTestAlg.cxx.

28 : Algorithm(name, pSvcLocator)
29{
30 declareProperty("OutputFileName", m_fout="histo.root");
31}

◆ NavigationTestAlg() [2/2]

NavigationTestAlg::NavigationTestAlg ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Member Function Documentation

◆ execute() [1/2]

StatusCode NavigationTestAlg::execute ( )

Definition at line 70 of file NavigationTestAlg.cxx.

70 {
71 MsgStream log(msgSvc(), name());
72
73 // Get EventNavigator from the TDS
74 SmartDataPtr<EventNavigator> navigator (eventSvc(),"/Event/Navigator");
75 if( ! navigator )
76 {
77 log << MSG::ERROR << " Unable to retrieve EventNavigator" << endreq;
78 return StatusCode::FAILURE;
79 }
80
81 log << MSG::INFO << "EventNavigator object" << endreq;
82 navigator->Print();
83
84
85 // ======== Analyse Monte-Carlo information using navigator
86 log << MSG::INFO << "=======================================================================" << endreq;
87 log << MSG::INFO << "MC Particles ===============================================================" << endreq;
88
89 // Get McParticle collection
90 SmartDataPtr<McParticleCol> mcParticles(eventSvc(),"/Event/MC/McParticleCol");
91 if( ! mcParticles )
92 {
93 log << MSG::ERROR << " Unable to retrieve McParticleCol" << endreq;
94 return StatusCode::FAILURE;
95 }
96
97 // For each McParticle
98 for( McParticleCol::iterator it = mcParticles->begin(); it != mcParticles->end(); it++ )
99 {
100 // get charge of McParticles
101 int pdg_code = (*it)->particleProperty();
102
103 // get true momentum value
104 double true_mom = (*it)->initialFourMomentum().vect().mag();
105
106 // fill true momentum
107 m_histo[5]->Fill(true_mom);
108
109 log<<MSG::INFO<<"Retrieved McParticle # "<<(*it)->trackIndex()
110 << " PDG " << pdg_code << " of true momentum "
111 << true_mom <<" GeV/c " <<endreq;
112
113 // get list of reconstructed tracks, which correspond to this particle, using EventNavigator
114 RecMdcTrackVector tracks = navigator->getMdcTracks(*it);
115 RecMdcKalTrackVector ktracks = navigator->getMdcKalTracks(*it);
116
117 // fill the distribution of number of rec. tracks corresponding to single particle
118 m_histo[2]->Fill(tracks.size());
119
120 log << MSG::INFO << " Found " << tracks.size() << " tracks and " << ktracks.size() << " Kalman tracks" << endreq;
121
122 // for each track
123 for(unsigned int i=0; i<ktracks.size(); i++)
124 {
125 // Hep3Vector p = momentum(tracks[i]);
126 double momentum = ktracks[i]->p();
127
128 log << MSG::INFO << "\t Track # " << i
129 << " id = " << ktracks[i]->trackId()
130 << " Pid " << ktracks[i]->getPidType()
131 << " Ptot " << momentum << " GeV/c" << endreq;
132
133 // fill difference of true momentum and rec momentum for the _same_ track
134 m_histo[0]->Fill(true_mom-momentum);
135
136 // fill delta_p/p wrt p for the _same_ track
137 m_histo2[0]->Fill(true_mom, fabs(true_mom-momentum)/true_mom);
138 }
139
140 // get list of reconstructed showers, which correspond to this particle, using EventNavigator
141 RecEmcShowerVector showers = navigator->getEmcRecShowers(*it);
142
143 m_histo[8]->Fill(showers.size());
144
145 log << MSG::INFO << " Found " << showers.size() << " showers" << endreq;
146
147 for(unsigned int i=0; i<showers.size(); i++)
148 {
149 double true_energy = (*it)->initialFourMomentum().e();
150 double rec_energy = showers[i]->energy()*1000; // MeV
151
152 log << MSG::INFO << "\t Shower # " << i
153 << " Id = " << showers[i]->getShowerId().get_value()
154 << " E = " << showers[i]->energy()*1000 << " MeV " << endreq;
155
156 m_histo[12]->Fill(true_energy);
157 m_histo[13]->Fill(rec_energy);
158 }
159 }
160
161 log << MSG::INFO << "MDC Tracks ==============================================================" << endreq;
162 // ======== Analyse reconstructed track information using navigator
163
164 SmartDataPtr<RecMdcKalTrackCol> mdcKalTracks(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
165 if( ! mdcKalTracks )
166 {
167 log << MSG::ERROR << " Unable to retrieve MdcKalTrackCol" << endreq;
168 return StatusCode::SUCCESS;
169 }
170
171 for( RecMdcKalTrackCol::iterator it = mdcKalTracks->begin(); it != mdcKalTracks->end(); it++ )
172 {
173 McParticleVector particles = navigator->getMcParticles(*it);
174
175 log << MSG::INFO << "Retrieved " << particles.size() << " McParticles for for MdcKalTrack # "
176 << (*it)->trackId() << " of reconstructed momentum " << (*it)->p() << " GeV/c (PID="
177 << (*it)->getPidType() << ")" << endreq;
178
179 // fill rec momentum
180 m_histo[6]->Fill((*it)->p());
181
182 // fill the distribution of number of parent particles corresponding to this track
183 m_histo[1]->Fill(particles.size());
184
185 // for each parent particle
186 for(unsigned int i=0; i<particles.size(); i++)
187 {
188 int pdg_code = particles[i]->particleProperty();
189
190 // get true momentum value
191 double true_mom = particles[i]->initialFourMomentum().vect().mag();
192
193 // get relevance
194 int relevance = navigator->getMcParticleRelevance(*it, particles[i]);
195
196 // dump particle names, momenta and their relevance
197 log << MSG::INFO << "\t" << pdg_code << " momentum "
198 << true_mom << " relevance " << relevance << endreq;
199 }
200 }
201
202 // ======== Analyse reconstructed shower information using navigator
203 log << MSG::INFO << "EMC Showers ==============================================================" << endreq;
204 SmartDataPtr<RecEmcShowerCol> emcShowers(eventSvc(),"/Event/Recon/RecEmcShowerCol");
205 if( ! emcShowers )
206 {
207 log << MSG::ERROR << " Unable to retrieve EmcRecShowerCol" << endreq;
208 return StatusCode::SUCCESS;
209 }
210
211 int ij = 0;
212 for( RecEmcShowerCol::iterator it = emcShowers->begin(); it != emcShowers->end(); it++ )
213 {
214 McParticleVector particles = navigator->getMcParticles(*it);
215 log << MSG::INFO << "Retrieved McParticles for EmcShower # " << ij++
216 << " Id " << (*it)->getShowerId().get_value()
217 << " Energy = " << (*it)->energy() << endreq;
218
219 for(unsigned int i=0; i<particles.size(); i++)
220 {
221 int pdg_code = particles[i]->particleProperty();
222
223 // get true energy
224 double true_e = particles[i]->initialFourMomentum().e();
225
226 // get relevance
227 int relevance = navigator->getMcParticleRelevance(*it, particles[i]);
228
229 log << "\t Particle " << i << " PDG " << pdg_code << " E=" << true_e
230 << " Relevance " << relevance << endreq;
231 }
232
233 // fill the distribution of number of parent particles corresponding to this track
234 m_histo[7]->Fill(particles.size());
235 }
236
237 log << MSG::INFO << "=======================================================================" << endreq;
238
239 return StatusCode::SUCCESS;
240}
**********INTEGER nmxhep !maximum number of particles DOUBLE PRECISION vhep INTEGER jdahep COMMON hepevt $ !serial number $ !number of particles $ !status code $ !particle ident KF $ !parent particles $ !childreen particles $ !four momentum
std::vector< const RecMdcTrack * > RecMdcTrackVector
std::vector< const Event::McParticle * > McParticleVector
std::vector< const RecMdcKalTrack * > RecMdcKalTrackVector
std::vector< const RecEmcShower * > RecEmcShowerVector

◆ execute() [2/2]

StatusCode NavigationTestAlg::execute ( )

◆ finalize() [1/2]

StatusCode NavigationTestAlg::finalize ( )

Definition at line 244 of file NavigationTestAlg.cxx.

244 {
245
246 // Save histograms to ROOT file
247 TFile of(m_fout.c_str(),"RECREATE");
248 of.cd();
249 for(vector<TH1F*>::iterator it=m_histo.begin(); it!=m_histo.end(); it++)
250 (*it)->Write();
251 for(vector<TH2F*>::iterator it=m_histo2.begin(); it!=m_histo2.end(); it++)
252 (*it)->Write();
253 of.Close();
254
255 return StatusCode::SUCCESS;
256}

◆ finalize() [2/2]

StatusCode NavigationTestAlg::finalize ( )

◆ initialize() [1/2]

StatusCode NavigationTestAlg::initialize ( )

Definition at line 35 of file NavigationTestAlg.cxx.

35 {
36 MsgStream log(msgSvc(), name());
37
38 log << MSG::INFO << " Creating histograms" << endreq;
39
40 // create histograms
41 m_histo.push_back(new TH1F("deltap","#delta p, true momentum - rec momentum, MeV/c",600,-100,500)); //0
42
43 m_histo.push_back(new TH1F("n_parents","Number of parents for reconstructed track",20,0,20)); // 1
44 m_histo.push_back(new TH1F("n_tracks","Number of tracks for one MC particle",20,0,20)); // 2
45
46 m_histo.push_back(new TH1F("hit_parent","Number of parent MC hits for reconstructed hit",20,-10,10)); // 3
47 m_histo.push_back(new TH1F("hit_track","Number of reconstructed hits for MC hit",20,-10,10)); // 4
48
49 m_histo.push_back(new TH1F("true_p","true momentum, MeV/c",1500,0,1500)); // 5
50 m_histo.push_back(new TH1F("rec_p","reconstructed momentum, MeV/c",1500,0,1500)); // 6
51
52 m_histo.push_back(new TH1F("n_emcparents","# of parents for reconstructed shower",20,0,20)); // 7
53 m_histo.push_back(new TH1F("n_showers","# of showers for one true particle",20,0,20)); // 8
54
55 m_histo.push_back(new TH1F("pdg_miss","Pdg of particles with # trk==0 and #showers > 0",6000,-3000,3000)); //9
56
57 m_histo.push_back(new TH1F("dE_true","True energy - shower energy",40,-2000,2000)); //10
58 m_histo.push_back(new TH1F("dE_rec","shower energy - sqrt(Prec**2+m_true**2)",40,-2000,2000)); //11
59
60 m_histo.push_back(new TH1F("E_true","True energy, MeV",100,0,1000)); //12
61 m_histo.push_back(new TH1F("E_rec","Rec energy, MeV",100,0,1000)); //13
62
63 m_histo2.push_back(new TH2F("deltap_p","#delta p/p vs. p, MeV",100,0,1000,20,0,2));
64
65 return StatusCode::SUCCESS;
66}

◆ initialize() [2/2]

StatusCode NavigationTestAlg::initialize ( )

The documentation for this class was generated from the following files: