Skimming/src/BTagSkimMC.cc

0001 #include "RecoBTag/Skimming/interface/BTagSkimMC.h"
0002 #include "DataFormats/Common/interface/Handle.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Utilities/interface/InputTag.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006
0007 using namespace edm;
0008 using namespace std;
0009
0010 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0011 using namespace reco;
0012
0013 BTagSkimMC::BTagSkimMC(const ParameterSet& p, const BTagSkimMCCount::Counters* count) : nEvents_(0), nAccepted_(0) {
0014   verbose = p.getUntrackedParameter<bool>("verbose", false);
0015   pthatMin = p.getParameter<double>("pthat_min");
0016   pthatMax = p.getParameter<double>("pthat_max");
0017   process_ = p.getParameter<string>("mcProcess");
0018   if (verbose)
0019     cout << " Requested:  " << process_ << endl;
0020 }
0021
0022 bool BTagSkimMC::filter(Event& evt, const EventSetup& es) {
0023   nEvents_++;
0024
0025   Handle<int> genProcessID;
0026   evt.getByLabel("genEventProcID", genProcessID);
0027   double processID = *genProcessID;
0028
0029   Handle<double> genEventScale;
0030   evt.getByLabel("genEventScale", genEventScale);
0031   double pthat = *genEventScale;
0032
0033   if (verbose)
0034     cout << "processID: " << processID << " - pthat: " << pthat;
0035
0036   if ((processID != 4) && (process_ == "QCD")) {  // the Pythia events (for ALPGEN see below)
0037
0038     Handle<double> genFilterEff;
0039     evt.getByLabel("genEventRunInfo", "FilterEfficiency", genFilterEff);
0040     double filter_eff = *genFilterEff;
0041     if (verbose)
0042       cout << " Is QCD ";
0043     // qcd (including min bias HS)
0044     if ((filter_eff == 1. || filter_eff == 0.964) && (processID == 11 || processID == 12 || processID == 13 ||
0045                                                       processID == 28 || processID == 68 || processID == 53)) {
0046       if (pthat > pthatMin && pthat < pthatMax) {
0047         if (verbose)
0048           cout << " ACCEPTED " << endl;
0049         nAccepted_++;
0050         return true;
0051       }
0052     }
0053
0054   }                           // ALPGEN
0055   else if (processID == 4) {  // this is the number for external ALPGEN events
0056
0057     Handle<GenParticleCollection> genParticles;
0058     evt.getByLabel("genParticles", genParticles);
0059
0060     for (size_t i = 0; i < genParticles->size(); ++i) {
0061       const Candidate& p = (*genParticles)[i];
0062       int id = p.pdgId();
0063       int st = p.status();
0064
0065       // tt+jets
0066       if (st == 3 && (id == 6 || id == -6)) {
0067         if (verbose)
0068           cout << "We have a ttbar event" << endl;
0069         nAccepted_++;
0070         return true;
0071       }
0072     }
0073   }
0074   if (verbose)
0075     cout << " REJECTED " << endl;
0076
0077   return false;
0078 }
0079
0080 void BTagSkimMC::endStream() {
0081   globalCache()->nEvents_ += nEvents_;
0082   globalCache()->nAccepted_ += nAccepted_;
0083 }
0084
0085 void BTagSkimMC::globalEndJob(const BTagSkimMCCount::Counters* count) {
0086   edm::LogVerbatim("BTagSkimMC") << "=============================================================================\n"
0087                                  << " Events read: " << count->nEvents_
0088                                  << "\n Events accepted by BTagSkimMC: " << count->nAccepted_
0089                                  << "\n Efficiency: " << (double)(count->nAccepted_) / (double)(count->nEvents_)
0090                                  << "\n==========================================================================="
0091                                  << std::endl;
0092 }
0093
0094 #include "FWCore/Framework/interface/MakerMacros.h"
0095
0096 DEFINE_FWK_MODULE(BTagSkimMC);