GenFilters/plugins/BCToEFilterAlgo.cc

0001 #include "GeneratorInterface/GenFilters/plugins/BCToEFilterAlgo.h"
0002 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0003 #include "FWCore/Framework/interface/ConsumesCollector.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/Utilities/interface/InputTag.h"
0007
0008 #include <cmath>
0009 #include <cstdlib>
0010
0011 BCToEFilterAlgo::BCToEFilterAlgo(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& iC)
0012     :  //set constants
0013       maxAbsEta_((float)iConfig.getParameter<double>("maxAbsEta")),
0014       eTThreshold_((float)iConfig.getParameter<double>("eTThreshold")),
0015       genParSource_(iC.consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("genParSource"))) {}
0016
0017 BCToEFilterAlgo::~BCToEFilterAlgo() {}
0018
0019 //look for status==1 e coming from b or c hadron
0020 //there is an eT threshold on the electron (configurable)
0021 bool BCToEFilterAlgo::filter(const edm::Event& iEvent) const {
0022   bool result = false;
0023
0024   edm::Handle<reco::GenParticleCollection> genParsHandle;
0025   iEvent.getByToken(genParSource_, genParsHandle);
0026   reco::GenParticleCollection genPars = *genParsHandle;
0027
0028   for (uint32_t ig = 0; ig < genPars.size(); ig++) {
0029     reco::GenParticle gp = genPars.at(ig);
0030     if (gp.status() == 1 && std::abs(gp.pdgId()) == 11 && gp.et() > eTThreshold_ && std::fabs(gp.eta()) < maxAbsEta_) {
0031       if (hasBCAncestors(gp)) {
0032         result = true;
0033       }
0034     }
0035   }
0036   return result;
0037 }
0038
0039 //does this particle have an ancestor (mother, mother of mother, etc.) that is a b or c hadron?
0040 //attention: the GenParticle argument must have the motherRef correctly filled for this
0041 //to work.  That is, you had better have got it out of the genParticles collection
0042 bool BCToEFilterAlgo::hasBCAncestors(const reco::GenParticle& gp) const {
0043   //stopping condition: this particle is a b or c hadron
0044   if (isBCHadron(gp))
0045     return true;
0046   //stopping condition: this particle has no mothers
0047   if (gp.numberOfMothers() == 0)
0048     return false;
0049   //otherwise continue
0050   bool retval = false;
0051   for (uint32_t im = 0; im < gp.numberOfMothers(); im++) {
0052     retval = retval || hasBCAncestors(*gp.motherRef(im));
0053   }
0054   return retval;
0055 }
0056
0057 bool BCToEFilterAlgo::isBCHadron(const reco::GenParticle& gp) const { return isBCMeson(gp) || isBCBaryon(gp); }
0058
0059 bool BCToEFilterAlgo::isBCMeson(const reco::GenParticle& gp) const {
0060   uint32_t pdgid = std::abs(gp.pdgId());
0061   uint32_t hundreds = pdgid % 1000;
0062   if (hundreds >= 400 && hundreds < 600) {
0063     return true;
0064   } else {
0065     return false;
0066   }
0067 }
0068
0069 bool BCToEFilterAlgo::isBCBaryon(const reco::GenParticle& gp) const {
0070   uint32_t pdgid = std::abs(gp.pdgId());
0071   uint32_t thousands = pdgid % 10000;
0072   if (thousands >= 4000 && thousands < 6000) {
0073     return true;
0074   } else {
0075     return false;
0076   }
0077 }