File indexing completed on 2024-04-06 12:27:37
0001 #include "RecoParticleFlow/PFTracking/interface/PFCheckHitPattern.h"
0002
0003
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/Utilities/interface/Exception.h"
0006 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0007 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0008 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0009
0010
0011 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0012 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0013 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0014
0015 #include <map>
0016
0017 using namespace reco;
0018 using namespace std;
0019
0020 void PFCheckHitPattern::init(const TrackerTopology* tkerTopo, const TrackerGeometry* tkerGeom) {
0021
0022
0023
0024
0025 geomInitDone_ = true;
0026
0027
0028 const TrackingGeometry::DetContainer& dets = tkerGeom->dets();
0029
0030
0031 for (unsigned int i = 0; i < dets.size(); i++) {
0032
0033 auto detId = dets[i]->geographicalId();
0034 auto detInfo = DetInfo(detId.subdetId(), tkerTopo->layer(detId));
0035 uint32_t subDet = detInfo.first;
0036
0037
0038 double r_or_z;
0039 if (this->barrel(subDet)) {
0040 r_or_z = dets[i]->position().perp();
0041 } else {
0042 r_or_z = fabs(dets[i]->position().z());
0043 }
0044
0045
0046 double minRZ = 999.;
0047 double maxRZ = 0.;
0048 if (rangeRorZ_.find(detInfo) != rangeRorZ_.end()) {
0049 minRZ = rangeRorZ_[detInfo].first;
0050 maxRZ = rangeRorZ_[detInfo].second;
0051 }
0052
0053
0054 if (minRZ > r_or_z)
0055 minRZ = r_or_z;
0056 if (maxRZ < r_or_z)
0057 maxRZ = r_or_z;
0058 rangeRorZ_[detInfo] = pair<double, double>(minRZ, maxRZ);
0059 }
0060
0061 #if 0
0062
0063 RZrangeMap::const_iterator d;
0064 for (d = rangeRorZ_.begin(); d != rangeRorZ_.end(); d++) {
0065 DetInfo detInfo = d->first;
0066 pair<double, double> rangeRZ = d->second;
0067 }
0068 #endif
0069 }
0070
0071 bool PFCheckHitPattern::barrel(uint32_t subDet) {
0072
0073 return (subDet == StripSubdetector::TIB || subDet == StripSubdetector::TOB ||
0074 subDet == PixelSubdetector::PixelBarrel);
0075 }
0076
0077 pair<PFCheckHitPattern::PFTrackHitInfo, PFCheckHitPattern::PFTrackHitInfo> PFCheckHitPattern::analyze(
0078 const TrackerTopology* tkerTopo,
0079 const TrackerGeometry* tkerGeom,
0080 const TrackBaseRef track,
0081 const TransientVertex& vert) {
0082
0083
0084
0085
0086
0087
0088 if (!geomInitDone_)
0089 this->init(tkerTopo, tkerGeom);
0090
0091
0092 const reco::HitPattern& hp = track.get()->hitPattern();
0093
0094
0095
0096 unsigned int nHitBefore = 0;
0097 unsigned int nHitAfter = 0;
0098
0099 for (int i = 0; i < hp.numberOfAllHits(HitPattern::TRACK_HITS); i++) {
0100 uint32_t hit = hp.getHitPattern(HitPattern::TRACK_HITS, i);
0101 if (hp.trackerHitFilter(hit) && hp.validHitFilter(hit)) {
0102 uint32_t subDet = hp.getSubStructure(hit);
0103 uint32_t layer = hp.getLayer(hit);
0104 DetInfo detInfo(subDet, layer);
0105 double maxRZ = rangeRorZ_[detInfo].second;
0106
0107 if (this->barrel(subDet)) {
0108 if (vert.position().perp() > maxRZ)
0109 nHitBefore++;
0110 else
0111 nHitAfter++;
0112 } else {
0113 if (fabs(vert.position().z()) > maxRZ)
0114 nHitBefore++;
0115 else
0116 nHitAfter++;
0117 }
0118 }
0119 }
0120
0121
0122
0123 unsigned int nMissHitAfter = 0;
0124 unsigned int nMissHitBefore = 0;
0125
0126 for (int i = 0; i < hp.numberOfAllHits(HitPattern::MISSING_INNER_HITS); i++) {
0127 uint32_t hit = hp.getHitPattern(HitPattern::MISSING_INNER_HITS, i);
0128 if (reco::HitPattern::trackerHitFilter(hit) && reco::HitPattern::missingHitFilter(hit)) {
0129 uint32_t subDet = reco::HitPattern::getSubStructure(hit);
0130 uint32_t layer = reco::HitPattern::getLayer(hit);
0131 DetInfo detInfo(subDet, layer);
0132 double minRZ = rangeRorZ_[detInfo].first;
0133
0134
0135
0136 if (this->barrel(subDet)) {
0137 if (vert.position().perp() < minRZ)
0138 nMissHitAfter++;
0139 else
0140 nMissHitBefore++;
0141 } else {
0142 if (fabs(vert.position().z()) < minRZ)
0143 nMissHitAfter++;
0144 else
0145 nMissHitBefore++;
0146 }
0147 }
0148 }
0149
0150 PFTrackHitInfo trackToVertex(nHitBefore, nMissHitBefore);
0151 PFTrackHitInfo trackFromVertex(nHitAfter, nMissHitAfter);
0152
0153 return pair<PFTrackHitInfo, PFTrackHitInfo>(trackToVertex, trackFromVertex);
0154 }
0155
0156 void PFCheckHitPattern::print(const TrackBaseRef track) const {
0157
0158 const reco::HitPattern& hp = track.get()->hitPattern();
0159
0160 cout << "=== Hits on Track ===" << endl;
0161 this->print(reco::HitPattern::TRACK_HITS, hp);
0162 cout << "=== Hits before track ===" << endl;
0163 this->print(reco::HitPattern::MISSING_INNER_HITS, hp);
0164 }
0165
0166 void PFCheckHitPattern::print(const reco::HitPattern::HitCategory category, const reco::HitPattern& hp) const {
0167 for (int i = 0; i < hp.numberOfAllHits(category); i++) {
0168 uint32_t hit = hp.getHitPattern(category, i);
0169 if (hp.trackerHitFilter(hit)) {
0170 uint32_t subdet = hp.getSubStructure(hit);
0171 uint32_t layer = hp.getLayer(hit);
0172 cout << "hit " << i << " subdet=" << subdet << " layer=" << layer << " type " << hp.getHitType(hit) << endl;
0173 }
0174 }
0175 }