Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 10:55:45

0001 /*
0002  *  See header file for a description of this class.
0003  *
0004  *  \author S. Bolognesi, Eric - CERN
0005  */
0006 
0007 #include "DQM/Physics/src/QcdHighPtDQM.h"
0008 
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "DataFormats/Common/interface/Handle.h"
0013 
0014 #include "DQMServices/Core/interface/DQMStore.h"
0015 #include "FWCore/ServiceRegistry/interface/Service.h"
0016 
0017 #include "DataFormats/METReco/interface/CaloMET.h"
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 #include <vector>
0020 
0021 #include <string>
0022 #include <cmath>
0023 
0024 using namespace std;
0025 using namespace edm;
0026 using namespace reco;
0027 using namespace math;
0028 
0029 // Get Jets and MET (no MET plots yet pending converging w/JetMET group)
0030 
0031 QcdHighPtDQM::QcdHighPtDQM(const ParameterSet &iConfig)
0032     : jetToken_(consumes<CaloJetCollection>(iConfig.getUntrackedParameter<edm::InputTag>("jetTag"))),
0033       metToken1_(consumes<CaloMETCollection>(iConfig.getUntrackedParameter<edm::InputTag>("metTag1"))),
0034       metToken2_(consumes<CaloMETCollection>(iConfig.getUntrackedParameter<edm::InputTag>("metTag2"))),
0035       metToken3_(consumes<CaloMETCollection>(iConfig.getUntrackedParameter<edm::InputTag>("metTag3"))),
0036       metToken4_(consumes<CaloMETCollection>(iConfig.getUntrackedParameter<edm::InputTag>("metTag4"))) {}
0037 
0038 QcdHighPtDQM::~QcdHighPtDQM() {}
0039 
0040 void QcdHighPtDQM::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &, edm::EventSetup const &) {
0041   iBooker.setCurrentFolder("Physics/QcdHighPt");
0042 
0043   MEcontainer_["dijet_mass"] =
0044       iBooker.book1D("dijet_mass", "dijet resonance invariant mass, barrel region", 100, 0, 1000);
0045   MEcontainer_["njets"] = iBooker.book1D("njets", "jet multiplicity", 10, 0, 10);
0046   MEcontainer_["etaphi"] = iBooker.book2D("etaphi", "eta/phi distribution", 83, -42, 42, 72, -M_PI, M_PI);
0047   MEcontainer_["njets30"] = iBooker.book1D("njets30", "jet multiplicity, pt > 30 GeV", 10, 0, 10);
0048 
0049   // book histograms for inclusive jet quantities
0050   MEcontainer_["inclusive_jet_pt"] = iBooker.book1D("inclusive_jet_pt", "inclusive jet Pt spectrum", 100, 0, 1000);
0051   MEcontainer_["inclusive_jet_pt_barrel"] =
0052       iBooker.book1D("inclusive_jet_pt_barrel", "inclusive jet Pt, eta < 1.3", 100, 0, 1000);
0053   MEcontainer_["inclusive_jet_pt_forward"] =
0054       iBooker.book1D("inclusive_jet_pt_forward", "inclusive jet Pt, 3.0 < eta < 5.0", 100, 0, 1000);
0055   MEcontainer_["inclusive_jet_pt_endcap"] =
0056       iBooker.book1D("inclusive_jet_pt_endcap", "inclusive jet Pt, 1.3 < eta < 3.0", 100, 0, 1000);
0057 
0058   // book histograms for leading jet quantities
0059   MEcontainer_["leading_jet_pt"] = iBooker.book1D("leading_jet_pt", "leading jet Pt", 100, 0, 1000);
0060   MEcontainer_["leading_jet_pt_barrel"] =
0061       iBooker.book1D("leading_jet_pt_barrel", "leading jet Pt, eta < 1.3", 100, 0, 1000);
0062   MEcontainer_["leading_jet_pt_forward"] =
0063       iBooker.book1D("leading_jet_pt_forward", "leading jet Pt, 3.0 < eta < 5.0", 100, 0, 1000);
0064   MEcontainer_["leading_jet_pt_endcap"] =
0065       iBooker.book1D("leading_jet_pt_endcap", "leading jet Pt, 1.3 < eta < 3.0", 100, 0, 1000);
0066 
0067   // book histograms for met over sum et and met over leading jet pt for various
0068   // flavors of MET
0069   MEcontainer_["movers_met"] = iBooker.book1D("movers_met", "MET over Sum ET for basic MET collection", 50, 0, 1);
0070   MEcontainer_["moverl_met"] =
0071       iBooker.book1D("moverl_met", "MET over leading jet Pt for basic MET collection", 50, 0, 2);
0072 
0073   MEcontainer_["movers_metho"] = iBooker.book1D("movers_metho", "MET over Sum ET for MET HO collection", 50, 0, 1);
0074   MEcontainer_["moverl_metho"] =
0075       iBooker.book1D("moverl_metho", "MET over leading jet Pt for MET HO collection", 50, 0, 2);
0076 
0077   MEcontainer_["movers_metnohf"] =
0078       iBooker.book1D("movers_metnohf", "MET over Sum ET for MET no HF collection", 50, 0, 1);
0079   MEcontainer_["moverl_metnohf"] =
0080       iBooker.book1D("moverl_metnohf", "MET over leading jet Pt for MET no HF collection", 50, 0, 2);
0081 
0082   MEcontainer_["movers_metnohfho"] =
0083       iBooker.book1D("movers_metnohfho", "MET over Sum ET for MET no HF HO collection", 50, 0, 1);
0084   MEcontainer_["moverl_metnohfho"] =
0085       iBooker.book1D("moverl_metnohfho", "MET over leading jet Pt for MET no HF HO collection", 50, 0, 2);
0086 
0087   // book histograms for EMF fraction for all jets and first 3 jets
0088   MEcontainer_["inclusive_jet_EMF"] = iBooker.book1D("inclusive_jet_EMF", "inclusive jet EMF", 50, -1, 1);
0089   MEcontainer_["leading_jet_EMF"] = iBooker.book1D("leading_jet_EMF", "leading jet EMF", 50, -1, 1);
0090   MEcontainer_["second_jet_EMF"] = iBooker.book1D("second_jet_EMF", "second jet EMF", 50, -1, 1);
0091   MEcontainer_["third_jet_EMF"] = iBooker.book1D("third_jet_EMF", "third jet EMF", 50, -1, 1);
0092 }
0093 
0094 // method to calculate MET over Sum ET from a particular MET collection
0095 float QcdHighPtDQM::movers(const CaloMETCollection &metcollection) {
0096   float metovers = 0;
0097   CaloMETCollection::const_iterator met_iter;
0098   for (met_iter = metcollection.begin(); met_iter != metcollection.end(); ++met_iter) {
0099     float mex = met_iter->momentum().x();
0100     float mey = met_iter->momentum().y();
0101     float met = sqrt(mex * mex + mey * mey);
0102     float sumet = met_iter->sumEt();
0103     metovers = (met / sumet);
0104   }
0105   return metovers;
0106 }
0107 
0108 // method to calculate MET over Leading jet PT for a particular MET collection
0109 float QcdHighPtDQM::moverl(const CaloMETCollection &metcollection, float &ljpt) {
0110   float metoverl = 0;
0111   CaloMETCollection::const_iterator met_iter;
0112   for (met_iter = metcollection.begin(); met_iter != metcollection.end(); ++met_iter) {
0113     float mex = met_iter->momentum().x();
0114     float mey = met_iter->momentum().y();
0115     float met = sqrt(mex * mex + mey * mey);
0116     metoverl = (met / ljpt);
0117   }
0118   return metoverl;
0119 }
0120 
0121 void QcdHighPtDQM::analyze(const Event &iEvent, const EventSetup &iSetup) {
0122   // Get Jets
0123   edm::Handle<CaloJetCollection> jetHandle;
0124   iEvent.getByToken(jetToken_, jetHandle);
0125   const CaloJetCollection &jets = *jetHandle;
0126   CaloJetCollection::const_iterator jet_iter;
0127 
0128   // Get MET collections
0129   edm::Handle<CaloMETCollection> metHandle;
0130   iEvent.getByToken(metToken1_, metHandle);
0131   const CaloMETCollection &met = *metHandle;
0132 
0133   edm::Handle<CaloMETCollection> metHOHandle;
0134   iEvent.getByToken(metToken2_, metHOHandle);
0135   const CaloMETCollection &metHO = *metHOHandle;
0136 
0137   edm::Handle<CaloMETCollection> metNoHFHandle;
0138   iEvent.getByToken(metToken3_, metNoHFHandle);
0139   const CaloMETCollection &metNoHF = *metNoHFHandle;
0140 
0141   edm::Handle<CaloMETCollection> metNoHFHOHandle;
0142   iEvent.getByToken(metToken4_, metNoHFHOHandle);
0143   const CaloMETCollection &metNoHFHO = *metNoHFHOHandle;
0144 
0145   // initialize leading jet value and jet multiplicity counter
0146   int njets = 0;
0147   int njets30 = 0;
0148   float leading_jetpt = 0;
0149   float leading_jeteta = 0;
0150 
0151   // initialize variables for picking out leading 2 barrel jets
0152   reco::CaloJet leadingbarreljet;
0153   reco::CaloJet secondbarreljet;
0154   int nbarreljets = 0;
0155 
0156   // get bins in eta.
0157   // Bins correspond to calotower regions.
0158 
0159   const float etabins[83] = {
0160       -5.191, -4.889, -4.716, -4.538, -4.363, -4.191, -4.013, -3.839, -3.664, -3.489, -3.314, -3.139, -2.964, -2.853,
0161       -2.650, -2.500, -2.322, -2.172, -2.043, -1.930, -1.830, -1.740, -1.653, -1.566, -1.479, -1.392, -1.305, -1.218,
0162       -1.131, -1.044, -.957,  -.879,  -.783,  -.696,  -.609,  -.522,  -.435,  -.348,  -.261,  -.174,  -.087,  0,
0163       .087,   .174,   .261,   .348,   .435,   .522,   .609,   .696,   .783,   .879,   .957,   1.044,  1.131,  1.218,
0164       1.305,  1.392,  1.479,  1.566,  1.653,  1.740,  1.830,  1.930,  2.043,  2.172,  2.322,  2.500,  2.650,  2.853,
0165       2.964,  3.139,  3.314,  3.489,  3.664,  3.839,  4.013,  4.191,  4.363,  4.538,  4.889,  5.191};
0166 
0167   for (jet_iter = jets.begin(); jet_iter != jets.end(); ++jet_iter) {
0168     njets++;
0169 
0170     // get Jet stats
0171     float jet_pt = jet_iter->pt();
0172     float jet_eta = jet_iter->eta();
0173     float jet_phi = jet_iter->phi();
0174 
0175     // fill jet Pt and jet EMF
0176     MEcontainer_["inclusive_jet_pt"]->Fill(jet_pt);
0177     MEcontainer_["inclusive_jet_EMF"]->Fill(jet_iter->emEnergyFraction());
0178 
0179     // pick out up to the first 2 leading barrel jets
0180     // for use in calculating dijet mass in barrel region
0181     // also fill jet Pt histogram for barrel
0182 
0183     if (jet_eta <= 1.3) {
0184       MEcontainer_["inclusive_jet_pt_barrel"]->Fill(jet_pt);
0185       if (nbarreljets == 0) {
0186         leadingbarreljet = jets[(njets - 1)];
0187         nbarreljets++;
0188       } else if (nbarreljets == 1) {
0189         secondbarreljet = jets[(njets - 1)];
0190         nbarreljets++;
0191       }
0192 
0193     }
0194 
0195     // fill jet Pt for endcap and forward regions
0196     else if (jet_eta <= 3.0 && jet_eta > 1.3) {
0197       MEcontainer_["inclusive_jet_pt_endcap"]->Fill(jet_pt);
0198     } else if (jet_eta <= 5.0 && jet_eta > 3.0) {
0199       MEcontainer_["inclusive_jet_pt_forward"]->Fill(jet_pt);
0200     }
0201 
0202     // count jet multiplicity for jets with Pt > 30
0203     if ((jet_pt) > 30)
0204       njets30++;
0205 
0206     // check leading jet quantities
0207     if (jet_pt > leading_jetpt) {
0208       leading_jetpt = jet_pt;
0209       leading_jeteta = jet_eta;
0210     }
0211 
0212     // fill eta-phi plot
0213     for (int eit = 0; eit < 81; eit++) {
0214       for (int pit = 0; pit < 72; pit++) {
0215         float low_eta = etabins[eit];
0216         float high_eta = etabins[eit + 1];
0217         float low_phi = (-M_PI) + pit * (M_PI / 36);
0218         float high_phi = low_phi + (M_PI / 36);
0219         if (jet_eta > low_eta && jet_eta < high_eta && jet_phi > low_phi && jet_phi < high_phi) {
0220           MEcontainer_["etaphi"]->Fill((eit - 41), jet_phi);
0221         }
0222       }
0223     }
0224   }
0225 
0226   // after iterating over all jets, fill leading jet quantity histograms
0227   // and jet multiplicity histograms
0228 
0229   MEcontainer_["leading_jet_pt"]->Fill(leading_jetpt);
0230 
0231   if (leading_jeteta <= 1.3) {
0232     MEcontainer_["leading_jet_pt_barrel"]->Fill(leading_jetpt);
0233   } else if (leading_jeteta <= 3.0 && leading_jeteta > 1.3) {
0234     MEcontainer_["leading_jet_pt_endcap"]->Fill(leading_jetpt);
0235   } else if (leading_jeteta <= 5.0 && leading_jeteta > 3.0) {
0236     MEcontainer_["leading_jet_pt_forward"]->Fill(leading_jetpt);
0237   }
0238 
0239   MEcontainer_["njets"]->Fill(njets);
0240 
0241   MEcontainer_["njets30"]->Fill(njets30);
0242 
0243   // fill MET over Sum ET and Leading jet PT for all MET flavors
0244   MEcontainer_["movers_met"]->Fill(movers(met));
0245   MEcontainer_["moverl_met"]->Fill(movers(met), leading_jetpt);
0246   MEcontainer_["movers_metho"]->Fill(movers(metHO));
0247   MEcontainer_["moverl_metho"]->Fill(movers(metHO), leading_jetpt);
0248   MEcontainer_["movers_metnohf"]->Fill(movers(metNoHF));
0249   MEcontainer_["moverl_metnohf"]->Fill(movers(metNoHF), leading_jetpt);
0250   MEcontainer_["movers_metnohfho"]->Fill(movers(metNoHFHO));
0251   MEcontainer_["moverl_metnohfho"]->Fill(movers(metNoHFHO), leading_jetpt);
0252 
0253   // fetch first 3 jet EMF
0254 
0255   if (!jets.empty()) {
0256     MEcontainer_["leading_jet_EMF"]->Fill(jets[0].emEnergyFraction());
0257     if (jets.size() >= 2) {
0258       MEcontainer_["second_jet_EMF"]->Fill(jets[1].emEnergyFraction());
0259       if (jets.size() >= 3) {
0260         MEcontainer_["third_jet_EMF"]->Fill(jets[2].emEnergyFraction());
0261       }
0262     }
0263   }
0264 
0265   // if 2 nontrivial barrel jets, reconstruct dijet mass
0266 
0267   if (nbarreljets == 2) {
0268     if (leadingbarreljet.energy() > 0 && secondbarreljet.energy() > 0) {
0269       math::XYZTLorentzVector DiJet = leadingbarreljet.p4() + secondbarreljet.p4();
0270       float dijet_mass = DiJet.mass();
0271       MEcontainer_["dijet_mass"]->Fill(dijet_mass);
0272     }
0273   }
0274 }
0275 
0276 // Local Variables:
0277 // show-trailing-whitespace: t
0278 // truncate-lines: t
0279 // End: