Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:21:24

0001 ////////////////////////////////////////////////////////////////////////////////
0002 //
0003 // L6SLBCorrector
0004 // --------------
0005 //
0006 //           25/10/2009   Hauke Held             <hauke.held@cern.ch>
0007 //                        Philipp Schieferdecker <philipp.schieferdecker@cern.ch
0008 ////////////////////////////////////////////////////////////////////////////////
0009 
0010 #include "JetMETCorrections/Algorithms/interface/L6SLBCorrector.h"
0011 
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/ParameterSet/interface/FileInPath.h"
0015 #include "FWCore/Framework/interface/Event.h"
0016 
0017 #include "DataFormats/TrackReco/interface/Track.h"
0018 
0019 #include <string>
0020 
0021 using namespace std;
0022 
0023 ////////////////////////////////////////////////////////////////////////////////
0024 // construction / destruction
0025 ////////////////////////////////////////////////////////////////////////////////
0026 
0027 //______________________________________________________________________________
0028 L6SLBCorrector::L6SLBCorrector(const JetCorrectorParameters& fParam, const edm::ParameterSet& fConfig)
0029     : addMuonToJet_(fConfig.getParameter<bool>("addMuonToJet")),
0030       srcBTagInfoElec_(fConfig.getParameter<edm::InputTag>("srcBTagInfoElectron")),
0031       srcBTagInfoMuon_(fConfig.getParameter<edm::InputTag>("srcBTagInfoMuon")),
0032       corrector_(nullptr) {
0033   vector<JetCorrectorParameters> vParam;
0034   vParam.push_back(fParam);
0035   corrector_ = new FactorizedJetCorrectorCalculator(vParam);
0036 }
0037 
0038 //______________________________________________________________________________
0039 L6SLBCorrector::~L6SLBCorrector() { delete corrector_; }
0040 
0041 ////////////////////////////////////////////////////////////////////////////////
0042 // implementation of member functions
0043 ////////////////////////////////////////////////////////////////////////////////
0044 
0045 //______________________________________________________________________________
0046 double L6SLBCorrector::correction(const LorentzVector& fJet) const {
0047   throw cms::Exception("EventRequired") << "Wrong interface correction(LorentzVector), event required!";
0048   return 1.0;
0049 }
0050 
0051 //______________________________________________________________________________
0052 double L6SLBCorrector::correction(const reco::Jet& fJet) const {
0053   throw cms::Exception("EventRequired") << "Wrong interface correction(reco::Jet), event required!";
0054   return 1.0;
0055 }
0056 
0057 //______________________________________________________________________________
0058 double L6SLBCorrector::correction(const reco::Jet& fJet,
0059                                   const edm::RefToBase<reco::Jet>& refToRawJet,
0060                                   const edm::Event& fEvent,
0061                                   const edm::EventSetup& fSetup) const {
0062   FactorizedJetCorrectorCalculator::VariableValues values;
0063   values.setJetPt(fJet.pt());
0064   values.setJetEta(fJet.eta());
0065   values.setJetPhi(fJet.phi());
0066   values.setJetE(fJet.energy());
0067 
0068   edm::Handle<vector<reco::SoftLeptonTagInfo> > muoninfos;
0069   fEvent.getByLabel(srcBTagInfoMuon_, muoninfos);
0070 
0071   const reco::SoftLeptonTagInfo& sltMuon = (*muoninfos)[getBTagInfoIndex(refToRawJet, *muoninfos)];
0072   if (sltMuon.leptons() > 0) {
0073     const edm::RefToBase<reco::Track>& trackRef = sltMuon.lepton(0);
0074     values.setLepPx(trackRef->px());
0075     values.setLepPy(trackRef->py());
0076     values.setLepPz(trackRef->pz());
0077     values.setAddLepToJet(addMuonToJet_);
0078     return corrector_->getCorrection(values);
0079   } else {
0080     edm::Handle<vector<reco::SoftLeptonTagInfo> > elecinfos;
0081     fEvent.getByLabel(srcBTagInfoElec_, elecinfos);
0082     const reco::SoftLeptonTagInfo& sltElec = (*elecinfos)[getBTagInfoIndex(refToRawJet, *elecinfos)];
0083     if (sltElec.leptons() > 0) {
0084       const edm::RefToBase<reco::Track>& trackRef = sltElec.lepton(0);
0085       values.setLepPx(trackRef->px());
0086       values.setLepPy(trackRef->py());
0087       values.setLepPz(trackRef->pz());
0088       values.setAddLepToJet(false);
0089       return corrector_->getCorrection(values);
0090     }
0091   }
0092   return 1.0;
0093 }
0094 
0095 ////////////////////////////////////////////////////////////////////////////////
0096 // implementation of private member functions
0097 ////////////////////////////////////////////////////////////////////////////////
0098 
0099 //______________________________________________________________________________
0100 int L6SLBCorrector::getBTagInfoIndex(const edm::RefToBase<reco::Jet>& refToRawJet,
0101                                      const vector<reco::SoftLeptonTagInfo>& tags) const {
0102   for (unsigned int i = 0; i < tags.size(); i++)
0103     if (tags[i].jet().get() == refToRawJet.get())
0104       return i;
0105   return -1;
0106 }