LArSoft  v07_13_02
Liquid Argon Software toolkit - http://larsoft.org/
HitCreationBaseTool.cc
Go to the documentation of this file.
1 
9 #include "Pandora/AlgorithmHeaders.h"
10 
12 
14 
15 using namespace pandora;
16 
17 namespace lar_content
18 {
19 
20 HitCreationBaseTool::HitCreationBaseTool() :
21  m_sigmaX2(1.),
22  m_chiSquaredCut(1.)
23 {
24 }
25 
26 //------------------------------------------------------------------------------------------------------------------------------------------
27 
29 {
30 }
31 
32 //------------------------------------------------------------------------------------------------------------------------------------------
33 
34 void HitCreationBaseTool::GetBestPosition3D(const HitType hitType1, const HitType hitType2, const CartesianPointVector &fitPositionList1,
35  const CartesianPointVector &fitPositionList2, ProtoHit &protoHit) const
36 {
37  if (fitPositionList1.empty() && fitPositionList2.empty())
38  {
39  throw StatusCodeException(STATUS_CODE_NOT_FOUND);
40  }
41  else if (fitPositionList1.empty())
42  {
43  if (fitPositionList2.size() != 1)
44  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
45 
46  this->GetBestPosition3D(hitType2, fitPositionList2.front(), protoHit);
47  }
48  else if (fitPositionList2.empty())
49  {
50  if (fitPositionList1.size() != 1)
51  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
52 
53  this->GetBestPosition3D(hitType1, fitPositionList1.front(), protoHit);
54  }
55  else
56  {
57  for (const CartesianVector &fitPosition1 : fitPositionList1)
58  {
59  for (const CartesianVector &fitPosition2 : fitPositionList2)
60  {
61  ProtoHit thisProtoHit(protoHit.GetParentCaloHit2D());
62  this->GetBestPosition3D(hitType1, hitType2, fitPosition1, fitPosition2, thisProtoHit);
63 
64  if (!protoHit.IsPositionSet() || (thisProtoHit.GetChi2() < protoHit.GetChi2()))
65  protoHit = thisProtoHit;
66  }
67  }
68  }
69 
70  // TODO Add additional chi-squared at edge of detector
71 }
72 
73 //------------------------------------------------------------------------------------------------------------------------------------------
74 
75 void HitCreationBaseTool::GetBestPosition3D(const HitType hitType1, const HitType hitType2, const CartesianVector &fitPosition1,
76  const CartesianVector &fitPosition2, ProtoHit &protoHit) const
77 {
78  // TODO Input better uncertainties into this method (sigmaHit, sigmaFit, sigmaX)
79  const CaloHit *const pCaloHit2D(protoHit.GetParentCaloHit2D());
80  const HitType hitType(pCaloHit2D->GetHitType());
81 
82  const double sigmaFit(LArGeometryHelper::GetSigmaUVW(this->GetPandora()));
83  const double sigmaHit(sigmaFit);
84 
85  CartesianVector position3D(0.f, 0.f, 0.f);
86  double chi2(std::numeric_limits<double>::max());
87 
88  const double u((TPC_VIEW_U == hitType) ? pCaloHit2D->GetPositionVector().GetZ() : (TPC_VIEW_U == hitType1) ? fitPosition1.GetZ() : fitPosition2.GetZ());
89  const double v((TPC_VIEW_V == hitType) ? pCaloHit2D->GetPositionVector().GetZ() : (TPC_VIEW_V == hitType1) ? fitPosition1.GetZ() : fitPosition2.GetZ());
90  const double w((TPC_VIEW_W == hitType) ? pCaloHit2D->GetPositionVector().GetZ() : (TPC_VIEW_W == hitType1) ? fitPosition1.GetZ() : fitPosition2.GetZ());
91 
92  const double sigmaU((TPC_VIEW_U == hitType) ? sigmaHit : sigmaFit);
93  const double sigmaV((TPC_VIEW_V == hitType) ? sigmaHit : sigmaFit);
94  const double sigmaW((TPC_VIEW_W == hitType) ? sigmaHit : sigmaFit);
95 
97  this->GetPandora().GetPlugins()->GetLArTransformationPlugin()->GetMinChiSquaredYZ(u, v, w, sigmaU, sigmaV, sigmaW, bestY, bestZ, chi2);
98  position3D.SetValues(pCaloHit2D->GetPositionVector().GetX(), static_cast<float>(bestY), static_cast<float>(bestZ));
99 
100  const double deltaX1(pCaloHit2D->GetPositionVector().GetX() - fitPosition1.GetX());
101  const double deltaX2(pCaloHit2D->GetPositionVector().GetX() - fitPosition2.GetX());
102  const double chi2X(((deltaX1 * deltaX1) / m_sigmaX2) + ((deltaX2 * deltaX2) / m_sigmaX2));
103 
104  protoHit.SetPosition3D(position3D, chi2 + chi2X);
105  protoHit.AddTrajectorySample(TrajectorySample(fitPosition1, hitType1, sigmaFit));
106  protoHit.AddTrajectorySample(TrajectorySample(fitPosition2, hitType2, sigmaFit));
107 }
108 
109 //------------------------------------------------------------------------------------------------------------------------------------------
110 
111 void HitCreationBaseTool::GetBestPosition3D(const HitType hitType, const CartesianVector &fitPosition, ProtoHit &protoHit) const
112 {
113  // TODO Input better uncertainties into this method (sigmaHit, sigmaFit, sigmaX)
114  const CaloHit *const pCaloHit2D(protoHit.GetParentCaloHit2D());
115  const double sigmaFit(LArGeometryHelper::GetSigmaUVW(this->GetPandora()));
116 
117  if (pCaloHit2D->GetHitType() == hitType)
118  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
119 
120  CartesianVector position3D(0.f, 0.f, 0.f);
121  float chi2(std::numeric_limits<float>::max());
122  LArGeometryHelper::MergeTwoPositions3D(this->GetPandora(), pCaloHit2D->GetHitType(), hitType, pCaloHit2D->GetPositionVector(),
123  fitPosition, position3D, chi2);
124 
125  // ATTN Replace chi2 from LArGeometryHelper for consistency with three-view treatment (purely a measure of delta x)
126  const double deltaX(pCaloHit2D->GetPositionVector().GetX() - fitPosition.GetX());
127  const double chi2X((deltaX * deltaX) / m_sigmaX2);
128 
129  protoHit.SetPosition3D(position3D, chi2X);
130  protoHit.AddTrajectorySample(TrajectorySample(fitPosition, hitType, sigmaFit));
131 }
132 
133 //------------------------------------------------------------------------------------------------------------------------------------------
134 
135 StatusCode HitCreationBaseTool::ReadSettings(const pandora::TiXmlHandle xmlHandle)
136 {
137  double sigmaX(std::sqrt(m_sigmaX2));
138 
139  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle,
140  "SigmaX", sigmaX));
141 
142  m_sigmaX2 = sigmaX * sigmaX;
143 
144  if (m_sigmaX2 < std::numeric_limits<double>::epsilon())
145  {
146  std::cout << "HitCreationBaseTool - Invalid parameter, SigmaX: " << sigmaX << std::endl;
147  return STATUS_CODE_INVALID_PARAMETER;
148  }
149 
150  PANDORA_RETURN_RESULT_IF_AND_IF(STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle,
151  "ChiSquaredCut", m_chiSquaredCut));
152 
153  return STATUS_CODE_SUCCESS;
154 }
155 
156 } // namespace lar_content
Proto hits are temporary constructs to be used during iterative 3D hit procedure. ...
static float GetSigmaUVW(const pandora::Pandora &pandora, const float maxSigmaDiscrepancy=0.01)
Find the sigmaUVW value for the detector geometry.
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.
double m_sigmaX2
The sigmaX squared value, for calculation of chi2 deltaX term.
const pandora::CaloHit * GetParentCaloHit2D() const
Get the address of the parent 2D calo hit.
ThreeDHitCreationAlgorithm::TrajectorySample TrajectorySample
bool IsPositionSet() const
Whether the proto hit position is set.
double m_chiSquaredCut
The chi squared cut (accept only values below the cut value)
void AddTrajectorySample(const TrajectorySample &trajectorySample)
Add a trajectory sample.
TFile f
Definition: plotHisto.C:6
Header file for the geometry helper class.
Int_t max
Definition: plot.C:27
double GetChi2() const
Get the chi squared value.
Header file for the hit creation base tool.
virtual void GetBestPosition3D(const pandora::HitType hitType1, const pandora::HitType hitType2, const pandora::CartesianPointVector &fitPositionList1, const pandora::CartesianPointVector &fitPositionList2, ProtoHit &protoHit) const
Get the three dimensional position using a provided two dimensional calo hit and candidate fit positi...
virtual ~HitCreationBaseTool()
Destructor.
void SetPosition3D(const pandora::CartesianVector &position3D, const double chi2)
Set position 3D.
virtual pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
Float_t w
Definition: plot.C:23