PFTau/interface/PFJetMonitor.h

0001 #ifndef DQMOffline_PFTau_PFJetMonitor_h
0002 #define DQMOffline_PFTau_PFJetMonitor_h
0003
0004 #include "DQMOffline/PFTau/interface/Benchmark.h"
0005 #include "DQMOffline/PFTau/interface/CandidateBenchmark.h"
0006 #include "DQMOffline/PFTau/interface/MatchCandidateBenchmark.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008
0009 #include "DataFormats/JetReco/interface/BasicJetCollection.h"
0010
0011 #include <vector>
0012 #include <numeric>    // std::iota
0013 #include <algorithm>  // std::sort
0014
0015 #include <TH1.h>  //needed by the deltaR->Fill() call
0016
0017 class PFJetMonitor : public Benchmark {
0018 public:
0019   PFJetMonitor(float dRMax = 0.3, bool matchCharge = true, Benchmark::Mode mode = Benchmark::DEFAULT);
0020
0021   ~PFJetMonitor() override;
0022
0023   /// set the parameters accessing them from ParameterSet
0024   void setParameters(const edm::ParameterSet &parameterSet);
0025
0026   /// set directory (to use in ROOT)
0027   void setDirectory(TDirectory *dir) override;
0028
0029   /// book histograms
0030   void setup(DQMStore::IBooker &b);
0031   void setup(DQMStore::IBooker &b, const edm::ParameterSet &parameterSet);
0032
0033   /// fill histograms with all particle
0034   template <class T, class C>
0035   void fill(const T &candidateCollection,
0036             const C &matchedCandCollection,
0037             float &minVal,
0038             float &maxVal,
0039             float &jetpT,
0040             const edm::ParameterSet &parameterSet);
0041
0042   void fillOne(const reco::Jet &jet, const reco::Jet &matchedJet);
0043
0044 protected:
0045   CandidateBenchmark candBench_;
0046   MatchCandidateBenchmark matchCandBench_;
0047
0048   TH2F *delta_frac_VS_frac_muon_;
0049   TH2F *delta_frac_VS_frac_photon_;
0050   TH2F *delta_frac_VS_frac_electron_;
0051   TH2F *delta_frac_VS_frac_charged_hadron_;
0052   TH2F *delta_frac_VS_frac_neutral_hadron_;
0053
0054   TH1F *deltaR_;
0055   float dRMax_;
0056   bool onlyTwoJets_;
0057   bool matchCharge_;
0058   bool createPFractionHistos_;
0059   bool histogramBooked_;
0060 };
0061
0062 #include "DQMOffline/PFTau/interface/Matchers.h"
0063
0064 template <class T, class C>
0065 void PFJetMonitor::fill(const T &jetCollection,
0066                         const C &matchedJetCollection,
0067                         float &minVal,
0068                         float &maxVal,
0069                         float &jetpT,
0070                         const edm::ParameterSet &parameterSet) {
0071   std::vector<int> matchIndices;
0072   PFB::match(jetCollection, matchedJetCollection, matchIndices, matchCharge_, dRMax_);
0073   // now matchIndices[i] stores the j-th closest matched jet
0074
0075   std::vector<uint32_t> sorted_pt_indices(jetCollection.size());
0076   std::iota(std::begin(sorted_pt_indices), std::end(sorted_pt_indices), 0);
0077   // Sort the vector of indices using the pt() as ordering variable
0078   std::sort(std::begin(sorted_pt_indices), std::end(sorted_pt_indices), [&](uint32_t i, uint32_t j) {
0079     return jetCollection[i].pt() < jetCollection[j].pt();
0080   });
0081   for (uint32_t i = 0; i < sorted_pt_indices.size(); ++i) {
0082     // If we want only the 2 pt-leading jets, now that they are orderd, simply
0083     // check if the index is either in the first or second location of the
0084     // sorted indices: if not, bail out.
0085     if (onlyTwoJets_ && i > 1)
0086       break;
0087
0088     const reco::Jet &jet = jetCollection[i];
0089
0090     if (!isInRange(jet.pt(), jet.eta(), jet.phi()))
0091       continue;
0092
0093     int iMatch = matchIndices[i];
0094     assert(iMatch < static_cast<int>(matchedJetCollection.size()));
0095
0096     if (iMatch != -1) {
0097       const reco::Jet &matchedJet = matchedJetCollection[iMatch];
0098       if (!isInRange(matchedJet.pt(), matchedJet.eta(), matchedJet.phi()))
0099         continue;
0100
0101       float ptRes = (jet.pt() - matchedJet.pt()) / matchedJet.pt();
0102
0103       jetpT = jet.pt();
0104       if (ptRes > maxVal)
0105         maxVal = ptRes;
0106       if (ptRes < minVal)
0107         minVal = ptRes;
0108
0109       candBench_.fillOne(jet);  // fill pt eta phi and charge histos for MATCHED candidate jet
0110       matchCandBench_.fillOne(jet, matchedJet, parameterSet);  // fill delta_x_VS_y histos for matched couple
0111       if (createPFractionHistos_ && histogramBooked_)
0112         fillOne(jet, matchedJet);  // book and fill delta_frac_VS_frac histos for matched couple
0113     }
0114
0115     for (unsigned j = 0; j < matchedJetCollection.size(); ++j)  // for DeltaR spectrum
0116       if (deltaR_)
0117         deltaR_->Fill(reco::deltaR(jetCollection[i], matchedJetCollection[j]));
0118   }  // end loop on jetCollection
0119 }
0120 #endif