Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:09:34 

0001 #ifndef DQMOFFLINE_L1TRIGGER_L1TPHASE2CORRELATOROFFLINE_H
0002 #define DQMOFFLINE_L1TRIGGER_L1TPHASE2CORRELATOROFFLINE_H
0003 
0004 // DataFormats
0005 #include "DataFormats/Common/interface/View.h"
0006 #include "DataFormats/Candidate/interface/Candidate.h"
0007 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0008 #include "DataFormats/JetReco/interface/GenJet.h"
0009 #include "DataFormats/Math/interface/deltaR.h"
0010 #include "DataFormats/L1TParticleFlow/interface/PFCandidate.h"
0011 
0012 // FWCore
0013 #include "FWCore/Common/interface/Provenance.h"
0014 #include "FWCore/Framework/interface/Frameworkfwd.h"
0015 #include "FWCore/Framework/interface/Event.h"
0016 #include "FWCore/Framework/interface/LuminosityBlock.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0019 #include "FWCore/Framework/interface/ConsumesCollector.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 
0022 // DQMServices
0023 #include "DQMServices/Core/interface/DQMOneEDAnalyzer.h"
0024 #include "DQMServices/Core/interface/DQMStore.h"
0025 #include "DQMOffline/L1Trigger/interface/HistDefinition.h"
0026 
0027 // MagneticField
0028 #include "MagneticField/Engine/interface/MagneticField.h"
0029 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0030 
0031 // HLTrigger
0032 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0033 
0034 // L1Trigger
0035 #include "L1Trigger/Phase2L1ParticleFlow/interface/L1TPFUtils.h"
0036 
0037 // CommonTools
0038 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0039 
0040 #include <iostream>
0041 #include <fstream>
0042 #include <utility>
0043 #include <vector>
0044 #include <TRandom3.h>
0045 #include <TTree.h>
0046 #include <cstdint>
0047 
0048 class L1TPhase2CorrelatorOffline : public DQMOneEDAnalyzer<> {
0049 public:
0050   L1TPhase2CorrelatorOffline(const edm::ParameterSet& ps);
0051   ~L1TPhase2CorrelatorOffline() override;
0052 
0053   enum PlotConfig { resVsPt, resVsEta, ptDist, etaDist };
0054 
0055   static const std::map<std::string, unsigned int> PlotConfigNames;
0056 
0057 protected:
0058   void dqmBeginRun(const edm::Run& run, const edm::EventSetup& iSetup) override;
0059   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0060   void analyze(edm::Event const& e, edm::EventSetup const& eSetup) override;
0061   void dqmEndRun(const edm::Run& run, const edm::EventSetup& iSetup) override;
0062 
0063   // Cut and Matching
0064 
0065 private:
0066   void bookPhase2CorrelatorHistos(DQMStore::IBooker&);
0067 
0068   // other functions
0069   double Distance(const reco::Candidate& c1, const reco::Candidate& c2);
0070   double DistancePhi(const reco::Candidate& c1, const reco::Candidate& c2);
0071   double calcDeltaPhi(double phi1, double phi2);
0072 
0073   void computeResponseResolution();
0074 
0075   std::vector<float> getQuantile(float quant, TH2F* hist);
0076   void medianResponseCorrResolution(MonitorElement* in2D, MonitorElement* response, MonitorElement* resolution);
0077 
0078   struct SimpleObject {
0079     float pt, eta, phi;
0080     SimpleObject(float apt, float aneta, float aphi) : pt(apt), eta(aneta), phi(aphi) {}
0081     bool operator<(const SimpleObject& other) const { return eta < other.eta; }
0082     bool operator<(const float& other) const { return eta < other; }
0083   };
0084   class MultiCollection {
0085   public:
0086     MultiCollection(const edm::ParameterSet& iConfig, const std::string& name, edm::ConsumesCollector&& coll)
0087         : name_(name), prop_(false), sel_("") {
0088       if (name.find("Ecal") != std::string::npos)
0089         prop_ = true;
0090       else if (name.find("Hcal") != std::string::npos)
0091         prop_ = true;
0092       else if (name.find("Calo") != std::string::npos)
0093         prop_ = true;
0094       const std::vector<edm::InputTag>& tags = iConfig.getParameter<std::vector<edm::InputTag>>(name);
0095       for (const auto& tag : tags)
0096         tokens_.push_back(coll.consumes<reco::CandidateView>(tag));
0097       if (iConfig.existsAs<std::string>(name + "_sel")) {
0098         sel_ = StringCutObjectSelector<reco::Candidate>(iConfig.getParameter<std::string>(name + "_sel"), true);
0099       }
0100     }
0101     const std::string& name() const { return name_; }
0102     bool prop() const { return prop_; }
0103     void get(const edm::Event& iEvent) {
0104       edm::Handle<reco::CandidateView> handle;
0105       for (const auto& token : tokens_) {
0106         iEvent.getByToken(token, handle);
0107         for (const reco::Candidate& c : *handle) {
0108           if (sel_(c))
0109             objects_.emplace_back(c.pt(), c.eta(), c.phi());
0110         }
0111       }
0112       std::sort(objects_.begin(), objects_.end());
0113     }
0114     const std::vector<SimpleObject>& objects() const { return objects_; }
0115     void clear() { objects_.clear(); }
0116 
0117   private:
0118     std::string name_;
0119     bool prop_;
0120     std::vector<edm::EDGetTokenT<reco::CandidateView>> tokens_;
0121     StringCutObjectSelector<reco::Candidate> sel_;
0122     std::vector<SimpleObject> objects_;
0123   };
0124   class InCone {
0125   public:
0126     InCone(const std::vector<SimpleObject>& objects, float eta, float phi, float dr) {
0127       auto first =
0128           std::lower_bound(objects.begin(), objects.end(), eta - dr - 0.01f);  // small offset to avoid dealing with ==
0129       auto end = std::lower_bound(objects.begin(), objects.end(), eta + dr + 0.01f);
0130       float dr2 = dr * dr;
0131       sum04_ = 0;
0132       for (auto it = first; it < end; ++it) {
0133         float mydr2 = ::deltaR2(eta, phi, it->eta, it->phi);
0134         if (mydr2 < dr2)
0135           ptdr2_.emplace_back(it->pt, mydr2);
0136         if (mydr2 < 0.16f)
0137           sum04_ += it->pt;
0138       }
0139     }
0140     float sum(float dr = 0.4) const {
0141       if (dr == 0.4f)
0142         return sum04_;
0143       float dr2 = dr * dr;
0144       float mysum = 0;
0145       for (const auto& part : ptdr2_) {
0146         if (part.second < dr2)
0147           mysum += part.first;
0148       }
0149       return mysum;
0150     }
0151     int number(float dr, float threshold) const {
0152       float dr2 = dr * dr, absthreshold = sum() * threshold;
0153       int mysum = 0;
0154       for (const auto& part : ptdr2_) {
0155         if (part.second < dr2 && part.first > absthreshold)
0156           mysum++;
0157       }
0158       return mysum;
0159     }
0160     float mindr(float threshold) const {
0161       float best = 9999, absthreshold = sum() * threshold;
0162       for (const auto& part : ptdr2_) {
0163         if (part.second < best && part.first > absthreshold)
0164           best = part.second;
0165       }
0166       return std::sqrt(best);
0167     }
0168     float nearest() const {
0169       std::pair<float, float> best(0, 9999);
0170       for (const auto& part : ptdr2_) {
0171         if (part.second < best.second)
0172           best = part;
0173       }
0174       return best.first;
0175     }
0176     float max(float dr = 0.4) const {
0177       float best = 0, dr2 = dr * dr;
0178       for (const auto& part : ptdr2_) {
0179         if (part.first > best && part.second < dr2)
0180           best = part.first;
0181       }
0182       return best;
0183     }
0184 
0185   private:
0186     std::vector<std::pair<float, float>> ptdr2_;
0187     float sum04_;
0188   };
0189 
0190   // variables from config file
0191   edm::EDGetTokenT<std::vector<reco::GenJet>> genJetToken_;
0192   edm::EDGetTokenT<std::vector<reco::GenParticle>> genParticleToken_;
0193   edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> BFieldTag_;
0194   edm::ParameterSet objs_;
0195   bool isParticleGun_;
0196   std::string histFolder_;
0197   std::string respresolFolder_;
0198   dqmoffline::l1t::HistDefinitions histDefinitions_;
0199 
0200   std::vector<edm::EDGetTokenT<std::vector<l1t::PFCandidate>>> phase2PFToken_;
0201   std::vector<edm::EDGetTokenT<std::vector<l1t::PFCandidate>>> phase2PuppiToken_;
0202 
0203   // config params
0204   struct McVars {
0205     float pt, pt02, eta, phi, iso02, iso04, iso08;
0206     int charge;
0207     float caloeta, calophi;
0208     int id;
0209     void fillP4(const reco::Candidate& cand) {
0210       pt = cand.pt();
0211       eta = cand.eta();
0212       phi = cand.phi();
0213       caloeta = eta;
0214       calophi = phi;
0215       charge = 0;
0216     }
0217     void fillPropagated(const reco::Candidate& cand, float bz) {
0218       if (cand.charge() != 0) {
0219         math::XYZTLorentzVector vertex(cand.vx(), cand.vy(), cand.vz(), 0.);
0220         auto caloetaphi = l1tpf::propagateToCalo(cand.p4(), vertex, cand.charge(), bz);
0221         caloeta = caloetaphi.first;
0222         calophi = caloetaphi.second;
0223       }
0224     }
0225 
0226   } mc_;
0227   struct RecoVars {
0228     float pt, pt02, pt08, ptbest, pthighest, mindr025;
0229     int n025, n010;
0230     void fill(const std::vector<SimpleObject>& objects, float eta, float phi) {
0231       InCone incone(objects, eta, phi, 0.8);
0232       pt = incone.sum();
0233       pt02 = incone.sum(0.2);
0234       pt08 = incone.sum(0.8);
0235       ptbest = incone.nearest();
0236       pthighest = incone.max();
0237       mindr025 = incone.mindr(0.25);
0238       n025 = incone.number(0.2, 0.25);
0239       n010 = incone.number(0.2, 0.10);
0240     }
0241     void clear() {
0242       pt = 0.;
0243       pt02 = 0.;
0244       pt08 = 0.;
0245       ptbest = 0.;
0246       pthighest = 0.;
0247       mindr025 = 0.;
0248       n025 = -1;
0249       n010 = -1;
0250     }
0251   };
0252   std::vector<std::pair<MultiCollection, RecoVars>> reco_;
0253   float bZ_;
0254 
0255   // Histograms
0256   MonitorElement* h_L1PF_pt_;
0257   MonitorElement* h_L1PF_eta_;
0258   MonitorElement* h_L1Puppi_pt_;
0259   MonitorElement* h_L1Puppi_eta_;
0260 
0261   MonitorElement* h_L1PF_pt_mu_;
0262   MonitorElement* h_L1PF_eta_mu_;
0263   MonitorElement* h_L1Puppi_pt_mu_;
0264   MonitorElement* h_L1Puppi_eta_mu_;
0265 
0266   MonitorElement* h_L1PF_pt_el_;
0267   MonitorElement* h_L1PF_eta_el_;
0268   MonitorElement* h_L1Puppi_pt_el_;
0269   MonitorElement* h_L1Puppi_eta_el_;
0270 
0271   MonitorElement* h_L1PF_pt_pho_;
0272   MonitorElement* h_L1PF_eta_pho_;
0273   MonitorElement* h_L1Puppi_pt_pho_;
0274   MonitorElement* h_L1Puppi_eta_pho_;
0275 
0276   MonitorElement* h_L1PF_pt_ch_;
0277   MonitorElement* h_L1PF_eta_ch_;
0278   MonitorElement* h_L1Puppi_pt_ch_;
0279   MonitorElement* h_L1Puppi_eta_ch_;
0280 
0281   MonitorElement* h_L1PF_pt_nh_;
0282   MonitorElement* h_L1PF_eta_nh_;
0283   MonitorElement* h_L1Puppi_pt_nh_;
0284   MonitorElement* h_L1Puppi_eta_nh_;
0285 
0286   MonitorElement* h_L1PF_part_ptratio_0p2_vs_pt_barrel_;
0287   MonitorElement* h_L1PF_part_ptratio_0p2_vs_pt_endcap_;
0288   MonitorElement* h_L1PF_part_ptratio_0p2_vs_pt_ecnotk_;
0289   MonitorElement* h_L1PF_part_ptratio_0p2_vs_pt_hf_;
0290   MonitorElement* h_L1PF_part_ptratio_0p2_vs_eta_;
0291   MonitorElement* h_L1Puppi_part_ptratio_0p2_vs_pt_barrel_;
0292   MonitorElement* h_L1Puppi_part_ptratio_0p2_vs_pt_endcap_;
0293   MonitorElement* h_L1Puppi_part_ptratio_0p2_vs_pt_ecnotk_;
0294   MonitorElement* h_L1Puppi_part_ptratio_0p2_vs_pt_hf_;
0295   MonitorElement* h_L1Puppi_part_ptratio_0p2_vs_eta_;
0296   MonitorElement* h_L1PF_jet_ptratio_vs_pt_barrel_;
0297   MonitorElement* h_L1PF_jet_ptratio_vs_pt_endcap_;
0298   MonitorElement* h_L1PF_jet_ptratio_vs_pt_ecnotk_;
0299   MonitorElement* h_L1PF_jet_ptratio_vs_pt_hf_;
0300   MonitorElement* h_L1PF_jet_ptratio_vs_eta_;
0301   MonitorElement* h_L1Puppi_jet_ptratio_vs_pt_barrel_;
0302   MonitorElement* h_L1Puppi_jet_ptratio_vs_pt_endcap_;
0303   MonitorElement* h_L1Puppi_jet_ptratio_vs_pt_ecnotk_;
0304   MonitorElement* h_L1Puppi_jet_ptratio_vs_pt_hf_;
0305   MonitorElement* h_L1Puppi_jet_ptratio_vs_eta_;
0306 
0307   MonitorElement* h_L1PF_part_response_0p2_pt_barrel_;
0308   MonitorElement* h_L1PF_part_response_0p2_pt_endcap_;
0309   MonitorElement* h_L1PF_part_response_0p2_pt_ecnotk_;
0310   MonitorElement* h_L1PF_part_response_0p2_pt_hf_;
0311   MonitorElement* h_L1PF_part_response_0p2_eta_;
0312   MonitorElement* h_L1Puppi_part_response_0p2_pt_barrel_;
0313   MonitorElement* h_L1Puppi_part_response_0p2_pt_endcap_;
0314   MonitorElement* h_L1Puppi_part_response_0p2_pt_ecnotk_;
0315   MonitorElement* h_L1Puppi_part_response_0p2_pt_hf_;
0316   MonitorElement* h_L1Puppi_part_response_0p2_eta_;
0317   MonitorElement* h_L1PF_jet_response_pt_barrel_;
0318   MonitorElement* h_L1PF_jet_response_pt_endcap_;
0319   MonitorElement* h_L1PF_jet_response_pt_ecnotk_;
0320   MonitorElement* h_L1PF_jet_response_pt_hf_;
0321   MonitorElement* h_L1PF_jet_response_eta_;
0322   MonitorElement* h_L1Puppi_jet_response_pt_barrel_;
0323   MonitorElement* h_L1Puppi_jet_response_pt_endcap_;
0324   MonitorElement* h_L1Puppi_jet_response_pt_ecnotk_;
0325   MonitorElement* h_L1Puppi_jet_response_pt_hf_;
0326   MonitorElement* h_L1Puppi_jet_response_eta_;
0327 
0328   MonitorElement* h_L1PF_part_resolution_0p2_pt_barrel_;
0329   MonitorElement* h_L1PF_part_resolution_0p2_pt_endcap_;
0330   MonitorElement* h_L1PF_part_resolution_0p2_pt_ecnotk_;
0331   MonitorElement* h_L1PF_part_resolution_0p2_pt_hf_;
0332   MonitorElement* h_L1Puppi_part_resolution_0p2_pt_barrel_;
0333   MonitorElement* h_L1Puppi_part_resolution_0p2_pt_endcap_;
0334   MonitorElement* h_L1Puppi_part_resolution_0p2_pt_ecnotk_;
0335   MonitorElement* h_L1Puppi_part_resolution_0p2_pt_hf_;
0336   MonitorElement* h_L1PF_jet_resolution_pt_barrel_;
0337   MonitorElement* h_L1PF_jet_resolution_pt_endcap_;
0338   MonitorElement* h_L1PF_jet_resolution_pt_ecnotk_;
0339   MonitorElement* h_L1PF_jet_resolution_pt_hf_;
0340   MonitorElement* h_L1Puppi_jet_resolution_pt_barrel_;
0341   MonitorElement* h_L1Puppi_jet_resolution_pt_endcap_;
0342   MonitorElement* h_L1Puppi_jet_resolution_pt_ecnotk_;
0343   MonitorElement* h_L1Puppi_jet_resolution_pt_hf_;
0344 };
0345 
0346 #endif
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!