Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-02-05 23:51:09

0001 #include "CommonTools/CandAlgos/interface/GenJetParticleSelector.h"
0002 #include "DataFormats/Candidate/interface/Candidate.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0005 #include "SimGeneral/HepPDTRecord/interface/PdtEntry.h"
0006 #include "FWCore/Utilities/interface/EDMException.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include <algorithm>
0009 using namespace std;
0010 using namespace edm;
0011 
0012 GenJetParticleSelector::GenJetParticleSelector(const ParameterSet& cfg, edm::ConsumesCollector& iC)
0013     : stableOnly_(cfg.getParameter<bool>("stableOnly")),
0014       partons_(false),
0015       bInclude_(false),
0016       tableToken_(iC.esConsumes()) {
0017   const string excludeString("excludeList");
0018   const string includeString("includeList");
0019   vpdt includeList, excludeList;
0020   vector<string> vPdtParams = cfg.getParameterNamesForType<vpdt>();
0021   bool found = std::find(vPdtParams.begin(), vPdtParams.end(), includeString) != vPdtParams.end();
0022   if (found)
0023     includeList = cfg.getParameter<vpdt>(includeString);
0024   found = find(vPdtParams.begin(), vPdtParams.end(), excludeString) != vPdtParams.end();
0025   if (found)
0026     excludeList = cfg.getParameter<vpdt>(excludeString);
0027   const string partonsString("partons");
0028   vector<string> vBoolParams = cfg.getParameterNamesForType<bool>();
0029   found = find(vBoolParams.begin(), vBoolParams.end(), partonsString) != vBoolParams.end();
0030   if (found)
0031     partons_ = cfg.getParameter<bool>(partonsString);
0032   bool bExclude = false;
0033   if (!includeList.empty())
0034     bInclude_ = true;
0035   if (!excludeList.empty())
0036     bExclude = true;
0037 
0038   if (bInclude_ && bExclude) {
0039     throw cms::Exception("ConfigError", "not allowed to use both includeList and excludeList at the same time\n");
0040   } else if (bInclude_) {
0041     pdtList_ = includeList;
0042   } else {
0043     pdtList_ = excludeList;
0044   }
0045   if (stableOnly_ && partons_) {
0046     throw cms::Exception("ConfigError", "not allowed to have both stableOnly and partons true at the same time\n");
0047   }
0048 }
0049 
0050 void GenJetParticleSelector::fillPSetDescription(edm::ParameterSetDescription& desc) {
0051   desc.add<bool>("stableOnly", true);
0052 }
0053 
0054 bool GenJetParticleSelector::operator()(const reco::Candidate& p) {
0055   int status = p.status();
0056   int id = abs(p.pdgId());
0057   if ((!stableOnly_ || status == 1) && !partons_ && ((pIds_.find(id) == pIds_.end()) ^ bInclude_))
0058     return true;
0059   else if (partons_ && (p.numberOfDaughters() > 0 && (p.daughter(0)->pdgId() == 91 || p.daughter(0)->pdgId() == 92)) &&
0060            (((pIds_.find(id) == pIds_.end()) ^ bInclude_)))
0061     return true;
0062   else
0063     return false;
0064 }
0065 
0066 void GenJetParticleSelector::init(const edm::EventSetup& es) {
0067   auto const& pdt = es.getData(tableToken_);
0068 
0069   for (vpdt::iterator i = pdtList_.begin(); i != pdtList_.end(); ++i)
0070     i->setup(pdt);
0071 }