Merge pull request #318 from peiffer/53X_20130814

filters for initial state partons

GeneratorInterface/GenFilters/interface/MCPdgIndexFilter.h

@@ -0,0 +1,36 @@
+#ifndef MCPdgIndexFilter_h
+#define MCPdgIndexFilter_h
+/*
+ Description: filter events based on the particle PDG ID at a given
+ index in the HepMC::GenEvent record.
+
+ Original Author: Burt Betchart, 2013/08/09
+*/
+
+#include <memory>
+
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+#include "FWCore/Framework/interface/EDFilter.h"
+
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+
+
+class MCPdgIndexFilter : public edm::EDFilter {
+   public:
+      explicit MCPdgIndexFilter(const edm::ParameterSet&);
+      ~MCPdgIndexFilter() {};
+
+      virtual bool filter(edm::Event&, const edm::EventSetup&);
+   private:
+      bool pass(const edm::Event&);
+      const std::string label_;
+      const std::vector<int> pdgID;
+      const std::vector<unsigned> index;
+      const unsigned maxIndex;
+      const bool taggingMode;
+      const std::string tag;
+};
+#endif

GeneratorInterface/GenFilters/src/MCPdgIndexFilter.cc

@@ -0,0 +1,52 @@
+#include "GeneratorInterface/GenFilters/interface/MCPdgIndexFilter.h"
+#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
+
+MCPdgIndexFilter::MCPdgIndexFilter(const edm::ParameterSet& cfg) :
+  label_(cfg.getUntrackedParameter("moduleLabel",std::string("generator"))),
+  pdgID(cfg.getParameter<std::vector<int> >("PdgId")),
+  index(cfg.getParameter<std::vector<unsigned> >("Index")),
+  maxIndex(*std::max_element(index.begin(),index.end())),
+  taggingMode(cfg.getUntrackedParameter<bool>("TagMode",false)),
+  tag(cfg.getUntrackedParameter<std::string>("Tag",""))
+{
+  if (pdgID.size() != index.size())
+    edm::LogWarning("MCPdgIndexFilter")
+      << "Configuration Error :"
+      << "Sizes of array parameters 'PdgId' and 'Index' differ.";
+
+  if (taggingMode) {
+    produces<bool>(tag);
+    edm::LogInfo("TagMode") << "Filter result in '" << tag << "', filtering disabled.";
+  }
+}
+
+
+bool MCPdgIndexFilter::filter(edm::Event& evt, const edm::EventSetup&) {
+  bool result = pass(evt);
+  LogDebug("FilterResult") << (result?"Pass":"Fail");
+  if (!taggingMode) return result;
+  evt.put( std::auto_ptr<bool>(new bool(result)), tag);
+  return true;
+}
+
+
+bool MCPdgIndexFilter::pass(const edm::Event& evt) {
+  edm::Handle<edm::HepMCProduct> hepmc;
+  evt.getByLabel(label_, hepmc);
+
+  const HepMC::GenEvent * genEvent = hepmc->GetEvent();
+
+  HepMC::GenEvent::particle_const_iterator
+    p(genEvent->particles_begin()),
+    p_end(genEvent->particles_end());
+
+  for ( unsigned i=0; p!=p_end && i<=maxIndex; ++p, i++ ) {
+    LogDebug("Particle") << "index: " << i << "   pdgID: " << (*p)->pdg_id();
+    for (unsigned j = 0; j < pdgID.size(); j++) {
+      if (i==index[j] && pdgID[j] != (*p)->pdg_id())
+	return false;
+    }
+  }
+  return true;
+}

GeneratorInterface/GenFilters/src/SealModule.cc

@@ -14,6 +14,7 @@
 #include "GeneratorInterface/GenFilters/interface/MCDijetResonance.h"
 #include "GeneratorInterface/GenFilters/interface/MCProcessFilter.h"
 #include "GeneratorInterface/GenFilters/interface/MCProcessRangeFilter.h"
+#include "GeneratorInterface/GenFilters/interface/MCPdgIndexFilter.h"
 #include "GeneratorInterface/GenFilters/interface/MCSingleParticleFilter.h"
 #include "GeneratorInterface/GenFilters/interface/MCSmartSingleParticleFilter.h"
 #include "GeneratorInterface/GenFilters/interface/MCZll.h"
@@ -59,6 +60,7 @@
   DEFINE_FWK_MODULE(MCDijetResonance);
   DEFINE_FWK_MODULE(MCProcessFilter);
   DEFINE_FWK_MODULE(MCProcessRangeFilter);
+  DEFINE_FWK_MODULE(MCPdgIndexFilter);
   DEFINE_FWK_MODULE(MCSingleParticleFilter);
   DEFINE_FWK_MODULE(MCSmartSingleParticleFilter);
   DEFINE_FWK_MODULE(MCZll);