11 inline T
sqr(
const T&
t) {
return t*
t;}
23 const double sraightLineCutoff = 1.e-7;
24 if (fabs(rho)*R < sraightLineCutoff &&
25 fabs(rho)*startingPos.
perp() < sraightLineCutoff) {
28 cyl.
toLocal(startingDir), propDir);
29 std::pair<bool,double> pl = slc.
pathLength( cyl);
41 double pt = startingDir.
perp();
42 Point center( startingPos.
x()-startingDir.
y()/(pt*
rho),
43 startingPos.
y()+startingDir.
x()/(pt*
rho));
44 double p2 = startingPos.
perp2();
47 if (fabs(center.
x()) > fabs(center.
y())) {
49 E = (R2cyl -
p2) / (2.*center.
x());
50 F = center.
y()/center.
x();
51 B = 2.*( startingPos.
y() - F*startingPos.
x() - E*F);
52 C = 2.*E*startingPos.
x() + E*E + p2 - R2cyl;
56 E = (R2cyl -
p2) / (2.*center.
y());
57 F = center.
x()/center.
y();
58 B = 2.*( startingPos.
x() - F*startingPos.
y() - E*F);
59 C = 2.*E*startingPos.
y() + E*E + p2 - R2cyl;
81 double ipabs = 1./startingDir.
mag();
82 double sinTheta = pt * ipabs;
83 double cosTheta = startingDir.
z() * ipabs;
86 double tmp = 0.5 * dMag *
rho;
87 if (
std::abs(tmp)>1.) tmp = ::copysign(1.,tmp);
90 startingPos.
y() +
theD.
y(),
91 startingPos.
z() +
theS*cosTheta);
94 double sinPhi = 2.*tmp*
sqrt(1.-tmp*tmp);
95 double cosPhi = 1.-2.*tmp*
tmp;
97 startingDir.
x()*sinPhi+startingDir.
y()*cosPhi,
106 double momProj1 = startingDir.
x()*d1.
x() + startingDir.
y()*d1.
y();
107 double momProj2 = startingDir.
x()*d2.
x() + startingDir.
y()*d2.
y();
122 if (momProj1*momProj2 < 0) {
125 theD = (momProj1*propSign > 0) ? d1 : d2;
128 else if (momProj1*propSign > 0) {
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
void chooseSolution(const Point &p1, const Point &p2, const PositionType &startingPos, const DirectionType &startingDir, PropagationDirection propDir)
GlobalVector DirectionType
std::pair< bool, double > pathLength(const Cylinder &cyl) const
Global3DPoint GlobalPoint
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Scalar radius() const
Radius of the cylinder.
Basic2DVector< TmpType > Point
LocalPoint toLocal(const GlobalPoint &gp) const
HelixBarrelCylinderCrossing(const GlobalPoint &startingPos, const GlobalVector &startingDir, double rho, PropagationDirection propDir, const Cylinder &cyl)
T y() const
Cartesian y coordinate.
T mag2() const
The vector magnitude squared. Equivalent to vec.dot(vec)
std::vector< std::vector< double > > tmp
Square< F >::type sqr(const F &f)
T x() const
Cartesian x coordinate.