PixelTrackFitting/src/ConformalMappingFit.cc

0001 #include "ConformalMappingFit.h"
0002
0003 using namespace std;
0004
0005 template <class T>
0006 T sqr(T t) {
0007   return t * t;
0008 }
0009
0010 ConformalMappingFit::ConformalMappingFit(const std::vector<PointXY>& hits,
0011                                          const std::vector<float>& errRPhi2,
0012                                          const Rotation* rot)
0013     : theRotation(rot), myRotation(rot == nullptr) {
0014   typedef ConformalMappingFit::MappedPoint<double> PointUV;
0015   int hits_size = hits.size();
0016   for (int i = 0; i < hits_size; i++) {
0017     if (!theRotation)
0018       findRot(hits[i]);
0019     PointUV point(hits[i], 1. / errRPhi2[i], theRotation);
0020     theFit.addPoint(point.u(), point.v(), point.weight());
0021   }
0022 }
0023
0024 void ConformalMappingFit::findRot(const PointXY& p) {
0025   myRotation = true;
0026   typedef Rotation::GlobalVector GlobalVector;  // ::GlobalVector is float!
0027   GlobalVector aX = GlobalVector(p.x(), p.y(), 0.).unit();
0028   GlobalVector aY(-aX.y(), aX.x(), 0.);
0029   GlobalVector aZ(0., 0., 1.);
0030   theRotation = new Rotation(aX, aY, aZ);
0031 }
0032
0033 ConformalMappingFit::~ConformalMappingFit() {
0034   if (myRotation)
0035     delete theRotation;
0036 }
0037
0038 double ConformalMappingFit::phiRot() const { return atan2(theRotation->xy(), theRotation->xx()); }
0039
0040 Measurement1D ConformalMappingFit::curvature() const {
0041   double val = fabs(2. * theFit.parA());
0042   double err = 2. * sqrt(theFit.varAA());
0043   return Measurement1D(val, err);
0044 }
0045
0046 Measurement1D ConformalMappingFit::directionPhi() const {
0047   double val = phiRot() + atan(theFit.parB());
0048   double err = sqrt(theFit.varBB());
0049   return Measurement1D(val, err);
0050 }
0051
0052 Measurement1D ConformalMappingFit::impactParameter() const {
0053   double val = -theFit.parC();
0054   double err = sqrt(theFit.varCC());
0055   return Measurement1D(val, err);
0056 }
0057
0058 int ConformalMappingFit::charge() const { return (theFit.parA() > 0.) ? -1 : 1; }