File indexing completed on 2023-04-15 01:47:32
0001 #ifndef RecoTracker_PixelLowPtUtilities_ClusterShapeHitFilter_h
0002 #define RecoTracker_PixelLowPtUtilities_ClusterShapeHitFilter_h
0003
0004 #include "TrackingTools/TrajectoryFiltering/interface/TrajectoryFilter.h"
0005
0006 #include "DataFormats/GeometryVector/interface/LocalVector.h"
0007 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0008
0009 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0010 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0011 #include "DataFormats/SiPixelCluster/interface/SiPixelClusterShapeCache.h"
0012
0013 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0014
0015 #include "RecoTracker/PixelLowPtUtilities/interface/ClusterData.h"
0016
0017 #include <utility>
0018 #include <unordered_map>
0019 #include <cstring>
0020
0021
0022 class PixelKeys {
0023 public:
0024 PixelKeys(int part, int dx, int dy) : key((part == 0) ? barrelPacking(dx, dy) : endcapPacking(dx, dy)) {}
0025
0026 operator unsigned int() const { return key; }
0027
0028 static unsigned char endcapPacking(int dx, int dy) {
0029 if (dx < 0 || dy < 0)
0030 return N;
0031 if (dx > 10 || dy > 4)
0032 return N;
0033 return N_barrel + dx * offset_endcap_dy + dy;
0034 }
0035
0036 static unsigned char barrelPacking(int dx, int dy) {
0037 if (dx < 0 || dy < 0)
0038 return N;
0039 if (dx > 10 || dy > 15)
0040 return N;
0041 if (dx < 8)
0042 return dx * 16 + dy;
0043 if (dy > 2)
0044 return N;
0045 return 128 + (dx - 8) * 3 + dy;
0046 }
0047
0048 bool isValid() const { return key < N; }
0049
0050 static const int offset_endcap_dy = 5;
0051 static const int offset_endcap_dx = 10;
0052 static const int N_endcap = 55;
0053 static const int N_barrel = 137;
0054 static const int N = N_barrel + N_endcap;
0055
0056 bool operator<(const PixelKeys& right) const { return key < right.key; }
0057
0058 private:
0059 unsigned char key;
0060 };
0061
0062 class StripKeys {
0063 public:
0064 static const int N = 40;
0065
0066 StripKeys(int width) : key(width > 0 ? width - 1 : N) {}
0067
0068 operator unsigned int() const { return key; }
0069
0070 bool isValid() const { return key < N; }
0071
0072 bool operator<(const StripKeys& right) const { return key < right.key; }
0073
0074 private:
0075 unsigned char key;
0076 };
0077
0078 struct PixelLimits {
0079 PixelLimits() {
0080
0081 auto limit = data[0];
0082 limit[0][0] = -10e10;
0083 limit[0][1] = 10e10;
0084 limit[1][0] = -10e10;
0085 limit[1][1] = 10e10;
0086 limit = data[1];
0087 limit[0][0] = -10e10;
0088 limit[0][1] = 10e10;
0089 limit[1][0] = -10e10;
0090 limit[1][1] = 10e10;
0091 }
0092
0093 float data[2][2][2];
0094
0095 bool isInside(const std::pair<float, float>& pred) const {
0096 auto limit = data[0];
0097 bool one = (pred.first > limit[0][0]) && (pred.first < limit[0][1]) && (pred.second > limit[1][0]) &&
0098 (pred.second < limit[1][1]);
0099
0100 limit = data[1];
0101 bool two = (pred.first > limit[0][0]) && (pred.first < limit[0][1]) && (pred.second > limit[1][0]) &&
0102 (pred.second < limit[1][1]);
0103
0104 return one || two;
0105 }
0106 };
0107
0108 struct StripLimits {
0109 StripLimits() {
0110 data[0][0] = -10e10;
0111 data[0][1] = 10e10;
0112 data[1][0] = -10e10;
0113 data[1][1] = 10e10;
0114 }
0115
0116 float data[2][2];
0117
0118 bool isInside(float pred) const {
0119 float const* limit = data[0];
0120 bool one = pred > limit[0] && pred < limit[1];
0121 limit = data[1];
0122 bool two = pred > limit[0] && pred < limit[1];
0123
0124 return one || two;
0125 }
0126 };
0127
0128
0129 namespace edm {
0130 class EventSetup;
0131 }
0132
0133 class TrackerGeometry;
0134 class TrackerTopology;
0135 class MagneticField;
0136 class SiPixelLorentzAngle;
0137 class SiStripLorentzAngle;
0138 class PixelGeomDetUnit;
0139 class StripGeomDetUnit;
0140 class StripTopology;
0141
0142
0143 namespace test {
0144 namespace ClusterShapeHitFilterTest {
0145 int test();
0146 }
0147 }
0148
0149 class ClusterShapeHitFilter {
0150
0151 friend int test::ClusterShapeHitFilterTest::test();
0152
0153 public:
0154 struct PixelData {
0155 const PixelGeomDetUnit* det;
0156 unsigned short part;
0157 unsigned short layer;
0158 std::pair<float, float> drift;
0159 std::pair<float, float> cotangent;
0160 };
0161
0162 struct StripData {
0163 const StripGeomDetUnit* det;
0164 StripTopology const* topology;
0165 float drift;
0166 float thickness;
0167 int nstrips;
0168 };
0169
0170 typedef TrajectoryFilter::Record Record;
0171
0172 ClusterShapeHitFilter(const TrackerGeometry* theTracker_,
0173 const TrackerTopology* theTkTopol_,
0174 const MagneticField* theMagneticField_,
0175 const SiPixelLorentzAngle* theSiPixelLorentzAngle_,
0176 const SiStripLorentzAngle* theSiStripLorentzAngle_,
0177 const std::string& pixelShapeFile_,
0178 const std::string& pixelShapeFileL1_);
0179
0180 ~ClusterShapeHitFilter();
0181
0182 void setShapeCuts(bool cutOnPixelShape, bool cutOnStripShape) {
0183 cutOnPixelShape_ = cutOnPixelShape;
0184 cutOnStripShape_ = cutOnStripShape;
0185 }
0186
0187 void setChargeCuts(bool cutOnPixelCharge, float minGoodPixelCharge, bool cutOnStripCharge, float minGoodStripCharge) {
0188 cutOnPixelCharge_ = cutOnPixelCharge;
0189 minGoodPixelCharge_ = minGoodPixelCharge;
0190 cutOnStripCharge_ = cutOnStripCharge;
0191 minGoodStripCharge_ = minGoodStripCharge;
0192 }
0193
0194 bool getSizes(const SiPixelRecHit& recHit,
0195 const LocalVector& ldir,
0196 const SiPixelClusterShapeCache& clusterShapeCache,
0197 int& part,
0198 ClusterData::ArrayType& meas,
0199 std::pair<float, float>& predr,
0200 PixelData const* pd = nullptr) const;
0201 bool isCompatible(const SiPixelRecHit& recHit,
0202 const LocalVector& ldir,
0203 const SiPixelClusterShapeCache& clusterShapeCache,
0204 PixelData const* pd = nullptr) const;
0205 bool isCompatible(const SiPixelRecHit& recHit,
0206 const GlobalVector& gdir,
0207 const SiPixelClusterShapeCache& clusterShapeCache,
0208 PixelData const* pd = nullptr) const;
0209
0210 bool getSizes(DetId detId,
0211 const SiStripCluster& cluster,
0212 const LocalPoint& lpos,
0213 const LocalVector& ldir,
0214 int& meas,
0215 float& pred) const;
0216 bool getSizes(
0217 const SiStripRecHit2D& recHit, const LocalPoint& lpos, const LocalVector& ldir, int& meas, float& pred) const {
0218 return getSizes(recHit.geographicalId(), recHit.stripCluster(), lpos, ldir, meas, pred);
0219 }
0220 bool isCompatible(DetId detId, const SiStripCluster& cluster, const LocalPoint& lpos, const LocalVector& ldir) const;
0221 bool isCompatible(DetId detId, const SiStripCluster& cluster, const LocalVector& ldir) const {
0222 return isCompatible(detId, cluster, LocalPoint(0, 0, 0), ldir);
0223 }
0224
0225 bool isCompatible(DetId detId, const SiStripCluster& cluster, const GlobalPoint& gpos, const GlobalVector& gdir) const;
0226 bool isCompatible(DetId detId, const SiStripCluster& cluster, const GlobalVector& gdir) const;
0227
0228 bool isCompatible(const SiStripRecHit2D& recHit, const LocalPoint& lpos, const LocalVector& ldir) const {
0229 return isCompatible(recHit.geographicalId(), recHit.stripCluster(), lpos, ldir);
0230 }
0231 bool isCompatible(const SiStripRecHit2D& recHit, const LocalVector& ldir) const {
0232 return isCompatible(recHit.geographicalId(), recHit.stripCluster(), ldir);
0233 }
0234 bool isCompatible(const SiStripRecHit2D& recHit, const GlobalPoint& gpos, const GlobalVector& gdir) const {
0235 return isCompatible(recHit.geographicalId(), recHit.stripCluster(), gpos, gdir);
0236 }
0237 bool isCompatible(const SiStripRecHit2D& recHit, const GlobalVector& gdir) const {
0238 return isCompatible(recHit.geographicalId(), recHit.stripCluster(), gdir);
0239 }
0240
0241 private:
0242
0243 ClusterShapeHitFilter(std::string const& f1, std::string const& f2) {
0244 loadPixelLimits(f1, pixelLimits);
0245 loadPixelLimits(f2, pixelLimitsL1);
0246 loadStripLimits();
0247 }
0248
0249 const PixelData& getpd(const SiPixelRecHit& recHit, PixelData const* pd = nullptr) const {
0250 if (pd)
0251 return *pd;
0252
0253 DetId id = recHit.geographicalId();
0254 auto p = pixelData.find(id);
0255 return (*p).second;
0256 }
0257
0258 const StripData& getsd(DetId id) const { return stripData.find(id)->second; }
0259
0260 void loadPixelLimits(std::string const& file, PixelLimits* plim);
0261 void loadStripLimits();
0262 void fillPixelData();
0263 void fillStripData();
0264
0265 std::pair<float, float> getCotangent(const PixelGeomDetUnit* pixelDet) const;
0266 float getCotangent(const ClusterShapeHitFilter::StripData& sd, const LocalPoint& p) const;
0267
0268 std::pair<float, float> getDrift(const PixelGeomDetUnit* pixelDet) const;
0269 float getDrift(const StripGeomDetUnit* stripDet) const;
0270
0271 bool isNormalOriented(const GeomDetUnit* geomDet) const;
0272
0273 const TrackerGeometry* theTracker;
0274 const TrackerTopology* theTkTopol;
0275 const MagneticField* theMagneticField;
0276
0277 const SiPixelLorentzAngle* theSiPixelLorentzAngle;
0278 const SiStripLorentzAngle* theSiStripLorentzAngle;
0279
0280 std::unordered_map<unsigned int, PixelData> pixelData;
0281 std::unordered_map<unsigned int, StripData> stripData;
0282
0283 PixelLimits pixelLimits[PixelKeys::N + 1];
0284 PixelLimits pixelLimitsL1[PixelKeys::N + 1];
0285
0286 StripLimits stripLimits[StripKeys::N + 1];
0287
0288 float theAngle[6];
0289 bool cutOnPixelCharge_, cutOnStripCharge_;
0290 float minGoodPixelCharge_, minGoodStripCharge_;
0291 bool cutOnPixelShape_, cutOnStripShape_;
0292 bool checkClusterCharge(DetId detId, const SiStripCluster& cluster, const LocalVector& ldir) const;
0293 bool checkClusterCharge(DetId detId, const SiPixelCluster& cluster, const LocalVector& ldir) const;
0294 };
0295
0296 #endif