CMS 3D CMS Logo

Public Member Functions | Private Member Functions | Private Attributes

BCToEFilterAlgo Class Reference

#include <BCToEFilterAlgo.h>

List of all members.

Public Member Functions

 BCToEFilterAlgo (const edm::ParameterSet &)
bool filter (const edm::Event &iEvent)
bool hasBCAncestors (reco::GenParticle gp)
 ~BCToEFilterAlgo ()

Private Member Functions

bool isBCBaryon (reco::GenParticle gp)
bool isBCHadron (reco::GenParticle gp)
bool isBCMeson (reco::GenParticle gp)

Private Attributes

float eTThreshold_
float FILTER_ETA_MAX_
edm::InputTag genParSource_

Detailed Description

BCToEFilterAlgo returns true for events that have an electron, above configurable eT threshold and within |eta|<2.5, that has an ancestor of a b or c quark

Author:
J Lamb, UCSB

Definition at line 28 of file BCToEFilterAlgo.h.


Constructor & Destructor Documentation

BCToEFilterAlgo::BCToEFilterAlgo ( const edm::ParameterSet iConfig)

Definition at line 7 of file BCToEFilterAlgo.cc.

References edm::ParameterSet::getParameter().

                                                               { 

  //set constants
  FILTER_ETA_MAX_=2.5;
  eTThreshold_=(float)iConfig.getParameter<double>("eTThreshold");
  genParSource_=iConfig.getParameter<edm::InputTag>("genParSource");

}
BCToEFilterAlgo::~BCToEFilterAlgo ( )

Definition at line 16 of file BCToEFilterAlgo.cc.

                                  {
}

Member Function Documentation

bool BCToEFilterAlgo::filter ( const edm::Event iEvent)

Definition at line 22 of file BCToEFilterAlgo.cc.

References abs, reco::LeafCandidate::et(), reco::LeafCandidate::eta(), edm::Event::getByLabel(), reco::LeafCandidate::pdgId(), query::result, and reco::LeafCandidate::status().

                                                    {

  bool result=false;

  
  
  Handle<reco::GenParticleCollection> genParsHandle;
  iEvent.getByLabel(genParSource_,genParsHandle);
  reco::GenParticleCollection genPars=*genParsHandle;

  for (uint32_t ig=0;ig<genPars.size();ig++) {
    reco::GenParticle gp=genPars.at(ig);
    if (gp.status()==1 && abs(gp.pdgId())==11 && gp.et()>eTThreshold_ && fabs(gp.eta())<FILTER_ETA_MAX_) {
      if (hasBCAncestors(gp)) {
        result=true;
      }
    }
  }
  return result;

}
bool BCToEFilterAlgo::hasBCAncestors ( reco::GenParticle  gp)

Definition at line 49 of file BCToEFilterAlgo.cc.

References reco::CompositeRefCandidateT< D >::motherRef(), and reco::CompositeRefCandidateT< D >::numberOfMothers().

                                                       {

  //stopping condition: this particle is a b or c hadron
  if (isBCHadron(gp)) return true;
  //stopping condition: this particle has no mothers
  if (gp.numberOfMothers()==0) return false;
  //otherwise continue
  bool retval=false;
  for (uint32_t im=0;im<gp.numberOfMothers();im++) {
    retval=retval || hasBCAncestors(*gp.motherRef(im));
  }
  return retval;
  
}
bool BCToEFilterAlgo::isBCBaryon ( reco::GenParticle  gp) [private]

Definition at line 80 of file BCToEFilterAlgo.cc.

References abs, and reco::LeafCandidate::pdgId().

                                                   {
  
  uint32_t pdgid=abs(gp.pdgId());
  uint32_t thousands=pdgid%10000;
  if (thousands>=4000 && thousands <6000) {
    return true;
  } else {
    return false;
  }

}
bool BCToEFilterAlgo::isBCHadron ( reco::GenParticle  gp) [private]

Definition at line 64 of file BCToEFilterAlgo.cc.

                                                   {
  return isBCMeson(gp) || isBCBaryon(gp);
}
bool BCToEFilterAlgo::isBCMeson ( reco::GenParticle  gp) [private]

Definition at line 68 of file BCToEFilterAlgo.cc.

References abs, and reco::LeafCandidate::pdgId().

                                                  {
  
  uint32_t pdgid=abs(gp.pdgId());
  uint32_t hundreds=pdgid%1000;
  if (hundreds>=400 && hundreds<600) {
    return true;
  } else {
    return false;
  }

}

Member Data Documentation

Definition at line 49 of file BCToEFilterAlgo.h.

Definition at line 47 of file BCToEFilterAlgo.h.

Definition at line 50 of file BCToEFilterAlgo.h.