File indexing completed on 2024-04-06 12:20:25
0001
0002
0003
0004
0005
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 }
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
0029 merging(clusters, towers, taus);
0030
0031
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
0039 l1t::CaloStage2Nav caloNav;
0040
0041
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
0053
0054 if (mainCluster.isValid()) {
0055 if (abs(mainCluster.hwEta()) > params_->isoTauEtaMax())
0056 continue;
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
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;
0086
0087
0088
0089
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
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
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
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
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
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
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
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
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);
0211 TTPos.at(0) = make_pair(-1, 1);
0212 TTPos.at(1) = make_pair(0, 1);
0213 TTPos.at(2) = make_pair(1, 1);
0214 TTPos.at(3) = make_pair(1, 0);
0215 TTPos.at(4) = make_pair(1, -1);
0216 TTPos.at(5) = make_pair(0, -1);
0217 TTPos.at(6) = make_pair(-1, -1);
0218 TTPos.at(7) = make_pair(-1, 0);
0219 TTPos.at(8) = make_pair(0, 2);
0220 TTPos.at(9) = make_pair(0, -2);
0221
0222 vector<CaloCluster::ClusterFlag> TTPosRemap(10);
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
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
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
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
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
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
0311 tau.setRawEt(tau.hwPt() + secondaryClusterHwPt);
0312 tau.setIsMerged(secondaryClusterHwPt > 0);
0313
0314
0315
0316
0317
0318 int calibPt = calibratedPt(mainCluster, towers, tau.rawEt(), tau.isMerged());
0319 tau.setHwPt(calibPt);
0320
0321
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;
0344
0345 isolBit = (((hwIsoEnergy < (params_->tauIsolationLUT()->data(LUTaddress))) ||
0346 (params_->tauIsolationLUT()->data(LUTaddress) > 255))
0347 ? 1
0348 : 0);
0349
0350
0351 tau.setHwIso(isolBit);
0352
0353
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
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;
0373 else if (mainCluster.fgEta() == 2)
0374 eta = seedEta + seedEtaSize * 0.251;
0375 else if (mainCluster.fgEta() == 1)
0376 eta = seedEta - seedEtaSize * 0.251;
0377
0378
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;
0407 else if (fgPhi == 2)
0408 phi = seedPhi + seedPhiSize * 0.251;
0409 else if (fgPhi == 1)
0410 phi = seedPhi - seedPhiSize * 0.251;
0411
0412
0413 math::PtEtaPhiMLorentzVector calibP4((double)calibPt * params_->tauLsb(), eta, phi, 0.);
0414 tau.setP4(calibP4);
0415
0416
0417 secClusters.clear();
0418 }
0419
0420 }
0421 }
0422
0423
0424 void l1t::Stage2Layer2TauAlgorithmFirmwareImp1::dosorting(std::vector<l1t::Tau>& taus) {
0425
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
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
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
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
0483
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);
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
0515 if (TT1.hwPt() < TT2.hwPt())
0516 return true;
0517 if (TT1.hwPt() > TT2.hwPt())
0518 return false;
0519
0520
0521 if (abs(TT1.hwEta()) > abs(TT2.hwEta()))
0522 return true;
0523 if (abs(TT1.hwEta()) < abs(TT2.hwEta()))
0524 return false;
0525
0526
0527 return (TT1.hwPhi() < TT2.hwPhi());
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
0537
0538 int mask[3][3] = {
0539 {1, 2, 2},
0540 {1, 0, 2},
0541 {1, 1, 2},
0542 };
0543
0544 bool vetoTT = false;
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);
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
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
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
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
0738 secClusters.emplace_back(std::move(secondaryCluster));
0739 }
0740 return secClusters;
0741 }
0742
0743
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;
0756
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
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
0776
0777 int rawPt = hwPt;
0778 if (rawPt > 8191)
0779 rawPt = 8191;
0780
0781 int corrXrawPt = corr * rawPt;
0782 int calibPt = (corrXrawPt >> 9);
0783 if (calibPt > 4095)
0784 calibPt = 4095;
0785
0786 return calibPt;
0787 }
0788
0789 unsigned int l1t::Stage2Layer2TauAlgorithmFirmwareImp1::isoLutIndex(int Et, int hweta, unsigned int nrTowers) {
0790
0791 int aeta = abs(hweta);
0792
0793
0794
0795
0796
0797 if (Et >= 255)
0798 Et = 255;
0799 if (aeta >= 31)
0800 aeta = 31;
0801 if (nrTowers >= 1023)
0802 nrTowers = 1023;
0803
0804
0805
0806
0807
0808
0809
0810
0811 int etaCmpr = params_->tauCompressLUT()->data(aeta);
0812 int etCmpr = params_->tauCompressLUT()->data((0x1 << 5) + Et);
0813 int nTTCmpr = params_->tauCompressLUT()->data((0x1 << 5) + (0x1 << 8) +
0814 nrTowers);
0815
0816
0817
0818 unsigned int address =
0819 ((etaCmpr << 10) | (etCmpr << 5) | nTTCmpr);
0820
0821 return address;
0822 }