HepMCCandAlgos/plugins/ModelFilter.cc

0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/global/EDFilter.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005
0006 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0007
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011 #include <iostream>
0012 #include <memory>
0013 #include <sstream>
0014 #include <cstdlib>
0015
0016 /**
0017     The ModelFilter class will select events in a "soup" MC
0018     (like the SUSY signal MC) from the comments of LHEEventProduct
0019     that match "modelTag". The user can require the value of that
0020     parameter to lie between a min and max value.
0021  */
0022
0023 namespace edm {
0024
0025   class ModelFilter : public edm::global::EDFilter<> {
0026   public:
0027     explicit ModelFilter(const edm::ParameterSet&);
0028
0029     static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0030     typedef std::vector<std::string>::const_iterator comments_const_iterator;
0031
0032   private:
0033     bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0034     static std::vector<std::string> split(std::string const& fstring, std::string const& splitter);
0035
0036     edm::EDGetTokenT<LHEEventProduct> tokenSource_;
0037     std::string modelTag_;
0038     std::vector<double> parameterMins_;
0039     std::vector<double> parameterMaxs_;
0040   };
0041
0042 }  // namespace edm
0043
0044 using namespace std;
0045 using namespace edm;
0046
0047 ModelFilter::ModelFilter(const edm::ParameterSet& iConfig) {
0048   tokenSource_ = consumes<LHEEventProduct>(iConfig.getParameter<InputTag>("source"));
0049   modelTag_ = iConfig.getParameter<string>("modelTag");
0050   parameterMins_ = iConfig.getParameter<vector<double> >("parameterMins");
0051   parameterMaxs_ = iConfig.getParameter<vector<double> >("parameterMaxs");
0052 }
0053
0054 bool ModelFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0055   Handle<LHEEventProduct> product;
0056   iEvent.getByToken(tokenSource_, product);
0057   comments_const_iterator comment;
0058
0059   string tempString;
0060   vector<string> parameters;
0061
0062   for (comment = product->comments_begin(); comment != product->comments_end(); comment++) {
0063     if (comment->find(modelTag_) != string::npos) {
0064       tempString = comment->substr(comment->find(modelTag_), comment->size());
0065       tempString = tempString.substr(0, tempString.find(' '));
0066       parameters = split(tempString, "_");
0067
0068       if (parameters.size() - 1 != parameterMins_.size()) {
0069         edm::LogError("ModelFilter") << "number of modeParameters does not match number of parameters in file";
0070         return false;
0071       } else if (parameterMins_.size() != parameterMaxs_.size()) {
0072         edm::LogError("ModelFilter") << "Error: umber of parameter mins != number parameter maxes";
0073       } else {
0074         for (unsigned i = 0; i < parameterMins_.size(); i++) {
0075           if (parameterMins_[i] > atof(parameters[i + 1].c_str()) ||
0076               parameterMaxs_[i] < atof(parameters[i + 1].c_str())) {
0077             return false;
0078           }
0079         }
0080         return true;
0081       }
0082     }
0083   }
0084   edm::LogInfo("ModelFilter") << "FAILED: " << *comment;
0085   return false;
0086 }
0087 void ModelFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0088   edm::ParameterSetDescription desc;
0089   desc.add<InputTag>("source");
0090   desc.add<string>("modelTag");
0091   desc.add<vector<double> >("parameterMins");
0092   desc.add<vector<double> >("parameterMaxs");
0093
0094   descriptions.addDefault(desc);
0095 }
0096 vector<string> ModelFilter::split(string const& fstring, string const& splitter) {
0097   vector<string> returnVector;
0098   size_t cursor;
0099   string beforeSplitter;
0100   string afterSplitter = fstring;
0101   if (fstring.find(splitter) == string::npos) {
0102     edm::LogInfo("ModelFilter") << "No " << splitter << " found";
0103     returnVector.push_back(fstring);
0104     return returnVector;
0105   } else {
0106     while (afterSplitter.find(splitter) != string::npos) {
0107       cursor = afterSplitter.find(splitter);
0108
0109       beforeSplitter = afterSplitter.substr(0, cursor);
0110       afterSplitter = afterSplitter.substr(cursor + 1, afterSplitter.size());
0111
0112       returnVector.push_back(beforeSplitter);
0113
0114       if (afterSplitter.find(splitter) == string::npos)
0115         returnVector.push_back(afterSplitter);
0116     }
0117     return returnVector;
0118   }
0119 }
0120 DEFINE_FWK_MODULE(ModelFilter);