Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:02:26

0001 #include "MeasurementTrackerImpl.h"
0002 
0003 #include "FWCore/Framework/interface/ESHandle.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 
0006 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0007 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0008 #include "Geometry/CommonDetUnit/interface/GluedGeomDet.h"
0009 #include "Geometry/CommonDetUnit/interface/StackGeomDet.h"
0010 #include "Geometry/CommonDetUnit/interface/DoubleSensGeomDet.h"
0011 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0012 
0013 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0014 
0015 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0016 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0017 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0018 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
0019 #include "DataFormats/Common/interface/ContainerMask.h"
0020 
0021 #include "TrackingTools/MeasurementDet/interface/MeasurementDetException.h"
0022 
0023 #include "RecoLocalTracker/ClusterParameterEstimator/interface/PixelClusterParameterEstimator.h"
0024 #include "RecoLocalTracker/Records/interface/TrackerCPERecord.h"
0025 #include "RecoLocalTracker/SiStripRecHitConverter/interface/SiStripRecHitMatcher.h"
0026 #include "RecoLocalTracker/SiStripRecHitConverter/interface/StripCPE.h"
0027 
0028 #include "RecoTracker/TkDetLayers/interface/GeometricSearchTracker.h"
0029 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
0030 #include "TkStripMeasurementDet.h"
0031 #include "TkPixelMeasurementDet.h"
0032 #include "TkPhase2OTMeasurementDet.h"
0033 #include "TkGluedMeasurementDet.h"
0034 #include "TkStackMeasurementDet.h"
0035 #include "TkDoubleSensMeasurementDet.h"
0036 
0037 #include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
0038 #include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"
0039 
0040 #include <iostream>
0041 #include <typeinfo>
0042 #include <map>
0043 #include <algorithm>
0044 
0045 //
0046 
0047 using namespace std;
0048 
0049 namespace {
0050 
0051   class StrictWeakOrdering {
0052   public:
0053     bool operator()(uint32_t p, const uint32_t& i) const { return p < i; }
0054   };
0055 
0056   struct CmpTKD {
0057     bool operator()(MeasurementDet const* rh, MeasurementDet const* lh) {
0058       return rh->fastGeomDet().geographicalId().rawId() < lh->fastGeomDet().geographicalId().rawId();
0059     }
0060     bool operator()(MeasurementDet const& rh, MeasurementDet const& lh) {
0061       return rh.fastGeomDet().geographicalId().rawId() < lh.fastGeomDet().geographicalId().rawId();
0062     }
0063   };
0064 
0065   template <typename TKD>
0066   void sortTKD(std::vector<TKD*>& det) {
0067     std::sort(det.begin(), det.end(), CmpTKD());
0068   }
0069   template <typename TKD>
0070   void sortTKD(std::vector<TKD>& det) {
0071     std::sort(det.begin(), det.end(), CmpTKD());
0072   }
0073 
0074 }  // namespace
0075 
0076 MeasurementTrackerImpl::MeasurementTrackerImpl(const BadStripCutsDet& badStripCuts,
0077                                                const PixelClusterParameterEstimator* pixelCPE,
0078                                                const StripClusterParameterEstimator* stripCPE,
0079                                                const SiStripRecHitMatcher* hitMatcher,
0080                                                const TrackerTopology* trackerTopology,
0081                                                const TrackerGeometry* trackerGeom,
0082                                                const GeometricSearchTracker* geometricSearchTracker,
0083                                                const SiStripQuality* stripQuality,
0084                                                int stripQualityFlags,
0085                                                int stripQualityDebugFlags,
0086                                                const SiPixelQuality* pixelQuality,
0087                                                const SiPixelFedCabling* pixelCabling,
0088                                                int pixelQualityFlags,
0089                                                int pixelQualityDebugFlags,
0090                                                const ClusterParameterEstimator<Phase2TrackerCluster1D>* phase2OTCPE)
0091     : MeasurementTracker(trackerGeom, geometricSearchTracker),
0092       theStDetConditions(hitMatcher, stripCPE),
0093       thePxDetConditions(pixelCPE),
0094       thePhase2DetConditions(phase2OTCPE) {
0095   this->initialize(trackerTopology);
0096   this->initializeStripStatus(badStripCuts, stripQuality, stripQualityFlags, stripQualityDebugFlags);
0097   this->initializePixelStatus(pixelQuality, pixelCabling, pixelQualityFlags, pixelQualityDebugFlags);
0098 }
0099 
0100 MeasurementTrackerImpl::~MeasurementTrackerImpl() {}
0101 
0102 void MeasurementTrackerImpl::initialize(const TrackerTopology* trackerTopology) {
0103   bool subIsPixel = false;
0104   //FIXME:just temporary solution for phase2 :
0105   //the OT is defined as PixelSubDetector!
0106   bool subIsOT = false;
0107 
0108   //if the TkGeometry has the subDet vector filled, the theDetMap is filled, otherwise nothing should happen
0109   if (!theTrackerGeom->detsPXB().empty()) {
0110     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0111         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsPXB().front()->geographicalId().subdetId()));
0112     addDets(theTrackerGeom->detsPXB(), subIsPixel, subIsOT);
0113   }
0114 
0115   if (!theTrackerGeom->detsPXF().empty()) {
0116     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0117         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsPXF().front()->geographicalId().subdetId()));
0118     addDets(theTrackerGeom->detsPXF(), subIsPixel, subIsOT);
0119   }
0120 
0121   subIsOT = true;
0122 
0123   if (!theTrackerGeom->detsTIB().empty()) {
0124     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0125         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsTIB().front()->geographicalId().subdetId()));
0126     addDets(theTrackerGeom->detsTIB(), subIsPixel, subIsOT);
0127   }
0128 
0129   if (!theTrackerGeom->detsTID().empty()) {
0130     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0131         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsTID().front()->geographicalId().subdetId()));
0132     addDets(theTrackerGeom->detsTID(), subIsPixel, subIsOT);
0133   }
0134 
0135   if (!theTrackerGeom->detsTOB().empty()) {
0136     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0137         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsTOB().front()->geographicalId().subdetId()));
0138     addDets(theTrackerGeom->detsTOB(), subIsPixel, subIsOT);
0139   }
0140 
0141   if (!theTrackerGeom->detsTEC().empty()) {
0142     subIsPixel = GeomDetEnumerators::isTrackerPixel(
0143         theTrackerGeom->geomDetSubDetector(theTrackerGeom->detsTEC().front()->geographicalId().subdetId()));
0144     addDets(theTrackerGeom->detsTEC(), subIsPixel, subIsOT);
0145   }
0146 
0147   // fist all stripdets
0148   sortTKD(theStripDets);
0149   initStMeasurementConditionSet(theStripDets);
0150   for (unsigned int i = 0; i != theStripDets.size(); ++i)
0151     theDetMap[theStDetConditions.id(i)] = &theStripDets[i];
0152 
0153   // now the glued dets
0154   sortTKD(theGluedDets);
0155   for (unsigned int i = 0; i != theGluedDets.size(); ++i)
0156     initGluedDet(theGluedDets[i], trackerTopology);
0157 
0158   // then the pixels
0159   sortTKD(thePixelDets);
0160   initPxMeasurementConditionSet(thePixelDets);
0161   for (unsigned int i = 0; i != thePixelDets.size(); ++i)
0162     theDetMap[thePxDetConditions.id(i)] = &thePixelDets[i];
0163 
0164   // then the phase2 dets
0165   sortTKD(thePhase2Dets);
0166   initPhase2OTMeasurementConditionSet(thePhase2Dets);
0167   for (unsigned int i = 0; i != thePhase2Dets.size(); ++i)
0168     theDetMap[thePhase2DetConditions.id(i)] = &thePhase2Dets[i];
0169 
0170   // and then the stack dets, at last
0171   sortTKD(theStackDets);
0172   for (unsigned int i = 0; i != theStackDets.size(); ++i)
0173     initStackDet(theStackDets[i]);
0174 
0175   // and then the double sensor dets
0176   sortTKD(theDoubleSensGeomDets);
0177   for (unsigned int i = 0; i != theDoubleSensGeomDets.size(); ++i)
0178     initDoubleSensDet(theDoubleSensGeomDets[i]);
0179 
0180   if (!checkDets())
0181     throw MeasurementDetException("Number of dets in MeasurementTracker not consistent with TrackerGeometry!");
0182 }
0183 
0184 void MeasurementTrackerImpl::initStMeasurementConditionSet(std::vector<TkStripMeasurementDet>& stripDets) {
0185   // assume vector is full and ordered!
0186   int size = stripDets.size();
0187   theStDetConditions.init(size);
0188   for (int i = 0; i != size; ++i) {
0189     auto& mdet = stripDets[i];
0190     mdet.setIndex(i);
0191     //intialize the detId !
0192     theStDetConditions.id_[i] = mdet.specificGeomDet().geographicalId().rawId();
0193     theStDetConditions.subId_[i] = SiStripDetId(theStDetConditions.id_[i]).subdetId() - 3;
0194     //initalize the total number of strips
0195     theStDetConditions.totalStrips_[i] = mdet.specificGeomDet().specificTopology().nstrips();
0196   }
0197 }
0198 
0199 void MeasurementTrackerImpl::initPxMeasurementConditionSet(std::vector<TkPixelMeasurementDet>& pixelDets) {
0200   // assume vector is full and ordered!
0201   int size = pixelDets.size();
0202   thePxDetConditions.init(size);
0203 
0204   for (int i = 0; i != size; ++i) {
0205     auto& mdet = pixelDets[i];
0206     mdet.setIndex(i);
0207     thePxDetConditions.id_[i] = mdet.specificGeomDet().geographicalId().rawId();
0208   }
0209 }
0210 
0211 void MeasurementTrackerImpl::initPhase2OTMeasurementConditionSet(std::vector<TkPhase2OTMeasurementDet>& phase2Dets) {
0212   // assume vector is full and ordered!
0213   int size = phase2Dets.size();
0214   thePhase2DetConditions.init(size);
0215 
0216   for (int i = 0; i != size; ++i) {
0217     auto& mdet = phase2Dets[i];
0218     mdet.setIndex(i);
0219     thePhase2DetConditions.id_[i] = mdet.specificGeomDet().geographicalId().rawId();
0220   }
0221 }
0222 
0223 void MeasurementTrackerImpl::addDets(const TrackingGeometry::DetContainer& dets, bool subIsPixel, bool subIsOT) {
0224   //in phase2, we can have composed subDetector made by Pixel or Strip
0225   for (TrackerGeometry::DetContainer::const_iterator gd = dets.begin(); gd != dets.end(); gd++) {
0226     const GeomDetUnit* gdu = dynamic_cast<const GeomDetUnit*>(*gd);
0227 
0228     //Pixel or Strip GeomDetUnit
0229     if (gdu->isLeaf()) {
0230       if (subIsPixel) {
0231         if (!subIsOT) {
0232           addPixelDet(*gd);
0233         } else {
0234           addPhase2Det(*gd);
0235         }
0236       } else {
0237         addStripDet(*gd);
0238       }
0239     } else {
0240       //Glued or Stack GeomDet
0241       const GluedGeomDet* gluedDet = dynamic_cast<const GluedGeomDet*>(*gd);
0242       const StackGeomDet* stackDet = dynamic_cast<const StackGeomDet*>(*gd);
0243       const DoubleSensGeomDet* doubleSensGeomDet = dynamic_cast<const DoubleSensGeomDet*>(*gd);
0244 
0245       if ((gluedDet == nullptr && stackDet == nullptr && doubleSensGeomDet == nullptr) ||
0246           (gluedDet != nullptr && stackDet != nullptr && doubleSensGeomDet != nullptr)) {
0247         throw MeasurementDetException(
0248             "MeasurementTracker ERROR: GeomDet neither DetUnit nor GluedDet nor StackDet nor DoubleSensGeomDet");
0249       }
0250       if (gluedDet != nullptr)
0251         addGluedDet(gluedDet);
0252       else if (stackDet != nullptr)
0253         addStackDet(stackDet);
0254       else
0255         addDoubleSensGeomDet(doubleSensGeomDet);
0256     }
0257   }
0258 }
0259 
0260 bool MeasurementTrackerImpl::checkDets() {
0261   if (theTrackerGeom->dets().size() == theDetMap.size())
0262     return true;
0263   return false;
0264 }
0265 
0266 void MeasurementTrackerImpl::addStripDet(const GeomDet* gd) {
0267   try {
0268     theStripDets.push_back(TkStripMeasurementDet(gd, theStDetConditions));
0269   } catch (MeasurementDetException& err) {
0270     edm::LogError("MeasurementDet") << "Oops, got a MeasurementDetException: " << err.what();
0271   }
0272 }
0273 
0274 void MeasurementTrackerImpl::addPixelDet(const GeomDet* gd) {
0275   try {
0276     thePixelDets.push_back(TkPixelMeasurementDet(gd, thePxDetConditions));
0277   } catch (MeasurementDetException& err) {
0278     edm::LogError("MeasurementDet") << "Oops, got a MeasurementDetException: " << err.what();
0279   }
0280 }
0281 
0282 void MeasurementTrackerImpl::addPhase2Det(const GeomDet* gd) {
0283   try {
0284     thePhase2Dets.push_back(TkPhase2OTMeasurementDet(gd, thePhase2DetConditions));
0285   } catch (MeasurementDetException& err) {
0286     edm::LogError("MeasurementDet") << "Oops, got a MeasurementDetException: " << err.what();
0287   }
0288 }
0289 
0290 void MeasurementTrackerImpl::addGluedDet(const GluedGeomDet* gd) {
0291   theGluedDets.push_back(TkGluedMeasurementDet(gd, theStDetConditions.matcher(), theStDetConditions.stripCPE()));
0292 }
0293 
0294 void MeasurementTrackerImpl::addStackDet(const StackGeomDet* gd) {
0295   //since the Stack will be composed by PS or 2S,
0296   //both cluster parameter estimators are needed? - right now just the thePixelCPE is used.
0297   theStackDets.push_back(TkStackMeasurementDet(gd, thePxDetConditions.pixelCPE()));
0298 }
0299 
0300 void MeasurementTrackerImpl::addDoubleSensGeomDet(const DoubleSensGeomDet* gd) {
0301   theDoubleSensGeomDets.push_back(TkDoubleSensMeasurementDet(gd, thePxDetConditions.pixelCPE()));
0302 }
0303 
0304 void MeasurementTrackerImpl::initGluedDet(TkGluedMeasurementDet& det, const TrackerTopology* trackerTopology) {
0305   const GluedGeomDet& gd = det.specificGeomDet();
0306   const MeasurementDet* monoDet = findDet(gd.monoDet()->geographicalId());
0307   const MeasurementDet* stereoDet = findDet(gd.stereoDet()->geographicalId());
0308   if (monoDet == nullptr || stereoDet == nullptr) {
0309     edm::LogError("MeasurementDet") << "MeasurementTracker ERROR: GluedDet components not found as MeasurementDets ";
0310     throw MeasurementDetException("MeasurementTracker ERROR: GluedDet components not found as MeasurementDets");
0311   }
0312   det.init(monoDet, stereoDet, trackerTopology);
0313   theDetMap[gd.geographicalId()] = &det;
0314 }
0315 
0316 void MeasurementTrackerImpl::initStackDet(TkStackMeasurementDet& det) {
0317   const StackGeomDet& gd = det.specificGeomDet();
0318   const MeasurementDet* lowerDet = findDet(gd.lowerDet()->geographicalId());
0319   const MeasurementDet* upperDet = findDet(gd.upperDet()->geographicalId());
0320   if (lowerDet == nullptr || upperDet == nullptr) {
0321     edm::LogError("MeasurementDet") << "MeasurementTracker ERROR: StackDet components not found as MeasurementDets ";
0322     throw MeasurementDetException("MeasurementTracker ERROR: StackDet components not found as MeasurementDets");
0323   }
0324   det.init(lowerDet, upperDet);
0325   theDetMap[gd.geographicalId()] = &det;
0326 }
0327 
0328 void MeasurementTrackerImpl::initDoubleSensDet(TkDoubleSensMeasurementDet& det) {
0329   const DoubleSensGeomDet& gd = det.specificGeomDet();
0330   const MeasurementDet* firstDet = findDet(gd.firstDet()->geographicalId());
0331   const MeasurementDet* secondDet = findDet(gd.secondDet()->geographicalId());
0332   if (firstDet == nullptr || secondDet == nullptr) {
0333     edm::LogError("MeasurementDet")
0334         << "MeasurementTracker ERROR: DoubleSensDet components not found as MeasurementDets ";
0335     throw MeasurementDetException("MeasurementTracker ERROR: DoubleSensDet components not found as MeasurementDets");
0336   }
0337   det.init(firstDet, secondDet);
0338   theDetMap[gd.geographicalId()] = &det;
0339 }
0340 
0341 void MeasurementTrackerImpl::initializeStripStatus(const BadStripCutsDet& badStripCuts,
0342                                                    const SiStripQuality* quality,
0343                                                    int qualityFlags,
0344                                                    int qualityDebugFlags) {
0345   if (qualityFlags & BadStrips) {
0346     theStDetConditions.badStripCuts_[SiStripDetId::TIB - 3] = badStripCuts.tib;
0347     theStDetConditions.badStripCuts_[SiStripDetId::TOB - 3] = badStripCuts.tob;
0348     theStDetConditions.badStripCuts_[SiStripDetId::TID - 3] = badStripCuts.tid;
0349     theStDetConditions.badStripCuts_[SiStripDetId::TEC - 3] = badStripCuts.tec;
0350   }
0351   theStDetConditions.setMaskBad128StripBlocks((qualityFlags & MaskBad128StripBlocks) != 0);
0352 
0353   if ((quality != nullptr) && (qualityFlags != 0)) {
0354     edm::LogInfo("MeasurementTracker") << "qualityFlags = " << qualityFlags;
0355     unsigned int on = 0, tot = 0;
0356     unsigned int foff = 0, ftot = 0, aoff = 0, atot = 0;
0357     for (int i = 0; i != theStDetConditions.nDet(); i++) {
0358       uint32_t detid = theStDetConditions.id(i);
0359       if (qualityFlags & BadModules) {
0360         bool isOn = quality->IsModuleUsable(detid);
0361         theStDetConditions.setActive(i, isOn);
0362         tot++;
0363         on += (unsigned int)isOn;
0364         if (qualityDebugFlags & BadModules) {
0365           edm::LogInfo("MeasurementTracker")
0366               << "MeasurementTrackerImpl::initializeStripStatus : detid " << detid << " is " << (isOn ? "on" : "off");
0367         }
0368       } else {
0369         theStDetConditions.setActive(i, true);
0370       }
0371       // first turn all APVs and fibers ON
0372       theStDetConditions.set128StripStatus(i, true);
0373       if (qualityFlags & BadAPVFibers) {
0374         short badApvs = quality->getBadApvs(detid);
0375         short badFibers = quality->getBadFibers(detid);
0376         for (int j = 0; j < 6; j++) {
0377           atot++;
0378           if (badApvs & (1 << j)) {
0379             theStDetConditions.set128StripStatus(i, false, j);
0380             aoff++;
0381           }
0382         }
0383         for (int j = 0; j < 3; j++) {
0384           ftot++;
0385           if (badFibers & (1 << j)) {
0386             theStDetConditions.set128StripStatus(i, false, 2 * j);
0387             theStDetConditions.set128StripStatus(i, false, 2 * j + 1);
0388             foff++;
0389           }
0390         }
0391       }
0392       auto& badStrips = theStDetConditions.getBadStripBlocks(i);
0393       badStrips.clear();
0394       if (qualityFlags & BadStrips) {
0395         SiStripBadStrip::Range range = quality->getRange(detid);
0396         for (SiStripBadStrip::ContainerIterator bit = range.first; bit != range.second; ++bit) {
0397           badStrips.push_back(quality->decode(*bit));
0398         }
0399       }
0400     }
0401     if (qualityDebugFlags & BadModules) {
0402       edm::LogInfo("MeasurementTracker StripModuleStatus")
0403           << " Total modules: " << tot << ", active " << on << ", inactive " << (tot - on);
0404     }
0405     if (qualityDebugFlags & BadAPVFibers) {
0406       edm::LogInfo("MeasurementTracker StripAPVStatus")
0407           << " Total APVs: " << atot << ", active " << (atot - aoff) << ", inactive " << (aoff);
0408       edm::LogInfo("MeasurementTracker StripFiberStatus")
0409           << " Total Fibers: " << ftot << ", active " << (ftot - foff) << ", inactive " << (foff);
0410     }
0411   } else {
0412     for (int i = 0; i != theStDetConditions.nDet(); i++) {
0413       theStDetConditions.setActive(i, true);          // module ON
0414       theStDetConditions.set128StripStatus(i, true);  // all APVs and fibers ON
0415     }
0416   }
0417 }
0418 
0419 void MeasurementTrackerImpl::initializePixelStatus(const SiPixelQuality* quality,
0420                                                    const SiPixelFedCabling* pixelCabling,
0421                                                    int qualityFlags,
0422                                                    int qualityDebugFlags) {
0423   if ((quality != nullptr) && (qualityFlags != 0)) {
0424     edm::LogInfo("MeasurementTracker") << "qualityFlags = " << qualityFlags;
0425     unsigned int on = 0, tot = 0, badrocs = 0;
0426     for (std::vector<TkPixelMeasurementDet>::iterator i = thePixelDets.begin(); i != thePixelDets.end(); i++) {
0427       uint32_t detid = ((*i).geomDet().geographicalId()).rawId();
0428       if (qualityFlags & BadModules) {
0429         bool isOn = quality->IsModuleUsable(detid);
0430         (i)->setActive(isOn);
0431         tot++;
0432         on += (unsigned int)isOn;
0433         if (qualityDebugFlags & BadModules) {
0434           edm::LogInfo("MeasurementTracker")
0435               << "MeasurementTrackerImpl::initializePixelStatus : detid " << detid << " is " << (isOn ? "on" : "off");
0436         }
0437       } else {
0438         (i)->setActive(true);
0439       }
0440       if ((qualityFlags & BadROCs) && (quality->getBadRocs(detid) != 0)) {
0441         std::vector<LocalPoint> badROCs = quality->getBadRocPositions(detid, *theTrackerGeom, pixelCabling);
0442         badrocs += badROCs.size();
0443         (i)->setBadRocPositions(badROCs);
0444       } else {
0445         (i)->clearBadRocPositions();
0446       }
0447     }
0448     if (qualityDebugFlags & BadModules) {
0449       edm::LogInfo("MeasurementTracker PixelModuleStatus")
0450           << " Total modules: " << tot << ", active " << on << ", inactive " << (tot - on);
0451     }
0452     if (qualityDebugFlags & BadROCs) {
0453       edm::LogInfo("MeasurementTracker PixelROCStatus") << " Total of bad ROCs: " << badrocs;
0454     }
0455   } else {
0456     for (std::vector<TkPixelMeasurementDet>::iterator i = thePixelDets.begin(); i != thePixelDets.end(); i++) {
0457       (i)->setActive(true);  // module ON
0458     }
0459   }
0460 }