All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
LArGeometryHelper.cc
Go to the documentation of this file.
1 /**
2  * @file larpandoracontent/LArHelpers/LArGeometryHelper.cc
3  *
4  * @brief Implementation of the geometry helper class.
5  *
6  * $Log: $
7  */
8 
9 #include "Managers/GeometryManager.h"
10 #include "Managers/PluginManager.h"
11 
12 #include "Geometry/DetectorGap.h"
13 #include "Geometry/LArTPC.h"
14 
15 #include "Objects/CartesianVector.h"
16 
17 #include "Pandora/Pandora.h"
18 
19 #include "larpandoracontent/LArHelpers/LArClusterHelper.h"
20 #include "larpandoracontent/LArHelpers/LArGeometryHelper.h"
21 
22 #include "larpandoracontent/LArObjects/LArTwoDSlidingFitResult.h"
23 
24 #include "Plugins/LArTransformationPlugin.h"
25 
26 using namespace pandora;
27 
28 namespace lar_content
29 {
30 
31 float LArGeometryHelper::MergeTwoPositions(const Pandora &pandora, const HitType view1, const HitType view2, const float position1, const float position2)
32 {
33  if (view1 == view2)
34  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
35 
36  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_V))
37  {
38  return pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoW(position1, position2);
39  }
40 
41  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_U))
42  {
43  return pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoW(position2, position1);
44  }
45 
46  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_U))
47  {
48  return pandora.GetPlugins()->GetLArTransformationPlugin()->WUtoV(position1, position2);
49  }
50 
51  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_W))
52  {
53  return pandora.GetPlugins()->GetLArTransformationPlugin()->WUtoV(position2, position1);
54  }
55 
56  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_W))
57  {
58  return pandora.GetPlugins()->GetLArTransformationPlugin()->VWtoU(position1, position2);
59  }
60 
61  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_V))
62  {
63  return pandora.GetPlugins()->GetLArTransformationPlugin()->VWtoU(position2, position1);
64  }
65 
66  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
67 }
68 
69 //------------------------------------------------------------------------------------------------------------------------------------------
70 
71 CartesianVector LArGeometryHelper::MergeTwoDirections(
72  const Pandora &pandora, const HitType view1, const HitType view2, const CartesianVector &direction1, const CartesianVector &direction2)
73 {
74  // x components must be consistent
75  if (direction1.GetX() * direction2.GetX() < 0.f)
76  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
77 
78  // scale x components to a common value
79  const float s1((std::fabs(direction1.GetX()) > std::numeric_limits<float>::epsilon()) ? 100.f * std::fabs(direction2.GetX()) : 1.f);
80  const float s2((std::fabs(direction2.GetX()) > std::numeric_limits<float>::epsilon()) ? 100.f * std::fabs(direction1.GetX()) : 1.f);
81 
82  float pX(s1 * direction1.GetX()), pU(0.f), pV(0.f), pW(0.f);
83  HitType newView(HIT_CUSTOM);
84 
85  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_V))
86  {
87  pU = s1 * direction1.GetZ();
88  pV = s2 * direction2.GetZ();
89  newView = TPC_VIEW_W;
90  }
91 
92  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_U))
93  {
94  pV = s1 * direction1.GetZ();
95  pU = s2 * direction2.GetZ();
96  newView = TPC_VIEW_W;
97  }
98 
99  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_U))
100  {
101  pW = s1 * direction1.GetZ();
102  pU = s2 * direction2.GetZ();
103  newView = TPC_VIEW_V;
104  }
105 
106  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_W))
107  {
108  pU = s1 * direction1.GetZ();
109  pW = s2 * direction2.GetZ();
110  newView = TPC_VIEW_V;
111  }
112 
113  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_W))
114  {
115  pV = s1 * direction1.GetZ();
116  pW = s2 * direction2.GetZ();
117  newView = TPC_VIEW_U;
118  }
119 
120  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_V))
121  {
122  pW = s1 * direction1.GetZ();
123  pV = s2 * direction2.GetZ();
124  newView = TPC_VIEW_U;
125  }
126 
127  if (newView == TPC_VIEW_W)
128  return CartesianVector(pX, 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoW(pU, pV)).GetUnitVector();
129 
130  if (newView == TPC_VIEW_U)
131  return CartesianVector(pX, 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->VWtoU(pV, pW)).GetUnitVector();
132 
133  if (newView == TPC_VIEW_V)
134  return CartesianVector(pX, 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->WUtoV(pW, pU)).GetUnitVector();
135 
136  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
137 }
138 
139 //------------------------------------------------------------------------------------------------------------------------------------------
140 
141 void LArGeometryHelper::MergeTwoPositions(const Pandora &pandora, const HitType view1, const HitType view2,
142  const CartesianVector &position1, const CartesianVector &position2, CartesianVector &position3, float &chiSquared)
143 {
144  if (view1 == view2)
145  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
146 
147  const float X3((position1.GetX() + position2.GetX()) / 2.f);
148  const float Z1(position1.GetZ());
149  const float Z2(position2.GetZ());
150  const float Z3(LArGeometryHelper::MergeTwoPositions(pandora, view1, view2, Z1, Z2));
151 
152  position3.SetValues(X3, 0.f, Z3);
153  const float sigmaUVW(LArGeometryHelper::GetSigmaUVW(pandora));
154  chiSquared = ((X3 - position1.GetX()) * (X3 - position1.GetX()) + (X3 - position2.GetX()) * (X3 - position2.GetX())) / (sigmaUVW * sigmaUVW);
155 }
156 
157 //------------------------------------------------------------------------------------------------------------------------------------------
158 
159 void LArGeometryHelper::MergeTwoPositions(const Pandora &pandora, const HitType view1, const HitType view2, const CartesianVector &position1,
160  const CartesianVector &position2, CartesianVector &outputU, CartesianVector &outputV, CartesianVector &outputW, float &chiSquared)
161 {
162  if (view1 == view2)
163  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
164 
165  CartesianVector position3(0.f, 0.f, 0.f);
166  LArGeometryHelper::MergeTwoPositions(pandora, view1, view2, position1, position2, position3, chiSquared);
167 
168  float aveU(0.f), aveV(0.f), aveW(0.f);
169  const float Z1(position1.GetZ()), Z2(position2.GetZ()), Z3(position3.GetZ()), aveX(position3.GetX());
170 
171  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_V))
172  {
173  aveU = Z1;
174  aveV = Z2;
175  aveW = Z3;
176  }
177 
178  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_W))
179  {
180  aveV = Z1;
181  aveW = Z2;
182  aveU = Z3;
183  }
184 
185  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_U))
186  {
187  aveW = Z1;
188  aveU = Z2;
189  aveV = Z3;
190  }
191 
192  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_U))
193  {
194  aveV = Z1;
195  aveU = Z2;
196  aveW = Z3;
197  }
198 
199  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_V))
200  {
201  aveW = Z1;
202  aveV = Z2;
203  aveU = Z3;
204  }
205 
206  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_W))
207  {
208  aveU = Z1;
209  aveW = Z2;
210  aveV = Z3;
211  }
212 
213  outputU.SetValues(aveX, 0.f, aveU);
214  outputV.SetValues(aveX, 0.f, aveV);
215  outputW.SetValues(aveX, 0.f, aveW);
216 }
217 
218 //------------------------------------------------------------------------------------------------------------------------------------------
219 
220 void LArGeometryHelper::MergeThreePositions(const Pandora &pandora, const HitType view1, const HitType view2, const HitType view3,
221  const CartesianVector &position1, const CartesianVector &position2, const CartesianVector &position3, CartesianVector &outputU,
222  CartesianVector &outputV, CartesianVector &outputW, float &chiSquared)
223 {
224  if ((view1 == view2) || (view2 == view3) || (view3 == view1))
225  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
226 
227  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_V))
228  {
229  return LArGeometryHelper::MergeThreePositions(pandora, position1, position2, position3, outputU, outputV, outputW, chiSquared);
230  }
231 
232  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_W))
233  {
234  return LArGeometryHelper::MergeThreePositions(pandora, position3, position1, position2, outputU, outputV, outputW, chiSquared);
235  }
236 
237  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_U))
238  {
239  return LArGeometryHelper::MergeThreePositions(pandora, position2, position3, position1, outputU, outputV, outputW, chiSquared);
240  }
241 
242  if ((view1 == TPC_VIEW_V) && (view2 == TPC_VIEW_U))
243  {
244  return LArGeometryHelper::MergeThreePositions(pandora, position2, position1, position3, outputU, outputV, outputW, chiSquared);
245  }
246 
247  if ((view1 == TPC_VIEW_W) && (view2 == TPC_VIEW_V))
248  {
249  return LArGeometryHelper::MergeThreePositions(pandora, position3, position2, position1, outputU, outputV, outputW, chiSquared);
250  }
251 
252  if ((view1 == TPC_VIEW_U) && (view2 == TPC_VIEW_W))
253  {
254  return LArGeometryHelper::MergeThreePositions(pandora, position1, position3, position2, outputU, outputV, outputW, chiSquared);
255  }
256 
257  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
258 }
259 
260 //------------------------------------------------------------------------------------------------------------------------------------------
261 
262 void LArGeometryHelper::MergeThreePositions(const Pandora &pandora, const CartesianVector &positionU, const CartesianVector &positionV,
263  const CartesianVector &positionW, CartesianVector &outputU, CartesianVector &outputV, CartesianVector &outputW, float &chiSquared)
264 {
265  const float YfromUV(pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoY(positionU.GetZ(), positionV.GetZ()));
266  const float YfromUW(pandora.GetPlugins()->GetLArTransformationPlugin()->UWtoY(positionU.GetZ(), positionW.GetZ()));
267  const float YfromVW(pandora.GetPlugins()->GetLArTransformationPlugin()->VWtoY(positionV.GetZ(), positionW.GetZ()));
268 
269  const float ZfromUV(pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoZ(positionU.GetZ(), positionV.GetZ()));
270  const float ZfromUW(pandora.GetPlugins()->GetLArTransformationPlugin()->UWtoZ(positionU.GetZ(), positionW.GetZ()));
271  const float ZfromVW(pandora.GetPlugins()->GetLArTransformationPlugin()->VWtoZ(positionV.GetZ(), positionW.GetZ()));
272 
273  // ATTN For detectors where w and z are equivalent, remain consistent with original treatment. TODO Use new treatment always.
274  const bool useOldWZEquivalentTreatment(std::fabs(ZfromUW - ZfromVW) < std::numeric_limits<float>::epsilon());
275  const float aveX((positionU.GetX() + positionV.GetX() + positionW.GetX()) / 3.f);
276  const float aveY(useOldWZEquivalentTreatment ? YfromUV : (YfromUV + YfromUW + YfromVW) / 3.f);
277  const float aveZ(useOldWZEquivalentTreatment ? (positionW.GetZ() + 2.f * ZfromUV) / 3.f : (ZfromUV + ZfromUW + ZfromVW) / 3.f);
278 
279  const float aveU(pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoU(aveY, aveZ));
280  const float aveV(pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoV(aveY, aveZ));
281  const float aveW(pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoW(aveY, aveZ));
282 
283  outputU.SetValues(aveX, 0.f, aveU);
284  outputV.SetValues(aveX, 0.f, aveV);
285  outputW.SetValues(aveX, 0.f, aveW);
286 
287  const float sigmaUVW(LArGeometryHelper::GetSigmaUVW(pandora));
288  chiSquared = ((outputU.GetX() - positionU.GetX()) * (outputU.GetX() - positionU.GetX()) +
289  (outputV.GetX() - positionV.GetX()) * (outputV.GetX() - positionV.GetX()) +
290  (outputW.GetX() - positionW.GetX()) * (outputW.GetX() - positionW.GetX()) +
291  (outputU.GetZ() - positionU.GetZ()) * (outputU.GetZ() - positionU.GetZ()) +
292  (outputV.GetZ() - positionV.GetZ()) * (outputV.GetZ() - positionV.GetZ()) +
293  (outputW.GetZ() - positionW.GetZ()) * (outputW.GetZ() - positionW.GetZ())) /
294  (sigmaUVW * sigmaUVW);
295 }
296 
297 //------------------------------------------------------------------------------------------------------------------------------------------
298 
299 void LArGeometryHelper::MergeTwoPositions3D(const Pandora &pandora, const HitType view1, const HitType view2,
300  const CartesianVector &position1, const CartesianVector &position2, CartesianVector &position3D, float &chiSquared)
301 {
302  CartesianVector positionU(0.f, 0.f, 0.f), positionV(0.f, 0.f, 0.f), positionW(0.f, 0.f, 0.f);
303  LArGeometryHelper::MergeTwoPositions(pandora, view1, view2, position1, position2, positionU, positionV, positionW, chiSquared);
304 
305  position3D.SetValues(positionW.GetX(), pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoY(positionU.GetZ(), positionV.GetZ()),
306  pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoZ(positionU.GetZ(), positionV.GetZ()));
307 }
308 
309 //------------------------------------------------------------------------------------------------------------------------------------------
310 
311 void LArGeometryHelper::MergeThreePositions3D(const Pandora &pandora, const HitType view1, const HitType view2, const HitType view3,
312  const CartesianVector &position1, const CartesianVector &position2, const CartesianVector &position3, CartesianVector &position3D, float &chiSquared)
313 {
314  CartesianVector positionU(0.f, 0.f, 0.f), positionV(0.f, 0.f, 0.f), positionW(0.f, 0.f, 0.f);
315  LArGeometryHelper::MergeThreePositions(pandora, view1, view2, view3, position1, position2, position3, positionU, positionV, positionW, chiSquared);
316 
317  position3D.SetValues(positionW.GetX(), pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoY(positionU.GetZ(), positionV.GetZ()),
318  pandora.GetPlugins()->GetLArTransformationPlugin()->UVtoZ(positionU.GetZ(), positionV.GetZ()));
319 }
320 
321 //------------------------------------------------------------------------------------------------------------------------------------------
322 
323 CartesianVector LArGeometryHelper::ProjectPosition(const Pandora &pandora, const CartesianVector &position3D, const HitType view)
324 {
325  if (view == TPC_VIEW_U)
326  {
327  return CartesianVector(
328  position3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoU(position3D.GetY(), position3D.GetZ()));
329  }
330 
331  else if (view == TPC_VIEW_V)
332  {
333  return CartesianVector(
334  position3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoV(position3D.GetY(), position3D.GetZ()));
335  }
336 
337  else if (view == TPC_VIEW_W)
338  {
339  return CartesianVector(
340  position3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoW(position3D.GetY(), position3D.GetZ()));
341  }
342 
343  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
344 }
345 
346 //------------------------------------------------------------------------------------------------------------------------------------------
347 
348 CartesianVector LArGeometryHelper::ProjectDirection(const Pandora &pandora, const CartesianVector &direction3D, const HitType view)
349 {
350  if (view == TPC_VIEW_U)
351  {
352  return CartesianVector(
353  direction3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoU(direction3D.GetY(), direction3D.GetZ()))
354  .GetUnitVector();
355  }
356 
357  else if (view == TPC_VIEW_V)
358  {
359  return CartesianVector(
360  direction3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoV(direction3D.GetY(), direction3D.GetZ()))
361  .GetUnitVector();
362  }
363 
364  else if (view == TPC_VIEW_W)
365  {
366  return CartesianVector(
367  direction3D.GetX(), 0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoW(direction3D.GetY(), direction3D.GetZ()))
368  .GetUnitVector();
369  }
370 
371  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
372 }
373 
374 //------------------------------------------------------------------------------------------------------------------------------------------
375 
376 float LArGeometryHelper::GetWirePitch(const Pandora &pandora, const HitType view, const float maxWirePitchDiscrepancy)
377 {
378  if (view != TPC_VIEW_U && view != TPC_VIEW_V && view != TPC_VIEW_W)
379  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
380 
381  const LArTPCMap &larTPCMap(pandora.GetGeometry()->GetLArTPCMap());
382 
383  if (larTPCMap.empty())
384  {
385  std::cout << "LArGeometryHelper::GetWirePitch - LArTPC description not registered with Pandora as required " << std::endl;
386  throw StatusCodeException(STATUS_CODE_NOT_INITIALIZED);
387  }
388 
389  const LArTPC *const pFirstLArTPC(larTPCMap.begin()->second);
390  const float wirePitch(view == TPC_VIEW_U ? pFirstLArTPC->GetWirePitchU()
391  : (view == TPC_VIEW_V ? pFirstLArTPC->GetWirePitchV() : pFirstLArTPC->GetWirePitchW()));
392 
393  for (const LArTPCMap::value_type &mapEntry : larTPCMap)
394  {
395  const LArTPC *const pLArTPC(mapEntry.second);
396  const float alternateWirePitch(
397  view == TPC_VIEW_U ? pLArTPC->GetWirePitchU() : (view == TPC_VIEW_V ? pLArTPC->GetWirePitchV() : pLArTPC->GetWirePitchW()));
398 
399  if (std::fabs(wirePitch - alternateWirePitch) > maxWirePitchDiscrepancy)
400  {
401  std::cout << "LArGeometryHelper::GetWirePitch - LArTPC configuration not supported" << std::endl;
402  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
403  }
404  }
405 
406  return wirePitch;
407 }
408 
409 //------------------------------------------------------------------------------------------------------------------------------------------
410 
411 CartesianVector LArGeometryHelper::GetWireAxis(const Pandora &pandora, const HitType view)
412 {
413  if (view == TPC_VIEW_U)
414  {
415  return CartesianVector(0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoU(1.f, 0.f),
416  pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoU(0.f, 1.f));
417  }
418 
419  else if (view == TPC_VIEW_V)
420  {
421  return CartesianVector(0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoV(1.f, 0.f),
422  pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoV(0.f, 1.f));
423  }
424 
425  else if (view == TPC_VIEW_W)
426  {
427  return CartesianVector(0.f, pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoW(1.f, 0.f),
428  pandora.GetPlugins()->GetLArTransformationPlugin()->YZtoW(0.f, 1.f));
429  }
430 
431  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
432 }
433 
434 //------------------------------------------------------------------------------------------------------------------------------------------
435 
436 void LArGeometryHelper::GetCommonDaughterVolumes(const Cluster *const pCluster1, const Cluster *const pCluster2, UIntSet &intersect)
437 {
438  UIntSet daughterVolumeIds1, daughterVolumeIds2;
439 
440  LArClusterHelper::GetDaughterVolumeIDs(pCluster1, daughterVolumeIds1);
441  LArClusterHelper::GetDaughterVolumeIDs(pCluster2, daughterVolumeIds2);
442 
443  std::set_intersection(daughterVolumeIds1.begin(), daughterVolumeIds1.end(), daughterVolumeIds2.begin(), daughterVolumeIds2.end(),
444  std::inserter(intersect, intersect.begin()));
445 }
446 
447 //------------------------------------------------------------------------------------------------------------------------------------------
448 
449 bool LArGeometryHelper::IsInGap(const Pandora &pandora, const CartesianVector &testPoint2D, const HitType hitType, const float gapTolerance)
450 {
451  // ATTN: input test point MUST be a 2D position vector
452  for (const DetectorGap *const pDetectorGap : pandora.GetGeometry()->GetDetectorGapList())
453  {
454  if (pDetectorGap->IsInGap(testPoint2D, hitType, gapTolerance))
455  return true;
456  }
457 
458  return false;
459 }
460 
461 //------------------------------------------------------------------------------------------------------------------------------------------
462 
463 bool LArGeometryHelper::IsInGap3D(const Pandora &pandora, const CartesianVector &testPoint3D, const HitType hitType, const float gapTolerance)
464 {
465  const CartesianVector testPoint2D(LArGeometryHelper::ProjectPosition(pandora, testPoint3D, hitType));
466  return LArGeometryHelper::IsInGap(pandora, testPoint2D, hitType, gapTolerance);
467 }
468 
469 //------------------------------------------------------------------------------------------------------------------------------------------
470 
471 bool LArGeometryHelper::IsXSamplingPointInGap(const Pandora &pandora, const float xSample, const TwoDSlidingFitResult &slidingFitResult, const float gapTolerance)
472 {
473  const HitType hitType(LArClusterHelper::GetClusterHitType(slidingFitResult.GetCluster()));
474  const CartesianVector minLayerPosition(slidingFitResult.GetGlobalMinLayerPosition());
475  const CartesianVector maxLayerPosition(slidingFitResult.GetGlobalMaxLayerPosition());
476 
477  const bool minLayerIsAtLowX(minLayerPosition.GetX() < maxLayerPosition.GetX());
478  const CartesianVector &lowXCoordinate(minLayerIsAtLowX ? minLayerPosition : maxLayerPosition);
479  const CartesianVector &highXCoordinate(minLayerIsAtLowX ? maxLayerPosition : minLayerPosition);
480 
481  if ((xSample > lowXCoordinate.GetX()) && (xSample < highXCoordinate.GetX()))
482  {
483  CartesianVector slidingFitPosition(0.f, 0.f, 0.f);
484 
485  if (STATUS_CODE_SUCCESS == slidingFitResult.GetGlobalFitPositionAtX(xSample, slidingFitPosition))
486  return (LArGeometryHelper::IsInGap(pandora, slidingFitPosition, hitType, gapTolerance));
487  }
488 
489  const CartesianVector lowXDirection(
490  minLayerIsAtLowX ? slidingFitResult.GetGlobalMinLayerDirection() : slidingFitResult.GetGlobalMaxLayerDirection());
491  const CartesianVector highXDirection(
492  minLayerIsAtLowX ? slidingFitResult.GetGlobalMaxLayerDirection() : slidingFitResult.GetGlobalMinLayerDirection());
493 
494  const bool sampleIsNearerToLowX(std::fabs(xSample - lowXCoordinate.GetX()) < std::fabs(xSample - highXCoordinate.GetX()));
495  const CartesianVector &startPosition(sampleIsNearerToLowX ? lowXCoordinate : highXCoordinate);
496  const CartesianVector &startDirection(sampleIsNearerToLowX ? lowXDirection : highXDirection);
497 
498  if (std::fabs(startDirection.GetX()) < std::numeric_limits<float>::epsilon())
499  throw StatusCodeException(STATUS_CODE_NOT_FOUND);
500 
501  const float pathLength((xSample - startPosition.GetX()) / startDirection.GetX());
502  const CartesianVector samplingPoint(startPosition + startDirection * pathLength);
503 
504  return (LArGeometryHelper::IsInGap(pandora, samplingPoint, hitType, gapTolerance));
505 }
506 
507 //------------------------------------------------------------------------------------------------------------------------------------------
508 
509 float LArGeometryHelper::CalculateGapDeltaZ(const Pandora &pandora, const float minZ, const float maxZ, const HitType hitType)
510 {
511  if (maxZ - minZ < std::numeric_limits<float>::epsilon())
512  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
513 
514  float gapDeltaZ(0.f);
515 
516  for (const DetectorGap *const pDetectorGap : pandora.GetGeometry()->GetDetectorGapList())
517  {
518  const LineGap *const pLineGap = dynamic_cast<const LineGap *>(pDetectorGap);
519 
520  if (!pLineGap)
521  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
522 
523  const LineGapType lineGapType(pLineGap->GetLineGapType());
524 
525  if (!(((TPC_VIEW_U == hitType) && (TPC_WIRE_GAP_VIEW_U == lineGapType)) || ((TPC_VIEW_V == hitType) && (TPC_WIRE_GAP_VIEW_V == lineGapType)) ||
526  ((TPC_VIEW_W == hitType) && (TPC_WIRE_GAP_VIEW_W == lineGapType))))
527  {
528  continue;
529  }
530 
531  if ((pLineGap->GetLineStartZ() > maxZ) || (pLineGap->GetLineEndZ() < minZ))
532  continue;
533 
534  const float gapMinZ(std::max(minZ, pLineGap->GetLineStartZ()));
535  const float gapMaxZ(std::min(maxZ, pLineGap->GetLineEndZ()));
536 
537  if ((gapMaxZ - gapMinZ) > std::numeric_limits<float>::epsilon())
538  gapDeltaZ += (gapMaxZ - gapMinZ);
539  }
540 
541  return gapDeltaZ;
542 }
543 
544 //------------------------------------------------------------------------------------------------------------------------------------------
545 
546 float LArGeometryHelper::GetSigmaUVW(const Pandora &pandora, const float maxSigmaDiscrepancy)
547 {
548  const LArTPCMap &larTPCMap(pandora.GetGeometry()->GetLArTPCMap());
549 
550  if (larTPCMap.empty())
551  {
552  std::cout << "LArGeometryHelper::GetSigmaUVW - LArTPC description not registered with Pandora as required " << std::endl;
553  throw StatusCodeException(STATUS_CODE_NOT_INITIALIZED);
554  }
555 
556  const LArTPC *const pFirstLArTPC(larTPCMap.begin()->second);
557  const float sigmaUVW(pFirstLArTPC->GetSigmaUVW());
558 
559  for (const LArTPCMap::value_type &mapEntry : larTPCMap)
560  {
561  const LArTPC *const pLArTPC(mapEntry.second);
562 
563  if (std::fabs(sigmaUVW - pLArTPC->GetSigmaUVW()) > maxSigmaDiscrepancy)
564  {
565  std::cout << "LArGeometryHelper::GetSigmaUVW - Plugin does not support provided LArTPC configurations " << std::endl;
566  throw StatusCodeException(STATUS_CODE_INVALID_PARAMETER);
567  }
568  }
569 
570  return sigmaUVW;
571 }
572 
573 } // namespace lar_content
lar_content::TwoDSlidingFitResult::GetGlobalMinLayerDirection
pandora::CartesianVector GetGlobalMinLayerDirection() const
Get global direction corresponding to the fit result in minimum fit layer.
Definition: LArTwoDSlidingFitResult.cc:265
lar_content::LArGeometryHelper::UIntSet
std::set< unsigned int > UIntSet
Definition: LArGeometryHelper.h:36
LArGeometryHelper.h
Header file for the geometry helper class.
LArClusterHelper.h
Header file for the cluster helper class.
LArTwoDSlidingFitResult.h
Header file for the lar two dimensional sliding fit result class.
lar_content::TwoDSlidingFitResult::GetGlobalMinLayerPosition
pandora::CartesianVector GetGlobalMinLayerPosition() const
Get global position corresponding to the fit result in minimum fit layer.
Definition: LArTwoDSlidingFitResult.cc:239
util::pathLength
double pathLength(BIter begin, EIter end)
Returns the total length of the specified path.
Definition: TrajectoryUtils.h:420
Z1
BEGIN_PROLOG supported so bottom corner of box bottom corner of box bottom corner of box top corner of box top corner of box Z1
Definition: radiological_gen.fcl:16
lar_content::TwoDSlidingFitResult::GetCluster
const pandora::Cluster * GetCluster() const
Get the address of the cluster, if originally provided.
Definition: LArTwoDSlidingFitResult.cc:143
pandora
required by fuzzyCluster table::sbnd_g4_services gaushitTruthMatch pandora
Definition: reco_sbnd.fcl:182
lar_content::TwoDSlidingFitResult::GetGlobalFitPositionAtX
pandora::StatusCode GetGlobalFitPositionAtX(const float x, pandora::CartesianVector &position) const
Get global fit position for a given input x coordinate.
Definition: LArTwoDSlidingFitResult.cc:375
lar_content::TwoDSlidingFitResult::GetGlobalMaxLayerDirection
pandora::CartesianVector GetGlobalMaxLayerDirection() const
Get global direction corresponding to the fit result in maximum fit layer.
Definition: LArTwoDSlidingFitResult.cc:278
lar_content::TwoDSlidingFitResult::GetGlobalMaxLayerPosition
pandora::CartesianVector GetGlobalMaxLayerPosition() const
Get global position corresponding to the fit result in maximum fit layer.
Definition: LArTwoDSlidingFitResult.cc:252
lar_content::TwoDSlidingFitResult
TwoDSlidingFitResult class.
Definition: LArTwoDSlidingFitResult.h:23
cout
BEGIN_PROLOG could also be cout
Definition: messageservice.fcl:10
Generated by   doxygen 1.8.5