Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:32:42

0001 // -*- C++ -*-
0002 //
0003 // Package:    L1T
0004 // Class:      L1Validator
0005 //
0006 /**
0007  * \class L1T L1Validator.cc Validation/L1T/plugins/L1Validator.cc
0008  *
0009  * Description: [one line class summary]
0010  *
0011  * Implementation:
0012  *    [Notes on implementation]
0013  */
0014 //
0015 // Original Author:  Scott Wilbur
0016 //         Created:  Wed, 28 Aug 2013 09:42:55 GMT
0017 // $Id$
0018 //
0019 //
0020 
0021 #include <string>
0022 
0023 #include <Validation/L1T/interface/L1Validator.h>
0024 
0025 #include "FWCore/ServiceRegistry/interface/Service.h"
0026 
0027 #include "DataFormats/Math/interface/deltaR.h"
0028 
0029 #include "TFile.h"
0030 
0031 // defining as a macro instead of a function because inheritance doesn't work:
0032 #define FINDRECOPART(TYPE, COLLECTION1, COLLECTION2)                                                    \
0033   const TYPE *RecoPart = NULL;                                                                          \
0034   double BestDist = 999.;                                                                               \
0035   for (uint i = 0; i < COLLECTION1->size(); i++) {                                                      \
0036     const TYPE *ThisPart = &COLLECTION1->at(i);                                                         \
0037     double ThisDist = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());   \
0038     if (ThisDist < 1.0 && ThisDist < BestDist) {                                                        \
0039       BestDist = ThisDist;                                                                              \
0040       RecoPart = ThisPart;                                                                              \
0041     }                                                                                                   \
0042   }                                                                                                     \
0043   if (COLLECTION1.product() != COLLECTION2.product()) {                                                 \
0044     for (uint i = 0; i < COLLECTION2->size(); i++) {                                                    \
0045       const TYPE *ThisPart = &COLLECTION2->at(i);                                                       \
0046       double ThisDist = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi()); \
0047       if (ThisDist < 1.0 && ThisDist < BestDist) {                                                      \
0048         BestDist = ThisDist;                                                                            \
0049         RecoPart = ThisPart;                                                                            \
0050       }                                                                                                 \
0051     }                                                                                                   \
0052   }
0053 
0054 L1Validator::L1Validator(const edm::ParameterSet &iConfig) {
0055   _dirName = iConfig.getParameter<std::string>("dirName");
0056   _GenSource = consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("GenSource"));
0057 
0058   _L1MuonBXSource = consumes<l1t::MuonBxCollection>(iConfig.getParameter<edm::InputTag>("L1MuonBXSource"));
0059   _L1EGammaBXSource = consumes<l1t::EGammaBxCollection>(iConfig.getParameter<edm::InputTag>("L1EGammaBXSource"));
0060   _L1TauBXSource = consumes<l1t::TauBxCollection>(iConfig.getParameter<edm::InputTag>("L1TauBXSource"));
0061   _L1JetBXSource = consumes<l1t::JetBxCollection>(iConfig.getParameter<edm::InputTag>("L1JetBXSource"));
0062   _srcToken = mayConsume<GenEventInfoProduct>(iConfig.getParameter<edm::InputTag>("srcToken"));
0063   _L1GenJetSource = consumes<reco::GenJetCollection>(iConfig.getParameter<edm::InputTag>("L1GenJetSource"));
0064 
0065   //_fileName = iConfig.getParameter<std::string>("fileName");
0066 }
0067 
0068 L1Validator::~L1Validator() {}
0069 
0070 void L1Validator::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &, edm::EventSetup const &) {
0071   // iBooker.setCurrentFolder(_dirName);
0072   _Hists.Book(iBooker, _dirName);
0073 };
0074 
0075 void L1Validator::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0076   using namespace edm;
0077   using namespace std;
0078   using namespace l1extra;
0079   using namespace reco;
0080 
0081   Handle<GenParticleCollection> GenParticles;
0082   Handle<l1t::MuonBxCollection> MuonsBX;
0083   Handle<l1t::EGammaBxCollection> EGammasBX;
0084   Handle<l1t::TauBxCollection> TausBX;
0085   Handle<l1t::JetBxCollection> JetsBX;
0086   Handle<GenEventInfoProduct> genEvtInfoProduct;
0087   Handle<reco::GenJetCollection> GenJets;
0088 
0089   bool GotEverything = true;
0090 
0091   if (!iEvent.getByToken(_GenSource, GenParticles))
0092     GotEverything = false;
0093   if (!iEvent.getByToken(_L1MuonBXSource, MuonsBX))
0094     GotEverything = false;
0095   if (!iEvent.getByToken(_L1EGammaBXSource, EGammasBX))
0096     GotEverything = false;
0097   if (!iEvent.getByToken(_L1TauBXSource, TausBX))
0098     GotEverything = false;
0099   if (!iEvent.getByToken(_L1JetBXSource, JetsBX))
0100     GotEverything = false;
0101   if (!iEvent.getByToken(_srcToken, genEvtInfoProduct))
0102     GotEverything = false;
0103   if (!iEvent.getByToken(_L1GenJetSource, GenJets))
0104     GotEverything = false;
0105 
0106   if (!GotEverything)
0107     return;
0108 
0109   /*
0110   std::string moduleName = "";
0111   if( genEvtInfoProduct.isValid() ) {
0112           const edm::Provenance& prov =
0113   iEvent.getProvenance(genEvtInfoProduct.id()); moduleName =
0114   edm::moduleName(prov);
0115           //cout<<" generator name: "<<moduleName<<endl;
0116   }
0117   */
0118 
0119   _Hists.NEvents++;
0120 
0121   int nL1Muons = 0, nL1EGammas = 0, nL1Taus = 0, nL1Jets = 0;
0122   if (MuonsBX->getFirstBX() >= 0)
0123     nL1Muons = MuonsBX->size(0);
0124   if (EGammasBX->getFirstBX() >= 0)
0125     nL1EGammas = EGammasBX->size(0);
0126   if (TausBX->getFirstBX() >= 0)
0127     nL1Taus = TausBX->size(0);
0128   if (JetsBX->getFirstBX() >= 0)
0129     nL1Jets = JetsBX->size(0);
0130 
0131   _Hists.FillNumber(L1ValidatorHists::Type::Muon, nL1Muons);
0132   _Hists.FillNumber(L1ValidatorHists::Type::Egamma, nL1EGammas);
0133   _Hists.FillNumber(L1ValidatorHists::Type::Tau, nL1Taus);
0134   _Hists.FillNumber(L1ValidatorHists::Type::Jet, nL1Jets);
0135 
0136   // For gen jet
0137 
0138   for (auto &Genjet : *GenJets) {
0139     // eta within calorimeter acceptance 4.7
0140     if (fabs((&Genjet)->eta()) > 4.7)
0141       continue;
0142 
0143     // only consider the gen jet with pt greater than 10 GeV
0144     if ((&Genjet)->pt() < 10.0)
0145       continue;
0146 
0147     double minDR = 999.0;
0148 
0149     // match L1T object
0150     const l1t::Jet *L1Part = nullptr;
0151     for (int iBx = JetsBX->getFirstBX(); iBx <= JetsBX->getLastBX(); ++iBx) {
0152       if (iBx > 0)
0153         continue;
0154       for (std::vector<l1t::Jet>::const_iterator jet = JetsBX->begin(iBx); jet != JetsBX->end(iBx); ++jet) {
0155         double idR = reco::deltaR((&Genjet)->eta(), (&Genjet)->phi(), jet->eta(), jet->phi());
0156         if (idR < minDR) {
0157           minDR = idR;
0158           L1Part = &(*jet);
0159         }
0160       }
0161     }
0162     _Hists.Fill(L1ValidatorHists::Type::Jet, &Genjet, L1Part);
0163   }
0164 
0165   for (uint i = 0; i < GenParticles->size(); i++) {
0166     const GenParticle *GenPart = &GenParticles->at(i);
0167 
0168     int pdg = GenPart->pdgId(), status = GenPart->status();
0169 
0170     double minDR = 999.0;
0171 
0172     // only consider the gen particle with pt greater than 10 GeV
0173     if (GenPart->pt() < 10.0)
0174       continue;
0175 
0176     /// select the final state (i.e status==1) muons (pdg==+/-13)
0177     if (status == 1 && abs(pdg) == 13) {  // Muon
0178 
0179       // eta within tracker acceptance 2.4
0180       if (fabs(GenPart->eta()) > 2.4)
0181         continue;
0182 
0183       // match L1T object
0184       const l1t::Muon *L1Part = nullptr;
0185       for (int iBx = MuonsBX->getFirstBX(); iBx <= MuonsBX->getLastBX(); ++iBx) {
0186         if (iBx > 0)
0187           continue;
0188         for (std::vector<l1t::Muon>::const_iterator mu = MuonsBX->begin(iBx); mu != MuonsBX->end(iBx); ++mu) {
0189           double idR = reco::deltaR(GenPart->eta(), GenPart->phi(), mu->eta(), mu->phi());
0190           if (idR < minDR) {
0191             minDR = idR;
0192             L1Part = &(*mu);
0193           }
0194         }
0195         _Hists.Fill(L1ValidatorHists::Type::Muon, GenPart, L1Part);
0196       }
0197 
0198       /// select the final state (i.e status==1) electrons (pdg==+/-11) and
0199       /// photons (pdg==22)
0200     } else if (status == 1 && (abs(pdg) == 11 || pdg == 22)) {  // Egamma
0201 
0202       // eta within EM calorimeter acceptance 2.5
0203       if (fabs(GenPart->eta()) > 2.5)
0204         continue;
0205 
0206       // exclude the calorimeter barrel and endcap overlap region
0207       if (fabs(GenPart->eta()) > 1.4442 && fabs(GenPart->eta()) < 1.5660)
0208         continue;
0209 
0210       // match L1T object
0211       const l1t::EGamma *L1Part = nullptr;
0212       for (int iBx = EGammasBX->getFirstBX(); iBx <= EGammasBX->getLastBX(); ++iBx) {
0213         if (iBx > 0)
0214           continue;
0215         for (std::vector<l1t::EGamma>::const_iterator eg = EGammasBX->begin(iBx); eg != EGammasBX->end(iBx); ++eg) {
0216           double idR = reco::deltaR(GenPart->eta(), GenPart->phi(), eg->eta(), eg->phi());
0217           if (idR < minDR) {
0218             minDR = idR;
0219             L1Part = &(*eg);
0220           }
0221         }
0222       }
0223       _Hists.Fill(L1ValidatorHists::Type::Egamma, GenPart, L1Part);
0224 
0225       /// select the matrix element (i.e status==2) taus (pdg==+/-15) before
0226       /// decay
0227     } else if (status == 2 && abs(pdg) == 15) {  // Tau
0228 
0229       // eta within tracker acceptance 2.4
0230       if (fabs(GenPart->eta()) > 2.4)
0231         continue;
0232 
0233       // match L1T object
0234       const l1t::Tau *L1Part = nullptr;
0235       for (int iBx = TausBX->getFirstBX(); iBx <= TausBX->getLastBX(); ++iBx) {
0236         if (iBx > 0)
0237           continue;
0238         for (std::vector<l1t::Tau>::const_iterator tau = TausBX->begin(iBx); tau != TausBX->end(iBx); ++tau) {
0239           double idR = reco::deltaR(GenPart->eta(), GenPart->phi(), tau->eta(), tau->phi());
0240           if (idR < minDR) {
0241             minDR = idR;
0242             L1Part = &(*tau);
0243           }
0244         }
0245       }
0246       _Hists.Fill(L1ValidatorHists::Type::Tau, GenPart, L1Part);
0247     }
0248   }
0249 }
0250 
0251 // The next three are exactly the same, but apparently inheritance doesn't work
0252 // like I thought it did.
0253 const reco::LeafCandidate *L1Validator::FindBest(const reco::GenParticle *GenPart,
0254                                                  const std::vector<l1extra::L1EmParticle> *Collection1,
0255                                                  const std::vector<l1extra::L1EmParticle> *Collection2 = nullptr) {
0256   const reco::LeafCandidate *BestPart = nullptr;
0257   double BestDR = 999.;
0258 
0259   for (uint i = 0; i < Collection1->size(); i++) {
0260     const reco::LeafCandidate *ThisPart = &Collection1->at(i);
0261     double ThisDR = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());
0262     if (ThisDR < BestDR) {
0263       BestDR = ThisDR;
0264       BestPart = ThisPart;
0265     }
0266   }
0267 
0268   if (Collection2 == nullptr)
0269     return BestPart;
0270 
0271   for (uint i = 0; i < Collection2->size(); i++) {
0272     const reco::LeafCandidate *ThisPart = &Collection2->at(i);
0273     double ThisDR = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());
0274     if (ThisDR < BestDR) {
0275       BestDR = ThisDR;
0276       BestPart = ThisPart;
0277     }
0278   }
0279 
0280   return BestPart;
0281 }
0282 
0283 const reco::LeafCandidate *L1Validator::FindBest(const reco::GenParticle *GenPart,
0284                                                  const std::vector<l1extra::L1JetParticle> *Collection1,
0285                                                  const std::vector<l1extra::L1JetParticle> *Collection2 = nullptr) {
0286   const reco::LeafCandidate *BestPart = nullptr;
0287   double BestDR = 999.;
0288 
0289   for (uint i = 0; i < Collection1->size(); i++) {
0290     const reco::LeafCandidate *ThisPart = &Collection1->at(i);
0291     double ThisDR = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());
0292     if (ThisDR < BestDR) {
0293       BestDR = ThisDR;
0294       BestPart = ThisPart;
0295     }
0296   }
0297 
0298   if (Collection2 == nullptr)
0299     return BestPart;
0300 
0301   for (uint i = 0; i < Collection2->size(); i++) {
0302     const reco::LeafCandidate *ThisPart = &Collection2->at(i);
0303     double ThisDR = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());
0304     if (ThisDR < BestDR) {
0305       BestDR = ThisDR;
0306       BestPart = ThisPart;
0307     }
0308   }
0309 
0310   return BestPart;
0311 }
0312 
0313 const reco::LeafCandidate *L1Validator::FindBest(const reco::GenParticle *GenPart,
0314                                                  const std::vector<l1extra::L1MuonParticle> *Collection1) {
0315   const reco::LeafCandidate *BestPart = nullptr;
0316   double BestDR = 999.;
0317 
0318   for (uint i = 0; i < Collection1->size(); i++) {
0319     const reco::LeafCandidate *ThisPart = &Collection1->at(i);
0320     double ThisDR = reco::deltaR(GenPart->eta(), GenPart->phi(), ThisPart->eta(), ThisPart->phi());
0321     if (ThisDR < BestDR) {
0322       BestDR = ThisDR;
0323       BestPart = ThisPart;
0324     }
0325   }
0326 
0327   return BestPart;
0328 }
0329 
0330 // ------------ method fills 'descriptions' with the allowed parameters for the
0331 // module  ------------
0332 void L1Validator::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0333   // The following says we do not know what parameters are allowed so do no
0334   // validation
0335   // Please change this to state exactly what you do use, even if it is no
0336   // parameters
0337   edm::ParameterSetDescription desc;
0338   desc.setUnknown();
0339   descriptions.addDefault(desc);
0340 }
0341 
0342 // define this as a plug-in
0343 DEFINE_FWK_MODULE(L1Validator);