Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-18 03:42:17

0001 #include "FWCore/Framework/interface/global/EDProducer.h"
0002 
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 
0007 #include "CondFormats/SiPixelObjects/interface/SiPixelQuality.h"
0008 #include "CondFormats/DataRecord/interface/SiPixelQualityRcd.h"
0009 
0010 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0011 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0012 
0013 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0014 #include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h"
0015 #include "DataFormats/TrackerCommon/interface/TrackerDetSide.h"
0016 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0017 
0018 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0019 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0020 
0021 #include "RecoTracker/MkFit/interface/MkFitEventOfHits.h"
0022 #include "RecoTracker/MkFit/interface/MkFitGeometry.h"
0023 #include "RecoTracker/MkFit/interface/MkFitHitWrapper.h"
0024 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
0025 
0026 // mkFit includes
0027 #include "RecoTracker/MkFitCore/interface/HitStructures.h"
0028 #include "RecoTracker/MkFitCMS/interface/MkStdSeqs.h"
0029 #include "RecoTracker/MkFitCMS/interface/LayerNumberConverter.h"
0030 
0031 class MkFitEventOfHitsProducer : public edm::global::EDProducer<> {
0032 public:
0033   explicit MkFitEventOfHitsProducer(edm::ParameterSet const& iConfig);
0034   ~MkFitEventOfHitsProducer() override = default;
0035 
0036   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0037 
0038 private:
0039   void produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const override;
0040 
0041   void fill(const std::vector<int>& hits, mkfit::EventOfHits& eventOfHits, const MkFitGeometry& mkFitGeom) const;
0042 
0043   const edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;
0044   const edm::EDGetTokenT<MkFitHitWrapper> pixelHitsToken_;
0045   const edm::EDGetTokenT<MkFitHitWrapper> stripHitsToken_;
0046   const edm::EDGetTokenT<std::vector<int>> pixelLayerIndexToHitToken_;
0047   const edm::EDGetTokenT<std::vector<int>> stripLayerIndexToHitToken_;
0048   const edm::ESGetToken<MkFitGeometry, TrackerRecoGeometryRecord> mkFitGeomToken_;
0049   edm::ESGetToken<SiPixelQuality, SiPixelQualityRcd> pixelQualityToken_;
0050   edm::ESGetToken<SiStripQuality, SiStripQualityRcd> stripQualityToken_;
0051   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0052   const edm::EDPutTokenT<MkFitEventOfHits> putToken_;
0053   const bool usePixelQualityDB_;
0054   const bool useStripStripQualityDB_;
0055 };
0056 
0057 MkFitEventOfHitsProducer::MkFitEventOfHitsProducer(edm::ParameterSet const& iConfig)
0058     : beamSpotToken_{consumes(iConfig.getParameter<edm::InputTag>("beamSpot"))},
0059       pixelHitsToken_{consumes(iConfig.getParameter<edm::InputTag>("pixelHits"))},
0060       stripHitsToken_{consumes(iConfig.getParameter<edm::InputTag>("stripHits"))},
0061       pixelLayerIndexToHitToken_{consumes(iConfig.getParameter<edm::InputTag>("pixelHits"))},
0062       stripLayerIndexToHitToken_{consumes(iConfig.getParameter<edm::InputTag>("stripHits"))},
0063       mkFitGeomToken_{esConsumes()},
0064       putToken_{produces<MkFitEventOfHits>()},
0065       usePixelQualityDB_{iConfig.getParameter<bool>("usePixelQualityDB")},
0066       useStripStripQualityDB_{iConfig.getParameter<bool>("useStripStripQualityDB")} {
0067   if (useStripStripQualityDB_ || usePixelQualityDB_)
0068     geomToken_ = esConsumes();
0069   if (usePixelQualityDB_) {
0070     pixelQualityToken_ = esConsumes();
0071   }
0072   if (useStripStripQualityDB_) {
0073     stripQualityToken_ = esConsumes();
0074   }
0075 }
0076 
0077 void MkFitEventOfHitsProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0078   edm::ParameterSetDescription desc;
0079 
0080   desc.add("beamSpot", edm::InputTag{"offlineBeamSpot"});
0081   desc.add("pixelHits", edm::InputTag{"mkFitSiPixelHits"});
0082   desc.add("stripHits", edm::InputTag{"mkFitSiStripHits"});
0083   desc.add("usePixelQualityDB", true)->setComment("Use SiPixelQuality DB information");
0084   desc.add("useStripStripQualityDB", true)->setComment("Use SiStrip quality DB information");
0085 
0086   descriptions.addWithDefaultLabel(desc);
0087 }
0088 
0089 void MkFitEventOfHitsProducer::produce(edm::StreamID iID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0090   const auto& pixelHits = iEvent.get(pixelHitsToken_);
0091   const auto& stripHits = iEvent.get(stripHitsToken_);
0092   const auto& mkFitGeom = iSetup.getData(mkFitGeomToken_);
0093 
0094   auto eventOfHits = std::make_unique<mkfit::EventOfHits>(mkFitGeom.trackerInfo());
0095   mkfit::StdSeq::cmssw_LoadHits_Begin(*eventOfHits, {&pixelHits.hits(), &stripHits.hits()});
0096 
0097   if (usePixelQualityDB_ || useStripStripQualityDB_) {
0098     std::vector<mkfit::DeadVec> deadvectors(mkFitGeom.layerNumberConverter().nLayers());
0099     const auto& trackerGeom = iSetup.getData(geomToken_);
0100 
0101     if (usePixelQualityDB_) {
0102       const auto& pixelQuality = iSetup.getData(pixelQualityToken_);
0103       const auto& badPixels = pixelQuality.getBadComponentList();
0104       for (const auto& bp : badPixels) {
0105         const DetId detid(bp.DetID);
0106         const auto& surf = trackerGeom.idToDet(detid)->surface();
0107         bool isBarrel = (mkFitGeom.topology()->side(detid) == static_cast<unsigned>(TrackerDetSide::Barrel));
0108         const auto ilay = mkFitGeom.mkFitLayerNumber(detid);
0109         const auto q1 = isBarrel ? surf.zSpan().first : surf.rSpan().first;
0110         const auto q2 = isBarrel ? surf.zSpan().second : surf.rSpan().second;
0111         if (bp.errorType == 0)
0112           deadvectors[ilay].push_back({surf.phiSpan().first, surf.phiSpan().second, q1, q2});
0113       }
0114     }
0115     // For Phase-2, disable (momentarily?) SiStrip quality check
0116     if (useStripStripQualityDB_) {
0117       const auto& siStripQuality = iSetup.getData(stripQualityToken_);
0118       const auto& badStrips = siStripQuality.getBadComponentList();
0119       for (const auto& bs : badStrips) {
0120         const DetId detid(bs.detid);
0121         const auto& surf = trackerGeom.idToDet(detid)->surface();
0122         bool isBarrel = (mkFitGeom.topology()->side(detid) == static_cast<unsigned>(TrackerDetSide::Barrel));
0123         const auto ilay = mkFitGeom.mkFitLayerNumber(detid);
0124         const auto q1 = isBarrel ? surf.zSpan().first : surf.rSpan().first;
0125         const auto q2 = isBarrel ? surf.zSpan().second : surf.rSpan().second;
0126         if (bs.BadModule)
0127           deadvectors[ilay].push_back({surf.phiSpan().first, surf.phiSpan().second, q1, q2});
0128         else {  //assume that BadApvs are filled in sync with BadFibers
0129           auto const& topo = dynamic_cast<const StripTopology&>(trackerGeom.idToDet(detid)->topology());
0130           int firstApv = -1;
0131           int lastApv = -1;
0132 
0133           auto addRangeAPV = [&topo, &surf, &q1, &q2](int first, int last, mkfit::DeadVec& dv) {
0134             auto const firstPoint = surf.toGlobal(topo.localPosition(first * sistrip::STRIPS_PER_APV));
0135             auto const lastPoint = surf.toGlobal(topo.localPosition((last + 1) * sistrip::STRIPS_PER_APV));
0136             float phi1 = firstPoint.phi();
0137             float phi2 = lastPoint.phi();
0138             if (reco::deltaPhi(phi1, phi2) > 0)
0139               std::swap(phi1, phi2);
0140             LogTrace("SiStripBadComponents")
0141                 << "insert bad range " << first << " to " << last << " " << phi1 << " " << phi2;
0142             dv.push_back({phi1, phi2, q1, q2});
0143           };
0144 
0145           const int nApvs = topo.nstrips() / sistrip::STRIPS_PER_APV;
0146           for (int apv = 0; apv < nApvs; ++apv) {
0147             const bool isBad = bs.BadApvs & (1 << apv);
0148             if (isBad) {
0149               LogTrace("SiStripBadComponents") << "bad apv " << apv << " on " << bs.detid;
0150               if (lastApv == -1) {
0151                 firstApv = apv;
0152                 lastApv = apv;
0153               } else if (lastApv + 1 == apv)
0154                 lastApv++;
0155 
0156               if (apv + 1 == nApvs)
0157                 addRangeAPV(firstApv, lastApv, deadvectors[ilay]);
0158             } else if (firstApv != -1) {
0159               addRangeAPV(firstApv, lastApv, deadvectors[ilay]);
0160               //and reset
0161               firstApv = -1;
0162               lastApv = -1;
0163             }
0164           }
0165         }
0166       }
0167     }
0168     mkfit::StdSeq::loadDeads(*eventOfHits, deadvectors);
0169   }
0170 
0171   fill(iEvent.get(pixelLayerIndexToHitToken_), *eventOfHits, mkFitGeom);
0172   fill(iEvent.get(stripLayerIndexToHitToken_), *eventOfHits, mkFitGeom);
0173 
0174   mkfit::StdSeq::cmssw_LoadHits_End(*eventOfHits);
0175 
0176   auto const bs = iEvent.get(beamSpotToken_);
0177   eventOfHits->setBeamSpot(
0178       mkfit::BeamSpot(bs.x0(), bs.y0(), bs.z0(), bs.sigmaZ(), bs.BeamWidthX(), bs.BeamWidthY(), bs.dxdz(), bs.dydz()));
0179 
0180   iEvent.emplace(putToken_, std::move(eventOfHits));
0181 }
0182 
0183 void MkFitEventOfHitsProducer::fill(const std::vector<int>& layerToHits,
0184                                     mkfit::EventOfHits& eventOfHits,
0185                                     const MkFitGeometry& mkFitGeom) const {
0186   for (unsigned int i = 0, end = layerToHits.size(); i < end; ++i) {
0187     const auto ilay = layerToHits[i];
0188     if (ilay >= 0) {
0189       eventOfHits[ilay].registerHit(i);
0190     }
0191   }
0192 }
0193 
0194 DEFINE_FWK_MODULE(MkFitEventOfHitsProducer);