File indexing completed on 2024-04-06 12:26:28
0001 #include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h"
0002
0003 #include <cmath>
0004
0005 SiStripClusterInfo::SiStripClusterInfo(edm::ConsumesCollector&& iC, const std::string& qualityLabel)
0006 : siStripNoisesToken_{iC.esConsumes<SiStripNoises, SiStripNoisesRcd>()},
0007 siStripGainToken_{iC.esConsumes<SiStripGain, SiStripGainRcd>()},
0008 siStripQualityToken_{iC.esConsumes<SiStripQuality, SiStripQualityRcd>(edm::ESInputTag("", qualityLabel))} {}
0009
0010 void SiStripClusterInfo::initEvent(const edm::EventSetup& iSetup) {
0011 siStripNoises_ = &iSetup.getData(siStripNoisesToken_);
0012 siStripGain_ = &iSetup.getData(siStripGainToken_);
0013 siStripQuality_ = &iSetup.getData(siStripQualityToken_);
0014 }
0015
0016 void SiStripClusterInfo::setCluster(const SiStripCluster& cluster, int detId) {
0017 cluster_ptr = &cluster;
0018 detId_ = detId;
0019 }
0020
0021 std::pair<uint16_t, uint16_t> SiStripClusterInfo::chargeLR() const {
0022 std::vector<uint8_t>::const_iterator begin(stripCharges().begin()), end(stripCharges().end()), max;
0023 max = max_element(begin, end);
0024 return std::make_pair(accumulate(begin, max, uint16_t(0)), accumulate(max + 1, end, uint16_t(0)));
0025 }
0026
0027 float SiStripClusterInfo::variance() const {
0028 float q(0), x1(0), x2(0);
0029 for (auto begin(stripCharges().begin()), end(stripCharges().end()), it(begin); it != end; ++it) {
0030 unsigned i = it - begin;
0031 q += (*it);
0032 x1 += (*it) * (i + 0.5);
0033 x2 += (*it) * (i * i + i + 1. / 3);
0034 }
0035 return (x2 - x1 * x1 / q) / q;
0036 }
0037
0038 std::vector<float> SiStripClusterInfo::stripNoisesRescaledByGain() const {
0039 SiStripNoises::Range detNoiseRange = siStripNoises_->getRange(detId_);
0040 SiStripApvGain::Range detGainRange = siStripGain_->getRange(detId_);
0041
0042 std::vector<float> results;
0043 results.reserve(width());
0044 for (size_t i = 0, e = width(); i < e; i++) {
0045 results.push_back(siStripNoises_->getNoise(firstStrip() + i, detNoiseRange) /
0046 siStripGain_->getStripGain(firstStrip() + i, detGainRange));
0047 }
0048 return results;
0049 }
0050
0051 std::vector<float> SiStripClusterInfo::stripNoises() const {
0052 SiStripNoises::Range detNoiseRange = siStripNoises_->getRange(detId_);
0053
0054 std::vector<float> noises;
0055 noises.reserve(width());
0056 for (size_t i = 0; i < width(); i++) {
0057 noises.push_back(siStripNoises_->getNoise(firstStrip() + i, detNoiseRange));
0058 }
0059 return noises;
0060 }
0061
0062 std::vector<float> SiStripClusterInfo::stripGains() const {
0063 SiStripApvGain::Range detGainRange = siStripGain_->getRange(detId_);
0064
0065 std::vector<float> gains;
0066 gains.reserve(width());
0067 for (size_t i = 0; i < width(); i++) {
0068 gains.push_back(siStripGain_->getStripGain(firstStrip() + i, detGainRange));
0069 }
0070 return gains;
0071 }
0072
0073 std::vector<bool> SiStripClusterInfo::stripQualitiesBad() const {
0074 std::vector<bool> isBad;
0075 isBad.reserve(width());
0076 for (int i = 0; i < width(); i++) {
0077 isBad.push_back(siStripQuality_->IsStripBad(detId_, firstStrip() + i));
0078 }
0079 return isBad;
0080 }
0081
0082 float SiStripClusterInfo::calculate_noise(const std::vector<float>& noise) const {
0083 float noiseSumInQuadrature = 0;
0084 int numberStripsOverThreshold = 0;
0085 for (int i = 0; i < width(); i++) {
0086 if (stripCharges()[i] != 0) {
0087 noiseSumInQuadrature += noise.at(i) * noise.at(i);
0088 numberStripsOverThreshold++;
0089 }
0090 }
0091 return std::sqrt(noiseSumInQuadrature / numberStripsOverThreshold);
0092 }
0093
0094 bool SiStripClusterInfo::IsAnythingBad() const {
0095 std::vector<bool> stripBad = stripQualitiesBad();
0096 return IsApvBad() || IsFiberBad() || IsModuleBad() ||
0097 accumulate(stripBad.begin(), stripBad.end(), false, std::logical_or<bool>());
0098 }
0099
0100 bool SiStripClusterInfo::IsApvBad() const {
0101 return siStripQuality_->IsApvBad(detId_, firstStrip() / 128) ||
0102 siStripQuality_->IsApvBad(detId_, (firstStrip() + width()) / 128);
0103 }
0104
0105 bool SiStripClusterInfo::IsFiberBad() const {
0106 return siStripQuality_->IsFiberBad(detId_, firstStrip() / 256) ||
0107 siStripQuality_->IsFiberBad(detId_, (firstStrip() + width()) / 256);
0108 }
0109
0110 bool SiStripClusterInfo::IsModuleBad() const { return siStripQuality_->IsModuleBad(detId_); }
0111
0112 bool SiStripClusterInfo::IsModuleUsable() const { return siStripQuality_->IsModuleUsable(detId_); }