LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
TrackMergeRefinementAlgorithm.cc
Go to the documentation of this file.
1 
8 #include "Pandora/AlgorithmHeaders.h"
9 
11 
15 
16 using namespace pandora;
17 
18 namespace lar_content
19 {
20 
21 TrackMergeRefinementAlgorithm::TrackMergeRefinementAlgorithm() :
22  m_maxLoopIterations(10),
23  m_minClusterLengthSum(75.f),
24  m_minSeparationDistance(0.f),
25  m_minDirectionDeviationCosAngle(0.99f),
26  m_maxPredictedMergePointOffset(5.f),
27  m_distanceToLine(0.35f),
28  m_boundaryTolerance(2.f)
29 {
30 }
31 
32 //------------------------------------------------------------------------------------------------------------------------------------------
33 
35 {
36  const ClusterList *pClusterList(nullptr);
37  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pClusterList));
38 
39  const CaloHitList *pCaloHitList(nullptr);
40  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::GetCurrentList(*this, pCaloHitList));
41 
42  ClusterVector clusterVector;
43  TwoDSlidingFitResultMap microSlidingFitResultMap, macroSlidingFitResultMap;
44  SlidingFitResultMapPair slidingFitResultMapPair({&microSlidingFitResultMap, &macroSlidingFitResultMap});
45 
46  this->InitialiseContainers(pClusterList, LArClusterHelper::SortByNHits, clusterVector, slidingFitResultMapPair);
47 
48  // ATTN: Keep track of created main track clusters so their hits can be protected in future iterations
49  unsigned int loopIterations(0);
50  ClusterList createdMainTrackClusters;
51  while (loopIterations < m_maxLoopIterations)
52  {
53  ++loopIterations;
54 
55  ClusterPairAssociation clusterAssociation;
56  if (!this->FindBestClusterAssociation(clusterVector, slidingFitResultMapPair, clusterAssociation))
57  break;
58 
59  ClusterList unavailableProtectedClusters;
60  this->GetUnavailableProtectedClusters(clusterAssociation, createdMainTrackClusters, unavailableProtectedClusters);
61 
62  ClusterToCaloHitListMap clusterToCaloHitListMap;
63  this->GetHitsInBoundingBox(clusterAssociation.GetUpstreamMergePoint(), clusterAssociation.GetDownstreamMergePoint(), pClusterList,
64  clusterToCaloHitListMap, unavailableProtectedClusters, m_distanceToLine);
65 
66  if (!this->AreExtrapolatedHitsGood(clusterToCaloHitListMap, clusterAssociation))
67  {
68  this->ConsiderClusterAssociation(clusterAssociation, clusterVector, slidingFitResultMapPair);
69  continue;
70  }
71 
72  const ClusterList::const_iterator upstreamIter(
73  std::find(createdMainTrackClusters.begin(), createdMainTrackClusters.end(), clusterAssociation.GetUpstreamCluster()));
74  if (upstreamIter != createdMainTrackClusters.end())
75  createdMainTrackClusters.erase(upstreamIter);
76 
77  const ClusterList::const_iterator downstreamIter(
78  std::find(createdMainTrackClusters.begin(), createdMainTrackClusters.end(), clusterAssociation.GetDownstreamCluster()));
79  if (downstreamIter != createdMainTrackClusters.end())
80  createdMainTrackClusters.erase(downstreamIter);
81 
82  createdMainTrackClusters.push_back(
83  this->CreateMainTrack(clusterAssociation, clusterToCaloHitListMap, pClusterList, clusterVector, slidingFitResultMapPair));
84  }
85 
86  return STATUS_CODE_SUCCESS;
87 }
88 
89 //------------------------------------------------------------------------------------------------------------------------------------------
90 
92  const SlidingFitResultMapPair &slidingFitResultMapPair, ClusterPairAssociation &clusterAssociation) const
93 {
94  bool foundAssociation(false);
95  float maxLength(0.f);
96 
97  // ATTN: Find the associated cluster pair with the longest length sum
98  for (ClusterVector::const_iterator currentIter = clusterVector.begin(); currentIter != clusterVector.end(); ++currentIter)
99  {
100  const Cluster *const pCurrentCluster(*currentIter);
101 
102  const TwoDSlidingFitResultMap::const_iterator currentMicroFitIter(slidingFitResultMapPair.first->find(pCurrentCluster));
103  if (currentMicroFitIter == slidingFitResultMapPair.first->end())
104  return false;
105 
106  const TwoDSlidingFitResultMap::const_iterator currentMacroFitIter(slidingFitResultMapPair.second->find(pCurrentCluster));
107  if (currentMacroFitIter == slidingFitResultMapPair.second->end())
108  return false;
109 
110  for (ClusterVector::const_iterator testIter = std::next(currentIter); testIter != clusterVector.end(); ++testIter)
111  {
112  const Cluster *const pTestCluster(*testIter);
113 
114  const float lengthSum(LArClusterHelper::GetLength(pCurrentCluster) + LArClusterHelper::GetLength(pTestCluster));
115  if ((lengthSum < maxLength) || (lengthSum < m_minClusterLengthSum))
116  continue;
117 
118  const TwoDSlidingFitResultMap::const_iterator testMicroFitIter(slidingFitResultMapPair.first->find(pTestCluster));
119  if (testMicroFitIter == slidingFitResultMapPair.first->end())
120  continue;
121 
122  const TwoDSlidingFitResultMap::const_iterator testMacroFitIter(slidingFitResultMapPair.second->find(pTestCluster));
123  if (testMacroFitIter == slidingFitResultMapPair.second->end())
124  continue;
125 
126  const bool isCurrentUpstream(LArClusterHelper::SortByPosition(pCurrentCluster, pTestCluster));
127 
128  // ATTN: Ensure that clusters are not contained within one another
129  const float currentMinLayerZ(currentMacroFitIter->second.GetGlobalMinLayerPosition().GetZ()),
130  currentMaxLayerZ(currentMacroFitIter->second.GetGlobalMaxLayerPosition().GetZ());
131  const float testMinLayerZ(testMacroFitIter->second.GetGlobalMinLayerPosition().GetZ()),
132  testMaxLayerZ(testMacroFitIter->second.GetGlobalMaxLayerPosition().GetZ());
133 
134  if (((currentMinLayerZ > testMinLayerZ) && (currentMaxLayerZ < testMaxLayerZ)) ||
135  ((testMinLayerZ > currentMinLayerZ) && (testMaxLayerZ < currentMaxLayerZ)))
136  continue;
137 
138  CartesianVector currentMergePoint(0.f, 0.f, 0.f), testMergePoint(0.f, 0.f, 0.f), currentMergeDirection(0.f, 0.f, 0.f),
139  testMergeDirection(0.f, 0.f, 0.f);
140  if (!this->GetClusterMergingCoordinates(currentMicroFitIter->second, currentMacroFitIter->second, testMacroFitIter->second,
141  !isCurrentUpstream, currentMergePoint, currentMergeDirection) ||
142  !this->GetClusterMergingCoordinates(testMicroFitIter->second, testMacroFitIter->second, currentMacroFitIter->second,
143  isCurrentUpstream, testMergePoint, testMergeDirection))
144  {
145  continue;
146  }
147 
148  if ((isCurrentUpstream && !this->AreClustersAssociated(currentMergePoint, currentMergeDirection, testMergePoint, testMergeDirection)) ||
149  (!isCurrentUpstream && !this->AreClustersAssociated(testMergePoint, testMergeDirection, currentMergePoint, currentMergeDirection)))
150  {
151  continue;
152  }
153 
154  if (isCurrentUpstream)
155  {
156  clusterAssociation = ClusterPairAssociation(
157  currentMergePoint, currentMergeDirection, testMergePoint, testMergeDirection * (-1.f), pCurrentCluster, pTestCluster);
158  }
159  else
160  {
161  clusterAssociation = ClusterPairAssociation(
162  testMergePoint, testMergeDirection, currentMergePoint, currentMergeDirection * (-1.f), pTestCluster, pCurrentCluster);
163  }
164 
165  foundAssociation = true;
166  maxLength = lengthSum;
167  }
168  }
169 
170  return foundAssociation;
171 }
172 
173 //------------------------------------------------------------------------------------------------------------------------------------------
174 
175 bool TrackMergeRefinementAlgorithm::AreClustersAssociated(const CartesianVector &upstreamPoint, const CartesianVector &upstreamDirection,
176  const CartesianVector &downstreamPoint, const CartesianVector &downstreamDirection) const
177 {
178  if (downstreamPoint.GetZ() < upstreamPoint.GetZ())
179  return false;
180 
181  const float separationDistance(std::sqrt(upstreamPoint.GetDistanceSquared(downstreamPoint)));
182  if (separationDistance < m_minSeparationDistance)
183  return false;
184 
185  if (upstreamDirection.GetCosOpeningAngle(downstreamDirection) < m_minDirectionDeviationCosAngle)
186  return false;
187 
188  const CartesianVector predictedDownstreamPoint(upstreamPoint + (upstreamDirection * separationDistance));
189  const float predictedDownstreamOffset((predictedDownstreamPoint - downstreamPoint).GetMagnitude());
190 
191  if (predictedDownstreamOffset > m_maxPredictedMergePointOffset)
192  return false;
193 
194  const CartesianVector predictedUpstreamPoint(downstreamPoint - (downstreamDirection * separationDistance));
195  const float predictedUpstreamOffset((predictedUpstreamPoint - upstreamPoint).GetMagnitude());
196 
197  if (predictedUpstreamOffset > m_maxPredictedMergePointOffset)
198  return false;
199 
200  return true;
201 }
202 
203 //------------------------------------------------------------------------------------------------------------------------------------------
204 
206  const ClusterPairAssociation &clusterAssociation, const ClusterList &createdMainTrackClusters, ClusterList &protectedClusters) const
207 {
208  for (const Cluster *const pMainTrackCluster : createdMainTrackClusters)
209  {
210  if ((pMainTrackCluster != clusterAssociation.GetUpstreamCluster()) && (pMainTrackCluster != clusterAssociation.GetDownstreamCluster()))
211  protectedClusters.push_back(pMainTrackCluster);
212  }
213 }
214 
215 //------------------------------------------------------------------------------------------------------------------------------------------
216 
217 bool TrackMergeRefinementAlgorithm::AreExtrapolatedHitsNearBoundaries(const CaloHitVector &extrapolatedHitVector, ClusterAssociation &clusterAssociation) const
218 {
219  if (extrapolatedHitVector.empty())
220  {
221  const float separationDistance((clusterAssociation.GetUpstreamMergePoint() - clusterAssociation.GetDownstreamMergePoint()).GetMagnitude());
222  return (separationDistance < m_boundaryTolerance);
223  }
224 
225  if (!this->IsNearBoundary(extrapolatedHitVector.front(), clusterAssociation.GetUpstreamMergePoint(), m_boundaryTolerance))
226  return false;
227 
228  if (!this->IsNearBoundary(extrapolatedHitVector.back(), clusterAssociation.GetDownstreamMergePoint(), m_boundaryTolerance))
229  return false;
230 
231  return true;
232 }
233 
234 //------------------------------------------------------------------------------------------------------------------------------------------
235 
237  pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
238 {
239  const Cluster *const upstreamCluster(clusterAssociation.GetUpstreamCluster()), *const downstreamCluster(clusterAssociation.GetDownstreamCluster());
240  const Cluster *const pConsideredCluster(upstreamCluster->GetNCaloHits() > downstreamCluster->GetNCaloHits() ? downstreamCluster : upstreamCluster);
241  RemoveClusterFromContainers(pConsideredCluster, clusterVector, slidingFitResultMapPair);
242 }
243 
244 //------------------------------------------------------------------------------------------------------------------------------------------
245 
247  const ClusterToCaloHitListMap &clusterToCaloHitListMap, const ClusterList *const pClusterList, ClusterVector &clusterVector,
248  SlidingFitResultMapPair &slidingFitResultMapPair) const
249 {
250  // Determine the shower clusters which contain hits that belong to the main track
251  ClusterVector showerClustersToFragment;
252  for (const auto &entry : clusterToCaloHitListMap)
253  {
254  if ((entry.first != clusterAssociation.GetUpstreamCluster()) && (entry.first != clusterAssociation.GetDownstreamCluster()))
255  showerClustersToFragment.push_back(entry.first);
256  }
257 
258  std::sort(showerClustersToFragment.begin(), showerClustersToFragment.end(), LArClusterHelper::SortByNHits);
259 
260  ClusterList remnantClusterList;
261  const Cluster *const pMainTrackCluster(
262  this->RemoveOffAxisHitsFromTrack(clusterAssociation.GetUpstreamCluster(), clusterAssociation.GetUpstreamMergePoint(), false,
263  clusterToCaloHitListMap, remnantClusterList, *slidingFitResultMapPair.first, *slidingFitResultMapPair.second));
264  const Cluster *const pClusterToDelete(
265  this->RemoveOffAxisHitsFromTrack(clusterAssociation.GetDownstreamCluster(), clusterAssociation.GetDownstreamMergePoint(), true,
266  clusterToCaloHitListMap, remnantClusterList, *slidingFitResultMapPair.first, *slidingFitResultMapPair.second));
267 
268  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::MergeAndDeleteClusters(*this, pMainTrackCluster, pClusterToDelete));
269 
270  for (const Cluster *const pShowerCluster : showerClustersToFragment)
271  {
272  const CaloHitList &caloHitsToMerge(clusterToCaloHitListMap.at(pShowerCluster));
273  this->AddHitsToMainTrack(pMainTrackCluster, pShowerCluster, caloHitsToMerge, clusterAssociation, remnantClusterList);
274  }
275 
276  ClusterList createdClusters;
277  this->ProcessRemnantClusters(remnantClusterList, pMainTrackCluster, pClusterList, createdClusters);
278 
279  // ATTN: Cleanup containers - choose to add created clusters back into containers
280  ClusterList modifiedClusters(showerClustersToFragment.begin(), showerClustersToFragment.end());
281  modifiedClusters.push_back(clusterAssociation.GetUpstreamCluster());
282  modifiedClusters.push_back(clusterAssociation.GetDownstreamCluster());
283  createdClusters.push_back(pMainTrackCluster);
284  this->UpdateContainers(createdClusters, modifiedClusters, LArClusterHelper::SortByNHits, clusterVector, slidingFitResultMapPair);
285 
286  return pMainTrackCluster;
287 }
288 
289 //------------------------------------------------------------------------------------------------------------------------------------------
290 
291 StatusCode TrackMergeRefinementAlgorithm::ReadSettings(const pandora::TiXmlHandle xmlHandle)
292 {
293  PANDORA_RETURN_RESULT_IF_AND_IF(
294  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxLoopIterations", m_maxLoopIterations));
295 
296  PANDORA_RETURN_RESULT_IF_AND_IF(
297  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinClusterLengthSum", m_minClusterLengthSum));
298 
299  PANDORA_RETURN_RESULT_IF_AND_IF(
300  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinSeparationDistance", m_minSeparationDistance));
301 
302  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
303  XmlHelper::ReadValue(xmlHandle, "MinDirectionDeviationCosAngle", m_minDirectionDeviationCosAngle));
304 
305  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
306  XmlHelper::ReadValue(xmlHandle, "MaxPredictedMergePointOffset", m_maxPredictedMergePointOffset));
307 
308  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "DistanceToLine", m_distanceToLine));
309 
310  PANDORA_RETURN_RESULT_IF_AND_IF(
311  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "BoundaryTolerance", m_boundaryTolerance));
312 
314 }
315 
316 } // namespace lar_content
static bool SortByNHits(const pandora::Cluster *const pLhs, const pandora::Cluster *const pRhs)
Sort clusters by number of hits, then layer span, then inner layer, then position, then pulse-height.
ClusterAssociation class.
Header file for the track merge refinement class.
Header file for the lar hit width helper class.
bool AreExtrapolatedHitsNearBoundaries(const pandora::CaloHitVector &extrapolatedHitVector, ClusterAssociation &clusterAssociation) const
Check the separation of the extremal extrapolated hits with the cluster merge points or...
void AddHitsToMainTrack(const pandora::Cluster *const pMainTrackCluster, const pandora::Cluster *const pShowerTrackCluster, const pandora::CaloHitList &caloHitsToMerge, const ClusterAssociation &clusterAssociation, pandora::ClusterList &remnantClusterList) const
Remove the hits from a shower cluster that belong to the main track and add them into the main track ...
std::pair< TwoDSlidingFitResultMap *, TwoDSlidingFitResultMap * > SlidingFitResultMapPair
const pandora::CartesianVector GetUpstreamMergePoint() const
Returns the upstream cluster merge point.
void ConsiderClusterAssociation(const ClusterPairAssociation &clusterAssociation, pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
Remove the cluster association from the cluster vector so that the same cluster pair is not considere...
intermediate_table::const_iterator const_iterator
std::unordered_map< const pandora::Cluster *, pandora::CaloHitList > ClusterToCaloHitListMap
virtual pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)=0
float m_maxPredictedMergePointOffset
The threshold separation distance between the predicted and true cluster merge points.
bool GetClusterMergingCoordinates(const TwoDSlidingFitResult &clusterMicroFitResult, const TwoDSlidingFitResult &clusterMacroFitResult, const TwoDSlidingFitResult &associatedMacroFitResult, const bool isEndUpstream, pandora::CartesianVector &clusterMergePosition, pandora::CartesianVector &clusterMergeDirection) const
Get the merging coordinate and direction for an input cluster with respect to an associated cluster...
TFile f
Definition: plotHisto.C:6
const pandora::Cluster * GetUpstreamCluster() const
Returns the address of the upstream cluster.
Header file for the geometry helper class.
static bool SortByPosition(const pandora::Cluster *const pLhs, const pandora::Cluster *const pRhs)
Sort clusters by position, then pulse-height.
float m_minClusterLengthSum
The threshold cluster and associated cluster length sum.
bool AreExtrapolatedHitsGood(const ClusterToCaloHitListMap &clusterToCaloHitListMap, ClusterAssociation &clusterAssociation) const
Perform topological checks on the collected hits to ensure no gaps are present.
void UpdateContainers(const pandora::ClusterList &clustersToAdd, const pandora::ClusterList &clustersToDelete, const T sortFunction, pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
Remove deleted clusters from the cluster vector and sliding fit maps and add in created clusters that...
Header file for the cluster helper class.
void ProcessRemnantClusters(const pandora::ClusterList &remnantClusterList, const pandora::Cluster *const pMainTrackCluster, const pandora::ClusterList *const pClusterList, pandora::ClusterList &createdClusters) const
Process the remnant clusters separating those that stradle the main track.
bool AreClustersAssociated(const pandora::CartesianVector &upstreamPoint, const pandora::CartesianVector &upstreamDirection, const pandora::CartesianVector &downstreamPoint, const pandora::CartesianVector &downstreamDirection) const
Whether two clusters are assoicated to one another.
bool FindBestClusterAssociation(const pandora::ClusterVector &clusterVector, const SlidingFitResultMapPair &slidingFitResultMapPair, ClusterPairAssociation &clusterAssociation) const
Find the best cluster association.
ClusterPairAssociation class.
static float GetLength(const pandora::Cluster *const pCluster)
Get length of cluster.
const pandora::CartesianVector GetDownstreamMergePoint() const
Returns the downstream cluster merge point.
std::unordered_map< const pandora::Cluster *, TwoDSlidingFitResult > TwoDSlidingFitResultMap
void InitialiseContainers(const pandora::ClusterList *pClusterList, const T sortFunction, pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
Fill the cluster vector and sliding fit maps with clusters that are determined to be track-like...
void RemoveClusterFromContainers(const pandora::Cluster *const pClustertoRemove, pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
Remove a cluster from the cluster vector and sliding fit maps.
const pandora::Cluster * CreateMainTrack(const ClusterPairAssociation &clusterAssociation, const ClusterToCaloHitListMap &clusterToCaloHitListMap, const pandora::ClusterList *pClusterList, pandora::ClusterVector &clusterVector, SlidingFitResultMapPair &slidingFitResultMapPair) const
Refine the cluster endpoints and merge together the associated clusters alongside any extrapolated hi...
float m_minSeparationDistance
The threshold separation distance between associated clusters.
unsigned int m_maxLoopIterations
The maximum number of main loop iterations.
float m_boundaryTolerance
The maximum allowed distance of an extremal extrapolate hit to a cluster merge point.
const pandora::Cluster * GetDownstreamCluster() const
Returns the address of the downstream cluster.
void GetHitsInBoundingBox(const pandora::CartesianVector &firstCorner, const pandora::CartesianVector &secondCorner, const pandora::ClusterList *const pClusterList, ClusterToCaloHitListMap &clusterToCaloHitListMap, const pandora::ClusterList &unavailableProtectedClusters=pandora::ClusterList(), const float distanceToLine=-1.f) const
Find the unprotected hits that are contained within a defined box with the option to apply a cut on t...
std::vector< art::Ptr< recob::Cluster > > ClusterVector
bool IsNearBoundary(const pandora::CaloHit *const pCaloHit, const pandora::CartesianVector &boundaryPosition2D, const float boundaryTolerance) const
Check whether a hit is close to a boundary point.
const pandora::Cluster * RemoveOffAxisHitsFromTrack(const pandora::Cluster *const pCluster, const pandora::CartesianVector &splitPosition, const bool isEndUpstream, const ClusterToCaloHitListMap &clusterToCaloHitListMap, pandora::ClusterList &remnantClusterList, TwoDSlidingFitResultMap &microSlidingFitResultMap, TwoDSlidingFitResultMap &macroSlidingFitResultMap) const
Remove any hits in the upstream/downstream cluster that lie off of the main track axis (i...
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
void GetUnavailableProtectedClusters(const ClusterPairAssociation &clusterAssociation, const pandora::ClusterList &createdMainTrackClusters, pandora::ClusterList &unavailableProtectedClusters) const
Obtain a list of clusters whos hits are protected and cannot be reassigned.
float m_distanceToLine
The threshold hit distance of an extrapolated hit from the segment connecting line.
float m_minDirectionDeviationCosAngle
The threshold cos opening angle of the associated cluster directions.