SiStripZeroSuppression/src/TT6CMNSubtractor.cc

0001 #include "RecoLocalTracker/SiStripZeroSuppression/interface/TT6CMNSubtractor.h"
0002
0003 #include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
0004 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0005 #include <cmath>
0006
0007 void TT6CMNSubtractor::init(const edm::EventSetup& es) {
0008   if (noiseWatcher_.check(es)) {
0009     noiseHandle = &es.getData(noiseToken_);
0010   }
0011   if (qualityWatcher_.check(es)) {
0012     qualityHandle = &es.getData(qualityToken_);
0013   }
0014 }
0015
0016 void TT6CMNSubtractor::subtract(uint32_t detId, uint16_t firstAPV, std::vector<int16_t>& digis) {
0017   subtract_(detId, firstAPV, digis);
0018 }
0019 void TT6CMNSubtractor::subtract(uint32_t detId, uint16_t firstAPV, std::vector<float>& digis) {
0020   subtract_(detId, firstAPV, digis);
0021 }
0022
0023 template <typename T>
0024 inline void TT6CMNSubtractor::subtract_(uint32_t detId, uint16_t firstAPV, std::vector<T>& digis) {
0025   short FixedBias = 128;
0026   SiStripNoises::Range detNoiseRange = noiseHandle->getRange(detId);
0027   SiStripQuality::Range detQualityRange = qualityHandle->getRange(detId);
0028
0029   typename std::vector<T>::iterator fs, ls;
0030   double sumVal = 0.0;
0031   double sumWt = 0.0;
0032
0033   for (uint16_t istrip = 0; istrip < digis.size(); ++istrip) {
0034     if (!qualityHandle->IsStripBad(detQualityRange, istrip + firstAPV * 128)) {
0035       float stripNoise = noiseHandle->getNoiseFast(istrip + firstAPV * 128, detNoiseRange);
0036
0037       if (std::abs(digis[istrip] - FixedBias) < cut_to_avoid_signal_ * stripNoise) {
0038         double nWeight = 1 / (stripNoise * stripNoise);
0039         sumVal += (digis[istrip] - FixedBias) * nWeight;
0040         sumWt += nWeight;
0041       }
0042     }
0043
0044     if (istrip % 128 == 127) {
0045       double CM = (sumWt) ? sumVal / sumWt : 0.0;
0046
0047       fs = digis.begin() + istrip - 127;
0048       ls = digis.begin() + istrip + 1;
0049
0050       while (fs < ls) {
0051         *fs = static_cast<T>(*fs - FixedBias - CM);
0052         fs++;
0053       }
0054
0055       sumVal = 0.0;
0056       sumWt = 0.0;
0057     }
0058   }
0059 }