LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
TrackUtils.cxx
Go to the documentation of this file.
1 
10 // our header
12 
13 // LArSoft libraries
22 #include "larcoreobj/SimpleTypesAndConstants/geo_vectors.h" // geo::Vector_t
24 
25 // framework libraries
26 #include "cetlib_except/exception.h"
27 
28 // ROOT libraries
29 
30 // C/C++ standard libraries
31 #include <cmath>
32 
33 //------------------------------------------------------------------------------
35 {
36 
37  if (view == geo::kUnknown) {
38  throw cet::exception("TrackProjectedLength") << "cannot provide projected length for "
39  << "unknown view\n";
40  }
41  double length = 0.;
42 
43  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout>()->Get();
44  double angleToVert = 0.;
45  for (auto const& plane : wireReadoutGeom.Iterate<geo::PlaneGeo>(geo::CryostatID{0})) {
46  if (plane.View() == view) {
47  angleToVert = plane.Wire(0).ThetaZ(false) - 0.5 * ::util::pi<>();
48  break;
49  }
50  }
51 
52  // now loop over all points in the trajectory and add the contribution to the
53  // to the desired view
54 
55  for (size_t p = 1; p < track.NumberTrajectoryPoints(); ++p) {
56  const auto& pos_cur = track.LocationAtPoint(p);
57  const auto& pos_prev = track.LocationAtPoint(p - 1);
58  double dist =
59  std::sqrt(std::pow(pos_cur.x() - pos_prev.x(), 2) + std::pow(pos_cur.y() - pos_prev.y(), 2) +
60  std::pow(pos_cur.z() - pos_prev.z(), 2));
61 
62  // (sin(angleToVert),cos(angleToVert)) is the direction perpendicular to wire
63  // fDir[p-1] is the direction between the two relevant points
64  const auto& dir_prev = track.DirectionAtPoint(p - 1);
65  double cosgamma =
66  std::abs(std::sin(angleToVert) * dir_prev.Y() + std::cos(angleToVert) * dir_prev.Z());
67 
69  length += dist / cosgamma;
70  } // end loop over distances between trajectory points
71 
72  return length;
73 } // lar::util::TrackProjectedLength()
74 
75 //------------------------------------------------------------------------------
77  geo::View_t view,
78  size_t trajectory_point /* = 0 */)
79 {
80  /*
81  * The plan:
82  * 1. find the wire plane we are talking about
83  * (in the right TPC and with the right view)
84  * 2. ask the plane the answer
85  *
86  */
87 
88  if (trajectory_point >= track.NumberTrajectoryPoints()) {
89  cet::exception("TrackPitchInView")
90  << "ERROR: Asking for trajectory point #" << trajectory_point
91  << " when trajectory vector size is of size " << track.NumberTrajectoryPoints() << ".\n";
92  }
93  recob::Track::TrajectoryPoint_t const& point = track.TrajectoryPoint(trajectory_point);
94 
95  // this throws if the position is not in any TPC,
96  // or if there is no view with specified plane
97  auto const& geom = *(lar::providerFrom<geo::Geometry>());
98  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout>()->Get();
99  geo::PlaneGeo const& plane = wireReadoutGeom.Plane(geom.PositionToTPCID(point.position), view);
100 
101  //
102  // 2. project the direction of the track on that plane
103  //
104  // this is the projection of the direction of the track on the wire plane;
105  // it is 2D and its second component ("Y()") is on wire coordinate direction;
106  // remember that the projection modulus is smaller than 1 because it is
107  // the 3D direction versor, deprived of its drift direction component
108  auto const& proj = plane.Projection(point.direction());
109 
110  if (lar::util::RealComparisons(1e-4).zero(proj.Y())) {
111  throw cet::exception("Track") << "track at point #" << trajectory_point
112  << " is almost parallel to the wires in view "
113  << geo::PlaneGeo::ViewName(view) << " (wire direction is "
114  << plane.GetWireDirection() << "; track direction is "
115  << point.direction() << ", its projection on plane " << plane.ID()
116  << " is " << proj << ").\n";
117  }
118 
119  //
120  // 3. scale that projection so that it covers a wire pitch worth in the wire
121  // coordinate direction;
122  // WirePitch() is what gives this vector a physical size [cm]
123  //
124  return proj.R() / std::abs(proj.Y()) * plane.WirePitch();
125 
126 } // lar::util::TrackPitchInView()
127 
128 //------------------------------------------------------------------------------
A point in the trajectory, with position and momentum.
Definition: TrackingTypes.h:85
Utilities related to art service access.
static std::string ViewName(View_t view)
Returns the name of the specified view.
Definition: PlaneGeo.cxx:584
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Unknown view.
Definition: geo_types.h:138
Provides simple real number checks.
constexpr bool zero(Value_t value) const
Returns whether the value is no farther from 0 than the threshold.
Point_t const & LocationAtPoint(size_t i) const
Access to track position at different points.
Definition: Track.h:160
double TrackProjectedLength(recob::Track const &track, geo::View_t view)
Returns the length of the projection of a track on a view.
Definition: TrackUtils.cxx:34
size_t NumberTrajectoryPoints() const
Various functions related to the presence and the number of (valid) points.
Definition: Track.h:136
constexpr auto abs(T v)
Returns the absolute value of the argument.
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
Class for approximate comparisons.
Access the description of the physical detector geometry.
Definitions of geometry vector data types.
Vector_t const & GetWireDirection() const
Returns the direction of the wires.
Definition: PlaneGeo.h:390
Geometry information for a single wire plane.The plane is represented in the geometry by a solid whic...
Definition: PlaneGeo.h:67
TrajectoryPoint_t TrajectoryPoint(size_t i) const
Access to i-th TrajectoryPoint or its Flags.
Definition: Track.h:151
Vector_t direction() const
Returns the direction of the trajectory (unit vector of the momentum).
Definition: TrackingTypes.h:97
Encapsulate the construction of a single detector plane .
double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Provides recob::Track data product.
Float_t proj
Definition: plot.C:35
double TrackPitchInView(recob::Track const &track, geo::View_t view, size_t trajectory_point=0U)
Returns the projected length of track on a wire pitch step [cm].
Definition: TrackUtils.cxx:76
Vector_t DirectionAtPoint(size_t i) const
Access to track direction at different points.
Definition: Track.h:168
WireCoordProjection_t Projection(Point_t const &point) const
Returns the projection of the specified point on the plane.
Definition: PlaneGeo.h:755
Collection of Physical constants used in LArSoft.
Float_t e
Definition: plot.C:35
Float_t track
Definition: plot.C:35
PlaneID const & ID() const
Returns the identifier of this plane.
Definition: PlaneGeo.h:173
Utility functions to extract information from recob::Track
art framework interface to geometry description
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
Definition: Track.h:49
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
double WirePitch() const
Return the wire pitch (in centimeters). It is assumed constant.
Definition: PlaneGeo.h:312
Encapsulate the construction of a single detector plane .
The data type to uniquely identify a cryostat.
Definition: geo_types.h:187
Point_t position
position in the trajectory [cm].
Definition: TrackingTypes.h:87