EgammaIsolationAlgos/src/ElectronTkIsolation.cc

0001 //*****************************************************************************
0002 // File:      ElectronTkIsolation.cc
0003 // ----------------------------------------------------------------------------
0004 // OrigAuth:  Matthias Mozer
0005 // Institute: IIHE-VUB
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 // unified acces to isolations
0046 std::pair<int, double> ElectronTkIsolation::getIso(const reco::Track* tmpTrack) const {
0047   int counter = 0;
0048   double ptSum = 0.;
0049   //Take the electron track
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   }  //end loop over tracks
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   //counter for the tracks in the isolation cone
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   //check also the originalAlgo in case the track is reconstructed by more than one
0111   algo = trk.originalAlgo();
0112   //reject only if both algo and originalAlgo are not acceptable
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 }