Higgs/plugins/HLTHiggsPlotter.cc

0001
0002 /** \file HLTHiggsPlotter.cc
0003 */
0004
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "DataFormats/Candidate/interface/CandMatchMap.h"
0007 #include "DataFormats/MuonReco/interface/Muon.h"
0008 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0009
0010 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0011
0012 #include "HLTHiggsPlotter.h"
0013 #include "HLTHiggsSubAnalysis.h"
0014 #include "EVTColContainer.h"
0015
0016 #include "TPRegexp.h"
0017
0018 #include <set>
0019 #include <cctype>
0020
0021 HLTHiggsPlotter::HLTHiggsPlotter(const edm::ParameterSet &pset,
0022                                  const std::string &hltPath,
0023                                  const std::vector<unsigned int> &objectsType,
0024                                  const unsigned int &NptPlots,
0025                                  const std::vector<double> &NminOneCuts)
0026     : _hltPath(hltPath),
0027       _hltProcessName(pset.getParameter<std::string>("hltProcessName")),
0028       _objectsType(std::set<unsigned int>(objectsType.begin(), objectsType.end())),
0029       _nObjects(objectsType.size()),
0030       _parametersEta(pset.getParameter<std::vector<double> >("parametersEta")),
0031       _parametersPhi(pset.getParameter<std::vector<double> >("parametersPhi")),
0032       _parametersTurnOn(pset.getParameter<std::vector<double> >("parametersTurnOn")),
0033       _NptPlots(NptPlots),
0034       _NminOneCuts(NminOneCuts) {
0035   for (std::set<unsigned int>::iterator it = _objectsType.begin(); it != _objectsType.end(); ++it) {
0036     // Some parameters extracted from the .py
0037     std::string objStr = EVTColContainer::getTypeString(*it);
0038     _cutMinPt[*it] = pset.getParameter<double>(std::string(objStr + "_cutMinPt").c_str());
0039     _cutMaxEta[*it] = pset.getParameter<double>(std::string(objStr + "_cutMaxEta").c_str());
0040   }
0041 }
0042
0043 HLTHiggsPlotter::~HLTHiggsPlotter() {}
0044
0045 void HLTHiggsPlotter::beginJob() {}
0046
0047 void HLTHiggsPlotter::beginRun(const edm::Run &iRun, const edm::EventSetup &iSetup) {}
0048
0049 void HLTHiggsPlotter::bookHistograms(DQMStore::IBooker &ibooker, const bool &useNminOneCuts) {
0050   for (std::set<unsigned int>::iterator it = _objectsType.begin(); it != _objectsType.end(); ++it) {
0051     std::vector<std::string> sources(2);
0052     sources[0] = "gen";
0053     sources[1] = "rec";
0054     TString maxPt;
0055
0056     const std::string objTypeStr = EVTColContainer::getTypeString(*it);
0057
0058     for (size_t i = 0; i < sources.size(); i++) {
0059       std::string source = sources[i];
0060
0061       if (useNminOneCuts && *it == EVTColContainer::PFJET) {
0062         if (source == "gen")
0063           continue;
0064         else {
0065           // N-1 jet plots (dEtaqq, mqq, dPhibb, CSV1, maxCSV_jets, maxCSV_E, PFMET, pt1, pt2, pt3, pt4)
0066           if (_NminOneCuts[0])
0067             bookHist(source, objTypeStr, "dEtaqq", ibooker);
0068           if (_NminOneCuts[1])
0069             bookHist(source, objTypeStr, "mqq", ibooker);
0070           if (_NminOneCuts[2])
0071             bookHist(source, objTypeStr, "dPhibb", ibooker);
0072           if (_NminOneCuts[3]) {
0073             if (_NminOneCuts[6])
0074               bookHist(source, objTypeStr, "maxCSV", ibooker);
0075             else
0076               bookHist(source, objTypeStr, "CSV1", ibooker);
0077           }
0078           if (_NminOneCuts[4])
0079             bookHist(source, objTypeStr, "CSV2", ibooker);
0080           if (_NminOneCuts[5])
0081             bookHist(source, objTypeStr, "CSV3", ibooker);
0082         }
0083       }
0084
0085       bookHist(source, objTypeStr, "Eta", ibooker);
0086       bookHist(source, objTypeStr, "Phi", ibooker);
0087       for (unsigned int i = 0; i < _NptPlots; i++) {
0088         maxPt = "MaxPt";
0089         maxPt += i + 1;
0090         bookHist(source, objTypeStr, maxPt.Data(), ibooker);
0091       }
0092     }
0093   }
0094 }
0095
0096 void HLTHiggsPlotter::analyze(const bool &isPassTrigger,
0097                               const std::string &source,
0098                               const std::vector<MatchStruct> &matches,
0099                               const unsigned int &minCandidates) {
0100   if (!isPassTrigger) {
0101     return;
0102   }
0103   std::map<unsigned int, int> countobjects;
0104   // Initializing the count of the used object
0105   for (std::set<unsigned int>::iterator co = _objectsType.begin(); co != _objectsType.end(); ++co) {
0106     countobjects[*co] = 0;
0107   }
0108
0109   int counttotal = 0;
0110   const int totalobjectssize2 = _NptPlots * countobjects.size();
0111   // Fill the histos if pass the trigger (just the two with higher pt)
0112   for (size_t j = 0; j < matches.size(); ++j) {
0113     // Is this object owned by this trigger? If not we are not interested...
0114     if (_objectsType.find(matches[j].objType) == _objectsType.end()) {
0115       continue;
0116     }
0117
0118     const unsigned int objType = matches[j].objType;
0119     const std::string objTypeStr = EVTColContainer::getTypeString(matches[j].objType);
0120
0121     float pt = matches[j].pt;
0122     float eta = matches[j].eta;
0123     float phi = matches[j].phi;
0124
0125     TString maxPt;
0126     if ((unsigned)countobjects[objType] < _NptPlots) {
0127       maxPt = "MaxPt";
0128       maxPt += (countobjects[objType] + 1);
0129       this->fillHist(isPassTrigger, source, objTypeStr, maxPt.Data(), pt);
0130       // Filled the high pt ...
0131       ++(countobjects[objType]);
0132       ++counttotal;
0133     } else {
0134       if ((unsigned)countobjects[objType] < minCandidates) {  // To get correct results for HZZ
0135         ++(countobjects[objType]);
0136         ++counttotal;
0137       } else
0138         continue;  //   Otherwise too many entries in Eta and Phi distributions
0139     }
0140
0141     this->fillHist(isPassTrigger, source, objTypeStr, "Eta", eta);
0142     this->fillHist(isPassTrigger, source, objTypeStr, "Phi", phi);
0143
0144     if (counttotal == totalobjectssize2) {
0145       break;
0146     }
0147   }
0148 }
0149
0150 void HLTHiggsPlotter::analyze(const bool &isPassTrigger,
0151                               const std::string &source,
0152                               const std::vector<MatchStruct> &matches,
0153                               std::map<std::string, bool> &nMinOne,
0154                               const float &dEtaqq,
0155                               const float &mqq,
0156                               const float &dPhibb,
0157                               const float &CSV1,
0158                               const float &CSV2,
0159                               const float &CSV3,
0160                               const bool &passAllCuts) {
0161   if (!isPassTrigger) {
0162     return;
0163   }
0164   std::map<unsigned int, int> countobjects;
0165   // Initializing the count of the used object
0166   for (std::set<unsigned int>::iterator co = _objectsType.begin(); co != _objectsType.end(); ++co) {
0167     if (!(*co == EVTColContainer::PFJET && source == "gen"))  // genJets are not there
0168       countobjects[*co] = 0;
0169   }
0170
0171   int counttotal = 0;
0172   const int totalobjectssize2 = _NptPlots * countobjects.size();
0173   // Fill the histos if pass the trigger (just the two with higher pt)
0174   for (size_t j = 0; j < matches.size(); ++j) {
0175     // Is this object owned by this trigger? If not we are not interested...
0176     if (_objectsType.find(matches[j].objType) == _objectsType.end()) {
0177       continue;
0178     }
0179
0180     const unsigned int objType = matches[j].objType;
0181     const std::string objTypeStr = EVTColContainer::getTypeString(matches[j].objType);
0182
0183     float pt = matches[j].pt;
0184     float eta = matches[j].eta;
0185     float phi = matches[j].phi;
0186
0187     // PFMET N-1 cut
0188     if (objType == EVTColContainer::PFMET && _NminOneCuts[8] && !nMinOne["PFMET"])
0189       continue;
0190
0191     TString maxPt;
0192     if ((unsigned)(countobjects)[objType] < _NptPlots) {
0193       maxPt = "MaxPt";
0194       maxPt += (countobjects[objType] + 1);
0195       if (objType != EVTColContainer::PFJET || nMinOne[maxPt.Data()]) {
0196         this->fillHist(isPassTrigger, source, objTypeStr, maxPt.Data(), pt);
0197       }
0198       ++(countobjects[objType]);
0199       ++counttotal;
0200     } else
0201       continue;  // if not needed (minCandidates == _NptPlots if _useNminOneCuts
0202     if (objType != EVTColContainer::PFJET || passAllCuts) {
0203       this->fillHist(isPassTrigger, source, objTypeStr, "Eta", eta);
0204       this->fillHist(isPassTrigger, source, objTypeStr, "Phi", phi);
0205     }
0206
0207     if (counttotal == totalobjectssize2) {
0208       break;
0209     }
0210   }
0211   if (source == "rec" && _objectsType.find(EVTColContainer::PFJET) != _objectsType.end()) {
0212     if (_NminOneCuts[0] && nMinOne["dEtaqq"]) {
0213       this->fillHist(isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), "dEtaqq", dEtaqq);
0214     }
0215     if (_NminOneCuts[1] && nMinOne["mqq"]) {
0216       this->fillHist(isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), "mqq", mqq);
0217     }
0218     if (_NminOneCuts[2] && nMinOne["dPhibb"]) {
0219       this->fillHist(isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), "dPhibb", dPhibb);
0220     }
0221     if (_NminOneCuts[3]) {
0222       std::string nameCSVplot = "CSV1";
0223       if (_NminOneCuts[6])
0224         nameCSVplot = "maxCSV";
0225       if (nMinOne[nameCSVplot])
0226         this->fillHist(
0227             isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), nameCSVplot, CSV1);
0228     }
0229     if (_NminOneCuts[4] && nMinOne["CSV2"]) {
0230       this->fillHist(isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), "CSV2", CSV2);
0231     }
0232     if (_NminOneCuts[5] && nMinOne["CSV3"]) {
0233       this->fillHist(isPassTrigger, source, EVTColContainer::getTypeString(EVTColContainer::PFJET), "CSV3", CSV3);
0234     }
0235   }
0236 }
0237
0238 void HLTHiggsPlotter::bookHist(const std::string &source,
0239                                const std::string &objType,
0240                                const std::string &variable,
0241                                DQMStore::IBooker &ibooker) {
0242   std::string sourceUpper = source;
0243   sourceUpper[0] = std::toupper(sourceUpper[0]);
0244   std::string name = source + objType + variable + "_" + _hltPath;
0245   TH1F *h = nullptr;
0246
0247   if (variable.find("MaxPt") != std::string::npos) {
0248     std::string desc;
0249     if (variable == "MaxPt1")
0250       desc = "Leading";
0251     else if (variable == "MaxPt2")
0252       desc = "Next-to-Leading";
0253     else
0254       desc = variable.substr(5, 6) + "th Leading";
0255     std::string title = "pT of " + desc + " " + sourceUpper + " " + objType +
0256                         " "
0257                         "where event pass the " +
0258                         _hltPath;
0259     const size_t nBinsStandard = _parametersTurnOn.size() - 1;
0260     size_t nBins = nBinsStandard;
0261     float *edges = new float[nBinsStandard + 1];
0262     for (size_t i = 0; i < nBinsStandard + 1; i++) {
0263       edges[i] = _parametersTurnOn[i];
0264     }
0265
0266     std::string jetObj = EVTColContainer::getTypeString(EVTColContainer::PFJET);
0267     if (objType == jetObj) {
0268       const size_t nBinsJets = 25;
0269       nBins = nBinsJets;
0270       delete[] edges;
0271       edges = new float[nBinsJets + 1];
0272       for (size_t i = 0; i < nBinsJets + 1; i++) {
0273         edges[i] = i * 10;
0274       }
0275     }
0276     if (objType == EVTColContainer::getTypeString(EVTColContainer::PFMET)) {
0277       const size_t nBinsJets = 30;
0278       nBins = nBinsJets;
0279       delete[] edges;
0280       edges = new float[nBinsJets + 1];
0281       for (size_t i = 0; i < nBinsJets + 1; i++) {
0282         edges[i] = i * 10;
0283       }
0284     }
0285     h = new TH1F(name.c_str(), title.c_str(), nBins, edges);
0286     delete[] edges;
0287   } else {
0288     if (variable == "dEtaqq") {
0289       std::string title = "#Delta #eta_{qq} of " + sourceUpper + " " + objType;
0290       int nBins = 20;
0291       double min = 0;
0292       double max = 4.8;
0293       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0294     } else if (variable == "mqq") {
0295       std::string title = "m_{qq} of " + sourceUpper + " " + objType;
0296       int nBins = 20;
0297       double min = 0;
0298       double max = 1000;
0299       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0300     } else if (variable == "dPhibb") {
0301       std::string title = "#Delta #phi_{bb} of " + sourceUpper + " " + objType;
0302       int nBins = 10;
0303       double min = 0;
0304       double max = 3.1416;
0305       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0306     } else if (variable == "CSV1") {
0307       std::string title = "CSV1 of " + sourceUpper + " " + objType;
0308       int nBins = 20;
0309       double min = 0;
0310       double max = 1;
0311       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0312     } else if (variable == "CSV2") {
0313       std::string title = "CSV2 of " + sourceUpper + " " + objType;
0314       int nBins = 20;
0315       double min = 0;
0316       double max = 1;
0317       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0318     } else if (variable == "CSV3") {
0319       std::string title = "CSV3 of " + sourceUpper + " " + objType;
0320       int nBins = 20;
0321       double min = 0;
0322       double max = 1;
0323       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0324     } else if (variable == "maxCSV") {
0325       std::string title = "max CSV of " + sourceUpper + " " + objType;
0326       int nBins = 20;
0327       double min = 0;
0328       double max = 1;
0329       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0330     } else {
0331       std::string symbol = (variable == "Eta") ? "#eta" : "#phi";
0332       std::string title = symbol + " of " + sourceUpper + " " + objType + " " + "where event pass the " + _hltPath;
0333       std::vector<double> params = (variable == "Eta") ? _parametersEta : _parametersPhi;
0334
0335       int nBins = (int)params[0];
0336       double min = params[1];
0337       double max = params[2];
0338       h = new TH1F(name.c_str(), title.c_str(), nBins, min, max);
0339     }
0340   }
0341   h->Sumw2();
0342   _elements[name] = ibooker.book1D(name.c_str(), h);
0343   delete h;
0344 }
0345
0346 void HLTHiggsPlotter::fillHist(const bool &passTrigger,
0347                                const std::string &source,
0348                                const std::string &objType,
0349                                const std::string &variable,
0350                                const float &value) {
0351   std::string sourceUpper = source;
0352   sourceUpper[0] = toupper(sourceUpper[0]);
0353   std::string name = source + objType + variable + "_" + _hltPath;
0354
0355   _elements[name]->Fill(value);
0356 }