File indexing completed on 2023-03-17 11:16:36
0001 #include "PhysicsTools/PatUtils/interface/LeptonJetIsolationAngle.h"
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "FWCore/Utilities/interface/Exception.h"
0005 #include "DataFormats/Math/interface/deltaR.h"
0006
0007 #include <vector>
0008
0009 using namespace pat;
0010
0011
0012 LeptonJetIsolationAngle::LeptonJetIsolationAngle(edm::ConsumesCollector&& iC)
0013 : jetToken_(iC.consumes<reco::CaloJetCollection>(edm::InputTag("iterativeCone5CaloJets"))),
0014 electronsToken_(iC.consumes<std::vector<reco::GsfElectron> >(edm::InputTag("pixelMatchGsfElectrons"))) {}
0015
0016
0017 LeptonJetIsolationAngle::~LeptonJetIsolationAngle() {}
0018
0019
0020 float LeptonJetIsolationAngle::calculate(const Electron& theElectron,
0021 const edm::Handle<edm::View<reco::Track> >& trackHandle,
0022 const edm::Event& iEvent) {
0023 CLHEP::HepLorentzVector theElectronHLV(theElectron.px(), theElectron.py(), theElectron.pz(), theElectron.energy());
0024 return this->calculate(theElectronHLV, trackHandle, iEvent);
0025 }
0026 float LeptonJetIsolationAngle::calculate(const Muon& theMuon,
0027 const edm::Handle<edm::View<reco::Track> >& trackHandle,
0028 const edm::Event& iEvent) {
0029 CLHEP::HepLorentzVector theMuonHLV(theMuon.px(), theMuon.py(), theMuon.pz(), theMuon.energy());
0030 return this->calculate(theMuonHLV, trackHandle, iEvent);
0031 }
0032
0033
0034 float LeptonJetIsolationAngle::calculate(const CLHEP::HepLorentzVector& aLepton,
0035 const edm::Handle<edm::View<reco::Track> >& trackHandle,
0036 const edm::Event& iEvent) {
0037
0038
0039 edm::Handle<reco::CaloJetCollection> jetHandle;
0040 iEvent.getByToken(jetToken_, jetHandle);
0041 reco::CaloJetCollection jetColl = *(jetHandle.product());
0042
0043 edm::Handle<std::vector<reco::GsfElectron> > electronsHandle;
0044 iEvent.getByToken(electronsToken_, electronsHandle);
0045 std::vector<reco::GsfElectron> electrons = *electronsHandle;
0046
0047 std::vector<Electron> isoElectrons;
0048 for (size_t ie = 0; ie < electrons.size(); ie++) {
0049 Electron anElectron(electrons[ie]);
0050 if (anElectron.pt() > 10 && trkIsolator_.calculate(anElectron, *trackHandle) < 3.0) {
0051 isoElectrons.push_back(electrons[ie]);
0052 }
0053 }
0054
0055 std::vector<reco::CaloJet> theJets;
0056 for (reco::CaloJetCollection::const_iterator itJet = jetColl.begin(); itJet != jetColl.end(); itJet++) {
0057 float mindr2 = 9999.;
0058 for (size_t ie = 0; ie < isoElectrons.size(); ie++) {
0059 float dr2 = ::deltaR2(*itJet, isoElectrons[ie]);
0060 if (dr2 < mindr2)
0061 mindr2 = dr2;
0062 }
0063 float mindr = std::sqrt(mindr2);
0064
0065 if (itJet->et() > 15 && mindr > 0.3)
0066 theJets.push_back(reco::CaloJet(*itJet));
0067 }
0068
0069 float isoAngle = 1000;
0070 for (std::vector<reco::CaloJet>::const_iterator itJet = theJets.begin(); itJet != theJets.end(); itJet++) {
0071 float curDR = this->spaceAngle(aLepton, *itJet);
0072 if (curDR < isoAngle)
0073 isoAngle = curDR;
0074 }
0075 return isoAngle;
0076 }
0077
0078
0079 float LeptonJetIsolationAngle::spaceAngle(const CLHEP::HepLorentzVector& aLepton, const reco::CaloJet& aJet) {
0080 return acos(sin(aJet.theta()) * cos(aJet.phi()) * sin(aLepton.theta()) * cos(aLepton.phi()) +
0081 sin(aJet.theta()) * sin(aJet.phi()) * sin(aLepton.theta()) * sin(aLepton.phi()) +
0082 cos(aJet.theta()) * cos(aLepton.theta()));
0083 }