File indexing completed on 2024-04-06 12:26:29
0001 #include "RecoLocalTracker/SiStripRecHitConverter/interface/SiStripRecHitConverterAlgorithm.h"
0002 #include "RecoLocalTracker/Records/interface/TkStripCPERecord.h"
0003 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0004
0005 #include "Geometry/CommonDetUnit/interface/GluedGeomDet.h"
0006 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0007
0008 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0009 #include "DataFormats/Common/interface/Ref.h"
0010
0011 #include "FWCore/Framework/interface/ConsumesCollector.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include <vector>
0017
0018 SiStripRecHitConverterAlgorithm::SiStripRecHitConverterAlgorithm(const edm::ParameterSet& conf,
0019 edm::ConsumesCollector iC)
0020 : useQuality(conf.getParameter<bool>("useSiStripQuality")),
0021 maskBad128StripBlocks(conf.getParameter<bool>("MaskBadAPVFibers")),
0022 doMatching(conf.getParameter<bool>("doMatching")),
0023 trackerToken(iC.esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0024 cpeToken(iC.esConsumes<StripClusterParameterEstimator, TkStripCPERecord>(
0025 conf.getParameter<edm::ESInputTag>("StripCPE"))) {
0026 if (doMatching) {
0027 matcherToken = iC.esConsumes<SiStripRecHitMatcher, TkStripCPERecord>(conf.getParameter<edm::ESInputTag>("Matcher"));
0028 }
0029 if (useQuality) {
0030 qualityToken =
0031 iC.esConsumes<SiStripQuality, SiStripQualityRcd>(conf.getParameter<edm::ESInputTag>("siStripQualityLabel"));
0032 }
0033 }
0034
0035 void SiStripRecHitConverterAlgorithm::fillPSetDescription(edm::ParameterSetDescription& desc) {
0036 desc.add<bool>("useSiStripQuality", false);
0037 desc.add<bool>("MaskBadAPVFibers", false);
0038 desc.add<bool>("doMatching", true);
0039 desc.add<edm::ESInputTag>("StripCPE", edm::ESInputTag("StripCPEfromTrackAngleESProducer", "StripCPEfromTrackAngle"));
0040 desc.add<edm::ESInputTag>("Matcher", edm::ESInputTag("SiStripRecHitMatcherESProducer", "StandardMatcher"));
0041 desc.add<edm::ESInputTag>("siStripQualityLabel", edm::ESInputTag());
0042 }
0043
0044 void SiStripRecHitConverterAlgorithm::initialize(const edm::EventSetup& es) {
0045 tracker = &es.getData(trackerToken);
0046 parameterestimator = &es.getData(cpeToken);
0047 if (doMatching) {
0048 matcher = &es.getData(matcherToken);
0049 }
0050 if (useQuality) {
0051 quality = &es.getData(qualityToken);
0052 }
0053 }
0054
0055 void SiStripRecHitConverterAlgorithm::run(edm::Handle<edmNew::DetSetVector<SiStripCluster>> input, products& output) {
0056 run(input, output, LocalVector(0., 0., 0.));
0057 }
0058
0059 void SiStripRecHitConverterAlgorithm::run(edm::Handle<edmNew::DetSetVector<SiStripCluster>> inputhandle,
0060 products& output,
0061 LocalVector trackdirection) {
0062 for (auto const& DS : *inputhandle) {
0063 auto id = DS.id();
0064 if (!useModule(id))
0065 continue;
0066
0067 Collector collector = StripSubdetector(id).stereo() ? Collector(*output.stereo, id) : Collector(*output.rphi, id);
0068
0069 bool bad128StripBlocks[6];
0070 fillBad128StripBlocks(id, bad128StripBlocks);
0071
0072 GeomDetUnit const& du = *(tracker->idToDetUnit(id));
0073 for (auto const& cluster : DS) {
0074 if (isMasked(cluster, bad128StripBlocks))
0075 continue;
0076
0077 StripClusterParameterEstimator::LocalValues parameters = parameterestimator->localParameters(cluster, du);
0078 collector.push_back(
0079 SiStripRecHit2D(parameters.first, parameters.second, du, DS.makeRefTo(inputhandle, &cluster)));
0080 }
0081
0082 if (collector.empty())
0083 collector.abort();
0084 }
0085 if (doMatching) {
0086 match(output, trackdirection);
0087 }
0088 }
0089
0090 namespace {
0091
0092 struct CollectorHelper {
0093 size_t nmatch;
0094
0095 using CollectorMatched = std::vector<std::unique_ptr<SiStripMatchedRecHit2D>>;
0096 typedef SiStripMatchedRecHit2DCollection::FastFiller Collector;
0097
0098 Collector& m_collector;
0099 CollectorMatched& m_collectorMatched;
0100 SiStripRecHit2DCollection::FastFiller& m_fillerRphiUnm;
0101 std::vector<SiStripRecHit2D::ClusterRef::key_type>& m_matchedSteroClusters;
0102
0103 static inline SiStripRecHit2D const& stereoHit(edmNew::DetSet<SiStripRecHit2D>::const_iterator iter) {
0104 return *iter;
0105 }
0106
0107 static inline SiStripRecHit2D const& monoHit(edmNew::DetSet<SiStripRecHit2D>::const_iterator iter) { return *iter; }
0108
0109 CollectorHelper& collector() { return *this; }
0110
0111 void operator()(SiStripMatchedRecHit2D const& rh) { m_collectorMatched.emplace_back(rh.clone()); }
0112
0113 CollectorHelper(Collector& i_collector,
0114 CollectorMatched& i_collectorMatched,
0115 SiStripRecHit2DCollection::FastFiller& i_fillerRphiUnm,
0116 std::vector<SiStripRecHit2D::ClusterRef::key_type>& i_matchedSteroClusters)
0117 : nmatch(0),
0118 m_collector(i_collector),
0119 m_collectorMatched(i_collectorMatched),
0120 m_fillerRphiUnm(i_fillerRphiUnm),
0121 m_matchedSteroClusters(i_matchedSteroClusters) {}
0122
0123 void closure(edmNew::DetSet<SiStripRecHit2D>::const_iterator it) {
0124 if (!m_collectorMatched.empty()) {
0125 nmatch += m_collectorMatched.size();
0126 for (auto const& itm : m_collectorMatched) {
0127 m_collector.push_back(*itm);
0128
0129 m_matchedSteroClusters.push_back(itm->stereoClusterRef().key());
0130 }
0131 m_collectorMatched.clear();
0132 } else {
0133
0134 m_fillerRphiUnm.push_back(*it);
0135 }
0136 }
0137 };
0138 }
0139
0140 void SiStripRecHitConverterAlgorithm::match(products& output, LocalVector trackdirection) const {
0141 int nmatch = 0;
0142 std::vector<std::unique_ptr<SiStripMatchedRecHit2D>> collectorMatched;
0143
0144
0145 SiStripRecHit2DCollection::const_iterator edStereoDet = output.stereo->end();
0146 SiStripRecHit2DCollection::const_iterator edRPhiDet = output.rphi->end();
0147
0148
0149 std::vector<SiStripRecHit2D::ClusterRef::key_type> matchedSteroClusters;
0150
0151 for (SiStripRecHit2DCollection::const_iterator itRPhiDet = output.rphi->begin(); itRPhiDet != edRPhiDet;
0152 ++itRPhiDet) {
0153 edmNew::DetSet<SiStripRecHit2D> rphiHits = *itRPhiDet;
0154 StripSubdetector specDetId(rphiHits.detId());
0155 uint32_t partnerId = specDetId.partnerDetId();
0156
0157
0158 if (partnerId == 0) {
0159
0160 if (!rphiHits.empty()) {
0161 SiStripRecHit2DCollection::FastFiller filler(*output.rphiUnmatched, rphiHits.detId());
0162 filler.resize(rphiHits.size());
0163 std::copy(rphiHits.begin(), rphiHits.end(), filler.begin());
0164 }
0165 continue;
0166 }
0167
0168 SiStripRecHit2DCollection::const_iterator itStereoDet = output.stereo->find(partnerId);
0169
0170
0171 if (itStereoDet == edStereoDet) {
0172
0173 if (!rphiHits.empty()) {
0174 SiStripRecHit2DCollection::FastFiller filler(*output.rphiUnmatched, rphiHits.detId());
0175 filler.resize(rphiHits.size());
0176 std::copy(rphiHits.begin(), rphiHits.end(), filler.begin());
0177 }
0178 continue;
0179 }
0180
0181 edmNew::DetSet<SiStripRecHit2D> stereoHits = *itStereoDet;
0182
0183
0184 const GluedGeomDet* gluedDet = (const GluedGeomDet*)tracker->idToDet(DetId(specDetId.glued()));
0185 typedef SiStripMatchedRecHit2DCollection::FastFiller Collector;
0186 Collector collector(*output.matched, specDetId.glued());
0187
0188
0189 SiStripRecHit2DCollection::FastFiller fillerRphiUnm(*output.rphiUnmatched, rphiHits.detId());
0190
0191
0192 matchedSteroClusters.clear();
0193
0194 #ifdef DOUBLE_MATCH
0195 CollectorHelper chelper(collector, collectorMatched, fillerRphiUnm, matchedSteroClusters);
0196 matcher->doubleMatch(
0197 rphiHits.begin(), rphiHits.end(), stereoHits.begin(), stereoHits.end(), gluedDet, trackdirection, chelper);
0198 nmatch += chelper.nmatch;
0199 #else
0200
0201 SiStripRecHitMatcher::SimpleHitCollection stereoSimpleHits;
0202
0203 stereoSimpleHits.reserve(stereoHits.size());
0204 for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = stereoHits.begin(), ed = stereoHits.end(); it != ed;
0205 ++it) {
0206 stereoSimpleHits.push_back(&*it);
0207 }
0208
0209 for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = rphiHits.begin(), ed = rphiHits.end(); it != ed; ++it) {
0210 matcher->match(
0211 &(*it), stereoSimpleHits.begin(), stereoSimpleHits.end(), collectorMatched, gluedDet, trackdirection);
0212 if (collectorMatched.size() > 0) {
0213 nmatch += collectorMatched.size();
0214 for (auto const& itm : collectorMatched) {
0215 collector.push_back(*itm);
0216
0217 matchedSteroClusters.push_back(itm->stereoClusterRef().key());
0218 }
0219 collectorMatched.clear();
0220 } else {
0221
0222 fillerRphiUnm.push_back(*it);
0223 }
0224 }
0225
0226 #endif
0227
0228
0229 if (collector.empty())
0230 collector.abort();
0231
0232
0233 if (fillerRphiUnm.empty())
0234 fillerRphiUnm.abort();
0235
0236
0237 SiStripRecHit2DCollection::FastFiller fillerStereoUnm(*output.stereoUnmatched, stereoHits.detId());
0238 std::sort(matchedSteroClusters.begin(), matchedSteroClusters.end());
0239 for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = stereoHits.begin(), ed = stereoHits.end(); it != ed;
0240 ++it) {
0241 if (!std::binary_search(matchedSteroClusters.begin(), matchedSteroClusters.end(), it->cluster().key())) {
0242 fillerStereoUnm.push_back(*it);
0243 }
0244 }
0245 if (fillerStereoUnm.empty())
0246 fillerStereoUnm.abort();
0247 }
0248
0249 for (SiStripRecHit2DCollection::const_iterator itStereoDet = output.stereo->begin(); itStereoDet != edStereoDet;
0250 ++itStereoDet) {
0251 edmNew::DetSet<SiStripRecHit2D> stereoHits = *itStereoDet;
0252 StripSubdetector specDetId(stereoHits.detId());
0253 uint32_t partnerId = specDetId.partnerDetId();
0254 if (partnerId == 0)
0255 continue;
0256 SiStripRecHit2DCollection::const_iterator itRPhiDet = output.rphi->find(partnerId);
0257 if (itRPhiDet == edRPhiDet) {
0258 if (!stereoHits.empty()) {
0259 SiStripRecHit2DCollection::FastFiller filler(*output.stereoUnmatched, stereoHits.detId());
0260 filler.resize(stereoHits.size());
0261 std::copy(stereoHits.begin(), stereoHits.end(), filler.begin());
0262 }
0263 }
0264 }
0265
0266 edm::LogInfo("SiStripRecHitConverter") << "found\n" << nmatch << " matched RecHits\n";
0267 }
0268
0269 void SiStripRecHitConverterAlgorithm::fillBad128StripBlocks(const uint32_t detid, bool bad128StripBlocks[6]) const {
0270 if (maskBad128StripBlocks) {
0271 short badApvs = quality->getBadApvs(detid);
0272 short badFibers = quality->getBadFibers(detid);
0273 for (int j = 0; j < 6; j++) {
0274 bad128StripBlocks[j] = (badApvs & (1 << j));
0275 }
0276 for (int j = 0; j < 3; j++) {
0277 if (badFibers & (1 << j)) {
0278 bad128StripBlocks[2 * j + 0] = true;
0279 bad128StripBlocks[2 * j + 1] = true;
0280 }
0281 }
0282 }
0283 }
0284
0285 inline bool SiStripRecHitConverterAlgorithm::isMasked(const SiStripCluster& cluster, bool bad128StripBlocks[6]) const {
0286 if (maskBad128StripBlocks) {
0287 if (bad128StripBlocks[cluster.firstStrip() >> 7]) {
0288 if (bad128StripBlocks[(cluster.firstStrip() + cluster.amplitudes().size()) >> 7] ||
0289 bad128StripBlocks[static_cast<int32_t>(cluster.barycenter() - 0.499999) >> 7]) {
0290 return true;
0291 }
0292 } else {
0293 if (bad128StripBlocks[(cluster.firstStrip() + cluster.amplitudes().size()) >> 7] &&
0294 bad128StripBlocks[static_cast<int32_t>(cluster.barycenter() - 0.499999) >> 7]) {
0295 return true;
0296 }
0297 }
0298 }
0299 return false;
0300 }
0301
0302 inline bool SiStripRecHitConverterAlgorithm::useModule(const uint32_t id) const {
0303 const StripGeomDetUnit* stripdet = (const StripGeomDetUnit*)tracker->idToDetUnit(id);
0304 if (stripdet == nullptr)
0305 edm::LogWarning("SiStripRecHitConverter") << "Detid=" << id << " not found";
0306 return stripdet != nullptr && (!useQuality || quality->IsModuleUsable(id));
0307 }