File indexing completed on 2023-03-17 11:17:43
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "RecoEgamma/EgammaIsolationAlgos/interface/ElectronTkIsolation.h"
0010 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0011
0012 #include <Math/VectorUtil.h>
0013
0014 ElectronTkIsolation::ElectronTkIsolation(double extRadius,
0015 double intRadiusBarrel,
0016 double intRadiusEndcap,
0017 double stripBarrel,
0018 double stripEndcap,
0019 double ptLow,
0020 double lip,
0021 double drb,
0022 const reco::TrackCollection* trackCollection,
0023 reco::TrackBase::Point beamPoint,
0024 const std::string& dzOptionString)
0025 : extRadius_(extRadius),
0026 intRadiusBarrel_(intRadiusBarrel),
0027 intRadiusEndcap_(intRadiusEndcap),
0028 stripBarrel_(stripBarrel),
0029 stripEndcap_(stripEndcap),
0030 ptLow_(ptLow),
0031 lip_(lip),
0032 drb_(drb),
0033 trackCollection_(trackCollection),
0034 beamPoint_(beamPoint) {
0035 setAlgosToReject();
0036 setDzOption(dzOptionString);
0037 }
0038
0039 ElectronTkIsolation::~ElectronTkIsolation() {}
0040
0041 std::pair<int, double> ElectronTkIsolation::getIso(const reco::GsfElectron* electron) const {
0042 return getIso(&(*(electron->gsfTrack())));
0043 }
0044
0045
0046 std::pair<int, double> ElectronTkIsolation::getIso(const reco::Track* tmpTrack) const {
0047 int counter = 0;
0048 double ptSum = 0.;
0049
0050 math::XYZVector tmpElectronMomentumAtVtx = (*tmpTrack).momentum();
0051 double tmpElectronEtaAtVertex = (*tmpTrack).eta();
0052
0053 for (reco::TrackCollection::const_iterator itrTr = (*trackCollection_).begin(); itrTr != (*trackCollection_).end();
0054 ++itrTr) {
0055 double this_pt = (*itrTr).pt();
0056 if (this_pt < ptLow_)
0057 continue;
0058
0059 double dzCut = 0;
0060 switch (dzOption_) {
0061 case egammaisolation::EgammaTrackSelector::dz:
0062 dzCut = fabs((*itrTr).dz() - (*tmpTrack).dz());
0063 break;
0064 case egammaisolation::EgammaTrackSelector::vz:
0065 dzCut = fabs((*itrTr).vz() - (*tmpTrack).vz());
0066 break;
0067 case egammaisolation::EgammaTrackSelector::bs:
0068 dzCut = fabs((*itrTr).dz(beamPoint_) - (*tmpTrack).dz(beamPoint_));
0069 break;
0070 case egammaisolation::EgammaTrackSelector::vtx:
0071 dzCut = fabs((*itrTr).dz(tmpTrack->vertex()));
0072 break;
0073 default:
0074 dzCut = fabs((*itrTr).vz() - (*tmpTrack).vz());
0075 break;
0076 }
0077 if (dzCut > lip_)
0078 continue;
0079 if (fabs((*itrTr).dxy(beamPoint_)) > drb_)
0080 continue;
0081 double dr = ROOT::Math::VectorUtil::DeltaR(itrTr->momentum(), tmpElectronMomentumAtVtx);
0082 double deta = (*itrTr).eta() - tmpElectronEtaAtVertex;
0083 bool isBarrel = std::abs(tmpElectronEtaAtVertex) < 1.479;
0084 double intRadius = isBarrel ? intRadiusBarrel_ : intRadiusEndcap_;
0085 double strip = isBarrel ? stripBarrel_ : stripEndcap_;
0086 if (dr < extRadius_ && dr >= intRadius && std::abs(deta) >= strip && passAlgo(*itrTr)) {
0087 ++counter;
0088 ptSum += this_pt;
0089 }
0090
0091 }
0092
0093 std::pair<int, double> retval;
0094 retval.first = counter;
0095 retval.second = ptSum;
0096
0097 return retval;
0098 }
0099
0100 int ElectronTkIsolation::getNumberTracks(const reco::GsfElectron* electron) const {
0101
0102 return getIso(electron).first;
0103 }
0104
0105 double ElectronTkIsolation::getPtTracks(const reco::GsfElectron* electron) const { return getIso(electron).second; }
0106
0107 bool ElectronTkIsolation::passAlgo(const reco::TrackBase& trk) const {
0108 int algo = trk.algo();
0109 bool rejAlgo = std::binary_search(algosToReject_.begin(), algosToReject_.end(), algo);
0110
0111 algo = trk.originalAlgo();
0112
0113 rejAlgo &= std::binary_search(algosToReject_.begin(), algosToReject_.end(), algo);
0114 return !rejAlgo;
0115 }
0116
0117 void ElectronTkIsolation::setAlgosToReject() {
0118 algosToReject_ = {reco::TrackBase::jetCoreRegionalStep};
0119 std::sort(algosToReject_.begin(), algosToReject_.end());
0120 }