File indexing completed on 2023-10-25 09:59:06
0001 #include "DataFormats/EgammaCandidates/interface/GsfElectronCore.h"
0002 #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrack.h"
0003 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0004 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0005 #include "RecoEgamma/EgammaElectronAlgos/interface/GsfElectronTools.h"
0006 #include "DataFormats/Math/interface/deltaR.h"
0007
0008 namespace egamma {
0009
0010 using namespace reco;
0011
0012
0013
0014
0015
0016 std::pair<TrackRef, float> getClosestCtfToGsf(GsfTrackRef const& gsfTrackRef,
0017 edm::Handle<reco::TrackCollection> const& ctfTracksH,
0018 edm::soa::EtaPhiTableView trackTable) {
0019 float maxFracShared = 0;
0020 TrackRef ctfTrackRef = TrackRef();
0021 const TrackCollection* ctfTrackCollection = ctfTracksH.product();
0022
0023 float gsfEta = gsfTrackRef->eta();
0024 float gsfPhi = gsfTrackRef->phi();
0025 const HitPattern& gsfHitPattern = gsfTrackRef->hitPattern();
0026
0027 constexpr float dR2 = 0.3 * 0.3;
0028
0029 unsigned int counter = 0;
0030 for (auto ctfTkIter = ctfTrackCollection->begin(); ctfTkIter != ctfTrackCollection->end(); ctfTkIter++, counter++) {
0031
0032 using namespace edm::soa::col;
0033 if (reco::deltaR2(gsfEta, gsfPhi, trackTable.get<Eta>(counter), trackTable.get<Phi>(counter)) > dR2)
0034 continue;
0035
0036 unsigned int shared = 0;
0037 int gsfHitCounter = 0;
0038 int numGsfInnerHits = 0;
0039 int numCtfInnerHits = 0;
0040
0041 const HitPattern& ctfHitPattern = ctfTkIter->hitPattern();
0042
0043 for (auto elHitsIt = gsfTrackRef->recHitsBegin(); elHitsIt != gsfTrackRef->recHitsEnd();
0044 elHitsIt++, gsfHitCounter++) {
0045 if (!((**elHitsIt).isValid()))
0046 {
0047 continue;
0048 }
0049
0050
0051 uint32_t gsfHit = gsfHitPattern.getHitPattern(HitPattern::TRACK_HITS, gsfHitCounter);
0052 if (!(HitPattern::pixelHitFilter(gsfHit) || HitPattern::stripTIBHitFilter(gsfHit) ||
0053 HitPattern::stripTIDHitFilter(gsfHit))) {
0054 continue;
0055 }
0056
0057 numGsfInnerHits++;
0058
0059 int ctfHitsCounter = 0;
0060 numCtfInnerHits = 0;
0061 for (auto ctfHitsIt = ctfTkIter->recHitsBegin(); ctfHitsIt != ctfTkIter->recHitsEnd();
0062 ctfHitsIt++, ctfHitsCounter++) {
0063 if (!((**ctfHitsIt).isValid()))
0064 {
0065 continue;
0066 }
0067
0068 uint32_t ctfHit = ctfHitPattern.getHitPattern(HitPattern::TRACK_HITS, ctfHitsCounter);
0069 if (!(HitPattern::pixelHitFilter(ctfHit) || HitPattern::stripTIBHitFilter(ctfHit) ||
0070 HitPattern::stripTIDHitFilter(ctfHit))) {
0071 continue;
0072 }
0073
0074 numCtfInnerHits++;
0075
0076 if ((**elHitsIt).sharesInput(&(**ctfHitsIt), TrackingRecHit::all)) {
0077 shared++;
0078 break;
0079 }
0080
0081 }
0082
0083 }
0084
0085 if ((numGsfInnerHits == 0) || (numCtfInnerHits == 0)) {
0086 continue;
0087 }
0088
0089 if (static_cast<float>(shared) / std::min(numGsfInnerHits, numCtfInnerHits) > maxFracShared) {
0090 maxFracShared = static_cast<float>(shared) / std::min(numGsfInnerHits, numCtfInnerHits);
0091 ctfTrackRef = TrackRef(ctfTracksH, counter);
0092 }
0093
0094 }
0095
0096 return {ctfTrackRef, maxFracShared};
0097 }
0098
0099 }