Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-28 03:05:39

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