LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
NViewDeltaRayMatchingAlgorithm.cc
Go to the documentation of this file.
1 
9 #include "Pandora/AlgorithmHeaders.h"
10 
15 
19 
21 
24 
26 
27 using namespace pandora;
28 
29 namespace lar_content
30 {
31 
32 template <typename T>
34  m_pseudoChi2Cut(1.5f),
35  m_xOverlapWindow(1.f),
36  m_minMatchedFraction(0.5),
37  m_minMatchedPoints(3),
38  m_minProjectedPositions(3),
39  m_maxCosmicRayHitFraction(0.05f),
40  m_maxDistanceToCluster(0.5f),
41  m_maxDistanceToReferencePoint(5.f),
42  m_strayClusterSeparation(2.f)
43 {
44 }
45 
46 //------------------------------------------------------------------------------------------------------------------------------------------
47 
48 template <typename T>
49 void NViewDeltaRayMatchingAlgorithm<T>::SelectInputClusters(const ClusterList *const pInputClusterList, ClusterList &selectedClusterList) const
50 {
51  for (const Cluster *const pCluster : *pInputClusterList)
52  {
53  if ((pCluster->IsAvailable()) && (this->DoesClusterPassTensorThreshold(pCluster)))
54  selectedClusterList.push_back(pCluster);
55  }
56 }
57 
58 //------------------------------------------------------------------------------------------------------------------------------------------
59 
60 template <typename T>
61 void NViewDeltaRayMatchingAlgorithm<T>::PrepareInputClusters(ClusterList &preparedClusterList)
62 {
63  if (preparedClusterList.empty())
64  return;
65 
66  const PfoList *pMuonPfoList(nullptr);
67 
68  PANDORA_THROW_RESULT_IF_AND_IF(
69  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_INITIALIZED, !=, PandoraContentApi::GetList(*this, m_muonPfoListName, pMuonPfoList));
70 
71  if ((!pMuonPfoList) || pMuonPfoList->empty())
72  return;
73 
74  const HitType hitType(LArClusterHelper::GetClusterHitType(preparedClusterList.front()));
75 
77 
78  this->FillStrayClusterList(LArClusterHelper::GetClusterHitType(preparedClusterList.front()));
79 }
80 
81 //------------------------------------------------------------------------------------------------------------------------------------------
82 
83 template <typename T>
85 {
86  const ClusterList &inputClusterList(this->GetInputClusterList(hitType));
87  ClusterList &strayClusterList((hitType == TPC_VIEW_U) ? m_strayClusterListU
88  : (hitType == TPC_VIEW_V) ? m_strayClusterListV
90 
91  for (const Cluster *const pCluster : inputClusterList)
92  {
93  if ((pCluster->IsAvailable()) && (!this->DoesClusterPassTensorThreshold(pCluster)))
94  strayClusterList.push_back(pCluster);
95  }
96 }
97 
98 //------------------------------------------------------------------------------------------------------------------------------------------
99 
100 template <typename T>
101 StatusCode NViewDeltaRayMatchingAlgorithm<T>::GetMuonCluster(const PfoList &commonMuonPfoList, const HitType hitType, const Cluster *&pMuonCluster) const
102 {
103  if (commonMuonPfoList.size() != 1)
104  return STATUS_CODE_NOT_FOUND;
105 
106  ClusterList muonClusterList;
107  LArPfoHelper::GetClusters(commonMuonPfoList.front(), hitType, muonClusterList);
108 
109  if (muonClusterList.size() != 1)
110  return STATUS_CODE_NOT_FOUND;
111 
112  pMuonCluster = muonClusterList.front();
113 
114  return STATUS_CODE_SUCCESS;
115 }
116 
117 //------------------------------------------------------------------------------------------------------------------------------------------
118 
119 template <typename T>
120 void NViewDeltaRayMatchingAlgorithm<T>::GetNearbyMuonPfos(const Cluster *const pCluster, ClusterList &consideredClusters, PfoList &nearbyMuonPfos) const
121 {
122  const HitType hitType(LArClusterHelper::GetClusterHitType(pCluster));
125 
126  consideredClusters.push_back(pCluster);
127 
128  const DeltaRayMatchingContainers::ClusterProximityMap::const_iterator clusterProximityIter(clusterProximityMap.find(pCluster));
129 
130  if (clusterProximityIter == clusterProximityMap.end())
131  return;
132 
133  for (const Cluster *const pNearbyCluster : clusterProximityIter->second)
134  {
135  if (std::find(consideredClusters.begin(), consideredClusters.end(), pNearbyCluster) != consideredClusters.end())
136  continue;
137 
138  const DeltaRayMatchingContainers::ClusterToPfoMap::const_iterator pfoIter(clusterToPfoMap.find(pNearbyCluster));
139 
140  if (pfoIter != clusterToPfoMap.end())
141  {
142  if (std::find(nearbyMuonPfos.begin(), nearbyMuonPfos.end(), pfoIter->second) == nearbyMuonPfos.end())
143  nearbyMuonPfos.push_back(pfoIter->second);
144 
145  continue;
146  }
147 
148  this->GetNearbyMuonPfos(pNearbyCluster, consideredClusters, nearbyMuonPfos);
149  }
150 }
151 
152 //------------------------------------------------------------------------------------------------------------------------------------------
153 
154 template <typename T>
156  const Cluster *const pCluster1, const Cluster *const pCluster2, const Cluster *const pCluster3, float &reducedChiSquared) const
157 {
158  float chiSquaredSum(0.f);
159  unsigned int nSamplingPoints(0), nMatchedSamplingPoints(0);
160  XOverlap xThreeViewOverlapObject(0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f);
161 
162  if (this->PerformThreeViewMatching(pCluster1, pCluster2, pCluster3, chiSquaredSum, nSamplingPoints, nMatchedSamplingPoints,
163  xThreeViewOverlapObject) == STATUS_CODE_NOT_FOUND)
164  return STATUS_CODE_NOT_FOUND;
165 
166  if (nSamplingPoints == 0)
167  return STATUS_CODE_NOT_FOUND;
168 
169  reducedChiSquared = chiSquaredSum / nSamplingPoints;
170 
171  return STATUS_CODE_SUCCESS;
172 }
173 
174 //------------------------------------------------------------------------------------------------------------------------------------------
175 
176 template <typename T>
177 StatusCode NViewDeltaRayMatchingAlgorithm<T>::PerformThreeViewMatching(const Cluster *const pClusterU, const Cluster *const pClusterV,
178  const Cluster *const pClusterW, float &chiSquaredSum, unsigned int &nSamplingPoints, unsigned int &nMatchedSamplingPoints, XOverlap &xOverlapObject) const
179 {
180  float xMinU(-std::numeric_limits<float>::max()), xMaxU(+std::numeric_limits<float>::max());
181  float xMinV(-std::numeric_limits<float>::max()), xMaxV(+std::numeric_limits<float>::max());
182  float xMinW(-std::numeric_limits<float>::max()), xMaxW(+std::numeric_limits<float>::max());
183 
184  pClusterU->GetClusterSpanX(xMinU, xMaxU);
185  pClusterV->GetClusterSpanX(xMinV, xMaxV);
186  pClusterW->GetClusterSpanX(xMinW, xMaxW);
187 
188  // Need to remove the xPitch from calculations to be consistent with view xSpan calculated in the xOverlapObject
189  const float xMinCentre(std::max(xMinU, std::max(xMinV, xMinW)));
190  const float xMaxCentre(std::min(xMaxU, std::min(xMaxV, xMaxW)));
191  const float xCentreOverlap(xMaxCentre - xMinCentre);
192 
193  if (xCentreOverlap < std::numeric_limits<float>::epsilon())
194  return STATUS_CODE_NOT_FOUND;
195 
196  const float xPitch(0.5 * m_xOverlapWindow);
197  const float xMin(std::max(xMinU, std::max(xMinV, xMinW)) - xPitch);
198  const float xMax(std::min(xMaxU, std::min(xMaxV, xMaxW)) + xPitch);
199  const float xOverlap(xMax - xMin);
200 
201  const HitType hitTypeU(LArClusterHelper::GetClusterHitType(pClusterU));
202  const HitType hitTypeV(LArClusterHelper::GetClusterHitType(pClusterV));
203  const HitType hitTypeW(LArClusterHelper::GetClusterHitType(pClusterW));
204  const float pitchU{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeU)};
205  const float pitchV{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeV)};
206  const float pitchW{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeW)};
207  const float pitchMax{std::max({pitchU, pitchV, pitchW})};
208 
209  if (hitTypeU == hitTypeV || hitTypeU == hitTypeW || hitTypeV == hitTypeW)
210  return STATUS_CODE_FAILURE;
211 
212  const unsigned int nPoints(1 + static_cast<unsigned int>(xOverlap / xPitch));
213 
214  chiSquaredSum = 0.f;
215  nSamplingPoints = 0;
216  nMatchedSamplingPoints = 0;
217 
218  for (unsigned int n = 0; n < nPoints; ++n)
219  {
220  const float x(xMin + (xMax - xMin) * (static_cast<float>(n) + 0.5f) / static_cast<float>(nPoints));
221  const float xmin(x - xPitch);
222  const float xmax(x + xPitch);
223 
224  try
225  {
226  float zMinU(0.f), zMinV(0.f), zMinW(0.f), zMaxU(0.f), zMaxV(0.f), zMaxW(0.f);
227  pClusterU->GetClusterSpanZ(xmin, xmax, zMinU, zMaxU);
228  pClusterV->GetClusterSpanZ(xmin, xmax, zMinV, zMaxV);
229  pClusterW->GetClusterSpanZ(xmin, xmax, zMinW, zMaxW);
230 
231  const float zU(0.5f * (zMinU + zMaxU));
232  const float zV(0.5f * (zMinV + zMaxV));
233  const float zW(0.5f * (zMinW + zMaxW));
234 
235  const float dzU(zMaxU - zMinU);
236  const float dzV(zMaxV - zMinV);
237  const float dzW(zMaxW - zMinW);
238  const float dzPitch(pitchMax);
239 
240  const float zprojU(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeV, hitTypeW, zV, zW));
241  const float zprojV(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeW, hitTypeU, zW, zU));
242  const float zprojW(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeU, hitTypeV, zU, zV));
243 
244  ++nSamplingPoints;
245 
246  const float deltaSquared(((zU - zprojU) * (zU - zprojU) + (zV - zprojV) * (zV - zprojV) + (zW - zprojW) * (zW - zprojW)) / 3.f);
247  const float sigmaSquared(dzU * dzU + dzV * dzV + dzW * dzW + dzPitch * dzPitch);
248  const float pseudoChi2(deltaSquared / sigmaSquared);
249 
250  chiSquaredSum += pseudoChi2;
251 
252  if (pseudoChi2 < m_pseudoChi2Cut)
253  ++nMatchedSamplingPoints;
254  }
255  catch (StatusCodeException &statusCodeException)
256  {
257  if (statusCodeException.GetStatusCode() != STATUS_CODE_NOT_FOUND)
258  return statusCodeException.GetStatusCode();
259 
260  continue;
261  }
262  }
263 
264  // Apply tensor threshold cuts
265  if (nSamplingPoints == 0)
266  return STATUS_CODE_NOT_FOUND;
267 
268  const float matchedFraction(static_cast<float>(nMatchedSamplingPoints) / static_cast<float>(nSamplingPoints));
269 
270  if ((matchedFraction < m_minMatchedFraction) || (nMatchedSamplingPoints < m_minMatchedPoints))
271  return STATUS_CODE_NOT_FOUND;
272 
273  xOverlapObject = XOverlap(xMinU, xMaxU, xMinV, xMaxV, xMinW, xMaxW, xCentreOverlap);
274 
275  return STATUS_CODE_SUCCESS;
276 }
277 
278 //------------------------------------------------------------------------------------------------------------------------------------------
279 
280 template <typename T>
282  const CaloHitList &pCluster1, const CaloHitList &pCluster2, const CaloHitList &pCluster3, float &reducedChiSquared) const
283 {
284  float chiSquaredSum(0.f);
285  unsigned int nSamplingPoints(0), nMatchedSamplingPoints(0);
286  XOverlap xThreeViewOverlapObject(0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f);
287 
288  if (this->PerformThreeViewMatching(pCluster1, pCluster2, pCluster3, chiSquaredSum, nSamplingPoints, nMatchedSamplingPoints,
289  xThreeViewOverlapObject) == STATUS_CODE_NOT_FOUND)
290  return STATUS_CODE_NOT_FOUND;
291 
292  if (nSamplingPoints == 0)
293  return STATUS_CODE_NOT_FOUND;
294 
295  reducedChiSquared = chiSquaredSum / nSamplingPoints;
296 
297  return STATUS_CODE_SUCCESS;
298 }
299 
300 //------------------------------------------------------------------------------------------------------------------------------------------
301 
302 template <typename T>
303 StatusCode NViewDeltaRayMatchingAlgorithm<T>::PerformThreeViewMatching(const CaloHitList &clusterU, const CaloHitList &clusterV,
304  const CaloHitList &clusterW, float &chiSquaredSum, unsigned int &nSamplingPoints, unsigned int &nMatchedSamplingPoints, XOverlap &xOverlapObject) const
305 {
306  float xMinU(-std::numeric_limits<float>::max()), xMaxU(+std::numeric_limits<float>::max());
307  float xMinV(-std::numeric_limits<float>::max()), xMaxV(+std::numeric_limits<float>::max());
308  float xMinW(-std::numeric_limits<float>::max()), xMaxW(+std::numeric_limits<float>::max());
309 
310  this->GetClusterSpanX(clusterU, xMinU, xMaxU);
311  this->GetClusterSpanX(clusterV, xMinV, xMaxV);
312  this->GetClusterSpanX(clusterW, xMinW, xMaxW);
313 
314  // Need to remove the xPitch from calculations to be consistent with view xSpan calculated in the xOverlapObject
315  const float xMinCentre(std::max(xMinU, std::max(xMinV, xMinW)));
316  const float xMaxCentre(std::min(xMaxU, std::min(xMaxV, xMaxW)));
317  const float xCentreOverlap(xMaxCentre - xMinCentre);
318 
319  if (xCentreOverlap < std::numeric_limits<float>::epsilon())
320  return STATUS_CODE_NOT_FOUND;
321 
322  const float xPitch(0.5 * m_xOverlapWindow);
323  const float xMin(std::max(xMinU, std::max(xMinV, xMinW)) - xPitch);
324  const float xMax(std::min(xMaxU, std::min(xMaxV, xMaxW)) + xPitch);
325  const float xOverlap(xMax - xMin);
326 
327  const HitType hitTypeU(clusterU.front()->GetHitType());
328  const HitType hitTypeV(clusterV.front()->GetHitType());
329  const HitType hitTypeW(clusterW.front()->GetHitType());
330  const float pitchU{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeU)};
331  const float pitchV{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeV)};
332  const float pitchW{LArGeometryHelper::GetWirePitch(this->GetPandora(), hitTypeW)};
333  const float pitchMax{std::max({pitchU, pitchV, pitchW})};
334 
335  if (hitTypeU == hitTypeV || hitTypeU == hitTypeW || hitTypeV == hitTypeW)
336  return STATUS_CODE_FAILURE;
337 
338  const unsigned int nPoints(1 + static_cast<unsigned int>(xOverlap / xPitch));
339 
340  chiSquaredSum = 0.f;
341  nSamplingPoints = 0;
342  nMatchedSamplingPoints = 0;
343 
344  for (unsigned int n = 0; n < nPoints; ++n)
345  {
346  const float x(xMin + (xMax - xMin) * (static_cast<float>(n) + 0.5f) / static_cast<float>(nPoints));
347  const float xmin(x - xPitch);
348  const float xmax(x + xPitch);
349 
350  try
351  {
352  float zMinU(0.f), zMinV(0.f), zMinW(0.f), zMaxU(0.f), zMaxV(0.f), zMaxW(0.f);
353  this->GetClusterSpanZ(clusterU, xmin, xmax, zMinU, zMaxU);
354  this->GetClusterSpanZ(clusterV, xmin, xmax, zMinV, zMaxV);
355  this->GetClusterSpanZ(clusterW, xmin, xmax, zMinW, zMaxW);
356 
357  const float zU(0.5f * (zMinU + zMaxU));
358  const float zV(0.5f * (zMinV + zMaxV));
359  const float zW(0.5f * (zMinW + zMaxW));
360 
361  const float dzU(zMaxU - zMinU);
362  const float dzV(zMaxV - zMinV);
363  const float dzW(zMaxW - zMinW);
364  const float dzPitch(pitchMax);
365 
366  const float zprojU(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeV, hitTypeW, zV, zW));
367  const float zprojV(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeW, hitTypeU, zW, zU));
368  const float zprojW(LArGeometryHelper::MergeTwoPositions(this->GetPandora(), hitTypeU, hitTypeV, zU, zV));
369 
370  ++nSamplingPoints;
371 
372  const float deltaSquared(((zU - zprojU) * (zU - zprojU) + (zV - zprojV) * (zV - zprojV) + (zW - zprojW) * (zW - zprojW)) / 3.f);
373  const float sigmaSquared(dzU * dzU + dzV * dzV + dzW * dzW + dzPitch * dzPitch);
374  const float pseudoChi2(deltaSquared / sigmaSquared);
375 
376  chiSquaredSum += pseudoChi2;
377 
378  if (pseudoChi2 < m_pseudoChi2Cut)
379  ++nMatchedSamplingPoints;
380  }
381  catch (StatusCodeException &statusCodeException)
382  {
383  if (statusCodeException.GetStatusCode() != STATUS_CODE_NOT_FOUND)
384  return statusCodeException.GetStatusCode();
385 
386  continue;
387  }
388  }
389 
390  // Apply tensor threshold cuts
391  if (nSamplingPoints == 0)
392  return STATUS_CODE_NOT_FOUND;
393 
394  const float matchedFraction(static_cast<float>(nMatchedSamplingPoints) / static_cast<float>(nSamplingPoints));
395 
396  if ((matchedFraction < m_minMatchedFraction) || (nMatchedSamplingPoints < m_minMatchedPoints))
397  return STATUS_CODE_NOT_FOUND;
398 
399  xOverlapObject = XOverlap(xMinU, xMaxU, xMinV, xMaxV, xMinW, xMaxW, xCentreOverlap);
400 
401  return STATUS_CODE_SUCCESS;
402 }
403 
404 //------------------------------------------------------------------------------------------------------------------------------------------
405 
406 template <typename T>
407 void NViewDeltaRayMatchingAlgorithm<T>::GetClusterSpanX(const CaloHitList &caloHitList, float &xMin, float &xMax) const
408 {
409  xMin = std::numeric_limits<float>::max();
410  xMax = -std::numeric_limits<float>::max();
411 
412  for (const CaloHit *const pCaloHit : caloHitList)
413  {
414  const float xCoordinate(pCaloHit->GetPositionVector().GetX());
415 
416  if (xCoordinate < xMin)
417  xMin = xCoordinate;
418 
419  if (xCoordinate > xMax)
420  xMax = xCoordinate;
421  }
422 }
423 
424 //------------------------------------------------------------------------------------------------------------------------------------------
425 
426 template <typename T>
428  const CaloHitList &caloHitList, const float xMin, const float xMax, float &zMin, float &zMax) const
429 {
430  zMin = std::numeric_limits<float>::max();
431  zMax = -std::numeric_limits<float>::max();
432 
433  bool found(false);
434  for (const CaloHit *const pCaloHit : caloHitList)
435  {
436  const float xCoordinate(pCaloHit->GetPositionVector().GetX());
437  const float zCoordinate(pCaloHit->GetPositionVector().GetZ());
438 
439  if ((xCoordinate < xMin) || (xCoordinate > xMax))
440  continue;
441 
442  found = true;
443 
444  if (zCoordinate < zMin)
445  zMin = zCoordinate;
446 
447  if (zCoordinate > zMax)
448  zMax = zCoordinate;
449  }
450 
451  if (!found)
452  throw StatusCodeException(STATUS_CODE_NOT_FOUND);
453 
454  return STATUS_CODE_SUCCESS;
455 }
456 
457 //------------------------------------------------------------------------------------------------------------------------------------------
458 
459 template <typename T>
461  const HitType &thirdViewHitType, const ParticleFlowObject *const pParentMuon, CartesianPointVector &projectedPositions) const
462 {
463  ClusterList muonClusterList1, muonClusterList2;
464 
465  HitTypeVector hitTypes({TPC_VIEW_U, TPC_VIEW_V, TPC_VIEW_W});
466 
467  hitTypes.erase(std::find(hitTypes.begin(), hitTypes.end(), thirdViewHitType));
468 
469  LArPfoHelper::GetClusters(pParentMuon, hitTypes[0], muonClusterList1);
470  LArPfoHelper::GetClusters(pParentMuon, hitTypes[1], muonClusterList2);
471 
472  if ((muonClusterList1.size() != 1) || (muonClusterList2.size() != 1))
473  return STATUS_CODE_NOT_FOUND;
474 
475  return (this->GetProjectedPositions(muonClusterList1.front(), muonClusterList2.front(), projectedPositions));
476 }
477 
478 //------------------------------------------------------------------------------------------------------------------------------------------
479 
480 template <typename T>
482  const Cluster *const pCluster1, const Cluster *const pCluster2, CartesianPointVector &projectedPositions) const
483 {
484  float xMin1(0.f), xMax1(0.f), xMin2(0.f), xMax2(0.f);
485  pCluster1->GetClusterSpanX(xMin1, xMax1);
486  pCluster2->GetClusterSpanX(xMin2, xMax2);
487 
488  const float xPitch(0.5 * m_xOverlapWindow);
489  const float xMin(std::max(xMin1, xMin2) - xPitch);
490  const float xMax(std::min(xMax1, xMax2) + xPitch);
491  const float xOverlap(xMax - xMin);
492 
493  if (xOverlap < std::numeric_limits<float>::epsilon())
494  return STATUS_CODE_NOT_FOUND;
495 
496  const HitType hitType1(LArClusterHelper::GetClusterHitType(pCluster1));
497  const HitType hitType2(LArClusterHelper::GetClusterHitType(pCluster2));
498 
499  if (hitType1 == hitType2)
500  throw StatusCodeException(STATUS_CODE_FAILURE);
501 
502  const unsigned int nPoints(1 + static_cast<unsigned int>(xOverlap / xPitch));
503 
504  // Projection into third view
505  for (unsigned int n = 0; n < nPoints; ++n)
506  {
507  const float x(xMin + (xMax - xMin) * (static_cast<float>(n) + 0.5f) / static_cast<float>(nPoints));
508  const float xmin(x - xPitch);
509  const float xmax(x + xPitch);
510 
511  try
512  {
513  float zMin1(0.f), zMin2(0.f), zMax1(0.f), zMax2(0.f);
514  pCluster1->GetClusterSpanZ(xmin, xmax, zMin1, zMax1);
515  pCluster2->GetClusterSpanZ(xmin, xmax, zMin2, zMax2);
516 
517  const float z1(0.5f * (zMin1 + zMax1));
518  const float z2(0.5f * (zMin2 + zMax2));
519 
520  float chi2;
521  CartesianVector projection(0.f, 0.f, 0.f);
523  this->GetPandora(), hitType1, hitType2, CartesianVector(x, 0.f, z1), CartesianVector(x, 0.f, z2), projection, chi2);
524 
525  projectedPositions.push_back(projection);
526  }
527  catch (StatusCodeException &statusCodeException)
528  {
529  if (statusCodeException.GetStatusCode() != STATUS_CODE_NOT_FOUND)
530  throw statusCodeException.GetStatusCode();
531 
532  continue;
533  }
534  }
535 
536  // Reject if projection is not good
537  if (projectedPositions.size() < m_minProjectedPositions)
538  return STATUS_CODE_NOT_FOUND;
539 
540  return STATUS_CODE_SUCCESS;
541 }
542 
543 //------------------------------------------------------------------------------------------------------------------------------------------
544 
545 template <typename T>
546 StatusCode NViewDeltaRayMatchingAlgorithm<T>::CollectHitsFromMuon(const Cluster *const pCluster1, const Cluster *const pCluster2,
547  const Cluster *const pThirdViewCluster, const ParticleFlowObject *const pParentMuon, const float minDistanceFromMuon,
548  const float maxDistanceToCollected, CaloHitList &collectedHits) const
549 {
550  HitType thirdViewHitType(TPC_VIEW_U);
551  CartesianPointVector deltaRayProjectedPositions;
552 
553  if (!pThirdViewCluster)
554  {
555  if (this->GetProjectedPositions(pCluster1, pCluster2, deltaRayProjectedPositions) != STATUS_CODE_SUCCESS)
556  return STATUS_CODE_NOT_FOUND;
557 
558  for (const HitType hitType : {TPC_VIEW_U, TPC_VIEW_V, TPC_VIEW_W})
559  {
560  if ((LArClusterHelper::GetClusterHitType(pCluster1) != hitType) && (LArClusterHelper::GetClusterHitType(pCluster2) != hitType))
561  {
562  thirdViewHitType = hitType;
563  break;
564  }
565  }
566  }
567  else
568  {
569  CaloHitList deltaRayCaloHitList;
570  pThirdViewCluster->GetOrderedCaloHitList().FillCaloHitList(deltaRayCaloHitList);
571 
572  for (const CaloHit *const pCaloHit : deltaRayCaloHitList)
573  deltaRayProjectedPositions.push_back(pCaloHit->GetPositionVector());
574 
575  thirdViewHitType = LArClusterHelper::GetClusterHitType(pThirdViewCluster);
576  }
577 
578  ClusterList muonClusterList;
579  LArPfoHelper::GetClusters(pParentMuon, thirdViewHitType, muonClusterList);
580 
581  if (muonClusterList.size() != 1)
582  return STATUS_CODE_NOT_FOUND;
583 
584  const Cluster *const pMuonCluster(muonClusterList.front());
585 
586  // To avoid fluctuatuions, parameterise the muon track
587  CartesianVector positionOnMuon(0.f, 0.f, 0.f), muonDirection(0.f, 0.f, 0.f);
588  if (this->ParameteriseMuon(pParentMuon, deltaRayProjectedPositions, thirdViewHitType, positionOnMuon, muonDirection) != STATUS_CODE_SUCCESS)
589  return STATUS_CODE_NOT_FOUND;
590 
591  this->CollectHitsFromMuon(
592  positionOnMuon, muonDirection, pMuonCluster, deltaRayProjectedPositions, minDistanceFromMuon, maxDistanceToCollected, collectedHits);
593 
594  if (collectedHits.empty())
595  return STATUS_CODE_NOT_FOUND;
596 
597  // Catch if delta ray has travelled along muon
598  if ((static_cast<float>(collectedHits.size()) / pMuonCluster->GetNCaloHits()) > m_maxCosmicRayHitFraction)
599  return STATUS_CODE_NOT_FOUND;
600 
601  return STATUS_CODE_SUCCESS;
602 }
603 
604 //------------------------------------------------------------------------------------------------------------------------------------------
605 
606 template <typename T>
607 void NViewDeltaRayMatchingAlgorithm<T>::CollectHitsFromMuon(const CartesianVector &positionOnMuon, const CartesianVector &muonDirection,
608  const Cluster *const pMuonCluster, const CartesianPointVector &deltaRayProjectedPositions, const float &minDistanceFromMuon,
609  const float maxDistanceToCollected, CaloHitList &collectedHits) const
610 {
611  CaloHitList cosmicRayHitList;
612  pMuonCluster->GetOrderedCaloHitList().FillCaloHitList(cosmicRayHitList);
613 
614  bool hitsAdded(true);
615  while (hitsAdded)
616  {
617  hitsAdded = false;
618 
619  for (const CaloHit *const pCaloHit : cosmicRayHitList)
620  {
621  if (std::find(collectedHits.begin(), collectedHits.end(), pCaloHit) != collectedHits.end())
622  continue;
623 
624  const float distanceToCollectedHits(std::min(LArMuonLeadingHelper::GetClosestDistance(pCaloHit, deltaRayProjectedPositions),
625  collectedHits.empty() ? std::numeric_limits<float>::max()
626  : LArClusterHelper::GetClosestDistance(pCaloHit->GetPositionVector(), collectedHits)));
627  const float distanceToMuonHits(muonDirection.GetCrossProduct(pCaloHit->GetPositionVector() - positionOnMuon).GetMagnitude());
628 
629  if ((std::fabs(distanceToMuonHits - distanceToCollectedHits) < std::numeric_limits<float>::epsilon()) ||
630  (distanceToMuonHits < minDistanceFromMuon) || (distanceToCollectedHits > distanceToMuonHits) ||
631  (distanceToCollectedHits > maxDistanceToCollected))
632  {
633  continue;
634  }
635 
636  collectedHits.push_back(pCaloHit);
637  hitsAdded = true;
638  }
639  }
640 }
641 
642 //------------------------------------------------------------------------------------------------------------------------------------------
643 
644 template <typename T>
645 StatusCode NViewDeltaRayMatchingAlgorithm<T>::ParameteriseMuon(const ParticleFlowObject *const pParentMuon,
646  const Cluster *const pDeltaRayCluster, CartesianVector &positionOnMuon, CartesianVector &muonDirection) const
647 {
648  CaloHitList deltaRayHitList;
649  pDeltaRayCluster->GetOrderedCaloHitList().FillCaloHitList(deltaRayHitList);
650 
651  CartesianPointVector deltaRayProjectedPositions;
652 
653  for (const CaloHit *const pCaloHit : deltaRayHitList)
654  deltaRayProjectedPositions.push_back(pCaloHit->GetPositionVector());
655 
656  return this->ParameteriseMuon(
657  pParentMuon, deltaRayProjectedPositions, LArClusterHelper::GetClusterHitType(pDeltaRayCluster), positionOnMuon, muonDirection);
658 }
659 
660 //------------------------------------------------------------------------------------------------------------------------------------------
661 
662 template <typename T>
663 StatusCode NViewDeltaRayMatchingAlgorithm<T>::ParameteriseMuon(const ParticleFlowObject *const pParentMuon,
664  const CartesianPointVector &deltaRayProjectedPositions, const HitType hitType, CartesianVector &positionOnMuon, CartesianVector &muonDirection) const
665 {
666  ClusterList muonClusterList;
667  LArPfoHelper::GetClusters(pParentMuon, hitType, muonClusterList);
668 
669  if (muonClusterList.size() != 1)
670  return STATUS_CODE_NOT_FOUND;
671 
672  CartesianPointVector muonProjectedPositions;
673  if (this->ProjectMuonPositions(hitType, pParentMuon, muonProjectedPositions) != STATUS_CODE_SUCCESS)
674  return STATUS_CODE_NOT_FOUND;
675 
676  const Cluster *const pMuonCluster(muonClusterList.front());
677  const float slidingFitPitch(LArGeometryHelper::GetWirePitch(this->GetPandora(), hitType));
678  const TwoDSlidingFitResult slidingFitResult(pMuonCluster, 40, slidingFitPitch);
679 
680  CartesianVector deltaRayVertex(0.f, 0.f, 0.f), muonVertex(0.f, 0.f, 0.f);
681  LArMuonLeadingHelper::GetClosestPositions(deltaRayProjectedPositions, pMuonCluster, deltaRayVertex, muonVertex);
682 
683  const StatusCode status(LArMuonLeadingHelper::GetClosestPosition(
684  muonVertex, muonProjectedPositions, pMuonCluster, m_maxDistanceToCluster, m_maxDistanceToReferencePoint, positionOnMuon));
685 
686  if (status != STATUS_CODE_SUCCESS)
687  return STATUS_CODE_NOT_FOUND;
688 
689  float rL(0.f), rT(0.f);
690  slidingFitResult.GetLocalPosition(positionOnMuon, rL, rT);
691 
692  if (slidingFitResult.GetGlobalFitDirection(rL, muonDirection) != STATUS_CODE_SUCCESS)
693  return STATUS_CODE_NOT_FOUND;
694 
695  return STATUS_CODE_SUCCESS;
696 }
697 
698 //------------------------------------------------------------------------------------------------------------------------------------------
699 
700 template <typename T>
701 void NViewDeltaRayMatchingAlgorithm<T>::SplitMuonCluster(const std::string &clusterListName, const Cluster *const pMuonCluster,
702  const CaloHitList &collectedHits, const Cluster *&pDeltaRayCluster) const
703 {
704  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::ReplaceCurrentList<Cluster>(*this, clusterListName));
705 
706  CaloHitList muonCaloHitList;
707  pMuonCluster->GetOrderedCaloHitList().FillCaloHitList(muonCaloHitList);
708 
709  for (const CaloHit *const pCaloHit : muonCaloHitList)
710  {
711  if (std::find(collectedHits.begin(), collectedHits.end(), pCaloHit) != collectedHits.end())
712  {
713  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::RemoveFromCluster(*this, pMuonCluster, pCaloHit));
714 
715  if (!pDeltaRayCluster)
716  {
717  const ClusterList *pTemporaryList(nullptr);
718  std::string temporaryListName, currentListName;
719  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentListName<Cluster>(*this, currentListName));
720  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
721  PandoraContentApi::CreateTemporaryListAndSetCurrent<ClusterList>(*this, pTemporaryList, temporaryListName));
722 
723  PandoraContentApi::Cluster::Parameters parameters;
724  parameters.m_caloHitList.push_back(pCaloHit);
725 
726  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::Cluster::Create(*this, parameters, pDeltaRayCluster));
727  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::SaveList<Cluster>(*this, temporaryListName, currentListName));
728  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::ReplaceCurrentList<Cluster>(*this, currentListName));
729  }
730  else
731  {
732  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::AddToCluster(*this, pDeltaRayCluster, pCaloHit));
733  }
734  }
735  }
736 }
737 
738 //------------------------------------------------------------------------------------------------------------------------------------------
739 
740 template <typename T>
742 {
743  std::sort(protoParticleVector.begin(), protoParticleVector.end(),
744  [](const ProtoParticle &a, const ProtoParticle &b) -> bool
745  {
746  unsigned int aHitTotal(0);
747  for (const Cluster *const pCluster : a.m_clusterList)
748  aHitTotal += pCluster->GetNCaloHits();
749 
750  unsigned int bHitTotal(0);
751  for (const Cluster *const pCluster : b.m_clusterList)
752  bHitTotal += pCluster->GetNCaloHits();
753 
754  return (aHitTotal > bHitTotal);
755  });
756 
757  for (ProtoParticle protoParticle : protoParticleVector)
758  {
759  float longestSpan(-std::numeric_limits<float>::max()), spanMinX(0.f), spanMaxX(0.f);
760 
761  for (const Cluster *const pCluster : protoParticle.m_clusterList)
762  {
763  float minX(0.f), maxX(0.f);
764  pCluster->GetClusterSpanX(minX, maxX);
765 
766  const float span(maxX - minX);
767 
768  if (span > longestSpan)
769  {
770  longestSpan = span;
771  spanMinX = minX;
772  spanMaxX = maxX;
773  }
774  }
775 
776  for (const Cluster *const pCluster : protoParticle.m_clusterList)
777  {
778  ClusterList collectedClusters;
779  this->CollectStrayClusters(pCluster, spanMinX, spanMaxX, collectedClusters);
780 
781  if (!collectedClusters.empty())
782  this->AddInStrayClusters(pCluster, collectedClusters);
783  }
784  }
785 
786  return (this->CreateThreeDParticles(protoParticleVector));
787 }
788 
789 //------------------------------------------------------------------------------------------------------------------------------------------
790 
791 template <typename T>
793  const Cluster *const pClusterToEnlarge, const float rangeMinX, const float rangeMaxX, ClusterList &collectedClusters)
794 {
795  const HitType hitType(LArClusterHelper::GetClusterHitType(pClusterToEnlarge));
796  const ClusterList &strayClusterList((hitType == TPC_VIEW_U) ? m_strayClusterListU
797  : (hitType == TPC_VIEW_V) ? m_strayClusterListV
800  const DeltaRayMatchingContainers::ClusterProximityMap::const_iterator clusterProximityIter(clusterProximityMap.find(pClusterToEnlarge));
801 
802  if (clusterProximityIter == clusterProximityMap.end())
803  return;
804 
805  for (const Cluster *const pNearbyCluster : clusterProximityIter->second)
806  {
807  if (std::find(strayClusterList.begin(), strayClusterList.end(), pNearbyCluster) == strayClusterList.end())
808  continue;
809 
810  float xMin(-std::numeric_limits<float>::max()), xMax(+std::numeric_limits<float>::max());
811  pNearbyCluster->GetClusterSpanX(xMin, xMax);
812 
813  if (!(((xMin > rangeMinX) && (xMin < rangeMaxX)) || ((xMax > rangeMinX) && (xMax < rangeMaxX))))
814  continue;
815 
816  if (LArClusterHelper::GetClosestDistance(pClusterToEnlarge, pNearbyCluster) > m_strayClusterSeparation)
817  continue;
818 
819  if (std::find(collectedClusters.begin(), collectedClusters.end(), pNearbyCluster) == collectedClusters.end())
820  collectedClusters.push_back(pNearbyCluster);
821  }
822 }
823 
824 //------------------------------------------------------------------------------------------------------------------------------------------
825 
826 template <typename T>
827 void NViewDeltaRayMatchingAlgorithm<T>::AddInStrayClusters(const Cluster *const pClusterToEnlarge, const ClusterList &collectedClusters)
828 {
829  this->UpdateUponDeletion(pClusterToEnlarge);
830 
831  for (const Cluster *const pCollectedCluster : collectedClusters)
832  {
833  this->UpdateUponDeletion(pCollectedCluster);
834 
835  std::string clusterListName(this->GetClusterListName(LArClusterHelper::GetClusterHitType(pClusterToEnlarge)));
836  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
837  PandoraContentApi::MergeAndDeleteClusters(*this, pClusterToEnlarge, pCollectedCluster, clusterListName, clusterListName));
838  }
839 
840  m_deltaRayMatchingContainers.AddClustersToContainers({pClusterToEnlarge}, {nullptr});
841 }
842 
843 //------------------------------------------------------------------------------------------------------------------------------------------
844 
845 template <typename T>
846 void NViewDeltaRayMatchingAlgorithm<T>::UpdateUponDeletion(const Cluster *const pDeletedCluster)
847 {
848  const HitType hitType(LArClusterHelper::GetClusterHitType(pDeletedCluster));
849  ClusterList &strayClusterList((hitType == TPC_VIEW_U) ? m_strayClusterListU
850  : (hitType == TPC_VIEW_V) ? m_strayClusterListV
852  const ClusterList::const_iterator strayClusterIter(std::find(strayClusterList.begin(), strayClusterList.end(), pDeletedCluster));
853 
854  if (strayClusterIter != strayClusterList.end())
855  strayClusterList.erase(strayClusterIter);
856 
858  const bool isMuonCluster(clusterToPfoMap.find(pDeletedCluster) != clusterToPfoMap.end());
859 
861 
862  if (!isMuonCluster)
864 }
865 
866 //------------------------------------------------------------------------------------------------------------------------------------------
867 
868 template <typename T>
869 void NViewDeltaRayMatchingAlgorithm<T>::UpdateForNewClusters(const ClusterVector &newClusterVector, const PfoVector &pfoVector)
870 {
871  m_deltaRayMatchingContainers.AddClustersToContainers(newClusterVector, pfoVector);
872 
873  for (unsigned int i = 0; i < newClusterVector.size(); i++)
874  {
875  const Cluster *const pNewCluster(newClusterVector.at(i));
876  const ParticleFlowObject *const pMuonPfo(pfoVector.at(i));
877 
878  // ATTN: Only add delta ray clusters into the tensor
879  if (!pMuonPfo)
881  }
882 }
883 
884 //------------------------------------------------------------------------------------------------------------------------------------------
885 
886 template <typename T>
888 {
890 
891  m_strayClusterListU.clear();
892  m_strayClusterListV.clear();
893  m_strayClusterListW.clear();
894 
896 }
897 
898 //------------------------------------------------------------------------------------------------------------------------------------------
899 
900 template <typename T>
901 StatusCode NViewDeltaRayMatchingAlgorithm<T>::ReadSettings(const TiXmlHandle xmlHandle)
902 {
903  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "MuonPfoListName", m_muonPfoListName));
904 
905  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
906  XmlHelper::ReadValue(xmlHandle, "SearchRegion1D", m_deltaRayMatchingContainers.m_searchRegion1D));
907 
908  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "PseudoChi2Cut", m_pseudoChi2Cut));
909 
910  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "OverlapWindow", m_xOverlapWindow));
911 
912  PANDORA_RETURN_RESULT_IF_AND_IF(
913  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinMatchedFraction", m_minMatchedFraction));
914 
915  PANDORA_RETURN_RESULT_IF_AND_IF(
916  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinMatchedPoints", m_minMatchedPoints));
917 
918  PANDORA_RETURN_RESULT_IF_AND_IF(
919  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinProjectedPositions", m_minProjectedPositions));
920 
921  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
922  XmlHelper::ReadValue(xmlHandle, "MaxCosmicRayHitFraction", m_maxCosmicRayHitFraction));
923 
924  PANDORA_RETURN_RESULT_IF_AND_IF(
925  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxDistanceToCluster", m_maxDistanceToCluster));
926 
927  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
928  XmlHelper::ReadValue(xmlHandle, "MaxDistanceToReferencePoint", m_maxDistanceToReferencePoint));
929 
930  PANDORA_RETURN_RESULT_IF_AND_IF(
931  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "StrayClusterSeparation", m_strayClusterSeparation));
932 
933  return NViewMatchingAlgorithm<T>::ReadSettings(xmlHandle);
934 }
935 
936 //------------------------------------------------------------------------------------------------------------------------------------------
937 
940 
941 } // namespace lar_content
Float_t x
Definition: compare.C:6
void SelectInputClusters(const pandora::ClusterList *const pInputClusterList, pandora::ClusterList &selectedClusterList) const
Select a subset of input clusters for processing in this algorithm.
std::vector< ProtoParticle > ProtoParticleVector
static float GetClosestDistance(const pandora::Cluster *const pCluster, const pandora::CartesianPointVector &cartesianPointVector)
Get closest distance between a specified cluster and list of positions.
Header file for the kd tree linker algo template class.
std::map< const pandora::Cluster *, const pandora::ParticleFlowObject * > ClusterToPfoMap
Header file for the pfo helper class.
pandora::StatusCode GetProjectedPositions(const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2, pandora::CartesianPointVector &projectedPositions) const
Use two clusters from different views to calculate projected positions in the remaining third view...
const ClusterProximityMap & GetClusterProximityMap(const pandora::HitType hitType) const
Get the mapping of clusters to to their neighbouring clusters.
pandora::StatusCode GetClusterSpanZ(const pandora::CaloHitList &caloHitList, const float xMin, const float xMax, float &zMin, float &zMax) const
Calculate the zSpan of a list of CaloHits in a specified x range.
unsigned int m_minProjectedPositions
The threshold number of projected points for a good projection.
static void GetClusters(const pandora::PfoList &pfoList, const pandora::HitType &hitType, pandora::ClusterList &clusterList)
Get a list of clusters of a particular hit type from a list of pfos.
pandora::StatusCode CollectHitsFromMuon(const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2, const pandora::Cluster *const pThirdViewCluster, const pandora::ParticleFlowObject *const pParentMuon, const float minDistanceFromMuon, const float maxDistanceToCollected, pandora::CaloHitList &collectedHits) const
In one view, pull out any hits from a cosmic ray cluster that belong to the child delta ray cluster...
pandora::ClusterList m_strayClusterListU
The list of U clusters that do not pass the tensor threshold requirement.
virtual void TidyUp()
Tidy member variables in derived class.
void UpdateForNewClusters(const pandora::ClusterVector &newClusterVector, const pandora::PfoVector &pfoVector)
Add a new cluster to algorithm ownership maps and, if it a delta ray cluster, to the underlying match...
void TidyUp()
Tidy member variables in derived class.
float m_pseudoChi2Cut
Pseudo chi2 cut for three view matching.
intermediate_table::const_iterator const_iterator
static pandora::HitType GetClusterHitType(const pandora::Cluster *const pCluster)
Get the hit type associated with a two dimensional cluster.
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
void FillContainers(const pandora::PfoList &inputPfoList, const pandora::ClusterList &inputClusterList1, const pandora::ClusterList &inputClusterList2=pandora::ClusterList(), const pandora::ClusterList &inputClusterList3=pandora::ClusterList())
Fill the HitToClusterMap, the ClusterProximityMap and the ClusterToPfoMap in all input views...
void UpdateUponDeletion(const pandora::Cluster *const pDeletedCluster)
Update to reflect cluster deletion.
void PrepareInputClusters(pandora::ClusterList &preparedClusterList)
Perform any preparatory steps required on the input clusters, e.g. caching expensive fit results...
void CollectStrayClusters(const pandora::Cluster *const pClusterToEnlarge, const float rangeMinX, const float rangeMaxX, pandora::ClusterList &collectedClusters)
Collect the stray clusters that are close to a specified cluster and that lie within a given x range...
void GetClusterSpanX(const pandora::CaloHitList &caloHitList, float &xMin, float &xMax) const
Calculate the xSpan of a list of CaloHits.
TFile f
Definition: plotHisto.C:6
pandora::StatusCode PerformThreeViewMatching(const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2, const pandora::Cluster *const pCluster3, float &reducedChiSquared) const
To determine how well three clusters (one in each view) map onto one another expressing this in terms...
void UpdateUponDeletion(const pandora::Cluster *const pDeletedCluster)
Update to reflect cluster deletion.
Header file for the geometry helper class.
void UpdateForNewCluster(const pandora::Cluster *const pNewCluster)
Update to reflect addition of a new cluster to the problem space.
float m_xOverlapWindow
The maximum allowed displacement in x position.
std::string m_muonPfoListName
The list of reconstructed cosmic ray pfos.
virtual bool DoesClusterPassTensorThreshold(const pandora::Cluster *const pCluster) const =0
To check whether a given cluster meets the requirements to be added into the matching container (tens...
Header file for the cluster helper class.
virtual pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
Header file for the muon leading helper class.
pandora::StatusCode ParameteriseMuon(const pandora::ParticleFlowObject *const pParentMuon, const pandora::Cluster *const pDeltaRayCluster, pandora::CartesianVector &positionOnMuon, pandora::CartesianVector &muonDirection) const
Parameterise the projection of a cosmic ray track in order to avoid poor/sparse projections.
void AddClustersToContainers(const pandora::ClusterVector &newClusterVector, const pandora::PfoVector &pfoVector)
Add a list of clusters to the hit to cluster and cluster proximity maps and, if appropriate, to the cluster to pfo map.
Header file for the lar two dimensional sliding fit result class.
static void GetClosestPositions(const pandora::CartesianPointVector &cartesianPointVector1, const pandora::Cluster *const pCluster2, pandora::CartesianVector &outputPosition1, pandora::CartesianVector &outputPosition2)
Get the closest positions between a list of positions and a cluster.
void FillStrayClusterList(const pandora::HitType hitType)
Fill the stray cluster list with clusters that do not pass the tensor threshold requirement.
pandora::ClusterList m_strayClusterListV
The list of V clusters that do not pass the tensor threshold requirement.
pandora::ClusterList m_clusterList
List of 2D clusters in a 3D proto particle.
pandora::ClusterList m_strayClusterListW
The list of W clusters that do not pass the tensor threshold requirement.
const pandora::ClusterList & GetInputClusterList(const pandora::HitType hitType) const
Get the input cluster list corresponding to a specified hit type.
const ClusterToPfoMap & GetClusterToPfoMap(const pandora::HitType hitType) const
Get the mapping of clusters to the pfos to which they belong.
static float MergeTwoPositions(const pandora::Pandora &pandora, const pandora::HitType view1, const pandora::HitType view2, const float position1, const float position2)
Merge two views (U,V) to give a third view (Z).
void SplitMuonCluster(const std::string &clusterListName, const pandora::Cluster *const pMuonCluster, const pandora::CaloHitList &collectedHits, const pandora::Cluster *&pDeltaRayCluster) const
Move a list of hits from a cosmic ray cluster into the given child delta ray cluster.
static float GetWirePitch(const pandora::Pandora &pandora, const pandora::HitType view, const float maxWirePitchDiscrepancy=0.01)
Return the wire pitch.
void AddInStrayClusters(const pandora::Cluster *const pClusterToEnlarge, const pandora::ClusterList &collectedClusters)
Merge in the collected stray clusters of a given delta ray cluster.
float m_strayClusterSeparation
The maximum allowed separation of a stray cluster and a delta ray cluster for merge.
void GetNearbyMuonPfos(const pandora::Cluster *const pCluster, pandora::ClusterList &consideredClusters, pandora::PfoList &nearbyMuonPfos) const
Use the cluster proximity map to travel along paths of nearby clusters finding the cosmic ray cluster...
static pandora::StatusCode GetClosestPosition(const pandora::CartesianVector &referencePoint, const pandora::CartesianPointVector &cartesianPointVector, const pandora::Cluster *const pCluster, const float maxDistanceToCluster, const float maxDistanceToReferencePoint, pandora::CartesianVector &closestPosition)
Get the closest position from an input list of projected positions that lies close to both a referenc...
float m_maxDistanceToReferencePoint
the maximum distance of a projected point to the cosmic ray vertex used when parameterising the cosmi...
pandora::StatusCode GetGlobalFitDirection(const float rL, pandora::CartesianVector &direction) const
Get global fit direction for a given longitudinal coordinate.
std::map< const pandora::Cluster *, pandora::ClusterList > ClusterProximityMap
Header file for the lar track overlap result class.
HitType
Definition: HitType.h:12
span(IterB &&b, IterE &&e, Adaptor &&adaptor) -> span< decltype(adaptor(std::forward< IterB >(b))), decltype(adaptor(std::forward< IterE >(e))) >
Header file for the three view matching control class.
const std::string & GetClusterListName(const pandora::HitType hitType) const
Get the cluster list name corresponding to a specified hit type.
Char_t n[5]
virtual bool CreateThreeDParticles(const ProtoParticleVector &protoParticleVector)
Create particles using findings from recent algorithm processing.
std::vector< art::Ptr< recob::Cluster > > ClusterVector
float m_maxCosmicRayHitFraction
The maximum allowed fraction of hits to be removed from the cosmic ray track.
bool CreatePfos(ProtoParticleVector &protoParticleVector)
Create delta ray pfos maxmising completeness by searching for and merging in any stray clusters...
pandora::StatusCode ProjectMuonPositions(const pandora::HitType &thirdViewHitType, const pandora::ParticleFlowObject *const pParentMuon, pandora::CartesianPointVector &projectedPositions) const
Use two views of a cosmic ray pfo to calculate projected positions in a given the third view...
DeltaRayMatchingContainers m_deltaRayMatchingContainers
The class of hit, cluster and pfo ownership and proximity maps.
XOverlap class.
Definition: LArXOverlap.h:17
void ClearContainers()
Empty all algorithm containers.
Header file for the two view matching control class.
void GetLocalPosition(const pandora::CartesianVector &position, float &rL, float &rT) const
Get local sliding fit coordinates for a given global position.
Header file for the lar track two view overlap result class.
pandora::StatusCode GetMuonCluster(const pandora::PfoList &commonMuonPfoList, const pandora::HitType hitType, const pandora::Cluster *&pMuonCluster) const
Return the cluster of the common cosmic ray pfo in a given view (function demands there to be only on...
void RemoveClusterFromContainers(const pandora::Cluster *const pDeletedCluster)
Remove an input cluster&#39;s hits from the hit to cluster and cluster proximity maps and...
float m_minMatchedFraction
The threshold matched fraction of sampling points for a good match.
unsigned int m_minMatchedPoints
The threshold number of matched sampling points for a good match.
float m_maxDistanceToCluster
the maximum distance of a projected point to the cosmic ray cluster used when parameterising the cosm...
float m_searchRegion1D
Search region, applied to each dimension, for look-up from kd-tree.
static float GetClosestDistance(const pandora::ClusterList &clusterList1, const pandora::ClusterList &clusterList2)
Get closest distance between clusters in a pair of cluster lists.