CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_3/src/RecoPixelVertexing/PixelTrackFitting/src/PixelTrackCleanerBySharedHits.cc

Go to the documentation of this file.
00001 #include "RecoPixelVertexing/PixelTrackFitting/interface/PixelTrackCleanerBySharedHits.h"
00002 
00003 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
00004 #include "DataFormats/TrackReco/interface/Track.h"
00005 #include "DataFormats/SiPixelDetId/interface/PXBDetId.h"
00006 #include "DataFormats/SiPixelDetId/interface/PXFDetId.h"
00007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00008 
00009 using namespace std;
00010 using namespace reco;
00011 using namespace pixeltrackfitting;
00012 
00013 PixelTrackCleanerBySharedHits::PixelTrackCleanerBySharedHits( const edm::ParameterSet& cfg)
00014 {}
00015 
00016 PixelTrackCleanerBySharedHits::~PixelTrackCleanerBySharedHits()
00017 {}
00018 
00019 TracksWithRecHits PixelTrackCleanerBySharedHits::cleanTracks(const TracksWithRecHits & trackHitPairs)
00020 {
00021   typedef std::vector<const TrackingRecHit *> RecHits;
00022   trackOk.clear();
00023 
00024   LogDebug("PixelTrackCleanerBySharedHits") << "Cleanering tracks" << "\n";
00025   int size = trackHitPairs.size();
00026   for (int i = 0; i < size; i++) trackOk.push_back(true);
00027 
00028   for (iTrack1 = 0; iTrack1 < size; iTrack1++)
00029   {
00030     track1 = trackHitPairs.at(iTrack1).first;
00031     const RecHits& recHits1 = trackHitPairs.at(iTrack1).second;
00032 
00033     if (!trackOk.at(iTrack1)) continue;
00034 
00035     for (iTrack2 = iTrack1 + 1; iTrack2 < size; iTrack2++)
00036     {
00037       if (!trackOk.at(iTrack1) || !trackOk.at(iTrack2)) continue;
00038 
00039       track2 = trackHitPairs.at(iTrack2).first;
00040       const RecHits& recHits2 = trackHitPairs.at(iTrack2).second;
00041 
00042       int commonRecHits = 0;
00043       for (int iRecHit1 = 0; iRecHit1 < (int)recHits1.size(); iRecHit1++)
00044       {
00045         for (int iRecHit2 = 0; iRecHit2 < (int)recHits2.size(); iRecHit2++)
00046         {
00047           if (recHitsAreEqual(recHits1.at(iRecHit1), recHits2.at(iRecHit2))) commonRecHits++;
00048         }
00049       }
00050       if (commonRecHits > 1) cleanTrack();
00051     }
00052   }
00053 
00054   vector<TrackWithRecHits> cleanedTracks;
00055 
00056   for (int i = 0; i < size; i++)
00057   {
00058     if (trackOk.at(i)) cleanedTracks.push_back(trackHitPairs.at(i));
00059     else delete trackHitPairs.at(i).first;
00060   }
00061   return cleanedTracks;
00062 }
00063 
00064 
00065 void PixelTrackCleanerBySharedHits::cleanTrack()
00066 {
00067   if (track1->pt() > track2->pt()) trackOk.at(iTrack2) = false;
00068   else trackOk.at(iTrack1) = false;
00069 }
00070 
00071 
00072 bool PixelTrackCleanerBySharedHits::recHitsAreEqual(const TrackingRecHit *recHit1, const TrackingRecHit *recHit2)
00073 {
00074   if (recHit1->geographicalId() != recHit2->geographicalId()) return false;
00075   LocalPoint pos1 = recHit1->localPosition();
00076   LocalPoint pos2 = recHit2->localPosition();
00077   return ((pos1.x() == pos2.x()) && (pos1.y() == pos2.y()));
00078 }