PixelTrackFitting/src/PixelFitterByConformalMappingAndLine.cc

0001 #include "RecoTracker/PixelTrackFitting/interface/PixelFitterByConformalMappingAndLine.h"
0002
0003 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0004 #include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h"
0005
0006 #include "MagneticField/Engine/interface/MagneticField.h"
0007 #include "CommonTools/Statistics/interface/LinearFit.h"
0008 #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h"
0009
0010 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHitBuilder.h"
0011 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0012 #include "TrackingTools/DetLayers/interface/DetLayer.h"
0013 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0014
0015 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0016 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0017
0018 #include "ConformalMappingFit.h"
0019 #include "RecoTracker/PixelTrackFitting/interface/PixelTrackBuilder.h"
0020 #include "RecoTracker/PixelTrackFitting/interface/RZLine.h"
0021
0022 #include "RecoTracker/TkMSParametrization/interface/PixelRecoUtilities.h"
0023 #include "RecoTracker/TkMSParametrization/interface/LongitudinalBendingCorrection.h"
0024
0025 using namespace std;
0026
0027 template <class T>
0028 T sqr(T t) {
0029   return t * t;
0030 }
0031
0032 PixelFitterByConformalMappingAndLine::PixelFitterByConformalMappingAndLine(
0033     const TransientTrackingRecHitBuilder *ttrhBuilder,
0034     const TrackerGeometry *tracker,
0035     const MagneticField *field,
0036     double fixImpactParameter,
0037     bool useFixImpactParameter)
0038     : theTTRHBuilder(ttrhBuilder),
0039       theTracker(tracker),
0040       theField(field),
0041       theFixImpactParameter(fixImpactParameter),
0042       theUseFixImpactParameter(useFixImpactParameter) {}
0043
0044 std::unique_ptr<reco::Track> PixelFitterByConformalMappingAndLine::run(const std::vector<const TrackingRecHit *> &hits,
0045                                                                        const TrackingRegion &region) const {
0046   int nhits = hits.size();
0047
0048   vector<GlobalPoint> points;
0049   vector<GlobalError> errors;
0050   vector<bool> isBarrel;
0051
0052   for (vector<const TrackingRecHit *>::const_iterator ih = hits.begin(); ih != hits.end(); ih++) {
0053     TransientTrackingRecHit::RecHitPointer recHit = theTTRHBuilder->build(*ih);
0054     points.push_back(recHit->globalPosition());
0055     errors.push_back(recHit->globalPositionError());
0056     isBarrel.push_back(recHit->detUnit()->type().isBarrel());
0057   }
0058
0059   //    if (useMultScatt) {
0060   //      MultipleScatteringParametrisation ms(hits[i].layer());
0061   //      float cotTheta = (p.z()-zVtx)/p.perp();
0062   //      err += sqr( ms( pt, cotTheta, PixelRecoPointRZ(0.,zVtx) ) );
0063   //    }
0064
0065   //
0066   // simple fit to get pt, phi0 used for precise calcul.
0067   //
0068   typedef ConformalMappingFit::PointXY PointXY;
0069   vector<PointXY> xy;
0070   vector<float> errRPhi2;
0071   for (int i = 0; i < nhits; ++i) {
0072     const GlobalPoint &point = points[i];
0073     xy.push_back(PointXY(point.x() - region.origin().x(), point.y() - region.origin().y()));
0074     float phiErr2 = errors[i].phierr(point);
0075     errRPhi2.push_back(point.perp2() * phiErr2);
0076   }
0077   ConformalMappingFit parabola(xy, errRPhi2);
0078   if (theUseFixImpactParameter)
0079     parabola.fixImpactParmaeter(theFixImpactParameter);
0080   else if (nhits < 3)
0081     parabola.fixImpactParmaeter(0.);
0082
0083   Measurement1D curv = parabola.curvature();
0084   float invPt = PixelRecoUtilities::inversePt(curv.value(), *theField);
0085   float valPt = (invPt > 1.e-4) ? 1. / invPt : 1.e4;
0086   float errPt = PixelRecoUtilities::inversePt(curv.error(), *theField) * sqr(valPt);
0087   Measurement1D pt(valPt, errPt);
0088   Measurement1D phi = parabola.directionPhi();
0089   Measurement1D tip = parabola.impactParameter();
0090
0091   //
0092   // precalculate theta to correct errors:
0093   //
0094   vector<float> r(nhits), z(nhits), errZ(nhits);
0095   float simpleCot = (points.back().z() - points.front().z()) / (points.back().perp() - points.front().perp());
0096   for (int i = 0; i < nhits; ++i) {
0097     const GlobalPoint &point = points[i];
0098     const GlobalError &error = errors[i];
0099     r[i] = sqrt(sqr(point.x() - region.origin().x()) + sqr(point.y() - region.origin().y()));
0100     r[i] += pixelrecoutilities::LongitudinalBendingCorrection(pt.value(), *theField)(r[i]);
0101     z[i] = point.z() - region.origin().z();
0102     errZ[i] = (isBarrel[i]) ? sqrt(error.czz()) : sqrt(error.rerr(point)) * simpleCot;
0103   }
0104
0105   //
0106   // line fit (R-Z plane)
0107   //
0108   RZLine rzLine(r, z, errZ);
0109
0110   //
0111   // parameters for track builder
0112   //
0113   Measurement1D zip(rzLine.intercept(), sqrt(rzLine.covii()));
0114   Measurement1D cotTheta(rzLine.cotTheta(), sqrt(rzLine.covss()));
0115   float chi2 = parabola.chi2() + rzLine.chi2();
0116   int charge = parabola.charge();
0117
0118   PixelTrackBuilder builder;
0119   return std::unique_ptr<reco::Track>(
0120       builder.build(pt, phi, cotTheta, tip, zip, chi2, charge, hits, theField, region.origin()));
0121 }