Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-03-09 00:51:15

0001 //
0002 // ** class l1t::Stage2Layer2TauAlgorithmFirmwareImp1
0003 // ** authors: J. Brooke, L. Cadamuro, L. Mastrolorenzo, J.B. Sauvan, T. Strebler, O. Davignon, etc.
0004 // ** date:   2 Oct 2015
0005 // ** Description: version of tau algorithm matching the jet-eg-tau merged implementation
0006 
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "L1Trigger/L1TCalorimeter/interface/Stage2Layer2TauAlgorithmFirmware.h"
0009 
0010 #include "L1Trigger/L1TCalorimeter/interface/CaloTools.h"
0011 #include "L1Trigger/L1TCalorimeter/interface/BitonicSort.h"
0012 #include "L1Trigger/L1TCalorimeter/interface/AccumulatingSort.h"
0013 
0014 namespace l1t {
0015   bool operator>(const l1t::Tau& a, const l1t::Tau& b) { return a.pt() > b.pt(); }
0016 }  // namespace l1t
0017 
0018 l1t::Stage2Layer2TauAlgorithmFirmwareImp1::Stage2Layer2TauAlgorithmFirmwareImp1(CaloParamsHelper const* params)
0019     : params_(params) {
0020   loadCalibrationLuts();
0021 }
0022 
0023 l1t::Stage2Layer2TauAlgorithmFirmwareImp1::~Stage2Layer2TauAlgorithmFirmwareImp1() {}
0024 
0025 void l1t::Stage2Layer2TauAlgorithmFirmwareImp1::processEvent(const std::vector<l1t::CaloCluster>& clusters,
0026                                                              const std::vector<l1t::CaloTower>& towers,
0027                                                              std::vector<l1t::Tau>& taus) {
0028   // fill L1 candidates collections from clusters, merging neighbour clusters
0029   merging(clusters, towers, taus);
0030   //FIXME: TO DO
0031   // isolation   (taus);
0032   dosorting(taus);
0033 }
0034 
0035 void l1t::Stage2Layer2TauAlgorithmFirmwareImp1::merging(const std::vector<l1t::CaloCluster>& clusters,
0036                                                         const std::vector<l1t::CaloTower>& towers,
0037                                                         std::vector<l1t::Tau>& taus) {
0038   // navigator
0039   l1t::CaloStage2Nav caloNav;
0040 
0041   // this is common to all taus in this event
0042   const int nrTowers = CaloTools::calNrTowers(-1 * params_->tauPUSParam(1),
0043                                               params_->tauPUSParam(1),
0044                                               1,
0045                                               72,
0046                                               towers,
0047                                               1 + params_->pileUpTowerThreshold(),
0048                                               999,
0049                                               CaloTools::CALO);
0050 
0051   for (const auto& mainCluster : clusters) {
0052     // loop only on valid clusters
0053     // by construction of the clustering, they are local maxima in the 9x3 jet window
0054     if (mainCluster.isValid()) {
0055       if (abs(mainCluster.hwEta()) > params_->isoTauEtaMax())
0056         continue;  // limit in main seed position in firmware
0057 
0058       taus.emplace_back(mainCluster);
0059       l1t::Tau& tau = taus.back();
0060 
0061       int iEta = mainCluster.hwEta();
0062       int iPhi = mainCluster.hwPhi();
0063       int iEtaP = caloNav.offsetIEta(iEta, 1);
0064       int iEtaM = caloNav.offsetIEta(iEta, -1);
0065       int iPhiP = caloNav.offsetIPhi(iPhi, 1);
0066       int iPhiM = caloNav.offsetIPhi(iPhi, -1);
0067       int iPhiP2 = caloNav.offsetIPhi(iPhi, 2);
0068       int iPhiP3 = caloNav.offsetIPhi(iPhi, 3);
0069       int iPhiM2 = caloNav.offsetIPhi(iPhi, -2);
0070       int iPhiM3 = caloNav.offsetIPhi(iPhi, -3);
0071 
0072       // get list of neighbor seeds and determine the highest E one
0073       std::vector<l1t::CaloTower> satellites;
0074       const l1t::CaloTower& towerN2 = l1t::CaloTools::getTower(towers, iEta, iPhiM2);
0075       const l1t::CaloTower& towerN3 = l1t::CaloTools::getTower(towers, iEta, iPhiM3);
0076       const l1t::CaloTower& towerN2W = l1t::CaloTools::getTower(towers, iEtaM, iPhiM2);
0077       const l1t::CaloTower& towerN2E = l1t::CaloTools::getTower(towers, iEtaP, iPhiM2);
0078       const l1t::CaloTower& towerS2 = l1t::CaloTools::getTower(towers, iEta, iPhiP2);
0079       const l1t::CaloTower& towerS3 = l1t::CaloTools::getTower(towers, iEta, iPhiP3);
0080       const l1t::CaloTower& towerS2W = l1t::CaloTools::getTower(towers, iEtaM, iPhiP2);
0081       const l1t::CaloTower& towerS2E = l1t::CaloTools::getTower(towers, iEtaP, iPhiP2);
0082 
0083       int seedThreshold = floor(params_->egSeedThreshold() / params_->towerLsbSum());
0084 
0085       std::vector<int> sites;  // numbering of the secondary cluster sites (seed positions)
0086       // get only local max, also ask that they are above seed threshold
0087       // FIXME : in firmware N --> larger phi but apparently only for these secondaries ... sigh ...
0088       // might need to revert ; ALSO check EAST / WEST
0089       // or at least check that everything is coherent here
0090       if (is3x3Maximum(towerN2, towers, caloNav) && towerN2.hwPt() >= seedThreshold &&
0091           !mainCluster.checkClusterFlag(CaloCluster::INCLUDE_NN))
0092         sites.push_back(5);
0093       if (is3x3Maximum(towerN3, towers, caloNav) && towerN3.hwPt() >= seedThreshold)
0094         sites.push_back(7);
0095       if (is3x3Maximum(towerN2W, towers, caloNav) && towerN2W.hwPt() >= seedThreshold)
0096         sites.push_back(4);
0097       if (is3x3Maximum(towerN2E, towers, caloNav) && towerN2E.hwPt() >= seedThreshold)
0098         sites.push_back(6);
0099       if (is3x3Maximum(towerS2, towers, caloNav) && towerS2.hwPt() >= seedThreshold &&
0100           !mainCluster.checkClusterFlag(CaloCluster::INCLUDE_SS))
0101         sites.push_back(2);
0102       if (is3x3Maximum(towerS3, towers, caloNav) && towerS3.hwPt() >= seedThreshold)
0103         sites.push_back(0);
0104       if (is3x3Maximum(towerS2W, towers, caloNav) && towerS2W.hwPt() >= seedThreshold)
0105         sites.push_back(1);
0106       if (is3x3Maximum(towerS2E, towers, caloNav) && towerS2E.hwPt() >= seedThreshold)
0107         sites.push_back(3);
0108 
0109       // find neighbor with highest energy, with some preference that is defined as in the firmware
0110       // Remember: the maxima requirement is already applied
0111       // For the four towers in a T
0112       // If cluster 0 is on a maxima, use that
0113       // Else if cluster 2 is on a maxima, use that
0114       // Else if both clusters 1 & 3 are both on maxima: Use the highest energy cluster
0115       // Else if cluster 1 is on a maxima, use that
0116       // Else if cluster 3 is on a maxima, use that
0117       // Else there is no candidate
0118       // Similarly for south
0119       // Then if N>S, use N
0120       // else S
0121       auto secClusters = makeSecClusters(towers, sites, mainCluster, caloNav);
0122       l1t::CaloCluster* secMaxN = nullptr;
0123       l1t::CaloCluster* secMaxS = nullptr;
0124       l1t::CaloCluster* secondaryCluster = nullptr;
0125 
0126       std::vector<int>::iterator isNeigh0 = find(sites.begin(), sites.end(), 0);
0127       std::vector<int>::iterator isNeigh1 = find(sites.begin(), sites.end(), 1);
0128       std::vector<int>::iterator isNeigh2 = find(sites.begin(), sites.end(), 2);
0129       std::vector<int>::iterator isNeigh3 = find(sites.begin(), sites.end(), 3);
0130       std::vector<int>::iterator isNeigh4 = find(sites.begin(), sites.end(), 4);
0131       std::vector<int>::iterator isNeigh5 = find(sites.begin(), sites.end(), 5);
0132       std::vector<int>::iterator isNeigh6 = find(sites.begin(), sites.end(), 6);
0133       std::vector<int>::iterator isNeigh7 = find(sites.begin(), sites.end(), 7);
0134 
0135       // N neighbor --------------------------------------------------
0136       if (isNeigh0 != sites.end())
0137         secMaxN = secClusters.at(isNeigh0 - sites.begin()).get();
0138       else if (isNeigh2 != sites.end())
0139         secMaxN = secClusters.at(isNeigh2 - sites.begin()).get();
0140 
0141       else if (isNeigh1 != sites.end() && isNeigh3 != sites.end()) {
0142         if ((secClusters.at(isNeigh1 - sites.begin()))->hwPt() == (secClusters.at(isNeigh3 - sites.begin()))->hwPt()) {
0143           // same E --> take 1
0144           if (mainCluster.hwEta() > 0)
0145             secMaxN = secClusters.at(isNeigh1 - sites.begin()).get();
0146           else
0147             secMaxN = secClusters.at(isNeigh3 - sites.begin()).get();
0148         }
0149 
0150         else {
0151           if ((secClusters.at(isNeigh1 - sites.begin()))->hwPt() > (secClusters.at(isNeigh3 - sites.begin()))->hwPt()) {
0152             secMaxN = secClusters.at(isNeigh1 - sites.begin()).get();
0153           } else
0154             secMaxN = secClusters.at(isNeigh3 - sites.begin()).get();
0155         }
0156 
0157       }
0158 
0159       else if (isNeigh1 != sites.end())
0160         secMaxN = secClusters.at(isNeigh1 - sites.begin()).get();
0161       else if (isNeigh3 != sites.end())
0162         secMaxN = secClusters.at(isNeigh3 - sites.begin()).get();
0163 
0164       // S neighbor --------------------------------------------------
0165       if (isNeigh7 != sites.end())
0166         secMaxS = secClusters.at(isNeigh7 - sites.begin()).get();
0167       else if (isNeigh5 != sites.end())
0168         secMaxS = secClusters.at(isNeigh5 - sites.begin()).get();
0169 
0170       else if (isNeigh4 != sites.end() && isNeigh6 != sites.end()) {
0171         if ((secClusters.at(isNeigh4 - sites.begin()))->hwPt() == (secClusters.at(isNeigh6 - sites.begin()))->hwPt()) {
0172           // same E --> take 1
0173           if (mainCluster.hwEta() > 0)
0174             secMaxS = secClusters.at(isNeigh4 - sites.begin()).get();
0175           else
0176             secMaxS = secClusters.at(isNeigh6 - sites.begin()).get();
0177         }
0178 
0179         else {
0180           if ((secClusters.at(isNeigh4 - sites.begin()))->hwPt() > (secClusters.at(isNeigh6 - sites.begin()))->hwPt())
0181             secMaxS = secClusters.at(isNeigh4 - sites.begin()).get();
0182           else
0183             secMaxS = secClusters.at(isNeigh6 - sites.begin()).get();
0184         }
0185 
0186       }
0187 
0188       else if (isNeigh4 != sites.end())
0189         secMaxS = secClusters.at(isNeigh4 - sites.begin()).get();
0190       else if (isNeigh6 != sites.end())
0191         secMaxS = secClusters.at(isNeigh6 - sites.begin()).get();
0192 
0193       // N vs S neighbor --------------------------------------------------
0194       if (secMaxN != nullptr && secMaxS != nullptr) {
0195         if (secMaxN->hwPt() > secMaxS->hwPt())
0196           secondaryCluster = secMaxN;
0197         else
0198           secondaryCluster = secMaxS;
0199       } else {
0200         if (secMaxN != nullptr)
0201           secondaryCluster = secMaxN;
0202         else if (secMaxS != nullptr)
0203           secondaryCluster = secMaxS;
0204       }
0205 
0206       // trim primary cluster to remove overlap of TT
0207       int neigEta[8] = {0, -1, 0, 1, -1, 0, 1, 0};
0208       int neigPhi[8] = {3, 2, 2, 2, -2, -2, -2, -3};
0209 
0210       vector<pair<int, int>> TTPos(10);  // numbering of TT in cluster; <iEta, iPhi>
0211       TTPos.at(0) = make_pair(-1, 1);    // SW
0212       TTPos.at(1) = make_pair(0, 1);     // S
0213       TTPos.at(2) = make_pair(1, 1);     // SE
0214       TTPos.at(3) = make_pair(1, 0);     // E
0215       TTPos.at(4) = make_pair(1, -1);    // NE
0216       TTPos.at(5) = make_pair(0, -1);    // N
0217       TTPos.at(6) = make_pair(-1, -1);   // NW
0218       TTPos.at(7) = make_pair(-1, 0);    // W
0219       TTPos.at(8) = make_pair(0, 2);     // SS
0220       TTPos.at(9) = make_pair(0, -2);    // NN
0221 
0222       vector<CaloCluster::ClusterFlag> TTPosRemap(10);  // using geographical notation
0223       TTPosRemap.at(0) = CaloCluster::INCLUDE_SW;
0224       TTPosRemap.at(1) = CaloCluster::INCLUDE_S;
0225       TTPosRemap.at(2) = CaloCluster::INCLUDE_SE;
0226       TTPosRemap.at(3) = CaloCluster::INCLUDE_E;
0227       TTPosRemap.at(4) = CaloCluster::INCLUDE_NE;
0228       TTPosRemap.at(5) = CaloCluster::INCLUDE_N;
0229       TTPosRemap.at(6) = CaloCluster::INCLUDE_NW;
0230       TTPosRemap.at(7) = CaloCluster::INCLUDE_W;
0231       TTPosRemap.at(8) = CaloCluster::INCLUDE_SS;
0232       TTPosRemap.at(9) = CaloCluster::INCLUDE_NN;
0233 
0234       // loop over TT of secondary cluster, if there is overlap remove this towers from main cluster
0235       l1t::CaloCluster mainClusterTrim = mainCluster;
0236       int secondaryClusterHwPt = 0;
0237 
0238       if (!secClusters.empty()) {
0239         secondaryClusterHwPt = secondaryCluster->hwPt();
0240 
0241         int iSecIdxPosition = find_if(secClusters.begin(),
0242                                       secClusters.end(),
0243                                       [&](auto const& element) { return element.get() == secondaryCluster; }) -
0244                               secClusters.begin();
0245         int secondaryClusterSite = sites.at(iSecIdxPosition);
0246 
0247         for (unsigned int iTT = 0; iTT < TTPos.size(); iTT++) {
0248           //get this TT in the "frame" of the main
0249           int thisTTinMainEta = neigEta[secondaryClusterSite] + TTPos.at(iTT).first;
0250           int thisTTinMainPhi = neigPhi[secondaryClusterSite] + TTPos.at(iTT).second;
0251           pair<int, int> thisTT = make_pair(thisTTinMainEta, thisTTinMainPhi);
0252           // check if main cluster has this tower included; if true, switch it off
0253           auto thisTTItr = find(TTPos.begin(), TTPos.end(), thisTT);
0254           if (thisTTItr != TTPos.end()) {
0255             int idx = thisTTItr - TTPos.begin();
0256             if (secondaryCluster->checkClusterFlag(TTPosRemap.at(iTT)))
0257               mainClusterTrim.setClusterFlag(TTPosRemap.at(idx), false);
0258           }
0259         }
0260       }
0261 
0262       // re-compute main cluster energy + position
0263       const l1t::CaloTower& seed = l1t::CaloTools::getTower(towers, iEta, iPhi);
0264       const l1t::CaloTower& towerNW = l1t::CaloTools::getTower(towers, iEtaM, iPhiM);
0265       const l1t::CaloTower& towerN = l1t::CaloTools::getTower(towers, iEta, iPhiM);
0266       const l1t::CaloTower& towerNE = l1t::CaloTools::getTower(towers, iEtaP, iPhiM);
0267       const l1t::CaloTower& towerE = l1t::CaloTools::getTower(towers, iEtaP, iPhi);
0268       const l1t::CaloTower& towerSE = l1t::CaloTools::getTower(towers, iEtaP, iPhiP);
0269       const l1t::CaloTower& towerS = l1t::CaloTools::getTower(towers, iEta, iPhiP);
0270       const l1t::CaloTower& towerSW = l1t::CaloTools::getTower(towers, iEtaM, iPhiP);
0271       const l1t::CaloTower& towerW = l1t::CaloTools::getTower(towers, iEtaM, iPhi);
0272       const l1t::CaloTower& towerNN = l1t::CaloTools::getTower(towers, iEta, iPhiM2);
0273       const l1t::CaloTower& towerSS = l1t::CaloTools::getTower(towers, iEta, iPhiP2);
0274 
0275       int seedEt = seed.hwPt();
0276       int towerEtNW = towerNW.hwPt();
0277       int towerEtN = towerN.hwPt();
0278       int towerEtNE = towerNE.hwPt();
0279       int towerEtE = towerE.hwPt();
0280       int towerEtSE = towerSE.hwPt();
0281       int towerEtS = towerS.hwPt();
0282       int towerEtSW = towerSW.hwPt();
0283       int towerEtW = towerW.hwPt();
0284       int towerEtNN = towerNN.hwPt();
0285       int towerEtSS = towerSS.hwPt();
0286 
0287       // Recompute hw energy from towers
0288       tau.setHwPt(seedEt);
0289       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NW))
0290         tau.setHwPt(tau.hwPt() + towerEtNW);
0291       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_N))
0292         tau.setHwPt(tau.hwPt() + towerEtN);
0293       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NE))
0294         tau.setHwPt(tau.hwPt() + towerEtNE);
0295       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_E))
0296         tau.setHwPt(tau.hwPt() + towerEtE);
0297       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SE))
0298         tau.setHwPt(tau.hwPt() + towerEtSE);
0299       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_S))
0300         tau.setHwPt(tau.hwPt() + towerEtS);
0301       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SW))
0302         tau.setHwPt(tau.hwPt() + towerEtSW);
0303       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_W))
0304         tau.setHwPt(tau.hwPt() + towerEtW);
0305       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NN))
0306         tau.setHwPt(tau.hwPt() + towerEtNN);
0307       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SS))
0308         tau.setHwPt(tau.hwPt() + towerEtSS);
0309 
0310       //Merging
0311       tau.setRawEt(tau.hwPt() + secondaryClusterHwPt);
0312       tau.setIsMerged(secondaryClusterHwPt > 0);
0313 
0314       // ==================================================================
0315       // Energy calibration
0316       // ==================================================================
0317       // Corrections function of ieta, ET, and cluster shape
0318       int calibPt = calibratedPt(mainCluster, towers, tau.rawEt(), tau.isMerged());
0319       tau.setHwPt(calibPt);
0320 
0321       // isolation
0322       int isoLeftExtension = params_->tauIsoAreaNrTowersEta();
0323       int isoRightExtension = params_->tauIsoAreaNrTowersEta();
0324       if (mainCluster.checkClusterFlag(CaloCluster::TRIM_LEFT))
0325         isoRightExtension++;
0326       else
0327         isoLeftExtension++;
0328 
0329       int isolBit = 0;
0330       int tauHwFootprint = tau.rawEt();
0331       unsigned int LUTaddress = isoLutIndex(tauHwFootprint, mainCluster.hwEta(), nrTowers);
0332       int hwEtSum = CaloTools::calHwEtSum(iEta,
0333                                           iPhi,
0334                                           towers,
0335                                           -isoLeftExtension,
0336                                           isoRightExtension,
0337                                           -1 * params_->tauIsoAreaNrTowersPhi(),
0338                                           params_->tauIsoAreaNrTowersPhi(),
0339                                           params_->tauPUSParam(2),
0340                                           CaloTools::CALO);
0341       int hwIsoEnergy = hwEtSum - tauHwFootprint;
0342       if (hwIsoEnergy < 0)
0343         hwIsoEnergy = 0;  // just in case the cluster is outside the window? should be very rare
0344 
0345       isolBit = (((hwIsoEnergy < (params_->tauIsolationLUT()->data(LUTaddress))) ||
0346                   (params_->tauIsolationLUT()->data(LUTaddress) > 255))
0347                      ? 1
0348                      : 0);
0349       //int isolBit2 = (((hwIsoEnergy < (params_->tauIsolationLUT2()->data(LUTaddress))) || (params_->tauIsolationLUT2()->data(LUTaddress)>255)) ? 1 : 0);
0350       //isolBit += (isolBit2 << 1);
0351       tau.setHwIso(isolBit);
0352 
0353       // development vars
0354       tau.setTowerIPhi(iPhi);
0355       tau.setTowerIEta(iEta);
0356       int qual = seed.hwQual();
0357       bool denomZeroFlag = ((qual & 0x1) > 0);
0358       bool eOverHFlag = ((qual & 0x2) > 0);
0359       int hasEM = (eOverHFlag || !denomZeroFlag);
0360       tau.setHasEM(hasEM > 0);
0361       tau.setIsoEt((short int)hwIsoEnergy);
0362       tau.setNTT((short int)nrTowers);
0363 
0364       // Physical eta/phi. Computed from ieta/iphi of the seed tower and the fine-grain position within the seed
0365       double eta = 0.;
0366       double phi = 0.;
0367       double seedEta = CaloTools::towerEta(iEta);
0368       double seedEtaSize = CaloTools::towerEtaSize(iEta);
0369       double seedPhi = CaloTools::towerPhi(iEta, iPhi);
0370       double seedPhiSize = CaloTools::towerPhiSize(iEta);
0371       if (mainCluster.fgEta() == 0)
0372         eta = seedEta;  // center
0373       else if (mainCluster.fgEta() == 2)
0374         eta = seedEta + seedEtaSize * 0.251;  // center + 1/4
0375       else if (mainCluster.fgEta() == 1)
0376         eta = seedEta - seedEtaSize * 0.251;  // center - 1/4
0377 
0378       //fgPhi is recomputed after trimming
0379       int fgPhi = 0;
0380 
0381       int EtUp = 0;
0382       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NE))
0383         EtUp += towerEtNE;
0384       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_N))
0385         EtUp += towerEtN;
0386       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NW))
0387         EtUp += towerEtNW;
0388       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_NN))
0389         EtUp += towerEtNN;
0390       int EtDown = 0;
0391       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SE))
0392         EtDown += towerEtSE;
0393       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_S))
0394         EtDown += towerEtS;
0395       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SW))
0396         EtDown += towerEtSW;
0397       if (mainClusterTrim.checkClusterFlag(CaloCluster::INCLUDE_SS))
0398         EtDown += towerEtSS;
0399       //
0400       if (EtDown > EtUp)
0401         fgPhi = 2;
0402       else if (EtUp > EtDown)
0403         fgPhi = 1;
0404 
0405       if (fgPhi == 0)
0406         phi = seedPhi;  // center
0407       else if (fgPhi == 2)
0408         phi = seedPhi + seedPhiSize * 0.251;  // center + 1/4
0409       else if (fgPhi == 1)
0410         phi = seedPhi - seedPhiSize * 0.251;  // center - 1/4
0411 
0412       // Set 4-vector
0413       math::PtEtaPhiMLorentzVector calibP4((double)calibPt * params_->tauLsb(), eta, phi, 0.);
0414       tau.setP4(calibP4);
0415 
0416       // delete all sec clusters that were allocated with new
0417       secClusters.clear();
0418     }
0419 
0420   }  // end loop on clusters
0421 }
0422 
0423 // -----------------------------------------------------------------------------------
0424 void l1t::Stage2Layer2TauAlgorithmFirmwareImp1::dosorting(std::vector<l1t::Tau>& taus) {
0425   // prepare content to be sorted -- each phi ring contains 18 elements, with Et = 0 if no candidate exists
0426   math::PtEtaPhiMLorentzVector emptyP4;
0427   l1t::Tau tempTau(emptyP4, 0, 0, 0, 0);
0428   std::vector<std::vector<l1t::Tau>> tauEtaPos(params_->isoTauEtaMax(), std::vector<l1t::Tau>(18, tempTau));
0429   std::vector<std::vector<l1t::Tau>> tauEtaNeg(params_->isoTauEtaMax(), std::vector<l1t::Tau>(18, tempTau));
0430 
0431   for (unsigned int iTau = 0; iTau < taus.size(); iTau++) {
0432     if (taus.at(iTau).hwEta() > 0)
0433       tauEtaPos.at(taus.at(iTau).hwEta() - 1).at((72 - taus.at(iTau).hwPhi()) / 4) = taus.at(iTau);
0434     else
0435       tauEtaNeg.at(-(taus.at(iTau).hwEta() + 1)).at((72 - taus.at(iTau).hwPhi()) / 4) = taus.at(iTau);
0436   }
0437 
0438   AccumulatingSort<l1t::Tau> etaPosSorter(6);
0439   AccumulatingSort<l1t::Tau> etaNegSorter(6);
0440   std::vector<l1t::Tau> accumEtaPos;
0441   std::vector<l1t::Tau> accumEtaNeg;
0442 
0443   for (int ieta = 0; ieta < params_->isoTauEtaMax(); ++ieta) {
0444     // eta +
0445     std::vector<l1t::Tau>::iterator start_, end_;
0446     start_ = tauEtaPos.at(ieta).begin();
0447     end_ = tauEtaPos.at(ieta).end();
0448     BitonicSort<l1t::Tau>(down, start_, end_);
0449     etaPosSorter.Merge(tauEtaPos.at(ieta), accumEtaPos);
0450 
0451     // eta -
0452     start_ = tauEtaNeg.at(ieta).begin();
0453     end_ = tauEtaNeg.at(ieta).end();
0454     BitonicSort<l1t::Tau>(down, start_, end_);
0455     etaNegSorter.Merge(tauEtaNeg.at(ieta), accumEtaNeg);
0456   }
0457 
0458   // put all 12 candidates in the original tau vector, removing zero energy ones
0459   taus.clear();
0460   for (const l1t::Tau& acctau : accumEtaPos) {
0461     if (acctau.hwPt() > 0)
0462       taus.push_back(acctau);
0463   }
0464   for (const l1t::Tau& acctau : accumEtaNeg) {
0465     if (acctau.hwPt() > 0)
0466       taus.push_back(acctau);
0467   }
0468 }
0469 
0470 // -----------------------------------------------------------------------------------
0471 
0472 void l1t::Stage2Layer2TauAlgorithmFirmwareImp1::loadCalibrationLuts() {
0473   float minScale = 0.;
0474   float maxScale = 2.;
0475   float minScaleEta = 0.5;
0476   float maxScaleEta = 1.5;
0477   offsetBarrelEH_ = 0.5;
0478   offsetBarrelH_ = 1.5;
0479   offsetEndcapsEH_ = 0.;
0480   offsetEndcapsH_ = 1.5;
0481 
0482   // In the combined calibration LUT, upper 3-bits are used as LUT index:
0483   // (0=BarrelA, 1=BarrelB, 2=BarrelC, 3=EndCapA, 4=EndCapA, 5=EndCapA, 6=Eta)
0484   enum { LUT_UPPER = 3 };
0485   enum { LUT_OFFSET = 0x80 };
0486   l1t::LUT const* lut = params_->tauCalibrationLUT();
0487   unsigned int size = (1 << lut->nrBitsData());
0488   unsigned int nBins = (1 << (lut->nrBitsAddress() - LUT_UPPER));
0489 
0490   std::vector<float> emptyCoeff;
0491   emptyCoeff.resize(nBins, 0.);
0492   float binSize = (maxScale - minScale) / (float)size;
0493   for (unsigned iLut = 0; iLut < 6; ++iLut) {
0494     coefficients_.push_back(emptyCoeff);
0495     for (unsigned addr = 0; addr < nBins; addr++) {
0496       float y = (float)lut->data(iLut * LUT_OFFSET + addr);
0497       coefficients_[iLut][addr] = minScale + binSize * y;
0498     }
0499   }
0500 
0501   size = (1 << lut->nrBitsData());
0502   nBins = (1 << 6);  // can't auto-extract this now due to combined LUT.
0503 
0504   emptyCoeff.resize(nBins, 0.);
0505   binSize = (maxScaleEta - minScaleEta) / (float)size;
0506   coefficients_.push_back(emptyCoeff);
0507   for (unsigned addr = 0; addr < nBins; addr++) {
0508     float y = (float)lut->data(6 * LUT_OFFSET + addr);
0509     coefficients_.back()[addr] = minScaleEta + binSize * y;
0510   }
0511 }
0512 
0513 bool l1t::Stage2Layer2TauAlgorithmFirmwareImp1::compareTowers(l1t::CaloTower TT1, l1t::CaloTower TT2) {
0514   // 1. compare hwPt (for the moment no switch with E and H only, always use E+H)
0515   if (TT1.hwPt() < TT2.hwPt())
0516     return true;
0517   if (TT1.hwPt() > TT2.hwPt())
0518     return false;
0519 
0520   // 2. if equal pT, most central -- eta is in the range -32, 32 with ieta = 0 skipped
0521   if (abs(TT1.hwEta()) > abs(TT2.hwEta()))
0522     return true;
0523   if (abs(TT1.hwEta()) < abs(TT2.hwEta()))
0524     return false;
0525 
0526   // 3. if equal eta, compare phi (arbitrary)
0527   return (TT1.hwPhi() < TT2.hwPhi());  // towers towards S are favored (remember: N --> smaller phi, S --> larger phi)
0528 }
0529 
0530 bool l1t::Stage2Layer2TauAlgorithmFirmwareImp1::is3x3Maximum(const l1t::CaloTower& tower,
0531                                                              const std::vector<CaloTower>& towers,
0532                                                              l1t::CaloStage2Nav& caloNav) {
0533   int iEta = tower.hwEta();
0534   int iPhi = tower.hwPhi();
0535 
0536   // 1 : >
0537   // 2 : >=
0538   int mask[3][3] = {
0539       {1, 2, 2},
0540       {1, 0, 2},
0541       {1, 1, 2},
0542   };
0543 
0544   bool vetoTT = false;  // false if it is a local max i.e. no veto
0545   for (int deta = -1; deta < 2; deta++) {
0546     for (int dphi = -1; dphi < 2; dphi++) {
0547       int iEtaNeigh = caloNav.offsetIEta(iEta, deta);
0548       int iPhiNeigh = caloNav.offsetIPhi(iPhi, dphi);
0549       const l1t::CaloTower& towerNeigh = l1t::CaloTools::getTower(towers, iEtaNeigh, iPhiNeigh);
0550       if (mask[2 - (dphi + 1)][deta + 1] == 0)
0551         continue;
0552       if (mask[2 - (dphi + 1)][deta + 1] == 1)
0553         vetoTT = (tower.hwPt() < towerNeigh.hwPt());
0554       if (mask[2 - (dphi + 1)][deta + 1] == 2)
0555         vetoTT = (tower.hwPt() <= towerNeigh.hwPt());
0556 
0557       if (vetoTT)
0558         break;
0559     }
0560     if (vetoTT)
0561       break;
0562   }
0563 
0564   return (!vetoTT);  // negate because I ask if is a local maxima
0565 }
0566 
0567 std::vector<std::unique_ptr<l1t::CaloCluster>> l1t::Stage2Layer2TauAlgorithmFirmwareImp1::makeSecClusters(
0568     const std::vector<l1t::CaloTower>& towers,
0569     std::vector<int>& sites,
0570     const l1t::CaloCluster& mainCluster,
0571     l1t::CaloStage2Nav& caloNav) {
0572   constexpr int neigEta[8] = {0, -1, 0, 1, -1, 0, 1, 0};
0573   constexpr int neigPhi[8] = {3, 2, 2, 2, -2, -2, -2, -3};
0574   int clusterThreshold = floor(params_->egNeighbourThreshold() / params_->towerLsbSum());
0575 
0576   int iEtamain = mainCluster.hwEta();
0577   int iPhimain = mainCluster.hwPhi();
0578 
0579   std::vector<unique_ptr<CaloCluster>> secClusters;
0580   for (unsigned int isite = 0; isite < sites.size(); isite++) {
0581     // build full cluster at this site
0582     const int siteNumber = sites.at(isite);
0583     int iSecEta = caloNav.offsetIEta(iEtamain, neigEta[siteNumber]);
0584     int iSecPhi = caloNav.offsetIPhi(iPhimain, neigPhi[siteNumber]);
0585 
0586     const l1t::CaloTower& towerSec = l1t::CaloTools::getTower(towers, iSecEta, iSecPhi);
0587 
0588     math::XYZTLorentzVector emptyP4;
0589     auto secondaryCluster =
0590         std::make_unique<l1t::CaloCluster>(emptyP4, towerSec.hwPt(), towerSec.hwEta(), towerSec.hwPhi());
0591 
0592     secondaryCluster->setHwPtEm(towerSec.hwEtEm());
0593     secondaryCluster->setHwPtHad(towerSec.hwEtHad());
0594     secondaryCluster->setHwSeedPt(towerSec.hwPt());
0595     secondaryCluster->setHwPt(towerSec.hwPt());
0596 
0597     int iSecEtaP = caloNav.offsetIEta(iSecEta, 1);
0598     int iSecEtaM = caloNav.offsetIEta(iSecEta, -1);
0599     int iSecPhiP = caloNav.offsetIPhi(iSecPhi, 1);
0600     int iSecPhiP2 = caloNav.offsetIPhi(iSecPhi, 2);
0601     int iSecPhiM = caloNav.offsetIPhi(iSecPhi, -1);
0602     int iSecPhiM2 = caloNav.offsetIPhi(iSecPhi, -2);
0603     const l1t::CaloTower& towerNW = l1t::CaloTools::getTower(towers, iSecEtaM, iSecPhiM);
0604     const l1t::CaloTower& towerN = l1t::CaloTools::getTower(towers, iSecEta, iSecPhiM);
0605     const l1t::CaloTower& towerNE = l1t::CaloTools::getTower(towers, iSecEtaP, iSecPhiM);
0606     const l1t::CaloTower& towerE = l1t::CaloTools::getTower(towers, iSecEtaP, iSecPhi);
0607     const l1t::CaloTower& towerSE = l1t::CaloTools::getTower(towers, iSecEtaP, iSecPhiP);
0608     const l1t::CaloTower& towerS = l1t::CaloTools::getTower(towers, iSecEta, iSecPhiP);
0609     const l1t::CaloTower& towerSW = l1t::CaloTools::getTower(towers, iSecEtaM, iSecPhiP);
0610     const l1t::CaloTower& towerW = l1t::CaloTools::getTower(towers, iSecEtaM, iSecPhi);
0611     const l1t::CaloTower& towerNN = l1t::CaloTools::getTower(towers, iSecEta, iSecPhiM2);
0612     const l1t::CaloTower& towerSS = l1t::CaloTools::getTower(towers, iSecEta, iSecPhiP2);
0613 
0614     // just use E+H for clustering
0615     int towerEtNW = towerNW.hwPt();
0616     int towerEtN = towerN.hwPt();
0617     int towerEtNE = towerNE.hwPt();
0618     int towerEtE = towerE.hwPt();
0619     int towerEtSE = towerSE.hwPt();
0620     int towerEtS = towerS.hwPt();
0621     int towerEtSW = towerSW.hwPt();
0622     int towerEtW = towerW.hwPt();
0623     int towerEtNN = towerNN.hwPt();
0624     int towerEtSS = towerSS.hwPt();
0625 
0626     int towerEtEmNW = towerNW.hwEtEm();
0627     int towerEtEmN = towerN.hwEtEm();
0628     int towerEtEmNE = towerNE.hwEtEm();
0629     int towerEtEmE = towerE.hwEtEm();
0630     int towerEtEmSE = towerSE.hwEtEm();
0631     int towerEtEmS = towerS.hwEtEm();
0632     int towerEtEmSW = towerSW.hwEtEm();
0633     int towerEtEmW = towerW.hwEtEm();
0634     int towerEtEmNN = towerNN.hwEtEm();
0635     int towerEtEmSS = towerSS.hwEtEm();
0636     //
0637     int towerEtHadNW = towerNW.hwEtHad();
0638     int towerEtHadN = towerN.hwEtHad();
0639     int towerEtHadNE = towerNE.hwEtHad();
0640     int towerEtHadE = towerE.hwEtHad();
0641     int towerEtHadSE = towerSE.hwEtHad();
0642     int towerEtHadS = towerS.hwEtHad();
0643     int towerEtHadSW = towerSW.hwEtHad();
0644     int towerEtHadW = towerW.hwEtHad();
0645     int towerEtHadNN = towerNN.hwEtHad();
0646     int towerEtHadSS = towerSS.hwEtHad();
0647 
0648     if (towerEtNW < clusterThreshold)
0649       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_NW, false);
0650     if (towerEtN < clusterThreshold) {
0651       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_N, false);
0652       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_NN, false);
0653     }
0654     if (towerEtNE < clusterThreshold)
0655       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_NE, false);
0656     if (towerEtE < clusterThreshold)
0657       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_E, false);
0658     if (towerEtSE < clusterThreshold)
0659       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_SE, false);
0660     if (towerEtS < clusterThreshold) {
0661       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_S, false);
0662       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_SS, false);
0663     }
0664     if (towerEtSW < clusterThreshold)
0665       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_SW, false);
0666     if (towerEtW < clusterThreshold)
0667       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_W, false);
0668     if (towerEtNN < clusterThreshold)
0669       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_NN, false);
0670     if (towerEtSS < clusterThreshold)
0671       secondaryCluster->setClusterFlag(CaloCluster::INCLUDE_SS, false);
0672 
0673     // compute cluster energy according to cluster flags
0674     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NW))
0675       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtNW);
0676     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_N))
0677       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtN);
0678     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NE))
0679       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtNE);
0680     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_E))
0681       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtE);
0682     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SE))
0683       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtSE);
0684     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_S))
0685       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtS);
0686     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SW))
0687       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtSW);
0688     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_W))
0689       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtW);
0690     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NN))
0691       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtNN);
0692     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SS))
0693       secondaryCluster->setHwPt(secondaryCluster->hwPt() + towerEtSS);
0694     //
0695     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NW))
0696       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmNW);
0697     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_N))
0698       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmN);
0699     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NE))
0700       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmNE);
0701     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_E))
0702       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmE);
0703     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SE))
0704       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmSE);
0705     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_S))
0706       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmS);
0707     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SW))
0708       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmSW);
0709     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_W))
0710       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmW);
0711     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NN))
0712       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmNN);
0713     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SS))
0714       secondaryCluster->setHwPtEm(secondaryCluster->hwPtEm() + towerEtEmSS);
0715     //
0716     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NW))
0717       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadNW);
0718     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_N))
0719       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadN);
0720     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NE))
0721       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadNE);
0722     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_E))
0723       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadE);
0724     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SE))
0725       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadSE);
0726     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_S))
0727       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadS);
0728     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SW))
0729       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadSW);
0730     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_W))
0731       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadW);
0732     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_NN))
0733       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadNN);
0734     if (secondaryCluster->checkClusterFlag(CaloCluster::INCLUDE_SS))
0735       secondaryCluster->setHwPtHad(secondaryCluster->hwPtHad() + towerEtHadSS);
0736 
0737     // save this cluster in the vector
0738     secClusters.emplace_back(std::move(secondaryCluster));
0739   }
0740   return secClusters;
0741 }
0742 
0743 // isMerged=0,1 ; hasEM=0,1
0744 unsigned int l1t::Stage2Layer2TauAlgorithmFirmwareImp1::calibLutIndex(int ieta, int Et, int hasEM, int isMerged) {
0745   int absieta = abs(ieta);
0746   if (absieta > 28)
0747     absieta = 28;
0748   if (Et > 255)
0749     Et = 255;
0750 
0751   unsigned int compressedEta = params_->tauCompressLUT()->data(absieta);
0752   unsigned int compressedEt = params_->tauCompressLUT()->data((0x1 << 5) + Et);
0753 
0754   unsigned int address = (compressedEta << 7) + (compressedEt << 2) + (hasEM << 1) +
0755                          isMerged;  //2 bits eta, 5 bits et, 1 bit hasEM, 1 bis isMerged
0756   // unsigned int address =  (compressedEt<<4)+(compressedEta<<2)+(hasEM<<1)+isMerged;
0757   return address;
0758 }
0759 
0760 int l1t::Stage2Layer2TauAlgorithmFirmwareImp1::calibratedPt(const l1t::CaloCluster& clus,
0761                                                             const std::vector<l1t::CaloTower>& towers,
0762                                                             int hwPt,
0763                                                             bool isMerged) {
0764   // get seed tower
0765   const l1t::CaloTower& seedTT = l1t::CaloTools::getTower(towers, CaloTools::caloEta(clus.hwEta()), clus.hwPhi());
0766   int qual = seedTT.hwQual();
0767   bool denomZeroFlag = ((qual & 0x1) > 0);
0768   bool eOverHFlag = ((qual & 0x2) > 0);
0769   int hasEM = (eOverHFlag || !denomZeroFlag);
0770   int isMergedI = (isMerged ? 1 : 0);
0771 
0772   unsigned int idx = calibLutIndex(clus.hwEta(), hwPt, hasEM, isMergedI);
0773   unsigned int corr = params_->tauCalibrationLUT()->data(idx);
0774 
0775   // now apply calibration factor: corrPt = rawPt * (corr[LUT] + 0.5)
0776   // where corr[LUT] is an integer mapped to the range [0, 2]
0777   int rawPt = hwPt;
0778   if (rawPt > 8191)
0779     rawPt = 8191;  // 13 bits for uncalibrated E
0780 
0781   int corrXrawPt = corr * rawPt;    // 17 bits
0782   int calibPt = (corrXrawPt >> 9);  // (10 bits) = (7 bits) + (9 bits)
0783   if (calibPt > 4095)
0784     calibPt = 4095;  // 12 bit in output
0785 
0786   return calibPt;
0787 }
0788 
0789 unsigned int l1t::Stage2Layer2TauAlgorithmFirmwareImp1::isoLutIndex(int Et, int hweta, unsigned int nrTowers) {
0790   // normalize to limits
0791   int aeta = abs(hweta);
0792 
0793   // input bits (NB: must be THE SAME in the input LUT for the compression)
0794   // int etaBits = 6  --> 64
0795   // int etBits  = 13 --> 8192
0796   // int nTTBits = 10 --> 1024
0797   if (Et >= 255)
0798     Et = 255;  // 8 bit for the input of compression LUT
0799   if (aeta >= 31)
0800     aeta = 31;
0801   if (nrTowers >= 1023)
0802     nrTowers = 1023;
0803 
0804   // get compressed value
0805   // NB: these also must MATCH the values in the LUT --> fix when new compression scheme is used
0806   // ultimately, the same compresison LUT as calibration will be used
0807   // etaCmprBits = 2;
0808   // EtCmprBits  = 4;//changed from 3, transparent to user
0809   // nTTCmprBits = 3;
0810 
0811   int etaCmpr = params_->tauCompressLUT()->data(aeta);
0812   int etCmpr = params_->tauCompressLUT()->data((0x1 << 5) + Et);  //offset: 5 bits from ieta
0813   int nTTCmpr = params_->tauCompressLUT()->data((0x1 << 5) + (0x1 << 8) +
0814                                                 nrTowers);  //offset non-compressed: 5 bits from ieta, 8 bits from iEt
0815 
0816   // get the address -- NOTE: this also depends on the compression scheme!
0817 
0818   unsigned int address =
0819       ((etaCmpr << 10) | (etCmpr << 5) | nTTCmpr);  //ordering compressed: 2 bits iEta, 5 bits iEt, 5 bits nTT
0820 
0821   return address;
0822 }