text: #include "DQMOffline/RecoB/plugins/BTagPerformanceAnalyzerOnData.h"
  #include "FWCore/Framework/interface/MakerMacros.h"
  #include "FWCore/Framework/interface/Event.h"
  #include "DQMOffline/RecoB/interface/TagInfoPlotterFactory.h"
  
  using namespace reco;
  using namespace edm;
  using namespace std;
  using namespace RecoBTag;
  
  BTagPerformanceAnalyzerOnData::BTagPerformanceAnalyzerOnData(const edm::ParameterSet& pSet) :
    jetSelector(
      pSet.getParameter<double>("etaMin"),
      pSet.getParameter<double>("etaMax"),
      pSet.getParameter<double>("ptRecJetMin"),
      pSet.getParameter<double>("ptRecJetMax"),
      0.0, 99999.0,
      pSet.getParameter<double>("ratioMin"),
      pSet.getParameter<double>("ratioMax"),
      pSet.getParameter<bool>( "doJetID" )
    ),
    etaRanges(pSet.getParameter< vector<double> >("etaRanges")),
    ptRanges(pSet.getParameter< vector<double> >("ptRanges")),
    doJEC(pSet.getParameter<bool>( "doJEC" )),
    moduleConfig(pSet.getParameter< vector<edm::ParameterSet> >("tagConfig"))
  {
    genToken = mayConsume<GenEventInfoProduct>(edm::InputTag("generator"));
    slInfoToken = consumes<SoftLeptonTagInfoCollection>(pSet.getParameter<InputTag>("softLeptonInfo"));
    jecMCToken = mayConsume<JetCorrector>(pSet.getParameter<edm::InputTag>( "JECsourceMC" ));
    jecDataToken = consumes<JetCorrector>(pSet.getParameter<edm::InputTag>( "JECsourceData" ));
    
    if (etaRanges.size() <= 1)
        etaRanges = { pSet.getParameter<double>("etaMin"), pSet.getParameter<double>("etaMax") };
    if (ptRanges.size() <= 1)
        ptRanges = { pSet.getParameter<double>("ptRecJetMin"), pSet.getParameter<double>("ptRecJetMax") };
    
    for (vector<edm::ParameterSet>::const_iterator iModule = moduleConfig.begin();
         iModule != moduleConfig.end(); ++iModule) {
  
      const string& dataFormatType = iModule->exists("type") ?
        iModule->getParameter<string>("type") :
        "JetTag";
      if (dataFormatType == "JetTag") {
        const InputTag& moduleLabel = iModule->getParameter<InputTag>("label");
        jetTagInputTags.push_back(moduleLabel);
        binJetTagPlotters.push_back(vector<std::unique_ptr<JetTagPlotter>>()) ;
        jetTagToken.push_back(consumes<JetTagCollection>(moduleLabel));
      }
      else if(dataFormatType == "TagCorrelation") {
        const InputTag& label1 = iModule->getParameter<InputTag>("label1");
        const InputTag& label2 = iModule->getParameter<InputTag>("label2");
        tagCorrelationInputTags.push_back(std::pair<edm::InputTag, edm::InputTag>(label1, label2));
        binTagCorrelationPlotters.push_back(vector<std::unique_ptr<TagCorrelationPlotter>>());
        tagCorrelationToken.push_back(std::pair< edm::EDGetTokenT<reco::JetTagCollection>, edm::EDGetTokenT<reco::JetTagCollection> >(consumes<JetTagCollection>(label1), consumes<JetTagCollection>(label2)));
      }
      else {
        vector<edm::InputTag> vIP;
        tiDataFormatType.push_back(dataFormatType);
        binTagInfoPlotters.push_back(vector<std::unique_ptr<BaseTagInfoPlotter>>()) ;
        std::vector< edm::EDGetTokenT<edm::View<reco::BaseTagInfo>> > tokens;
        if(dataFormatType == "GenericMVA") {
          const std::vector<InputTag> listInfo = iModule->getParameter<vector<InputTag>>("listTagInfos"); 
          for(unsigned int ITi=0; ITi<listInfo.size(); ITi++) { 
            tokens.push_back(consumes< View<BaseTagInfo> >(listInfo[ITi]));
            vIP.push_back(listInfo[ITi]);
          }
        }
        else {
          const InputTag& moduleLabel = iModule->getParameter<InputTag>("label");
          tokens.push_back(consumes< View<BaseTagInfo> >(moduleLabel));
          vIP.push_back(moduleLabel);
        }
        tagInfoToken.push_back(tokens);
        tagInfoInputTags.push_back(vIP);
      }
    }
  }
  
  void BTagPerformanceAnalyzerOnData::bookHistograms(DQMStore::IBooker & ibook, edm::Run const & run, edm::EventSetup const & es)
  {
    // Book all histograms.
  
    // iterate over ranges:
    const int iEtaStart = -1                   ;  // this will be the inactive one
    const int iEtaEnd   = etaRanges.size() > 2 ? etaRanges.size() - 1 : 0; // if there is only one bin defined, leave it as the inactive one
    const int iPtStart  = -1                   ;  // this will be the inactive one
    const int iPtEnd    = ptRanges.size() > 2 ? ptRanges.size() - 1 : 0; // if there is only one bin defined, leave it as the inactive one
    setTDRStyle();
  
    TagInfoPlotterFactory theFactory;
    int iTag = -1; int iTagCorr = -1; int iInfoTag = -1;
    for (vector<edm::ParameterSet>::const_iterator iModule = moduleConfig.begin(); 
         iModule != moduleConfig.end(); ++iModule) {
  
      const string& dataFormatType = iModule->exists("type") ? 
                                     iModule->getParameter<string>("type") :
                                     "JetTag";
      if (dataFormatType == "JetTag") {
        iTag++;
        const string& folderName    = iModule->getParameter<string>("folder");
        
        bool doDifferentialPlots = false;
        double discrCut = -999.;
        if (iModule->exists("differentialPlots") && iModule->getParameter<bool>("differentialPlots") == true) {
            doDifferentialPlots = true;
            discrCut = iModule->getParameter<double>("discrCut");
        }
        
        // eta loop
        for (int iEta = iEtaStart ; iEta < iEtaEnd ; ++iEta) {
          // pt loop
          for (int iPt = iPtStart ; iPt < iPtEnd ; ++iPt) {
            
            const EtaPtBin& etaPtBin = getEtaPtBin(iEta, iPt);
            
            // Instantiate the generic b tag plotter
            binJetTagPlotters.at(iTag).push_back(std::make_unique<JetTagPlotter>(folderName, etaPtBin, iModule->getParameter<edm::ParameterSet>("parameters"), 0, false, ibook, false, doDifferentialPlots, discrCut));
            
          }
        }
        
      } else if (dataFormatType == "TagCorrelation") { 
          iTagCorr++;
          const InputTag& label1 = iModule->getParameter<InputTag>("label1");
          const InputTag& label2 = iModule->getParameter<InputTag>("label2");
  
          // eta loop
          for (int iEta = iEtaStart ; iEta != iEtaEnd ; ++iEta) {
            // pt loop
            for(int iPt = iPtStart ; iPt != iPtEnd ; ++iPt) {
              const EtaPtBin& etaPtBin = getEtaPtBin(iEta, iPt);
              // Instantiate the generic b tag correlation plotter
              binTagCorrelationPlotters.at(iTagCorr).push_back(
                      std::make_unique<TagCorrelationPlotter>(label1.label(), label2.label(), etaPtBin,
                                                              iModule->getParameter<edm::ParameterSet>("parameters"),
                                                              0, false, false, ibook)
                          );
            }
          }
      } else {
        iInfoTag++;
        // tag info retrievel is deferred (needs availability of EventSetup)
        const InputTag& moduleLabel = iModule->getParameter<InputTag>("label");
        const string& folderName    = iModule->getParameter<string>("folder");
        
        // eta loop
        for (int iEta = iEtaStart ; iEta < iEtaEnd ; ++iEta) {
          // pt loop
          for (int iPt = iPtStart ; iPt < iPtEnd ; ++iPt) {
            const EtaPtBin& etaPtBin = getEtaPtBin(iEta, iPt);
            
            // Instantiate the tagInfo plotter
            binTagInfoPlotters.at(iInfoTag).push_back(theFactory.buildPlotter(dataFormatType, moduleLabel.label(),
                                            etaPtBin, iModule->getParameter<edm::ParameterSet>("parameters"), folderName,
                                            0, false, ibook)
                    );
          }
        }
      }
    }
  }
  
  EtaPtBin BTagPerformanceAnalyzerOnData::getEtaPtBin(const int& iEta, const int& iPt)
  {
    // DEFINE BTagBin:
    bool    etaActive_ , ptActive_;
    double  etaMin_, etaMax_, ptMin_, ptMax_;
    
    if ( iEta != -1 ) {
      etaActive_ = true;
      etaMin_    = etaRanges[iEta];
      etaMax_    = etaRanges[iEta+1];
    }
    else {
      etaActive_ = false;
      etaMin_    = etaRanges[0];
      etaMax_    = etaRanges[etaRanges.size() - 1];
    }
  
    if ( iPt != -1 ) {
      ptActive_ = true;
      ptMin_    = ptRanges[iPt];
      ptMax_    = ptRanges[iPt+1];
    }
    else {
      ptActive_ = false;
      ptMin_    = ptRanges[0];
      ptMax_    = ptRanges[ptRanges.size() - 1];
    }
    return EtaPtBin(etaActive_, etaMin_, etaMax_, ptActive_, ptMin_ , ptMax_);
  }
  
  BTagPerformanceAnalyzerOnData::~BTagPerformanceAnalyzerOnData()
  {
    /*for (unsigned int iJetLabel = 0; iJetLabel != binJetTagPlotters.size(); ++iJetLabel) {
      int plotterSize =  binJetTagPlotters[iJetLabel].size();
      for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
        delete binJetTagPlotters[iJetLabel][iPlotter];
      }
    }
  
    for (vector<vector<TagCorrelationPlotter*> >::iterator iJetLabel = binTagCorrelationPlotters.begin();
         iJetLabel != binTagCorrelationPlotters.end(); ++iJetLabel) 
      for(vector<TagCorrelationPlotter*>::iterator iPlotter = iJetLabel->begin(); iPlotter != iJetLabel->end(); ++iPlotter)
        delete *iPlotter;
  
    for (unsigned int iJetLabel = 0; iJetLabel != binTagInfoPlotters.size(); ++iJetLabel) {
      int plotterSize =  binTagInfoPlotters[iJetLabel].size();
      for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
        delete binTagInfoPlotters[iJetLabel][iPlotter];
      }
    }*/
  }
  
  void BTagPerformanceAnalyzerOnData::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
  {
  
    edm::Handle<reco::SoftLeptonTagInfoCollection> infoHandle;
    iEvent.getByToken(slInfoToken, infoHandle);
  
    //Get JEC
    const JetCorrector* corrector = nullptr;
    if(doJEC) {
      edm::Handle<GenEventInfoProduct> genInfoHandle; //check if data or MC                                                                                                                                 
      iEvent.getByToken(genToken, genInfoHandle);
      edm::Handle<JetCorrector> corrHandle;
      if( !genInfoHandle.isValid() ) iEvent.getByToken(jecDataToken, corrHandle);
      else iEvent.getByToken(jecMCToken, corrHandle);
      corrector = corrHandle.product();
    }
  
    // Look first at the jetTags
  
    for (unsigned int iJetLabel = 0; iJetLabel != jetTagInputTags.size(); ++iJetLabel) {
      edm::Handle<reco::JetTagCollection> tagHandle;
      iEvent.getByToken(jetTagToken[iJetLabel], tagHandle);
      //
      // insert check on the presence of the collections
      //
  
      if (!tagHandle.isValid()) {
        edm::LogWarning("BTagPerformanceAnalyzerOnData")<<" Collection "<<jetTagInputTags[iJetLabel]<<" not present. Skipping it for this event.";
        continue;
      }
      
      const reco::JetTagCollection & tagColl = *(tagHandle.product());
      LogDebug("Info") << "Found " << tagColl.size() << " B candidates in collection " << jetTagInputTags[iJetLabel];
  
      int plotterSize =  binJetTagPlotters[iJetLabel].size();
      for (JetTagCollection::const_iterator tagI = tagColl.begin(); tagI != tagColl.end(); ++tagI) {
        
        //JEC
        double jec = 1.0;
        if(doJEC && corrector) {
          jec = corrector->correction(*(tagI->first));
        }
  
        if (!jetSelector(*(tagI->first), -1, infoHandle, jec)) 
          continue;
        
        for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
          bool inBin = binJetTagPlotters[iJetLabel][iPlotter]->etaPtBin().inBin(*tagI->first, jec);
          // Fill histograms if in desired pt/rapidity bin.
          if (inBin)
            binJetTagPlotters[iJetLabel][iPlotter]->analyzeTag(*tagI, jec, -1);
        }
      }
      for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
        binJetTagPlotters[iJetLabel][iPlotter]->analyzeTag();
      }
    }
  
    // Now look at Tag Correlations
    for (unsigned int iJetLabel = 0; iJetLabel != tagCorrelationInputTags.size(); ++iJetLabel) {
      const std::pair< edm::EDGetTokenT<reco::JetTagCollection>, edm::EDGetTokenT<reco::JetTagCollection> >& inputTokens = tagCorrelationToken[iJetLabel];
      edm::Handle<reco::JetTagCollection> tagHandle1;
      iEvent.getByToken(inputTokens.first, tagHandle1);
      const reco::JetTagCollection& tagColl1 = *(tagHandle1.product());
  
      edm::Handle<reco::JetTagCollection> tagHandle2;
      iEvent.getByToken(inputTokens.second, tagHandle2);
      const reco::JetTagCollection& tagColl2 = *(tagHandle2.product());
  
      int plotterSize = binTagCorrelationPlotters[iJetLabel].size();
      for (JetTagCollection::const_iterator tagI = tagColl1.begin(); tagI != tagColl1.end(); ++tagI) {
        //JEC
        double jec = 1.0;
        if(doJEC && corrector) {
          jec = corrector->correction(*(tagI->first));
        }
  
        if (!jetSelector(*(tagI->first), -1, infoHandle, jec))
          continue;
        
        for(int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
          bool inBin = binTagCorrelationPlotters[iJetLabel][iPlotter]->etaPtBin().inBin(*(tagI->first), jec);
  
          if(inBin)
          {
            double discr2 = tagColl2[tagI->first];
            binTagCorrelationPlotters[iJetLabel][iPlotter]->analyzeTags(tagI->second, discr2, -1);
          }
        }
      }
    }
  
    // Now look at the TagInfos
    for (unsigned int iJetLabel = 0; iJetLabel != tiDataFormatType.size(); ++iJetLabel) {
      int plotterSize = binTagInfoPlotters[iJetLabel].size();
      for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter)
        binTagInfoPlotters[iJetLabel][iPlotter]->setEventSetup(iSetup);
  
      vector<edm::EDGetTokenT<edm::View<reco::BaseTagInfo>> > & tokens = tagInfoToken[iJetLabel];
      //check number of tag infos = expected number of tag infos
      vector<string> labels = binTagInfoPlotters[iJetLabel][0]->tagInfoRequirements();
      if (labels.empty())
        labels.push_back("label");
      if (labels.size() != tokens.size()) throw cms::Exception("Configuration") << "Different number of Tag Infos than expected" << labels.size() << tokens.size() << endl;
  
      unsigned int nInputTags = tokens.size();
      vector< edm::Handle< View<BaseTagInfo> > > tagInfoHandles(nInputTags);
      edm::ProductID jetProductID;
      unsigned int nTagInfos = 0;
      for (unsigned int iInputTags = 0; iInputTags < tokens.size(); ++iInputTags) {
        edm::Handle< View<BaseTagInfo> > & tagInfoHandle = tagInfoHandles[iInputTags];
        iEvent.getByToken(tokens[iInputTags], tagInfoHandle); 
        //
        // protect against missing products
        //
        if (tagInfoHandle.isValid() == false) {
          edm::LogWarning("BTagPerformanceAnalyzerOnData")<<" Collection "<<tagInfoInputTags[iJetLabel][iInputTags]<<" not present. Skipping it for this event.";
          continue;
        }
  
        unsigned int size = tagInfoHandle->size();
        LogDebug("Info") << "Found " << size << " B candidates in collection " << tagInfoInputTags[iJetLabel][iInputTags];
        edm::ProductID thisProductID = (size > 0) ? (*tagInfoHandle)[0].jet().id() : edm::ProductID();
        if (iInputTags == 0) {
          jetProductID = thisProductID;
          nTagInfos = size;
        } else if (jetProductID != thisProductID)
          throw cms::Exception("Configuration") << "TagInfos are referencing a different jet collection." << endl;
        else if (nTagInfos != size)
          throw cms::Exception("Configuration") << "TagInfo collections are having a different size." << endl;
      }
  
      for (unsigned int iTagInfos = 0; iTagInfos < nTagInfos; ++iTagInfos) {
        vector<const BaseTagInfo*> baseTagInfos(nInputTags);
        edm::RefToBase<Jet> jetRef;
        for (unsigned int iTagInfo = 0; iTagInfo < nInputTags; iTagInfo++) {
          const BaseTagInfo &baseTagInfo = (*tagInfoHandles[iTagInfo])[iTagInfos];
          if (iTagInfo == 0)
            jetRef = baseTagInfo.jet();
          else if (baseTagInfo.jet() != jetRef)
            throw cms::Exception("Configuration") << "TagInfos pointing to different jets." << endl;
          baseTagInfos[iTagInfo] = &baseTagInfo;
        }
  
        //JEC
        double jec = 1.0;
        if(doJEC && corrector) {
          jec = corrector->correction(*(jetRef));
        }
  
        if (!jetSelector(*jetRef, -1, infoHandle, jec))
          continue;
        
        for (int iPlotter = 0; iPlotter != plotterSize; ++iPlotter) {
          bool inBin = binTagInfoPlotters[iJetLabel][iPlotter]->etaPtBin().inBin(*jetRef, jec);
          // Fill histograms if in desired pt/rapidity bin.
          if (inBin)
            binTagInfoPlotters[iJetLabel][iPlotter]->analyzeTag(baseTagInfos, jec, -1);
        }
      }
    }
  }
  
  //define this as a plug-in
  DEFINE_FWK_MODULE(BTagPerformanceAnalyzerOnData);
  