CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_9_patch3/src/CommonTools/ParticleFlow/src/PFPileUpAlgo.cc

Go to the documentation of this file.
00001 #include "CommonTools/ParticleFlow/interface/PFPileUpAlgo.h"
00002 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
00003 #include "DataFormats/VertexReco/interface/Vertex.h"
00004 
00005 void PFPileUpAlgo::process(const reco::PFCandidateCollection & pfCandidates, 
00006                            const reco::VertexCollection & vertices, 
00007                            const edm::Handle<reco::PFCandidateCollection> * handle)  {
00008 
00009   pfCandidatesFromVtx_.clear();
00010   pfCandidatesFromPU_.clear();
00011 
00012   for( unsigned i=0; i<pfCandidates.size(); i++ ) {
00013     
00014     const reco::PFCandidate& cand = pfCandidates[i];
00015     
00016     int ivertex;
00017 
00018     switch( cand.particleId() ) {
00019     case reco::PFCandidate::h:
00020       ivertex = chargedHadronVertex( vertices, cand );
00021       break;
00022     default:
00023       continue;
00024     } 
00025     
00026     // no associated vertex, or primary vertex
00027     // not pile-up
00028     if( ivertex == -1  || 
00029         ivertex == 0 ) {
00030       pfCandidatesFromVtx_.push_back(cand);
00031       pfCandidatesFromVtx_.back().setSourceCandidatePtr( reco::PFCandidatePtr(*handle,i) );
00032     } else {
00033       // associated to a vertex
00034       pfCandidatesFromPU_.push_back(cand);
00035       if(handle) {
00036         pfCandidatesFromPU_.back().setSourceCandidatePtr( reco::PFCandidatePtr(*handle,i) );
00037       }
00038     }
00039   }
00040 }
00041 
00042 
00043 int 
00044 PFPileUpAlgo::chargedHadronVertex( const reco::VertexCollection& vertices, const reco::PFCandidate& pfcand ) const {
00045 
00046   
00047   reco::TrackBaseRef trackBaseRef( pfcand.trackRef() );
00048   
00049   size_t  iVertex = 0;
00050   unsigned index=0;
00051   unsigned nFoundVertex = 0;
00052   typedef reco::VertexCollection::const_iterator IV;
00053   typedef reco::Vertex::trackRef_iterator IT;
00054   float bestweight=0;
00055   for(IV iv=vertices.begin(); iv!=vertices.end(); ++iv, ++index) {
00056 
00057     const reco::Vertex& vtx = *iv;
00058     
00059     // loop on tracks in vertices
00060     for(IT iTrack=vtx.tracks_begin(); 
00061         iTrack!=vtx.tracks_end(); ++iTrack) {
00062          
00063       const reco::TrackBaseRef& baseRef = *iTrack;
00064 
00065       // one of the tracks in the vertex is the same as 
00066       // the track considered in the function
00067       if(baseRef == trackBaseRef ) {
00068         float w = vtx.trackWeight(baseRef);
00069         //select the vertex for which the track has the highest weight
00070         if (w > bestweight){
00071           bestweight=w;
00072           iVertex=index;
00073           nFoundVertex++;
00074         }               
00075       }
00076     }
00077   }
00078 
00079   if (nFoundVertex>0){
00080     if (nFoundVertex!=1)
00081       edm::LogWarning("TrackOnTwoVertex")<<"a track is shared by at least two verteces. Used to be an assert";
00082     return iVertex;
00083   }
00084   // no vertex found with this track. 
00085 
00086   // optional: as a secondary solution, associate the closest vertex in z
00087   if ( checkClosestZVertex_ ) {
00088 
00089     double dzmin = 10000;
00090     double ztrack = pfcand.vertex().z();
00091     bool foundVertex = false;
00092     index = 0;
00093     for(IV iv=vertices.begin(); iv!=vertices.end(); ++iv, ++index) {
00094 
00095       double dz = fabs(ztrack - iv->z());
00096       if(dz<dzmin) {
00097         dzmin = dz; 
00098         iVertex = index;
00099         foundVertex = true;
00100       }
00101     }
00102 
00103     if( foundVertex ) 
00104       return iVertex;  
00105 
00106   }
00107 
00108 
00109   return -1 ;
00110 }
00111