CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/RecoPixelVertexing/PixelLowPtUtilities/src/PixelTripletLowPtGenerator.cc

Go to the documentation of this file.
00001 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/PixelTripletLowPtGenerator.h"
00002 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/ThirdHitPrediction.h"
00003 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/TripletFilter.h"
00004 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/HitInfo.h"
00005 
00006 #include "RecoTracker/TkMSParametrization/interface/PixelRecoPointRZ.h"
00007 
00008 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00009 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00010 #include "FWCore/Framework/interface/ESHandle.h"
00011 
00012 #include "DataFormats/SiPixelDetId/interface/PXBDetId.h"
00013 #include "DataFormats/SiPixelDetId/interface/PXFDetId.h"
00014 
00015 #undef Debug
00016 
00017 using namespace std;
00018 using namespace ctfseeding;
00019 
00020 /*****************************************************************************/
00021 void PixelTripletLowPtGenerator::init(const HitPairGenerator & pairs,
00022       const vector<SeedingLayer> & layers,
00023       LayerCacheType* layerCache)
00024 {
00025   thePairGenerator = pairs.clone();
00026   theLayers        = layers;
00027   theLayerCache    = layerCache;
00028 
00029   checkMultipleScattering = ps.getParameter<bool>("checkMultipleScattering");
00030   nSigMultipleScattering  = ps.getParameter<double>("nSigMultipleScattering");
00031   checkClusterShape       = ps.getParameter<bool>("checkClusterShape"); 
00032   rzTolerance             = ps.getParameter<double>("rzTolerance");
00033   maxAngleRatio           = ps.getParameter<double>("maxAngleRatio");
00034   builderName             = ps.getParameter<string>("TTRHBuilder");
00035 }
00036 
00037 /*****************************************************************************/
00038 void PixelTripletLowPtGenerator::getTracker
00039   (const edm::EventSetup& es)
00040 {
00041   if(theTracker == 0)
00042   {
00043     // Get tracker geometry
00044     edm::ESHandle<TrackerGeometry> tracker;
00045     es.get<TrackerDigiGeometryRecord>().get(tracker);
00046 
00047     theTracker = tracker.product();
00048   }
00049 
00050   if(theFilter == 0)
00051   {
00052     theFilter = new TripletFilter(es); 
00053   }
00054 }
00055 
00056 /*****************************************************************************/
00057 GlobalPoint PixelTripletLowPtGenerator::getGlobalPosition
00058   (const TrackingRecHit* recHit)
00059 {
00060   DetId detId = recHit->geographicalId();
00061 
00062   return
00063     theTracker->idToDet(detId)->toGlobal(recHit->localPosition());
00064 }
00065 
00066 /*****************************************************************************/
00067 void PixelTripletLowPtGenerator::hitTriplets(
00068     const TrackingRegion& region,
00069     OrderedHitTriplets & result,
00070     const edm::Event & ev,
00071     const edm::EventSetup& es) 
00072 {
00073   // Generate pairs
00074   OrderedHitPairs pairs; pairs.reserve(30000);
00075   thePairGenerator->hitPairs(region,pairs,ev,es);
00076 
00077   if (pairs.size() == 0) return;
00078 
00079   int size = theLayers.size(); 
00080 
00081   // Set aliases
00082   const RecHitsSortedInPhi **thirdHitMap = new const RecHitsSortedInPhi*[size]; 
00083   for(int il=0; il<size; il++)
00084     thirdHitMap[il] = &(*theLayerCache)(&theLayers[il], region, ev, es);
00085 
00086   // Get tracker
00087   getTracker(es);
00088 
00089   // Look at all generated pairs
00090   for(OrderedHitPairs::const_iterator ip = pairs.begin();
00091                                       ip!= pairs.end(); ip++)
00092   {
00093     // Fill rechits and points
00094     vector<const TrackingRecHit*> recHits(3);
00095     vector<GlobalPoint> points(3);
00096 
00097     recHits[0] = (*ip).inner()->hit();
00098     recHits[1] = (*ip).outer()->hit();
00099 
00100 #ifdef Debug
00101     cerr << " RecHits " + HitInfo::getInfo(*recHits[0]) +
00102                           HitInfo::getInfo(*recHits[1]) << endl;
00103 #endif
00104 
00105     for(int i=0; i<2; i++)
00106       points[i] = getGlobalPosition(recHits[i]);
00107 
00108     // Initialize helix prediction
00109     ThirdHitPrediction
00110       thePrediction(region,
00111                     points[0],points[1], es,
00112                     nSigMultipleScattering,maxAngleRatio,builderName);
00113 
00114     // Look at all layers
00115     for(int il=0; il<size; il++)
00116     {
00117       const SeedingLayer & layerwithhits = theLayers[il];
00118       const DetLayer * layer = layerwithhits.detLayer();
00119 
00120 #ifdef Debug
00121       cerr << "  check layer " << layer->subDetector()
00122                         << " " << layer->location() << endl;
00123 #endif
00124 
00125       // Get ranges for the third hit
00126       float phi[2],rz[2];
00127       thePrediction.getRanges(layer, phi,rz);
00128 
00129       PixelRecoRange<float> phiRange(phi[0]              , phi[1]             );
00130       PixelRecoRange<float>  rzRange( rz[0] - rzTolerance, rz[1] + rzTolerance);
00131 
00132       // Get third hit candidates from cache
00133       typedef RecHitsSortedInPhi::Hit Hit;
00134       vector<Hit> thirdHits = thirdHitMap[il]->hits(phiRange.min(),phiRange.max());
00135       typedef vector<Hit>::const_iterator IH;
00136 
00137       for (IH th=thirdHits.begin(), eh=thirdHits.end(); th < eh; ++th) 
00138       {
00139         // Fill rechit and point
00140         recHits[2] = (*th)->hit();
00141         points[2]  = getGlobalPosition(recHits[2]);
00142 
00143 #ifdef Debug
00144         cerr << "  third hit " + HitInfo::getInfo(*recHits[2]) << endl;
00145 #endif
00146 
00147         // Check if third hit is compatible with multiple scattering
00148         vector<GlobalVector> globalDirs;
00149         if(thePrediction.isCompatibleWithMultipleScattering
00150              (points[2], recHits, globalDirs, es) == false)
00151         {
00152 #ifdef Debug
00153           cerr << "  not compatible: multiple scattering" << endl;
00154 #endif
00155           if(checkMultipleScattering) continue;
00156         }
00157 
00158         // Convert to localDirs
00159 /*
00160         vector<LocalVector> localDirs;
00161         vector<GlobalVector>::const_iterator globalDir = globalDirs.begin();
00162         for(vector<const TrackingRecHit *>::const_iterator
00163                                             recHit  = recHits.begin();
00164                                             recHit != recHits.end(); recHit++)
00165         {
00166           localDirs.push_back(theTracker->idToDet(
00167                              (*recHit)->geographicalId())->toLocal(*globalDir));
00168           globalDir++;
00169         }
00170 */
00171 
00172         // Check if the cluster shapes are compatible with thrusts
00173         if(checkClusterShape)
00174         {
00175           if(! theFilter->checkTrack(recHits,globalDirs))
00176           {
00177 #ifdef Debug
00178             cerr << "  not compatible: cluster shape" << endl;
00179 #endif
00180             continue;
00181           }
00182         }
00183 
00184         // All checks passed, put triplet back
00185         result.push_back(OrderedHitTriplet((*ip).inner(),(*ip).outer(),*th));
00186       }
00187     }
00188   } 
00189   delete [] thirdHitMap;
00190 
00191   return;
00192 }
00193 
00194