Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DMCThrownMatching_factory.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: DMCThrownMatching_factory.cc
4 // Created: Tue Aug 9 14:29:24 EST 2011
5 // Creator: pmatt (on Linux ifarml6 2.6.18-128.el5 x86_64)
6 //
7 
8 #ifdef VTRACE
9 #include "vt_user.h"
10 #endif
11 
13 #include "TMath.h"
14 
15 //------------------
16 // init
17 //------------------
19 {
20  dDebugLevel = 0;
21  dMaximumTOFMatchDistance = 10.0; //cm
22  dMaximumFCALMatchDistance = 10.0; //cm
24  dTargetCenter = 0.0; //cm
27 
28  return NOERROR;
29 }
30 
31 //------------------
32 // brun
33 //------------------
34 jerror_t DMCThrownMatching_factory::brun(jana::JEventLoop* locEventLoop, int32_t runnumber)
35 {
36  gPARMS->SetDefaultParameter("MCMATCH:DEBUG_LEVEL", dDebugLevel);
37  gPARMS->SetDefaultParameter("MCMATCH:MAX_TOF_DISTANCE", dMaximumTOFMatchDistance);
38  gPARMS->SetDefaultParameter("MCMATCH:MAX_FCAL_DISTANCE", dMaximumFCALMatchDistance);
39  gPARMS->SetDefaultParameter("MCMATCH:MAX_BCAL_ANGLE", dMaximumBCALMatchAngleDegrees);
40  gPARMS->SetDefaultParameter("MCMATCH:MAX_TOTAL_ERROR", dMaxTotalParticleErrorForMatch);
41  gPARMS->SetDefaultParameter("MCMATCH:MIN_TRACK_MATCH", dMinTrackMatchHitFraction);
42 
43  DApplication* locApplication = dynamic_cast<DApplication*>(locEventLoop->GetJApplication());
44  DGeometry* locGeometry = locApplication->GetDGeometry(locEventLoop->GetJEvent().GetRunNumber());
45  locGeometry->GetTargetZ(dTargetCenter);
46 
47  return NOERROR;
48 }
49 
50 //------------------
51 // evnt
52 //------------------
53 jerror_t DMCThrownMatching_factory::evnt(jana::JEventLoop* locEventLoop, uint64_t eventnumber)
54 {
55 #ifdef VTRACE
56  VT_TRACER("DMCThrownMatching_factory::evnt()");
57 #endif
58 
59  vector<const DMCThrown*> locMCThrowns;
60  locEventLoop->Get(locMCThrowns);
61 
62  if(locMCThrowns.empty())
63  return NOERROR;
64 
65  locMCThrowns.clear();
66  locEventLoop->Get(locMCThrowns, "FinalState");
67 
68  vector<const DMCThrown*> locMCThrowns_Charged;
69  vector<const DMCThrown*> locMCThrowns_Neutral;
70 
71  for(size_t loc_i = 0; loc_i < locMCThrowns.size(); ++loc_i)
72  {
73  if(ParticleCharge((Particle_t)locMCThrowns[loc_i]->type) == 0)
74  locMCThrowns_Neutral.push_back(locMCThrowns[loc_i]);
75  else
76  locMCThrowns_Charged.push_back(locMCThrowns[loc_i]);
77  }
78  if(dDebugLevel > 0)
79  cout << "input #thrown, ok charged # thrown, ok neutral # thrown = " << locMCThrowns.size() << ", " << locMCThrowns_Charged.size() << ", " << locMCThrowns_Neutral.size() << endl;
80 
81  vector<const DChargedTrack*> locChargedTracks;
82  locEventLoop->Get(locChargedTracks);
83 
84  vector<const DNeutralParticle*> locNeutralParticles;
85  locEventLoop->Get(locNeutralParticles);
86 
87  vector<const DChargedTrackHypothesis*> locChargedTrackHypotheses;
88  locEventLoop->Get(locChargedTrackHypotheses);
89 
90  vector<const DNeutralParticleHypothesis*> locNeutralParticleHypotheses;
91  locEventLoop->Get(locNeutralParticleHypotheses);
92 
93  DMCThrownMatching* locMCThrownMatching = new DMCThrownMatching();
94 
95  Find_GenReconMatches_BeamPhotons(locEventLoop, locMCThrownMatching);
96 
97  Find_GenReconMatches_ChargedHypo(locMCThrowns_Charged, locChargedTrackHypotheses, locMCThrownMatching);
98  Find_GenReconMatches_ChargedTrack(locChargedTracks, locMCThrownMatching);
99 
100  Find_GenReconMatches_NeutralHypo(locMCThrowns_Neutral, locNeutralParticleHypotheses, locMCThrownMatching);
101  Find_GenReconMatches_NeutralParticle(locNeutralParticles, locMCThrownMatching);
102 
103  Find_GenReconMatches_TOFPoints(locEventLoop, locMCThrownMatching);
104  Find_GenReconMatches_BCALShowers(locEventLoop, locMCThrownMatching);
105  Find_GenReconMatches_FCALShowers(locEventLoop, locMCThrownMatching);
106 
107  if(dDebugLevel > 0)
108  {
109  double locMatchFOM = 0.0;
110  cout << "Charged Track Matching Summary:" << endl;
111  for(size_t loc_i = 0; loc_i < locChargedTrackHypotheses.size(); ++loc_i)
112  {
113  double locP = locChargedTrackHypotheses[loc_i]->momentum().Mag();
114  double locTheta = locChargedTrackHypotheses[loc_i]->momentum().Theta()*180.0/TMath::Pi();
115  double locPhi = locChargedTrackHypotheses[loc_i]->momentum().Phi()*180.0/TMath::Pi();
116  Particle_t locPID = locChargedTrackHypotheses[loc_i]->PID();
117  cout << "charged info: " << locChargedTrackHypotheses[loc_i]->Get_TrackTimeBased()->candidateid << ", " << ParticleType(locPID) << ", " << locP << ", " << locTheta << ", " << locPhi << endl;
118  cout << "matched info: ";
119  const DMCThrown* locMCThrown = locMCThrownMatching->Get_MatchingMCThrown(locChargedTrackHypotheses[loc_i], locMatchFOM);
120  if(locMCThrown == NULL)
121  {
122  cout << "NO MATCH." << endl;
123  continue;
124  }
125  locP = locMCThrown->momentum().Mag();
126  locTheta = locMCThrown->momentum().Theta()*180.0/TMath::Pi();
127  locPhi = locMCThrown->momentum().Phi()*180.0/TMath::Pi();
128  locPID = locMCThrown->PID();
129  cout << ParticleType(locPID) << ", " << locP << ", " << locTheta << ", " << locPhi << ", " << locMatchFOM << endl;
130  }
131  cout << "Neutral Particle Matching Summary:" << endl;
132  for(size_t loc_i = 0; loc_i < locNeutralParticleHypotheses.size(); ++loc_i)
133  {
134  double locP = locNeutralParticleHypotheses[loc_i]->momentum().Mag();
135  double locTheta = locNeutralParticleHypotheses[loc_i]->momentum().Theta()*180.0/TMath::Pi();
136  double locPhi = locNeutralParticleHypotheses[loc_i]->momentum().Phi()*180.0/TMath::Pi();
137  Particle_t locPID = locNeutralParticleHypotheses[loc_i]->PID();
138  cout << "neutral info: " << locNeutralParticleHypotheses[loc_i] << ", " << ParticleType(locPID) << ", " << locP << ", " << locTheta << ", " << locPhi << endl;
139  cout << "matched info: ";
140  const DMCThrown* locMCThrown = locMCThrownMatching->Get_MatchingMCThrown(locNeutralParticleHypotheses[loc_i], locMatchFOM);
141  if(locMCThrown == NULL)
142  {
143  cout << "NO MATCH." << endl;
144  continue;
145  }
146  locP = locMCThrown->momentum().Mag();
147  locTheta = locMCThrown->momentum().Theta()*180.0/TMath::Pi();
148  locPhi = locMCThrown->momentum().Phi()*180.0/TMath::Pi();
149  locPID = locMCThrown->PID();
150  cout << ParticleType(locPID) << ", " << locP << ", " << locTheta << ", " << locPhi << ", " << locMatchFOM << endl;
151  }
152  cout << "Unmatched Charged DMCThrowns:" << endl;
153  for(size_t loc_i = 0; loc_i < locMCThrowns.size(); ++loc_i)
154  {
155  if(ParticleCharge(locMCThrowns[loc_i]->PID()) == 0)
156  continue;
157  if(locMCThrownMatching->Get_MatchingChargedTrack(locMCThrowns[loc_i], locMatchFOM) != NULL)
158  continue;
159  double locP = locMCThrowns[loc_i]->momentum().Mag();
160  double locTheta = locMCThrowns[loc_i]->momentum().Theta()*180.0/TMath::Pi();
161  double locPhi = locMCThrowns[loc_i]->momentum().Phi()*180.0/TMath::Pi();
162  Particle_t locPID = locMCThrowns[loc_i]->PID();
163  cout << "thrown info: " << ParticleType(locPID) << ", " << locP << ", " << locTheta << ", " << locPhi << endl;
164  }
165  }
166 
167  _data.push_back(locMCThrownMatching);
168 
169  return NOERROR;
170 }
171 
172 void DMCThrownMatching_factory::Find_GenReconMatches_BeamPhotons(JEventLoop* locEventLoop, DMCThrownMatching* locMCThrownMatching) const
173 {
174  vector<const DBeamPhoton*> locBeamPhotons;
175  locEventLoop->Get(locBeamPhotons);
176 
177  vector<const DBeamPhoton*> locBeamPhotons_MCGEN;
178  locEventLoop->Get(locBeamPhotons_MCGEN, "MCGEN");
179 
180  vector<const DBeamPhoton*> locBeamPhotons_TAGGEDMCGEN;
181  locEventLoop->Get(locBeamPhotons_TAGGEDMCGEN, "TAGGEDMCGEN");
182 
183  vector<const DBeamPhoton*> locBeamPhotons_TRUTH;
184  locEventLoop->Get(locBeamPhotons_TRUTH, "TRUTH");
185 
186  //first set MCGEN & tagged MCGEN info
187  locMCThrownMatching->Set_MCGENBeamPhoton(locBeamPhotons_MCGEN.empty() ? nullptr : locBeamPhotons_MCGEN[0]);
188  locMCThrownMatching->Set_TaggedMCGENBeamPhoton(locBeamPhotons_TAGGEDMCGEN.empty() ? nullptr : locBeamPhotons_TAGGEDMCGEN[0]);
189 
190  //recon photon, truth photon, delta-t
191  map<const DBeamPhoton*, deque<pair<const DBeamPhoton*, double> > > locPossibleMatches;
192  for(size_t loc_i = 0; loc_i < locBeamPhotons.size(); ++loc_i)
193  {
194  for(size_t loc_j = 0; loc_j < locBeamPhotons_TRUTH.size(); ++loc_j)
195  {
196  const DTAGMHit* locTAGMHit = NULL;
197  locBeamPhotons[loc_i]->GetSingle(locTAGMHit);
198  const DTAGHHit* locTAGHHit = NULL;
199  locBeamPhotons[loc_i]->GetSingle(locTAGHHit);
200 
201  const DTAGMHit* locTAGMHit_TRUTH = NULL;
202  locBeamPhotons_TRUTH[loc_j]->GetSingle(locTAGMHit_TRUTH);
203  const DTAGHHit* locTAGHHit_TRUTH = NULL;
204  locBeamPhotons_TRUTH[loc_j]->GetSingle(locTAGHHit_TRUTH);
205 
206  if(locTAGMHit != NULL)
207  {
208  if(locTAGMHit_TRUTH == NULL)
209  continue;
210  if((locTAGMHit->row != locTAGMHit_TRUTH->row) || (locTAGMHit->column != locTAGMHit_TRUTH->column))
211  continue;
212  }
213  else
214  {
215  if(locTAGHHit_TRUTH == NULL)
216  continue;
217  if(locTAGHHit->counter_id != locTAGHHit_TRUTH->counter_id)
218  continue;
219  }
220 
221  double locDeltaT = fabs(locBeamPhotons[loc_i]->time() - locBeamPhotons_TRUTH[loc_j]->time());
222  pair<const DBeamPhoton*, double> locMatchPair(locBeamPhotons_TRUTH[loc_j], locDeltaT);
223  locPossibleMatches[locBeamPhotons[loc_i]].push_back(locMatchPair);
224  }
225  }
226 
227  map<const DBeamPhoton*, const DBeamPhoton*> locBeamPhotonToTruthMap;
228  map<const DBeamPhoton*, const DBeamPhoton*> locBeamTruthToPhotonMap;
229 
230  map<const DBeamPhoton*, deque<pair<const DBeamPhoton*, double> > >::iterator locIterator = locPossibleMatches.begin();
231  for(; locIterator != locPossibleMatches.end(); ++locIterator)
232  {
233  deque<pair<const DBeamPhoton*, double> > locPairs = locIterator->second;
234  double locBestDeltaT = 9.9E9;
235  const DBeamPhoton* locBestBeamPhoton = NULL;
236  for(size_t loc_i = 0; loc_i < locPairs.size(); ++loc_i)
237  {
238  double locDeltaT = locPairs[loc_i].second;
239  if(locDeltaT >= locBestDeltaT)
240  continue;
241  locBestDeltaT = locDeltaT;
242  locBestBeamPhoton = locPairs[loc_i].first;
243  }
244 
245  if((dDebugLevel > 20) && (locBestBeamPhoton != NULL))
246  cout << "Photon match: delta-E, delta-t = " << fabs(locIterator->first->energy() - locBestBeamPhoton->energy()) << ", " << locBestDeltaT << endl;
247  locBeamPhotonToTruthMap[locIterator->first] = locBestBeamPhoton;
248  locBeamTruthToPhotonMap[locBestBeamPhoton] = locIterator->first;
249  }
250 
251  locMCThrownMatching->Set_BeamPhotonToTruthMap(locBeamPhotonToTruthMap);
252  locMCThrownMatching->Set_BeamTruthToPhotonMap(locBeamTruthToPhotonMap);
253 }
254 
255 void DMCThrownMatching_factory::Find_GenReconMatches_FCALShowers(JEventLoop* locEventLoop, DMCThrownMatching* locMCThrownMatching) const
256 {
257  vector<const DFCALTruthShower*> locFCALTruthShowers;
258  const DFCALTruthShower* locFCALTruthShower;
259  locEventLoop->Get(locFCALTruthShowers);
260 
261  vector<const DFCALShower*> locFCALShowers;
262  const DFCALShower* locFCALShower;
263  locEventLoop->Get(locFCALShowers);
264 
265  map<const DFCALShower*, pair<const DFCALTruthShower*, double> > locFCALShowerToTruthMap;
266  map<const DFCALTruthShower*, pair<const DFCALShower*, double> > locFCALTruthToShowerMap;
267 
268  size_t locBestFCALShowerIndex = 0, locBestFCALTruthShowerIndex = 0;
269  double locBestMatchDistance = 0.0, locMatchDistance;
270  while((!locFCALShowers.empty()) && (!locFCALTruthShowers.empty()))
271  {
272  bool locMatchFoundFlag = false;
273  double locLargestEnergy = 0.0;
274  for(size_t loc_i = 0; loc_i < locFCALTruthShowers.size(); ++loc_i)
275  {
276  locFCALTruthShower = locFCALTruthShowers[loc_i];
277  for(size_t loc_j = 0; loc_j < locFCALShowers.size(); ++loc_j)
278  {
279  locFCALShower = locFCALShowers[loc_j];
280 
281  DVector3 locFCALShowerPosition = locFCALShower->getPosition();
282 
283  //propagate truth shower from the FCAL face to the depth of the measured shower
284  double locDeltaZ = locFCALShowerPosition.Z() - locFCALTruthShower->z();
285  DVector3 locTrueMomentum(locFCALTruthShower->px(), locFCALTruthShower->py(), locFCALTruthShower->pz());
286  double locDeltaPathLength = locDeltaZ/cos(locTrueMomentum.Theta());
287  double locTrueX = locFCALTruthShower->x() + locDeltaPathLength*sin(locTrueMomentum.Theta())*cos(locTrueMomentum.Phi());
288  double locTrueY = locFCALTruthShower->y() + locDeltaPathLength*sin(locTrueMomentum.Theta())*sin(locTrueMomentum.Phi());
289  DVector3 locFCALTruthShowerPosition(locTrueX, locTrueY, locFCALShowerPosition.Z());
290 
291  locMatchDistance = (locFCALShowerPosition - locFCALTruthShowerPosition).Mag();
292  if(locMatchDistance > dMaximumFCALMatchDistance)
293  continue;
294 
295  //keep the shower with the largest energy
296  if(locFCALShower->getEnergy() > locLargestEnergy)
297  {
298  locBestMatchDistance = locMatchDistance;
299  locMatchFoundFlag = true;
300  locLargestEnergy = locFCALShower->getEnergy();
301  locBestFCALTruthShowerIndex = loc_i;
302  locBestFCALShowerIndex = loc_j;
303  }
304  }
305  }
306 
307  if(!locMatchFoundFlag) //no more good matches!1
308  break;
309 
310  locFCALTruthShower = locFCALTruthShowers[locBestFCALTruthShowerIndex];
311  locFCALShower = locFCALShowers[locBestFCALShowerIndex];
312 
313  locFCALShowerToTruthMap[locFCALShower] = pair<const DFCALTruthShower*, double>(locFCALTruthShower, locBestMatchDistance);
314  locFCALTruthToShowerMap[locFCALTruthShower] = pair<const DFCALShower*, double>(locFCALShower, locBestMatchDistance);
315 
316  locFCALShowers.erase(locFCALShowers.begin() + locBestFCALShowerIndex);
317  locFCALTruthShowers.erase(locFCALTruthShowers.begin() + locBestFCALTruthShowerIndex);
318  }
319 
320  locMCThrownMatching->Set_FCALShowerToTruthMap(locFCALShowerToTruthMap);
321  locMCThrownMatching->Set_FCALTruthToShowerMap(locFCALTruthToShowerMap);
322 }
323 
324 void DMCThrownMatching_factory::Find_GenReconMatches_BCALShowers(JEventLoop* locEventLoop, DMCThrownMatching* locMCThrownMatching) const
325 {
326  vector<const DBCALTruthShower*> locBCALTruthShowers;
327  const DBCALTruthShower* locBCALTruthShower;
328  locEventLoop->Get(locBCALTruthShowers);
329 
330  vector<const DBCALShower*> locBCALShowers;
331  const DBCALShower* locBCALShower;
332  locEventLoop->Get(locBCALShowers);
333 
334  map<const DBCALShower*, pair<const DBCALTruthShower*, double> > locBCALShowerToTruthMap;
335  map<const DBCALTruthShower*, pair<const DBCALShower*, double> > locBCALTruthToShowerMap;
336 
337  vector<const DMCThrown*> locMCThrowns;
338  locEventLoop->Get(locMCThrowns);
339 
340  //map of truth shower to dmcthrown
341  map<const DBCALTruthShower*, const DMCThrown*> locBCALTruthToMCThrownMap;
342  for(size_t loc_i = 0; loc_i < locBCALTruthShowers.size(); ++loc_i)
343  {
344  for(size_t loc_j = 0; loc_j < locMCThrowns.size(); ++loc_j)
345  {
346  if(locBCALTruthShowers[loc_i]->track != locMCThrowns[loc_j]->myid)
347  continue;
348  locBCALTruthToMCThrownMap[locBCALTruthShowers[loc_i]] = locMCThrowns[loc_j];
349  break;
350  }
351  }
352 
353  size_t locBestBCALShowerIndex = 0, locBestBCALTruthShowerIndex = 0;
354  while((!locBCALShowers.empty()) && (!locBCALTruthShowers.empty()))
355  {
356  bool locMatchFoundFlag = false;
357  double locLargestEnergy = 0.0;
358  double locBestUnitCircleArcLength = 0.0;
359  for(size_t loc_i = 0; loc_i < locBCALTruthShowers.size(); ++loc_i)
360  {
361  locBCALTruthShower = locBCALTruthShowers[loc_i];
362  for(size_t loc_j = 0; loc_j < locBCALShowers.size(); ++loc_j)
363  {
364  locBCALShower = locBCALShowers[loc_j];
365 
366  DVector3 locBCALShowerPosition(locBCALShower->x, locBCALShower->y, locBCALShower->z);
367  DVector3 locBCALTruthShowerPosition(locBCALTruthShower->r*cos(locBCALTruthShower->phi), locBCALTruthShower->r*sin(locBCALTruthShower->phi), locBCALTruthShower->z);
368 
369  //compare theta & phi
370  map<const DBCALTruthShower*, const DMCThrown*>::const_iterator locIterator = locBCALTruthToMCThrownMap.find(locBCALTruthShower);
371  DVector3 locProductionVertex(0.0, 0.0, dTargetCenter); //target center
372  if(locIterator != locBCALTruthToMCThrownMap.end())
373  locProductionVertex = locIterator->second->position();
374  double locTrueTheta = (locBCALTruthShowerPosition - locProductionVertex).Theta();
375  double locReconstructedTheta = (locBCALShowerPosition - locProductionVertex).Theta();
376 
377  double locDeltaPhi = locBCALShowerPosition.Phi() - locBCALTruthShower->phi;
378  while(locDeltaPhi > TMath::Pi())
379  locDeltaPhi -= 2.0*TMath::Pi();
380  while(locDeltaPhi < -1.0*TMath::Pi())
381  locDeltaPhi += 2.0*TMath::Pi();
382 
383  double locUnitCircleArcLength = acos(sin(locTrueTheta)*sin(locReconstructedTheta)*cos(locDeltaPhi) + cos(locTrueTheta)*cos(locReconstructedTheta));
384 //x = r*cos(phi)*sin(theta)
385 //y = r*sin(phi)*sin(theta)
386 //z = r*cos(theta)
387 
388 //unit arclength = acos(unit dot product)
389 //unit dot product
390  //cos(phi1)*sin(theta1)*cos(phi2)*sin(theta2) + sin(phi1)*sin(theta1)*sin(phi2)*sin(theta2) + cos(theta1)*cos(theta2)
391  //sin(theta1)*sin(theta2)*(cos(phi1)*cos(phi2) + sin(phi1)*sin(phi2)) + cos(theta1)*cos(theta2)
392 //unit arclength = acos(sin(theta1)*sin(theta2)*cos(phi1 - phi2) + cos(theta1)*cos(theta2));
393 
394  if((locUnitCircleArcLength*180.0/TMath::Pi()) > dMaximumBCALMatchAngleDegrees)
395  continue;
396 
397  //keep the shower with the largest energy
398  if(locBCALShower->E > locLargestEnergy)
399  {
400  locBestUnitCircleArcLength = locUnitCircleArcLength;
401  locMatchFoundFlag = true;
402  locLargestEnergy = locBCALShower->E;
403  locBestBCALTruthShowerIndex = loc_i;
404  locBestBCALShowerIndex = loc_j;
405  }
406  }
407  }
408 
409  if(!locMatchFoundFlag) //no more good matches!
410  break;
411 
412  locBCALTruthShower = locBCALTruthShowers[locBestBCALTruthShowerIndex];
413  locBCALShower = locBCALShowers[locBestBCALShowerIndex];
414 
415  locBCALShowerToTruthMap[locBCALShower] = pair<const DBCALTruthShower*, double>(locBCALTruthShower, locBestUnitCircleArcLength);
416  locBCALTruthToShowerMap[locBCALTruthShower] = pair<const DBCALShower*, double>(locBCALShower, locBestUnitCircleArcLength);
417 
418  locBCALShowers.erase(locBCALShowers.begin() + locBestBCALShowerIndex);
419  locBCALTruthShowers.erase(locBCALTruthShowers.begin() + locBestBCALTruthShowerIndex);
420  }
421 
422  locMCThrownMatching->Set_BCALShowerToTruthMap(locBCALShowerToTruthMap);
423  locMCThrownMatching->Set_BCALTruthToShowerMap(locBCALTruthToShowerMap);
424 }
425 
426 void DMCThrownMatching_factory::Find_GenReconMatches_TOFPoints(JEventLoop* locEventLoop, DMCThrownMatching* locMCThrownMatching) const
427 {
428  vector<const DTOFTruth*> locTOFTruths;
429  const DTOFTruth* locTOFTruth;
430  locEventLoop->Get(locTOFTruths);
431 
432  vector<const DTOFPoint*> locTOFPoints;
433  const DTOFPoint* locTOFPoint;
434  locEventLoop->Get(locTOFPoints);
435 
436  map<const DTOFPoint*, pair<const DTOFTruth*, double> > locTOFPointToTruthMap;
437  map<const DTOFTruth*, pair<const DTOFPoint*, double> > locTOFTruthToPointMap;
438 
439  size_t locBestTOFPointIndex = 0, locBestTOFTruthIndex = 0;
440  double locBestMatchDistance = 0.0, locMatchDistance;
441  while((!locTOFPoints.empty()) && (!locTOFTruths.empty()))
442  {
443  bool locMatchFoundFlag = false;
444  double locLargestEnergy = 0.0;
445  for(size_t loc_i = 0; loc_i < locTOFTruths.size(); ++loc_i)
446  {
447  locTOFTruth = locTOFTruths[loc_i];
448  for(size_t loc_j = 0; loc_j < locTOFPoints.size(); ++loc_j)
449  {
450  locTOFPoint = locTOFPoints[loc_j];
451  DVector3 locTOFPointPosition = locTOFPoint->pos;
452 
453  //DTOFPoint and DTOFTruth reported at different z's (I think center vs. detector face): propagate truth information to the reconstructed z
454  double locDeltaZ = locTOFPointPosition.Z() - locTOFTruth->z;
455  DVector3 locTrueMomentum(locTOFTruth->px, locTOFTruth->py, locTOFTruth->pz);
456  double locDeltaPathLength = locDeltaZ/cos(locTrueMomentum.Theta());
457  double locTrueX = locTOFTruth->x + locDeltaPathLength*sin(locTrueMomentum.Theta())*cos(locTrueMomentum.Phi());
458  double locTrueY = locTOFTruth->y + locDeltaPathLength*sin(locTrueMomentum.Theta())*sin(locTrueMomentum.Phi());
459  DVector3 locTOFTruthPosition(locTrueX, locTrueY, locTOFPointPosition.Z());
460 
461  locMatchDistance = (locTOFTruthPosition - locTOFPointPosition).Mag();
462  if(locMatchDistance > dMaximumTOFMatchDistance)
463  continue;
464 
465  //keep the hit with the largest energy
466  if(locTOFTruth->E > locLargestEnergy)
467  {
468  locBestMatchDistance = locMatchDistance;
469  locMatchFoundFlag = true;
470  locLargestEnergy = locTOFTruth->E;
471  locBestTOFTruthIndex = loc_i;
472  locBestTOFPointIndex = loc_j;
473  }
474  }
475  }
476 
477  if(!locMatchFoundFlag) //no more good matches!
478  break;
479 
480  locTOFTruth = locTOFTruths[locBestTOFTruthIndex];
481  locTOFPoint = locTOFPoints[locBestTOFPointIndex];
482 
483  locTOFPointToTruthMap[locTOFPoint] = pair<const DTOFTruth*, double>(locTOFTruth, locBestMatchDistance);
484  locTOFTruthToPointMap[locTOFTruth] = pair<const DTOFPoint*, double>(locTOFPoint, locBestMatchDistance);
485 
486  locTOFPoints.erase(locTOFPoints.begin() + locBestTOFPointIndex);
487  locTOFTruths.erase(locTOFTruths.begin() + locBestTOFTruthIndex);
488  }
489 
490  locMCThrownMatching->Set_TOFPointToTruthMap(locTOFPointToTruthMap);
491  locMCThrownMatching->Set_TOFTruthToPointMap(locTOFTruthToPointMap);
492 }
493 
494 void DMCThrownMatching_factory::Find_GenReconMatches_ChargedTrack(const vector<const DChargedTrack*>& locChargedTracks, DMCThrownMatching* locMCThrownMatching) const
495 {
496  //assumes Find_GenReconMatches_ChargedHypo has been called first!
497  map<const DChargedTrackHypothesis*, pair<const DMCThrown*, double> > locChargedHypoToThrownMap;
498  map<const DMCThrown*, pair<deque<const DChargedTrackHypothesis*>, double> > locThrownToChargedHypoMap;
499  locMCThrownMatching->Get_ChargedHypoToThrownMap(locChargedHypoToThrownMap);
500  locMCThrownMatching->Get_ThrownToChargedHypoMap(locThrownToChargedHypoMap);
501 
502  map<const DChargedTrack*, pair<const DMCThrown*, double> > locChargedToThrownMap;
503  map<const DMCThrown*, pair<const DChargedTrack*, double> > locThrownToChargedMap;
504 
505  map<const DMCThrown*, pair<deque<const DChargedTrackHypothesis*>, double> >::iterator locIterator;
506  for(locIterator = locThrownToChargedHypoMap.begin(); locIterator != locThrownToChargedHypoMap.end(); ++locIterator)
507  {
508  for(size_t loc_j = 0; loc_j < locChargedTracks.size(); ++loc_j)
509  {
510  if(locChargedTracks[loc_j]->candidateid == (locIterator->second.first)[0]->Get_TrackTimeBased()->candidateid) //match
511  {
512  locChargedToThrownMap[locChargedTracks[loc_j]] = pair<const DMCThrown*, double>(locIterator->first, locIterator->second.second);
513  locThrownToChargedMap[locIterator->first] = pair<const DChargedTrack*, double>(locChargedTracks[loc_j], locIterator->second.second);
514  }
515  }
516  }
517 
518  locMCThrownMatching->Set_ChargedToThrownMap(locChargedToThrownMap);
519  locMCThrownMatching->Set_ThrownToChargedMap(locThrownToChargedMap);
520 }
521 
522 
523 void DMCThrownMatching_factory::Find_GenReconMatches_ChargedHypo(const vector<const DMCThrown*>& locMCThrownVector, const vector<const DChargedTrackHypothesis*>& locChargedTrackHypothesisVector, DMCThrownMatching* locMCThrownMatching) const
524 {
525  map<const DChargedTrackHypothesis*, pair<const DMCThrown*, double> > locChargedToThrownMap;
526  map<const DMCThrown*, pair<deque<const DChargedTrackHypothesis*>, double> > locThrownToChargedMap;
527 
528  map<int, const DMCThrown*> locMyIDToThrownMap;
529  for(size_t loc_i = 0; loc_i < locMCThrownVector.size(); ++loc_i)
530  locMyIDToThrownMap[locMCThrownVector[loc_i]->myid] = locMCThrownVector[loc_i];
531 
532  if(dDebugLevel > 0)
533  cout << "START IT!" << endl;
534 
535  //calculate weighted num-matched hits
536  set<pair<double, pair<const DMCThrown*, const DChargedTrackHypothesis*> > > locParticleMatches;
537  for(size_t loc_j = 0; loc_j < locChargedTrackHypothesisVector.size(); ++loc_j)
538  {
539  const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTrackHypothesisVector[loc_j];
540  auto locTrackTimeBased = locChargedTrackHypothesis->Get_TrackTimeBased();
541 
542  if(locMyIDToThrownMap.find(locTrackTimeBased->dMCThrownMatchMyID) == locMyIDToThrownMap.end())
543  continue;
544 
545  const DMCThrown* locMCThrown = locMyIDToThrownMap[locTrackTimeBased->dMCThrownMatchMyID];
546  double locNumTrackHits = Get_NumTrackHits(locTrackTimeBased);
547  double locHitFraction = locTrackTimeBased->dNumHitsMatchedToThrown/locNumTrackHits;
548  if(locHitFraction < dMinTrackMatchHitFraction)
549  continue; //not good enough!
550 
551  //#-matched-hits * hit_fraction
552  double locMatchFOM = locTrackTimeBased->dNumHitsMatchedToThrown * locHitFraction;
553 
554  if(dDebugLevel > 0)
555  {
556  cout << "MATCHING: MCTHROWN: ";
557  cout << ParticleType((Particle_t)(locMCThrown->type)) << ", " << locMCThrown->momentum().Mag() << ", " << locMCThrown->momentum().Theta()*180.0/TMath::Pi() << ", " << locMCThrown->momentum().Phi()*180.0/TMath::Pi() << endl;
558  cout << "MATCHING: CHARGEDHYPO: ";
559  cout << ParticleType(locChargedTrackHypothesis->PID()) << ", " << locChargedTrackHypothesis->momentum().Mag() << ", " << locChargedTrackHypothesis->momentum().Theta()*180.0/TMath::Pi() << ", " << locChargedTrackHypothesis->momentum().Phi()*180.0/TMath::Pi() << endl;
560  cout << "MATCHING: FOM, candidate id: " << locMatchFOM << ", " << locChargedTrackHypothesis->Get_TrackTimeBased()->candidateid << endl;
561  }
562 
563  pair<const DMCThrown*, const DChargedTrackHypothesis*> locTrackPair(locMCThrown, locChargedTrackHypothesis);
564  pair<double, pair<const DMCThrown*, const DChargedTrackHypothesis*> > locMatchPair(locMatchFOM, locTrackPair);
565  locParticleMatches.insert(locMatchPair);
566  }
567  if(locParticleMatches.empty())
568  return; //nothing to set
569 
570  //loop over sets, save the best matches //sorted from least to greatest
571  set<pair<double, pair<const DMCThrown*, const DChargedTrackHypothesis*> > >::iterator locIterator = locParticleMatches.end();
572  set<const DMCThrown*> locMatchedThrowns;
573  set<const DChargedTrackHypothesis*> locMatchedHypotheses;
574  for(--locIterator; locIterator != locParticleMatches.begin(); --locIterator)
575  {
576  double locMatchFOM = locIterator->first;
577  const DMCThrown* locMCThrown = locIterator->second.first;
578  const DChargedTrackHypothesis* locChargedTrackHypothesis = locIterator->second.second;
579 
580  if(locMatchedThrowns.find(locMCThrown) != locMatchedThrowns.end())
581  continue; //track match already saved
582  if(locMatchedHypotheses.find(locChargedTrackHypothesis) != locMatchedHypotheses.end())
583  continue; //track match already saved
584 
585  locMatchedThrowns.insert(locMCThrown);
586 
587  //automatically add all other DChargedTrackHypothesis objects from the same DChargedTrack to this match.
588  deque<const DChargedTrackHypothesis*> locMatchedChargedHypos;
589  for(size_t loc_i = 0; loc_i < locChargedTrackHypothesisVector.size(); ++loc_i)
590  {
591  if(dDebugLevel > 0)
592  {
593  cout << "CHECKING: CHARGEDHYPO: ";
594  cout << ParticleType(locChargedTrackHypothesisVector[loc_i]->PID()) << ", " << locChargedTrackHypothesisVector[loc_i]->momentum().Mag() << ", " << locChargedTrackHypothesisVector[loc_i]->momentum().Theta()*180.0/TMath::Pi() << ", " << locChargedTrackHypothesisVector[loc_i]->momentum().Phi()*180.0/TMath::Pi() << endl;
595  cout << "best id, test id = " << locChargedTrackHypothesis->Get_TrackTimeBased()->candidateid << ", " << locChargedTrackHypothesisVector[loc_i]->Get_TrackTimeBased()->candidateid << endl;
596  }
597  if(locChargedTrackHypothesisVector[loc_i]->Get_TrackTimeBased()->candidateid == locChargedTrackHypothesis->Get_TrackTimeBased()->candidateid)
598  {
599  if(dDebugLevel > 0)
600  cout << "save!" << endl;
601  locChargedToThrownMap[locChargedTrackHypothesisVector[loc_i]] = pair<const DMCThrown*, double>(locMCThrown, locMatchFOM);
602  locMatchedChargedHypos.push_back(locChargedTrackHypothesisVector[loc_i]);
603  locMatchedHypotheses.insert(locChargedTrackHypothesisVector[loc_i]);
604  }
605  }
606  locThrownToChargedMap[locMCThrown] = pair<deque<const DChargedTrackHypothesis*>, double>(locMatchedChargedHypos, locMatchFOM);
607  }
608 
609  locMCThrownMatching->Set_ChargedHypoToThrownMap(locChargedToThrownMap);
610  locMCThrownMatching->Set_ThrownToChargedHypoMap(locThrownToChargedMap);
611 }
612 
613 void DMCThrownMatching_factory::Find_GenReconMatches_NeutralParticle(const vector<const DNeutralParticle*>& locNeutralParticles, DMCThrownMatching* locMCThrownMatching) const
614 {
615  //assumes Find_GenReconMatches_NeutralHypo has been called first!
616  map<const DNeutralParticleHypothesis*, pair<const DMCThrown*, double> > locNeutralHypoToThrownMap;
617  map<const DMCThrown*, pair<deque<const DNeutralParticleHypothesis*>, double> > locThrownToNeutralHypoMap;
618  locMCThrownMatching->Get_NeutralHypoToThrownMap(locNeutralHypoToThrownMap);
619  locMCThrownMatching->Get_ThrownToNeutralHypoMap(locThrownToNeutralHypoMap);
620 
621  map<const DNeutralParticle*, pair<const DMCThrown*, double> > locNeutralToThrownMap;
622  map<const DMCThrown*, pair<const DNeutralParticle*, double> > locThrownToNeutralMap;
623 
624  map<const DMCThrown*, pair<deque<const DNeutralParticleHypothesis*>, double> >::iterator locIterator;
625  for(locIterator = locThrownToNeutralHypoMap.begin(); locIterator != locThrownToNeutralHypoMap.end(); ++locIterator)
626  {
627  for(size_t loc_j = 0; loc_j < locNeutralParticles.size(); ++loc_j)
628  {
629  if((locIterator->second.first)[0]->Get_NeutralShower() == locNeutralParticles[loc_j]->dNeutralShower)
630  {
631  locNeutralToThrownMap[locNeutralParticles[loc_j]] = pair<const DMCThrown*, double>(locIterator->first, locIterator->second.second);
632  locThrownToNeutralMap[locIterator->first] = pair<const DNeutralParticle*, double>(locNeutralParticles[loc_j], locIterator->second.second);
633  }
634  }
635  }
636 
637  locMCThrownMatching->Set_NeutralToThrownMap(locNeutralToThrownMap);
638  locMCThrownMatching->Set_ThrownToNeutralMap(locThrownToNeutralMap);
639 }
640 
641 void DMCThrownMatching_factory::Find_GenReconMatches_NeutralHypo(const vector<const DMCThrown*>& locMCThrownVector, const vector<const DNeutralParticleHypothesis*>& locNeutralParticleHypothesisVector, DMCThrownMatching* locMCThrownMatching) const
642 {
643  map<const DNeutralParticleHypothesis*, pair<const DMCThrown*, double> > locNeutralToThrownMap;
644  map<const DMCThrown*, pair<deque<const DNeutralParticleHypothesis*>, double> > locThrownToNeutralMap;
645 
646  if(dDebugLevel > 0)
647  cout << "START IT!" << endl;
648 
649  //build inverse covariance matrix map
650  map<const DNeutralParticleHypothesis*, TMatrixDSym> locInverseCovMatrixMap;
651  for(size_t loc_i = 0; loc_i < locNeutralParticleHypothesisVector.size(); ++loc_i)
652  {
653  const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticleHypothesisVector[loc_i];
654  TMatrixDSym locInverse3x3Matrix(3);
655  if(Calc_InverseMatrix(*(locNeutralParticleHypothesis->errorMatrix()), locInverse3x3Matrix))
656  locInverseCovMatrixMap.insert(pair<const DNeutralParticleHypothesis*, TMatrixFSym>(locNeutralParticleHypothesis, locInverse3x3Matrix));
657  }
658 
659  //calculate match FOMs
660  set<pair<double, pair<const DMCThrown*, const DNeutralParticleHypothesis*> > > locParticleMatches;
661  for(size_t loc_i = 0; loc_i < locMCThrownVector.size(); ++loc_i)
662  {
663  const DMCThrown* locMCThrown = locMCThrownVector[loc_i];
664  for(size_t loc_j = 0; loc_j < locNeutralParticleHypothesisVector.size(); ++loc_j)
665  {
666  const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticleHypothesisVector[loc_j];
667  if(locInverseCovMatrixMap.find(locNeutralParticleHypothesis) == locInverseCovMatrixMap.end())
668  continue;
669  TMatrixDSym& locInverse3x3Matrix = locInverseCovMatrixMap[locNeutralParticleHypothesis];
670 
671  double locMatchFOM = Calc_MatchFOM(locMCThrown->momentum(), locNeutralParticleHypothesis->momentum(), locInverse3x3Matrix);
672 
673  if(dDebugLevel > 0)
674  {
675  cout << "MATCHING: MCTHROWN: ";
676  cout << ParticleType((Particle_t)(locMCThrown->type)) << ", " << locMCThrown->momentum().Mag() << ", " << locMCThrown->momentum().Theta()*180.0/TMath::Pi() << ", " << locMCThrown->momentum().Phi()*180.0/TMath::Pi() << endl;
677  cout << "MATCHING: NEUTRALHYPO: ";
678  cout << ParticleType(locNeutralParticleHypothesis->PID()) << ", " << locNeutralParticleHypothesis->momentum().Mag() << ", " << locNeutralParticleHypothesis->momentum().Theta()*180.0/TMath::Pi() << ", " << locNeutralParticleHypothesis->momentum().Phi()*180.0/TMath::Pi() << endl;
679  cout << "MATCHING: FOM: " << locMatchFOM << endl;
680  }
681 
682  pair<const DMCThrown*, const DNeutralParticleHypothesis*> locTrackPair(locMCThrown, locNeutralParticleHypothesis);
683  pair<double, pair<const DMCThrown*, const DNeutralParticleHypothesis*> > locMatchPair(locMatchFOM, locTrackPair);
684  locParticleMatches.insert(locMatchPair);
685  }
686  }
687 
688  if(locParticleMatches.empty())
689  return; //nothing to set
690 
691  //loop over sets, save the best matches //sorted from least to greatest
692  set<pair<double, pair<const DMCThrown*, const DNeutralParticleHypothesis*> > >::iterator locIterator = locParticleMatches.end();
693  set<const DMCThrown*> locMatchedThrowns;
694  set<const DNeutralParticleHypothesis*> locMatchedHypotheses;
695  for(--locIterator; locIterator != locParticleMatches.begin(); --locIterator)
696  {
697  double locMatchFOM = locIterator->first;
698  const DMCThrown* locMCThrown = locIterator->second.first;
699  const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locIterator->second.second;
700 
701  if(locMatchedThrowns.find(locMCThrown) != locMatchedThrowns.end())
702  continue; //track match already saved
703  if(locMatchedHypotheses.find(locNeutralParticleHypothesis) != locMatchedHypotheses.end())
704  continue; //track match already saved
705 
706  locMatchedThrowns.insert(locMCThrown);
707 
708  //automatically add all other DNeutralParticleHypothesis objects from the same DNeutralShower to this match.
709  deque<const DNeutralParticleHypothesis*> locMatchedNeutralHypos(1, locNeutralParticleHypothesis);
710  for(size_t loc_i = 0; loc_i < locNeutralParticleHypothesisVector.size(); ++loc_i)
711  {
712  if(locNeutralParticleHypothesis->Get_NeutralShower() == locNeutralParticleHypothesisVector[loc_i]->Get_NeutralShower())
713  {
714  locNeutralToThrownMap[locNeutralParticleHypothesisVector[loc_i]] = pair<const DMCThrown*, double>(locMCThrown, locMatchFOM);
715  locMatchedNeutralHypos.push_back(locNeutralParticleHypothesisVector[loc_i]);
716  locMatchedHypotheses.insert(locNeutralParticleHypothesisVector[loc_i]);
717  }
718  }
719  locThrownToNeutralMap[locMCThrown] = pair<deque<const DNeutralParticleHypothesis*>, double>(locMatchedNeutralHypos, locMatchFOM);
720  }
721 
722  locMCThrownMatching->Set_NeutralHypoToThrownMap(locNeutralToThrownMap);
723  locMCThrownMatching->Set_ThrownToNeutralHypoMap(locThrownToNeutralMap);
724 }
725 
726 bool DMCThrownMatching_factory::Calc_InverseMatrix(const TMatrixFSym& locInputCovarianceMatrix, TMatrixDSym& locInverse3x3Matrix) const
727 {
728  double locTotalError = sqrt(locInputCovarianceMatrix(0, 0) + locInputCovarianceMatrix(1, 1) + locInputCovarianceMatrix(2, 2));
729  if(locTotalError >= dMaxTotalParticleErrorForMatch)
730  return false;
731 
732  locInverse3x3Matrix.ResizeTo(3, 3);
733  for(unsigned int loc_i = 0; loc_i < 3; ++loc_i)
734  {
735  for(unsigned int loc_j = 0; loc_j < 3; ++loc_j)
736  locInverse3x3Matrix(loc_i, loc_j) = locInputCovarianceMatrix(loc_i, loc_j);
737  }
738 
739  //invert matrix
740  TDecompLU locDecompLU(locInverse3x3Matrix);
741  //check to make sure that the matrix is decomposable and has a non-zero determinant
742  if((!locDecompLU.Decompose()) || (fabs(locInverse3x3Matrix.Determinant()) < 1.0E-300))
743  return false; // matrix is not invertible
744  locInverse3x3Matrix.Invert();
745  return true;
746 }
747 
748 double DMCThrownMatching_factory::Calc_MatchFOM(const DVector3& locMomentum_Thrown, const DVector3& locMomentum_Detected, TMatrixDSym locInverse3x3Matrix) const
749 {
750  DVector3 locDeltaP3 = locMomentum_Detected - locMomentum_Thrown;
751  TMatrix locDeltas(3, 1);
752  locDeltas(0, 0) = locDeltaP3.Px();
753  locDeltas(1, 0) = locDeltaP3.Py();
754  locDeltas(2, 0) = locDeltaP3.Pz();
755 
756  double locChiSq = (locInverse3x3Matrix.SimilarityT(locDeltas))(0, 0);
757  double locFOM = TMath::Prob(locChiSq, 3);
758  if(dDebugLevel > 10)
759  cout << "delta pxyz, chisq, FOM = " << locDeltaP3.Px() << ", " << locDeltaP3.Py() << ", " << locDeltaP3.Pz() << ", " << locChiSq << ", " << locFOM << endl;
760 
761  return locFOM;
762 }
763 
764 //------------------
765 // erun
766 //------------------
768 {
769  return NOERROR;
770 }
771 
772 //------------------
773 // fini
774 //------------------
776 {
777  return NOERROR;
778 }
779 
780 
float py() const
void Set_ThrownToChargedMap(const map< const DMCThrown *, pair< const DChargedTrack *, double > > &locThrownToChargedMap)
Definition: track.h:16
void Set_ThrownToNeutralMap(const map< const DMCThrown *, pair< const DNeutralParticle *, double > > &locThrownToNeutralMap)
jerror_t brun(jana::JEventLoop *locEventLoop, int32_t runnumber)
Called everytime a new run number is detected.
double getEnergy() const
Definition: DFCALShower.h:156
jerror_t init(void)
Called once at program start.
void Find_GenReconMatches_NeutralHypo(const vector< const DMCThrown * > &locInputMCThrownVector, const vector< const DNeutralParticleHypothesis * > &locInputNeutralParticleHypothesisVector, DMCThrownMatching *locMCThrownMatching) const
bool Calc_InverseMatrix(const TMatrixFSym &locInputCovarianceMatrix, TMatrixDSym &locInverse3x3Matrix) const
void Set_FCALTruthToShowerMap(map< const DFCALTruthShower *, pair< const DFCALShower *, double > > &locFCALTruthToShowerMap)
void Set_FCALShowerToTruthMap(map< const DFCALShower *, pair< const DFCALTruthShower *, double > > &locFCALShowerToTruthMap)
float z
Definition: DTOFTruth.h:23
float x
Definition: DTOFTruth.h:23
void Set_NeutralHypoToThrownMap(const map< const DNeutralParticleHypothesis *, pair< const DMCThrown *, double > > &locNeutralHypoToThrownMap)
TVector3 DVector3
Definition: DVector3.h:14
float y() const
double energy(void) const
void Find_GenReconMatches_TOFPoints(JEventLoop *locEventLoop, DMCThrownMatching *locMCThrownMatching) const
int row
Definition: DTAGMHit.h:20
float E
Definition: DTOFTruth.h:26
void Set_BCALTruthToShowerMap(map< const DBCALTruthShower *, pair< const DBCALShower *, double > > &locBCALTruthToShowerMap)
const DTrackTimeBased * Get_TrackTimeBased(void) const
void Find_GenReconMatches_NeutralParticle(const vector< const DNeutralParticle * > &locNeutralParticles, DMCThrownMatching *locMCThrownMatching) const
static char * ParticleType(Particle_t p)
Definition: particleType.h:142
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
void Set_TOFTruthToPointMap(const map< const DTOFTruth *, pair< const DTOFPoint *, double > > &locTOFTruthToPointMap)
static int ParticleCharge(Particle_t p)
oid_t candidateid
id of DTrackCandidate corresponding to this track
const DNeutralShower * Get_NeutralShower(void) const
void Set_TaggedMCGENBeamPhoton(const DBeamPhoton *locBeamPhoton)
float y
Definition: DTOFTruth.h:23
void Set_ThrownToChargedHypoMap(const map< const DMCThrown *, pair< deque< const DChargedTrackHypothesis * >, double > > &locThrownToChargedHypoMap)
void Set_ThrownToNeutralHypoMap(const map< const DMCThrown *, pair< deque< const DNeutralParticleHypothesis * >, double > > &locThrownToNeutralHypoMap)
size_t Get_NumTrackHits(const DTrackTimeBased *locTrackTimeBased)
int counter_id
Definition: DTAGHHit.h:20
jerror_t evnt(jana::JEventLoop *locEventLoop, uint64_t eventnumber)
Called every event.
const DChargedTrack * Get_MatchingChargedTrack(const DMCThrown *locInputMCThrown, double &locMatchFOM) const
void Find_GenReconMatches_BeamPhotons(JEventLoop *locEventLoop, DMCThrownMatching *locMCThrownMatching) const
DVector3 pos
Definition: DTOFPoint.h:33
void Get_ChargedHypoToThrownMap(map< const DChargedTrackHypothesis *, pair< const DMCThrown *, double > > &locChargedHypoToThrownMap) const
void Find_GenReconMatches_ChargedTrack(const vector< const DChargedTrack * > &locChargedTracks, DMCThrownMatching *locMCThrownMatching) const
void Set_BeamTruthToPhotonMap(map< const DBeamPhoton *, const DBeamPhoton * > &locBeamTruthToPhotonMap)
double Calc_MatchFOM(const DVector3 &locMomentum_Thrown, const DVector3 &locMomentum_Detected, TMatrixDSym locInverse3x3Matrix) const
void Set_TOFPointToTruthMap(const map< const DTOFPoint *, pair< const DTOFTruth *, double > > &locTOFPointToTruthMap)
DGeometry * GetDGeometry(unsigned int run_number)
float px
Definition: DTOFTruth.h:24
float px() const
jerror_t fini(void)
Called after last event of last event source has been processed.
void Set_ChargedToThrownMap(const map< const DChargedTrack *, pair< const DMCThrown *, double > > &locChargedToThrownMap)
float py
Definition: DTOFTruth.h:24
void Find_GenReconMatches_ChargedHypo(const vector< const DMCThrown * > &locInputMCThrownVector, const vector< const DChargedTrackHypothesis * > &locInputChargedTrackHypothesisVector, DMCThrownMatching *locMCThrownMatching) const
void Set_NeutralToThrownMap(const map< const DNeutralParticle *, pair< const DMCThrown *, double > > &locNeutralToThrownMap)
void Find_GenReconMatches_FCALShowers(JEventLoop *locEventLoop, DMCThrownMatching *locMCThrownMatching) const
double sqrt(double)
double sin(double)
void Get_NeutralHypoToThrownMap(map< const DNeutralParticleHypothesis *, pair< const DMCThrown *, double > > &locNeutralHypoToThrownMap) const
float pz
Definition: DTOFTruth.h:24
const DVector3 & momentum(void) const
void Set_MCGENBeamPhoton(const DBeamPhoton *locBeamPhoton)
void Set_ChargedHypoToThrownMap(const map< const DChargedTrackHypothesis *, pair< const DMCThrown *, double > > &locChargedHypoToThrownMap)
void Set_BeamPhotonToTruthMap(map< const DBeamPhoton *, const DBeamPhoton * > &locBeamPhotonToTruthMap)
void Find_GenReconMatches_BCALShowers(JEventLoop *locEventLoop, DMCThrownMatching *locMCThrownMatching) const
void Get_ThrownToChargedHypoMap(map< const DMCThrown *, pair< deque< const DChargedTrackHypothesis * >, double > > &locThrownToChargedHypoMap) const
float x() const
int column
Definition: DTAGMHit.h:21
void Set_BCALShowerToTruthMap(map< const DBCALShower *, pair< const DBCALTruthShower *, double > > &locBCALShowerToTruthMap)
shared_ptr< const TMatrixFSym > errorMatrix(void) const
float pz() const
float z() const
int type
GEANT particle ID.
Definition: DMCThrown.h:20
DVector3 getPosition() const
Definition: DFCALShower.h:151
bool GetTargetZ(double &z_target) const
z-location of center of target
Definition: DGeometry.cc:1933
Particle_t PID(void) const
const DMCThrown * Get_MatchingMCThrown(const DChargedTrackHypothesis *locChargedTrackHypothesis, double &locMatchFOM) const
Particle_t
Definition: particleType.h:12
void Get_ThrownToNeutralHypoMap(map< const DMCThrown *, pair< deque< const DNeutralParticleHypothesis * >, double > > &locThrownToNeutralHypoMap) const