CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch12/src/RecoTracker/TkSeedGenerator/src/SeedFromConsecutiveHitsCreator.cc

Go to the documentation of this file.
00001 #include "RecoTracker/TkSeedGenerator/interface/SeedFromConsecutiveHitsCreator.h"
00002 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegion.h"
00003 
00004 #include "TrackingTools/KalmanUpdators/interface/KFUpdator.h"
00005 #include "RecoTracker/TkSeedGenerator/interface/FastHelix.h"
00006 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00007 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00009 #include "MagneticField/Engine/interface/MagneticField.h"
00010 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00011 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h" 
00012 #include "TrackingTools/Records/interface/TransientRecHitRecord.h" 
00013 #include "TrackingTools/MaterialEffects/interface/PropagatorWithMaterial.h"
00014 #include "FWCore/Framework/interface/ESHandle.h"
00015 #include "TrackingTools/GeomPropagators/interface/PropagationExceptions.h"
00016 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
00017 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00018 
00019 template <class T> T sqr( T t) {return t*t;}
00020 
00021 const TrajectorySeed * SeedFromConsecutiveHitsCreator::trajectorySeed(
00022     TrajectorySeedCollection & seedCollection,
00023     const SeedingHitSet & hits,
00024     const TrackingRegion & region,
00025     const edm::EventSetup& es)
00026 {
00027   if ( hits.size() < 2) return 0;
00028 
00029   GlobalTrajectoryParameters kine = initialKinematic(hits, region, es);
00030   float sinTheta = sin(kine.momentum().theta());
00031 
00032   CurvilinearTrajectoryError error = initialError(region,  sinTheta);
00033   FreeTrajectoryState fts(kine, error);
00034 
00035   return buildSeed(seedCollection,hits,fts,es); 
00036 }
00037 
00038 
00039 GlobalTrajectoryParameters SeedFromConsecutiveHitsCreator::initialKinematic(
00040       const SeedingHitSet & hits, 
00041       const TrackingRegion & region, 
00042       const edm::EventSetup& es) const
00043 {
00044   edm::ESHandle<MagneticField> bfield;
00045   es.get<IdealMagneticFieldRecord>().get(bfield);
00046 
00047   GlobalTrajectoryParameters kine;
00048 
00049   TransientTrackingRecHit::ConstRecHitPointer tth1 = hits[0];
00050   TransientTrackingRecHit::ConstRecHitPointer tth2 = hits[1];
00051   const GlobalPoint& vertexPos = region.origin();
00052 
00053   FastHelix helix(tth2->globalPosition(), tth1->globalPosition(), vertexPos, es);
00054   if (helix.isValid()) {
00055     kine = helix.stateAtVertex().parameters();
00056   } else {
00057     GlobalVector initMomentum(tth2->globalPosition() - vertexPos);
00058     initMomentum *= (100./initMomentum.perp()); 
00059     kine = GlobalTrajectoryParameters(vertexPos, initMomentum, 1, &*bfield);
00060   } 
00061 
00062   bool isBOFF = ( std::abs(bfield->inTesla(GlobalPoint(0,0,0)).z()) < 1e-3 );
00063   if (isBOFF && (theBOFFMomentum > 0)) {
00064     kine = GlobalTrajectoryParameters(kine.position(),
00065                               kine.momentum().unit() * theBOFFMomentum,
00066                               kine.charge(),
00067                               &*bfield);
00068   }
00069   return kine;
00070 }
00071 
00072 
00073 
00074 CurvilinearTrajectoryError SeedFromConsecutiveHitsCreator::
00075    initialError( const TrackingRegion& region, float sinTheta) const
00076 {
00077   // Set initial uncertainty on track parameters, using only P.V. constraint and no hit
00078   // information.
00079   GlobalError vertexErr( sqr(region.originRBound()), 0, sqr(region.originRBound()),
00080                    0, 0, sqr(region.originZBound()));
00081 
00082   float ptMin = region.ptMin();
00083 
00084 
00085   AlgebraicSymMatrix C(5,1);
00086 
00087 // FIXME: minC00. Prevent apriori uncertainty in 1/P from being too small, 
00088 // to avoid instabilities.
00089 // N.B. This parameter needs optimising ...
00090   float sin2th = sqr(sinTheta);
00091   float minC00 = 1.0;
00092   C[0][0] = std::max(sin2th/sqr(ptMin), minC00);
00093   float zErr = vertexErr.czz();
00094   float transverseErr = vertexErr.cxx(); // assume equal cxx cyy
00095   C[3][3] = transverseErr;
00096   C[4][4] = zErr*sin2th + transverseErr*(1-sin2th);
00097 
00098   return CurvilinearTrajectoryError(C);
00099 }
00100 
00101 const TrajectorySeed * SeedFromConsecutiveHitsCreator::buildSeed(
00102     TrajectorySeedCollection & seedCollection,
00103     const SeedingHitSet & hits,
00104     const FreeTrajectoryState & fts,
00105     const edm::EventSetup& es) const
00106 {
00107   // get tracker
00108   edm::ESHandle<TrackerGeometry> tracker;
00109   es.get<TrackerDigiGeometryRecord>().get(tracker);
00110   
00111   // get propagator
00112   edm::ESHandle<Propagator>  propagatorHandle;
00113   es.get<TrackingComponentsRecord>().get(thePropagatorLabel, propagatorHandle);
00114   const Propagator*  propagator = &(*propagatorHandle);
00115   
00116   // get updator
00117   KFUpdator  updator;
00118   
00119   // Now update initial state track using information from seed hits.
00120   
00121   TrajectoryStateOnSurface updatedState;
00122   edm::OwnVector<TrackingRecHit> seedHits;
00123   
00124   const TrackingRecHit* hit = 0;
00125   for ( unsigned int iHit = 0; iHit < hits.size(); iHit++) {
00126     hit = hits[iHit]->hit();
00127     TrajectoryStateOnSurface state = (iHit==0) ? 
00128       propagator->propagate(fts,tracker->idToDet(hit->geographicalId())->surface())
00129       : propagator->propagate(updatedState, tracker->idToDet(hit->geographicalId())->surface());
00130     if (!state.isValid()) return 0;
00131     
00132     TransientTrackingRecHit::ConstRecHitPointer tth = hits[iHit]; 
00133     
00134     TransientTrackingRecHit::RecHitPointer newtth = refitHit( tth, state);
00135     
00136     if (!checkHit(state,newtth,es)) return 0;
00137 
00138     updatedState =  updator.update(state, *newtth);
00139     if (!updatedState.isValid()) return 0;
00140     
00141     seedHits.push_back(newtth->hit()->clone());
00142   } 
00143 
00144   TrajectoryStateTransform transformer;
00145   boost::shared_ptr<PTrajectoryStateOnDet> PTraj(
00146       transformer.persistentState(updatedState, hit->geographicalId().rawId()));
00147   
00148   seedCollection.push_back( TrajectorySeed(*PTraj,seedHits,alongMomentum));
00149   return &seedCollection.back();
00150 }
00151 
00152 TransientTrackingRecHit::RecHitPointer SeedFromConsecutiveHitsCreator::refitHit(
00153       const TransientTrackingRecHit::ConstRecHitPointer &hit, 
00154       const TrajectoryStateOnSurface &state) const
00155 {
00156   return hit->clone(state);
00157 }