CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_6/src/PhysicsTools/HepMCCandAlgos/src/GenParticlesHelper.cc

Go to the documentation of this file.
00001 #include "PhysicsTools/HepMCCandAlgos/interface/GenParticlesHelper.h"
00002 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00003 
00004 
00005 using namespace reco;
00006 
00007 namespace GenParticlesHelper {
00008 
00009   void 
00010   findParticles(const reco::GenParticleCollection& sourceParticles, reco::GenParticleRefVector& particleRefs, int pdgId, int status ) {
00011 
00012     unsigned index = 0;
00013     for(IG ig = sourceParticles.begin(); 
00014         ig!= sourceParticles.end(); ++ig, ++index) {
00015  
00016       const GenParticle& gen = *ig;
00017     
00018       // status has been specified, and this one does not have the correct
00019       // status
00020       if(status && gen.status()!=status ) continue;
00021     
00022       if( std::abs(gen.pdgId()) == pdgId ) {
00023         GenParticleRef genref( &sourceParticles, index );
00024         particleRefs.push_back( genref );
00025       }
00026     }
00027   }
00028 
00029 
00030   void 
00031   findDescendents(const reco::GenParticleRef& base, 
00032                   reco::GenParticleRefVector& descendents, 
00033                   int status, int pdgId ) {
00034 
00035 
00036     const GenParticleRefVector& daughterRefs = base->daughterRefVector();
00037   
00038     for(IGR idr = daughterRefs.begin(); 
00039         idr!= daughterRefs.end(); ++idr ) {
00040     
00041       if( (*idr)->status() == status && 
00042           (!pdgId || std::abs((*idr)->pdgId()) == pdgId) ) {
00043       
00044         // cout<<"adding "<<(*idr)<<endl;
00045         descendents.push_back(*idr);
00046       }
00047       else 
00048         findDescendents( *idr, descendents, status, pdgId );
00049     }
00050   }
00051 
00052 
00053 
00054   void 
00055   findSisters(const reco::GenParticleRef& baseSister, 
00056               GenParticleRefVector& sisterRefs) {
00057   
00058     assert( baseSister->numberOfMothers() > 0 );
00059   
00060     // get first mother 
00061     const GenParticleRefVector& mothers = baseSister->motherRefVector();
00062   
00063     // get sisters 
00064     const GenParticleRefVector allRefs 
00065       = mothers[0]->daughterRefVector();
00066   
00067     typedef GenParticleRefVector::const_iterator IT;
00068     for(IT id = allRefs.begin();
00069         id != allRefs.end();
00070         ++id ) {
00071     
00072       if( *id == baseSister ) { 
00073         continue; // this is myself
00074       }
00075       else 
00076         sisterRefs.push_back( *id );
00077     }
00078   }
00079 
00080   bool 
00081   isDirect(const reco::GenParticleRef& particle) {
00082     assert( (particle->status() != 0) && (particle->status() < 4 ) );
00083     if( particle->status() == 3 )
00084       return true;
00085     else {
00086       assert( particle->numberOfMothers() > 0 );
00087   
00088       // get first mother 
00089       const GenParticleRefVector& mothers = particle->motherRefVector();
00090       if( mothers[0]->status() == 3 )
00091         return true;
00092       else
00093         return false;
00094     }
00095   }
00096 
00097 
00098   bool hasAncestor( const reco::GenParticle* particle,
00099                     int pdgId, int status )  {    
00100  
00101     if( particle->pdgId() == pdgId && 
00102         particle->status() == status )
00103       return true;
00104 
00105     const GenParticleRefVector& mothers = particle->motherRefVector();
00106     
00107     for( IGR im = mothers.begin(); im!=mothers.end(); ++im) {
00108       const GenParticle& part = **im;
00109       if( hasAncestor( &part, pdgId, status) )
00110         return true;
00111     } 
00112 
00113     return false;
00114   }
00115 
00116 
00117   std::ostream& operator<<( std::ostream& out, 
00118                             const reco::GenParticleRef& genRef ) {
00119     
00120     if(!out) return out;
00121     
00122     out<<genRef.key()<<" "<<genRef->pt();
00123     
00124     return out;
00125   }
00126 
00127 }