File indexing completed on 2024-04-06 12:21:52
0001 #include "L1Trigger/TrackFindingTMTT/interface/StubKiller.h"
0002 #include "DataFormats/Math/interface/deltaPhi.h"
0003 #include "FWCore/Utilities/interface/Exception.h"
0004
0005 using namespace std;
0006
0007 namespace tmtt {
0008
0009 StubKiller::StubKiller(StubKiller::KillOptions killScenario,
0010 const TrackerTopology* trackerTopology,
0011 const TrackerGeometry* trackerGeometry,
0012 const edm::Event& event)
0013 : killScenario_(killScenario),
0014 trackerTopology_(trackerTopology),
0015 trackerGeometry_(trackerGeometry),
0016 minPhiToKill_(0),
0017 maxPhiToKill_(0),
0018 minZToKill_(0),
0019 maxZToKill_(0),
0020 minRToKill_(0),
0021 maxRToKill_(0),
0022 fractionOfStubsToKillInLayers_(0),
0023 fractionOfStubsToKillEverywhere_(0),
0024 fractionOfModulesToKillEverywhere_(0) {
0025 if (rndmService_.isAvailable()) {
0026 rndmEngine_ = &(rndmService_->getEngine(event.streamID()));
0027 } else {
0028 throw cms::Exception("BadConfig")
0029 << "StubKiller: requires RandomNumberGeneratorService, not present in cfg file, namely:" << endl
0030 << "process.RandomNumberGeneratorService=cms.Service('RandomNumberGeneratorService',TMTrackProducer=cms.PSet("
0031 "initialSeed=cms.untracked.uint32(12345)))";
0032 }
0033
0034
0035
0036
0037
0038 if (killScenario_ == KillOptions::layer5) {
0039 layersToKill_ = {5};
0040 minPhiToKill_ = 0;
0041 maxPhiToKill_ = 0.5 * M_PI;
0042 minZToKill_ = -1000;
0043 maxZToKill_ = 0;
0044 minRToKill_ = 0;
0045 maxRToKill_ = 1000;
0046 fractionOfStubsToKillInLayers_ = 1;
0047 fractionOfStubsToKillEverywhere_ = 0;
0048 fractionOfModulesToKillEverywhere_ = 0.05;
0049 }
0050
0051
0052 else if (killScenario_ == KillOptions::layer1) {
0053 layersToKill_ = {1};
0054 minPhiToKill_ = 0;
0055 maxPhiToKill_ = 0.5 * M_PI;
0056 minZToKill_ = -1000;
0057 maxZToKill_ = 0;
0058 minRToKill_ = 0;
0059 maxRToKill_ = 1000;
0060 fractionOfStubsToKillInLayers_ = 1;
0061 fractionOfStubsToKillEverywhere_ = 0;
0062 fractionOfModulesToKillEverywhere_ = 0.05;
0063 }
0064
0065
0066 else if (killScenario_ == KillOptions::layer1layer2) {
0067 layersToKill_ = {1, 2};
0068 minPhiToKill_ = 0;
0069 maxPhiToKill_ = 0.5 * M_PI;
0070 minZToKill_ = -1000;
0071 maxZToKill_ = 0;
0072 minRToKill_ = 0;
0073 maxRToKill_ = 1000;
0074 fractionOfStubsToKillInLayers_ = 1;
0075 fractionOfStubsToKillEverywhere_ = 0;
0076 fractionOfModulesToKillEverywhere_ = 0;
0077 }
0078
0079
0080 else if (killScenario_ == KillOptions::layer1disk1) {
0081 layersToKill_ = {1, 11};
0082 minPhiToKill_ = 0;
0083 maxPhiToKill_ = 0.5 * M_PI;
0084 minZToKill_ = -1000;
0085 maxZToKill_ = 0;
0086 minRToKill_ = 0;
0087 maxRToKill_ = 66.5;
0088 fractionOfStubsToKillInLayers_ = 1;
0089 fractionOfStubsToKillEverywhere_ = 0;
0090 fractionOfModulesToKillEverywhere_ = 0;
0091 }
0092
0093
0094 else if (killScenario_ == KillOptions::random) {
0095 layersToKill_ = {};
0096 fractionOfStubsToKillInLayers_ = 0;
0097 fractionOfStubsToKillEverywhere_ = 0.;
0098 fractionOfModulesToKillEverywhere_ = 0.05;
0099 }
0100
0101 deadModules_.clear();
0102 if (fractionOfModulesToKillEverywhere_ > 0) {
0103 this->chooseModulesToKill();
0104 }
0105 this->addDeadLayerModulesToDeadModuleList();
0106 }
0107
0108
0109
0110 bool StubKiller::killStub(const TTStub<Ref_Phase2TrackerDigi_>* stub) const {
0111 if (killScenario_ == KillOptions::none)
0112 return false;
0113 else {
0114
0115
0116 bool killStubRandomly = killStub(stub,
0117 layersToKill_,
0118 minPhiToKill_,
0119 maxPhiToKill_,
0120 minZToKill_,
0121 maxZToKill_,
0122 minRToKill_,
0123 maxRToKill_,
0124 fractionOfStubsToKillInLayers_,
0125 fractionOfStubsToKillEverywhere_);
0126
0127
0128 bool killStubInDeadModules = killStubInDeadModule(stub);
0129 return killStubRandomly || killStubInDeadModules;
0130 }
0131 }
0132
0133
0134
0135
0136
0137
0138
0139
0140
0141 bool StubKiller::killStub(const TTStub<Ref_Phase2TrackerDigi_>* stub,
0142 const vector<int>& layersToKill,
0143 const double minPhiToKill,
0144 const double maxPhiToKill,
0145 const double minZToKill,
0146 const double maxZToKill,
0147 const double minRToKill,
0148 const double maxRToKill,
0149 const double fractionOfStubsToKillInLayers,
0150 const double fractionOfStubsToKillEverywhere) const {
0151
0152 if (not layersToKill.empty()) {
0153
0154 DetId stackDetid = stub->getDetId();
0155 DetId geoDetId(stackDetid.rawId() + 1);
0156
0157
0158 if (deadModules_.empty() || deadModules_.find(geoDetId) == deadModules_.end()) {
0159 bool isInBarrel = geoDetId.subdetId() == StripSubdetector::TOB || geoDetId.subdetId() == StripSubdetector::TIB;
0160
0161 int layerID = 0;
0162 if (isInBarrel) {
0163 layerID = trackerTopology_->layer(geoDetId);
0164 } else {
0165 layerID = 10 * trackerTopology_->side(geoDetId) + trackerTopology_->tidWheel(geoDetId);
0166 }
0167
0168 if (find(layersToKill.begin(), layersToKill.end(), layerID) != layersToKill.end()) {
0169
0170 const GeomDetUnit* det0 = trackerGeometry_->idToDetUnit(geoDetId);
0171 const PixelGeomDetUnit* theGeomDet = dynamic_cast<const PixelGeomDetUnit*>(det0);
0172 const PixelTopology* topol = dynamic_cast<const PixelTopology*>(&(theGeomDet->specificTopology()));
0173 MeasurementPoint measurementPoint = stub->clusterRef(0)->findAverageLocalCoordinatesCentered();
0174 LocalPoint clustlp = topol->localPosition(measurementPoint);
0175 GlobalPoint pos = theGeomDet->surface().toGlobal(clustlp);
0176
0177 double stubPhi = reco::deltaPhi(pos.phi(), 0.);
0178
0179 if (stubPhi > minPhiToKill && stubPhi < maxPhiToKill && pos.z() > minZToKill && pos.z() < maxZToKill &&
0180 pos.perp() > minRToKill && pos.perp() < maxRToKill) {
0181
0182 if (fractionOfStubsToKillInLayers == 1) {
0183 return true;
0184 } else {
0185 if (rndmEngine_->flat() < fractionOfStubsToKillInLayers) {
0186 return true;
0187 }
0188 }
0189 }
0190 }
0191 }
0192 }
0193
0194
0195 if (fractionOfStubsToKillEverywhere > 0) {
0196 if (rndmEngine_->flat() < fractionOfStubsToKillEverywhere) {
0197 return true;
0198 }
0199 }
0200
0201 return false;
0202 }
0203
0204
0205
0206 bool StubKiller::killStubInDeadModule(const TTStub<Ref_Phase2TrackerDigi_>* stub) const {
0207 if (not deadModules_.empty()) {
0208 DetId stackDetid = stub->getDetId();
0209 DetId geoDetId(stackDetid.rawId() + 1);
0210 auto deadModule = deadModules_.find(geoDetId);
0211 if (deadModule != deadModules_.end()) {
0212 if (deadModule->second == 1) {
0213 return true;
0214 } else {
0215 if (rndmEngine_->flat() < deadModule->second) {
0216 return true;
0217 }
0218 }
0219 }
0220 }
0221
0222 return false;
0223 }
0224
0225
0226
0227 void StubKiller::chooseModulesToKill() {
0228 for (const GeomDetUnit* gd : trackerGeometry_->detUnits()) {
0229 if (!trackerTopology_->isLower(gd->geographicalId()))
0230 continue;
0231 if (rndmEngine_->flat() < fractionOfModulesToKillEverywhere_) {
0232 deadModules_[gd->geographicalId()] = 1;
0233 }
0234 }
0235 }
0236
0237
0238
0239 void StubKiller::addDeadLayerModulesToDeadModuleList() {
0240 for (const GeomDetUnit* gd : trackerGeometry_->detUnits()) {
0241 float moduleR = gd->position().perp();
0242 float moduleZ = gd->position().z();
0243 float modulePhi = reco::deltaPhi(gd->position().phi(), 0.);
0244 DetId geoDetId = gd->geographicalId();
0245 bool isInBarrel = geoDetId.subdetId() == StripSubdetector::TOB || geoDetId.subdetId() == StripSubdetector::TIB;
0246
0247 int layerID = 0;
0248 if (isInBarrel) {
0249 layerID = trackerTopology_->layer(geoDetId);
0250 } else {
0251 layerID = 10 * trackerTopology_->side(geoDetId) + trackerTopology_->tidWheel(geoDetId);
0252 }
0253 if (find(layersToKill_.begin(), layersToKill_.end(), layerID) != layersToKill_.end()) {
0254 if (modulePhi > minPhiToKill_ && modulePhi < maxPhiToKill_ && moduleZ > minZToKill_ && moduleZ < maxZToKill_ &&
0255 moduleR > minRToKill_ && moduleR < maxRToKill_) {
0256 if (deadModules_.find(gd->geographicalId()) == deadModules_.end()) {
0257 deadModules_[gd->geographicalId()] = fractionOfStubsToKillInLayers_;
0258 }
0259 }
0260 }
0261 }
0262 }
0263 };