LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
CandidateVertexCreationAlgorithm.cc
Go to the documentation of this file.
1 
9 #include "Pandora/AlgorithmHeaders.h"
10 
13 
15 
16 #include <utility>
17 
18 using namespace pandora;
19 
20 namespace lar_content
21 {
22 
23 CandidateVertexCreationAlgorithm::CandidateVertexCreationAlgorithm() :
24  m_replaceCurrentVertexList(true),
25  m_slidingFitWindow(20),
26  m_minClusterCaloHits(5),
27  m_minClusterLengthSquared(3.f * 3.f),
28  m_chiSquaredCut(2.f),
29  m_enableEndpointCandidates(true),
30  m_maxEndpointXDiscrepancy(4.f),
31  m_enableCrossingCandidates(false),
32  m_nMaxCrossingCandidates(500),
33  m_maxCrossingXDiscrepancy(0.5f),
34  m_extrapolationNSteps(200),
35  m_extrapolationStepSize(0.1f),
36  m_maxCrossingSeparationSquared(2.f * 2.f),
37  m_minNearbyCrossingDistanceSquared(0.5f * 0.5f),
38  m_reducedCandidates(false),
39  m_selectionCutFactorMax(2.f),
40  m_nClustersPassingMaxCutsPar(26.f)
41 {
42 }
43 
44 //------------------------------------------------------------------------------------------------------------------------------------------
45 
47 {
48  try
49  {
50  ClusterVector clusterVectorU, clusterVectorV, clusterVectorW;
51  this->SelectClusters(clusterVectorU, clusterVectorV, clusterVectorW);
52 
53  const VertexList *pVertexList(NULL);
54  std::string temporaryListName;
55  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::CreateTemporaryListAndSetCurrent(*this, pVertexList, temporaryListName));
56 
58  {
59  this->CreateEndpointCandidates(clusterVectorU, clusterVectorV);
60  this->CreateEndpointCandidates(clusterVectorU, clusterVectorW);
61  this->CreateEndpointCandidates(clusterVectorV, clusterVectorW);
62  }
63 
65  this->CreateCrossingCandidates(clusterVectorU, clusterVectorV, clusterVectorW);
66 
67  if (!m_inputVertexListName.empty())
68  this->AddInputVertices();
69 
70  if (!pVertexList->empty())
71  {
72  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::SaveList<Vertex>(*this, m_outputVertexListName));
73 
75  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::ReplaceCurrentList<Vertex>(*this, m_outputVertexListName));
76  }
77  }
78  catch (StatusCodeException &statusCodeException)
79  {
80  this->TidyUp();
81  throw statusCodeException;
82  }
83 
84  this->TidyUp();
85 
86  return STATUS_CODE_SUCCESS;
87 }
88 
89 //------------------------------------------------------------------------------------------------------------------------------------------
90 
91 void CandidateVertexCreationAlgorithm::SelectClusters(ClusterVector &clusterVectorU, ClusterVector &clusterVectorV, ClusterVector &clusterVectorW)
92 {
93  for (const std::string &clusterListName : m_inputClusterListNames)
94  {
95  const ClusterList *pClusterList(NULL);
96  PANDORA_THROW_RESULT_IF_AND_IF(
97  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_INITIALIZED, !=, PandoraContentApi::GetList(*this, clusterListName, pClusterList));
98 
99  if (!pClusterList || pClusterList->empty())
100  {
101  if (PandoraContentApi::GetSettings(*this)->ShouldDisplayAlgorithmInfo())
102  std::cout << "CandidateVertexCreationAlgorithm: unable to find cluster list " << clusterListName << std::endl;
103 
104  continue;
105  }
106 
107  const HitType hitType(LArClusterHelper::GetClusterHitType(*(pClusterList->begin())));
108 
109  if ((TPC_VIEW_U != hitType) && (TPC_VIEW_V != hitType) && (TPC_VIEW_W != hitType))
110  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
111 
112  ClusterVector &selectedClusterVector((TPC_VIEW_U == hitType) ? clusterVectorU
113  : (TPC_VIEW_V == hitType) ? clusterVectorV
114  : clusterVectorW);
115 
116  if (!selectedClusterVector.empty())
117  throw StatusCodeException(STATUS_CODE_FAILURE);
118 
119  ClusterVector sortedClusters(pClusterList->begin(), pClusterList->end());
120  std::sort(sortedClusters.begin(), sortedClusters.end(), LArClusterHelper::SortByNHits);
121 
122  unsigned int nClustersPassingMaxCuts(0);
124  {
125  for (const Cluster *const pCluster : sortedClusters)
126  {
127  float selectionCutFactor(1.f);
128 
129  if (pCluster->GetParticleId() == E_MINUS)
130  selectionCutFactor = m_selectionCutFactorMax;
131 
132  if (pCluster->GetNCaloHits() < m_minClusterCaloHits * selectionCutFactor)
133  continue;
134 
135  if (LArClusterHelper::GetLengthSquared(pCluster) < m_minClusterLengthSquared * selectionCutFactor * selectionCutFactor)
136  continue;
137 
138  nClustersPassingMaxCuts++;
139  }
140  }
141 
142  for (const Cluster *const pCluster : sortedClusters)
143  {
144  float selectionCutFactor(1.f);
145 
146  if (pCluster->GetParticleId() == E_MINUS && m_reducedCandidates)
147  {
148  selectionCutFactor = (m_selectionCutFactorMax + 1.f) * 0.5f +
149  (m_selectionCutFactorMax - 1.f) * 0.5f * std::tanh(static_cast<float>(nClustersPassingMaxCuts) - m_nClustersPassingMaxCutsPar);
150  }
151 
152  if (pCluster->GetNCaloHits() < m_minClusterCaloHits * selectionCutFactor)
153  continue;
154 
155  if (LArClusterHelper::GetLengthSquared(pCluster) < m_minClusterLengthSquared * selectionCutFactor * selectionCutFactor)
156  continue;
157 
158  try
159  {
160  this->AddToSlidingFitCache(pCluster);
161  selectedClusterVector.push_back(pCluster);
162  }
163  catch (StatusCodeException &statusCodeException)
164  {
165  if (STATUS_CODE_FAILURE == statusCodeException.GetStatusCode())
166  throw statusCodeException;
167  }
168  }
169  }
170 }
171 
172 //------------------------------------------------------------------------------------------------------------------------------------------
173 
174 void CandidateVertexCreationAlgorithm::CreateEndpointCandidates(const ClusterVector &clusterVector1, const ClusterVector &clusterVector2) const
175 {
176  for (const Cluster *const pCluster1 : clusterVector1)
177  {
178  const HitType hitType1(LArClusterHelper::GetClusterHitType(pCluster1));
179 
180  const TwoDSlidingFitResult &fitResult1(this->GetCachedSlidingFitResult(pCluster1));
181  const CartesianVector minLayerPosition1(fitResult1.GetGlobalMinLayerPosition());
182  const CartesianVector maxLayerPosition1(fitResult1.GetGlobalMaxLayerPosition());
183 
184  for (const Cluster *const pCluster2 : clusterVector2)
185  {
186  const HitType hitType2(LArClusterHelper::GetClusterHitType(pCluster2));
187 
188  const TwoDSlidingFitResult &fitResult2(this->GetCachedSlidingFitResult(pCluster2));
189  const CartesianVector minLayerPosition2(fitResult2.GetGlobalMinLayerPosition());
190  const CartesianVector maxLayerPosition2(fitResult2.GetGlobalMaxLayerPosition());
191 
192  this->CreateEndpointVertex(maxLayerPosition1, hitType1, fitResult2);
193  this->CreateEndpointVertex(minLayerPosition1, hitType1, fitResult2);
194  this->CreateEndpointVertex(maxLayerPosition2, hitType2, fitResult1);
195  this->CreateEndpointVertex(minLayerPosition2, hitType2, fitResult1);
196  }
197  }
198 }
199 
200 //------------------------------------------------------------------------------------------------------------------------------------------
201 
203  const CartesianVector &position1, const HitType hitType1, const TwoDSlidingFitResult &fitResult2) const
204 {
205  const CartesianVector minLayerPosition2(fitResult2.GetGlobalMinLayerPosition());
206  const CartesianVector maxLayerPosition2(fitResult2.GetGlobalMaxLayerPosition());
207 
208  if ((((position1.GetX() < minLayerPosition2.GetX()) && (position1.GetX() < maxLayerPosition2.GetX())) ||
209  ((position1.GetX() > minLayerPosition2.GetX()) && (position1.GetX() > maxLayerPosition2.GetX()))) &&
210  (std::fabs(position1.GetX() - minLayerPosition2.GetX()) > m_maxEndpointXDiscrepancy) &&
211  (std::fabs(position1.GetX() - maxLayerPosition2.GetX()) > m_maxEndpointXDiscrepancy))
212  {
213  return;
214  }
215 
216  CartesianVector position2(0.f, 0.f, 0.f);
217  if (STATUS_CODE_SUCCESS != fitResult2.GetExtrapolatedPositionAtX(position1.GetX(), position2))
218  return;
219 
220  const HitType hitType2(LArClusterHelper::GetClusterHitType(fitResult2.GetCluster()));
221 
222  float chiSquared(0.f);
223  CartesianVector position3D(0.f, 0.f, 0.f);
224  LArGeometryHelper::MergeTwoPositions3D(this->GetPandora(), hitType1, hitType2, position1, position2, position3D, chiSquared);
225 
226  if (chiSquared > m_chiSquaredCut)
227  return;
228 
229  PandoraContentApi::Vertex::Parameters parameters;
230  parameters.m_position = position3D;
231  parameters.m_vertexLabel = VERTEX_INTERACTION;
232  parameters.m_vertexType = VERTEX_3D;
233 
234  const Vertex *pVertex(NULL);
235  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::Vertex::Create(*this, parameters, pVertex));
236 }
237 
238 //------------------------------------------------------------------------------------------------------------------------------------------
239 
241  const ClusterVector &clusterVectorU, const ClusterVector &clusterVectorV, const ClusterVector &clusterVectorW) const
242 {
243  CartesianPointVector crossingsU, crossingsV, crossingsW;
244  this->FindCrossingPoints(clusterVectorU, crossingsU);
245  this->FindCrossingPoints(clusterVectorV, crossingsV);
246  this->FindCrossingPoints(clusterVectorW, crossingsW);
247 
248  unsigned int nCrossingCandidates(0);
249  this->CreateCrossingVertices(crossingsU, crossingsV, TPC_VIEW_U, TPC_VIEW_V, nCrossingCandidates);
250  this->CreateCrossingVertices(crossingsU, crossingsW, TPC_VIEW_U, TPC_VIEW_W, nCrossingCandidates);
251  this->CreateCrossingVertices(crossingsV, crossingsW, TPC_VIEW_V, TPC_VIEW_W, nCrossingCandidates);
252 }
253 
254 //------------------------------------------------------------------------------------------------------------------------------------------
255 
256 void CandidateVertexCreationAlgorithm::FindCrossingPoints(const ClusterVector &clusterVector, CartesianPointVector &crossingPoints) const
257 {
258  ClusterToSpacepointsMap clusterToSpacepointsMap;
259 
260  for (const Cluster *const pCluster : clusterVector)
261  {
262  ClusterToSpacepointsMap::iterator mapIter(clusterToSpacepointsMap.emplace(pCluster, CartesianPointVector()).first);
263  this->GetSpacepoints(pCluster, mapIter->second);
264  }
265 
266  for (const Cluster *const pCluster1 : clusterVector)
267  {
268  for (const Cluster *const pCluster2 : clusterVector)
269  {
270  if (pCluster1 == pCluster2)
271  continue;
272 
273  this->FindCrossingPoints(clusterToSpacepointsMap.at(pCluster1), clusterToSpacepointsMap.at(pCluster2), crossingPoints);
274  }
275  }
276 }
277 
278 //------------------------------------------------------------------------------------------------------------------------------------------
279 
280 void CandidateVertexCreationAlgorithm::GetSpacepoints(const Cluster *const pCluster, CartesianPointVector &spacepoints) const
281 {
282  LArClusterHelper::GetCoordinateVector(pCluster, spacepoints);
283 
284  const TwoDSlidingFitResult &fitResult(this->GetCachedSlidingFitResult(pCluster));
285  const float minLayerRL(fitResult.GetL(fitResult.GetMinLayer()));
286  const float maxLayerRL(fitResult.GetL(fitResult.GetMaxLayer()));
287 
288  for (unsigned int iStep = 0; iStep < m_extrapolationNSteps; ++iStep)
289  {
290  const float deltaRL(static_cast<float>(iStep) * m_extrapolationStepSize);
291 
292  CartesianVector positionPositive(0.f, 0.f, 0.f), positionNegative(0.f, 0.f, 0.f);
293  fitResult.GetExtrapolatedPosition(maxLayerRL + deltaRL, positionPositive);
294  fitResult.GetExtrapolatedPosition(minLayerRL - deltaRL, positionNegative);
295 
296  spacepoints.push_back(positionPositive);
297  spacepoints.push_back(positionNegative);
298  }
299 
300  std::sort(spacepoints.begin(), spacepoints.end(), LArClusterHelper::SortCoordinatesByPosition);
301 }
302 
303 //------------------------------------------------------------------------------------------------------------------------------------------
304 
306  const CartesianPointVector &spacepoints1, const CartesianPointVector &spacepoints2, CartesianPointVector &crossingPoints) const
307 {
308  bool bestCrossingFound(false);
309  float bestSeparationSquared(m_maxCrossingSeparationSquared);
310  CartesianVector bestPosition1(0.f, 0.f, 0.f), bestPosition2(0.f, 0.f, 0.f);
311 
312  for (const CartesianVector &position1 : spacepoints1)
313  {
314  for (const CartesianVector &position2 : spacepoints2)
315  {
316  const float separationSquared((position1 - position2).GetMagnitudeSquared());
317 
318  if (separationSquared < bestSeparationSquared)
319  {
320  bestCrossingFound = true;
321  bestSeparationSquared = separationSquared;
322  bestPosition1 = position1;
323  bestPosition2 = position2;
324  }
325  }
326  }
327 
328  if (bestCrossingFound)
329  {
330  bool alreadyPopulated(false);
331 
332  for (const CartesianVector &existingPosition : crossingPoints)
333  {
334  if (((existingPosition - bestPosition1).GetMagnitudeSquared() < m_minNearbyCrossingDistanceSquared) ||
335  ((existingPosition - bestPosition2).GetMagnitudeSquared() < m_minNearbyCrossingDistanceSquared))
336  {
337  alreadyPopulated = true;
338  break;
339  }
340  }
341 
342  if (!alreadyPopulated)
343  {
344  crossingPoints.push_back(bestPosition1);
345  crossingPoints.push_back(bestPosition2);
346  }
347  }
348 }
349 
350 //------------------------------------------------------------------------------------------------------------------------------------------
351 
352 void CandidateVertexCreationAlgorithm::CreateCrossingVertices(const CartesianPointVector &crossingPoints1,
353  const CartesianPointVector &crossingPoints2, const HitType hitType1, const HitType hitType2, unsigned int &nCrossingCandidates) const
354 {
355 
356  for (const CartesianVector &position1 : crossingPoints1)
357  {
358  for (const CartesianVector &position2 : crossingPoints2)
359  {
360  if (nCrossingCandidates > m_nMaxCrossingCandidates)
361  return;
362 
363  if (std::fabs(position1.GetX() - position2.GetX()) > m_maxCrossingXDiscrepancy)
364  continue;
365 
366  float chiSquared(0.f);
367  CartesianVector position3D(0.f, 0.f, 0.f);
368  LArGeometryHelper::MergeTwoPositions3D(this->GetPandora(), hitType1, hitType2, position1, position2, position3D, chiSquared);
369 
370  if (chiSquared > m_chiSquaredCut)
371  continue;
372 
373  PandoraContentApi::Vertex::Parameters parameters;
374  parameters.m_position = position3D;
375  parameters.m_vertexLabel = VERTEX_INTERACTION;
376  parameters.m_vertexType = VERTEX_3D;
377 
378  const Vertex *pVertex(NULL);
379  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::Vertex::Create(*this, parameters, pVertex));
380  ++nCrossingCandidates;
381  }
382  }
383 }
384 
385 //------------------------------------------------------------------------------------------------------------------------------------------
386 
388 {
389  const VertexList *pInputVertexList{nullptr};
390  try
391  { // ATTN - No guarantee the list has been initialised, but silent failure here is ok
392  PandoraContentApi::GetList(*this, m_inputVertexListName, pInputVertexList);
393  if (!pInputVertexList)
394  return;
395 
396  for (const Vertex *pInputVertex : *pInputVertexList)
397  {
398  PandoraContentApi::Vertex::Parameters parameters;
399  parameters.m_position = pInputVertex->GetPosition();
400  parameters.m_vertexLabel = VERTEX_INTERACTION;
401  parameters.m_vertexType = VERTEX_3D;
402 
403  const Vertex *pVertex(nullptr);
404  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::Vertex::Create(*this, parameters, pVertex));
405  }
406  }
407  catch (const StatusCodeException &)
408  {
409  return;
410  }
411 }
412 
413 //------------------------------------------------------------------------------------------------------------------------------------------
414 
416 {
417  const float slidingFitPitch(LArGeometryHelper::GetWireZPitch(this->GetPandora()));
418  const TwoDSlidingFitResult slidingFitResult(pCluster, m_slidingFitWindow, slidingFitPitch);
419 
420  if (!m_slidingFitResultMap.insert(TwoDSlidingFitResultMap::value_type(pCluster, slidingFitResult)).second)
421  throw StatusCodeException(STATUS_CODE_FAILURE);
422 }
423 
424 //------------------------------------------------------------------------------------------------------------------------------------------
425 
427 {
429 
430  if (m_slidingFitResultMap.end() == iter)
431  throw StatusCodeException(STATUS_CODE_NOT_FOUND);
432 
433  return iter->second;
434 }
435 
436 //------------------------------------------------------------------------------------------------------------------------------------------
437 
439 {
440  m_slidingFitResultMap.clear();
441 }
442 
443 //------------------------------------------------------------------------------------------------------------------------------------------
444 
445 StatusCode CandidateVertexCreationAlgorithm::ReadSettings(const TiXmlHandle xmlHandle)
446 {
447  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadVectorOfValues(xmlHandle, "InputClusterListNames", m_inputClusterListNames));
448 
449  PANDORA_RETURN_RESULT_IF_AND_IF(
450  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "InputVertexListName", m_inputVertexListName));
451 
452  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "OutputVertexListName", m_outputVertexListName));
453 
454  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
455  XmlHelper::ReadValue(xmlHandle, "ReplaceCurrentVertexList", m_replaceCurrentVertexList));
456 
457  PANDORA_RETURN_RESULT_IF_AND_IF(
458  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "SlidingFitWindow", m_slidingFitWindow));
459 
460  PANDORA_RETURN_RESULT_IF_AND_IF(
461  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinClusterCaloHits", m_minClusterCaloHits));
462 
463  float minClusterLength = std::sqrt(m_minClusterLengthSquared);
464  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MinClusterLength", minClusterLength));
465  m_minClusterLengthSquared = minClusterLength * minClusterLength;
466 
467  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ChiSquaredCut", m_chiSquaredCut));
468 
469  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
470  XmlHelper::ReadValue(xmlHandle, "EnableEndpointCandidates", m_enableEndpointCandidates));
471 
472  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
473  XmlHelper::ReadValue(xmlHandle, "MaxEndpointXDiscrepancy", m_maxEndpointXDiscrepancy));
474 
475  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
476  XmlHelper::ReadValue(xmlHandle, "EnableCrossingCandidates", m_enableCrossingCandidates));
477 
478  PANDORA_RETURN_RESULT_IF_AND_IF(
479  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "NMaxCrossingCandidates", m_nMaxCrossingCandidates));
480 
481  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
482  XmlHelper::ReadValue(xmlHandle, "MaxCrossingXDiscrepancy", m_maxCrossingXDiscrepancy));
483 
484  PANDORA_RETURN_RESULT_IF_AND_IF(
485  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ExtrapolationNSteps", m_extrapolationNSteps));
486 
487  PANDORA_RETURN_RESULT_IF_AND_IF(
488  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ExtrapolationStepSize", m_extrapolationStepSize));
489 
490  PANDORA_RETURN_RESULT_IF_AND_IF(
491  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "ReducedCandidates", m_reducedCandidates));
492 
493  PANDORA_RETURN_RESULT_IF_AND_IF(
494  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "SelectionCutFactorMax", m_selectionCutFactorMax));
495 
496  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
497  XmlHelper::ReadValue(xmlHandle, "NClustersPassingMaxCutsPar", m_nClustersPassingMaxCutsPar));
498 
499  float maxCrossingSeparation = std::sqrt(m_maxCrossingSeparationSquared);
500  PANDORA_RETURN_RESULT_IF_AND_IF(
501  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxCrossingSeparation", maxCrossingSeparation));
502  m_maxCrossingSeparationSquared = maxCrossingSeparation * maxCrossingSeparation;
503 
504  float minNearbyCrossingDistance = std::sqrt(m_minNearbyCrossingDistanceSquared);
505  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
506  XmlHelper::ReadValue(xmlHandle, "MinNearbyCrossingDistance", minNearbyCrossingDistance));
507  m_minNearbyCrossingDistanceSquared = minNearbyCrossingDistance * minNearbyCrossingDistance;
508 
509  return STATUS_CODE_SUCCESS;
510 }
511 
512 } // namespace lar_content
bool m_reducedCandidates
Whether to reduce the number of candidates.
intermediate_table::iterator iterator
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.
float m_minClusterLengthSquared
The min length (squared) in base cluster selection method.
void CreateEndpointVertex(const pandora::CartesianVector &position1, const pandora::HitType hitType1, const TwoDSlidingFitResult &fitResult2) const
Create a candidate vertex position, using an end-point position from one cluster and sliding fit to a...
void CreateCrossingCandidates(const pandora::ClusterVector &clusterVectorU, const pandora::ClusterVector &clusterVectorV, const pandora::ClusterVector &clusterVectorW) const
Extrapolate 2D clusters, find where they cross, and match crossing points between views to create ver...
static void MergeTwoPositions3D(const pandora::Pandora &pandora, const pandora::HitType view1, const pandora::HitType view2, const pandora::CartesianVector &position1, const pandora::CartesianVector &position2, pandora::CartesianVector &position3D, float &chiSquared)
Merge 2D positions from two views to give unified 3D position.
void GetSpacepoints(const pandora::Cluster *const pCluster, pandora::CartesianPointVector &spacePoints) const
Get a list of spacepoints representing cluster 2D hit positions and extrapolated positions.
float m_extrapolationStepSize
The extrapolation step size in cm.
float m_maxCrossingSeparationSquared
The separation (squared) between spacepoints below which a crossing can be identified.
pandora::StatusCode GetExtrapolatedPositionAtX(const float x, pandora::CartesianVector &position) const
Get extrapolated position (beyond span) for a given input x coordinate.
bool m_enableEndpointCandidates
Whether to create endpoint-based candidates.
void AddToSlidingFitCache(const pandora::Cluster *const pCluster)
Creates a 2D sliding fit of a cluster and stores it for later use.
intermediate_table::const_iterator const_iterator
bool m_replaceCurrentVertexList
Whether to replace the current vertex list with the output list.
float m_maxEndpointXDiscrepancy
The max cluster endpoint discrepancy.
static pandora::HitType GetClusterHitType(const pandora::Cluster *const pCluster)
Get the hit type associated with a two dimensional cluster.
void FindCrossingPoints(const pandora::ClusterVector &clusterVector, pandora::CartesianPointVector &crossingPoints) const
Identify where (extrapolated) clusters plausibly cross in 2D.
static float GetWireZPitch(const pandora::Pandora &pandora, const float maxWirePitchDiscrepancy=0.01)
Return the wire pitch.
TFile f
Definition: plotHisto.C:6
static bool SortCoordinatesByPosition(const pandora::CartesianVector &lhs, const pandora::CartesianVector &rhs)
Sort cartesian vectors by their position (use Z, followed by X, followed by Y)
Header file for the geometry helper class.
Header file for the candidate vertex creation algorithm class.
int GetMaxLayer() const
Get the maximum occupied layer in the sliding fit.
int GetMinLayer() const
Get the minimum occupied layer in the sliding fit.
std::string m_inputVertexListName
The list name for existing candidate vertices.
Header file for the cluster helper class.
bool m_enableCrossingCandidates
Whether to create crossing vertex candidates.
float m_maxCrossingXDiscrepancy
The max cluster endpoint discrepancy.
TwoDSlidingFitResultMap m_slidingFitResultMap
The sliding fit result map.
void CreateEndpointCandidates(const pandora::ClusterVector &clusterVector1, const pandora::ClusterVector &clusterVector2) const
Create candidate vertex positions by comparing pairs of cluster end positions.
float m_nClustersPassingMaxCutsPar
Parameter for number of clusters passing the max base cluster selection cuts.
pandora::StatusCode GetExtrapolatedPosition(const float rL, pandora::CartesianVector &position) const
Get extrapolated position (beyond span) for a given input coordinate.
float m_minNearbyCrossingDistanceSquared
The minimum allowed distance between identified crossing positions.
unsigned int m_minClusterCaloHits
The min number of hits in base cluster selection method.
void TidyUp()
Clear relevant algorithm member variables between events.
pandora::CartesianVector GetGlobalMinLayerPosition() const
Get global position corresponding to the fit result in minimum fit layer.
unsigned int m_extrapolationNSteps
Number of extrapolation steps, at each end of cluster, of specified size.
std::string m_outputVertexListName
The name under which to save the output vertex list.
const TwoDSlidingFitResult & GetCachedSlidingFitResult(const pandora::Cluster *const pCluster) const
Get a sliding fit result from the algorithm cache.
pandora::StringVector m_inputClusterListNames
The list of cluster list names.
unsigned int m_nMaxCrossingCandidates
The max number of crossing candidates to create.
const pandora::Cluster * GetCluster() const
Get the address of the cluster, if originally provided.
std::unordered_map< const pandora::Cluster *, pandora::CartesianPointVector > ClusterToSpacepointsMap
float m_selectionCutFactorMax
Maximum factor to multiply the base cluster selection cuts.
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
HitType
Definition: HitType.h:12
void SelectClusters(pandora::ClusterVector &clusterVectorU, pandora::ClusterVector &clusterVectorV, pandora::ClusterVector &clusterVectorW)
Select a subset of input clusters (contained in the input list names) for processing in this algorith...
float m_chiSquaredCut
The chi squared cut (accept only 3D vertex positions with values below cut)
float GetL(const int layer) const
Get longitudinal coordinate for a given sliding linear fit layer number.
static float GetLengthSquared(const pandora::Cluster *const pCluster)
Get length squared of cluster.
static void GetCoordinateVector(const pandora::Cluster *const pCluster, pandora::CartesianPointVector &coordinateVector)
Get vector of hit coordinates from an input cluster.
boost::graph_traits< ModuleGraph >::vertex_descriptor Vertex
Definition: ModuleGraph.h:25
std::vector< art::Ptr< recob::Cluster > > ClusterVector
void AddInputVertices() const
Add candidate vertices from any input vertices.
void CreateCrossingVertices(const pandora::CartesianPointVector &crossingPoints1, const pandora::CartesianPointVector &crossingPoints2, const pandora::HitType hitType1, const pandora::HitType hitType2, unsigned int &nCrossingCandidates) const
Attempt to create candidate vertex positions, using 2D crossing points in 2 views.
std::list< Vertex > VertexList
Definition: DCEL.h:169
pandora::CartesianVector GetGlobalMaxLayerPosition() const
Get global position corresponding to the fit result in maximum fit layer.
unsigned int m_slidingFitWindow
The layer window for the sliding linear fits.