GenFilters/plugins/LHEGenericMassFilter.cc

0001 // system include files
0002 #include <memory>
0003 #include <iostream>
0004
0005 // user include files
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/global/EDFilter.h"
0008
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/PluginManager/interface/ModuleDef.h"
0012
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0015
0016 using namespace edm;
0017 using namespace std;
0018
0019 //
0020 // class declaration
0021 //
0022
0023 class LHEGenericMassFilter : public edm::global::EDFilter<> {
0024 public:
0025   explicit LHEGenericMassFilter(const edm::ParameterSet&);
0026   ~LHEGenericMassFilter() override = default;
0027   static void fillDescriptions(edm::ConfigurationDescriptions&);
0028
0029 private:
0030   bool filter(edm::StreamID, edm::Event&, edm::EventSetup const&) const override;
0031
0032   // ----------member data ---------------------------
0033
0034   const edm::EDGetTokenT<LHEEventProduct> src_;
0035   const int numRequired_;              // number of particles required to pass filter
0036   const std::vector<int> particleID_;  // vector of particle IDs to look for
0037   const double minMass_;
0038   const double maxMass_;
0039   const bool requiredOutgoingStatus_;  // Whether particles required to pass filter must have outgoing status
0040 };
0041
0042 LHEGenericMassFilter::LHEGenericMassFilter(const edm::ParameterSet& iConfig)
0043     : src_(consumes<LHEEventProduct>(iConfig.getParameter<edm::InputTag>("src"))),
0044       numRequired_(iConfig.getParameter<int>("NumRequired")),
0045       particleID_(iConfig.getParameter<std::vector<int>>("ParticleID")),
0046       minMass_(iConfig.getParameter<double>("MinMass")),
0047       maxMass_(iConfig.getParameter<double>("MaxMass")),
0048       requiredOutgoingStatus_(iConfig.getParameter<bool>("RequiredOutgoingStatus")) {}
0049
0050 // ------------ method called to skim the data  ------------
0051 bool LHEGenericMassFilter::filter(edm::StreamID iID, edm::Event& iEvent, edm::EventSetup const& iSetup) const {
0052   edm::Handle<LHEEventProduct> EvtHandle;
0053   iEvent.getByToken(src_, EvtHandle);
0054
0055   int nFound = 0;
0056
0057   double Px = 0.;
0058   double Py = 0.;
0059   double Pz = 0.;
0060   double E = 0.;
0061
0062   for (int i = 0; i < EvtHandle->hepeup().NUP; ++i) {
0063     // if requiredOutgoingStatus_ keep only outgoing particles, otherwise keep them all
0064     if (requiredOutgoingStatus_ && EvtHandle->hepeup().ISTUP[i] != 1) {
0065       continue;
0066     }
0067     for (unsigned int j = 0; j < particleID_.size(); ++j) {
0068       if (abs(particleID_[j]) == abs(EvtHandle->hepeup().IDUP[i])) {
0069         nFound++;
0070         Px = Px + EvtHandle->hepeup().PUP[i][0];
0071         Py = Py + EvtHandle->hepeup().PUP[i][1];
0072         Pz = Pz + EvtHandle->hepeup().PUP[i][2];
0073         E = E + EvtHandle->hepeup().PUP[i][3];
0074
0075         break;  // only match a given particle once!
0076       }
0077     }  // loop over targets
0078
0079   }  // loop over particles
0080
0081   // event accept/reject logic
0082   if (nFound == numRequired_) {
0083     double sqrdMass = E * E - (Px * Px + Py * Py + Pz * Pz);
0084     if (sqrdMass > minMass_ * minMass_ && sqrdMass < maxMass_ * maxMass_) {
0085       return true;
0086     }
0087   }
0088   return false;
0089 }
0090
0091 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0092 void LHEGenericMassFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0093   edm::ParameterSetDescription desc;
0094   desc.add<edm::InputTag>("src", edm::InputTag("externalLHEProducer"));
0095   desc.add<int>("NumRequired", 1);
0096   desc.add<vector<int>>("ParticleID", std::vector<int>{1});
0097   desc.add<double>("MinMass", 0.0);
0098   desc.add<double>("MaxMass", 1.0);
0099   desc.add<bool>("RequiredOutgoingStatus", true);
0100   descriptions.addDefault(desc);
0101 }
0102
0103 //define this as a plug-in
0104 DEFINE_FWK_MODULE(LHEGenericMassFilter);