LArSoft  v06_85_00
Liquid Argon Software toolkit - http://larsoft.org/
BoxBoundedGeo.cxx
Go to the documentation of this file.
1 
9 
10 // LArSoft libraries
11 #include "larcorealg/Geometry/geo_vectors_utils.h" // geo::vect namespace
12 
13 
14 namespace geo
15 {
16  //----------------------------------------------------------------------------
18  (TVector3 const& point, double wiggle /* = 1.0 */) const
19  { return ContainsPosition(geo::vect::toPoint(point), wiggle); }
20 
21  //----------------------------------------------------------------------------
23  (double const* point, double wiggle /* = 1.0 */) const
24  { return ContainsPosition(geo::vect::makePointFromCoords(point), wiggle); }
25 
26  //----------------------------------------------------------------------------
27  std::vector<geo::Point_t> BoxBoundedGeo::GetIntersections(
28  geo::Point_t const& TrajectoryStart,
29  geo::Vector_t const& TrajectoryDirect
30  ) const {
31 
32  std::vector<geo::Point_t> IntersectionPoints;
33  std::vector<double> LineParameters;
34 
35  // Generate normal vectors and offsets for every plane of the box
36  // All normal vectors are headed outwards
37  // BUG the double brace syntax is required to work around clang bug 21629
38  // (https://bugs.llvm.org/show_bug.cgi?id=21629)
39  static std::array<geo::Vector_t, 6U> const NormalVectors = {{
40  -geo::Xaxis<geo::Vector_t>(), geo::Xaxis<geo::Vector_t>(), // anode, cathode,
41  -geo::Yaxis<geo::Vector_t>(), geo::Yaxis<geo::Vector_t>(), // bottom, top,
42  -geo::Zaxis<geo::Vector_t>(), geo::Zaxis<geo::Vector_t>() // upstream, downstream
43  }};
44  // BUG the double brace syntax is required to work around clang bug 21629
45  // (https://bugs.llvm.org/show_bug.cgi?id=21629)
46  std::array<geo::Point_t, 6U> const NormalVectorOffsets = {{
47  geo::Point_t{ Min().X(), Min().Y(), Min().Z() }, // Anode offset
48  geo::Point_t{ Max().X(), Min().Y(), Min().Z() }, // Cathode
49  geo::Point_t{ Min().X(), Min().Y(), Min().Z() }, // Bottom
50  geo::Point_t{ Min().X(), Max().Y(), Min().Z() }, // Top
51  geo::Point_t{ Min().X(), Min().Y(), Min().Z() }, // upstream
52  geo::Point_t{ Min().X(), Min().Y(), Max().Z() } // downstream
53  }};
54 
55  // Loop over all surfaces of the box
56  for(unsigned int face_no = 0; face_no < NormalVectors.size(); face_no++)
57  {
58  // Check if trajectory and surface are not parallel
59  if(NormalVectors[face_no].Dot(TrajectoryDirect))
60  {
61  // Calculate the line parameter for the intersection points
62  LineParameters.push_back( NormalVectors[face_no].Dot(NormalVectorOffsets.at(face_no) - TrajectoryStart)
63  / NormalVectors[face_no].Dot(TrajectoryDirect) );
64  }
65  else continue;
66 
67  // Calculate intersection point using the line parameter
68  IntersectionPoints.push_back( TrajectoryStart + LineParameters.back()*TrajectoryDirect );
69 
70  // Coordinate which should be ignored when checking for limits added by Christoph Rudolf von Rohr 05/21/2016
71  unsigned int NoCheckCoord;
72 
73  // Calculate NoCheckCoord out of the face_no
74  if(face_no % 2)
75  {
76  // Convert odd face number to coordinate
77  NoCheckCoord = (face_no - 1)/2;
78  }
79  else
80  {
81  // Convert even face number to coordinate
82  NoCheckCoord = face_no/2;
83  }
84 
85  // Loop over all three space coordinates
86  unsigned int coord = 0;
87  for(auto extractCoord: geo::vect::coordReaders<geo::Point_t>())
88  {
89  auto const lastPointCoord = geo::vect::bindCoord(IntersectionPoints.back(), extractCoord);
90  auto const minCoord = geo::vect::bindCoord(c_min, extractCoord);
91  auto const maxCoord = geo::vect::bindCoord(c_max, extractCoord);
92 
93  // Changed by Christoph Rudolf von Rohr 05/21/2016
94  // Then check if point is not within the surface limits at this coordinate, without looking
95  // at the plane normal vector coordinate. We can assume, that our algorithm already found this coordinate correctily.
96  // In rare cases, looking for boundaries in this coordinate causes unexpected behavior due to floating point inaccuracies.
97  if( coord++ != NoCheckCoord && ((lastPointCoord() > maxCoord()) || (lastPointCoord() < minCoord)) )
98  {
99  // if off limits, get rid of the useless data and break the coordinate loop
100  LineParameters.pop_back();
101  IntersectionPoints.pop_back();
102  break;
103  }
104  }// coordinate loop
105  }// Surcaces loop
106 
107  // sort points according to their parameter value (first is entry, second is exit)
108  if(LineParameters.size() == 2 && LineParameters.front() > LineParameters.back())
109  {
110  std::swap(IntersectionPoints.front(),IntersectionPoints.back());
111  }
112 
113  return IntersectionPoints;
114  } // GetIntersections()
115 
116  //----------------------------------------------------------------------------
117  std::vector<TVector3> BoxBoundedGeo::GetIntersections
118  (TVector3 const& TrajectoryStart, TVector3 const& TrajectoryDirect) const
119  {
120  std::vector<TVector3> intersections;
121  for (auto const& point: GetIntersections(geo::Point_t(TrajectoryStart), geo::Vector_t(TrajectoryDirect)))
122  intersections.emplace_back(point.X(), point.Y(), point.Z());
123  return intersections;
124  } // GetIntersections(TVector3)
125 
126  //----------------------------------------------------------------------------
128 
129  for (auto coordMan: geo::vect::coordManagers<geo::Point_t>()) {
130  auto min = geo::vect::bindCoord(c_min, coordMan);
131  auto max = geo::vect::bindCoord(c_max, coordMan);
132  if (min() > max()) {
133  auto temp = min();
134  min = max();
135  max = temp;
136  }
137  } // for
138 
139  } // BoxBoundedGeo::SortCoordinates()
140 
141  //----------------------------------------------------------------------------
142 
143 
144 } // namespace geo
145 
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:167
auto coord(Vector &v, unsigned int n) noexcept
Returns an object to manage the coordinate n of a vector.
::geo::Point_t toPoint(Point const &p)
Convert the specified point into a geo::Point_t.
constexpr auto bindCoord(Vector const &v, CoordReader_t< Vector > helper)
Binds the specified constant vector to the coordinate reader.
Int_t max
Definition: plot.C:27
Utilities to extend the interface of geometry vectors.
geo::Point_t Min() const
Returns the corner point with the smallest coordinates.
Provides a base class aware of world box coordinates.
Coords_t c_max
maximum coordinates (x, y, z)
Coords_t c_min
minimum coordinates (x, y, z)
void SortCoordinates()
Makes sure each coordinate of the minimum point is smaller than maximum.
Int_t min
Definition: plot.C:26
std::vector< TVector3 > GetIntersections(TVector3 const &TrajectoryStart, TVector3 const &TrajectoryDirect) const
Calculates the entry and exit points of a trajectory on the box surface.
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:187
GENVECTOR_CONSTEXPR::geo::Point_t makePointFromCoords(Coords &&coords)
Creates a geo::Point_t from its coordinates (see makeFromCoords()).
Namespace collecting geometry-related classes utilities.
geo::Point_t Max() const
Returns the corner point with the largest coordinates.
bool ContainsPosition(geo::Point_t const &point, double wiggle=1.0) const
Returns whether this volume contains the specified point.