File indexing completed on 2025-07-01 00:07:49
0001 #include "SimCalorimetry/HcalTrigPrimAlgos/interface/HcalTriggerPrimitiveAlgo.h"
0002
0003 #include "CalibFormats/CaloObjects/interface/IntegerCaloSamples.h"
0004 #include "CondFormats/HcalObjects/interface/HcalTPParameters.h"
0005 #include "CondFormats/HcalObjects/interface/HcalTPChannelParameters.h"
0006
0007 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0008 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0009 #include "DataFormats/HcalDetId/interface/HcalElectronicsId.h"
0010
0011 #include "EventFilter/HcalRawToDigi/interface/HcalDCCHeader.h"
0012 #include "EventFilter/HcalRawToDigi/interface/HcalHTRData.h"
0013
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016
0017 #include "Geometry/HcalTowerAlgo/interface/HcalTrigTowerGeometry.h"
0018
0019 #include <iostream>
0020
0021 using namespace std;
0022
0023 HcalTriggerPrimitiveAlgo::HcalTriggerPrimitiveAlgo(bool pf,
0024 const std::vector<double>& w,
0025 int latency,
0026 uint32_t FG_threshold,
0027 const std::vector<uint32_t>& FG_HF_thresholds,
0028 uint32_t ZS_threshold,
0029 int numberOfSamples,
0030 int numberOfPresamples,
0031 int numberOfFilterPresamplesHBQIE11,
0032 int numberOfFilterPresamplesHEQIE11,
0033 int numberOfSamplesHF,
0034 int numberOfPresamplesHF,
0035 int numberOfSamplesZDC,
0036 int numberOfPresamplesZDC,
0037 bool useTDCInMinBiasBits,
0038 uint32_t minSignalThreshold,
0039 uint32_t PMT_NoiseThreshold)
0040 : incoder_(nullptr),
0041 outcoder_(nullptr),
0042 theThreshold(0),
0043 peakfind_(pf),
0044 weights_(w),
0045 latency_(latency),
0046 FG_threshold_(FG_threshold),
0047 FG_HF_thresholds_(FG_HF_thresholds),
0048 ZS_threshold_(ZS_threshold),
0049 numberOfSamples_(numberOfSamples),
0050 numberOfPresamples_(numberOfPresamples),
0051 numberOfFilterPresamplesHBQIE11_(numberOfFilterPresamplesHBQIE11),
0052 numberOfFilterPresamplesHEQIE11_(numberOfFilterPresamplesHEQIE11),
0053 numberOfSamplesHF_(numberOfSamplesHF),
0054 numberOfPresamplesHF_(numberOfPresamplesHF),
0055 numberOfSamplesZDC_(numberOfSamplesZDC),
0056 numberOfPresamplesZDC_(numberOfPresamplesZDC),
0057 useTDCInMinBiasBits_(useTDCInMinBiasBits),
0058 minSignalThreshold_(minSignalThreshold),
0059 PMT_NoiseThreshold_(PMT_NoiseThreshold),
0060 NCTScaleShift(0),
0061 RCTScaleShift(0),
0062 peak_finder_algorithm_(2),
0063 override_parameters_() {
0064
0065 if (!peakfind_) {
0066 numberOfSamples_ = 1;
0067 numberOfPresamples_ = 0;
0068 numberOfSamplesHF_ = 1;
0069 numberOfPresamplesHF_ = 0;
0070 numberOfSamplesZDC_ = 1;
0071 numberOfPresamplesZDC_ = 0;
0072 }
0073
0074 ZS_threshold_I_ = ZS_threshold_;
0075 }
0076
0077 HcalTriggerPrimitiveAlgo::~HcalTriggerPrimitiveAlgo() {}
0078
0079 void HcalTriggerPrimitiveAlgo::setUpgradeFlags(bool hb, bool he, bool hf) {
0080 upgrade_hb_ = hb;
0081 upgrade_he_ = he;
0082 upgrade_hf_ = hf;
0083 }
0084
0085 void HcalTriggerPrimitiveAlgo::setFixSaturationFlag(bool fix_saturation) { fix_saturation_ = fix_saturation; }
0086
0087 void HcalTriggerPrimitiveAlgo::overrideParameters(const edm::ParameterSet& ps) {
0088 override_parameters_ = ps;
0089
0090 if (override_parameters_.exists("ADCThresholdHF")) {
0091 override_adc_hf_ = true;
0092 override_adc_hf_value_ = override_parameters_.getParameter<uint32_t>("ADCThresholdHF");
0093 }
0094 if (override_parameters_.exists("TDCMaskHF")) {
0095 override_tdc_hf_ = true;
0096 override_tdc_hf_value_ = override_parameters_.getParameter<unsigned long long>("TDCMaskHF");
0097 }
0098 }
0099
0100 void HcalTriggerPrimitiveAlgo::addSignal(const HBHEDataFrame& frame) {
0101
0102
0103 if (frame.id().depth() == 5)
0104 return;
0105
0106 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(frame.id());
0107 assert(ids.size() == 1 || ids.size() == 2);
0108 IntegerCaloSamples samples1(ids[0], int(frame.size()));
0109
0110 samples1.setPresamples(frame.presamples());
0111 incoder_->adc2Linear(frame, samples1);
0112
0113 std::vector<bool> msb;
0114 incoder_->lookupMSB(frame, msb);
0115
0116 if (ids.size() == 2) {
0117
0118 IntegerCaloSamples samples2(ids[1], samples1.size());
0119 for (int i = 0; i < samples1.size(); ++i) {
0120 samples1[i] = uint32_t(samples1[i] * 0.5);
0121 samples2[i] = samples1[i];
0122 }
0123 samples2.setPresamples(frame.presamples());
0124 addSignal(samples2);
0125 addFG(ids[1], msb);
0126 }
0127 addSignal(samples1);
0128 addFG(ids[0], msb);
0129 }
0130
0131 void HcalTriggerPrimitiveAlgo::addSignal(const HFDataFrame& frame) {
0132 if (frame.id().depth() == 1 || frame.id().depth() == 2) {
0133 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(frame.id());
0134 std::vector<HcalTrigTowerDetId>::const_iterator it;
0135 for (it = ids.begin(); it != ids.end(); ++it) {
0136 HcalTrigTowerDetId trig_tower_id = *it;
0137 IntegerCaloSamples samples(trig_tower_id, frame.size());
0138 samples.setPresamples(frame.presamples());
0139 incoder_->adc2Linear(frame, samples);
0140
0141
0142
0143 IntegerCaloSamples zero_samples(trig_tower_id, frame.size());
0144 zero_samples.setPresamples(frame.presamples());
0145 addSignal(zero_samples);
0146
0147
0148 if (trig_tower_id.version() == 0) {
0149
0150 uint32_t fgid = (frame.id().maskDepth());
0151
0152 if (theTowerMapFGSum.find(trig_tower_id) == theTowerMapFGSum.end()) {
0153 SumFGContainer sumFG;
0154 theTowerMapFGSum.insert(std::pair<HcalTrigTowerDetId, SumFGContainer>(trig_tower_id, sumFG));
0155 }
0156
0157 SumFGContainer& sumFG = theTowerMapFGSum[trig_tower_id];
0158 SumFGContainer::iterator sumFGItr;
0159 for (sumFGItr = sumFG.begin(); sumFGItr != sumFG.end(); ++sumFGItr) {
0160 if (sumFGItr->id() == fgid) {
0161 break;
0162 }
0163 }
0164
0165 if (sumFGItr != sumFG.end()) {
0166 for (int i = 0; i < samples.size(); ++i) {
0167 (*sumFGItr)[i] += samples[i];
0168 }
0169 } else {
0170
0171 IntegerCaloSamples sumFGSamples(DetId(fgid), samples.size());
0172 sumFGSamples.setPresamples(samples.presamples());
0173 for (int i = 0; i < samples.size(); ++i) {
0174 sumFGSamples[i] = samples[i];
0175 }
0176 sumFG.push_back(sumFGSamples);
0177 }
0178
0179
0180 if (HF_Veto.find(fgid) == HF_Veto.end()) {
0181 vector<bool> vetoBits(samples.size(), false);
0182 HF_Veto[fgid] = vetoBits;
0183 }
0184 for (int i = 0; i < samples.size(); ++i) {
0185 if (samples[i] < minSignalThreshold_) {
0186 HF_Veto[fgid][i] = true;
0187 }
0188 }
0189 }
0190
0191 else if (trig_tower_id.version() == 1) {
0192 uint32_t fgid = (frame.id().maskDepth());
0193 HFDetails& details = theHFDetailMap[trig_tower_id][fgid];
0194
0195 if (frame.id().depth() == 1) {
0196 details.long_fiber = samples;
0197 details.LongDigi = frame;
0198 } else if (frame.id().depth() == 2) {
0199 details.short_fiber = samples;
0200 details.ShortDigi = frame;
0201 } else {
0202
0203 edm::LogWarning("HcalTPAlgo") << "Unable to figure out what to do with data frame for " << frame.id();
0204 return;
0205 }
0206 }
0207
0208 else {
0209 return;
0210 }
0211 }
0212 }
0213 }
0214
0215 void HcalTriggerPrimitiveAlgo::addSignal(const QIE10DataFrame& frame) {
0216
0217 DetId detId = DetId(frame.detid());
0218 if (detId.det() == DetId::Hcal) {
0219 HcalDetId detId = frame.detid();
0220
0221 if (detId.subdet() != HcalForward)
0222 return;
0223
0224 auto ids = theTrigTowerGeometry->towerIds(frame.id());
0225 for (const auto& id : ids) {
0226 if (id.version() == 0) {
0227 edm::LogError("HcalTPAlgo") << "Encountered QIE10 data frame mapped to TP version 0:" << id;
0228 continue;
0229 }
0230 int nsamples = frame.samples();
0231
0232 IntegerCaloSamples samples(id, nsamples);
0233 samples.setPresamples(frame.presamples());
0234 incoder_->adc2Linear(frame, samples, false);
0235
0236
0237
0238 IntegerCaloSamples zero_samples(id, nsamples);
0239 zero_samples.setPresamples(frame.presamples());
0240 addSignal(zero_samples);
0241
0242 auto fid = HcalDetId(frame.id());
0243 auto& details = theHFUpgradeDetailMap[id][fid.maskDepth()];
0244 auto& detail = details[fid.depth() - 1];
0245 detail.samples = samples;
0246 detail.digi = frame;
0247 detail.validity.resize(nsamples);
0248 detail.passTDC.resize(nsamples);
0249 incoder_->lookupMSB(frame, detail.fgbits);
0250 for (int idx = 0; idx < nsamples; ++idx) {
0251 detail.validity[idx] = validChannel(frame, idx);
0252 detail.passTDC[idx] = passTDC(frame, idx);
0253 }
0254 }
0255 } else if (detId.det() == DetId::Calo && detId.subdetId() == HcalZDCDetId::SubdetectorId) {
0256 HcalZDCDetId detId = frame.detid();
0257
0258 if (detId.section() != HcalZDCDetId::EM && detId.section() != HcalZDCDetId::HAD) {
0259 return;
0260 }
0261
0262 if (detId.section() == HcalZDCDetId::EM && detId.channel() > 5) {
0263 return;
0264 }
0265
0266 auto ids = theTrigTowerGeometry->towerIds_ZDC(frame.id());
0267 for (const auto& id : ids) {
0268 int nsamples = frame.samples();
0269
0270 IntegerCaloSamples samples(id, nsamples);
0271 IntegerCaloSamples samples_PUsub(id, nsamples);
0272
0273 samples.setPresamples(frame.presamples());
0274 samples_PUsub.setPresamples(frame.presamples());
0275
0276 incoder_->adc2Linear(frame, samples, false);
0277 incoder_->adc2Linear(frame, samples_PUsub, true);
0278
0279 for (int i = 1; i < samples.size(); ++i) {
0280 if (samples_PUsub[i - 1] > samples[i])
0281 samples[i] = 0;
0282 else
0283 samples[i] -= samples_PUsub[i - 1];
0284 }
0285
0286 addSignal(samples);
0287 }
0288 }
0289 }
0290
0291 void HcalTriggerPrimitiveAlgo::addSignal(const QIE11DataFrame& frame) {
0292 HcalDetId detId(frame.id());
0293
0294 if (detId.subdet() != HcalEndcap && detId.subdet() != HcalBarrel)
0295 return;
0296
0297 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(detId);
0298 assert(ids.size() == 1 || ids.size() == 2);
0299 IntegerCaloSamples samples1(ids[0], int(frame.samples()));
0300
0301 samples1.setPresamples(frame.presamples());
0302 incoder_->adc2Linear(frame, samples1);
0303
0304 std::vector<std::bitset<2>> msb(frame.samples(), 0);
0305 incoder_->lookupMSB(frame, msb);
0306
0307 if (ids.size() == 2) {
0308
0309 IntegerCaloSamples samples2(ids[1], samples1.size());
0310 for (int i = 0; i < samples1.size(); ++i) {
0311 samples1[i] = uint32_t(samples1[i]);
0312 samples2[i] = samples1[i];
0313 }
0314 samples2.setPresamples(frame.presamples());
0315 addSignal(samples2);
0316 addUpgradeFG(ids[1], detId.depth(), msb);
0317 addUpgradeTDCFG(ids[1], frame);
0318 }
0319 addSignal(samples1);
0320 addUpgradeFG(ids[0], detId.depth(), msb);
0321 addUpgradeTDCFG(ids[0], frame);
0322 }
0323
0324 void HcalTriggerPrimitiveAlgo::addSignal(const IntegerCaloSamples& samples) {
0325 HcalTrigTowerDetId id(samples.id());
0326 SumMap::iterator itr = theSumMap.find(id);
0327
0328 if (itr == theSumMap.end()) {
0329 theSumMap.insert(std::make_pair(id, samples));
0330 } else {
0331
0332 for (int i = 0; i < samples.size(); ++i) {
0333 (itr->second)[i] += samples[i];
0334 }
0335 }
0336
0337
0338 if (fix_saturation_) {
0339 SatMap::iterator itr_sat = theSatMap.find(id);
0340
0341 assert((itr == theSumMap.end()) == (itr_sat == theSatMap.end()));
0342
0343 if (itr_sat == theSatMap.end()) {
0344 vector<bool> check_sat;
0345 for (int i = 0; i < samples.size(); ++i) {
0346 if (!(samples[i] < QIE11_LINEARIZATION_ET)) {
0347 check_sat.push_back(true);
0348 } else
0349 check_sat.push_back(false);
0350 }
0351 theSatMap.insert(std::make_pair(id, check_sat));
0352 } else {
0353 for (int i = 0; i < samples.size(); ++i) {
0354 if (!(samples[i] < QIE11_LINEARIZATION_ET))
0355 (itr_sat->second)[i] = true;
0356 }
0357 }
0358 }
0359 }
0360
0361 void HcalTriggerPrimitiveAlgo::analyze(IntegerCaloSamples& samples, HcalTriggerPrimitiveDigi& result) {
0362 int shrink = weights_.size() - 1;
0363 std::vector<bool>& msb = fgMap_[samples.id()];
0364 IntegerCaloSamples sum(samples.id(), samples.size());
0365
0366
0367 for (int ibin = 0; ibin < int(samples.size()) - shrink; ++ibin) {
0368 int algosumvalue = 0;
0369 for (unsigned int i = 0; i < weights_.size(); i++) {
0370
0371 algosumvalue += int(samples[ibin + i] * weights_[i]);
0372 }
0373 if (algosumvalue < 0)
0374 sum[ibin] = 0;
0375
0376
0377 else
0378 sum[ibin] = algosumvalue;
0379 }
0380
0381
0382 int dgPresamples = samples.presamples();
0383 int tpPresamples = numberOfPresamples_;
0384 int shift = dgPresamples - tpPresamples;
0385 int dgSamples = samples.size();
0386 int tpSamples = numberOfSamples_;
0387 if (peakfind_) {
0388 if ((shift < shrink) || (shift + tpSamples + shrink > dgSamples - (peak_finder_algorithm_ - 1))) {
0389 edm::LogInfo("HcalTriggerPrimitiveAlgo::analyze")
0390 << "TP presample or size from the configuration file is out of the accessible range. Using digi values from "
0391 "data instead...";
0392 shift = shrink;
0393 tpPresamples = dgPresamples - shrink;
0394 tpSamples = dgSamples - (peak_finder_algorithm_ - 1) - shrink - shift;
0395 }
0396 }
0397
0398 std::vector<int> finegrain(tpSamples, false);
0399
0400 IntegerCaloSamples output(samples.id(), tpSamples);
0401 output.setPresamples(tpPresamples);
0402
0403 for (int ibin = 0; ibin < tpSamples; ++ibin) {
0404
0405
0406 int idx = ibin + shift;
0407
0408
0409 if (peakfind_) {
0410 bool isPeak = false;
0411 switch (peak_finder_algorithm_) {
0412 case 1:
0413 isPeak = (samples[idx] > samples[idx - 1] && samples[idx] >= samples[idx + 1] && samples[idx] > theThreshold);
0414 break;
0415 case 2:
0416 isPeak = (sum[idx] > sum[idx - 1] && sum[idx] >= sum[idx + 1] && sum[idx] > theThreshold);
0417 break;
0418 default:
0419 break;
0420 }
0421
0422 if (isPeak) {
0423 output[ibin] = std::min<unsigned int>(sum[idx], QIE8_LINEARIZATION_ET);
0424 finegrain[ibin] = msb[idx];
0425 }
0426
0427 else
0428 output[ibin] = 0;
0429 } else {
0430 output[ibin] = std::min<unsigned int>(sum[idx], QIE8_LINEARIZATION_ET);
0431 finegrain[ibin] = msb[idx];
0432 }
0433
0434
0435 if (peak_finder_algorithm_ == 1) {
0436 if (samples[idx] >= QIE8_LINEARIZATION_ET)
0437 output[ibin] = QIE8_LINEARIZATION_ET;
0438 }
0439 }
0440 outcoder_->compress(output, finegrain, result);
0441 }
0442
0443 void HcalTriggerPrimitiveAlgo::analyzeQIE11(IntegerCaloSamples& samples,
0444 vector<bool> sample_saturation,
0445 HcalTriggerPrimitiveDigi& result,
0446 const HcalFinegrainBit& fg_algo) {
0447 HcalDetId detId(samples.id());
0448
0449
0450 int theIeta = detId.ietaAbs();
0451
0452 unsigned int dgSamples = samples.size();
0453 unsigned int dgPresamples = samples.presamples();
0454
0455 unsigned int tpSamples = numberOfSamples_;
0456 unsigned int tpPresamples = numberOfPresamples_;
0457
0458 unsigned int filterSamples = weightsQIE11_[theIeta].size();
0459 unsigned int filterPresamples = theIeta > theTrigTowerGeometry->topology().lastHBRing()
0460 ? numberOfFilterPresamplesHEQIE11_
0461 : numberOfFilterPresamplesHBQIE11_;
0462
0463 unsigned int shift = dgPresamples - tpPresamples;
0464
0465
0466 unsigned int shrink = filterSamples - 1;
0467
0468 auto& msb = fgUpgradeMap_[samples.id()];
0469 auto& timingTDC = fgUpgradeTDCMap_[samples.id()];
0470 IntegerCaloSamples sum(samples.id(), samples.size());
0471
0472 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(detId);
0473
0474
0475 bool force_saturation[samples.size()];
0476 for (int i = 0; i < samples.size(); i++) {
0477 force_saturation[i] = false;
0478 }
0479
0480
0481 for (unsigned int ibin = 0; ibin < dgSamples - shrink; ++ibin) {
0482 int algosumvalue = 0;
0483 bool check_sat = false;
0484
0485 if (weightsQIE11_[theIeta][0] == 255) {
0486 for (unsigned int i = 0; i < filterSamples; i++) {
0487
0488
0489
0490 unsigned int sample = samples[ibin + i];
0491
0492 if (fix_saturation_ && (sample_saturation.size() > ibin + i))
0493 check_sat = (check_sat | sample_saturation[ibin + i] | (sample > QIE11_MAX_LINEARIZATION_ET));
0494
0495 if (sample > QIE11_MAX_LINEARIZATION_ET)
0496 sample = QIE11_MAX_LINEARIZATION_ET;
0497
0498
0499 int segmentationFactor = 1;
0500 if (ids.size() == 2) {
0501 segmentationFactor = 2;
0502 }
0503
0504 algosumvalue += int(sample / segmentationFactor);
0505 }
0506 if (algosumvalue < 0)
0507 sum[ibin] = 0;
0508
0509
0510 else
0511 sum[ibin] = algosumvalue;
0512
0513 if (check_sat)
0514 force_saturation[ibin] = true;
0515
0516 } else {
0517
0518
0519
0520 int sampleTS = samples[ibin + 1];
0521 int sampleTSminus1 = samples[ibin];
0522
0523 if (fix_saturation_ && (sample_saturation.size() > ibin + 1))
0524 check_sat |= sample_saturation[ibin + 1] | (sampleTS >= QIE11_MAX_LINEARIZATION_ET);
0525
0526 if (sampleTS > QIE11_MAX_LINEARIZATION_ET)
0527 sampleTS = QIE11_MAX_LINEARIZATION_ET;
0528
0529 if (sampleTSminus1 > QIE11_MAX_LINEARIZATION_ET || sample_saturation[ibin])
0530 sampleTSminus1 = QIE11_MAX_LINEARIZATION_ET;
0531
0532
0533 int segmentationFactor = 1;
0534 if (ids.size() == 2) {
0535 segmentationFactor = 2;
0536 }
0537
0538
0539 int theWeight = weightsQIE11_[theIeta][0];
0540
0541 algosumvalue = ((sampleTS << 8) - (sampleTSminus1 * theWeight)) / 256 / segmentationFactor;
0542
0543 if (algosumvalue < 0)
0544 sum[ibin] = 0;
0545
0546
0547 else
0548 sum[ibin] = algosumvalue;
0549
0550 if (check_sat)
0551 force_saturation[ibin] = true;
0552 }
0553 }
0554
0555 std::vector<int> finegrain(tpSamples, false);
0556
0557 IntegerCaloSamples output(samples.id(), tpSamples);
0558 output.setPresamples(tpPresamples);
0559
0560
0561
0562 unsigned int codedVetoThreshold = codedVetoThresholds_[theIeta];
0563
0564
0565 unsigned int actualVetoThreshold = codedVetoThreshold;
0566 bool applyVetoThreshold = codedVetoThreshold > 0 && codedVetoThreshold <= 2048;
0567
0568
0569 if (codedVetoThreshold > 0) {
0570 if (codedVetoThreshold <= 2048) {
0571
0572 if (codedVetoThreshold == 2048)
0573 actualVetoThreshold = 0;
0574 } else {
0575 edm::LogWarning("HcalTPAlgo") << "Specified veto threshold value " << codedVetoThreshold
0576 << " is not in range (1, 2048) ! Vetoing in PFA1' will not be enabled !";
0577 }
0578 }
0579
0580 for (unsigned int ibin = 0; ibin < tpSamples; ++ibin) {
0581
0582
0583 int idx = ibin + shift - filterPresamples;
0584
0585
0586 if (idx <= 0) {
0587 output[ibin] = 0;
0588 continue;
0589 }
0590
0591
0592 if (weightsQIE11_[theIeta][0] == 255) {
0593 bool isPeak = (sum[idx] > sum[idx - 1] && sum[idx] >= sum[idx + 1] && sum[idx] > theThreshold);
0594 if (isPeak) {
0595 output[ibin] = std::min<unsigned int>(sum[idx], QIE11_MAX_LINEARIZATION_ET);
0596
0597 if (fix_saturation_ && force_saturation[idx] && ids.size() == 2)
0598 output[ibin] = QIE11_MAX_LINEARIZATION_ET / 2;
0599 else if (fix_saturation_ && force_saturation[idx])
0600 output[ibin] = QIE11_MAX_LINEARIZATION_ET;
0601
0602 } else {
0603
0604 output[ibin] = 0;
0605 }
0606 } else {
0607
0608
0609
0610 if (applyVetoThreshold && sum[idx + 1] >= actualVetoThreshold &&
0611 (sum[idx] < sum[idx + 1] || force_saturation[idx + 1]) && !force_saturation[idx])
0612 output[ibin] = 0;
0613 else {
0614
0615
0616 output[ibin] = std::min<unsigned int>(sum[idx], QIE11_MAX_LINEARIZATION_ET);
0617 if (fix_saturation_ && force_saturation[idx]) {
0618 output[ibin] = QIE11_MAX_LINEARIZATION_ET;
0619 if (ids.size() == 2)
0620 output[ibin] /= 2;
0621 }
0622 }
0623 }
0624
0625
0626 finegrain[ibin] = fg_algo.compute(timingTDC[idx + filterPresamples], ids[0]).to_ulong() |
0627 fg_algo.compute(msb[idx + filterPresamples]).to_ulong() << 4;
0628 if (ibin == tpPresamples && (idx + filterPresamples) != dgPresamples)
0629 edm::LogError("HcalTriggerPritimveAlgo")
0630 << "TP SOI (tpPresamples = " << tpPresamples
0631 << ") is not aligned with digi SOI (dgPresamples = " << dgPresamples << ")";
0632 }
0633 outcoder_->compress(output, finegrain, result);
0634 }
0635
0636 void HcalTriggerPrimitiveAlgo::analyzeZDC(IntegerCaloSamples& samples, HcalTriggerPrimitiveDigi& result) {
0637 HcalTrigTowerDetId detId(samples.id());
0638
0639 unsigned int tpSamples;
0640 unsigned int tpPresamples;
0641
0642 tpSamples = samples.size();
0643 tpPresamples = samples.presamples();
0644 result.setSize(tpSamples);
0645 result.setPresamples(tpPresamples);
0646
0647 IntegerCaloSamples output(samples.id(), tpSamples);
0648 output.setPresamples(tpPresamples);
0649
0650 for (int i = 0; i < samples.size(); i++) {
0651 if (samples[i] > QIE10_ZDC_MAX_LINEARIZATION_ET)
0652 output[i] = QIE10_ZDC_MAX_LINEARIZATION_ET;
0653 else
0654 output[i] = samples[i];
0655 HcalTriggerPrimitiveSample zdcSample(output[i]);
0656 result.setSample(i, zdcSample);
0657 }
0658 }
0659
0660 void HcalTriggerPrimitiveAlgo::analyzeHF(IntegerCaloSamples& samples,
0661 HcalTriggerPrimitiveDigi& result,
0662 const int hf_lumi_shift) {
0663 HcalTrigTowerDetId detId(samples.id());
0664
0665
0666 int dgPresamples = samples.presamples();
0667 int tpPresamples = numberOfPresamplesHF_;
0668 int shift = dgPresamples - tpPresamples;
0669 int dgSamples = samples.size();
0670 int tpSamples = numberOfSamplesHF_;
0671 if (shift < 0 || shift + tpSamples > dgSamples) {
0672 edm::LogInfo("HcalTriggerPrimitiveAlgo::analyzeHF")
0673 << "TP presample or size from the configuration file is out of the accessible range. Using digi values from "
0674 "data instead...";
0675 tpPresamples = dgPresamples;
0676 shift = 0;
0677 tpSamples = dgSamples;
0678 }
0679
0680 std::vector<int> finegrain(tpSamples, false);
0681
0682 TowerMapFGSum::const_iterator tower2fg = theTowerMapFGSum.find(detId);
0683 assert(tower2fg != theTowerMapFGSum.end());
0684
0685 const SumFGContainer& sumFG = tower2fg->second;
0686
0687
0688 for (SumFGContainer::const_iterator sumFGItr = sumFG.begin(); sumFGItr != sumFG.end(); ++sumFGItr) {
0689 const std::vector<bool>& veto = HF_Veto[sumFGItr->id().rawId()];
0690 for (int ibin = 0; ibin < tpSamples; ++ibin) {
0691 int idx = ibin + shift;
0692
0693 bool vetoed = idx < int(veto.size()) && veto[idx];
0694 if (!(vetoed && (*sumFGItr)[idx] > PMT_NoiseThreshold_)) {
0695 samples[idx] += (*sumFGItr)[idx];
0696 finegrain[ibin] = (finegrain[ibin] || (*sumFGItr)[idx] >= FG_threshold_);
0697 }
0698 }
0699 }
0700
0701 IntegerCaloSamples output(samples.id(), tpSamples);
0702 output.setPresamples(tpPresamples);
0703
0704 for (int ibin = 0; ibin < tpSamples; ++ibin) {
0705 int idx = ibin + shift;
0706 output[ibin] = samples[idx] >> hf_lumi_shift;
0707 static const int MAX_OUTPUT = QIE8_LINEARIZATION_ET;
0708 if (output[ibin] > MAX_OUTPUT)
0709 output[ibin] = MAX_OUTPUT;
0710 }
0711 outcoder_->compress(output, finegrain, result);
0712 }
0713
0714 void HcalTriggerPrimitiveAlgo::analyzeHF2016(const IntegerCaloSamples& samples,
0715 HcalTriggerPrimitiveDigi& result,
0716 const int hf_lumi_shift,
0717 const HcalFeatureBit* embit) {
0718
0719 const int SHIFT = samples.presamples() - numberOfPresamplesHF_;
0720 assert(SHIFT >= 0);
0721 assert((SHIFT + numberOfSamplesHF_) <= samples.size());
0722
0723
0724 const HcalTrigTowerDetId detId(samples.id());
0725 HFDetailMap::const_iterator it = theHFDetailMap.find(detId);
0726
0727 if (it == theHFDetailMap.end()) {
0728 return;
0729 }
0730
0731 std::vector<std::bitset<2>> finegrain(numberOfSamplesHF_, false);
0732
0733
0734 IntegerCaloSamples output(samples.id(), numberOfSamplesHF_);
0735 output.setPresamples(numberOfPresamplesHF_);
0736
0737 for (const auto& item : it->second) {
0738 auto& details = item.second;
0739 for (int ibin = 0; ibin < numberOfSamplesHF_; ++ibin) {
0740 const int IDX = ibin + SHIFT;
0741 int long_fiber_val = 0;
0742 if (IDX < details.long_fiber.size()) {
0743 long_fiber_val = details.long_fiber[IDX];
0744 }
0745 int short_fiber_val = 0;
0746 if (IDX < details.short_fiber.size()) {
0747 short_fiber_val = details.short_fiber[IDX];
0748 }
0749 output[ibin] += (long_fiber_val + short_fiber_val);
0750
0751 uint32_t ADCLong = details.LongDigi[ibin].adc();
0752 uint32_t ADCShort = details.ShortDigi[ibin].adc();
0753
0754 if (details.LongDigi.id().ietaAbs() >= FIRST_FINEGRAIN_TOWER) {
0755 finegrain[ibin][1] = (ADCLong > FG_HF_thresholds_[0] || ADCShort > FG_HF_thresholds_[0]);
0756
0757 if (embit != nullptr)
0758 finegrain[ibin][0] = embit->fineGrainbit(details.ShortDigi, details.LongDigi, ibin);
0759 }
0760 }
0761 }
0762
0763 for (int bin = 0; bin < numberOfSamplesHF_; ++bin) {
0764 static const unsigned int MAX_OUTPUT = QIE8_LINEARIZATION_ET;
0765 output[bin] = min({MAX_OUTPUT, output[bin] >> hf_lumi_shift});
0766 }
0767
0768 std::vector<int> finegrain_converted;
0769 finegrain_converted.reserve(finegrain.size());
0770 for (const auto& fg : finegrain)
0771 finegrain_converted.push_back(fg.to_ulong());
0772 outcoder_->compress(output, finegrain_converted, result);
0773 }
0774
0775 bool HcalTriggerPrimitiveAlgo::passTDC(const QIE10DataFrame& digi, int ts) const {
0776 auto parameters = conditions_->getHcalTPParameters();
0777 auto adc_threshold = parameters->getADCThresholdHF();
0778 auto tdc_mask = parameters->getTDCMaskHF();
0779
0780 if (override_adc_hf_)
0781 adc_threshold = override_adc_hf_value_;
0782 if (override_tdc_hf_)
0783 tdc_mask = override_tdc_hf_value_;
0784
0785 if (digi[ts].adc() < adc_threshold)
0786 return true;
0787
0788 return (1ul << digi[ts].le_tdc()) & tdc_mask;
0789 }
0790
0791 bool HcalTriggerPrimitiveAlgo::validChannel(const QIE10DataFrame& digi, int ts) const {
0792
0793 if (digi.linkError() || ts >= digi.samples() || !digi[ts].ok())
0794 return false;
0795
0796 auto mask = conditions_->getHcalTPChannelParameter(HcalDetId(digi.id()))->getMask();
0797 if (mask)
0798 return false;
0799
0800 return true;
0801 }
0802
0803 void HcalTriggerPrimitiveAlgo::analyzeHFQIE10(const IntegerCaloSamples& samples,
0804 HcalTriggerPrimitiveDigi& result,
0805 const int hf_lumi_shift,
0806 const HcalFeatureBit* embit) {
0807
0808 const int shift = samples.presamples() - numberOfPresamplesHF_;
0809 assert(shift >= 0);
0810 assert((shift + numberOfSamplesHF_) <= samples.size());
0811 assert(hf_lumi_shift >= 2);
0812
0813
0814 const HcalTrigTowerDetId detId(samples.id());
0815 auto it = theHFUpgradeDetailMap.find(detId);
0816
0817 if (it == theHFUpgradeDetailMap.end()) {
0818 return;
0819 }
0820
0821 std::vector<std::bitset<2>> finegrain(numberOfSamplesHF_, false);
0822
0823
0824 IntegerCaloSamples output(samples.id(), numberOfSamplesHF_);
0825 output.setPresamples(numberOfPresamplesHF_);
0826
0827 for (const auto& item : it->second) {
0828 auto& details = item.second;
0829 for (int ibin = 0; ibin < numberOfSamplesHF_; ++ibin) {
0830 const int idx = ibin + shift;
0831
0832 int long_fiber_val = 0;
0833 int long_fiber_count = 0;
0834 int short_fiber_val = 0;
0835 int short_fiber_count = 0;
0836
0837 bool saturated = false;
0838
0839 for (auto i : {0, 2}) {
0840 if (idx < details[i].samples.size() and details[i].validity[idx] and details[i].passTDC[idx]) {
0841 long_fiber_val += details[i].samples[idx];
0842 saturated = saturated || (details[i].samples[idx] == QIE10_LINEARIZATION_ET);
0843 ++long_fiber_count;
0844 }
0845 }
0846 for (auto i : {1, 3}) {
0847 if (idx < details[i].samples.size() and details[i].validity[idx] and details[i].passTDC[idx]) {
0848 short_fiber_val += details[i].samples[idx];
0849 saturated = saturated || (details[i].samples[idx] == QIE10_LINEARIZATION_ET);
0850 ++short_fiber_count;
0851 }
0852 }
0853
0854 if (saturated) {
0855 output[ibin] = QIE10_MAX_LINEARIZATION_ET;
0856 } else {
0857
0858
0859
0860
0861 if (long_fiber_count == 2)
0862 long_fiber_val >>= 1;
0863 if (short_fiber_count == 2)
0864 short_fiber_val >>= 1;
0865
0866 auto sum = long_fiber_val + short_fiber_val;
0867
0868
0869
0870 if (long_fiber_count > 0 and short_fiber_count > 0)
0871 sum >>= 1;
0872
0873 output[ibin] += sum;
0874 }
0875
0876 for (const auto& detail : details) {
0877 if (idx < int(detail.digi.size()) and detail.validity[idx] and
0878 HcalDetId(detail.digi.id()).ietaAbs() >= FIRST_FINEGRAIN_TOWER) {
0879 if (useTDCInMinBiasBits_ && !detail.passTDC[idx])
0880 continue;
0881 finegrain[ibin][1] = finegrain[ibin][1] or detail.fgbits[idx][0];
0882
0883
0884
0885 finegrain[ibin][0] = finegrain[ibin][0] or detail.fgbits[idx][1];
0886 }
0887 }
0888
0889 if (embit != nullptr and FG_HF_thresholds_.at(1) != 255) {
0890 finegrain[ibin][0] = embit->fineGrainbit(details[1].digi,
0891 details[3].digi,
0892 details[0].digi,
0893 details[2].digi,
0894 details[1].validity[idx],
0895 details[3].validity[idx],
0896 details[0].validity[idx],
0897 details[2].validity[idx],
0898 idx);
0899 }
0900 }
0901 }
0902
0903 for (int bin = 0; bin < numberOfSamplesHF_; ++bin) {
0904 output[bin] = min({(unsigned int)QIE10_MAX_LINEARIZATION_ET, output[bin] >> (hf_lumi_shift - 2)});
0905 }
0906 std::vector<int> finegrain_converted;
0907 finegrain_converted.reserve(finegrain.size());
0908 for (const auto& fg : finegrain)
0909 finegrain_converted.push_back(fg.to_ulong());
0910 outcoder_->compress(output, finegrain_converted, result);
0911 }
0912
0913 void HcalTriggerPrimitiveAlgo::runZS(HcalTrigPrimDigiCollection& result) {
0914 for (HcalTrigPrimDigiCollection::iterator tp = result.begin(); tp != result.end(); ++tp) {
0915 bool ZS = true;
0916 for (int i = 0; i < tp->size(); ++i) {
0917 if (tp->sample(i).compressedEt() > ZS_threshold_I_) {
0918 ZS = false;
0919 break;
0920 }
0921 }
0922 if (ZS)
0923 tp->setZSInfo(false, true);
0924 else
0925 tp->setZSInfo(true, false);
0926 }
0927 }
0928
0929 void HcalTriggerPrimitiveAlgo::runFEFormatError(const FEDRawDataCollection* rawraw,
0930 const HcalElectronicsMap* emap,
0931 HcalTrigPrimDigiCollection& result) {
0932 std::set<uint32_t> FrontEndErrors;
0933
0934 for (int i = FEDNumbering::MINHCALFEDID; i <= FEDNumbering::MAXHCALFEDID; ++i) {
0935 const FEDRawData& raw = rawraw->FEDData(i);
0936 if (raw.size() < 12)
0937 continue;
0938 const HcalDCCHeader* dccHeader = (const HcalDCCHeader*)(raw.data());
0939 if (!dccHeader)
0940 continue;
0941 HcalHTRData htr;
0942 for (int spigot = 0; spigot < HcalDCCHeader::SPIGOT_COUNT; spigot++) {
0943 if (!dccHeader->getSpigotPresent(spigot))
0944 continue;
0945 dccHeader->getSpigotData(spigot, htr, raw.size());
0946 int dccid = dccHeader->getSourceId();
0947 int errWord = htr.getErrorsWord() & 0x1FFFF;
0948 bool HTRError = (!htr.check() || htr.isHistogramEvent() || (errWord & 0x800) != 0);
0949
0950 if (HTRError) {
0951 bool valid = false;
0952 for (int fchan = 0; fchan < 3 && !valid; fchan++) {
0953 for (int fib = 0; fib < 9 && !valid; fib++) {
0954 HcalElectronicsId eid(fchan, fib, spigot, dccid - FEDNumbering::MINHCALFEDID);
0955 eid.setHTR(htr.readoutVMECrateId(), htr.htrSlot(), htr.htrTopBottom());
0956 DetId detId = emap->lookup(eid);
0957 if (detId.null())
0958 continue;
0959 HcalSubdetector subdet = (HcalSubdetector(detId.subdetId()));
0960 if (detId.det() != 4 || (subdet != HcalBarrel && subdet != HcalEndcap && subdet != HcalForward))
0961 continue;
0962 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(detId);
0963 for (std::vector<HcalTrigTowerDetId>::const_iterator triggerId = ids.begin(); triggerId != ids.end();
0964 ++triggerId) {
0965 FrontEndErrors.insert(triggerId->rawId());
0966 }
0967
0968 }
0969 }
0970 }
0971 }
0972 }
0973
0974
0975
0976 HcalTriggerPrimitiveSample zeroSample(0);
0977 for (HcalTrigPrimDigiCollection::iterator tp = result.begin(); tp != result.end(); ++tp) {
0978 if (FrontEndErrors.find(tp->id().rawId()) != FrontEndErrors.end()) {
0979 for (int i = 0; i < tp->size(); ++i)
0980 tp->setSample(i, zeroSample);
0981 }
0982 }
0983 }
0984
0985 void HcalTriggerPrimitiveAlgo::addFG(const HcalTrigTowerDetId& id, std::vector<bool>& msb) {
0986 FGbitMap::iterator itr = fgMap_.find(id);
0987 if (itr != fgMap_.end()) {
0988 std::vector<bool>& _msb = itr->second;
0989 for (size_t i = 0; i < msb.size(); ++i)
0990 _msb[i] = _msb[i] || msb[i];
0991 } else
0992 fgMap_[id] = msb;
0993 }
0994
0995 bool HcalTriggerPrimitiveAlgo::validUpgradeFG(const HcalTrigTowerDetId& id, int depth) const {
0996 if (depth > LAST_FINEGRAIN_DEPTH)
0997 return false;
0998 if (id.ietaAbs() > LAST_FINEGRAIN_TOWER)
0999 return false;
1000 if (id.ietaAbs() == HBHE_OVERLAP_TOWER and not upgrade_hb_)
1001 return false;
1002 return true;
1003 }
1004
1005 bool HcalTriggerPrimitiveAlgo::needLegacyFG(const HcalTrigTowerDetId& id) const {
1006
1007
1008
1009 if (id.ietaAbs() == HBHE_OVERLAP_TOWER and not upgrade_hb_)
1010 return true;
1011 return false;
1012 }
1013
1014 bool HcalTriggerPrimitiveAlgo::needUpgradeID(const HcalTrigTowerDetId& id, int depth) const {
1015
1016
1017
1018 unsigned int aieta = id.ietaAbs();
1019 if (aieta >= FIRST_DEPTH7_TOWER and aieta <= LAST_FINEGRAIN_TOWER and depth > LAST_FINEGRAIN_DEPTH)
1020 return true;
1021 return false;
1022 }
1023
1024 void HcalTriggerPrimitiveAlgo::addUpgradeFG(const HcalTrigTowerDetId& id,
1025 int depth,
1026 const std::vector<std::bitset<2>>& bits) {
1027 if (not validUpgradeFG(id, depth)) {
1028 if (needLegacyFG(id)) {
1029 std::vector<bool> pseudo(bits.size(), false);
1030 addFG(id, pseudo);
1031 } else if (needUpgradeID(id, depth)) {
1032
1033
1034
1035 auto it = fgUpgradeMap_.find(id);
1036 if (it == fgUpgradeMap_.end()) {
1037 FGUpgradeContainer element;
1038 element.resize(bits.size());
1039 fgUpgradeMap_.insert(std::make_pair(id, element));
1040 }
1041 }
1042
1043 return;
1044 }
1045
1046 auto it = fgUpgradeMap_.find(id);
1047 if (it == fgUpgradeMap_.end()) {
1048 FGUpgradeContainer element;
1049 element.resize(bits.size());
1050 it = fgUpgradeMap_.insert(std::make_pair(id, element)).first;
1051 }
1052 for (unsigned int i = 0; i < bits.size(); ++i) {
1053 it->second[i][0][depth - 1] = bits[i][0];
1054 it->second[i][1][depth - 1] = bits[i][1];
1055 }
1056 }
1057
1058 void HcalTriggerPrimitiveAlgo::addUpgradeTDCFG(const HcalTrigTowerDetId& id, const QIE11DataFrame& frame) {
1059 HcalDetId detId(frame.id());
1060 if (detId.subdet() != HcalEndcap && detId.subdet() != HcalBarrel)
1061 return;
1062
1063 std::vector<HcalTrigTowerDetId> ids = theTrigTowerGeometry->towerIds(detId);
1064 assert(ids.size() == 1 || ids.size() == 2);
1065 IntegerCaloSamples samples1(ids[0], int(frame.samples()));
1066 samples1.setPresamples(frame.presamples());
1067 incoder_->adc2Linear(frame, samples1);
1068 std::vector<unsigned short> bits12_15 = incoder_->group0FGbits(frame);
1069
1070 bool is_compressed = false;
1071 if (detId.subdet() == HcalBarrel) {
1072 is_compressed = (frame.flavor() == 3);
1073
1074 }
1075
1076 auto it = fgUpgradeTDCMap_.find(id);
1077 if (it == fgUpgradeTDCMap_.end()) {
1078 FGUpgradeTDCContainer element;
1079 element.resize(frame.samples());
1080 it = fgUpgradeTDCMap_.insert(std::make_pair(id, element)).first;
1081 }
1082 for (int i = 0; i < frame.samples(); i++) {
1083 it->second[i][detId.depth() - 1] =
1084 std::make_pair(std::make_pair(bits12_15[i], is_compressed), std::make_pair(frame[i].tdc(), samples1[i]));
1085 }
1086 }
1087
1088 void HcalTriggerPrimitiveAlgo::setWeightsQIE11(const edm::ParameterSet& weightsQIE11) {
1089
1090 std::vector<std::string> ietaStrs = weightsQIE11.getParameterNames();
1091 for (auto& ietaStr : ietaStrs) {
1092
1093 auto const& v = weightsQIE11.getParameter<std::vector<int>>(ietaStr);
1094 weightsQIE11_[std::stoi(ietaStr.substr(4))] = {{v[0], v[1]}};
1095 }
1096 }
1097
1098 void HcalTriggerPrimitiveAlgo::setWeightQIE11(int aieta, int weight) {
1099
1100
1101 weightsQIE11_[aieta] = {{weight, 255}};
1102 }
1103
1104 void HcalTriggerPrimitiveAlgo::setCodedVetoThresholds(const edm::ParameterSet& codedVetoThresholds) {
1105
1106 std::vector<std::string> ietaStrs = codedVetoThresholds.getParameterNames();
1107 for (auto& ietaStr : ietaStrs) {
1108
1109 auto const& v = codedVetoThresholds.getParameter<int>(ietaStr);
1110 codedVetoThresholds_[std::stoi(ietaStr.substr(4))] = {v};
1111 }
1112 }
1113
1114 void HcalTriggerPrimitiveAlgo::setCodedVetoThreshold(int aieta, int codedVetoThreshold) {
1115
1116 codedVetoThresholds_[aieta] = {codedVetoThreshold};
1117 }
1118
1119 void HcalTriggerPrimitiveAlgo::setPeakFinderAlgorithm(int algo) {
1120 if (algo <= 0 || algo > 2)
1121 throw cms::Exception("ERROR: Only algo 1 & 2 are supported.") << std::endl;
1122 peak_finder_algorithm_ = algo;
1123 }
1124
1125 void HcalTriggerPrimitiveAlgo::setNCTScaleShift(int shift) { NCTScaleShift = shift; }
1126
1127 void HcalTriggerPrimitiveAlgo::setRCTScaleShift(int shift) { RCTScaleShift = shift; }