LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
DeltaRayMergeTool.cc
Go to the documentation of this file.
1 
9 #include "Pandora/AlgorithmHeaders.h"
10 
14 
16 
17 using namespace pandora;
18 
19 namespace lar_content
20 {
21 
22 DeltaRayMergeTool::DeltaRayMergeTool() :
23  m_maxDRSeparationFromTrack(1.5f),
24  m_maxClusterSeparation(3.f),
25  m_maxVertexSeparation(10.f),
26  m_maxGoodMatchReducedChiSquared(1.f)
27 {
28 }
29 
30 //------------------------------------------------------------------------------------------------------------------------------------------
31 
33 {
34  m_pParentAlgorithm = pAlgorithm;
35 
36  if (PandoraContentApi::GetSettings(*m_pParentAlgorithm)->ShouldDisplayAlgorithmInfo())
37  std::cout << "----> Running Algorithm Tool: " << this->GetInstanceName() << ", " << this->GetType() << std::endl;
38 
39  return this->ExamineConnectedElements(overlapTensor);
40 }
41 
42 //------------------------------------------------------------------------------------------------------------------------------------------
43 
45 {
46  bool mergeMade(false), mergesMade(false), finishedTwoViewMerges(false);
47 
48  do
49  {
50  mergeMade = false;
51 
52  ClusterVector sortedKeyClusters;
53  overlapTensor.GetSortedKeyClusters(sortedKeyClusters);
54 
55  ClusterSet usedKeyClusters;
56  for (const Cluster *const pKeyCluster : sortedKeyClusters)
57  {
58  if (usedKeyClusters.count(pKeyCluster))
59  continue;
60 
61  TensorType::ElementList elementList;
62  overlapTensor.GetConnectedElements(pKeyCluster, true, elementList);
63 
64  for (const TensorType::Element &element : elementList)
65  usedKeyClusters.insert(element.GetClusterU());
66 
67  if (elementList.size() < 2)
68  continue;
69 
70  if (!finishedTwoViewMerges && this->MakeTwoCommonViewMerges(elementList))
71  {
72  mergeMade = true;
73  mergesMade = true;
74  break;
75  }
76 
77  finishedTwoViewMerges = true;
78 
79  if (this->MakeOneCommonViewMerges(elementList))
80  {
81  mergeMade = true;
82  mergesMade = true;
83  break;
84  }
85  }
86  } while (mergeMade);
87 
88  return mergesMade;
89 }
90 
91 //------------------------------------------------------------------------------------------------------------------------------------------
92 
94 {
95  for (auto iter1 = elementList.begin(); iter1 != elementList.end(); ++iter1)
96  {
97  const TensorType::Element &element1(*iter1);
98 
99  for (auto iter2 = std::next(iter1); iter2 != elementList.end(); ++iter2)
100  {
101  const TensorType::Element &element2(*iter2);
102 
103  for (const HitType hitType1 : {TPC_VIEW_U, TPC_VIEW_V})
104  {
105  if ((element1.GetCluster(hitType1) == element2.GetCluster(hitType1)))
106  {
107  for (const HitType hitType2 : {TPC_VIEW_V, TPC_VIEW_W})
108  {
109  if (hitType1 == hitType2)
110  continue;
111 
112  if ((element1.GetCluster(hitType2) == element2.GetCluster(hitType2)))
113  {
114  const HitType mergeHitType(hitType1 == TPC_VIEW_U ? (hitType2 == TPC_VIEW_V ? TPC_VIEW_W : TPC_VIEW_V) : TPC_VIEW_U);
115  const Cluster *pClusterToEnlarge(element1.GetCluster(mergeHitType)), *pClusterToDelete(element2.GetCluster(mergeHitType));
116 
117  if (this->AreAssociated(element1, element2, mergeHitType))
118  {
119  m_pParentAlgorithm->UpdateUponDeletion(pClusterToEnlarge);
120  m_pParentAlgorithm->UpdateUponDeletion(pClusterToDelete);
121 
122  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
123  PandoraContentApi::ReplaceCurrentList<Cluster>(
125 
126  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
127  PandoraContentApi::MergeAndDeleteClusters(*m_pParentAlgorithm, pClusterToEnlarge, pClusterToDelete));
128 
129  m_pParentAlgorithm->UpdateForNewClusters({pClusterToEnlarge}, {nullptr});
130 
131  return true;
132  }
133  }
134  }
135  }
136  }
137  }
138  }
139 
140  return false;
141 }
142 
143 //------------------------------------------------------------------------------------------------------------------------------------------
144 
145 bool DeltaRayMergeTool::AreAssociated(const TensorType::Element &element1, const TensorType::Element &element2, const HitType &mergeHitType) const
146 {
147  const PfoList &commonMuonPfoList1(element1.GetOverlapResult().GetCommonMuonPfoList());
148  const PfoList &commonMuonPfoList2(element2.GetOverlapResult().GetCommonMuonPfoList());
149 
150  // Demand the elements to have a shared common muon
151  PfoList commonMuonPfoList;
152  this->CombineCommonMuonPfoLists(commonMuonPfoList1, commonMuonPfoList2, commonMuonPfoList);
153 
154  if (commonMuonPfoList.empty())
155  return false;
156 
157  const Cluster *const pCluster1(element1.GetCluster(mergeHitType)), *const pCluster2(element2.GetCluster(mergeHitType));
158 
159  PfoList connectedMuonPfoList1, connectedMuonPfoList2;
160  this->GetConnectedMuons(pCluster1, commonMuonPfoList1, connectedMuonPfoList1);
161  this->GetConnectedMuons(pCluster2, commonMuonPfoList2, connectedMuonPfoList2);
162 
163  if (connectedMuonPfoList1.empty() || connectedMuonPfoList2.empty())
164  return (this->IsBrokenCluster(pCluster1, pCluster2));
165 
166  for (const ParticleFlowObject *const pConnectedMuon1 : connectedMuonPfoList1)
167  {
168  for (const ParticleFlowObject *const pConnectedMuon2 : connectedMuonPfoList2)
169  {
170  if ((pConnectedMuon1 == pConnectedMuon2) && (this->IsHiddenByTrack(pConnectedMuon1, pCluster1, pCluster2)))
171  return true;
172  }
173  }
174 
175  return false;
176 }
177 
178 //------------------------------------------------------------------------------------------------------------------------------------------
179 
180 void DeltaRayMergeTool::CombineCommonMuonPfoLists(const PfoList &commonMuonPfoList1, const PfoList &commonMuonPfoList2, PfoList &commonMuonPfoList) const
181 {
182  for (const ParticleFlowObject *const pCommonMuonPfo1 : commonMuonPfoList1)
183  {
184  for (const ParticleFlowObject *const pCommonMuonPfo2 : commonMuonPfoList2)
185  {
186  if (pCommonMuonPfo1 == pCommonMuonPfo2)
187  commonMuonPfoList.push_back(pCommonMuonPfo1);
188  }
189  }
190 }
191 
192 //------------------------------------------------------------------------------------------------------------------------------------------
193 
194 void DeltaRayMergeTool::GetConnectedMuons(const Cluster *const pDeltaRayCluster, const PfoList &commonMuonPfoList, PfoList &connectedMuonPfoList) const
195 {
196  for (const ParticleFlowObject *const pCommonMuonPfo : commonMuonPfoList)
197  {
198  if (this->IsConnected(pDeltaRayCluster, pCommonMuonPfo))
199  connectedMuonPfoList.push_back(pCommonMuonPfo);
200  }
201 }
202 
203 //------------------------------------------------------------------------------------------------------------------------------------------
204 
205 bool DeltaRayMergeTool::IsConnected(const Cluster *const pCluster, const Pfo *const pCommonMuonPfo) const
206 {
207  HitType hitType(LArClusterHelper::GetClusterHitType(pCluster));
208 
209  ClusterList muonClusterList;
210  LArPfoHelper::GetClusters(pCommonMuonPfo, hitType, muonClusterList);
211 
212  if (muonClusterList.size() != 1)
213  return false;
214 
215  const float separation(LArClusterHelper::GetClosestDistance(pCluster, muonClusterList));
216 
217  return separation < m_maxDRSeparationFromTrack;
218 }
219 
220 //------------------------------------------------------------------------------------------------------------------------------------------
221 
222 bool DeltaRayMergeTool::IsBrokenCluster(const Cluster *const pClusterToEnlarge, const Cluster *const pClusterToDelete) const
223 {
224  const float clusterSeparation(LArClusterHelper::GetClosestDistance(pClusterToEnlarge, pClusterToDelete));
225 
226  return clusterSeparation < m_maxClusterSeparation;
227 }
228 
229 //------------------------------------------------------------------------------------------------------------------------------------------
230 
231 bool DeltaRayMergeTool::IsHiddenByTrack(const ParticleFlowObject *const pMuonPfo, const Cluster *const pCluster1, const Cluster *const pCluster2) const
232 {
233  CaloHitList vertices1, vertices2;
234  this->FindVertices(pMuonPfo, pCluster1, vertices1);
235  this->FindVertices(pMuonPfo, pCluster2, vertices2);
236 
237  if (vertices1.empty() || vertices2.empty())
238  return false;
239 
240  for (const CaloHit *const pCaloHit : vertices1)
241  {
242  const float separation(LArClusterHelper::GetClosestDistance(pCaloHit->GetPositionVector(), vertices2));
243 
244  if (separation < m_maxVertexSeparation)
245  return true;
246  }
247 
248  return false;
249 }
250 
251 //------------------------------------------------------------------------------------------------------------------------------------------
252 
253 void DeltaRayMergeTool::FindVertices(const Pfo *const pCommonMuonPfo, const Cluster *const pCluster, CaloHitList &vertexList) const
254 {
255  HitType hitType(LArClusterHelper::GetClusterHitType(pCluster));
256 
257  ClusterList muonClusterList;
258  LArPfoHelper::GetClusters(pCommonMuonPfo, hitType, muonClusterList);
259 
260  if (muonClusterList.size() != 1)
261  return;
262 
263  CaloHitList caloHitList;
264  muonClusterList.front()->GetOrderedCaloHitList().FillCaloHitList(caloHitList);
265 
266  for (const CaloHit *const pCaloHit : caloHitList)
267  {
268  if (LArClusterHelper::GetClosestDistance(pCaloHit->GetPositionVector(), pCluster) < m_maxDRSeparationFromTrack)
269  vertexList.push_back(pCaloHit);
270  }
271 }
272 
273 //------------------------------------------------------------------------------------------------------------------------------------------
274 
276 {
277  for (auto iter1 = elementList.begin(); iter1 != elementList.end(); ++iter1)
278  {
279  const TensorType::Element &element1(*iter1);
280 
281  for (auto iter2 = std::next(iter1); iter2 != elementList.end(); ++iter2)
282  {
283  const TensorType::Element &element2(*iter2);
284 
285  for (const HitType hitType : {TPC_VIEW_U, TPC_VIEW_V, TPC_VIEW_W})
286  {
287  if (element1.GetCluster(hitType) == element2.GetCluster(hitType))
288  {
289  const HitType mergeHitType1(hitType == TPC_VIEW_U ? TPC_VIEW_V : hitType == TPC_VIEW_V ? TPC_VIEW_W : TPC_VIEW_U);
290  const HitType mergeHitType2(mergeHitType1 == TPC_VIEW_U ? TPC_VIEW_V
291  : mergeHitType1 == TPC_VIEW_V ? TPC_VIEW_W
292  : TPC_VIEW_U);
293 
294  const Cluster *pClusterToEnlarge1 = element1.GetCluster(mergeHitType1), *pClusterToDelete1 = element2.GetCluster(mergeHitType1);
295  const Cluster *pClusterToEnlarge2 = element1.GetCluster(mergeHitType2), *pClusterToDelete2 = element2.GetCluster(mergeHitType2);
296 
297  if ((pClusterToEnlarge1 == pClusterToDelete1) || (pClusterToEnlarge2 == pClusterToDelete2))
298  continue;
299 
300  if (!this->AreAssociated(element1, element2, mergeHitType1))
301  continue;
302 
303  if (!this->AreAssociated(element1, element2, mergeHitType2))
304  continue;
305 
306  CaloHitList caloHitList1, caloHitList2, caloHitList3;
307  pClusterToEnlarge1->GetOrderedCaloHitList().FillCaloHitList(caloHitList1);
308  pClusterToDelete1->GetOrderedCaloHitList().FillCaloHitList(caloHitList1);
309  pClusterToEnlarge2->GetOrderedCaloHitList().FillCaloHitList(caloHitList2);
310  pClusterToDelete2->GetOrderedCaloHitList().FillCaloHitList(caloHitList2);
311  element1.GetCluster(hitType)->GetOrderedCaloHitList().FillCaloHitList(caloHitList3);
312 
313  float reducedChiSquared(std::numeric_limits<float>::max());
314  StatusCode status(m_pParentAlgorithm->PerformThreeViewMatching(caloHitList1, caloHitList2, caloHitList3, reducedChiSquared));
315 
316  if (status == STATUS_CODE_NOT_FOUND)
317  continue;
318 
319  if (reducedChiSquared < m_maxGoodMatchReducedChiSquared)
320  {
321  m_pParentAlgorithm->UpdateUponDeletion(pClusterToEnlarge1);
322  m_pParentAlgorithm->UpdateUponDeletion(pClusterToDelete1);
323 
324  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
325  PandoraContentApi::ReplaceCurrentList<Cluster>(*m_pParentAlgorithm, m_pParentAlgorithm->GetClusterListName(mergeHitType1)));
326 
327  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
328  PandoraContentApi::MergeAndDeleteClusters(*m_pParentAlgorithm, pClusterToEnlarge1, pClusterToDelete1));
329 
330  m_pParentAlgorithm->UpdateForNewClusters({pClusterToEnlarge1}, {nullptr});
331 
332  m_pParentAlgorithm->UpdateUponDeletion(pClusterToEnlarge2);
333  m_pParentAlgorithm->UpdateUponDeletion(pClusterToDelete2);
334 
335  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
336  PandoraContentApi::ReplaceCurrentList<Cluster>(*m_pParentAlgorithm, m_pParentAlgorithm->GetClusterListName(mergeHitType2)));
337 
338  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=,
339  PandoraContentApi::MergeAndDeleteClusters(*m_pParentAlgorithm, pClusterToEnlarge2, pClusterToDelete2));
340 
341  m_pParentAlgorithm->UpdateForNewClusters({pClusterToEnlarge2}, {nullptr});
342 
343  return true;
344  }
345  }
346  }
347  }
348  }
349 
350  return false;
351 }
352 
353 //------------------------------------------------------------------------------------------------------------------------------------------
354 
355 StatusCode DeltaRayMergeTool::ReadSettings(const TiXmlHandle xmlHandle)
356 {
357  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
358  XmlHelper::ReadValue(xmlHandle, "MaxDRSeparationFromTrack", m_maxDRSeparationFromTrack));
359 
360  PANDORA_RETURN_RESULT_IF_AND_IF(
361  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxClusterSeparation", m_maxClusterSeparation));
362 
363  PANDORA_RETURN_RESULT_IF_AND_IF(
364  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxVertexSeparation", m_maxVertexSeparation));
365 
366  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
367  XmlHelper::ReadValue(xmlHandle, "MaxGoodMatchReducedChiSquared", m_maxGoodMatchReducedChiSquared));
368 
369  return STATUS_CODE_SUCCESS;
370 }
371 
372 } // namespace lar_content
float m_maxVertexSeparation
The maximum separation of the connection points of two delta ray clusters that are hidden by a CR tra...
float m_maxClusterSeparation
The maximum separation of two broken clusters that should be merged.
Header file for the pfo helper class.
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.
void GetConnectedElements(const pandora::Cluster *const pCluster, const bool ignoreUnavailable, ElementList &elementList) const
Get a list of elements connected to a specified cluster.
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
Header file for the delta ray merge tool 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...
static pandora::HitType GetClusterHitType(const pandora::Cluster *const pCluster)
Get the hit type associated with a two dimensional cluster.
ThreeViewDeltaRayMatchingAlgorithm * m_pParentAlgorithm
Address of the parent matching algorithm.
void GetConnectedMuons(const pandora::Cluster *const pDeltaRayCluster, const pandora::PfoList &commonMuonPfoList, pandora::PfoList &connectedMuonPfoList) const
Return the list of muon pfos that a specified delta ray cluster is directly connected to...
void UpdateUponDeletion(const pandora::Cluster *const pDeletedCluster)
Update to reflect cluster deletion.
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...
bool IsHiddenByTrack(const pandora::ParticleFlowObject *const pMuonPfo, const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2) const
Determine whether two delta ray clusters are actually a single cluster that is hidden behind a cosmic...
bool MakeOneCommonViewMerges(const TensorType::ElementList &elementList) const
Search for two matches with a single common cluster and attempt to merge the clusters in the other tw...
Header file for the cluster helper class.
Header file for the muon leading helper class.
void CombineCommonMuonPfoLists(const pandora::PfoList &commonMuonPfoList1, const pandora::PfoList &commonMuonPfoList2, pandora::PfoList &commonMuonPfoList) const
Create a list of the shared common muon pfos of two elements.
void FindVertices(const pandora::Pfo *const pCommonMuonPfo, const pandora::Cluster *const pCluster, pandora::CaloHitList &vertexList) const
Find all connection points of a delta ray cluster and a cosmic ray pfo.
void GetSortedKeyClusters(pandora::ClusterVector &sortedKeyClusters) const
Get a sorted vector of key clusters (U clusters with current implementation)
bool Run(ThreeViewDeltaRayMatchingAlgorithm *const pAlgorithm, TensorType &overlapTensor)
Run the algorithm tool.
float m_maxGoodMatchReducedChiSquared
The threshold reduced chi squared value for a potential two view merge to go ahead.
bool ExamineConnectedElements(TensorType &overlapTensor) const
Identify ambiguous matches (e.g. 3:2:1) and attempt to merge clusters together.
HitType
Definition: HitType.h:12
bool IsConnected(const pandora::Cluster *const pCluster, const pandora::Pfo *const pCommonMuonPfo) const
Determine whether a given cluster is connected to a cosmic ray pfo.
bool IsBrokenCluster(const pandora::Cluster *const pClusterToEnlarge, const pandora::Cluster *const pClusterToDelete) const
Determine whether two delta ray clusters have been split.
const std::string & GetClusterListName(const pandora::HitType hitType) const
Get the cluster list name corresponding to a specified hit type.
float m_maxDRSeparationFromTrack
The maximum distance of a connected delta ray from a cosmic ray track.
std::vector< art::Ptr< recob::Cluster > > ClusterVector
bool MakeTwoCommonViewMerges(const TensorType::ElementList &elementList) const
Search for two matches with two common clusters and attempt to merge the clusters in the third view t...
bool AreAssociated(const TensorType::Element &element1, const TensorType::Element &element2, const pandora::HitType &mergeHitType) const
Determine, from a topological point of view, whether two delta ray clusters should be merged together...
static float GetClosestDistance(const pandora::ClusterList &clusterList1, const pandora::ClusterList &clusterList2)
Get closest distance between clusters in a pair of cluster lists.