LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
lar_content::ThreeDHitCreationAlgorithm Class Reference

ThreeDHitCreationAlgorithm::Algorithm class. More...

#include "ThreeDHitCreationAlgorithm.h"

Inheritance diagram for lar_content::ThreeDHitCreationAlgorithm:

Classes

class  ProtoHit
 Proto hits are temporary constructs to be used during iterative 3D hit procedure. More...
 
class  TrajectorySample
 Trajectory samples record the results of sampling a particles in a particular view. More...
 

Public Types

typedef std::vector< TrajectorySampleTrajectorySampleVector
 
typedef std::vector< ProtoHitProtoHitVector
 

Public Member Functions

 ThreeDHitCreationAlgorithm ()
 Default constructor. More...
 
void FilterCaloHitsByType (const pandora::CaloHitVector &inputCaloHitVector, const pandora::HitType hitType, pandora::CaloHitVector &outputCaloHitVector) const
 Get the subset of a provided calo hit vector corresponding to a specified hit type. More...
 

Private Types

typedef std::vector< HitCreationBaseTool * > HitCreationToolVector
 

Private Member Functions

pandora::StatusCode Run ()
 
void SeparateTwoDHits (const pandora::ParticleFlowObject *const pPfo, const ProtoHitVector &protoHitVector, pandora::CaloHitVector &remainingHitVector) const
 Get the list of 2D calo hits in a pfo for which 3D hits have and have not been created. More...
 
void IterativeTreatment (ProtoHitVector &protoHitVector) const
 Improve initial 3D hits by fitting proto hits and iteratively creating consisted 3D hit trajectory. More...
 
void ExtractResults (const ProtoHitVector &protoHitVector, double &chi2, pandora::CartesianPointVector &pointVector) const
 Extract key results from a provided proto hit vector. More...
 
double GetChi2WrtFit (const ThreeDSlidingFitResult &slidingFitResult, const ProtoHitVector &protoHitVector) const
 Receive a chi2 value indicating consistency of a list of proto hits with a provided 3D sliding fit trajectory. More...
 
double GetHitMovementChi2 (const ProtoHitVector &protoHitVector) const
 Receive a chi2 value indicating consistency of a list of proto hits with the original, input hit positions. More...
 
void RefineHitPositions (const ThreeDSlidingFitResult &slidingFitResult, ProtoHitVector &protoHitVector) const
 Refine the 3D hit positions (and chi2) for a list of proto hits, in accordance with a provided 3D sliding fit trajectory. More...
 
void CreateThreeDHits (const ProtoHitVector &protoHitVector, pandora::CaloHitList &newThreeDHits) const
 Create new three dimensional hits from two dimensional hits. More...
 
void CreateThreeDHit (const ProtoHit &protoHit, const pandora::CaloHit *&pCaloHit3D) const
 Create a new three dimensional hit from a two dimensional hit. More...
 
bool CheckThreeDHit (const ProtoHit &protoHit) const
 Check that a new three dimensional position is not unphysical. More...
 
void AddThreeDHitsToPfo (const pandora::ParticleFlowObject *const pPfo, const pandora::CaloHitList &caloHitList) const
 Add a specified list of three dimensional hits to a cluster in a pfo, creating the new cluster if required. More...
 
pandora::StatusCode ReadSettings (const pandora::TiXmlHandle xmlHandle)
 

Private Attributes

HitCreationToolVector m_algorithmToolVector
 The algorithm tool vector. More...
 
std::string m_inputPfoListName
 The name of the input pfo list. More...
 
std::string m_outputCaloHitListName
 The name of the output calo hit list. More...
 
std::string m_outputClusterListName
 The name of the output cluster list. More...
 
bool m_iterateTrackHits
 Whether to enable iterative improvement of 3D hits for track trajectories. More...
 
bool m_iterateShowerHits
 Whether to enable iterative improvement of 3D hits for showers. More...
 
unsigned int m_slidingFitHalfWindow
 The sliding linear fit half window. More...
 
unsigned int m_nHitRefinementIterations
 The maximum number of hit refinement iterations. More...
 
double m_sigma3DFitMultiplier
 Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit. More...
 
double m_iterationMaxChi2Ratio
 Max ratio between current and previous chi2 values to cease iterations. More...
 

Detailed Description

ThreeDHitCreationAlgorithm::Algorithm class.

Definition at line 27 of file ThreeDHitCreationAlgorithm.h.

Member Typedef Documentation

Constructor & Destructor Documentation

lar_content::ThreeDHitCreationAlgorithm::ThreeDHitCreationAlgorithm ( )

Default constructor.

Definition at line 27 of file ThreeDHitCreationAlgorithm.cc.

27  :
28  m_iterateTrackHits(true),
29  m_iterateShowerHits(false),
34 {
35 }
unsigned int m_slidingFitHalfWindow
The sliding linear fit half window.
unsigned int m_nHitRefinementIterations
The maximum number of hit refinement iterations.
bool m_iterateTrackHits
Whether to enable iterative improvement of 3D hits for track trajectories.
bool m_iterateShowerHits
Whether to enable iterative improvement of 3D hits for showers.
double m_sigma3DFitMultiplier
Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit.
double m_iterationMaxChi2Ratio
Max ratio between current and previous chi2 values to cease iterations.

Member Function Documentation

void lar_content::ThreeDHitCreationAlgorithm::AddThreeDHitsToPfo ( const pandora::ParticleFlowObject *const  pPfo,
const pandora::CaloHitList &  caloHitList 
) const
private

Add a specified list of three dimensional hits to a cluster in a pfo, creating the new cluster if required.

Parameters
pPfothe address of the pfo
caloHitListthe list of three dimensional hits

Definition at line 391 of file ThreeDHitCreationAlgorithm.cc.

References lar_content::LArPfoHelper::GetThreeDClusterList(), and m_outputClusterListName.

Referenced by Run().

392 {
393  if (caloHitList.empty())
394  throw StatusCodeException(STATUS_CODE_NOT_INITIALIZED);
395 
396  ClusterList threeDClusterList;
397  LArPfoHelper::GetThreeDClusterList(pPfo, threeDClusterList);
398 
399  if (!threeDClusterList.empty())
400  throw StatusCodeException(STATUS_CODE_FAILURE);
401 
402  const ClusterList *pClusterList(nullptr);
403  std::string clusterListName;
404  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::CreateTemporaryListAndSetCurrent(*this, pClusterList, clusterListName));
405 
406  PandoraContentApi::Cluster::Parameters parameters;
407  parameters.m_caloHitList.insert(parameters.m_caloHitList.end(), caloHitList.begin(), caloHitList.end());
408 
409  const Cluster *pCluster3D(nullptr);
410  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::Cluster::Create(*this, parameters, pCluster3D));
411 
412  if (!pCluster3D || !pClusterList || pClusterList->empty())
413  throw StatusCodeException(STATUS_CODE_FAILURE);
414 
415  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::SaveList<Cluster>(*this, m_outputClusterListName));
416  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::AddToPfo(*this, pPfo, pCluster3D));
417 }
std::string m_outputClusterListName
The name of the output cluster list.
static void GetThreeDClusterList(const pandora::ParticleFlowObject *const pPfo, pandora::ClusterList &clusterList)
Get the list of 3D clusters from an input pfo.
bool lar_content::ThreeDHitCreationAlgorithm::CheckThreeDHit ( const ProtoHit protoHit) const
private

Check that a new three dimensional position is not unphysical.

Parameters
protoHitthe proto hit
boolean

Definition at line 373 of file ThreeDHitCreationAlgorithm.cc.

References lar_content::ThreeDHitCreationAlgorithm::ProtoHit::GetPosition3D().

Referenced by CreateThreeDHit().

374 {
375  try
376  {
377  // Check that corresponding pseudo layer is within range - TODO use full LArTPC geometry here
378  (void)PandoraContentApi::GetPlugins(*this)->GetPseudoLayerPlugin()->GetPseudoLayer(protoHit.GetPosition3D());
379  }
380  catch (StatusCodeException &)
381  {
382  return false;
383  }
384 
385  // TODO Check against detector geometry
386  return true;
387 }
void lar_content::ThreeDHitCreationAlgorithm::CreateThreeDHit ( const ProtoHit protoHit,
const pandora::CaloHit *&  pCaloHit3D 
) const
private

Create a new three dimensional hit from a two dimensional hit.

Parameters
protoHitthe proto hit containing all required information
pCaloHit3Dto receive the address of the new three dimensional calo hit

Definition at line 338 of file ThreeDHitCreationAlgorithm.cc.

References CheckThreeDHit(), lar_content::ThreeDHitCreationAlgorithm::ProtoHit::GetParentCaloHit2D(), and lar_content::ThreeDHitCreationAlgorithm::ProtoHit::GetPosition3D().

Referenced by CreateThreeDHits().

339 {
340  if (!this->CheckThreeDHit(protoHit))
341  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
342 
343  PandoraContentApi::CaloHit::Parameters parameters;
344  parameters.m_positionVector = protoHit.GetPosition3D();
345  parameters.m_hitType = TPC_3D;
346 
347  const CaloHit *const pCaloHit2D(protoHit.GetParentCaloHit2D());
348  parameters.m_pParentAddress = static_cast<const void *>(pCaloHit2D);
349 
350  // TODO Check these parameters, especially new cell dimensions
351  parameters.m_cellThickness = pCaloHit2D->GetCellThickness();
352  parameters.m_cellGeometry = RECTANGULAR;
353  parameters.m_cellSize0 = pCaloHit2D->GetCellLengthScale();
354  parameters.m_cellSize1 = pCaloHit2D->GetCellLengthScale();
355  parameters.m_cellNormalVector = pCaloHit2D->GetCellNormalVector();
356  parameters.m_expectedDirection = pCaloHit2D->GetExpectedDirection();
357  parameters.m_nCellRadiationLengths = pCaloHit2D->GetNCellRadiationLengths();
358  parameters.m_nCellInteractionLengths = pCaloHit2D->GetNCellInteractionLengths();
359  parameters.m_time = pCaloHit2D->GetTime();
360  parameters.m_inputEnergy = pCaloHit2D->GetInputEnergy();
361  parameters.m_mipEquivalentEnergy = pCaloHit2D->GetMipEquivalentEnergy();
362  parameters.m_electromagneticEnergy = pCaloHit2D->GetElectromagneticEnergy();
363  parameters.m_hadronicEnergy = pCaloHit2D->GetHadronicEnergy();
364  parameters.m_isDigital = pCaloHit2D->IsDigital();
365  parameters.m_hitRegion = pCaloHit2D->GetHitRegion();
366  parameters.m_layer = pCaloHit2D->GetLayer();
367  parameters.m_isInOuterSamplingLayer = pCaloHit2D->IsInOuterSamplingLayer();
368  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::CaloHit::Create(*this, parameters, pCaloHit3D));
369 }
bool CheckThreeDHit(const ProtoHit &protoHit) const
Check that a new three dimensional position is not unphysical.
void lar_content::ThreeDHitCreationAlgorithm::CreateThreeDHits ( const ProtoHitVector protoHitVector,
pandora::CaloHitList &  newThreeDHits 
) const
private

Create new three dimensional hits from two dimensional hits.

Parameters
protoHitVectorthe input proto hit vector
newThreeDHitsto receive the addresses of the new three dimensional calo hits

Definition at line 322 of file ThreeDHitCreationAlgorithm.cc.

References CreateThreeDHit().

Referenced by Run().

323 {
324  for (const ProtoHit &protoHit : protoHitVector)
325  {
326  const CaloHit *pCaloHit3D(nullptr);
327  this->CreateThreeDHit(protoHit, pCaloHit3D);
328 
329  if (!pCaloHit3D)
330  throw StatusCodeException(STATUS_CODE_FAILURE);
331 
332  newThreeDHits.push_back(pCaloHit3D);
333  }
334 }
void CreateThreeDHit(const ProtoHit &protoHit, const pandora::CaloHit *&pCaloHit3D) const
Create a new three dimensional hit from a two dimensional hit.
void lar_content::ThreeDHitCreationAlgorithm::ExtractResults ( const ProtoHitVector protoHitVector,
double &  chi2,
pandora::CartesianPointVector &  pointVector 
) const
private

Extract key results from a provided proto hit vector.

Parameters
protoHitVectorthe proto hit vector
chi2to receive the sum of the proto hit chi2 values
pointVectorto receive a vector of proto hit 3D positions

Definition at line 183 of file ThreeDHitCreationAlgorithm.cc.

Referenced by IterativeTreatment().

184 {
185  chi2 = 0.;
186  pointVector.clear();
187 
188  for (const ProtoHit &protoHit : protoHitVector)
189  {
190  chi2 += protoHit.GetChi2();
191  pointVector.push_back(protoHit.GetPosition3D());
192  }
193 }
void lar_content::ThreeDHitCreationAlgorithm::FilterCaloHitsByType ( const pandora::CaloHitVector &  inputCaloHitVector,
const pandora::HitType  hitType,
pandora::CaloHitVector &  outputCaloHitVector 
) const

Get the subset of a provided calo hit vector corresponding to a specified hit type.

Parameters
inputCaloHitVectorthe input calo hit vector
hitTypethe hit type to filter upon
outputCaloHitVectorto receive the output calo hit vector

Definition at line 39 of file ThreeDHitCreationAlgorithm.cc.

Referenced by lar_content::ShowerHitsBaseTool::Run().

40 {
41  for (const CaloHit *const pCaloHit : inputCaloHitVector)
42  {
43  if (hitType == pCaloHit->GetHitType())
44  outputCaloHitVector.push_back(pCaloHit);
45  }
46 }
double lar_content::ThreeDHitCreationAlgorithm::GetChi2WrtFit ( const ThreeDSlidingFitResult slidingFitResult,
const ProtoHitVector protoHitVector 
) const
private

Receive a chi2 value indicating consistency of a list of proto hits with a provided 3D sliding fit trajectory.

Parameters
slidingFitResultthe 3D sliding fit result
protoHitVectorthe proto hit vector
Returns
the chi2 value

Definition at line 197 of file ThreeDHitCreationAlgorithm.cc.

References f, lar_content::ThreeDSlidingFitResult::GetGlobalFitPosition(), lar_content::ThreeDSlidingFitResult::GetLongitudinalDisplacement(), lar_content::LArGeometryHelper::GetSigmaUVW(), and m_sigma3DFitMultiplier.

Referenced by IterativeTreatment().

198 {
199  const double sigmaUVW(LArGeometryHelper::GetSigmaUVW(this->GetPandora()));
200  const double sigma3DFit(sigmaUVW * m_sigma3DFitMultiplier);
201 
202  double chi2WrtFit(0.);
203 
204  for (const ProtoHit &protoHit : protoHitVector)
205  {
206  CartesianVector pointOnFit(0.f, 0.f, 0.f);
207  const double rL(slidingFitResult.GetLongitudinalDisplacement(protoHit.GetPosition3D()));
208 
209  if (STATUS_CODE_SUCCESS != slidingFitResult.GetGlobalFitPosition(rL, pointOnFit))
210  continue;
211 
212  const double uFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoU(pointOnFit.GetY(), pointOnFit.GetZ()));
213  const double vFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoV(pointOnFit.GetY(), pointOnFit.GetZ()));
214  const double wFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoW(pointOnFit.GetY(), pointOnFit.GetZ()));
215 
216  const double outputU(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoU(
217  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ()));
218  const double outputV(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoV(
219  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ()));
220  const double outputW(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoW(
221  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ()));
222 
223  const double deltaUFit(uFit - outputU), deltaVFit(vFit - outputV), deltaWFit(wFit - outputW);
224  chi2WrtFit += ((deltaUFit * deltaUFit) / (sigma3DFit * sigma3DFit)) + ((deltaVFit * deltaVFit) / (sigma3DFit * sigma3DFit)) +
225  ((deltaWFit * deltaWFit) / (sigma3DFit * sigma3DFit));
226  }
227 
228  return chi2WrtFit;
229 }
static float GetSigmaUVW(const pandora::Pandora &pandora, const float maxSigmaDiscrepancy=0.01)
Find the sigmaUVW value for the detector geometry.
TFile f
Definition: plotHisto.C:6
double m_sigma3DFitMultiplier
Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit.
double lar_content::ThreeDHitCreationAlgorithm::GetHitMovementChi2 ( const ProtoHitVector protoHitVector) const
private

Receive a chi2 value indicating consistency of a list of proto hits with the original, input hit positions.

Parameters
protoHitVectorthe proto hit vector
Returns
the chi2 value

Definition at line 233 of file ThreeDHitCreationAlgorithm.cc.

References lar_content::LArGeometryHelper::GetSigmaUVW(), and lar_content::LArGeometryHelper::ProjectPosition().

234 {
235  const double sigmaUVW(LArGeometryHelper::GetSigmaUVW(this->GetPandora()));
236  double hitMovementChi2(0.);
237 
238  for (const ProtoHit &protoHit : protoHitVector)
239  {
240  const CaloHit *const pCaloHit2D(protoHit.GetParentCaloHit2D());
241  const HitType hitType(pCaloHit2D->GetHitType());
242 
243  const CartesianVector projectedPosition(LArGeometryHelper::ProjectPosition(this->GetPandora(), protoHit.GetPosition3D(), hitType));
244  const double delta(static_cast<double>(pCaloHit2D->GetPositionVector().GetZ() - projectedPosition.GetZ()));
245 
246  hitMovementChi2 += (delta * delta) / (sigmaUVW * sigmaUVW);
247  }
248 
249  return hitMovementChi2;
250 }
static float GetSigmaUVW(const pandora::Pandora &pandora, const float maxSigmaDiscrepancy=0.01)
Find the sigmaUVW value for the detector geometry.
static pandora::CartesianVector ProjectPosition(const pandora::Pandora &pandora, const pandora::CartesianVector &position3D, const pandora::HitType view)
Project 3D position into a given 2D view.
HitType
Definition: HitType.h:12
void lar_content::ThreeDHitCreationAlgorithm::IterativeTreatment ( ProtoHitVector protoHitVector) const
private

Improve initial 3D hits by fitting proto hits and iteratively creating consisted 3D hit trajectory.

Parameters
protoHitVectorthe vector of proto hits, describing current state of 3D hit construction

Definition at line 138 of file ThreeDHitCreationAlgorithm.cc.

References ExtractResults(), GetChi2WrtFit(), lar_content::LArGeometryHelper::GetWirePitch(), m_iterationMaxChi2Ratio, m_nHitRefinementIterations, m_slidingFitHalfWindow, and RefineHitPositions().

Referenced by Run().

139 {
140  const float pitchU{LArGeometryHelper::GetWirePitch(this->GetPandora(), TPC_VIEW_U)};
141  const float pitchV{LArGeometryHelper::GetWirePitch(this->GetPandora(), TPC_VIEW_V)};
142  const float pitchW{LArGeometryHelper::GetWirePitch(this->GetPandora(), TPC_VIEW_W)};
143  const float layerPitch(std::max({pitchU, pitchV, pitchW}));
144  const unsigned int layerWindow(m_slidingFitHalfWindow);
145 
146  double originalChi2(0.);
147  CartesianPointVector currentPoints3D;
148  this->ExtractResults(protoHitVector, originalChi2, currentPoints3D);
149 
150  try
151  {
152  const ThreeDSlidingFitResult originalSlidingFitResult(&currentPoints3D, layerWindow, layerPitch);
153  const double originalChi2WrtFit(this->GetChi2WrtFit(originalSlidingFitResult, protoHitVector));
154  double currentChi2(originalChi2 + originalChi2WrtFit);
155 
156  unsigned int nIterations(0);
157 
158  while (nIterations++ < m_nHitRefinementIterations)
159  {
160  ProtoHitVector newProtoHitVector(protoHitVector);
161  const ThreeDSlidingFitResult newSlidingFitResult(&currentPoints3D, layerWindow, layerPitch);
162  this->RefineHitPositions(newSlidingFitResult, newProtoHitVector);
163 
164  double newChi2(0.);
165  CartesianPointVector newPoints3D;
166  this->ExtractResults(newProtoHitVector, newChi2, newPoints3D);
167 
168  if (newChi2 > m_iterationMaxChi2Ratio * currentChi2)
169  break;
170 
171  currentChi2 = newChi2;
172  currentPoints3D = newPoints3D;
173  protoHitVector = newProtoHitVector;
174  }
175  }
176  catch (const StatusCodeException &)
177  {
178  }
179 }
unsigned int m_slidingFitHalfWindow
The sliding linear fit half window.
void RefineHitPositions(const ThreeDSlidingFitResult &slidingFitResult, ProtoHitVector &protoHitVector) const
Refine the 3D hit positions (and chi2) for a list of proto hits, in accordance with a provided 3D sli...
unsigned int m_nHitRefinementIterations
The maximum number of hit refinement iterations.
double GetChi2WrtFit(const ThreeDSlidingFitResult &slidingFitResult, const ProtoHitVector &protoHitVector) const
Receive a chi2 value indicating consistency of a list of proto hits with a provided 3D sliding fit tr...
static float GetWirePitch(const pandora::Pandora &pandora, const pandora::HitType view, const float maxWirePitchDiscrepancy=0.01)
Return the wire pitch.
void ExtractResults(const ProtoHitVector &protoHitVector, double &chi2, pandora::CartesianPointVector &pointVector) const
Extract key results from a provided proto hit vector.
double m_iterationMaxChi2Ratio
Max ratio between current and previous chi2 values to cease iterations.
StatusCode lar_content::ThreeDHitCreationAlgorithm::ReadSettings ( const pandora::TiXmlHandle  xmlHandle)
private

Definition at line 463 of file ThreeDHitCreationAlgorithm.cc.

References m_algorithmToolVector, m_inputPfoListName, m_iterateShowerHits, m_iterateTrackHits, m_iterationMaxChi2Ratio, m_nHitRefinementIterations, m_outputCaloHitListName, m_outputClusterListName, m_sigma3DFitMultiplier, and m_slidingFitHalfWindow.

464 {
465  AlgorithmToolVector algorithmToolVector;
466  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ProcessAlgorithmToolList(*this, xmlHandle, "HitCreationTools", algorithmToolVector));
467 
468  for (AlgorithmToolVector::const_iterator iter = algorithmToolVector.begin(), iterEnd = algorithmToolVector.end(); iter != iterEnd; ++iter)
469  {
470  HitCreationBaseTool *const pHitCreationTool(dynamic_cast<HitCreationBaseTool *>(*iter));
471 
472  if (!pHitCreationTool)
473  return STATUS_CODE_INVALID_PARAMETER;
474 
475  m_algorithmToolVector.push_back(pHitCreationTool);
476  }
477 
478  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "InputPfoListName", m_inputPfoListName));
479  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "OutputCaloHitListName", m_outputCaloHitListName));
480  PANDORA_RETURN_RESULT_IF(STATUS_CODE_SUCCESS, !=, XmlHelper::ReadValue(xmlHandle, "OutputClusterListName", m_outputClusterListName));
481 
482  PANDORA_RETURN_RESULT_IF_AND_IF(
483  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "IterateTrackHits", m_iterateTrackHits));
484 
485  PANDORA_RETURN_RESULT_IF_AND_IF(
486  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "IterateShowerHits", m_iterateShowerHits));
487 
488  PANDORA_RETURN_RESULT_IF_AND_IF(
489  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "SlidingFitHalfWindow", m_slidingFitHalfWindow));
490 
491  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=,
492  XmlHelper::ReadValue(xmlHandle, "NHitRefinementIterations", m_nHitRefinementIterations));
493 
494  PANDORA_RETURN_RESULT_IF_AND_IF(
495  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "Sigma3DFitMultiplier", m_sigma3DFitMultiplier));
496 
497  PANDORA_RETURN_RESULT_IF_AND_IF(
498  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "IterationMaxChi2Ratio", m_iterationMaxChi2Ratio));
499 
500  return STATUS_CODE_SUCCESS;
501 }
std::string m_inputPfoListName
The name of the input pfo list.
unsigned int m_slidingFitHalfWindow
The sliding linear fit half window.
std::string m_outputCaloHitListName
The name of the output calo hit list.
intermediate_table::const_iterator const_iterator
unsigned int m_nHitRefinementIterations
The maximum number of hit refinement iterations.
std::string m_outputClusterListName
The name of the output cluster list.
bool m_iterateTrackHits
Whether to enable iterative improvement of 3D hits for track trajectories.
bool m_iterateShowerHits
Whether to enable iterative improvement of 3D hits for showers.
HitCreationToolVector m_algorithmToolVector
The algorithm tool vector.
double m_sigma3DFitMultiplier
Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit.
double m_iterationMaxChi2Ratio
Max ratio between current and previous chi2 values to cease iterations.
void lar_content::ThreeDHitCreationAlgorithm::RefineHitPositions ( const ThreeDSlidingFitResult slidingFitResult,
ProtoHitVector protoHitVector 
) const
private

Refine the 3D hit positions (and chi2) for a list of proto hits, in accordance with a provided 3D sliding fit trajectory.

Parameters
slidingFitResultthe 3D sliding fit result
protoHitVectorthe proto hit vector, non const as proto hit properties will be updated

Definition at line 254 of file ThreeDHitCreationAlgorithm.cc.

References f, lar_content::ThreeDSlidingFitResult::GetGlobalFitPosition(), lar_content::ThreeDSlidingFitResult::GetLongitudinalDisplacement(), lar_content::LArGeometryHelper::GetSigmaUVW(), m_sigma3DFitMultiplier, and w.

Referenced by IterativeTreatment().

255 {
256  const double sigmaUVW(LArGeometryHelper::GetSigmaUVW(this->GetPandora()));
257  const double sigmaFit(sigmaUVW); // ATTN sigmaFit and sigmaHit here should agree with treatment in HitCreation tools
258  const double sigmaHit(sigmaUVW);
259  const double sigma3DFit(sigmaUVW * m_sigma3DFitMultiplier);
260 
261  for (ProtoHit &protoHit : protoHitVector)
262  {
263  CartesianVector pointOnFit(0.f, 0.f, 0.f);
264  const double rL(slidingFitResult.GetLongitudinalDisplacement(protoHit.GetPosition3D()));
265 
266  if (STATUS_CODE_SUCCESS != slidingFitResult.GetGlobalFitPosition(rL, pointOnFit))
267  continue;
268 
269  const CaloHit *const pCaloHit2D(protoHit.GetParentCaloHit2D());
270  const HitType hitType(pCaloHit2D->GetHitType());
271 
272  const double uFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoU(pointOnFit.GetY(), pointOnFit.GetZ()));
273  const double vFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoV(pointOnFit.GetY(), pointOnFit.GetZ()));
274  const double wFit(PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoW(pointOnFit.GetY(), pointOnFit.GetZ()));
275 
276  const double sigmaU((TPC_VIEW_U == hitType) ? sigmaHit : sigmaFit);
277  const double sigmaV((TPC_VIEW_V == hitType) ? sigmaHit : sigmaFit);
278  const double sigmaW((TPC_VIEW_W == hitType) ? sigmaHit : sigmaFit);
279 
280  CartesianVector position3D(0.f, 0.f, 0.f);
281  double chi2(std::numeric_limits<double>::max());
282  double u(std::numeric_limits<double>::max()), v(std::numeric_limits<double>::max()), w(std::numeric_limits<double>::max());
283 
284  if (protoHit.GetNTrajectorySamples() == 2)
285  {
286  u = (TPC_VIEW_U == hitType) ? pCaloHit2D->GetPositionVector().GetZ()
287  : (TPC_VIEW_U == protoHit.GetFirstTrajectorySample().GetHitType()) ? protoHit.GetFirstTrajectorySample().GetPosition().GetZ()
288  : protoHit.GetLastTrajectorySample().GetPosition().GetZ();
289  v = (TPC_VIEW_V == hitType) ? pCaloHit2D->GetPositionVector().GetZ()
290  : (TPC_VIEW_V == protoHit.GetFirstTrajectorySample().GetHitType()) ? protoHit.GetFirstTrajectorySample().GetPosition().GetZ()
291  : protoHit.GetLastTrajectorySample().GetPosition().GetZ();
292  w = (TPC_VIEW_W == hitType) ? pCaloHit2D->GetPositionVector().GetZ()
293  : (TPC_VIEW_W == protoHit.GetFirstTrajectorySample().GetHitType()) ? protoHit.GetFirstTrajectorySample().GetPosition().GetZ()
294  : protoHit.GetLastTrajectorySample().GetPosition().GetZ();
295  }
296  else if (protoHit.GetNTrajectorySamples() == 1)
297  {
298  u = PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoU(
299  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ());
300  v = PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoV(
301  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ());
302  w = PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->YZtoW(
303  protoHit.GetPosition3D().GetY(), protoHit.GetPosition3D().GetZ());
304  }
305  else
306  {
307  std::cout << "ThreeDHitCreationAlgorithm::IterativeTreatment - Unexpected number of trajectory samples" << std::endl;
308  throw StatusCodeException(STATUS_CODE_FAILURE);
309  }
310 
311  double bestY(std::numeric_limits<double>::max()), bestZ(std::numeric_limits<double>::max());
312  PandoraContentApi::GetPlugins(*this)->GetLArTransformationPlugin()->GetMinChiSquaredYZ(
313  u, v, w, sigmaU, sigmaV, sigmaW, uFit, vFit, wFit, sigma3DFit, bestY, bestZ, chi2);
314  position3D.SetValues(pCaloHit2D->GetPositionVector().GetX(), static_cast<float>(bestY), static_cast<float>(bestZ));
315 
316  protoHit.SetPosition3D(position3D, chi2);
317  }
318 }
static float GetSigmaUVW(const pandora::Pandora &pandora, const float maxSigmaDiscrepancy=0.01)
Find the sigmaUVW value for the detector geometry.
TFile f
Definition: plotHisto.C:6
HitType
Definition: HitType.h:12
double m_sigma3DFitMultiplier
Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit.
Float_t w
Definition: plot.C:20
StatusCode lar_content::ThreeDHitCreationAlgorithm::Run ( )
private

Definition at line 50 of file ThreeDHitCreationAlgorithm.cc.

References AddThreeDHitsToPfo(), CreateThreeDHits(), lar_content::LArPfoHelper::IsShower(), lar_content::LArPfoHelper::IsTrack(), IterativeTreatment(), m_algorithmToolVector, m_inputPfoListName, m_iterateShowerHits, m_iterateTrackHits, m_outputCaloHitListName, SeparateTwoDHits(), and lar_content::LArPfoHelper::SortByNHits().

51 {
52  const PfoList *pPfoList(nullptr);
53  PANDORA_RETURN_RESULT_IF_AND_IF(
54  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_INITIALIZED, !=, PandoraContentApi::GetList(*this, m_inputPfoListName, pPfoList));
55 
56  if (!pPfoList || pPfoList->empty())
57  {
58  if (PandoraContentApi::GetSettings(*this)->ShouldDisplayAlgorithmInfo())
59  std::cout << "ThreeDHitCreationAlgorithm: unable to find pfo list " << m_inputPfoListName << std::endl;
60 
61  return STATUS_CODE_SUCCESS;
62  }
63 
64  CaloHitList allNewThreeDHits;
65 
66  PfoVector pfoVector(pPfoList->begin(), pPfoList->end());
67  std::sort(pfoVector.begin(), pfoVector.end(), LArPfoHelper::SortByNHits);
68 
69  for (const ParticleFlowObject *const pPfo : pfoVector)
70  {
71  ProtoHitVector protoHitVector;
72 
73  for (HitCreationBaseTool *const pHitCreationTool : m_algorithmToolVector)
74  {
75  CaloHitVector remainingTwoDHits;
76  this->SeparateTwoDHits(pPfo, protoHitVector, remainingTwoDHits);
77 
78  if (remainingTwoDHits.empty())
79  break;
80 
81  pHitCreationTool->Run(this, pPfo, remainingTwoDHits, protoHitVector);
82  }
83 
85  this->IterativeTreatment(protoHitVector);
86 
87  if (protoHitVector.empty())
88  continue;
89 
90  CaloHitList newThreeDHits;
91  this->CreateThreeDHits(protoHitVector, newThreeDHits);
92  this->AddThreeDHitsToPfo(pPfo, newThreeDHits);
93 
94  allNewThreeDHits.insert(allNewThreeDHits.end(), newThreeDHits.begin(), newThreeDHits.end());
95  }
96 
97  if (!allNewThreeDHits.empty())
98  PANDORA_THROW_RESULT_IF(STATUS_CODE_SUCCESS, !=, PandoraContentApi::SaveList(*this, allNewThreeDHits, m_outputCaloHitListName));
99 
100  return STATUS_CODE_SUCCESS;
101 }
void IterativeTreatment(ProtoHitVector &protoHitVector) const
Improve initial 3D hits by fitting proto hits and iteratively creating consisted 3D hit trajectory...
static bool SortByNHits(const pandora::ParticleFlowObject *const pLhs, const pandora::ParticleFlowObject *const pRhs)
Sort pfos by number of constituent hits.
std::string m_inputPfoListName
The name of the input pfo list.
std::string m_outputCaloHitListName
The name of the output calo hit list.
void SeparateTwoDHits(const pandora::ParticleFlowObject *const pPfo, const ProtoHitVector &protoHitVector, pandora::CaloHitVector &remainingHitVector) const
Get the list of 2D calo hits in a pfo for which 3D hits have and have not been created.
static bool IsTrack(const pandora::ParticleFlowObject *const pPfo)
Return track flag based on Pfo Particle ID.
static bool IsShower(const pandora::ParticleFlowObject *const pPfo)
Return shower flag based on Pfo Particle ID.
bool m_iterateTrackHits
Whether to enable iterative improvement of 3D hits for track trajectories.
bool m_iterateShowerHits
Whether to enable iterative improvement of 3D hits for showers.
void CreateThreeDHits(const ProtoHitVector &protoHitVector, pandora::CaloHitList &newThreeDHits) const
Create new three dimensional hits from two dimensional hits.
HitCreationToolVector m_algorithmToolVector
The algorithm tool vector.
void AddThreeDHitsToPfo(const pandora::ParticleFlowObject *const pPfo, const pandora::CaloHitList &caloHitList) const
Add a specified list of three dimensional hits to a cluster in a pfo, creating the new cluster if req...
void lar_content::ThreeDHitCreationAlgorithm::SeparateTwoDHits ( const pandora::ParticleFlowObject *const  pPfo,
const ProtoHitVector protoHitVector,
pandora::CaloHitVector &  remainingHitVector 
) const
private

Get the list of 2D calo hits in a pfo for which 3D hits have and have not been created.

Parameters
pPfothe address of the pfo
protoHitVectorthe vector of proto hits, describing current state of 3D hit construction
remainingHitVectorto receive the vector of 2D calo hits for which 3D hits have not been created

Definition at line 105 of file ThreeDHitCreationAlgorithm.cc.

References lar_content::LArClusterHelper::GetClusterHitType(), lar_content::LArPfoHelper::GetTwoDClusterList(), and lar_content::LArClusterHelper::SortHitsByPosition().

Referenced by Run().

107 {
108  ClusterList twoDClusterList;
109  LArPfoHelper::GetTwoDClusterList(pPfo, twoDClusterList);
110  CaloHitList remainingHitList;
111 
112  for (const Cluster *const pCluster : twoDClusterList)
113  {
114  if (TPC_3D == LArClusterHelper::GetClusterHitType(pCluster))
115  throw StatusCodeException(STATUS_CODE_FAILURE);
116 
117  pCluster->GetOrderedCaloHitList().FillCaloHitList(remainingHitList);
118  }
119 
120  CaloHitSet remainingHitSet(remainingHitList.begin(), remainingHitList.end());
121 
122  for (const ProtoHit &protoHit : protoHitVector)
123  {
124  CaloHitSet::iterator eraseIter = remainingHitSet.find(protoHit.GetParentCaloHit2D());
125 
126  if (remainingHitSet.end() == eraseIter)
127  throw StatusCodeException(STATUS_CODE_FAILURE);
128 
129  remainingHitSet.erase(eraseIter);
130  }
131 
132  remainingHitVector.insert(remainingHitVector.end(), remainingHitSet.begin(), remainingHitSet.end());
133  std::sort(remainingHitVector.begin(), remainingHitVector.end(), LArClusterHelper::SortHitsByPosition);
134 }
intermediate_table::iterator iterator
static void GetTwoDClusterList(const pandora::ParticleFlowObject *const pPfo, pandora::ClusterList &clusterList)
Get the list of 2D clusters from an input pfo.
static pandora::HitType GetClusterHitType(const pandora::Cluster *const pCluster)
Get the hit type associated with a two dimensional cluster.
static bool SortHitsByPosition(const pandora::CaloHit *const pLhs, const pandora::CaloHit *const pRhs)
Sort calo hits by their position (use Z, followed by X, followed by Y)

Member Data Documentation

HitCreationToolVector lar_content::ThreeDHitCreationAlgorithm::m_algorithmToolVector
private

The algorithm tool vector.

Definition at line 272 of file ThreeDHitCreationAlgorithm.h.

Referenced by ReadSettings(), and Run().

std::string lar_content::ThreeDHitCreationAlgorithm::m_inputPfoListName
private

The name of the input pfo list.

Definition at line 274 of file ThreeDHitCreationAlgorithm.h.

Referenced by ReadSettings(), and Run().

bool lar_content::ThreeDHitCreationAlgorithm::m_iterateShowerHits
private

Whether to enable iterative improvement of 3D hits for showers.

Definition at line 279 of file ThreeDHitCreationAlgorithm.h.

Referenced by ReadSettings(), and Run().

bool lar_content::ThreeDHitCreationAlgorithm::m_iterateTrackHits
private

Whether to enable iterative improvement of 3D hits for track trajectories.

Definition at line 278 of file ThreeDHitCreationAlgorithm.h.

Referenced by ReadSettings(), and Run().

double lar_content::ThreeDHitCreationAlgorithm::m_iterationMaxChi2Ratio
private

Max ratio between current and previous chi2 values to cease iterations.

Definition at line 283 of file ThreeDHitCreationAlgorithm.h.

Referenced by IterativeTreatment(), and ReadSettings().

unsigned int lar_content::ThreeDHitCreationAlgorithm::m_nHitRefinementIterations
private

The maximum number of hit refinement iterations.

Definition at line 281 of file ThreeDHitCreationAlgorithm.h.

Referenced by IterativeTreatment(), and ReadSettings().

std::string lar_content::ThreeDHitCreationAlgorithm::m_outputCaloHitListName
private

The name of the output calo hit list.

Definition at line 275 of file ThreeDHitCreationAlgorithm.h.

Referenced by ReadSettings(), and Run().

std::string lar_content::ThreeDHitCreationAlgorithm::m_outputClusterListName
private

The name of the output cluster list.

Definition at line 276 of file ThreeDHitCreationAlgorithm.h.

Referenced by AddThreeDHitsToPfo(), and ReadSettings().

double lar_content::ThreeDHitCreationAlgorithm::m_sigma3DFitMultiplier
private

Multiplicative factor: sigmaUVW (same as sigmaHit and sigma2DFit) to sigma3DFit.

Definition at line 282 of file ThreeDHitCreationAlgorithm.h.

Referenced by GetChi2WrtFit(), ReadSettings(), and RefineHitPositions().

unsigned int lar_content::ThreeDHitCreationAlgorithm::m_slidingFitHalfWindow
private

The sliding linear fit half window.

Definition at line 280 of file ThreeDHitCreationAlgorithm.h.

Referenced by IterativeTreatment(), and ReadSettings().


The documentation for this class was generated from the following files: