28 m_minClusterCaloHits(20),
29 m_minClusterLengthSquared(5.f * 5.f),
30 m_minClusterXSpan(0.25f),
31 m_vertexClusterMode(false),
32 m_minVertexLongitudinalDistance(-2.5f),
33 m_maxVertexTransverseDistance(1.5f),
34 m_minXOverlapFraction(0.75f),
35 m_minXOverlapFractionGaps(0.75f),
36 m_sampleStepSize(0.5f),
37 m_slidingFitHalfWindow(10),
46 ClusterList inputClusterListU, inputClusterListV, inputClusterListW;
47 this->
GetInputClusters(inputClusterListU, inputClusterListV, inputClusterListW);
49 ClusterList selectedClusterListU, selectedClusterListV, selectedClusterListW;
55 this->
FindOverlaps(selectedClusterListU, selectedClusterListV, overlapTensor);
56 this->
FindOverlaps(selectedClusterListV, selectedClusterListW, overlapTensor);
57 this->
FindOverlaps(selectedClusterListW, selectedClusterListU, overlapTensor);
60 return STATUS_CODE_SUCCESS;
72 if (!pClusterList || pClusterList->empty())
75 std::cout <<
"ParticleRecoveryAlgorithm: unable to find cluster list " << *iter << std::endl;
80 for (ClusterList::const_iterator cIter = pClusterList->begin(), cIterEnd = pClusterList->end(); cIter != cIterEnd; ++cIter)
82 const Cluster *
const pCluster(*cIter);
89 clusterList.push_back(pCluster);
114 for (ClusterList::const_iterator iter = inputClusterList.begin(), iterEnd = inputClusterList.end(); iter != iterEnd; ++iter)
116 const Cluster *
const pCluster = *iter;
127 float xMin(0.f), xMax(0.f);
133 selectedClusterList.push_back(pCluster);
143 for (ClusterList::const_iterator iter = inputClusterList.begin(), iterEnd = inputClusterList.end(); iter != iterEnd; ++iter)
147 if (!(*iter)->IsAvailable())
159 for (ClusterList::const_iterator iter = inputClusterList.begin(), iterEnd = inputClusterList.end(); iter != iterEnd; ++iter)
163 const Cluster *
const pCluster = *iter;
170 for (CartesianPointVector::const_iterator vIter = vertexList.begin(), vIterEnd = vertexList.end(); vIter != vIterEnd; ++vIter)
175 selectedClusterList.push_back(pCluster);
190 for (ClusterList::const_iterator iter1 = clusterList1.begin(), iter1End = clusterList1.end(); iter1 != iter1End; ++iter1)
192 for (ClusterList::const_iterator iter2 = clusterList2.begin(), iter2End = clusterList2.end(); iter2 != iter2End; ++iter2)
210 float xMin1(0.f), xMax1(0.f), xMin2(0.f), xMax2(0.f);
214 const float xSpan1(xMax1 - xMin1), xSpan2(xMax2 - xMin2);
216 if ((xSpan1 < std::numeric_limits<float>::epsilon()) || (xSpan2 < std::numeric_limits<float>::epsilon()))
219 const float xOverlap(std::min(xMax1, xMax2) - std::max(xMin1, xMin2));
221 float xOverlapFraction1(xOverlap / xSpan1), xOverlapFraction2(xOverlap / xSpan2);
235 const Cluster *
const pCluster2,
const float xMin2,
const float xMax2,
float &xOverlapFraction1,
float &xOverlapFraction2)
const
240 const float xMin(std::min(xMin1, xMin2));
241 const float xMax(std::max(xMax1, xMax2));
242 float xMinEff1(xMin1), xMaxEff1(xMax1), xMinEff2(xMin2), xMaxEff2(xMax2);
247 const float effectiveXSpan1(xMaxEff1 - xMinEff1), effectiveXSpan2(xMaxEff2 - xMinEff2);
248 const float effectiveXOverlapSpan(std::min(xMaxEff1, xMaxEff2) - std::max(xMinEff1, xMinEff2));
250 if ((effectiveXSpan1 > std::numeric_limits<float>::epsilon()) && (effectiveXSpan2 > std::numeric_limits<float>::epsilon()))
252 xOverlapFraction1 = std::min(1.f, (effectiveXOverlapSpan / effectiveXSpan1));
253 xOverlapFraction2 = std::min(1.f, (effectiveXOverlapSpan / effectiveXSpan2));
260 const pandora::Cluster *
const pCluster,
const float xMin,
const float xMax,
float &xMinEff,
float &xMaxEff)
const
269 const int nSamplingPointsLeft(1 +
static_cast<int>((xMinEff - xMin) /
m_sampleStepSize));
270 const int nSamplingPointsRight(1 +
static_cast<int>((xMax - xMaxEff) /
m_sampleStepSize));
271 float dxMin(0.f), dxMax(0.f);
273 for (
int iSample = 1; iSample <= nSamplingPointsLeft; ++iSample)
275 const float xSample(std::max(xMin, xMinEff -
static_cast<float>(iSample) *
m_sampleStepSize));
280 dxMin = xMinEff - xSample;
283 for (
int iSample = 1; iSample <= nSamplingPointsRight; ++iSample)
285 const float xSample(std::min(xMax, xMaxEff +
static_cast<float>(iSample) *
m_sampleStepSize));
290 dxMax = xSample - xMaxEff;
293 xMinEff = xMinEff - dxMin;
294 xMaxEff = xMaxEff + dxMax;
308 for (
const Cluster *
const pKeyCluster : sortedKeyClusters)
310 ClusterList clusterListU, clusterListV, clusterListW;
312 overlapTensor.
GetConnectedElements(pKeyCluster,
true, clusterListU, clusterListV, clusterListW);
313 const unsigned int nU(clusterListU.size()), nV(clusterListV.size()), nW(clusterListW.size());
315 if ((0 == nU * nV) && (0 == nV * nW) && (0 == nW * nU))
319 clusterListAll.insert(clusterListAll.end(), clusterListU.begin(), clusterListU.end());
320 clusterListAll.insert(clusterListAll.end(), clusterListV.begin(), clusterListV.end());
321 clusterListAll.insert(clusterListAll.end(), clusterListW.begin(), clusterListW.end());
323 if ((1 == nU * nV * nW) && this->
CheckConsistency(*(clusterListU.begin()), *(clusterListV.begin()), *(clusterListW.begin())))
327 else if ((0 == nU * nV * nW) && ((1 == nU && 1 == nV) || (1 == nV && 1 == nW) || (1 == nW && 1 == nU)))
344 float xMinU(0.f), xMinV(0.f), xMinW(0.f), xMaxU(0.f), xMaxV(0.f), xMaxW(0.f);
349 const float xMin(std::max(xMinU, std::max(xMinV, xMinW)));
350 const float xMax(std::min(xMaxU, std::min(xMaxV, xMaxW)));
351 const float xOverlap(xMax - xMin);
353 if (xOverlap < std::numeric_limits<float>::epsilon())
357 float u(std::numeric_limits<float>::max()), v(std::numeric_limits<float>::max()), w(std::numeric_limits<float>::max());
370 const float pseudoChi2(((u - vw2u) * (u - vw2u) + (v - wu2v) * (v - wu2v) + (w - uv2w) * (w - uv2w)) / 3.f);
382 if (clusterList.size() < 2)
385 const PfoList *pPfoList = NULL;
386 std::string pfoListName;
390 PandoraContentApi::ParticleFlowObject::Parameters pfoParameters;
391 pfoParameters.m_particleId = MU_MINUS;
394 pfoParameters.m_energy = 0.f;
396 pfoParameters.m_clusterList = clusterList;
399 PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::ParticleFlowObject::Create(*
this, pfoParameters, pPfo));
401 if (!pPfoList->empty())
420 if (pClusterU && pClusterV && !pClusterW)
427 else if (!pClusterU && pClusterV && pClusterW)
434 else if (pClusterU && !pClusterV && pClusterW)
462 (
TPC_VIEW_U == hitType) ? m_clusterNavigationMapUV : (
TPC_VIEW_V == hitType) ? m_clusterNavigationMapVW : m_clusterNavigationMapWU);
464 if (clusterList.end() != std::find(clusterList.begin(), clusterList.end(), pCluster))
467 clusterList.push_back(pCluster);
469 ClusterNavigationMap::const_iterator iter = navigationMap.find(pCluster);
471 if (navigationMap.end() == iter)
474 for (ClusterList::const_iterator cIter = iter->second.begin(), cIterEnd = iter->second.end(); cIter != cIterEnd; ++cIter)
475 this->GetConnectedElements(*cIter, ignoreUnavailable, clusterListU, clusterListV, clusterListW);
516 std::cout <<
"ParticleRecoveryAlgorithm: Invalid value for SampleStepSize " <<
m_sampleStepSize << std::endl;
525 return STATUS_CODE_SUCCESS;
Header file for the cluster helper class.
Header file for the geometry helper class.
Header file for the lar pointing cluster class.
Header file for the track recovery algorithm class.
#define PANDORA_THROW_RESULT_IF_AND_IF(StatusCode1, StatusCode2, Operator, Command)
#define PANDORA_THROW_RESULT_IF(StatusCode1, Operator, Command)
#define PANDORA_RETURN_RESULT_IF_AND_IF(StatusCode1, StatusCode2, Operator, Command)
#define PANDORA_RETURN_RESULT_IF(StatusCode1, Operator, Command)
static const pandora::GeometryManager * GetGeometry(const pandora::Algorithm &algorithm)
Get the pandora geometry instance.
static pandora::StatusCode ReplaceCurrentList(const pandora::Algorithm &algorithm, const std::string &newListName)
Replace the current list with a pre-saved list; use this new list as a permanent replacement for the ...
static pandora::StatusCode CreateTemporaryListAndSetCurrent(const pandora::Algorithm &algorithm, const T *&pT, std::string &temporaryListName)
Create a temporary list and set it to be the current list, enabling object creation.
static pandora::StatusCode GetList(const pandora::Algorithm &algorithm, const std::string &listName, const T *&pT)
Get a named list.
static const pandora::PandoraSettings * GetSettings(const pandora::Algorithm &algorithm)
Get the pandora settings instance.
static pandora::HitType GetClusterHitType(const pandora::Cluster *const pCluster)
Get the hit type associated with a two dimensional cluster.
static pandora::StatusCode GetAverageZ(const pandora::Cluster *const pCluster, const float xmin, const float xmax, float &averageZ)
Get average Z positions of the calo hits in a cluster in range xmin to xmax.
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,...
static float GetLengthSquared(const pandora::Cluster *const pCluster)
Get length squared of cluster.
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).
static float GetWireZPitch(const pandora::Pandora &pandora, const float maxWirePitchDiscrepancy=0.01)
Return the wire pitch.
static bool IsXSamplingPointInGap(const pandora::Pandora &pandora, const float xSample, const TwoDSlidingFitResult &slidingFitResult, const float gapTolerance=0.f)
Whether there is a gap in a cluster (described via its sliding fit result) at a specified x sampling ...
const pandora::CartesianVector & GetPosition() const
Get the vertex position.
static bool IsNode(const pandora::CartesianVector &parentVertex, const LArPointingCluster::Vertex &daughterVertex, const float minLongitudinalDistance, const float maxTransverseDistance)
Whether pointing vertex is adjacent to a given position.
LArPointingCluster class.
const Vertex & GetInnerVertex() const
Get the inner vertex.
const Vertex & GetOuterVertex() const
Get the outer vertex.
SimpleOverlapTensor class.
ClusterNavigationMap m_clusterNavigationMapUV
The cluster navigation map U->V.
void GetConnectedElements(const pandora::Cluster *const pCluster, const bool ignoreUnavailable, pandora::ClusterList &clusterListU, pandora::ClusterList &clusterListV, pandora::ClusterList &clusterListW) const
Get elements connected to a specified cluster.
ClusterNavigationMap m_clusterNavigationMapWU
The cluster navigation map W->U.
std::unordered_map< const pandora::Cluster *, pandora::ClusterList > ClusterNavigationMap
void AddAssociation(const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2)
Add an association between two clusters to the simple overlap tensor.
ClusterNavigationMap m_clusterNavigationMapVW
The cluster navigation map V->W.
const pandora::ClusterList & GetKeyClusters() const
Get the list of key clusters.
pandora::ClusterList m_keyClusters
The list of key clusters.
float m_minClusterXSpan
The min x span required in order to consider a cluster.
std::string m_outputPfoListName
The output pfo list name.
unsigned int m_minClusterCaloHits
The min number of hits in base cluster selection method.
bool m_vertexClusterMode
Whether to demand clusters are associated with vertices of existing particles.
void CreateTrackParticle(const pandora::ClusterList &clusterList) const
Create and save a track particle containing the provided clusters.
pandora::StringVector m_inputClusterListNames
The list of cluster list names.
float m_minClusterLengthSquared
The min length (squared) in base cluster selection method.
void ExamineTensor(const SimpleOverlapTensor &overlapTensor) const
Identify unambiguous cluster overlaps and resolve ambiguous overlaps, creating new track particles.
unsigned int m_slidingFitHalfWindow
The half window for the fit sliding result constructor.
void GetInputClusters(pandora::ClusterList &inputClusterListU, pandora::ClusterList &inputClusterListV, pandora::ClusterList &inputClusterListW) const
Get the input cluster lists for processing in this algorithm.
void VertexClusterSelection(const pandora::ClusterList &inputClusterList, pandora::ClusterList &selectedClusterList) const
Select a subset of input clusters nodally associated with the vertices of existing particles.
float m_minXOverlapFractionGaps
The min x overlap fraction when there are gaps involved.
void StandardClusterSelection(const pandora::ClusterList &inputClusterList, pandora::ClusterList &selectedClusterList) const
Select a subset of input clusters for processing in this algorithm.
bool m_checkGaps
Whether to check for gaps in the calculation of the overlap.
float m_maxVertexTransverseDistance
Vertex association check: max transverse distance cut.
float m_sampleStepSize
The sampling step size used in association checks, units cm.
float m_pseudoChi2Cut
The selection cut on the matched chi2.
void CalculateEffectiveOverlapFractions(const pandora::Cluster *const pCluster1, const float xMin1, const float xMax1, const pandora::Cluster *const pCluster2, const float xMin2, const float xMax2, float &xOverlapFraction1, float &xOverlapFraction2) const
Calculate effective overlap fractions taking into account gaps.
ParticleRecoveryAlgorithm()
Default constructor.
pandora::StatusCode Run()
Run the algorithm.
bool CheckConsistency(const pandora::Cluster *const pClusterU, const pandora::Cluster *const pClusterV, const pandora::Cluster *const pClusterW) const
Whether a trio of clusters are consistent with representing projections of the same 3d trajectory.
void SelectInputClusters(const pandora::ClusterList &inputClusterList, pandora::ClusterList &selectedClusterList) const
Select a subset of input clusters for processing in this algorithm.
void FindOverlaps(const pandora::ClusterList &clusterList1, const pandora::ClusterList &clusterList2, SimpleOverlapTensor &overlapTensor) const
Find cluster overlaps and record these in the overlap tensor.
bool IsOverlap(const pandora::Cluster *const pCluster1, const pandora::Cluster *const pCluster2) const
Whether two clusters overlap convincingly in x.
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
Read the algorithm settings.
float m_minXOverlapFraction
The min x overlap fraction required in order to id overlapping clusters.
float m_minVertexLongitudinalDistance
Vertex association check: min longitudinal distance cut.
void CalculateEffectiveSpan(const pandora::Cluster *const pCluster, const float xMin, const float xMax, float &xMinEff, float &xMaxEff) const
Calculate effective span for a given clsuter taking gaps into account.
TwoDSlidingFitResult class.
bool IsAvailable() const
Whether the cluster is available to be added to a particle flow object.
unsigned int GetNCaloHits() const
Get the number of calo hits in the cluster.
void GetClusterSpanX(float &xmin, float &xmax) const
Get minimum and maximum X positions of the calo hits in this cluster.
const DetectorGapList & GetDetectorGapList() const
Get the list of gaps in the active detector volume.
bool ShouldDisplayAlgorithmInfo() const
Whether to display algorithm information during processing.
ParticleFlowObject class.
static float GetParticleMass(const int pdgCode)
Get the mass of a particle type.
static int GetParticleCharge(const int pdgCode)
Get the charge of a particle type.
const Pandora & GetPandora() const
Get the associated pandora instance.
StatusCodeException class.
static StatusCode ReadVectorOfValues(const TiXmlHandle &xmlHandle, const std::string &xmlElementName, std::vector< T > &vector)
Read a vector of values from a (space separated) list in an xml element.
static StatusCode ReadValue(const TiXmlHandle &xmlHandle, const std::string &xmlElementName, T &t)
Read a value from an xml element.
HitType
Calorimeter hit type enum.
std::vector< const Cluster * > ClusterVector
MANAGED_CONTAINER< const Cluster * > ClusterList
std::vector< CartesianVector > CartesianPointVector
StatusCode
The StatusCode enum.
MANAGED_CONTAINER< const ParticleFlowObject * > PfoList