SiStripRecHitConverter/src/SiStripRecHitConverterAlgorithm.cc

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 SiStripRecHitConverterAlgorithm::SiStripRecHitConverterAlgorithm(const edm::ParameterSet& conf,
0017                                                                  edm::ConsumesCollector iC)
0018     : useQuality(conf.getParameter<bool>("useSiStripQuality")),
0019       maskBad128StripBlocks(conf.getParameter<bool>("MaskBadAPVFibers")),
0020       doMatching(conf.getParameter<bool>("doMatching")),
0021       trackerToken(iC.esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0022       cpeToken(iC.esConsumes<StripClusterParameterEstimator, TkStripCPERecord>(
0023           conf.getParameter<edm::ESInputTag>("StripCPE"))) {
0024   if (doMatching) {
0025     matcherToken = iC.esConsumes<SiStripRecHitMatcher, TkStripCPERecord>(conf.getParameter<edm::ESInputTag>("Matcher"));
0026   }
0027   if (useQuality) {
0028     qualityToken =
0029         iC.esConsumes<SiStripQuality, SiStripQualityRcd>(conf.getParameter<edm::ESInputTag>("siStripQualityLabel"));
0030   }
0031 }
0032
0033 void SiStripRecHitConverterAlgorithm::fillPSetDescription(edm::ParameterSetDescription& desc) {
0034   desc.add<bool>("useSiStripQuality", false);
0035   desc.add<bool>("MaskBadAPVFibers", false);
0036   desc.add<bool>("doMatching", true);
0037   desc.add<edm::ESInputTag>("StripCPE", edm::ESInputTag("StripCPEfromTrackAngleESProducer", "StripCPEfromTrackAngle"));
0038   desc.add<edm::ESInputTag>("Matcher", edm::ESInputTag("SiStripRecHitMatcherESProducer", "StandardMatcher"));
0039   desc.add<edm::ESInputTag>("siStripQualityLabel", edm::ESInputTag());
0040 }
0041
0042 void SiStripRecHitConverterAlgorithm::initialize(const edm::EventSetup& es) {
0043   tracker = &es.getData(trackerToken);
0044   parameterestimator = &es.getData(cpeToken);
0045   if (doMatching) {
0046     matcher = &es.getData(matcherToken);
0047   }
0048   if (useQuality) {
0049     quality = &es.getData(qualityToken);
0050   }
0051 }
0052
0053 void SiStripRecHitConverterAlgorithm::run(edm::Handle<edmNew::DetSetVector<SiStripCluster> > input, products& output) {
0054   run(input, output, LocalVector(0., 0., 0.));
0055 }
0056
0057 void SiStripRecHitConverterAlgorithm::run(edm::Handle<edmNew::DetSetVector<SiStripCluster> > inputhandle,
0058                                           products& output,
0059                                           LocalVector trackdirection) {
0060   for (auto const& DS : *inputhandle) {
0061     auto id = DS.id();
0062     if (!useModule(id))
0063       continue;
0064
0065     Collector collector = StripSubdetector(id).stereo() ? Collector(*output.stereo, id) : Collector(*output.rphi, id);
0066
0067     bool bad128StripBlocks[6];
0068     fillBad128StripBlocks(id, bad128StripBlocks);
0069
0070     GeomDetUnit const& du = *(tracker->idToDetUnit(id));
0071     for (auto const& cluster : DS) {
0072       if (isMasked(cluster, bad128StripBlocks))
0073         continue;
0074
0075       StripClusterParameterEstimator::LocalValues parameters = parameterestimator->localParameters(cluster, du);
0076       collector.push_back(
0077           SiStripRecHit2D(parameters.first, parameters.second, du, DS.makeRefTo(inputhandle, &cluster)));
0078     }
0079
0080     if (collector.empty())
0081       collector.abort();
0082   }
0083   if (doMatching) {
0084     match(output, trackdirection);
0085   }
0086 }
0087
0088 namespace {
0089
0090   struct CollectorHelper {
0091     size_t nmatch;
0092
0093     typedef edm::OwnVector<SiStripMatchedRecHit2D> CollectorMatched;
0094     typedef SiStripMatchedRecHit2DCollection::FastFiller Collector;
0095
0096     Collector& m_collector;
0097     CollectorMatched& m_collectorMatched;
0098     SiStripRecHit2DCollection::FastFiller& m_fillerRphiUnm;
0099     std::vector<SiStripRecHit2D::ClusterRef::key_type>& m_matchedSteroClusters;
0100
0101     static inline SiStripRecHit2D const& stereoHit(edmNew::DetSet<SiStripRecHit2D>::const_iterator iter) {
0102       return *iter;
0103     }
0104
0105     static inline SiStripRecHit2D const& monoHit(edmNew::DetSet<SiStripRecHit2D>::const_iterator iter) { return *iter; }
0106
0107     struct Add {
0108       Add(CollectorHelper& ih) : h(ih) {}
0109       CollectorHelper& h;
0110       void operator()(SiStripMatchedRecHit2D const& rh) { h.m_collectorMatched.push_back(rh); }
0111     };
0112
0113     CollectorHelper& collector() { return *this; }
0114
0115     void operator()(SiStripMatchedRecHit2D const& rh) { m_collectorMatched.push_back(rh); }
0116
0117     CollectorHelper(Collector& i_collector,
0118                     CollectorMatched& i_collectorMatched,
0119                     SiStripRecHit2DCollection::FastFiller& i_fillerRphiUnm,
0120                     std::vector<SiStripRecHit2D::ClusterRef::key_type>& i_matchedSteroClusters)
0121         : nmatch(0),
0122           m_collector(i_collector),
0123           m_collectorMatched(i_collectorMatched),
0124           m_fillerRphiUnm(i_fillerRphiUnm),
0125           m_matchedSteroClusters(i_matchedSteroClusters) {}
0126
0127     void closure(edmNew::DetSet<SiStripRecHit2D>::const_iterator it) {
0128       if (!m_collectorMatched.empty()) {
0129         nmatch += m_collectorMatched.size();
0130         for (edm::OwnVector<SiStripMatchedRecHit2D>::const_iterator itm = m_collectorMatched.begin(),
0131                                                                     edm = m_collectorMatched.end();
0132              itm != edm;
0133              ++itm) {
0134           m_collector.push_back(*itm);
0135           // mark the stereo hit cluster as used, so that the hit won't go in the unmatched stereo ones
0136           m_matchedSteroClusters.push_back(itm->stereoClusterRef().key());
0137         }
0138         m_collectorMatched.clear();
0139       } else {
0140         // store a copy of this rphi hit as an unmatched rphi hit
0141         m_fillerRphiUnm.push_back(*it);
0142       }
0143     }
0144   };
0145 }  // namespace
0146
0147 void SiStripRecHitConverterAlgorithm::match(products& output, LocalVector trackdirection) const {
0148   int nmatch = 0;
0149   edm::OwnVector<SiStripMatchedRecHit2D> collectorMatched;  // gp/FIXME: avoid this
0150
0151   // Remember the ends of the collections, as we will use them a lot
0152   SiStripRecHit2DCollection::const_iterator edStereoDet = output.stereo->end();
0153   SiStripRecHit2DCollection::const_iterator edRPhiDet = output.rphi->end();
0154
0155   // two work vectors for bookeeping clusters used by the stereo part of the matched hits
0156   std::vector<SiStripRecHit2D::ClusterRef::key_type> matchedSteroClusters;
0157
0158   for (SiStripRecHit2DCollection::const_iterator itRPhiDet = output.rphi->begin(); itRPhiDet != edRPhiDet;
0159        ++itRPhiDet) {
0160     edmNew::DetSet<SiStripRecHit2D> rphiHits = *itRPhiDet;
0161     StripSubdetector specDetId(rphiHits.detId());
0162     uint32_t partnerId = specDetId.partnerDetId();
0163
0164     // if not part of a glued pair
0165     if (partnerId == 0) {
0166       // I must copy these as unmatched
0167       if (!rphiHits.empty()) {
0168         SiStripRecHit2DCollection::FastFiller filler(*output.rphiUnmatched, rphiHits.detId());
0169         filler.resize(rphiHits.size());
0170         std::copy(rphiHits.begin(), rphiHits.end(), filler.begin());
0171       }
0172       continue;
0173     }
0174
0175     SiStripRecHit2DCollection::const_iterator itStereoDet = output.stereo->find(partnerId);
0176
0177     // if the partner is not found (which probably can happen if it's empty)
0178     if (itStereoDet == edStereoDet) {
0179       // I must copy these as unmatched
0180       if (!rphiHits.empty()) {
0181         SiStripRecHit2DCollection::FastFiller filler(*output.rphiUnmatched, rphiHits.detId());
0182         filler.resize(rphiHits.size());
0183         std::copy(rphiHits.begin(), rphiHits.end(), filler.begin());
0184       }
0185       continue;
0186     }
0187
0188     edmNew::DetSet<SiStripRecHit2D> stereoHits = *itStereoDet;
0189
0190     // Get ready for making glued hits
0191     const GluedGeomDet* gluedDet = (const GluedGeomDet*)tracker->idToDet(DetId(specDetId.glued()));
0192     typedef SiStripMatchedRecHit2DCollection::FastFiller Collector;
0193     Collector collector(*output.matched, specDetId.glued());
0194
0195     // Prepare also the list for unmatched rphi hits
0196     SiStripRecHit2DCollection::FastFiller fillerRphiUnm(*output.rphiUnmatched, rphiHits.detId());
0197
0198     // a list of clusters used by the matched part of the stereo hits in this detector
0199     matchedSteroClusters.clear();  // at the beginning, empty
0200
0201 #ifdef DOUBLE_MATCH
0202     CollectorHelper chelper(collector, collectorMatched, fillerRphiUnm, matchedSteroClusters);
0203     matcher->doubleMatch(
0204         rphiHits.begin(), rphiHits.end(), stereoHits.begin(), stereoHits.end(), gluedDet, trackdirection, chelper);
0205     nmatch += chelper.nmatch;
0206 #else
0207     // Make simple collection of this (gp:FIXME: why do we need it?)
0208     SiStripRecHitMatcher::SimpleHitCollection stereoSimpleHits;
0209     // gp:FIXME: use std::transform
0210     stereoSimpleHits.reserve(stereoHits.size());
0211     for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = stereoHits.begin(), ed = stereoHits.end(); it != ed;
0212          ++it) {
0213       stereoSimpleHits.push_back(&*it);
0214     }
0215
0216     for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = rphiHits.begin(), ed = rphiHits.end(); it != ed; ++it) {
0217       matcher->match(
0218           &(*it), stereoSimpleHits.begin(), stereoSimpleHits.end(), collectorMatched, gluedDet, trackdirection);
0219       if (collectorMatched.size() > 0) {
0220         nmatch += collectorMatched.size();
0221         for (edm::OwnVector<SiStripMatchedRecHit2D>::const_iterator itm = collectorMatched.begin(),
0222                                                                     edm = collectorMatched.end();
0223              itm != edm;
0224              ++itm) {
0225           collector.push_back(*itm);
0226           // mark the stereo hit cluster as used, so that the hit won't go in the unmatched stereo ones
0227           matchedSteroClusters.push_back(itm->stereoClusterRef().key());
0228         }
0229         collectorMatched.clear();
0230       } else {
0231         // store a copy of this rphi hit as an unmatched rphi hit
0232         fillerRphiUnm.push_back(*it);
0233       }
0234     }
0235
0236 #endif
0237
0238     // discard matched hits if the collection is empty
0239     if (collector.empty())
0240       collector.abort();
0241
0242     // discard unmatched rphi hits if there are none
0243     if (fillerRphiUnm.empty())
0244       fillerRphiUnm.abort();
0245
0246     // now look for unmatched stereo hits
0247     SiStripRecHit2DCollection::FastFiller fillerStereoUnm(*output.stereoUnmatched, stereoHits.detId());
0248     std::sort(matchedSteroClusters.begin(), matchedSteroClusters.end());
0249     for (edmNew::DetSet<SiStripRecHit2D>::const_iterator it = stereoHits.begin(), ed = stereoHits.end(); it != ed;
0250          ++it) {
0251       if (!std::binary_search(matchedSteroClusters.begin(), matchedSteroClusters.end(), it->cluster().key())) {
0252         fillerStereoUnm.push_back(*it);
0253       }
0254     }
0255     if (fillerStereoUnm.empty())
0256       fillerStereoUnm.abort();
0257   }
0258
0259   for (SiStripRecHit2DCollection::const_iterator itStereoDet = output.stereo->begin(); itStereoDet != edStereoDet;
0260        ++itStereoDet) {
0261     edmNew::DetSet<SiStripRecHit2D> stereoHits = *itStereoDet;
0262     StripSubdetector specDetId(stereoHits.detId());
0263     uint32_t partnerId = specDetId.partnerDetId();
0264     if (partnerId == 0)
0265       continue;
0266     SiStripRecHit2DCollection::const_iterator itRPhiDet = output.rphi->find(partnerId);
0267     if (itRPhiDet == edRPhiDet) {
0268       if (!stereoHits.empty()) {
0269         SiStripRecHit2DCollection::FastFiller filler(*output.stereoUnmatched, stereoHits.detId());
0270         filler.resize(stereoHits.size());
0271         std::copy(stereoHits.begin(), stereoHits.end(), filler.begin());
0272       }
0273     }
0274   }
0275
0276   edm::LogInfo("SiStripRecHitConverter") << "found\n" << nmatch << "  matched RecHits\n";
0277 }
0278
0279 void SiStripRecHitConverterAlgorithm::fillBad128StripBlocks(const uint32_t detid, bool bad128StripBlocks[6]) const {
0280   if (maskBad128StripBlocks) {
0281     short badApvs = quality->getBadApvs(detid);
0282     short badFibers = quality->getBadFibers(detid);
0283     for (int j = 0; j < 6; j++) {
0284       bad128StripBlocks[j] = (badApvs & (1 << j));
0285     }
0286     for (int j = 0; j < 3; j++) {
0287       if (badFibers & (1 << j)) {
0288         bad128StripBlocks[2 * j + 0] = true;
0289         bad128StripBlocks[2 * j + 1] = true;
0290       }
0291     }
0292   }
0293 }
0294
0295 inline bool SiStripRecHitConverterAlgorithm::isMasked(const SiStripCluster& cluster, bool bad128StripBlocks[6]) const {
0296   if (maskBad128StripBlocks) {
0297     if (bad128StripBlocks[cluster.firstStrip() >> 7]) {
0298       if (bad128StripBlocks[(cluster.firstStrip() + cluster.amplitudes().size()) >> 7] ||
0299           bad128StripBlocks[static_cast<int32_t>(cluster.barycenter() - 0.499999) >> 7]) {
0300         return true;
0301       }
0302     } else {
0303       if (bad128StripBlocks[(cluster.firstStrip() + cluster.amplitudes().size()) >> 7] &&
0304           bad128StripBlocks[static_cast<int32_t>(cluster.barycenter() - 0.499999) >> 7]) {
0305         return true;
0306       }
0307     }
0308   }
0309   return false;
0310 }
0311
0312 inline bool SiStripRecHitConverterAlgorithm::useModule(const uint32_t id) const {
0313   const StripGeomDetUnit* stripdet = (const StripGeomDetUnit*)tracker->idToDetUnit(id);
0314   if (stripdet == nullptr)
0315     edm::LogWarning("SiStripRecHitConverter") << "Detid=" << id << " not found";
0316   return stripdet != nullptr && (!useQuality || quality->IsModuleUsable(id));
0317 }