File indexing completed on 2024-11-19 23:20:18
0001 #include "RecoLocalTracker/SiStripZeroSuppression/interface/SiStripRawProcessingAlgorithms.h"
0002
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/EventSetup.h"
0005 #include "DataFormats/Common/interface/Handle.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "DataFormats/Common/interface/DetSetVector.h"
0008 #include "DataFormats/SiStripDigi/interface/SiStripDigi.h"
0009 #include "DataFormats/SiStripDigi/interface/SiStripRawDigi.h"
0010 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
0011 #include "RecoLocalTracker/SiStripZeroSuppression/interface/SiStripRawProcessingFactory.h"
0012 #include <memory>
0013
0014 SiStripRawProcessingAlgorithms::SiStripRawProcessingAlgorithms(edm::ConsumesCollector iC,
0015 std::unique_ptr<SiStripPedestalsSubtractor> ped,
0016 std::unique_ptr<SiStripCommonModeNoiseSubtractor> cmn,
0017 std::unique_ptr<SiStripFedZeroSuppression> zs,
0018 std::unique_ptr<SiStripAPVRestorer> res,
0019 bool doAPVRest,
0020 bool useCMMap)
0021 : subtractorPed(std::move(ped)),
0022 subtractorCMN(std::move(cmn)),
0023 suppressor(std::move(zs)),
0024 restorer(std::move(res)),
0025 doAPVRestore(doAPVRest),
0026 useCMMeanMap(useCMMap),
0027 tkGeomToken_(iC.esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()) {}
0028
0029 void SiStripRawProcessingAlgorithms::initialize(const edm::EventSetup& es) {
0030 subtractorPed->init(es);
0031 subtractorCMN->init(es);
0032 suppressor->init(es);
0033 if (restorer.get())
0034 restorer->init(es);
0035
0036 trGeo = &es.getData(tkGeomToken_);
0037 }
0038
0039 void SiStripRawProcessingAlgorithms::initialize(const edm::EventSetup& es, const edm::Event& e) {
0040 initialize(es);
0041 if (restorer.get() && doAPVRestore && useCMMeanMap)
0042 restorer->loadMeanCMMap(e);
0043 }
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058 uint16_t SiStripRawProcessingAlgorithms::suppressHybridData(const uint16_t maxNStrips,
0059 const edm::DetSet<SiStripDigi>& hybridDigis,
0060 edm::DetSet<SiStripDigi>& suppressedDigis,
0061 uint16_t firstAPV) {
0062 uint16_t nAPVFlagged = 0;
0063 auto currentDigi = hybridDigis.begin();
0064 const auto endDigi = hybridDigis.end();
0065 auto currentAPV = firstAPV;
0066
0067
0068 while (currentDigi != endDigi) {
0069
0070 const auto nextAPVBoundary = SiStripDigi((currentAPV + 1) * 128, 0);
0071
0072
0073 if (nextAPVBoundary.strip() > maxNStrips) {
0074 edm::LogError("SiStripRawProcessingAlgorithms")
0075 << "In DetId " << suppressedDigis.id << " encountered APV boundary with strip number "
0076 << nextAPVBoundary.strip() << ", which exceeds the maximum allowed value for this module (" << maxNStrips
0077 << "). Exiting loop.";
0078 break;
0079 }
0080
0081 const auto nextAPVDigi = std::lower_bound(currentDigi, endDigi, nextAPVBoundary);
0082 const auto nDigisInAPV = std::distance(currentDigi, nextAPVDigi);
0083
0084
0085 if (nDigisInAPV > 64) {
0086
0087 digivector_t workDigis(128, -1024);
0088
0089
0090 for (auto it = currentDigi; it != nextAPVDigi; ++it) {
0091 workDigis[it->strip() - 128 * currentAPV] = it->adc() * 2 - 1024;
0092 }
0093
0094
0095 digivector_t workDigisPedSubtracted(workDigis);
0096 subtractorCMN->subtract(hybridDigis.id, currentAPV, workDigis);
0097
0098
0099 const auto apvFlagged = restorer->inspectAndRestore(
0100 hybridDigis.id, currentAPV, workDigisPedSubtracted, workDigis, subtractorCMN->getAPVsCM());
0101 nAPVFlagged += apvFlagged;
0102
0103
0104 if (getAPVFlags()[currentAPV]) {
0105
0106 suppressor->suppress(workDigis, currentAPV, suppressedDigis);
0107 } else {
0108
0109 for (uint16_t i = 0; i < 128; ++i) {
0110 const auto digi = workDigisPedSubtracted[i];
0111 if (digi > 0) {
0112 suppressedDigis.push_back(SiStripDigi(currentAPV * 128 + i, suppressor->truncate(digi)));
0113 }
0114 }
0115 }
0116 } else {
0117
0118 std::transform(currentDigi, nextAPVDigi, std::back_inserter(suppressedDigis), [this](const SiStripDigi& inDigi) {
0119 return SiStripDigi(inDigi.strip(), suppressor->truncate(inDigi.adc()));
0120 });
0121 }
0122
0123
0124 currentDigi = nextAPVDigi;
0125 ++currentAPV;
0126 }
0127
0128 return nAPVFlagged;
0129 }
0130
0131
0132
0133
0134
0135
0136
0137
0138
0139
0140
0141
0142
0143 uint16_t SiStripRawProcessingAlgorithms::suppressVirginRawData(uint32_t id,
0144 uint16_t firstAPV,
0145 digivector_t& procRawDigis,
0146 edm::DetSet<SiStripDigi>& output) {
0147 subtractorPed->subtract(id, firstAPV * 128, procRawDigis);
0148 return suppressProcessedRawData(id, firstAPV, procRawDigis, output);
0149 }
0150
0151
0152
0153
0154
0155
0156
0157
0158
0159
0160
0161 uint16_t SiStripRawProcessingAlgorithms::suppressVirginRawData(const edm::DetSet<SiStripRawDigi>& rawDigis,
0162 edm::DetSet<SiStripDigi>& output) {
0163 digivector_t rawdigis;
0164 rawdigis.reserve(rawDigis.size());
0165 std::transform(std::begin(rawDigis), std::end(rawDigis), std::back_inserter(rawdigis), [](SiStripRawDigi digi) {
0166 return digi.adc();
0167 });
0168 return suppressVirginRawData(rawDigis.id, 0, rawdigis, output);
0169 }
0170
0171
0172
0173
0174
0175
0176
0177
0178
0179
0180
0181
0182
0183 uint16_t SiStripRawProcessingAlgorithms::suppressProcessedRawData(uint32_t id,
0184 uint16_t firstAPV,
0185 digivector_t& procRawDigis,
0186 edm::DetSet<SiStripDigi>& output) {
0187 digivector_t procRawDigisPedSubtracted;
0188
0189 int16_t nAPVFlagged = 0;
0190 if (doAPVRestore)
0191 procRawDigisPedSubtracted.assign(procRawDigis.begin(), procRawDigis.end());
0192 subtractorCMN->subtract(id, firstAPV, procRawDigis);
0193 if (doAPVRestore)
0194 nAPVFlagged =
0195 restorer->inspectAndRestore(id, firstAPV, procRawDigisPedSubtracted, procRawDigis, subtractorCMN->getAPVsCM());
0196 suppressor->suppress(procRawDigis, firstAPV, output);
0197 return nAPVFlagged;
0198 }
0199
0200
0201
0202
0203
0204
0205
0206
0207
0208
0209
0210 uint16_t SiStripRawProcessingAlgorithms::suppressProcessedRawData(const edm::DetSet<SiStripRawDigi>& rawDigis,
0211 edm::DetSet<SiStripDigi>& output) {
0212 digivector_t rawdigis;
0213 rawdigis.reserve(rawDigis.size());
0214 std::transform(std::begin(rawDigis), std::end(rawDigis), std::back_inserter(rawdigis), [](SiStripRawDigi digi) {
0215 return digi.adc();
0216 });
0217 return suppressProcessedRawData(rawDigis.id, 0, rawdigis, output);
0218 }
0219
0220
0221
0222
0223
0224
0225
0226
0227
0228
0229
0230
0231
0232
0233
0234
0235 uint16_t SiStripRawProcessingAlgorithms::convertVirginRawToHybrid(uint32_t id,
0236 uint16_t firstAPV,
0237 digivector_t& procRawDigis,
0238 edm::DetSet<SiStripDigi>& output) {
0239 digivector_t procRawDigisPedSubtracted;
0240
0241 for (auto& digi : procRawDigis) {
0242 digi += 1024;
0243 }
0244
0245 subtractorPed->subtract(id, firstAPV * 128, procRawDigis);
0246
0247 for (auto& digi : procRawDigis) {
0248 digi /= 2;
0249 }
0250
0251 procRawDigisPedSubtracted.assign(procRawDigis.begin(), procRawDigis.end());
0252
0253 subtractorCMN->subtract(id, firstAPV, procRawDigis);
0254
0255 const auto nAPVFlagged = restorer->inspect(id, firstAPV, procRawDigis, subtractorCMN->getAPVsCM());
0256
0257 for (auto& digi : procRawDigis) {
0258 digi *= 2;
0259 }
0260
0261 const std::vector<bool>& apvf = getAPVFlags();
0262 const std::size_t nAPVs = procRawDigis.size() / 128;
0263 for (uint16_t iAPV = firstAPV; iAPV < nAPVs + firstAPV; ++iAPV) {
0264 if (apvf[iAPV]) {
0265
0266 for (uint16_t i = 0; i < 128; ++i) {
0267 const int16_t digi = procRawDigisPedSubtracted[128 * (iAPV - firstAPV) + i];
0268 output.push_back(SiStripDigi(128 * iAPV + i, (digi < 0 ? 0 : suppressor->truncate(digi))));
0269 }
0270 } else {
0271 const auto firstDigiIt = std::begin(procRawDigis) + 128 * (iAPV - firstAPV);
0272 std::vector<int16_t> singleAPVdigi(firstDigiIt, firstDigiIt + 128);
0273 suppressor->suppress(singleAPVdigi, iAPV, output);
0274 }
0275 }
0276
0277 return nAPVFlagged;
0278 }
0279
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289
0290
0291
0292 uint16_t SiStripRawProcessingAlgorithms::convertVirginRawToHybrid(const edm::DetSet<SiStripRawDigi>& rawDigis,
0293 edm::DetSet<SiStripDigi>& suppressedDigis) {
0294 digivector_t rawdigis;
0295 rawdigis.reserve(rawDigis.size());
0296 std::transform(std::begin(rawDigis), std::end(rawDigis), std::back_inserter(rawdigis), [](SiStripRawDigi digi) {
0297 return digi.adc();
0298 });
0299 return convertVirginRawToHybrid(rawDigis.id, 0, rawdigis, suppressedDigis);
0300 }