CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiStripZeroSuppression.cc
Go to the documentation of this file.
2 
14 #include <memory>
15 
17  : algorithms(
18  SiStripRawProcessingFactory::create(conf.getParameter<edm::ParameterSet>("Algorithms"), consumesCollector())),
19  produceRawDigis(conf.getParameter<bool>("produceRawDigis")),
20  storeCM(conf.getParameter<bool>("storeCM")),
21  fixCM(conf.getParameter<bool>("fixCM")),
22  produceCalculatedBaseline(conf.getParameter<bool>("produceCalculatedBaseline")),
23  produceBaselinePoints(conf.getParameter<bool>("produceBaselinePoints")),
24  storeInZScollBadAPV(conf.getParameter<bool>("storeInZScollBadAPV")),
25  produceHybridFormat(conf.getParameter<bool>("produceHybridFormat")) {
26  for (const auto& inputTag : conf.getParameter<std::vector<edm::InputTag>>("RawDigiProducersList")) {
27  const auto& tagName = inputTag.instance();
28  produces<edm::DetSetVector<SiStripDigi>>(tagName);
29  if (produceRawDigis)
30  produces<edm::DetSetVector<SiStripRawDigi>>(tagName);
31  if (storeCM)
32  produces<edm::DetSetVector<SiStripProcessedRawDigi>>("APVCM" + tagName);
34  produces<edm::DetSetVector<SiStripProcessedRawDigi>>("BADAPVBASELINE" + tagName);
36  produces<edm::DetSetVector<SiStripDigi>>("BADAPVBASELINEPOINTS" + tagName);
37 
38  RawType inputType = RawType::Unknown;
39  if (tagName == "ProcessedRaw") {
40  inputType = RawType::ProcessedRaw;
42  throw cms::Exception("Processed Raw Cannot be converted in hybrid Format");
43  } else if (tagName == "VirginRaw") {
44  inputType = RawType::VirginRaw;
45  } else if (tagName == "ScopeMode") {
46  inputType = RawType::ScopeMode;
48  throw cms::Exception("Scope Mode cannot be converted in hybrid Format");
49  }
50  if (RawType::Unknown != inputType) {
51  rawInputs.emplace_back(tagName, inputType, consumes<edm::DetSetVector<SiStripRawDigi>>(inputTag));
52  } else if (tagName == "ZeroSuppressed") {
54  } else {
55  throw cms::Exception("Unknown input type")
56  << tagName << " unknown. "
57  << "SiStripZeroZuppression can only process types \"VirginRaw\", \"ProcessedRaw\" and \"ZeroSuppressed\"";
58  }
59  }
60 
61  if (produceHybridFormat &&
62  ("HybridEmulation" !=
63  conf.getParameter<edm::ParameterSet>("Algorithms").getParameter<std::string>("APVInspectMode")))
64  throw cms::Exception("Invalid option") << "When producing data in the hybrid format, the APV restorer must be "
65  "configured with APVInspectMode='HybridEmulation'";
66 
67  if ((!hybridInputs.empty()) && produceRawDigis) {
68  edm::LogInfo("SiStripZeroSuppression") << "Raw digis will not be saved for hybrid inputs";
69  }
70 
71  if (!(rawInputs.empty() && hybridInputs.empty())) {
72  output_base.reserve(16000);
73  if (produceRawDigis && !rawInputs.empty())
74  output_base_raw.reserve(16000);
75  if (storeCM)
76  output_apvcm.reserve(16000);
78  output_baseline.reserve(16000);
80  output_baseline_points.reserve(16000);
81  }
82 }
83 
85  algorithms->initialize(es, e);
86 
87  for (const auto& input : rawInputs) {
88  clearOutputs();
90  e.getByToken(std::get<rawtoken_t>(input), inDigis);
91  if (!inDigis->empty())
92  processRaw(*inDigis, std::get<RawType>(input));
93  putOutputs(e, std::get<std::string>(input));
94  }
95  for (const auto& input : hybridInputs) {
96  clearOutputs();
98  e.getByToken(std::get<zstoken_t>(input), inDigis);
99  if (!inDigis->empty()) {
100  processHybrid(*inDigis);
101  }
102  putOutputs(e, std::get<std::string>(input));
103  }
104 }
105 
107  output_base.clear();
108  output_base_raw.clear();
109  output_baseline.clear();
110  output_baseline_points.clear();
111  output_apvcm.clear();
112 }
114  evt.put(std::make_unique<edm::DetSetVector<SiStripDigi>>(output_base), tagName);
115  if (produceRawDigis && !rawInputs.empty())
118  evt.put(std::make_unique<edm::DetSetVector<SiStripProcessedRawDigi>>(output_baseline), "BADAPVBASELINE" + tagName);
120  evt.put(std::make_unique<edm::DetSetVector<SiStripDigi>>(output_baseline_points), "BADAPVBASELINEPOINTS" + tagName);
121  if (storeCM)
122  evt.put(std::make_unique<edm::DetSetVector<SiStripProcessedRawDigi>>(output_apvcm), "APVCM" + tagName);
123 }
124 
126  for (const auto& rawDigis : input) {
127  edm::DetSet<SiStripDigi> suppressedDigis(rawDigis.id);
128 
129  int16_t nAPVflagged = 0;
130  if (RawType::ProcessedRaw == inType) {
131  nAPVflagged = algorithms->suppressProcessedRawData(rawDigis, suppressedDigis);
132  } else if (RawType::ScopeMode == inType) {
133  nAPVflagged = algorithms->suppressVirginRawData(rawDigis, suppressedDigis);
134  } else if (RawType::VirginRaw == inType) {
135  if (produceHybridFormat) {
136  nAPVflagged = algorithms->convertVirginRawToHybrid(rawDigis, suppressedDigis);
137  } else {
138  nAPVflagged = algorithms->suppressVirginRawData(rawDigis, suppressedDigis);
139  }
140  }
141 
142  storeExtraOutput(rawDigis.id, nAPVflagged);
143  if (!suppressedDigis.empty() && (storeInZScollBadAPV || nAPVflagged == 0))
144  output_base.push_back(std::move(suppressedDigis));
145 
146  if (produceRawDigis && nAPVflagged > 0) {
147  output_base_raw.push_back(formatRawDigis(rawDigis));
148  }
149  }
150 }
151 
153  for (const auto& inDigis : input) {
154  edm::DetSet<SiStripDigi> suppressedDigis(inDigis.id);
155 
156  unsigned int detId = inDigis.id;
157  uint16_t maxNStrips = SiStripDetId(detId).numberOfAPVs() * 128;
158 
159  uint16_t nAPVflagged = 0;
160  nAPVflagged = algorithms->suppressHybridData(maxNStrips, inDigis, suppressedDigis);
161 
162  storeExtraOutput(inDigis.id, nAPVflagged);
163  if (!suppressedDigis.empty())
164  output_base.push_back(std::move(suppressedDigis));
165  }
166 }
167 
169  edm::DetSet<SiStripRawDigi> outRawDigis(rawDigis.id);
170  outRawDigis.reserve(rawDigis.size());
171  const std::vector<bool>& apvf = algorithms->getAPVFlags();
172  uint32_t strip = 0;
173  for (const auto rawDigi : rawDigis) {
174  int16_t apvN = strip / 128;
175  if (apvf[apvN])
176  outRawDigis.push_back(rawDigi);
177  else
178  outRawDigis.push_back(SiStripRawDigi(0));
179  ++strip;
180  }
181  return outRawDigis;
182 }
183 
185  const std::vector<int16_t>& rawDigis) {
186  edm::DetSet<SiStripRawDigi> outRawDigis(detId);
187  outRawDigis.reserve(rawDigis.size());
188  const std::vector<bool>& apvf = algorithms->getAPVFlags();
189  uint32_t strip = 0;
190  for (const auto rawDigi : rawDigis) {
191  int16_t apvN = strip / 128;
192  if (apvf[apvN])
193  outRawDigis.push_back(SiStripRawDigi(rawDigi));
194  else
195  outRawDigis.push_back(SiStripRawDigi(0));
196  ++strip;
197  }
198  return outRawDigis;
199 }
200 
201 inline void SiStripZeroSuppression::storeExtraOutput(uint32_t id, int16_t nAPVflagged) {
202  const auto& vmedians = algorithms->getAPVsCM();
203  if (storeCM)
204  storeCMN(id, vmedians);
205  if (nAPVflagged > 0) {
207  storeBaseline(id, vmedians);
210  }
211 }
212 
213 inline void SiStripZeroSuppression::storeBaseline(uint32_t id, const medians_t& vmedians) {
214  const auto& baselinemap = algorithms->getBaselineMap();
215 
216  edm::DetSet<SiStripProcessedRawDigi> baselineDetSet(id);
217  baselineDetSet.reserve(vmedians.size() * 128);
218  for (const auto& vmed : vmedians) {
219  const uint16_t apvN = vmed.first;
220  const float median = vmed.second;
221  auto itBaselineMap = baselinemap.find(apvN);
222  if (baselinemap.end() == itBaselineMap) {
223  for (size_t strip = 0; strip < 128; ++strip)
224  baselineDetSet.push_back(SiStripProcessedRawDigi(median));
225  } else {
226  for (size_t strip = 0; strip < 128; ++strip)
227  baselineDetSet.push_back(SiStripProcessedRawDigi((itBaselineMap->second)[strip]));
228  }
229  }
230 
231  if (!baselineDetSet.empty())
232  output_baseline.push_back(baselineDetSet);
233 }
234 
236  edm::DetSet<SiStripDigi> baspointDetSet(id);
237  for (const auto& itBaselinePointVect : algorithms->getSmoothedPoints()) {
238  const uint16_t apvN = itBaselinePointVect.first;
239  for (const auto& itBaselinePointMap : itBaselinePointVect.second) {
240  const uint16_t bpstrip = itBaselinePointMap.first + apvN * 128;
241  const int16_t bp = itBaselinePointMap.second;
242  baspointDetSet.push_back(SiStripDigi(bpstrip, bp + 128));
243  }
244  }
245 
246  if (!baspointDetSet.empty())
247  output_baseline_points.push_back(std::move(baspointDetSet));
248 }
249 
250 inline void SiStripZeroSuppression::storeCMN(uint32_t id, const medians_t& vmedians) {
251  std::vector<bool> apvf(6, false);
252  if (fixCM) {
253  const auto& apvFlagged = algorithms->getAPVFlags();
254  std::copy(std::begin(apvFlagged), std::end(apvFlagged), std::begin(apvf));
255  }
256 
258  short apvNb = 0;
259  for (const auto& vmed : vmedians) {
260  if (vmed.first > apvNb) {
261  for (int i{0}; i < vmed.first - apvNb; ++i) {
262  apvDetSet.push_back(SiStripProcessedRawDigi(-999.));
263  apvNb++;
264  }
265  }
266 
267  if (fixCM && apvf[vmed.first]) {
268  apvDetSet.push_back(SiStripProcessedRawDigi(-999.));
269  } else {
270  apvDetSet.push_back(SiStripProcessedRawDigi(vmed.second));
271  }
272  apvNb++;
273  }
274 
275  if (!apvDetSet.empty())
276  output_apvcm.push_back(std::move(apvDetSet));
277 }
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:133
void push_back(const T &t)
Definition: DetSet.h:66
std::vector< edm::DetSet< SiStripDigi > > output_base
std::vector< edm::DetSet< SiStripProcessedRawDigi > > output_baseline
def create(alignables, pedeDump, additionalData, outputFile, config)
void produce(edm::Event &, const edm::EventSetup &) override
std::vector< std::tuple< std::string, zstoken_t > > hybridInputs
A signed Digi for the silicon strip detector, containing only adc information, and suitable for stori...
void processRaw(const edm::DetSetVector< SiStripRawDigi > &input, RawType inType)
std::vector< edm::DetSet< SiStripProcessedRawDigi > > output_apvcm
static std::string const input
Definition: EdmProvDump.cc:50
unsigned int numberOfAPVs() const
Definition: SiStripDetId.h:190
void reserve(size_t s)
Definition: DetSet.h:65
void storeExtraOutput(uint32_t, int16_t)
std::unique_ptr< SiStripRawProcessingAlgorithms > algorithms
std::vector< edm::DetSet< SiStripDigi > > output_baseline_points
std::vector< std::tuple< std::string, RawType, rawtoken_t > > rawInputs
A Digi for the silicon strip detector, containing both strip and adc information, and suitable for st...
Definition: SiStripDigi.h:12
Log< level::Info, false > LogInfo
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:19
void storeBaseline(uint32_t, const medians_t &)
std::vector< edm::DetSet< SiStripRawDigi > > output_base_raw
void storeCMN(uint32_t, const medians_t &)
bool empty() const
Definition: DetSet.h:62
std::vector< std::pair< short, float > > medians_t
void processHybrid(const edm::DetSetVector< SiStripDigi > &input)
SiStripZeroSuppression(const edm::ParameterSet &)
T median(std::vector< T > values)
Definition: median.h:16
edm::DetSet< SiStripRawDigi > formatRawDigis(const edm::DetSet< SiStripRawDigi > &rawDigis)
HLT enums.
det_id_type id
Definition: DetSet.h:79
void putOutputs(edm::Event &evt, const std::string &tagName)
A Digi for the silicon strip detector, containing only adc information, and suitable for storing raw ...
size_type size() const
Definition: DetSet.h:61
def move(src, dest)
Definition: eostools.py:511