CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_4_5_patch3/src/RecoPixelVertexing/PixelTrackFitting/src/PixelFitterByConformalMappingAndLine.cc

Go to the documentation of this file.
00001 #include "RecoPixelVertexing/PixelTrackFitting/interface/PixelFitterByConformalMappingAndLine.h"
00002 
00003 #include "FWCore/Framework/interface/EventSetup.h"
00004 #include "FWCore/Framework/interface/ESHandle.h"
00005 
00006 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
00007 #include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h"
00008 
00009 
00010 #include "MagneticField/Engine/interface/MagneticField.h"
00011 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00012 #include "CommonTools/Statistics/interface/LinearFit.h"
00013 #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h"
00014 
00015 
00016 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHitBuilder.h"
00017 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
00018 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
00019 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00020 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00021 #include "TrackingTools/DetLayers/interface/DetLayer.h"
00022 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
00023 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
00024 
00025 #include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
00026 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
00027 
00028 #include "ConformalMappingFit.h"
00029 #include "RecoPixelVertexing/PixelTrackFitting/interface/PixelTrackBuilder.h"
00030 #include "RZLine.h"
00031 
00032 #include "RecoTracker/TkMSParametrization/interface/PixelRecoUtilities.h"
00033 #include "RecoTracker/TkMSParametrization/interface/LongitudinalBendingCorrection.h"
00034 
00035 using namespace std;
00036 
00037 template <class T> T sqr( T t) {return t*t;}
00038 
00039 
00040 PixelFitterByConformalMappingAndLine::PixelFitterByConformalMappingAndLine(
00041     const edm::ParameterSet& cfg) : theConfig(cfg)
00042 { }
00043 
00044 PixelFitterByConformalMappingAndLine::PixelFitterByConformalMappingAndLine()
00045 { }
00046 
00047 reco::Track* PixelFitterByConformalMappingAndLine::run(
00048     const edm::EventSetup& es,
00049     const std::vector<const TrackingRecHit * > & hits,
00050     const TrackingRegion & region) const
00051 {
00052 
00053   int nhits = hits.size();
00054 
00055   vector<GlobalPoint> points;
00056   vector<GlobalError> errors;
00057   vector<bool> isBarrel;
00058   
00059 
00060   edm::ESHandle<TrackerGeometry> tracker;
00061   es.get<TrackerDigiGeometryRecord>().get(tracker);
00062 
00063   edm::ESHandle<MagneticField> field;
00064   es.get<IdealMagneticFieldRecord>().get(field);
00065 
00066   edm::ESHandle<TransientTrackingRecHitBuilder> ttrhBuilder;
00067   string ttrhBuilderName = theConfig.getParameter<std::string>("TTRHBuilder");
00068   es.get<TransientRecHitRecord>().get( ttrhBuilderName, ttrhBuilder);
00069 
00070   for (vector<const TrackingRecHit*>::const_iterator ih=hits.begin();  ih!=hits.end(); ih++) {
00071     TransientTrackingRecHit::RecHitPointer recHit = ttrhBuilder->build(*ih);
00072     points.push_back( recHit->globalPosition() );
00073     errors.push_back( recHit->globalPositionError() );
00074     isBarrel.push_back( recHit->detUnit()->type().isBarrel() );
00075   }
00076 
00077 //    if (useMultScatt) {
00078 //      MultipleScatteringParametrisation ms(hits[i].layer());
00079 //      float cotTheta = (p.z()-zVtx)/p.perp();
00080 //      err += sqr( ms( pt, cotTheta, PixelRecoPointRZ(0.,zVtx) ) );
00081 //    }
00082 
00083   //
00084   // simple fit to get pt, phi0 used for precise calcul.
00085   //
00086   typedef ConformalMappingFit::PointXY PointXY;
00087   vector<PointXY> xy; vector<float> errRPhi2;
00088   for (int i=0; i < nhits; ++i) {
00089     const GlobalPoint & point = points[i];
00090     xy.push_back(PointXY( point.x()-region.origin().x(), point.y()-region.origin().y()));
00091     float phiErr2 = errors[i].phierr(point);
00092     errRPhi2.push_back( point.perp2()*phiErr2);
00093   }
00094   ConformalMappingFit parabola(xy, errRPhi2);
00095   if (theConfig.exists("fixImpactParameter")) 
00096     parabola.fixImpactParmaeter(theConfig.getParameter<double>("fixImpactParameter"));
00097   else if (nhits < 3) parabola.fixImpactParmaeter(0.);
00098  
00099 
00100   Measurement1D curv = parabola.curvature();
00101   float invPt = PixelRecoUtilities::inversePt( curv.value(), es);
00102   float valPt =  (invPt > 1.e-4) ? 1./invPt : 1.e4;
00103   float errPt =PixelRecoUtilities::inversePt(curv.error(), es) * sqr(valPt);
00104   Measurement1D pt (valPt,errPt);
00105   Measurement1D phi = parabola.directionPhi();
00106   Measurement1D tip = parabola.impactParameter();
00107 
00108   //
00109   // precalculate theta to correct errors:
00110   //
00111   vector<float> r(nhits),z(nhits),errZ(nhits);
00112   float simpleCot = ( points.back().z()-points.front().z() )/ (points.back().perp() - points.front().perp() );
00113   for (int i=0; i< nhits; ++i) {
00114     const GlobalPoint & point = points[i]; 
00115     const GlobalError & error = errors[i];
00116     r[i] = sqrt( sqr(point.x()-region.origin().x()) + sqr(point.y()-region.origin().y()) );
00117     r[i] += pixelrecoutilities::LongitudinalBendingCorrection(pt.value(),es)(r[i]);
00118     z[i] = point.z()-region.origin().z();  
00119     errZ[i] =  (isBarrel[i]) ? sqrt(error.czz()) : sqrt( error.rerr(point) )*simpleCot;
00120   }
00121 
00122   //
00123   // line fit (R-Z plane)
00124   //
00125   RZLine rzLine(r,z,errZ);
00126   float      cottheta, intercept, covss, covii, covsi;
00127   rzLine.fit(cottheta, intercept, covss, covii, covsi);
00128   
00129 //
00130 // parameters for track builder 
00131 //
00132   Measurement1D zip(intercept, sqrt(covii));
00133   Measurement1D cotTheta(cottheta, sqrt(covss));  
00134   float chi2 = parabola.chi2() +  rzLine.chi2(cottheta, intercept);
00135   int charge = parabola.charge();
00136 
00137 
00138   PixelTrackBuilder builder;
00139   return builder.build(pt, phi, cotTheta, tip, zip, chi2, charge, hits,  field.product(), region.origin());
00140 }
00141 
00142