ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
KFParticle_truthAndDetTools.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file KFParticle_truthAndDetTools.cc
2 
3 #include <trackbase/TrkrHitSet.h>
5 
6 std::map<std::string, int> Use =
7  {
8  {"MVTX", 1},
9  {"INTT", 1},
10  {"TPC", 1},
11  {"EMCAL", 0},
12  {"OHCAL", 0},
13  {"IHCAL", 0}};
14 
16  : m_svtx_evalstack(nullptr)
17 {
18 } //Constructor
19 
21 
23 {
24  SvtxTrack *matched_track = NULL;
25 
26  for (SvtxTrackMap::Iter iter = trackmap->begin();
27  iter != trackmap->end();
28  ++iter)
29  {
30  if (iter->first == track_id) matched_track = iter->second;
31  }
32 
33  return matched_track;
34 }
35 
37 {
38  SvtxVertex *matched_vertex = NULL;
39 
40  for (SvtxVertexMap::Iter iter = vertexmap->begin();
41  iter != vertexmap->end();
42  ++iter)
43  {
44  if (iter->first == vertex_id) matched_vertex = iter->second;
45  }
46 
47  return matched_vertex;
48 }
49 
51 {
52  if (!m_svtx_evalstack)
53  {
54  m_svtx_evalstack = new SvtxEvalStack(topNode);
60  }
61 
62  m_svtx_evalstack->next_event(topNode);
63 
64  TrkrDefs::cluskey clusKey = *thisTrack->begin_cluster_keys();
66 
67  return particle;
68 }
69 
70 void KFParticle_truthAndDetTools::initializeTruthBranches(TTree *m_tree, int daughter_id, std::string daughter_number, bool m_constrain_to_vertex_truthMatch)
71 {
72  m_tree->Branch(TString(daughter_number) + "_true_ID", &m_true_daughter_id[daughter_id], TString(daughter_number) + "_true_ID/I");
73  if (m_constrain_to_vertex_truthMatch) m_tree->Branch(TString(daughter_number) + "_true_IP", &m_true_daughter_ip[daughter_id], TString(daughter_number) + "_true_IP/F");
74  if (m_constrain_to_vertex_truthMatch) m_tree->Branch(TString(daughter_number) + "_true_IP_xy", &m_true_daughter_ip_xy[daughter_id], TString(daughter_number) + "_true_IP_xy/F");
75  m_tree->Branch(TString(daughter_number) + "_true_px", &m_true_daughter_px[daughter_id], TString(daughter_number) + "_true_px/F");
76  m_tree->Branch(TString(daughter_number) + "_true_py", &m_true_daughter_py[daughter_id], TString(daughter_number) + "_true_py/F");
77  m_tree->Branch(TString(daughter_number) + "_true_pz", &m_true_daughter_pz[daughter_id], TString(daughter_number) + "_true_pz/F");
78  m_tree->Branch(TString(daughter_number) + "_true_p", &m_true_daughter_p[daughter_id], TString(daughter_number) + "_true_p/F");
79  m_tree->Branch(TString(daughter_number) + "_true_pT", &m_true_daughter_pt[daughter_id], TString(daughter_number) + "_true_pT/F");
80  m_tree->Branch(TString(daughter_number) + "_true_EV_x", &m_true_daughter_vertex_x[daughter_id], TString(daughter_number) + "_true_EV_x/F");
81  m_tree->Branch(TString(daughter_number) + "_true_EV_y", &m_true_daughter_vertex_y[daughter_id], TString(daughter_number) + "_true_EV_y/F");
82  m_tree->Branch(TString(daughter_number) + "_true_EV_z", &m_true_daughter_vertex_z[daughter_id], TString(daughter_number) + "_true_EV_z/F");
83  if (m_constrain_to_vertex_truthMatch)
84  {
85  m_tree->Branch(TString(daughter_number) + "_true_PV_x", &m_true_daughter_pv_x[daughter_id], TString(daughter_number) + "_true_pv_x/F");
86  m_tree->Branch(TString(daughter_number) + "_true_PV_y", &m_true_daughter_pv_y[daughter_id], TString(daughter_number) + "_true_pv_y/F");
87  m_tree->Branch(TString(daughter_number) + "_true_PV_z", &m_true_daughter_pv_z[daughter_id], TString(daughter_number) + "_true_pv_x/F");
88  }
89  m_tree->Branch(TString(daughter_number) + "_true_track_history_PDG_ID", &m_true_daughter_track_history_PDG_ID[daughter_id]);
90  m_tree->Branch(TString(daughter_number) + "_true_track_history_PDG_mass", &m_true_daughter_track_history_PDG_mass[daughter_id]);
91  m_tree->Branch(TString(daughter_number) + "_true_track_history_px", &m_true_daughter_track_history_px[daughter_id]);
92  m_tree->Branch(TString(daughter_number) + "_true_track_history_py", &m_true_daughter_track_history_py[daughter_id]);
93  m_tree->Branch(TString(daughter_number) + "_true_track_history_pz", &m_true_daughter_track_history_pz[daughter_id]);
94  m_tree->Branch(TString(daughter_number) + "_true_track_history_pE", &m_true_daughter_track_history_pE[daughter_id]);
95  m_tree->Branch(TString(daughter_number) + "_true_track_history_pT", &m_true_daughter_track_history_pT[daughter_id]);
96 }
97 
98 void KFParticle_truthAndDetTools::fillTruthBranch(PHCompositeNode *topNode, TTree */*m_tree*/, KFParticle daughter, int daughter_id, KFParticle vertex, bool m_constrain_to_vertex_truthMatch)
99 {
100  float true_px, true_py, true_pz, true_p, true_pt;
101 
102  PHNodeIterator nodeIter(topNode);
103  PHNode *findNode = dynamic_cast<PHNode*>(nodeIter.findFirst(m_trk_map_node_name_nTuple));
104  if (findNode)
105  {
106  dst_trackmap = findNode::getClass<SvtxTrackMap>(topNode, m_trk_map_node_name_nTuple);
107  }
108  else
109  {
110  std::cout << "KFParticle truth matching: " << m_trk_map_node_name_nTuple << " does not exist" << std::endl;
111  }
112  findNode = dynamic_cast<PHNode*>(nodeIter.findFirst(m_vtx_map_node_name_nTuple));
113  if (findNode)
114  {
115  dst_vertexmap = findNode::getClass<SvtxVertexMap>(topNode, m_vtx_map_node_name_nTuple);
116  }
117  else
118  {
119  std::cout << "KFParticle truth matching: " << m_vtx_map_node_name_nTuple << " does not exist" << std::endl;
120  }
121 
122  track = getTrack(daughter.Id(), dst_trackmap);
123  g4particle = getTruthTrack(track, topNode);
124 
125  bool isParticleValid = g4particle == nullptr ? 0 : 1;
126 
127  true_px = isParticleValid ? (Float_t) g4particle->get_px() : 0.;
128  true_py = isParticleValid ? (Float_t) g4particle->get_py() : 0.;
129  true_pz = isParticleValid ? (Float_t) g4particle->get_pz() : 0.;
130  true_p = sqrt(pow(true_px, 2) + pow(true_py, 2) + pow(true_pz, 2));
131  true_pt = sqrt(pow(true_px, 2) + pow(true_py, 2));
132 
133  m_true_daughter_px[daughter_id] = true_px;
134  m_true_daughter_py[daughter_id] = true_py;
135  m_true_daughter_pz[daughter_id] = true_pz;
136  m_true_daughter_p[daughter_id] = true_p;
137  m_true_daughter_pt[daughter_id] = true_pt;
138  m_true_daughter_id[daughter_id] = isParticleValid ? g4particle->get_pid() : 0;
139 
140  if (isParticleValid)
141  {
142  g4vertex_point = trutheval->get_vertex(g4particle);
143  }
144 
145  m_true_daughter_vertex_x[daughter_id] = isParticleValid ? g4vertex_point->get_x() : 0.;
146  m_true_daughter_vertex_y[daughter_id] = isParticleValid ? g4vertex_point->get_y() : 0.;
147  m_true_daughter_vertex_z[daughter_id] = isParticleValid ? g4vertex_point->get_z() : 0.;
148 
149  if (m_constrain_to_vertex_truthMatch)
150  {
151  //Calculate true DCA
152  SvtxVertex* recoVertex = getVertex(vertex.Id(), dst_vertexmap);
153  PHG4VtxPoint* truePoint = vertexeval->max_truth_point_by_ntracks(recoVertex);
154 
155  KFParticle trueKFParticleVertex;
156 
157  float f_vertexParameters[6] = {0};
158 
159  if (truePoint == NULL)
160  {
161  std::cout << "KFParticle truth matching: This event has no PHG4VtxPoint information!\n";
162  std::cout << "Your truth track DCA will be measured wrt a reconstructed vertex!" << std::endl;
163 
164  f_vertexParameters[0] = recoVertex->get_x();
165  f_vertexParameters[1] = recoVertex->get_y();
166  f_vertexParameters[2] = recoVertex->get_z();
167  }
168  else
169  {
170  f_vertexParameters[0] = truePoint->get_x();
171  f_vertexParameters[1] = truePoint->get_y();
172  f_vertexParameters[2] = truePoint->get_z();
173  }
174 
175  float f_vertexCovariance[21] = {0};
176 
177  trueKFParticleVertex.Create(f_vertexParameters, f_vertexCovariance, 0, -1);
178 
179  KFParticle trueKFParticle;
180 
181  float f_trackParameters[6] = {m_true_daughter_vertex_x[daughter_id],
182  m_true_daughter_vertex_y[daughter_id],
183  m_true_daughter_vertex_z[daughter_id],
184  true_px,
185  true_py,
186  true_pz};
187 
188  float f_trackCovariance[21] = {0};
189 
190  trueKFParticle.Create(f_trackParameters, f_trackCovariance, 1, -1);
191 
192  m_true_daughter_ip[daughter_id] = trueKFParticle.GetDistanceFromVertex(trueKFParticleVertex);
193  m_true_daughter_ip_xy[daughter_id] = trueKFParticle.GetDistanceFromVertexXY(trueKFParticleVertex);
194 
195  m_true_daughter_pv_x[daughter_id] = truePoint == NULL ? -99. : truePoint->get_x();
196  m_true_daughter_pv_y[daughter_id] = truePoint == NULL ? -99. : truePoint->get_y();
197  m_true_daughter_pv_z[daughter_id] = truePoint == NULL ? -99. : truePoint->get_z();
198 
199  }
200 }
201 
202 void KFParticle_truthAndDetTools::fillHepMCBranch(HepMC::GenParticle *particle, int daughter_id)
203 {
204  HepMC::FourVector myFourVector = particle->momentum();
205 
206  m_true_daughter_track_history_PDG_ID[daughter_id].push_back(particle->pdg_id());
207  m_true_daughter_track_history_PDG_mass[daughter_id].push_back((Float_t) particle->generatedMass());
208  m_true_daughter_track_history_px[daughter_id].push_back((Float_t) myFourVector.px());
209  m_true_daughter_track_history_py[daughter_id].push_back((Float_t) myFourVector.py());
210  m_true_daughter_track_history_pz[daughter_id].push_back((Float_t) myFourVector.pz());
211  m_true_daughter_track_history_pE[daughter_id].push_back((Float_t) myFourVector.e());
212  m_true_daughter_track_history_pT[daughter_id].push_back((Float_t) myFourVector.perp());
213 }
214 
215 int KFParticle_truthAndDetTools::getHepMCInfo(PHCompositeNode *topNode, TTree */*m_tree*/, KFParticle daughter, int daughter_id)
216 {
217  //Make dummy particle for null pointers and missing nodes
218  HepMC::GenParticle* dummyParticle = new HepMC::GenParticle();
219  HepMC::FourVector dummyFourVector(0,0,0,0);
220  dummyParticle->set_momentum(dummyFourVector);
221  dummyParticle->set_pdg_id(0);
222  dummyParticle->set_generated_mass(0.);
223 
224  PHNodeIterator nodeIter(topNode);
225  PHNode *findNode = dynamic_cast<PHNode*>(nodeIter.findFirst(m_trk_map_node_name_nTuple));
226  if (findNode)
227  {
228  dst_trackmap = findNode::getClass<SvtxTrackMap>(topNode, m_trk_map_node_name_nTuple);
229  }
230  else
231  {
232  std::cout << "KFParticle truth matching: " << m_trk_map_node_name_nTuple << " does not exist" << std::endl;
233  }
234 
235  track = getTrack(daughter.Id(), dst_trackmap);
236  g4particle = getTruthTrack(track, topNode);
237 
238  bool isParticleValid = g4particle == nullptr ? 0 : 1;
239 
240  if (!isParticleValid)
241  {
242  std::cout << "KFParticle truth matching: this track is a ghost" << std::endl;
243  fillHepMCBranch(dummyParticle, daughter_id);
244  return 0;
245  }
246 
247  m_geneventmap = findNode::getClass<PHHepMCGenEventMap>(topNode, "PHHepMCGenEventMap");
248  if (!m_geneventmap)
249  {
250  std::cout << "KFParticle truth matching: Missing node PHHepMCGenEventMap" << std::endl;
251  std::cout << "You will have no mother information" << std::endl;
252  fillHepMCBranch(dummyParticle, daughter_id);
253  return 0;
254  }
255 
256  m_genevt = m_geneventmap->get(1);
257  if (!m_genevt)
258  {
259  std::cout << "KFParticle truth matching: Missing node PHHepMCGenEvent" << std::endl;
260  std::cout << "You will have no mother information" << std::endl;
261  fillHepMCBranch(dummyParticle, daughter_id);
262  return 0;
263  }
264 
265  HepMC::GenEvent* theEvent = m_genevt->getEvent();
266  HepMC::GenParticle* prevParticle = nullptr;
267 
268  int forbiddenPDGIDs[] = {21,22}; //Stop tracing history when we reach quarks, gluons and photons
269 
270  for (HepMC::GenEvent::particle_const_iterator p = theEvent->particles_begin(); p != theEvent->particles_end(); ++p)
271  {
272  if (((*p)->barcode() == g4particle->get_barcode()))
273  {
274  prevParticle = *p;
275  while (!prevParticle->is_beam())
276  {
277  bool breakOut = false;
278  for (HepMC::GenVertex::particle_iterator mother = prevParticle->production_vertex()->particles_begin(HepMC::parents);
279  mother != prevParticle->production_vertex()->particles_end(HepMC::parents); ++mother)
280  {
281  if (std::find(std::begin(forbiddenPDGIDs), std::end(forbiddenPDGIDs),
282  abs((*mother)->pdg_id())) != std::end(forbiddenPDGIDs)) {breakOut = true; break;}
283 
284  fillHepMCBranch((*mother), daughter_id);
285  prevParticle = *mother;
286  }
287  if (breakOut) break;
288  }
289  }
290  }
291 
292  return 0;
293 }//End of function
294 
295 void KFParticle_truthAndDetTools::initializeCaloBranches(TTree *m_tree, int daughter_id, std::string daughter_number)
296 {
297  m_tree->Branch(TString(daughter_number) + "_EMCAL_DeltaPhi", &detector_emcal_deltaphi[daughter_id], TString(daughter_number) + "_EMCAL_DeltaPhi/F");
298  m_tree->Branch(TString(daughter_number) + "_EMCAL_DeltaEta", &detector_emcal_deltaeta[daughter_id], TString(daughter_number) + "_EMCAL_DeltaEta/F");
299  m_tree->Branch(TString(daughter_number) + "_EMCAL_energy_3x3", &detector_emcal_energy_3x3[daughter_id], TString(daughter_number) + "_EMCAL_energy_3x3/F");
300  m_tree->Branch(TString(daughter_number) + "_EMCAL_energy_5x5", &detector_emcal_energy_5x5[daughter_id], TString(daughter_number) + "_EMCAL_energy_5x5/F");
301  m_tree->Branch(TString(daughter_number) + "_EMCAL_energy_cluster", &detector_emcal_cluster_energy[daughter_id], TString(daughter_number) + "_EMCAL_energy_cluster/F");
302  m_tree->Branch(TString(daughter_number) + "_IHCAL_DeltaPhi", &detector_ihcal_deltaphi[daughter_id], TString(daughter_number) + "_IHCAL_DeltaPhi/F");
303  m_tree->Branch(TString(daughter_number) + "_IHCAL_DeltaEta", &detector_ihcal_deltaeta[daughter_id], TString(daughter_number) + "_IHCAL_DeltaEta/F");
304  m_tree->Branch(TString(daughter_number) + "_IHCAL_energy_3x3", &detector_ihcal_energy_3x3[daughter_id], TString(daughter_number) + "_IHCAL_energy_3x3/F");
305  m_tree->Branch(TString(daughter_number) + "_IHCAL_energy_5x5", &detector_ihcal_energy_5x5[daughter_id], TString(daughter_number) + "_IHCAL_energy_5x5/F");
306  m_tree->Branch(TString(daughter_number) + "_IHCAL_energy_cluster", &detector_ihcal_cluster_energy[daughter_id], TString(daughter_number) + "_IHCAL_energy_cluster/F");
307  m_tree->Branch(TString(daughter_number) + "_OHCAL_DeltaPhi", &detector_ohcal_deltaphi[daughter_id], TString(daughter_number) + "_OHCAL_DeltaEta/F");
308  m_tree->Branch(TString(daughter_number) + "_OHCAL_DeltaEta", &detector_ohcal_deltaeta[daughter_id], TString(daughter_number) + "_OHCAL_DeltaEta/F");
309  m_tree->Branch(TString(daughter_number) + "_OHCAL_energy_3x3", &detector_ohcal_energy_3x3[daughter_id], TString(daughter_number) + "_OHCAL_energy_3x3/F");
310  m_tree->Branch(TString(daughter_number) + "_OHCAL_energy_5x5", &detector_ohcal_energy_5x5[daughter_id], TString(daughter_number) + "_OHCAL_energy_5x5/F");
311  m_tree->Branch(TString(daughter_number) + "_OHCAL_energy_cluster", &detector_ohcal_cluster_energy[daughter_id], TString(daughter_number) + "_OHCAL_energy_cluster/F");
312 }
313 
315  TTree */*m_tree*/, KFParticle daughter, int daughter_id)
316 {
317  PHNodeIterator nodeIter(topNode);
318  PHNode *findNode = dynamic_cast<PHNode*>(nodeIter.findFirst(m_trk_map_node_name_nTuple));
319  if (findNode)
320  {
321  dst_trackmap = findNode::getClass<SvtxTrackMap>(topNode, m_trk_map_node_name_nTuple);
322  }
323  else
324  {
325  std::cout << "KFParticle truth matching: " << m_trk_map_node_name_nTuple << " does not exist" << std::endl;
326  }
327 
328  track = getTrack(daughter.Id(), dst_trackmap);
329 
335 
341 
347 }
348 
349 void KFParticle_truthAndDetTools::initializeDetectorBranches(TTree *m_tree, int daughter_id, std::string daughter_number)
350 {
351  m_tree->Branch(TString(daughter_number) + "_local_x", &detector_local_x[daughter_id]);
352  m_tree->Branch(TString(daughter_number) + "_local_y", &detector_local_y[daughter_id]);
353  m_tree->Branch(TString(daughter_number) + "_local_z", &detector_local_z[daughter_id]);
354  m_tree->Branch(TString(daughter_number) + "_layer", &detector_layer[daughter_id]);
355 
356  for (auto const &subdetector : Use)
357  {
358  if (subdetector.second) initializeSubDetectorBranches(m_tree, subdetector.first, daughter_id, daughter_number);
359  }
360 }
361 
362 void KFParticle_truthAndDetTools::initializeSubDetectorBranches(TTree *m_tree, std::string detectorName, int daughter_id, std::string daughter_number)
363 {
364  if (detectorName == "MVTX")
365  {
366  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_staveID", &mvtx_staveID[daughter_id]);
367  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_chipID", &mvtx_chipID[daughter_id]);
368  }
369  if (detectorName == "INTT")
370  {
371  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_ladderZID", &intt_ladderZID[daughter_id]);
372  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_ladderPhiID", &intt_ladderPhiID[daughter_id]);
373  }
374  if (detectorName == "TPC")
375  {
376  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_sectorID", &tpc_sectorID[daughter_id]);
377  m_tree->Branch(TString(daughter_number) + "_" + TString(detectorName) + "_side", &tpc_side[daughter_id]);
378  }
379 }
380 
382  TTree */*m_tree*/, KFParticle daughter, int daughter_id)
383 {
384  PHNodeIterator nodeIter(topNode);
385 
386  PHNode *findNode = dynamic_cast<PHNode*>(nodeIter.findFirst(m_trk_map_node_name_nTuple));
387  if (findNode)
388  {
389  dst_trackmap = findNode::getClass<SvtxTrackMap>(topNode, m_trk_map_node_name_nTuple);
390  }
391  else
392  {
393  std::cout << "KFParticle truth matching: " << m_trk_map_node_name_nTuple << " does not exist" << std::endl;
394  }
395 
396  findNode = dynamic_cast<PHNode*>(nodeIter.findFirst("TRKR_CLUSTER"));
397  if (findNode)
398  {
399  dst_clustermap = findNode::getClass<TrkrClusterContainer>(topNode, "TRKR_CLUSTER");
400  }
401  else
402  {
403  std::cout << "KFParticle detector info: TRKR_CLUSTER does not exist" << std::endl;
404  }
405 
406  track = getTrack(daughter.Id(), dst_trackmap);
407 
409  iter != track->end_cluster_keys();
410  ++iter)
411  {
412  TrkrDefs::cluskey clusKey = *iter;
413  TrkrCluster *cluster = dst_clustermap->findCluster(clusKey);
414  const unsigned int trkrId = TrkrDefs::getTrkrId(clusKey);
415 
416  detector_local_x[daughter_id].push_back(cluster->getX());
417  detector_local_y[daughter_id].push_back(cluster->getY());
418  detector_local_z[daughter_id].push_back(cluster->getZ());
419  detector_layer[daughter_id].push_back(TrkrDefs::getLayer(clusKey));
420  unsigned int staveId, chipId, ladderZId, ladderPhiId, sectorId, side;
421  staveId = chipId = ladderZId = ladderPhiId = sectorId = side = -99;
422 
423  if (Use["MVTX"] && trkrId == TrkrDefs::mvtxId)
424  {
425  staveId = MvtxDefs::getStaveId(clusKey);
426  chipId = MvtxDefs::getChipId(clusKey);
427  }
428  else if (Use["INTT"] && trkrId == TrkrDefs::inttId)
429  {
430  ladderZId = InttDefs::getLadderZId(clusKey);
431  ladderPhiId = InttDefs::getLadderPhiId(clusKey);
432  }
433  else if (Use["TPC"] && trkrId == TrkrDefs::tpcId)
434  {
435  sectorId = TpcDefs::getSectorId(clusKey);
436  side = TpcDefs::getSide(clusKey);
437  }
438 
439  mvtx_staveID[daughter_id].push_back(staveId);
440  mvtx_chipID[daughter_id].push_back(chipId);
441  intt_ladderZID[daughter_id].push_back(ladderZId);
442  intt_ladderPhiID[daughter_id].push_back(ladderPhiId);
443  tpc_sectorID[daughter_id].push_back(sectorId);
444  tpc_side[daughter_id].push_back(side);
445  }
446 }
447 
449  TTree */*m_tree*/,
450  KFParticle motherParticle,
451  std::vector<KFParticle> daughters,
452  std::vector<KFParticle> intermediates)
453 {
455  std::vector<KFParticle> primaryVertices = kfpTupleTools.makeAllPrimaryVertices(topNode, m_vtx_map_node_name_nTuple);
456 
457  for (unsigned int i = 0; i < primaryVertices.size(); ++i)
458  {
459  allPV_x.push_back(primaryVertices[i].GetX());
460  allPV_y.push_back(primaryVertices[i].GetY());
461  allPV_z.push_back(primaryVertices[i].GetZ());
462 
463  allPV_mother_IP.push_back(motherParticle.GetDistanceFromVertex(primaryVertices[i]));
464  allPV_mother_IPchi2.push_back(motherParticle.GetDeviationFromVertex(primaryVertices[i]));
465 
466  for (unsigned int j = 0; j < daughters.size(); ++j)
467  {
468  allPV_daughter_IP[j].push_back(daughters[j].GetDistanceFromVertex(primaryVertices[i]));
469  allPV_daughter_IPchi2[j].push_back(daughters[j].GetDeviationFromVertex(primaryVertices[i]));
470  }
471 
472  for (unsigned int j = 0; j < intermediates.size(); ++j)
473  {
474  allPV_intermediates_IP[j].push_back(intermediates[j].GetDistanceFromVertex(primaryVertices[i]));
475  allPV_intermediates_IPchi2[j].push_back(intermediates[j].GetDeviationFromVertex(primaryVertices[i]));
476  }
477  }
478 }
479 
481 {
482  for (int i = 0; i < m_num_tracks_nTuple; ++i)
483  {
484  //Truth vectors
492 
493  //Detector vectors
494  detector_local_x[i].clear();
495  detector_local_y[i].clear();
496  detector_local_z[i].clear();
497  detector_layer[i].clear();
498  mvtx_staveID[i].clear();
499  mvtx_chipID[i].clear();
500  intt_ladderZID[i].clear();
501  intt_ladderPhiID[i].clear();
502  tpc_sectorID[i].clear();
503  tpc_side[i].clear();
504 
505  //PV vectors
506  allPV_daughter_IP[i].clear();
507  allPV_daughter_IPchi2[i].clear();
508  }
509 
510  allPV_x.clear();
511  allPV_y.clear();
512  allPV_z.clear();
513  allPV_z.clear();
514 
515  allPV_mother_IP.clear();
516  allPV_mother_IPchi2.clear();
517 
518  for (int i = 0; i < m_num_intermediate_states_nTuple; ++i)
519  {
520  allPV_intermediates_IP[i].clear();
521  allPV_intermediates_IPchi2[i].clear();
522  }
523 }