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