File indexing completed on 2021-02-14 14:26:45
0001 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/TripletFilter.h"
0002
0003 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/ClusterShapeHitFilter.h"
0004 #include "RecoPixelVertexing/PixelLowPtUtilities/interface/HitInfo.h"
0005
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/Framework/interface/ESHandle.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009
0010 #include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h"
0011 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
0012
0013 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0014 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0015
0016 #include "RecoTracker/Record/interface/CkfComponentsRecord.h"
0017 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0018
0019 using namespace std;
0020
0021
0022 TripletFilter::TripletFilter(const edm::EventSetup& es) {
0023
0024 edm::ESHandle<ClusterShapeHitFilter> shape;
0025 es.get<CkfComponentsRecord>().get("ClusterShapeHitFilter", shape);
0026 theFilter = shape.product();
0027 }
0028
0029
0030 TripletFilter::~TripletFilter() {}
0031
0032
0033 bool TripletFilter::checkTrack(const vector<const TrackingRecHit*>& recHits,
0034 const vector<LocalVector>& localDirs,
0035 const TrackerTopology* tTopo,
0036 const SiPixelClusterShapeCache& clusterShapeCache) {
0037 bool ok = true;
0038
0039 vector<LocalVector>::const_iterator localDir = localDirs.begin();
0040 for (vector<const TrackingRecHit*>::const_iterator recHit = recHits.begin(); recHit != recHits.end(); recHit++) {
0041 const SiPixelRecHit* pixelRecHit = dynamic_cast<const SiPixelRecHit*>(*recHit);
0042
0043 if (!pixelRecHit->isValid()) {
0044 ok = false;
0045 break;
0046 }
0047
0048 if (!theFilter->isCompatible(*pixelRecHit, *localDir, clusterShapeCache)) {
0049 LogTrace("MinBiasTracking") << " [TripletFilter] clusShape problem" << HitInfo::getInfo(**recHit, tTopo);
0050
0051 ok = false;
0052 break;
0053 }
0054
0055 localDir++;
0056 }
0057
0058 return ok;
0059 }
0060
0061
0062 bool TripletFilter::checkTrack(const vector<const TrackingRecHit*>& recHits,
0063 const vector<GlobalVector>& globalDirs,
0064 const TrackerTopology* tTopo,
0065 const SiPixelClusterShapeCache& clusterShapeCache) {
0066 bool ok = true;
0067
0068 vector<GlobalVector>::const_iterator globalDir = globalDirs.begin();
0069 for (vector<const TrackingRecHit*>::const_iterator recHit = recHits.begin(); recHit != recHits.end(); recHit++) {
0070 const SiPixelRecHit* pixelRecHit = dynamic_cast<const SiPixelRecHit*>(*recHit);
0071
0072 if (!pixelRecHit->isValid()) {
0073 ok = false;
0074 break;
0075 }
0076
0077 if (!theFilter->isCompatible(*pixelRecHit, *globalDir, clusterShapeCache)) {
0078 LogTrace("MinBiasTracking") << " [TripletFilter] clusShape problem" << HitInfo::getInfo(**recHit, tTopo);
0079
0080 ok = false;
0081 break;
0082 }
0083
0084 globalDir++;
0085 }
0086
0087 return ok;
0088 }