File indexing completed on 2025-01-09 23:33:50
0001 #ifndef RecoLocalTracker_SiStripClusterizer_ThreeThresholdAlgorithm_h
0002 #define RecoLocalTracker_SiStripClusterizer_ThreeThresholdAlgorithm_h
0003
0004 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0005 #include "DataFormats/SiStripCluster/interface/SiStripClusterTools.h"
0006 #include "DataFormats/SiStripDigi/interface/SiStripDigi.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithm.h"
0009 #include "RecoLocalTracker/SiStripClusterizer/interface/SiStripApvShotCleaner.h"
0010
0011 #include <cmath>
0012 #include <numeric>
0013
0014 class ThreeThresholdAlgorithm final : public StripClusterizerAlgorithm {
0015 friend class StripClusterizerAlgorithmFactory;
0016
0017 public:
0018 using State = StripClusterizerAlgorithm::State;
0019 using Det = StripClusterizerAlgorithm::Det;
0020 void clusterizeDetUnit(const edm::DetSet<SiStripDigi>&, output_t::TSFastFiller&) const override;
0021 void clusterizeDetUnit(const edmNew::DetSet<SiStripDigi>&, output_t::TSFastFiller&) const override;
0022
0023
0024 void stripByStripAdd(State& state, uint16_t strip, uint8_t adc, std::vector<SiStripCluster>& out) const override;
0025 void stripByStripEnd(State& state, std::vector<SiStripCluster>& out) const override;
0026
0027 void stripByStripAdd(State& state, uint16_t strip, uint8_t adc, output_t::TSFastFiller& out) const override {
0028 if (candidateEnded(state, strip))
0029 endCandidate(state, out);
0030 addToCandidate(state, strip, adc);
0031 }
0032
0033 void stripByStripEnd(State& state, output_t::TSFastFiller& out) const override { endCandidate(state, out); }
0034
0035 private:
0036 template <class T>
0037 void clusterizeDetUnit_(const T&, output_t::TSFastFiller&) const;
0038
0039 ThreeThresholdAlgorithm(const edm::ESGetToken<SiStripClusterizerConditions, SiStripClusterizerConditionsRcd>&,
0040 float,
0041 float,
0042 float,
0043 unsigned,
0044 unsigned,
0045 unsigned,
0046 unsigned,
0047 bool removeApvShots,
0048 float minGoodCharge);
0049
0050
0051 uint16_t firstStrip(State const& state) const { return state.lastStrip - state.ADCs.size() + 1; }
0052 bool candidateEnded(State const& state, const uint16_t&) const;
0053 bool candidateAccepted(State const& state) const;
0054
0055
0056 template <class T>
0057 void endCandidate(State& state, T&) const;
0058 void clearCandidate(State& state) const {
0059 state.candidateLacksSeed = true;
0060 state.noiseSquared = 0;
0061 state.ADCs.clear();
0062 }
0063 void addToCandidate(State& state, const SiStripDigi& digi) const { addToCandidate(state, digi.strip(), digi.adc()); }
0064 void addToCandidate(State& state, uint16_t strip, uint8_t adc) const;
0065 void appendBadNeighbors(State& state) const;
0066 void applyGains(State& state) const;
0067
0068 float ChannelThreshold, SeedThreshold, ClusterThresholdSquared;
0069 uint8_t MaxSequentialHoles, MaxSequentialBad, MaxAdjacentBad;
0070 unsigned MaxClusterSize;
0071 bool RemoveApvShots;
0072 float minGoodCharge;
0073 };
0074
0075 template <class digiDetSet>
0076 inline void ThreeThresholdAlgorithm::clusterizeDetUnit_(const digiDetSet& digis, output_t::TSFastFiller& output) const {
0077 const auto& cond = conditions();
0078 if (cond.isModuleBad(digis.detId()))
0079 return;
0080
0081 auto const& det = cond.findDetId(digis.detId());
0082 if (!det.valid())
0083 return;
0084
0085 #ifdef EDM_ML_DEBUG
0086 if (!cond.isModuleUsable(digis.detId()))
0087 edm::LogWarning("ThreeThresholdAlgorithm") << " id " << digis.detId() << " not usable???" << std::endl;
0088 #endif
0089
0090 typename digiDetSet::const_iterator scan(digis.begin()), end(digis.end());
0091
0092 SiStripApvShotCleaner ApvCleaner;
0093 if (RemoveApvShots) {
0094 ApvCleaner.clean(digis, scan, end);
0095 }
0096
0097 output.reserve(16);
0098 State state(det);
0099 while (scan != end) {
0100 while (scan != end && !candidateEnded(state, scan->strip()))
0101 addToCandidate(state, *scan++);
0102 endCandidate(state, output);
0103 }
0104 }
0105
0106 inline bool ThreeThresholdAlgorithm::candidateEnded(State const& state, const uint16_t& testStrip) const {
0107 uint16_t holes = testStrip - state.lastStrip - 1;
0108 return (((!state.ADCs.empty()) &
0109 (holes > MaxSequentialHoles)
0110 ) &&
0111 (holes > MaxSequentialBad ||
0112 !state.det().allBadBetween(state.lastStrip, testStrip)
0113 ));
0114 }
0115
0116 inline void ThreeThresholdAlgorithm::addToCandidate(State& state, uint16_t strip, uint8_t adc) const {
0117 float Noise = state.det().noise(strip);
0118 if (adc < static_cast<uint8_t>(Noise * ChannelThreshold) || state.det().bad(strip))
0119 return;
0120
0121 if (state.candidateLacksSeed)
0122 state.candidateLacksSeed = adc < static_cast<uint8_t>(Noise * SeedThreshold);
0123 if (state.ADCs.empty())
0124 state.lastStrip = strip - 1;
0125 while (++state.lastStrip < strip)
0126 state.ADCs.push_back(0);
0127
0128 if (state.ADCs.size() <= MaxClusterSize)
0129 state.ADCs.push_back(adc);
0130 state.noiseSquared += Noise * Noise;
0131 }
0132
0133 inline void ThreeThresholdAlgorithm::clusterizeDetUnit(const edmNew::DetSet<SiStripDigi>& digis,
0134 output_t::TSFastFiller& output) const {
0135 clusterizeDetUnit_(digis, output);
0136 }
0137
0138 template <class T>
0139 inline void ThreeThresholdAlgorithm::endCandidate(State& state, T& out) const {
0140 if (candidateAccepted(state)) {
0141 applyGains(state);
0142 if (MaxAdjacentBad > 0)
0143 appendBadNeighbors(state);
0144 if (minGoodCharge <= 0 ||
0145 siStripClusterTools::chargePerCM(state.det().detId, state.ADCs.begin(), state.ADCs.end()) > minGoodCharge)
0146 out.push_back(std::move(SiStripCluster(firstStrip(state), state.ADCs.begin(), state.ADCs.end())));
0147 }
0148 clearCandidate(state);
0149 }
0150
0151 inline bool ThreeThresholdAlgorithm::candidateAccepted(State const& state) const {
0152 return (!state.candidateLacksSeed && state.ADCs.size() <= MaxClusterSize &&
0153 state.noiseSquared * ClusterThresholdSquared <=
0154 std::pow(float(std::accumulate(state.ADCs.begin(), state.ADCs.end(), int(0))), 2.f));
0155 }
0156
0157 inline void ThreeThresholdAlgorithm::applyGains(State& state) const {
0158 uint16_t strip = firstStrip(state);
0159 for (auto& adc : state.ADCs) {
0160 #ifdef EDM_ML_DEBUG
0161
0162 #endif
0163
0164 auto charge = int(float(adc) * state.det().weight(strip++) + 0.5f);
0165 if (adc < 254)
0166 adc = (charge > 1022 ? 255 : (charge > 253 ? 254 : charge));
0167 }
0168 }
0169
0170 inline void ThreeThresholdAlgorithm::appendBadNeighbors(State& state) const {
0171 uint8_t max = MaxAdjacentBad;
0172 while (0 < max--) {
0173 if (state.det().bad(firstStrip(state) - 1)) {
0174 state.ADCs.insert(state.ADCs.begin(), 0);
0175 }
0176 if (state.det().bad(state.lastStrip + 1)) {
0177 state.ADCs.push_back(0);
0178 state.lastStrip++;
0179 }
0180 }
0181 }
0182
0183 #endif