Changes to PartonShowerBsHepMCFilter and LHEGenericFilter

GeneratorInterface/Core/BuildFile.xml

@@ -3,6 +3,7 @@
 <use   name="FWCore/Utilities"/>
 <use   name="SimDataFormats/GeneratorProducts"/>
 <use   name="GeneratorInterface/LHEInterface"/>
+<use   name="heppdt"/>
 <use   name="boost"/>
 <use   name="clhep"/>
 <use   name="lhapdf"/>

GeneratorInterface/Core/interface/PartonShowerBsHepMCFilter.h

@@ -24,10 +24,7 @@ class PartonShowerBsHepMCFilter : public BaseHepMCFilter{
 
 private:
 
-  int particle_id;
-  int exclude_status_id;
-  int status_id;
 };
 
 
-#endif
+#endif

GeneratorInterface/Core/src/HepMCFilterDriver.cc

@@ -87,4 +87,4 @@ void HepMCFilterDriver::resetStatistics() {
   sumtotal_w_        = 0;
   sumtotal_w2_       = 0;
 
-}
+}

GeneratorInterface/Core/src/PartonShowerBsHepMCFilter.cc

@@ -1,22 +1,15 @@
 #include "GeneratorInterface/Core/interface/PartonShowerBsHepMCFilter.h"
 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
 #include <iostream>
+#include "HepPDT/ParticleID.hh"
 
 
 using namespace edm;
 using namespace std;
 
 
 //constructor
-PartonShowerBsHepMCFilter::PartonShowerBsHepMCFilter(const edm::ParameterSet& iConfig) :
-
-  // particle id of the gen particles that you want to filter
-  particle_id(iConfig.getParameter<int>("Particle_id")),
-  // status id of the particles that you want to exclude from the filter
-  exclude_status_id(iConfig.getUntrackedParameter<int>("Exclude_status_id",-1)),
-  // status id of the particles that you want to filetr on
-  status_id(iConfig.getUntrackedParameter<int>("Status_id",-1))
-
+PartonShowerBsHepMCFilter::PartonShowerBsHepMCFilter(const edm::ParameterSet& iConfig) 
 {
 
 }
@@ -35,26 +28,21 @@ PartonShowerBsHepMCFilter::~PartonShowerBsHepMCFilter()
 // ------------ method called to produce the data  ------------
 bool PartonShowerBsHepMCFilter::filter(const HepMC::GenEvent* evt)
 {
-
-  if( exclude_status_id > 0. && status_id > 0.){
-    std::cout << "ERROR: Skipping event: Configuration has both exclude and status id set to a value > 0. They can not be used simultaneously." << std::endl;
-    return false; // skip event
-  }  
 
-  for ( HepMC::GenEvent::particle_const_iterator p = evt->particles_begin();
-	p != evt->particles_end(); ++p ) {
+  // loop over gen particles
+  for ( HepMC::GenEvent::particle_const_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p ){
 
-    if( abs((*p)->pdg_id()) == particle_id ){	
-      if( exclude_status_id > 0. && (*p)->status() != exclude_status_id ) 
-        return true; // keep event
-      else if( status_id > 0. && (*p)->status() == status_id )
-        return true; // keep event
-      else 
-        return true; // keep event
+    // check only status 2 particles
+    if( (*p)->status()==2 ){
+      // if one of the status 2 particles is a B-hadron, accept the event
+      HepPDT::ParticleID pid((*p)->pdg_id());
+      if( pid.hasBottom() ){
+        return true; // accept event
+      }
     }
-
+    
   }
 
   return false; // skip event
 
-}
+}

GeneratorInterface/GenFilters/interface/LHEGenericFilter.h

@@ -49,8 +49,10 @@ class LHEGenericFilter : public edm::EDFilter {
 
   edm::EDGetTokenT<LHEEventProduct> src_;
   int numRequired_;                // number of particles required to pass filter
-  bool acceptMore_;                // if true (default), accept numRequired or more.
-                                   // if false, accept events with exactly equal to numRequired.
+  std::string acceptLogic_;        // LT  meaning <
+                                   // GT          >
+                                   // EQ          =
+                                   // NE          !=
   std::vector<int> particleID_;    // vector of particle IDs to look for
   int totalEvents_;                // counters
   int passedEvents_;

GeneratorInterface/GenFilters/python/LHEGenericFilter_cfi.py

@@ -0,0 +1,8 @@
+import FWCore.ParameterSet.Config as cms
+
+lheGenericFilter = cms.EDFilter("LHEGenericFilter",
+    src = cms.InputTag("source"),
+    NumRequired = cms.int32(2),
+    ParticleID = cms.vint32(5),
+    AcceptLogic = cms.string("LT") # LT meaning < NumRequired, GT >, EQ =, NE !=
+)                                

GeneratorInterface/GenFilters/src/LHEGenericFilter.cc

@@ -1,20 +1,22 @@
 #include "GeneratorInterface/GenFilters/interface/LHEGenericFilter.h"
 
+using namespace edm;
+using namespace std;
+
 LHEGenericFilter::LHEGenericFilter(const edm::ParameterSet& iConfig) :
-  numRequired_(iConfig.getParameter<int>("NumRequired")),
-  acceptMore_(iConfig.getParameter<bool>("AcceptMore")),
-  particleID_(iConfig.getParameter< std::vector<int> >("ParticleID")),
-  totalEvents_(0), passedEvents_(0)
+numRequired_(iConfig.getParameter<int>("NumRequired")),
+acceptLogic_(iConfig.getParameter<std::string>("AcceptLogic")),
+particleID_(iConfig.getParameter< std::vector<int> >("ParticleID")),
+totalEvents_(0), passedEvents_(0)
 {
   //here do whatever other initialization is needed
   src_ = consumes<LHEEventProduct>(iConfig.getParameter<edm::InputTag>("src"));
 }
 
 LHEGenericFilter::~LHEGenericFilter()
 {
-
-   // do anything here that needs to be done at destruction time
-   // (e.g. close files, deallocate resources etc.)
+  // do anything here that needs to be done at destruction time
+  // (e.g. close files, deallocate resources etc.)
 
 }
 
@@ -29,20 +31,25 @@ bool LHEGenericFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup)
   int nFound = 0;
 
   for (int i = 0; i < EvtHandle->hepeup().NUP; ++i) {
-     if (EvtHandle->hepeup().ISTUP[i] != 1) {
-       continue;
-     }    
+    if (EvtHandle->hepeup().ISTUP[i] != 1) { // keep only outgoing particles
+      continue;
+    }
     for (unsigned int j = 0; j < particleID_.size(); ++j) {
       if (particleID_[j] == 0 || abs(particleID_[j]) == abs(EvtHandle->hepeup().IDUP[i]) ) {
-	nFound++;
-	break; // only match a given particle once!
+        nFound++;
+        break; // only match a given particle once!
       }
     } // loop over targets
 
-    if (acceptMore_ && nFound == numRequired_) break; // stop looking if we don't mind having more
   } // loop over particles
 
-  if (nFound == numRequired_) {
+  // event accept/reject logic
+  if (
+         (acceptLogic_.compare("LT")==0 && nFound < numRequired_)
+      || (acceptLogic_.compare("GT")==0 && nFound > numRequired_)
+      || (acceptLogic_.compare("EQ")==0 && nFound == numRequired_)
+      || (acceptLogic_.compare("NE")==0 && nFound != numRequired_)
+      ) {
     passedEvents_++;
     return true;
   } else {
@@ -54,7 +61,7 @@ bool LHEGenericFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup)
 // ------------ method called once each job just after ending the event loop  ------------
 void LHEGenericFilter::endJob() {
   edm::LogInfo("LHEGenericFilter") << "=== Results of LHEGenericFilter: passed "
-                                        << passedEvents_ << "/" << totalEvents_ << " events" << std::endl;
+  << passedEvents_ << "/" << totalEvents_ << " events" << std::endl;
 }
 
 //define this as a plug-in