Trajectory.cxx

Go to the documentation of this file.

#include "lardataobj/RecoBase/Trajectory.h"
#include "lardataobj/RecoBase/TrackingPlane.h"

// LArSoft libraries
#include "larcoreobj/SimpleTypesAndConstants/PhysicalConstants.h" // util::pi()

// C/C++ standard libraries
#include <ostream>
#include <stdexcept> // std::runtime_error
#include <utility>   // std::move()

//------------------------------------------------------------------------------
recob::Trajectory::Trajectory(Positions_t&& positions, Momenta_t&& momenta, bool hasMomenta)
  : fPositions(std::move(positions)), fMomenta(std::move(momenta)), fHasMomentum(hasMomenta)
{
  // invariant check
  if (fPositions.size() != fMomenta.size()) {
    throw std::runtime_error("recob::Trajectory constructed with " +
                             std::to_string(fPositions.size()) + " positions and " +
                             std::to_string(fMomenta.size()) +
                             " momenta! it requires the same number for both.");
  }
  if (fPositions.size() < 2) {
    throw std::runtime_error("recob::Trajectory constructed with " +
                             std::to_string(fPositions.size()) +
                             " trajectory points! it requires at least 2.");
  }
} // recob::Trajectory::Trajectory()

//------------------------------------------------------------------------------
double recob::Trajectory::Length(size_t startAt /* = 0 */) const
{

  // sanity check
  if (startAt >= LastPoint()) return 0.;

  // just sum the distance between all locations in the trajectory
  Point_t const* curr = &(LocationAtPoint(startAt));
  Point_t const* next = curr;
  Point_t const* last = &(End());
  Coord_t length = 0.0;
  while (next++ != last) {
    length += (*next - *curr).R();
    curr = next;
  } // while
  return length;
} // recob::Trajectory::Length()

//------------------------------------------------------------------------------
double recob::Trajectory::ZenithAngle(size_t p /* = 0 */) const
{

  // The zenith angle is defined by the angle between the track starting
  // direction and the y axis.
  // The y component of the starting direction is in fact the cosine of that
  // angle (and std::acos() conveniently returns a value in [0;pi]).
  // Our convention has the zenith angle the supplemental of the standard one.

  return util::pi<Coord_t>() - std::acos(DirectionAtPoint(p).Y());

} // recob::Trajectory::ZenithAngle()

//------------------------------------------------------------------------------
double recob::Trajectory::AzimuthAngle(size_t p /* = 0 */) const
{
  //
  // std::atan2(y, x) returns the angle of a point (x,y) respect to x axis.
  // In our case, the definition of the angle (0 for z axis, pi/2 for x axis)
  // translates atan2's y into our x, and x into z.
  //
  decltype(auto) startDir = DirectionAtPoint(p);
  return std::atan2(startDir.X(), startDir.Z());
} // recob::Trajectory::AzimuthAngle()

//------------------------------------------------------------------------------
recob::Trajectory::Vector_t recob::Trajectory::DirectionAtPoint(size_t i) const
{

  auto const& mom = MomentumVectorAtPoint(i);
  return HasMomentum() ? (mom / mom.R()) : mom;

} // recob::Trajectory::DirectionAtPoint()

//------------------------------------------------------------------------------
recob::Trajectory::Rotation_t recob::Trajectory::GlobalToLocalRotationAtPoint(size_t p) const
{
  return recob::tracking::Plane::Global3DToLocal3DRotation(DirectionAtPoint(p));
} // recob::Trajectory::GlobalToLocalRotationAtPoint()

//------------------------------------------------------------------------------
recob::Trajectory::Rotation_t recob::Trajectory::LocalToGlobalRotationAtPoint(size_t p) const
{
  return recob::tracking::Plane::Local3DToGlobal3DRotation(DirectionAtPoint(p));
} // recob::Trajectory::GlobalToLocalRotationAtPoint()

//------------------------------------------------------------------------------
std::ostream& recob::operator<<(std::ostream& out, recob::Trajectory const& traj)
{
  traj.Dump(out);
  return out;
}

//------------------------------------------------------------------------------