21 return compute(chargeA, momentumA, positionA, chargeB, momentumB, positionB);
39 return compute(chargeA, momentumA, positionA, chargeB, momentumB, positionB);
43 pair<GlobalPoint, GlobalPoint> ClosestApproachInRPhi::points()
const
46 throw cms::Exception(
"TrackingTools/PatternTools",
"ClosestApproachInRPhi::could not compute track crossing. Check status before calling this method!");
47 return pair<GlobalPoint, GlobalPoint> (posA, posB);
52 ClosestApproachInRPhi::crossingPoint()
const
55 throw cms::Exception(
"TrackingTools/PatternTools",
"ClosestApproachInRPhi::could not compute track crossing. Check status before calling this method!");
56 return GlobalPoint(0.5*(posA.basicVector() + posB.basicVector()));
61 float ClosestApproachInRPhi::distance()
const
64 throw cms::Exception(
"TrackingTools/PatternTools",
"ClosestApproachInRPhi::could not compute track crossing. Check status before calling this method!");
65 return (posB - posA).mag();
80 circleParameters(chargeA, momentumA, positionA, xca, yca, ra, bz);
82 circleParameters(chargeB, momentumB, positionB, xcb, ycb, rb, bz);
85 double xg1, yg1, xg2, yg2;
86 int flag = transverseCoord(xca, yca, ra, xcb, ycb, rb, xg1, yg1, xg2, yg2);
92 double xga, yga, zga, xgb, ygb, zgb;
97 double za1 = zCoord(momentumA, positionA, ra, xca, yca, xg1, yg1);
98 double zb1 = zCoord(momentumB, positionB, rb, xcb, ycb, xg1, yg1);
99 double za2 = zCoord(momentumA, positionA, ra, xca, yca, xg2, yg2);
100 double zb2 = zCoord(momentumB, positionB, rb, xcb, ycb, xg2, yg2);
102 if (
abs(zb1 - za1) <
abs(zb2 - za2)) {
103 xga = xg1; yga = yg1; zga = za1; zgb = zb1;
106 xga = xg2; yga = yg2; zga = za2; zgb = zb2;
108 xgb = xga; ygb = yga;
112 xga = xg1; yga = yg1;
113 zga = zCoord(momentumA, positionA, ra, xca, yca, xga, yga);
114 xgb = xg2; ygb = yg2;
115 zgb = zCoord(momentumB, positionB, rb, xcb, ycb, xgb, ygb);
124 pair <GlobalTrajectoryParameters, GlobalTrajectoryParameters>
125 ClosestApproachInRPhi::trajectoryParameters ()
const
128 throw cms::Exception(
"TrackingTools/PatternTools",
"ClosestApproachInRPhi::could not compute track crossing. Check status before calling this method!");
129 pair <GlobalTrajectoryParameters, GlobalTrajectoryParameters>
130 ret ( newTrajectory( posA, paramA, bz),
131 newTrajectory( posB, paramB, bz) );
139 double qob = oldgtp.
charge()/bz;
144 double npx = (newpt.
y()-yc)*(bz/oldgtp.
charge());
145 double npy = (xc-newpt.
x())*(bz/oldgtp.
charge());
177 double& xc,
double& yc,
double&
r,
188 double qob = charge/bz;
198 xc = position.
x() + qob * momentum.
y();
199 yc = position.
y() - qob * momentum.
x();
205 ClosestApproachInRPhi::transverseCoord(
double cxa,
double cya,
double ra,
206 double cxb,
double cyb,
double rb,
207 double & xg1,
double & yg1,
208 double & xg2,
double & yg2)
211 double x1, y1, x2, y2;
216 double d_ab =
sqrt((cxb - cxa)*(cxb - cxa) + (cyb - cya)*(cyb - cya));
221 double u = (cxb - cxa) / d_ab;
222 double v = (cyb - cya) / d_ab;
225 if (d_ab <= ra + rb && d_ab >=
abs(rb - ra)) {
231 double cosphi = (ra*ra - rb*rb + d_ab*d_ab) / (2*ra*d_ab);
232 double sinphi2 = 1. - cosphi*cosphi;
233 if (sinphi2 < 0.) { sinphi2 = 0.; cosphi = 1.; }
236 double sinphi =
sqrt(sinphi2);
237 x1 = ra*cosphi; y1 = ra*sinphi; x2 = x1; y2 = -y1;
239 else if (d_ab > ra + rb) {
246 x1 = ra; y1 = 0; x2 = d_ab - rb; y2 = 0;
248 else if (d_ab <
abs(rb - ra)) {
256 if (ra <= rb) sign = -1.;
257 x1 = sign*ra; y1 = 0; x2 = d_ab + sign*rb; y2 = 0;
264 xg1 = u*x1 - v*y1 + cxa; yg1 = v*x1 + u*y1 + cya;
265 xg2 = u*x2 - v*y2 + cxa; yg2 = v*x2 + u*y2 + cya;
274 double r,
double xc,
double yc,
275 double xg,
double yg)
279 double x = pos.
x();
double y = pos.
y();
double z = pos.
z();
281 double px = mom.
x();
double py = mom.
y();
double pz = mom.
z();
287 double sinHalfPhi =
sqrt((x-xg)*(x-xg) + (y-yg)*(y-yg))/(2*
r);
288 if (sinHalfPhi < 0.383) {
289 phi = 2*asin(sinHalfPhi);
292 double cosPhi = ((x-xc)*(xg-xc) + (y-yc)*(yg-yc))/(r*r);
293 if (
std::abs(cosPhi) > 1) cosPhi = (cosPhi > 0 ? 1 : -1);
294 phi =
abs(acos(cosPhi));
297 double signPhi = ((x - xc)*(yg - yc) - (xg - xc)*(y - yc) > 0) ? 1. : -1.;
301 double signOmega = ((x - xc)*py - (y - yc)*px > 0) ? 1. : -1.;
306 double dz = signPhi*signOmega*(pz/mom.
transverse())*phi*r;
TrackCharge charge() const
const GlobalTrajectoryParameters & parameters() const
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
Global3DPoint GlobalPoint
GlobalPoint globalPosition() const
TrackCharge charge() const
static int position[TOTALCHAMBERS][3]
FreeTrajectoryState const * freeState(bool withErrors=true) const
GlobalVector momentum() const
Abs< T >::type abs(const T &t)
GlobalVector momentum() const
GlobalPoint position() const
GlobalPoint position() const
const GlobalTrajectoryParameters & globalParameters() const
GlobalVector globalMomentum() const
const MagneticField & magneticField() const
TrackCharge charge() const