Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:25:11

0001 #include "RecoLocalTracker/SiStripClusterizer/interface/ThreeThresholdAlgorithm.h"
0002 #include "DataFormats/SiStripDigi/interface/SiStripDigi.h"
0003 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0004 #include <cmath>
0005 #include <numeric>
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 
0008 #include "DataFormats/SiStripCluster/interface/SiStripClusterTools.h"
0009 
0010 ThreeThresholdAlgorithm::ThreeThresholdAlgorithm(
0011     const edm::ESGetToken<SiStripClusterizerConditions, SiStripClusterizerConditionsRcd>& conditionsToken,
0012     float chan,
0013     float seed,
0014     float cluster,
0015     unsigned holes,
0016     unsigned bad,
0017     unsigned adj,
0018     bool removeApvShots,
0019     float minGoodCharge)
0020     : StripClusterizerAlgorithm(conditionsToken),
0021       ChannelThreshold(chan),
0022       SeedThreshold(seed),
0023       ClusterThresholdSquared(cluster * cluster),
0024       MaxSequentialHoles(holes),
0025       MaxSequentialBad(bad),
0026       MaxAdjacentBad(adj),
0027       RemoveApvShots(removeApvShots),
0028       minGoodCharge(minGoodCharge) {}
0029 
0030 template <class digiDetSet>
0031 inline void ThreeThresholdAlgorithm::clusterizeDetUnit_(const digiDetSet& digis, output_t::TSFastFiller& output) const {
0032   const auto& cond = conditions();
0033   if (cond.isModuleBad(digis.detId()))
0034     return;
0035 
0036   auto const& det = cond.findDetId(digis.detId());
0037   if (!det.valid())
0038     return;
0039 
0040 #ifdef EDM_ML_DEBUG
0041   if (!cond.isModuleUsable(digis.detId()))
0042     edm::LogWarning("ThreeThresholdAlgorithm") << " id " << digis.detId() << " not usable???" << std::endl;
0043 #endif
0044 
0045   typename digiDetSet::const_iterator scan(digis.begin()), end(digis.end());
0046 
0047   SiStripApvShotCleaner ApvCleaner;
0048   if (RemoveApvShots) {
0049     ApvCleaner.clean(digis, scan, end);
0050   }
0051 
0052   output.reserve(16);
0053   State state(det);
0054   while (scan != end) {
0055     while (scan != end && !candidateEnded(state, scan->strip()))
0056       addToCandidate(state, *scan++);
0057     endCandidate(state, output);
0058   }
0059 }
0060 
0061 inline bool ThreeThresholdAlgorithm::candidateEnded(State const& state, const uint16_t& testStrip) const {
0062   uint16_t holes = testStrip - state.lastStrip - 1;
0063   return (((!state.ADCs.empty()) &       // a candidate exists, and
0064            (holes > MaxSequentialHoles)  // too many holes if not all are bad strips, and
0065            ) &&
0066           (holes > MaxSequentialBad ||                             // (too many bad strips anyway, or
0067            !state.det().allBadBetween(state.lastStrip, testStrip)  // not all holes are bad strips)
0068            ));
0069 }
0070 
0071 inline void ThreeThresholdAlgorithm::addToCandidate(State& state, uint16_t strip, uint8_t adc) const {
0072   float Noise = state.det().noise(strip);
0073   if (adc < static_cast<uint8_t>(Noise * ChannelThreshold) || state.det().bad(strip))
0074     return;
0075 
0076   if (state.candidateLacksSeed)
0077     state.candidateLacksSeed = adc < static_cast<uint8_t>(Noise * SeedThreshold);
0078   if (state.ADCs.empty())
0079     state.lastStrip = strip - 1;  // begin candidate
0080   while (++state.lastStrip < strip)
0081     state.ADCs.push_back(0);  // pad holes
0082 
0083   state.ADCs.push_back(adc);
0084   state.noiseSquared += Noise * Noise;
0085 }
0086 
0087 template <class T>
0088 inline void ThreeThresholdAlgorithm::endCandidate(State& state, T& out) const {
0089   if (candidateAccepted(state)) {
0090     applyGains(state);
0091     if (MaxAdjacentBad > 0)
0092       appendBadNeighbors(state);
0093     if (minGoodCharge <= 0 ||
0094         siStripClusterTools::chargePerCM(state.det().detId, state.ADCs.begin(), state.ADCs.end()) > minGoodCharge)
0095       out.push_back(std::move(SiStripCluster(firstStrip(state), state.ADCs.begin(), state.ADCs.end())));
0096   }
0097   clearCandidate(state);
0098 }
0099 
0100 inline bool ThreeThresholdAlgorithm::candidateAccepted(State const& state) const {
0101   return (!state.candidateLacksSeed &&
0102           state.noiseSquared * ClusterThresholdSquared <=
0103               std::pow(float(std::accumulate(state.ADCs.begin(), state.ADCs.end(), int(0))), 2.f));
0104 }
0105 
0106 inline void ThreeThresholdAlgorithm::applyGains(State& state) const {
0107   uint16_t strip = firstStrip(state);
0108   for (auto& adc : state.ADCs) {
0109 #ifdef EDM_ML_DEBUG
0110     // if(adc > 255) throw InvalidChargeException( SiStripDigi(strip,adc) );
0111 #endif
0112     // if(adc > 253) continue; //saturated, do not scale
0113     auto charge = int(float(adc) * state.det().weight(strip++) + 0.5f);  //adding 0.5 turns truncation into rounding
0114     if (adc < 254)
0115       adc = (charge > 1022 ? 255 : (charge > 253 ? 254 : charge));
0116   }
0117 }
0118 
0119 inline void ThreeThresholdAlgorithm::appendBadNeighbors(State& state) const {
0120   uint8_t max = MaxAdjacentBad;
0121   while (0 < max--) {
0122     if (state.det().bad(firstStrip(state) - 1)) {
0123       state.ADCs.insert(state.ADCs.begin(), 0);
0124     }
0125     if (state.det().bad(state.lastStrip + 1)) {
0126       state.ADCs.push_back(0);
0127       state.lastStrip++;
0128     }
0129   }
0130 }
0131 
0132 void ThreeThresholdAlgorithm::clusterizeDetUnit(const edm::DetSet<SiStripDigi>& digis,
0133                                                 output_t::TSFastFiller& output) const {
0134   clusterizeDetUnit_(digis, output);
0135 }
0136 void ThreeThresholdAlgorithm::clusterizeDetUnit(const edmNew::DetSet<SiStripDigi>& digis,
0137                                                 output_t::TSFastFiller& output) const {
0138   clusterizeDetUnit_(digis, output);
0139 }
0140 
0141 void ThreeThresholdAlgorithm::stripByStripAdd(State& state,
0142                                               uint16_t strip,
0143                                               uint8_t adc,
0144                                               std::vector<SiStripCluster>& out) const {
0145   if (candidateEnded(state, strip))
0146     endCandidate(state, out);
0147   addToCandidate(state, SiStripDigi(strip, adc));
0148 }
0149 
0150 void ThreeThresholdAlgorithm::stripByStripEnd(State& state, std::vector<SiStripCluster>& out) const {
0151   endCandidate(state, out);
0152 }