WireLine.h

/**************************************************************************
 * BASF2 (Belle Analysis Framework 2)                                     *
 * Copyright(C) 2012 - Belle II Collaboration                             *
 *                                                                        *
 * Author: The Belle II Collaboration                                     *
 * Contributors: Oliver Frost                                             *
 *                                                                        *
 * This software is provided "as is" without any warranty.                *
 **************************************************************************/
#pragma once
 
#include <tracking/trackFindingCDC/geometry/Vector3D.h>
#include <tracking/trackFindingCDC/geometry/Vector2D.h>
 
#include <cmath>
 
namespace Belle2 {
  namespace TrackFindingCDC {
    class WireLine  {
 
    public:
      WireLine()
      {}
 
      WireLine(const Vector3D& forward, const Vector3D& backward, double sagCoeff);
 
      WireLine movedBy(const Vector3D& offset) const;
 
      WireLine movedBy(const Vector2D& offset) const;
 
      Vector3D nominalPos3DAtZ(const double z) const
      { return Vector3D(nominalPos2DAtZ(z), z); }
 
      Vector2D nominalPos2DAtZ(const double z) const
      { return refPos2D() + nominalMovePerZ() * z; }
 
      Vector3D sagPos3DAtZ(const double z) const
      { return Vector3D(sagPos2DAtZ(z), z); }
 
      Vector2D sagPos2DAtZ(const double z) const
      {
        Vector2D pos2D = nominalPos2DAtZ(z);
        double xzFactorSquared = 1 + nominalMovePerZ().x() * nominalMovePerZ().x();
        double sagY = - sagCoeff() * xzFactorSquared * (forwardZ() - z) * (z - backwardZ());
        pos2D.setY(pos2D.y() + sagY);
        return pos2D;
      }
 
      const Vector2D& nominalMovePerZ() const
      { return m_nominalMovePerZ; }
 
      Vector2D sagMovePerZ(const double z) const
      {
        Vector2D movePerZ = nominalMovePerZ();
        double xzFactorSquared = 1 + nominalMovePerZ().x() * nominalMovePerZ().x();
        double sagDYDZ = - sagCoeff() * xzFactorSquared * (forwardZ() + backwardZ() - 2 * z);
        movePerZ.setY(movePerZ.y() + sagDYDZ);
        return movePerZ;
      }
 
      double nominalDistance(const Vector3D& pos3D) const
      { return (pos3D - refPos3D()).orthogonalComp(Vector3D(nominalMovePerZ(), 1)); }
 
      double sagDistance(const Vector3D& pos3D) const
      {
        Vector3D wirePos3D = sagPos3DAtZ(pos3D.z());
        Vector3D movePerZ(sagMovePerZ(pos3D.z()), 1);
        return (pos3D - wirePos3D).orthogonalComp(movePerZ);
      }
 
      Vector3D nominalClosest3D(const Vector3D& point) const
      { return  refPos3D() - (point - refPos3D()).parallelVector(Vector3D(nominalMovePerZ(), 1)); }
 
      Vector3D sagClosest3D(const Vector3D& point) const
      {
        Vector3D wirePos3D = sagPos3DAtZ(point.z());
        Vector3D movePerZ(sagMovePerZ(point.z()), 1);
        return  wirePos3D - (point - wirePos3D).parallelVector(movePerZ);
      }
 
      Vector3D forward3D() const
      { return nominalPos3DAtZ(forwardZ()); }
 
      Vector2D forward2D() const
      { return nominalPos2DAtZ(forwardZ()); }
 
      Vector3D backward3D() const
      { return nominalPos3DAtZ(backwardZ()); }
 
      Vector2D backward2D() const
      { return nominalPos2DAtZ(backwardZ()); }
 
      Vector3D wireVector() const
      { return Vector3D(nominalMovePerZ() * deltaZ(), deltaZ()); }
 
      double forwardZ() const
      { return m_forwardZ; }
 
      double backwardZ() const
      { return m_backwardZ; }
 
      double deltaZ() const
      { return forwardZ() - backwardZ(); }
 
      double outOfZBoundsFactor(double z) const
      { return std::fmax(backwardZ() - z, z - forwardZ()) / deltaZ(); }
 
      double forwardPhi() const
      { return forward2D().phi(); }
 
      double backwardPhi() const
      { return backward2D().phi(); }
 
      double forwardCylindricalR() const
      { return forward2D().cylindricalR(); }
 
      double backwardCylindricalR() const
      { return backward2D().cylindricalR(); }
 
      double forwardPhiToRef() const
      { return forward2D().angleWith(refPos2D()); }
 
      double backwardPhiToRef() const
      { return backward2D().angleWith(refPos2D()); }
 
      /* backwardToForwardAngle means how far the backward position has to be rotated in the xy projection
         in the mathematical positiv sense that it seems to be coaligned with the forward position. */
      double backwardToForwardAngle() const
      { return backward2D().angleWith(forward2D()) ; }
 
      double tanLambda() const
      { return 1 / nominalMovePerZ().norm(); }
 
      double lambda() const
      { return std::atan(tanLambda()); }
 
      double tanTheta() const
      { return std::atan(nominalMovePerZ().norm()); }
 
      double theta() const
      { return std::atan(nominalMovePerZ().norm()); }
 
      double nominalPerigeeZ() const
      { return -refPos2D().dot(nominalMovePerZ()) / nominalMovePerZ().normSquared(); }
 
      Vector3D nominalPerigee3D() const
      { return nominalPos3DAtZ(nominalPerigeeZ()); }
 
      Vector2D nominalPerigee2D() const
      { return refPos2D().orthogonalVector(nominalMovePerZ()); }
 
      double refX() const
      { return m_refPos3D.x(); }
 
      double refY() const
      { return m_refPos3D.y(); }
 
      double refZ() const
      { return m_refPos3D.z(); }
 
      double refCylindricalRSquared() const
      { return m_refPos3D.cylindricalRSquared(); }
 
      const Vector2D& refPos2D() const
      { return m_refPos3D.xy(); }
 
      const Vector3D& refPos3D() const
      { return m_refPos3D; }
 
      double sagCoeff() const
      { return m_sagCoeff; }
 
    private:
      Vector3D m_refPos3D;
 
      Vector2D m_nominalMovePerZ;
 
      double m_forwardZ = 0.0;
 
      double m_backwardZ = 0.0;
 
      double m_sagCoeff = 0.0;
 
    };
  }
}