Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-11-21 00:28:29

0001 #include <iostream>
0002 #include <fstream>
0003 #include <cmath>
0004 #include <string>
0005 #include "CalibCalorimetry/HcalTPGAlgos/interface/HcaluLUTTPGCoder.h"
0006 #include "CalibFormats/HcalObjects/interface/HcalCoderDb.h"
0007 #include "CalibFormats/HcalObjects/interface/HcalCalibrations.h"
0008 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
0009 #include "DataFormats/HcalDigi/interface/HcalQIENum.h"
0010 #include "DataFormats/HcalDigi/interface/QIE10DataFrame.h"
0011 #include "DataFormats/HcalDigi/interface/QIE11DataFrame.h"
0012 #include "DataFormats/HcalDetId/interface/HcalTrigTowerDetId.h"
0013 #include "Geometry/HcalTowerAlgo/interface/HcalTrigTowerGeometry.h"
0014 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0015 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0016 #include "FWCore/Utilities/interface/Exception.h"
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 #include "FWCore/Framework/interface/Event.h"
0019 #include "FWCore/Framework/interface/EventSetup.h"
0020 #include "FWCore/Framework/interface/Frameworkfwd.h"
0021 #include "FWCore/Framework/interface/MakerMacros.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "CondFormats/HcalObjects/interface/HcalL1TriggerObjects.h"
0024 #include "CondFormats/HcalObjects/interface/HcalL1TriggerObject.h"
0025 #include "CalibCalorimetry/HcalAlgos/interface/HcalDbASCIIIO.h"
0026 #include "CalibCalorimetry/HcalAlgos/interface/HcalSiPMnonlinearity.h"
0027 #include "CalibCalorimetry/HcalAlgos/interface/HcalPulseContainmentCorrection.h"
0028 #include "CalibCalorimetry/HcalTPGAlgos/interface/XMLProcessor.h"
0029 #include "CalibCalorimetry/HcalTPGAlgos/interface/LutXml.h"
0030 
0031 const float HcaluLUTTPGCoder::lsb_ = 1. / 16;
0032 const float HcaluLUTTPGCoder::zdc_lsb_ = 50.;
0033 
0034 const int HcaluLUTTPGCoder::QIE8_LUT_BITMASK;
0035 const int HcaluLUTTPGCoder::QIE10_LUT_BITMASK;
0036 const int HcaluLUTTPGCoder::QIE11_LUT_BITMASK;
0037 const int HcaluLUTTPGCoder::QIE10_ZDC_LUT_BITMASK;
0038 
0039 constexpr double MaximumFractionalError = 0.002;  // 0.2% error allowed from this source
0040 
0041 HcaluLUTTPGCoder::HcaluLUTTPGCoder()
0042     : topo_{},
0043       delay_{},
0044       LUTGenerationMode_{},
0045       FG_HF_thresholds_{},
0046       bitToMask_{},
0047       firstHBEta_{},
0048       lastHBEta_{},
0049       nHBEta_{},
0050       maxDepthHB_{},
0051       sizeHB_{},
0052       firstHEEta_{},
0053       lastHEEta_{},
0054       nHEEta_{},
0055       maxDepthHE_{},
0056       sizeHE_{},
0057       firstHFEta_{},
0058       lastHFEta_{},
0059       nHFEta_{},
0060       maxDepthHF_{},
0061       sizeHF_{},
0062       sizeZDC_{},
0063       cosh_ieta_28_HE_low_depths_{},
0064       cosh_ieta_28_HE_high_depths_{},
0065       cosh_ieta_29_HE_{},
0066       allLinear_{},
0067       contain1TSHB_{},
0068       contain1TSHE_{},
0069       applyFixPCC_{},
0070       linearLSB_QIE8_{},
0071       linearLSB_QIE11_{},
0072       linearLSB_QIE11Overlap_{} {}
0073 
0074 HcaluLUTTPGCoder::HcaluLUTTPGCoder(const HcalTopology* top, const HcalTimeSlew* delay) { init(top, delay); }
0075 
0076 void HcaluLUTTPGCoder::init(const HcalTopology* top, const HcalTimeSlew* delay) {
0077   topo_ = top;
0078   delay_ = delay;
0079   LUTGenerationMode_ = true;
0080   FG_HF_thresholds_ = {0, 0};
0081   bitToMask_ = 0;
0082   allLinear_ = false;
0083   contain1TSHB_ = false;
0084   contain1TSHE_ = false;
0085   applyFixPCC_ = false;
0086   linearLSB_QIE8_ = 1.;
0087   linearLSB_QIE11_ = 1.;
0088   linearLSB_QIE11Overlap_ = 1.;
0089   pulseCorr_ = std::make_unique<HcalPulseContainmentManager>(MaximumFractionalError, false);
0090   firstHBEta_ = topo_->firstHBRing();
0091   lastHBEta_ = topo_->lastHBRing();
0092   nHBEta_ = (lastHBEta_ - firstHBEta_ + 1);
0093   maxDepthHB_ = topo_->maxDepth(HcalBarrel);
0094   sizeHB_ = 2 * nHBEta_ * nFi_ * maxDepthHB_;
0095   firstHEEta_ = topo_->firstHERing();
0096   lastHEEta_ = topo_->lastHERing();
0097   nHEEta_ = (lastHEEta_ - firstHEEta_ + 1);
0098   maxDepthHE_ = topo_->maxDepth(HcalEndcap);
0099   sizeHE_ = 2 * nHEEta_ * nFi_ * maxDepthHE_;
0100   firstHFEta_ = topo_->firstHFRing();
0101   lastHFEta_ = topo_->lastHFRing();
0102   nHFEta_ = (lastHFEta_ - firstHFEta_ + 1);
0103   maxDepthHF_ = topo_->maxDepth(HcalForward);
0104   sizeHF_ = 2 * nHFEta_ * nFi_ * maxDepthHF_;
0105   sizeZDC_ = 2 * 9;
0106   size_t nluts = (size_t)(sizeHB_ + sizeHE_ + sizeHF_ + sizeZDC_ * 2 + 1);
0107   inputLUT_ = std::vector<HcaluLUTTPGCoder::Lut>(nluts);
0108   gain_ = std::vector<float>(nluts, 0.);
0109   ped_ = std::vector<float>(nluts, 0.);
0110   make_cosh_ieta_map();
0111 }
0112 
0113 void HcaluLUTTPGCoder::compress(const IntegerCaloSamples& ics,
0114                                 const std::vector<bool>& featureBits,
0115                                 HcalTriggerPrimitiveDigi& tp) const {
0116   throw cms::Exception("PROBLEM: This method should never be invoked!");
0117 }
0118 
0119 HcaluLUTTPGCoder::~HcaluLUTTPGCoder() {}
0120 
0121 int HcaluLUTTPGCoder::getLUTId(HcalSubdetector id, int ieta, int iphi, int depth) const {
0122   int retval(0);
0123   if (id == HcalBarrel) {
0124     retval = (depth - 1) + maxDepthHB_ * (iphi - 1);
0125     if (ieta > 0)
0126       retval += maxDepthHB_ * nFi_ * (ieta - firstHBEta_);
0127     else
0128       retval += maxDepthHB_ * nFi_ * (ieta + lastHBEta_ + nHBEta_);
0129   } else if (id == HcalEndcap) {
0130     retval = sizeHB_;
0131     retval += (depth - 1) + maxDepthHE_ * (iphi - 1);
0132     if (ieta > 0)
0133       retval += maxDepthHE_ * nFi_ * (ieta - firstHEEta_);
0134     else
0135       retval += maxDepthHE_ * nFi_ * (ieta + lastHEEta_ + nHEEta_);
0136   } else if (id == HcalForward) {
0137     retval = sizeHB_ + sizeHE_;
0138     retval += (depth - 1) + maxDepthHF_ * (iphi - 1);
0139     if (ieta > 0)
0140       retval += maxDepthHF_ * nFi_ * (ieta - firstHFEta_);
0141     else
0142       retval += maxDepthHF_ * nFi_ * (ieta + lastHFEta_ + nHFEta_);
0143   } else if (id == HcalOther) {
0144     retval = sizeHB_ + sizeHE_ + sizeHF_;
0145     retval += depth;
0146     if (iphi != 1)
0147       retval += 5;
0148     if (ieta <= 0)
0149       retval += 9;
0150   }
0151   return retval;
0152 }
0153 
0154 int HcaluLUTTPGCoder::getLUTId(uint32_t rawid) const {
0155   HcalDetId detid(rawid);
0156   return getLUTId(detid.subdet(), detid.ieta(), detid.iphi(), detid.depth());
0157 }
0158 
0159 int HcaluLUTTPGCoder::getLUTId(const HcalDetId& detid) const {
0160   return getLUTId(detid.subdet(), detid.ieta(), detid.iphi(), detid.depth());
0161 }
0162 
0163 int HcaluLUTTPGCoder::getLUTId(const HcalZDCDetId& detid) const {
0164   if (detid.section() == HcalZDCDetId::EM) {
0165     return getLUTId(HcalOther, detid.zside(), 1, detid.channel());
0166   } else if (detid.section() == HcalZDCDetId::HAD) {
0167     return getLUTId(HcalOther, detid.zside(), 2, detid.channel());
0168   } else if (detid.section() == HcalZDCDetId::LUM) {
0169     return getLUTId(HcalOther, detid.zside(), 3, detid.channel());
0170   } else {
0171     return getLUTId(HcalOther, detid.zside(), 4, detid.channel());
0172   }
0173 }
0174 
0175 void HcaluLUTTPGCoder::update(const char* filename, bool appendMSB) {
0176   std::ifstream file(filename, std::ios::in);
0177   assert(file.is_open());
0178 
0179   std::vector<HcalSubdetector> subdet;
0180   std::string buffer;
0181 
0182   // Drop first (comment) line
0183   std::getline(file, buffer);
0184   std::getline(file, buffer);
0185 
0186   unsigned int index = buffer.find('H', 0);
0187   while (index < buffer.length()) {
0188     std::string subdetStr = buffer.substr(index, 2);
0189     if (subdetStr == "HB")
0190       subdet.push_back(HcalBarrel);
0191     else if (subdetStr == "HE")
0192       subdet.push_back(HcalEndcap);
0193     else if (subdetStr == "HF")
0194       subdet.push_back(HcalForward);
0195     //TODO Check subdet
0196     //else exception
0197     index += 2;
0198     index = buffer.find('H', index);
0199   }
0200 
0201   // Get upper/lower ranges for ieta/iphi/depth
0202   size_t nCol = subdet.size();
0203   assert(nCol > 0);
0204 
0205   std::vector<int> ietaU;
0206   std::vector<int> ietaL;
0207   std::vector<int> iphiU;
0208   std::vector<int> iphiL;
0209   std::vector<int> depU;
0210   std::vector<int> depL;
0211   std::vector<Lut> lutFromFile(nCol);
0212   LutElement lutValue;
0213 
0214   for (size_t i = 0; i < nCol; ++i) {
0215     int ieta;
0216     file >> ieta;
0217     ietaL.push_back(ieta);
0218   }
0219 
0220   for (size_t i = 0; i < nCol; ++i) {
0221     int ieta;
0222     file >> ieta;
0223     ietaU.push_back(ieta);
0224   }
0225 
0226   for (size_t i = 0; i < nCol; ++i) {
0227     int iphi;
0228     file >> iphi;
0229     iphiL.push_back(iphi);
0230   }
0231 
0232   for (size_t i = 0; i < nCol; ++i) {
0233     int iphi;
0234     file >> iphi;
0235     iphiU.push_back(iphi);
0236   }
0237 
0238   for (size_t i = 0; i < nCol; ++i) {
0239     int dep;
0240     file >> dep;
0241     depL.push_back(dep);
0242   }
0243 
0244   for (size_t i = 0; i < nCol; ++i) {
0245     int dep;
0246     file >> dep;
0247     depU.push_back(dep);
0248   }
0249 
0250   // Read Lut Entry
0251   for (size_t i = 0; file >> lutValue; i = (i + 1) % nCol) {
0252     lutFromFile[i].push_back(lutValue);
0253   }
0254 
0255   // Check lut size
0256   for (size_t i = 0; i < nCol; ++i)
0257     assert(lutFromFile[i].size() == INPUT_LUT_SIZE);
0258 
0259   for (size_t i = 0; i < nCol; ++i) {
0260     for (int ieta = ietaL[i]; ieta <= ietaU[i]; ++ieta) {
0261       for (int iphi = iphiL[i]; iphi <= iphiU[i]; ++iphi) {
0262         for (int depth = depL[i]; depth <= depU[i]; ++depth) {
0263           HcalDetId id(subdet[i], ieta, iphi, depth);
0264           if (!topo_->valid(id))
0265             continue;
0266 
0267           int lutId = getLUTId(id);
0268           for (size_t adc = 0; adc < INPUT_LUT_SIZE; ++adc) {
0269             if (appendMSB) {
0270               // Append FG bit LUT to MSB
0271               // MSB = Most Significant Bit = bit 10
0272               // Overwrite bit 10
0273               LutElement msb = (lutFromFile[i][adc] != 0 ? QIE8_LUT_MSB : 0);
0274               inputLUT_[lutId][adc] = (msb | (inputLUT_[lutId][adc] & QIE8_LUT_BITMASK));
0275             } else
0276               inputLUT_[lutId][adc] = lutFromFile[i][adc];
0277           }  // for adc
0278         }    // for depth
0279       }      // for iphi
0280     }        // for ieta
0281   }          // for nCol
0282 }
0283 
0284 void HcaluLUTTPGCoder::updateXML(const char* filename) {
0285   LutXml* _xml = new LutXml(filename);
0286   _xml->create_lut_map();
0287   HcalSubdetector subdet[3] = {HcalBarrel, HcalEndcap, HcalForward};
0288   for (int ieta = -HcalDetId::kHcalEtaMask2; ieta <= (int)(HcalDetId::kHcalEtaMask2); ++ieta) {
0289     for (unsigned int iphi = 0; iphi <= HcalDetId::kHcalPhiMask2; ++iphi) {
0290       for (unsigned int depth = 1; depth < HcalDetId::kHcalDepthMask2; ++depth) {
0291         for (int isub = 0; isub < 3; ++isub) {
0292           HcalDetId detid(subdet[isub], ieta, iphi, depth);
0293           if (!topo_->valid(detid))
0294             continue;
0295           int id = getLUTId(subdet[isub], ieta, iphi, depth);
0296           std::vector<unsigned int>* lut = _xml->getLutFast(detid);
0297           if (lut == nullptr)
0298             throw cms::Exception("PROBLEM: No inputLUT_ in xml file for ") << detid << std::endl;
0299           if (lut->size() != UPGRADE_LUT_SIZE)
0300             throw cms::Exception("PROBLEM: Wrong inputLUT_ size in xml file for ") << detid << std::endl;
0301           Lut& lutRaw = inputLUT_[id];
0302           lutRaw.resize(UPGRADE_LUT_SIZE, 0);
0303           for (unsigned int i = 0; i < UPGRADE_LUT_SIZE; ++i)
0304             inputLUT_[id][i] = (LutElement)lut->at(i);
0305         }
0306       }
0307     }
0308   }
0309   delete _xml;
0310   XMLProcessor::getInstance()->terminate();
0311 }
0312 
0313 double HcaluLUTTPGCoder::cosh_ieta(int ieta, int depth, HcalSubdetector subdet) {
0314   // ieta = 28 and 29 are both associated with trigger tower 28
0315   // so special handling is required. HF ieta=29 channels included in TT30
0316   // are already handled correctly in cosh_ieta_
0317   if (abs(ieta) >= 28 && subdet == HcalEndcap && allLinear_) {
0318     if (abs(ieta) == 29)
0319       return cosh_ieta_29_HE_;
0320     if (abs(ieta) == 28) {
0321       if (depth <= 3)
0322         return cosh_ieta_28_HE_low_depths_;
0323       else
0324         return cosh_ieta_28_HE_high_depths_;
0325     }
0326   }
0327 
0328   return cosh_ieta_[ieta];
0329 }
0330 
0331 void HcaluLUTTPGCoder::make_cosh_ieta_map(void) {
0332   cosh_ieta_ = std::vector<double>(lastHFEta_ + 1, -1.0);
0333 
0334   HcalTrigTowerGeometry triggeo(topo_);
0335 
0336   for (int i = 1; i <= firstHFEta_; ++i) {
0337     double eta_low = 0., eta_high = 0.;
0338     triggeo.towerEtaBounds(i, 0, eta_low, eta_high);
0339     cosh_ieta_[i] = cosh((eta_low + eta_high) / 2.);
0340   }
0341   for (int i = firstHFEta_; i <= lastHFEta_; ++i) {
0342     std::pair<double, double> etas = topo_->etaRange(HcalForward, i);
0343     double eta1 = etas.first;
0344     double eta2 = etas.second;
0345     cosh_ieta_[i] = cosh((eta1 + eta2) / 2.);
0346   }
0347 
0348   // trigger tower 28 in HE has a more complicated geometry
0349   std::pair<double, double> eta28 = topo_->etaRange(HcalEndcap, 28);
0350   std::pair<double, double> eta29 = topo_->etaRange(HcalEndcap, 29);
0351   cosh_ieta_29_HE_ = cosh((eta29.first + eta29.second) / 2.);
0352   cosh_ieta_28_HE_low_depths_ = cosh((eta28.first + eta28.second) / 2.);
0353   // for higher depths in ieta = 28, the trigger tower extends past
0354   // the ieta = 29 channels
0355   cosh_ieta_28_HE_high_depths_ = cosh((eta28.first + eta29.second) / 2.);
0356 }
0357 
0358 void HcaluLUTTPGCoder::update(const HcalDbService& conditions) {
0359   const HcalLutMetadata* metadata = conditions.getHcalLutMetadata();
0360   assert(metadata != nullptr);
0361   float nominalgain_ = metadata->getNominalGain();
0362 
0363   pulseCorr_ = std::make_unique<HcalPulseContainmentManager>(MaximumFractionalError, applyFixPCC_);
0364   pulseCorr_->beginRun(&conditions, delay_);
0365 
0366   make_cosh_ieta_map();
0367 
0368   // Here we will determine if we are using new version of TPs (1TS)
0369   // i.e. are we using a new pulse filter scheme.
0370   const HcalElectronicsMap* emap = conditions.getHcalMapping();
0371 
0372   int lastHBRing = topo_->lastHBRing();
0373   int lastHERing = topo_->lastHERing();
0374 
0375   // First, determine if we should configure for the filter scheme
0376   // Check the tp version to make this determination
0377   bool foundHB = false;
0378   bool foundHE = false;
0379   bool newHBtp = false;
0380   bool newHEtp = false;
0381   std::vector<HcalElectronicsId> vIds = emap->allElectronicsIdTrigger();
0382   for (std::vector<HcalElectronicsId>::const_iterator eId = vIds.begin(); eId != vIds.end(); eId++) {
0383     // The first HB or HE id is enough to tell whether to use new scheme in HB or HE
0384     if (foundHB and foundHE)
0385       break;
0386 
0387     HcalTrigTowerDetId hcalTTDetId(emap->lookupTrigger(*eId));
0388     if (hcalTTDetId.null())
0389       continue;
0390 
0391     int aieta = abs(hcalTTDetId.ieta());
0392     // The absence of TT channels in the HcalTPChannelParameters
0393     // is intepreted as to not use the new filter
0394     int weight = -1.0;
0395     auto tpParam = conditions.getHcalTPChannelParameter(hcalTTDetId, false);
0396     if (tpParam)
0397       weight = tpParam->getauxi1();
0398 
0399     if (aieta <= lastHBRing) {
0400       foundHB = true;
0401       if (weight != -1.0)
0402         newHBtp = true;
0403     } else if (aieta > lastHBRing and aieta < lastHERing) {
0404       foundHE = true;
0405       if (weight != -1.0)
0406         newHEtp = true;
0407     }
0408   }
0409 
0410   for (const auto& id : metadata->getAllChannels()) {
0411     if (id.det() == DetId::Hcal and topo_->valid(id)) {
0412       HcalDetId cell(id);
0413       HcalSubdetector subdet = cell.subdet();
0414 
0415       if (subdet != HcalBarrel and subdet != HcalEndcap and subdet != HcalForward and subdet != HcalOther)
0416         continue;
0417 
0418       const HcalQIECoder* channelCoder = conditions.getHcalCoder(cell);
0419       const HcalQIEShape* shape = conditions.getHcalShape(cell);
0420       HcalCoderDb coder(*channelCoder, *shape);
0421       const HcalLutMetadatum* meta = metadata->getValues(cell);
0422 
0423       // defaults for energy requirement for bits 12-15 are high / low to avoid FG bit 0-4 being set when not intended
0424       unsigned int bit12_energy = 0;
0425       unsigned int bit13_energy = 0;
0426       unsigned int bit14_energy = 999;
0427       unsigned int bit15_energy = 999;
0428 
0429       bool is2018OrLater = topo_->triggerMode() >= HcalTopologyMode::TriggerMode_2018 or
0430                            topo_->triggerMode() == HcalTopologyMode::TriggerMode_2018legacy;
0431       if (is2018OrLater or topo_->dddConstants()->isPlan1(cell)) {
0432         bit12_energy = 16;  // depths 1,2 max energy
0433         bit13_energy = 80;  // depths 3+ min energy
0434         bit14_energy = 64;  // prompt min energy
0435         bit15_energy = 64;  // delayed min energy
0436       }
0437 
0438       int lutId = getLUTId(cell);
0439       Lut& lut = inputLUT_[lutId];
0440       float ped = 0;
0441       float gain = 0;
0442       uint32_t status = 0;
0443 
0444       if (LUTGenerationMode_) {
0445         const HcalCalibrations& calibrations = conditions.getHcalCalibrations(cell);
0446         for (auto capId : {0, 1, 2, 3}) {
0447           ped += calibrations.effpedestal(capId);
0448           gain += calibrations.LUTrespcorrgain(capId);
0449         }
0450         ped /= 4.0;
0451         gain /= 4.0;
0452 
0453         //Get Channel Quality
0454         const HcalChannelStatus* channelStatus = conditions.getHcalChannelStatus(cell);
0455         status = channelStatus->getValue();
0456 
0457       } else {
0458         const HcalL1TriggerObject* myL1TObj = conditions.getHcalL1TriggerObject(cell);
0459         ped = myL1TObj->getPedestal();
0460         gain = myL1TObj->getRespGain();
0461         status = myL1TObj->getFlag();
0462       }  // LUTGenerationMode_
0463 
0464       ped_[lutId] = ped;
0465       gain_[lutId] = gain;
0466       bool isMasked = ((status & bitToMask_) > 0);
0467       float rcalib = meta->getRCalib();
0468 
0469       auto adc2fC = [channelCoder, shape](unsigned int adc) {
0470         float fC = 0;
0471         for (auto capId : {0, 1, 2, 3})
0472           fC += channelCoder->charge(*shape, adc, capId);
0473         return fC / 4;
0474       };
0475 
0476       int qieType = conditions.getHcalQIEType(cell)->getValue();
0477 
0478       const size_t SIZE = qieType == QIE8 ? INPUT_LUT_SIZE : UPGRADE_LUT_SIZE;
0479       const int MASK = qieType == QIE8 ? QIE8_LUT_BITMASK : qieType == QIE10 ? QIE10_LUT_BITMASK : QIE11_LUT_BITMASK;
0480       double linearLSB = linearLSB_QIE8_;
0481       if (qieType == QIE11 and cell.ietaAbs() == topo_->lastHBRing())
0482         linearLSB = linearLSB_QIE11Overlap_;
0483       else if (qieType == QIE11)
0484         linearLSB = linearLSB_QIE11_;
0485 
0486       lut.resize(SIZE, 0);
0487 
0488       // Input LUT for HB/HE/HF
0489       if (subdet == HcalBarrel || subdet == HcalEndcap) {
0490         int granularity = meta->getLutGranularity();
0491 
0492         double correctionPhaseNS = conditions.getHcalRecoParam(cell)->correctionPhaseNS();
0493 
0494         if (qieType == QIE11) {
0495           if (overrideDBweightsAndFilterHB_ and cell.ietaAbs() <= lastHBRing)
0496             correctionPhaseNS = containPhaseNSHB_;
0497           else if (overrideDBweightsAndFilterHE_ and cell.ietaAbs() > lastHBRing)
0498             correctionPhaseNS = containPhaseNSHE_;
0499         }
0500         for (unsigned int adc = 0; adc < SIZE; ++adc) {
0501           if (isMasked)
0502             lut[adc] = 0;
0503           else {
0504             double nonlinearityCorrection = 1.0;
0505             double containmentCorrection = 1.0;
0506             // SiPM nonlinearity was not corrected in 2017
0507             // and containment corrections  were not
0508             // ET-dependent prior to 2018
0509             if (is2018OrLater) {
0510               double containmentCorrection1TS = pulseCorr_->correction(cell, 1, correctionPhaseNS, adc2fC(adc));
0511               // Use the 1-TS containment correction to estimate the charge of the pulse
0512               // from the individual samples
0513               double correctedCharge = containmentCorrection1TS * adc2fC(adc);
0514               double containmentCorrection2TSCorrected =
0515                   pulseCorr_->correction(cell, 2, correctionPhaseNS, correctedCharge);
0516               if (qieType == QIE11) {
0517                 // When contain1TS_ is set, it should still only apply for QIE11-related things
0518                 if ((((contain1TSHB_ and overrideDBweightsAndFilterHB_) or newHBtp) and cell.ietaAbs() <= lastHBRing) or
0519                     (((contain1TSHE_ and overrideDBweightsAndFilterHE_) or newHEtp) and cell.ietaAbs() > lastHBRing)) {
0520                   containmentCorrection = containmentCorrection1TS;
0521                 } else {
0522                   containmentCorrection = containmentCorrection2TSCorrected;
0523                 }
0524 
0525                 const HcalSiPMParameter& siPMParameter(*conditions.getHcalSiPMParameter(cell));
0526                 HcalSiPMnonlinearity corr(
0527                     conditions.getHcalSiPMCharacteristics()->getNonLinearities(siPMParameter.getType()));
0528                 const double fcByPE = siPMParameter.getFCByPE();
0529                 const double effectivePixelsFired = correctedCharge / fcByPE;
0530                 nonlinearityCorrection = corr.getRecoCorrectionFactor(effectivePixelsFired);
0531               } else {
0532                 containmentCorrection = containmentCorrection2TSCorrected;
0533               }
0534             }
0535             if (allLinear_)
0536               lut[adc] = (LutElement)std::min(
0537                   std::max(0,
0538                            int((adc2fC(adc) - ped) * gain * rcalib * nonlinearityCorrection * containmentCorrection /
0539                                linearLSB / cosh_ieta(cell.ietaAbs(), cell.depth(), HcalEndcap))),
0540                   MASK);
0541             else
0542               lut[adc] =
0543                   (LutElement)std::min(std::max(0,
0544                                                 int((adc2fC(adc) - ped) * gain * rcalib * nonlinearityCorrection *
0545                                                     containmentCorrection / nominalgain_ / granularity)),
0546                                        MASK);
0547 
0548             unsigned int linearizedADC =
0549                 lut[adc];  // used for bits 12, 13, 14, 15 for Group 0 LUT for LLP time and depth bits that rely on linearized energies
0550 
0551             if (qieType == QIE11) {
0552               if (subdet == HcalBarrel) {  // edit since bits 12-15 not supported in HE yet
0553                 if ((linearizedADC < bit12_energy and cell.depth() <= 2) or (cell.depth() >= 3))
0554                   lut[adc] |= 1 << 12;
0555                 if (linearizedADC >= bit13_energy and cell.depth() >= 3)
0556                   lut[adc] |= 1 << 13;
0557                 if (linearizedADC >= bit14_energy)
0558                   lut[adc] |= 1 << 14;
0559                 if (linearizedADC >= bit15_energy)
0560                   lut[adc] |= 1 << 15;
0561               }
0562             }
0563 
0564             //Zeroing the 4th depth in the trigger towers where |ieta| = 16 to match the behavior in the uHTR firmware in Run3, where the 4th depth is not included in the sum over depths when constructing the TP energy for this tower.
0565             if (abs(cell.ieta()) == 16 && cell.depth() == 4 &&
0566                 topo_->triggerMode() >= HcalTopologyMode::TriggerMode_2021) {
0567               lut[adc] = 0;
0568             }
0569           }
0570         }
0571       } else if (subdet == HcalForward) {
0572         for (unsigned int adc = 0; adc < SIZE; ++adc) {
0573           if (isMasked)
0574             lut[adc] = 0;
0575           else {
0576             lut[adc] = std::min(
0577                 std::max(0, int((adc2fC(adc) - ped) * gain * rcalib / lsb_ / cosh_ieta_[cell.ietaAbs()])), MASK);
0578             if (adc > FG_HF_thresholds_[0])
0579               lut[adc] |= QIE10_LUT_MSB0;
0580             if (adc > FG_HF_thresholds_[1])
0581               lut[adc] |= QIE10_LUT_MSB1;
0582           }
0583         }
0584       }
0585     } else if (id.det() == DetId::Calo && id.subdetId() == HcalZDCDetId::SubdetectorId) {
0586       HcalZDCDetId cell(id.rawId());
0587 
0588       if (cell.section() != HcalZDCDetId::EM && cell.section() != HcalZDCDetId::HAD)
0589         continue;
0590 
0591       const HcalQIECoder* channelCoder = conditions.getHcalCoder(cell);
0592       const HcalQIEShape* shape = conditions.getHcalShape(cell);
0593       const HcalPedestal* effPedestals = conditions.getEffectivePedestal(cell);
0594       HcalCoderDb coder(*channelCoder, *shape);
0595       const HcalLutMetadatum* meta = metadata->getValues(cell);
0596 
0597       auto tpParam = conditions.getHcalTPChannelParameter(cell, false);
0598       int weight = tpParam->getauxi1();
0599 
0600       int lutId = getLUTId(cell);
0601       int lutId_ootpu = lutId + sizeZDC_;
0602       Lut& lut = inputLUT_[lutId];
0603       Lut& lut_ootpu = inputLUT_[lutId_ootpu];
0604       float ped = 0;
0605       float gain = 0;
0606       float pedWidth = 0;
0607       uint32_t status = 0;
0608 
0609       if (LUTGenerationMode_) {
0610         const HcalCalibrations& calibrations = conditions.getHcalCalibrations(cell);
0611         for (auto capId : {0, 1, 2, 3}) {
0612           ped += calibrations.effpedestal(capId);
0613           gain += calibrations.LUTrespcorrgain(capId);
0614           pedWidth += effPedestals->getWidth(capId);
0615         }
0616         ped /= 4.0;
0617         gain /= 4.0;
0618         pedWidth /= 4.0;
0619 
0620         const HcalChannelStatus* channelStatus = conditions.getHcalChannelStatus(cell);
0621         status = channelStatus->getValue();
0622 
0623       } else {
0624         const HcalL1TriggerObject* myL1TObj = conditions.getHcalL1TriggerObject(cell);
0625         ped = myL1TObj->getPedestal();
0626         gain = myL1TObj->getRespGain();
0627         status = myL1TObj->getFlag();
0628       }  // LUTGenerationMode_
0629 
0630       ped_[lutId] = ped;
0631       gain_[lutId] = gain;
0632       bool isMasked = ((status & bitToMask_) > 0);
0633       float rcalib = meta->getRCalib();
0634 
0635       auto adc2fC = [channelCoder, shape](unsigned int adc) {
0636         float fC = 0;
0637         for (auto capId : {0, 1, 2, 3})
0638           fC += channelCoder->charge(*shape, adc, capId);
0639         return fC / 4;
0640       };
0641 
0642       int qieType = conditions.getHcalQIEType(cell)->getValue();
0643 
0644       const size_t SIZE = qieType == QIE8 ? INPUT_LUT_SIZE : UPGRADE_LUT_SIZE;
0645       const int MASK = QIE10_ZDC_LUT_BITMASK;
0646 
0647       lut.resize(SIZE, 0);
0648       lut_ootpu.resize(SIZE, 0);
0649 
0650       for (unsigned int adc = 0; adc < SIZE; ++adc) {
0651         if (isMasked) {
0652           lut[adc] = 0;
0653           lut_ootpu[adc] = 0;
0654         } else {
0655           if ((adc2fC(adc) - ped) < (pedWidth * 5.)) {
0656             lut[adc] = 0;
0657             lut_ootpu[adc] = 0;
0658           } else {
0659             lut[adc] = std::min(std::max(0, int((adc2fC(adc) - ped) * gain * rcalib / zdc_lsb_)), MASK);
0660             lut_ootpu[adc] =
0661                 std::min(std::max(0, int((adc2fC(adc) - ped) * gain * rcalib * weight / (zdc_lsb_ * 256))), MASK);
0662           }
0663         }
0664       }
0665     } else {
0666       continue;
0667     }
0668   }
0669 }
0670 
0671 void HcaluLUTTPGCoder::adc2Linear(const HBHEDataFrame& df, IntegerCaloSamples& ics) const {
0672   int lutId = getLUTId(df.id());
0673   const Lut& lut = inputLUT_.at(lutId);
0674   for (int i = 0; i < df.size(); i++) {
0675     ics[i] = (lut.at(df[i].adc()) & QIE8_LUT_BITMASK);
0676   }
0677 }
0678 
0679 void HcaluLUTTPGCoder::adc2Linear(const HFDataFrame& df, IntegerCaloSamples& ics) const {
0680   int lutId = getLUTId(df.id());
0681   const Lut& lut = inputLUT_.at(lutId);
0682   for (int i = 0; i < df.size(); i++) {
0683     ics[i] = (lut.at(df[i].adc()) & QIE8_LUT_BITMASK);
0684   }
0685 }
0686 
0687 void HcaluLUTTPGCoder::adc2Linear(const QIE10DataFrame& df, IntegerCaloSamples& ics, bool ootpu_lut) const {
0688   DetId detId = DetId(df.detid());
0689   if (detId.det() == DetId::Hcal) {
0690     int lutId = getLUTId(HcalDetId(df.id()));
0691     const Lut& lut = inputLUT_.at(lutId);
0692     for (int i = 0; i < df.samples(); i++) {
0693       ics[i] = (lut.at(df[i].adc()) & QIE10_LUT_BITMASK);
0694     }
0695   } else {
0696     int lutId = getLUTId(HcalZDCDetId(df.id()));
0697     if (ootpu_lut)
0698       lutId = lutId + sizeZDC_;
0699     const Lut& lut = inputLUT_.at(lutId);
0700     for (int i = 0; i < df.samples(); i++) {
0701       ics[i] = (lut.at(df[i].adc()) & QIE10_ZDC_LUT_BITMASK);
0702     }
0703   }
0704 }
0705 
0706 void HcaluLUTTPGCoder::adc2Linear(const QIE11DataFrame& df, IntegerCaloSamples& ics) const {
0707   int lutId = getLUTId(HcalDetId(df.id()));
0708   const Lut& lut = inputLUT_.at(lutId);
0709   for (int i = 0; i < df.samples(); i++) {
0710     ics[i] = (lut.at(df[i].adc()) & QIE11_LUT_BITMASK);
0711   }
0712 }
0713 
0714 unsigned short HcaluLUTTPGCoder::adc2Linear(HcalQIESample sample, HcalDetId id) const {
0715   int lutId = getLUTId(id);
0716   return ((inputLUT_.at(lutId)).at(sample.adc()) & QIE8_LUT_BITMASK);
0717 }
0718 
0719 std::vector<unsigned short> HcaluLUTTPGCoder::group0FGbits(const QIE11DataFrame& df) const {
0720   int lutId = getLUTId(HcalDetId(df.id()));
0721   const Lut& lut = inputLUT_.at(lutId);
0722   std::vector<unsigned short> group0LLPbits;
0723   group0LLPbits.reserve(df.samples());
0724   for (int i = 0; i < df.samples(); i++) {
0725     group0LLPbits.push_back((lut.at(df[i].adc()) >> 12) &
0726                             0xF);  // four bits (12-15) of LUT used to set 6 finegrain bits from uHTR
0727   }
0728   return group0LLPbits;
0729 }
0730 
0731 float HcaluLUTTPGCoder::getLUTPedestal(HcalDetId id) const {
0732   int lutId = getLUTId(id);
0733   return ped_.at(lutId);
0734 }
0735 
0736 float HcaluLUTTPGCoder::getLUTGain(HcalDetId id) const {
0737   int lutId = getLUTId(id);
0738   return gain_.at(lutId);
0739 }
0740 
0741 std::vector<unsigned short> HcaluLUTTPGCoder::getLinearizationLUT(HcalDetId id) const {
0742   int lutId = getLUTId(id);
0743   return inputLUT_.at(lutId);
0744 }
0745 
0746 std::vector<unsigned short> HcaluLUTTPGCoder::getLinearizationLUT(HcalZDCDetId id, bool ootpu_lut) const {
0747   int lutId = getLUTId(id);
0748   if (ootpu_lut)
0749     return inputLUT_.at(lutId + sizeZDC_);
0750   else
0751     return inputLUT_.at(lutId);
0752 }
0753 
0754 void HcaluLUTTPGCoder::lookupMSB(const HBHEDataFrame& df, std::vector<bool>& msb) const {
0755   msb.resize(df.size());
0756   for (int i = 0; i < df.size(); ++i)
0757     msb[i] = getMSB(df.id(), df.sample(i).adc());
0758 }
0759 
0760 bool HcaluLUTTPGCoder::getMSB(const HcalDetId& id, int adc) const {
0761   int lutId = getLUTId(id);
0762   const Lut& lut = inputLUT_.at(lutId);
0763   return (lut.at(adc) & QIE8_LUT_MSB);
0764 }
0765 
0766 void HcaluLUTTPGCoder::lookupMSB(const QIE10DataFrame& df, std::vector<std::bitset<2>>& msb) const {
0767   msb.resize(df.samples());
0768   int lutId = getLUTId(HcalDetId(df.id()));
0769   const Lut& lut = inputLUT_.at(lutId);
0770   for (int i = 0; i < df.samples(); ++i) {
0771     msb[i][0] = lut.at(df[i].adc()) & QIE10_LUT_MSB0;
0772     msb[i][1] = lut.at(df[i].adc()) & QIE10_LUT_MSB1;
0773   }
0774 }
0775 
0776 void HcaluLUTTPGCoder::lookupMSB(const QIE11DataFrame& df, std::vector<std::bitset<2>>& msb) const {
0777   int lutId = getLUTId(HcalDetId(df.id()));
0778   const Lut& lut = inputLUT_.at(lutId);
0779   for (int i = 0; i < df.samples(); ++i) {
0780     msb[i][0] = lut.at(df[i].adc()) & QIE11_LUT_MSB0;
0781     msb[i][1] = lut.at(df[i].adc()) & QIE11_LUT_MSB1;
0782   }
0783 }