PixelTrackFitting/src/ConformalMappingFit.h

0001 #ifndef RecoTracker_PixelTrackFitting_src_ConformalMappingFit_h
0002 #define RecoTracker_PixelTrackFitting_src_ConformalMappingFit_h
0003
0004 #include "DataFormats/GeometryVector/interface/Basic2DVector.h"
0005 #include "DataFormats/GeometryVector/interface/Basic3DVector.h"
0006 #include "DataFormats/GeometrySurface/interface/TkRotation.h"
0007 #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h"
0008
0009 #include "ParabolaFit.h"
0010 #include <vector>
0011
0012 namespace edm {
0013   class ParameterSet;
0014 }
0015
0016 class ConformalMappingFit {
0017 public:
0018   typedef TkRotation<double> Rotation;
0019   typedef Basic2DVector<double> PointXY;
0020
0021   ConformalMappingFit(const std::vector<PointXY> &hits,
0022                       const std::vector<float> &errRPhi2,
0023                       const Rotation *rotation = nullptr);
0024
0025   ~ConformalMappingFit();
0026
0027   Measurement1D curvature() const;
0028   Measurement1D directionPhi() const;
0029   Measurement1D impactParameter() const;
0030
0031   int charge() const;
0032   double chi2() const { return theFit.chi2(); }
0033
0034   const Rotation *rotation() const { return theRotation; }
0035
0036   void fixImpactParmaeter(double ip) { theFit.fixParC(ip); }
0037   void skipErrorCalculation() { theFit.skipErrorCalculationByDefault(); }
0038
0039 private:
0040   double phiRot() const;
0041   void findRot(const PointXY &);
0042
0043 private:
0044   const Rotation *theRotation;
0045   bool myRotation;
0046   ParabolaFit theFit;
0047
0048   template <class T>
0049   class MappedPoint {
0050   public:
0051     typedef Basic2DVector<T> PointXY;
0052     MappedPoint() : theU(0), theV(0), theW(0), pRot(0) {}
0053     MappedPoint(const T &aU, const T &aV, const T &aWeight, const TkRotation<T> *aRot)
0054         : theU(aU), theV(aV), theW(aWeight), pRot(aRot) {}
0055     MappedPoint(const PointXY &point, const T &weight, const TkRotation<T> *aRot) : pRot(aRot) {
0056       T radius2 = point.mag2();
0057       Basic3DVector<T> rotated = (*pRot) * point;
0058       theU = rotated.x() / radius2;
0059       theV = rotated.y() / radius2;
0060       theW = weight * radius2 * radius2;
0061     }
0062     T u() const { return theU; }
0063     T v() const { return theV; }
0064     T weight() const { return theW; }
0065     PointXY unmap() const {
0066       T invRadius2 = theU * theU + theV * theV;
0067       Basic3DVector<T> tmp = (*pRot).multiplyInverse(Basic2DVector<T>(theU, theV));
0068       return PointXY(tmp.x() / invRadius2, tmp.y() / invRadius2);
0069     }
0070     T unmappedWeight() const {
0071       T invRadius2 = theU * theU + theV * theV;
0072       return theW * invRadius2 * invRadius2;
0073     }
0074
0075   private:
0076     T theU, theV, theW;
0077     const TkRotation<T> *pRot;
0078   };
0079 };
0080
0081 #endif