TrackFindingTMTT/src/StubKiller.cc

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     // These scenarios correspond to slide 12 of  https://indico.cern.ch/event/719985/contributions/2970687/attachments/1634587/2607365/StressTestTF-Acosta-Apr18.pdf
0035     // Scenario 1
0036
0037     // kill layer 5 in one quadrant +5 % random module loss to connect to what was done before
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     // Scenario 2
0051     // kill layer 1 in one quadrant +5 % random module loss
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     // Scenario 3
0065     // kill layer 1 + layer 2, both in same quadrant
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     // Scenario 4
0079     // kill layer 1 and disk 1, both in same quadrant
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     // An extra scenario not listed in the slides
0093     // 5% random module loss throughout tracker
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   // Indicate if given stub was killed by dead tracker module, based on dead module scenario.
0109
0110   bool StubKiller::killStub(const TTStub<Ref_Phase2TrackerDigi_>* stub) const {
0111     if (killScenario_ == KillOptions::none)
0112       return false;
0113     else {
0114       // Check if stub is in dead region specified by *ToKill_
0115       // And randomly kill stubs throughout tracker (not just thos in specific regions/modules)
0116       bool killStubRandomly = killStub(stub,
0117                                        layersToKill_,
0118                                        minPhiToKill_,
0119                                        maxPhiToKill_,
0120                                        minZToKill_,
0121                                        maxZToKill_,
0122                                        minRToKill_,
0123                                        maxRToKill_,
0124                                        fractionOfStubsToKillInLayers_,
0125                                        fractionOfStubsToKillEverywhere_);
0126       // Kill modules in specifid modules
0127       // Random modules throughout the tracker, and those modules in specific regions (so may already have been killed by killStub above)
0128       bool killStubInDeadModules = killStubInDeadModule(stub);
0129       return killStubRandomly || killStubInDeadModules;
0130     }
0131   }
0132
0133   // Indicate if given stub was killed by dead tracker module, based on specified dead regions
0134   // rather than based on the dead module scenario.
0135   // layersToKill - a vector stating the layers we are killing stubs in.  Can be an empty vector.
0136   // Barrel layers are encoded as 1-6. The endcap layers are encoded as 11-15 (-z) and 21-25 (+z)
0137   // min/max Phi/Z/R - stubs within the region specified by these boundaries and layersToKill are flagged for killing
0138   // fractionOfStubsToKillInLayers - The fraction of stubs to kill in the specified layers/region.
0139   // fractionOfStubsToKillEverywhere - The fraction of stubs to kill throughout the tracker
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     // Only kill stubs in specified layers
0152     if (not layersToKill.empty()) {
0153       // Get the layer the stub is in, and check if it's in the layer you want to kill
0154       DetId stackDetid = stub->getDetId();
0155       DetId geoDetId(stackDetid.rawId() + 1);
0156
0157       // If this module is in the deadModule list, don't also try to kill the stub here
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           // Get the phi and z of stub, and check if it's in the region you want to kill
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             // Kill fraction of stubs
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     // Kill fraction of stubs throughout tracker
0195     if (fractionOfStubsToKillEverywhere > 0) {
0196       if (rndmEngine_->flat() < fractionOfStubsToKillEverywhere) {
0197         return true;
0198       }
0199     }
0200
0201     return false;
0202   }
0203
0204   // Indicate if given stub was in (partially) dead tracker module, based on dead module scenario.
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   // Identify modules to be killed, chosen randomly from those in the whole tracker.
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   //  Identify modules to be killed, chosen based on location in tracker.
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 };  // namespace tmtt