File indexing completed on 2025-01-14 23:17:00
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013 #include "L1Trigger/L1TGlobal/interface/CorrCondition.h"
0014 #include "L1Trigger/L1TGlobal/interface/CorrThreeBodyCondition.h"
0015
0016
0017 #include <iostream>
0018 #include <iomanip>
0019
0020 #include <string>
0021 #include <vector>
0022 #include <algorithm>
0023
0024
0025
0026 #include "L1Trigger/L1TGlobal/interface/CorrelationTemplate.h"
0027 #include "L1Trigger/L1TGlobal/interface/CorrelationThreeBodyTemplate.h"
0028 #include "L1Trigger/L1TGlobal/interface/ConditionEvaluation.h"
0029
0030 #include "L1Trigger/L1TGlobal/interface/MuCondition.h"
0031 #include "L1Trigger/L1TGlobal/interface/MuonTemplate.h"
0032 #include "L1Trigger/L1TGlobal/interface/GlobalScales.h"
0033 #include "L1Trigger/L1TGlobal/interface/GlobalBoard.h"
0034
0035 #include "DataFormats/L1Trigger/interface/L1Candidate.h"
0036 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0037 #include "FWCore/MessageLogger/interface/MessageDrop.h"
0038
0039
0040
0041 l1t::CorrThreeBodyCondition::CorrThreeBodyCondition() : ConditionEvaluation() {}
0042
0043
0044 l1t::CorrThreeBodyCondition::CorrThreeBodyCondition(const GlobalCondition* corrTemplate,
0045 const GlobalCondition* cond0Condition,
0046 const GlobalCondition* cond1Condition,
0047 const GlobalCondition* cond2Condition,
0048 const GlobalBoard* ptrGTB,
0049 const GlobalScales* ptrGS)
0050 : ConditionEvaluation(),
0051 m_gtCorrelationThreeBodyTemplate(static_cast<const CorrelationThreeBodyTemplate*>(corrTemplate)),
0052 m_gtCond0(cond0Condition),
0053 m_gtCond1(cond1Condition),
0054 m_gtCond2(cond2Condition),
0055 m_uGtB(ptrGTB),
0056 m_gtScales(ptrGS) {}
0057
0058
0059 void l1t::CorrThreeBodyCondition::copy(const l1t::CorrThreeBodyCondition& cp) {
0060 m_gtCorrelationThreeBodyTemplate = cp.gtCorrelationThreeBodyTemplate();
0061 m_uGtB = cp.getuGtB();
0062 m_gtScales = cp.getScales();
0063
0064 m_condMaxNumberObjects = cp.condMaxNumberObjects();
0065 m_condLastResult = cp.condLastResult();
0066 m_combinationsInCond = cp.getCombinationsInCond();
0067
0068 m_verbosity = cp.m_verbosity;
0069 }
0070
0071 l1t::CorrThreeBodyCondition::CorrThreeBodyCondition(const l1t::CorrThreeBodyCondition& cp) : ConditionEvaluation() {
0072 copy(cp);
0073 }
0074
0075
0076 l1t::CorrThreeBodyCondition::~CorrThreeBodyCondition() {
0077
0078 }
0079
0080
0081 l1t::CorrThreeBodyCondition& l1t::CorrThreeBodyCondition::operator=(const l1t::CorrThreeBodyCondition& cp) {
0082 copy(cp);
0083 return *this;
0084 }
0085
0086
0087 void l1t::CorrThreeBodyCondition::setuGtB(const GlobalBoard* ptrGTB) { m_uGtB = ptrGTB; }
0088
0089 void l1t::CorrThreeBodyCondition::setScales(const GlobalScales* sc) { m_gtScales = sc; }
0090
0091
0092 const bool l1t::CorrThreeBodyCondition::evaluateCondition(const int bxEval) const {
0093 if (m_verbosity) {
0094 std::ostringstream myCout;
0095 m_gtCorrelationThreeBodyTemplate->print(myCout);
0096 LogDebug("L1TGlobal") << "Three-body Correlation Condition Evaluation..." << std::endl;
0097 }
0098
0099 bool condResult = false;
0100
0101
0102 int nObjInCond = 3;
0103 std::vector<GlobalObject> cndObjTypeVec(nObjInCond);
0104
0105
0106 const GtConditionCategory cond0Categ = m_gtCorrelationThreeBodyTemplate->cond0Category();
0107 const GtConditionCategory cond1Categ = m_gtCorrelationThreeBodyTemplate->cond1Category();
0108 const GtConditionCategory cond2Categ = m_gtCorrelationThreeBodyTemplate->cond2Category();
0109
0110 const MuonTemplate* corrMuon = nullptr;
0111
0112 CombinationsWithBxInCond cond0Comb{};
0113 CombinationsWithBxInCond cond1Comb{};
0114 CombinationsWithBxInCond cond2Comb{};
0115
0116
0117 bool reqObjResult = false;
0118
0119 if (cond0Categ == CondMuon) {
0120 LogDebug("L1TGlobal") << "\n --------------------- First muon checks ---------------------" << std::endl;
0121 corrMuon = static_cast<const MuonTemplate*>(m_gtCond0);
0122 MuCondition muCondition(corrMuon, m_uGtB, 0, 0);
0123
0124 muCondition.evaluateConditionStoreResult(bxEval);
0125 reqObjResult = muCondition.condLastResult();
0126
0127 cond0Comb = (muCondition.getCombinationsInCond());
0128
0129 cndObjTypeVec[0] = (corrMuon->objectType())[0];
0130
0131 if (m_verbosity) {
0132 std::ostringstream myCout;
0133 muCondition.print(myCout);
0134 LogDebug("L1TGlobal") << myCout.str() << std::endl;
0135 }
0136 } else {
0137
0138 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 0" << std::endl;
0139 return false;
0140 }
0141
0142
0143 if (!reqObjResult) {
0144 LogDebug("L1TGlobal") << "\n First subcondition false, second subcondition not evaluated and not printed."
0145 << std::endl;
0146 return false;
0147 }
0148
0149
0150 if (cond1Categ == CondMuon) {
0151 LogDebug("L1TGlobal") << "\n --------------------- Second muon checks ---------------------" << std::endl;
0152 corrMuon = static_cast<const MuonTemplate*>(m_gtCond1);
0153 MuCondition muCondition(corrMuon, m_uGtB, 0, 0);
0154
0155 muCondition.evaluateConditionStoreResult(bxEval);
0156 reqObjResult = muCondition.condLastResult();
0157
0158 cond1Comb = (muCondition.getCombinationsInCond());
0159
0160 cndObjTypeVec[1] = (corrMuon->objectType())[0];
0161
0162 if (m_verbosity) {
0163 std::ostringstream myCout;
0164 muCondition.print(myCout);
0165 LogDebug("L1TGlobal") << myCout.str() << std::endl;
0166 }
0167 }
0168
0169 else {
0170
0171 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 1" << std::endl;
0172 return false;
0173 }
0174
0175
0176 if (!reqObjResult) {
0177 LogDebug("L1TGlobal") << "\n Second subcondition false, third subcondition not evaluated and not printed."
0178 << std::endl;
0179 return false;
0180 }
0181
0182
0183 if (cond2Categ == CondMuon) {
0184 LogDebug("L1TGlobal") << "\n --------------------- Third muon checks ---------------------" << std::endl;
0185 corrMuon = static_cast<const MuonTemplate*>(m_gtCond2);
0186 MuCondition muCondition(corrMuon, m_uGtB, 0, 0);
0187
0188 muCondition.evaluateConditionStoreResult(bxEval);
0189 reqObjResult = muCondition.condLastResult();
0190
0191 cond2Comb = (muCondition.getCombinationsInCond());
0192
0193 cndObjTypeVec[2] = (corrMuon->objectType())[0];
0194
0195 if (m_verbosity) {
0196 std::ostringstream myCout;
0197 muCondition.print(myCout);
0198 LogDebug("L1TGlobal") << myCout.str() << std::endl;
0199 }
0200 }
0201
0202 else {
0203
0204 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 2" << std::endl;
0205 return false;
0206 }
0207
0208
0209 if (!reqObjResult) {
0210 return false;
0211 } else {
0212 LogDebug("L1TGlobal")
0213 << "\n"
0214 << "Found three objects satisfying subconditions: evaluate three-body correlation requirements.\n"
0215 << std::endl;
0216 }
0217
0218
0219 CorrelationThreeBodyTemplate::CorrelationThreeBodyParameter corrPar =
0220 *(m_gtCorrelationThreeBodyTemplate->correlationThreeBodyParameter());
0221
0222
0223 SingleCombWithBxInCond objectsInComb;
0224 objectsInComb.reserve(nObjInCond);
0225
0226
0227
0228 (combinationsInCond()).clear();
0229
0230
0231 const BXVector<const l1t::Muon*>* candMuVec = nullptr;
0232
0233
0234 int phiIndex0 = 0;
0235 double phi0Phy = 0.;
0236 int phiIndex1 = 0;
0237 double phi1Phy = 0.;
0238 int phiIndex2 = 0;
0239 double phi2Phy = 0.;
0240
0241 int etaIndex0 = 0;
0242 double eta0Phy = 0.;
0243 int etaBin0 = 0;
0244 int etaIndex1 = 0;
0245 double eta1Phy = 0.;
0246 int etaBin1 = 0;
0247 int etaIndex2 = 0;
0248 double eta2Phy = 0.;
0249 int etaBin2 = 0;
0250
0251 int etIndex0 = 0;
0252 int etBin0 = 0;
0253 double et0Phy = 0.;
0254 int etIndex1 = 0;
0255 int etBin1 = 0;
0256 double et1Phy = 0.;
0257 int etIndex2 = 0;
0258 int etBin2 = 0;
0259 double et2Phy = 0.;
0260
0261
0262 int chrg0 = -1;
0263 int chrg1 = -1;
0264 int chrg2 = -1;
0265
0266
0267 int phiBound = 0;
0268 if (cond0Categ == CondMuon || cond1Categ == CondMuon || cond2Categ == CondMuon) {
0269 const GlobalScales::ScaleParameters& par = m_gtScales->getMUScales();
0270 phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0271 } else {
0272
0273 const GlobalScales::ScaleParameters& par = m_gtScales->getEGScales();
0274 phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0275 }
0276 LogDebug("L1TGlobal") << "Phi Bound = " << phiBound << std::endl;
0277
0278
0279 std::string lutObj0 = "NULL";
0280 std::string lutObj1 = "NULL";
0281 std::string lutObj2 = "NULL";
0282
0283 LogTrace("L1TGlobal") << " Number of objects satisfying the subcondition 0: " << (cond0Comb.size()) << std::endl;
0284 LogTrace("L1TGlobal") << " Number of objects satisfying the subcondition 1: " << (cond1Comb.size()) << std::endl;
0285 LogTrace("L1TGlobal") << " Number of objects satisfying the subcondition 2: " << (cond2Comb.size()) << std::endl;
0286
0287
0288
0289
0290 unsigned int preShift = 0;
0291
0292
0293 for (auto const& it0Comb : cond0Comb) {
0294
0295
0296 LogDebug("L1TGlobal") << "Looking at first subcondition" << std::endl;
0297
0298 if (it0Comb.empty()) {
0299 LogTrace("L1TGlobal") << "\n SingleCombWithBxInCond it0Comb.size() " << it0Comb.size();
0300 return false;
0301 }
0302
0303 auto const cond0bx = it0Comb[0].first;
0304 auto const obj0Index = it0Comb[0].second;
0305
0306
0307 if (cond0Categ == CondMuon) {
0308 lutObj0 = "MU";
0309 candMuVec = m_uGtB->getCandL1Mu();
0310 auto const* mu0 = candMuVec->at(cond0bx, obj0Index);
0311 phiIndex0 = mu0->hwPhiAtVtx();
0312 etaIndex0 = mu0->hwEtaAtVtx();
0313 etIndex0 = mu0->hwPt();
0314 chrg0 = mu0->hwCharge();
0315
0316 etaBin0 = etaIndex0;
0317 if (etaBin0 < 0)
0318 etaBin0 = m_gtScales->getMUScales().etaBins.size() + etaBin0;
0319
0320 etBin0 = etIndex0;
0321 int ssize = m_gtScales->getMUScales().etBins.size();
0322 if (etBin0 >= ssize) {
0323 etBin0 = ssize - 1;
0324 LogTrace("L1TGlobal") << "MU0 hw et" << etBin0 << " out of scale range. Setting to maximum.";
0325 }
0326
0327
0328 std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex0);
0329 phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0330 binEdges = m_gtScales->getMUScales().etaBins.at(etaBin0);
0331 eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0332 binEdges = m_gtScales->getMUScales().etBins.at(etBin0);
0333 et0Phy = 0.5 * (binEdges.second + binEdges.first);
0334 LogDebug("L1TGlobal") << "Found all quantities for MU0" << std::endl;
0335 } else {
0336
0337 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 0" << std::endl;
0338 return false;
0339 }
0340
0341
0342 for (auto const& it1Comb : cond1Comb) {
0343 LogDebug("L1TGlobal") << "Looking at second subcondition" << std::endl;
0344
0345 if (it1Comb.empty()) {
0346 LogTrace("L1TGlobal") << "\n SingleCombWithBxInCond it1Comb.size() " << it1Comb.size();
0347 return false;
0348 }
0349
0350 auto const cond1bx = it1Comb[0].first;
0351 auto const obj1Index = it1Comb[0].second;
0352
0353
0354 if (cndObjTypeVec[0] == cndObjTypeVec[1] && obj0Index == obj1Index && cond0bx == cond1bx) {
0355 LogDebug("L1TGlobal") << "Corr Condition looking at same leg...skip" << std::endl;
0356 continue;
0357 }
0358
0359 if (cond1Categ == CondMuon) {
0360 lutObj1 = "MU";
0361 candMuVec = m_uGtB->getCandL1Mu();
0362 auto const* mu1 = candMuVec->at(cond1bx, obj1Index);
0363 phiIndex1 = mu1->hwPhiAtVtx();
0364 etaIndex1 = mu1->hwEtaAtVtx();
0365 etIndex1 = mu1->hwPt();
0366 chrg1 = mu1->hwCharge();
0367
0368 etaBin1 = etaIndex1;
0369 if (etaBin1 < 0)
0370 etaBin1 = m_gtScales->getMUScales().etaBins.size() + etaBin1;
0371
0372 etBin1 = etIndex1;
0373 int ssize = m_gtScales->getMUScales().etBins.size();
0374 if (etBin1 >= ssize) {
0375 LogTrace("L1TGlobal") << "MU1 hw et" << etBin1 << " out of scale range. Setting to maximum.";
0376 etBin1 = ssize - 1;
0377 }
0378
0379
0380 std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex1);
0381 phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0382 binEdges = m_gtScales->getMUScales().etaBins.at(etaBin1);
0383 eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0384 binEdges = m_gtScales->getMUScales().etBins.at(etBin1);
0385 et1Phy = 0.5 * (binEdges.second + binEdges.first);
0386 LogDebug("L1TGlobal") << "Found all quantities for MU1" << std::endl;
0387 } else {
0388
0389 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 1" << std::endl;
0390 return false;
0391 }
0392
0393
0394 for (auto const& it2Comb : cond2Comb) {
0395 LogDebug("L1TGlobal") << "Looking at the third object for the three-body condition" << std::endl;
0396
0397 if (it2Comb.empty()) {
0398 LogTrace("L1TGlobal") << "\n SingleCombWithBxInCond it2Comb.size() " << it2Comb.size();
0399 return false;
0400 }
0401
0402 auto const cond2bx = it2Comb[0].first;
0403 auto const obj2Index = it2Comb[0].second;
0404
0405
0406
0407 if ((cndObjTypeVec[0] == cndObjTypeVec[2] && obj0Index == obj2Index && cond0bx == cond2bx) ||
0408 (cndObjTypeVec[1] == cndObjTypeVec[2] && obj1Index == obj2Index && cond1bx == cond2bx)) {
0409 LogDebug("L1TGlobal") << "Corr Condition looking at same leg...skip" << std::endl;
0410 continue;
0411 }
0412
0413 if (cond2Categ == CondMuon) {
0414 lutObj2 = "MU";
0415 candMuVec = m_uGtB->getCandL1Mu();
0416 auto const* mu2 = candMuVec->at(cond2bx, obj2Index);
0417 phiIndex2 = mu2->hwPhiAtVtx();
0418 etaIndex2 = mu2->hwEtaAtVtx();
0419 etIndex2 = mu2->hwPt();
0420 chrg2 = mu2->hwCharge();
0421
0422 etaBin2 = etaIndex2;
0423 if (etaBin2 < 0)
0424 etaBin2 = m_gtScales->getMUScales().etaBins.size() + etaBin2;
0425
0426 etBin2 = etIndex2;
0427 int ssize = m_gtScales->getMUScales().etBins.size();
0428 if (etBin2 >= ssize) {
0429 LogTrace("L1TGlobal") << "MU2 hw et" << etBin2 << " out of scale range. Setting to maximum.";
0430 etBin2 = ssize - 1;
0431 }
0432
0433
0434 std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex2);
0435 phi2Phy = 0.5 * (binEdges.second + binEdges.first);
0436 binEdges = m_gtScales->getMUScales().etaBins.at(etaBin2);
0437 eta2Phy = 0.5 * (binEdges.second + binEdges.first);
0438 binEdges = m_gtScales->getMUScales().etBins.at(etBin2);
0439 et2Phy = 0.5 * (binEdges.second + binEdges.first);
0440 LogDebug("L1TGlobal") << "Found all quantities for MU2" << std::endl;
0441 }
0442
0443 else {
0444
0445 LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 2" << std::endl;
0446 return false;
0447 };
0448
0449 if (m_verbosity) {
0450 LogDebug("L1TGlobal") << "\n >>>>>> THREE-MUON EVENT!" << std::endl;
0451 LogDebug("L1TGlobal") << ">>>>>> Object involved in the three-body correlation condition are ["
0452 << l1t::GlobalObjectEnumToString(cndObjTypeVec[0]) << ", "
0453 << l1t::GlobalObjectEnumToString(cndObjTypeVec[1]) << ", "
0454 << l1t::GlobalObjectEnumToString(cndObjTypeVec[2]) << "] with collection indices ["
0455 << obj0Index << ", " << obj1Index << obj2Index << "] "
0456 << " having: \n"
0457 << " Et values = [" << etIndex0 << ", " << etIndex1 << ", " << etIndex2 << "]\n"
0458 << " phi indices = [" << phiIndex0 << ", " << phiIndex1 << ", " << phiIndex2
0459 << "]\n"
0460 << " eta indices = [" << etaIndex0 << ", " << etaIndex1 << ", " << etaIndex2
0461 << "]\n"
0462 << " charge values = [" << chrg0 << ", " << chrg1 << ", " << chrg2 << "]\n";
0463 }
0464
0465
0466
0467 bool reqResult = false;
0468 bool chrgCorrel = true;
0469
0470
0471
0472 if (cond0Categ == CondMuon && cond1Categ == CondMuon && cond2Categ == CondMuon) {
0473
0474 if (corrPar.chargeCorrelation == 4 && ((chrg0 + chrg1 + chrg2) == 3 || (chrg0 + chrg1 + chrg2) == 0)) {
0475 chrgCorrel = false;
0476 }
0477
0478 if (corrPar.chargeCorrelation == 2 && ((chrg0 + chrg1 + chrg2) == 1 || (chrg0 + chrg1 + chrg2) == 2)) {
0479 chrgCorrel = false;
0480 }
0481
0482 if (corrPar.chargeCorrelation == 1) {
0483 chrgCorrel = true;
0484 }
0485 }
0486
0487
0488 objectsInComb.clear();
0489 objectsInComb.emplace_back(cond0bx, obj0Index);
0490 objectsInComb.emplace_back(cond1bx, obj1Index);
0491 objectsInComb.emplace_back(cond2bx, obj2Index);
0492
0493
0494 double deltaPhiPhy_01 = fabs(phi1Phy - phi0Phy);
0495 if (deltaPhiPhy_01 > M_PI)
0496 deltaPhiPhy_01 = 2. * M_PI - deltaPhiPhy_01;
0497 double deltaEtaPhy_01 = fabs(eta1Phy - eta0Phy);
0498
0499 double deltaPhiPhy_02 = fabs(phi2Phy - phi0Phy);
0500 if (deltaPhiPhy_02 > M_PI)
0501 deltaPhiPhy_02 = 2. * M_PI - deltaPhiPhy_02;
0502 double deltaEtaPhy_02 = fabs(eta2Phy - eta0Phy);
0503
0504 double deltaPhiPhy_12 = fabs(phi2Phy - phi1Phy);
0505 if (deltaPhiPhy_12 > M_PI)
0506 deltaPhiPhy_12 = 2. * M_PI - deltaPhiPhy_12;
0507 double deltaEtaPhy_12 = fabs(eta2Phy - eta1Phy);
0508
0509
0510 int deltaPhiFW_01 = abs(phiIndex0 - phiIndex1);
0511 if (deltaPhiFW_01 >= phiBound)
0512 deltaPhiFW_01 = 2 * phiBound - deltaPhiFW_01;
0513 std::string lutName_01 = lutObj0;
0514 lutName_01 += "-";
0515 lutName_01 += lutObj1;
0516 long long deltaPhiLUT_01 = m_gtScales->getLUT_DeltaPhi(lutName_01, deltaPhiFW_01);
0517 unsigned int precDeltaPhiLUT_01 = m_gtScales->getPrec_DeltaPhi(lutName_01);
0518
0519 int deltaEtaFW_01 = abs(etaIndex0 - etaIndex1);
0520 long long deltaEtaLUT_01 = 0;
0521 unsigned int precDeltaEtaLUT_01 = 0;
0522 deltaEtaLUT_01 = m_gtScales->getLUT_DeltaEta(lutName_01, deltaEtaFW_01);
0523 precDeltaEtaLUT_01 = m_gtScales->getPrec_DeltaEta(lutName_01);
0524
0525 int deltaPhiFW_02 = abs(phiIndex0 - phiIndex2);
0526 if (deltaPhiFW_02 >= phiBound)
0527 deltaPhiFW_02 = 2 * phiBound - deltaPhiFW_02;
0528 std::string lutName_02 = lutObj0;
0529 lutName_02 += "-";
0530 lutName_02 += lutObj2;
0531 long long deltaPhiLUT_02 = m_gtScales->getLUT_DeltaPhi(lutName_02, deltaPhiFW_02);
0532 unsigned int precDeltaPhiLUT_02 = m_gtScales->getPrec_DeltaPhi(lutName_02);
0533
0534 int deltaEtaFW_02 = abs(etaIndex0 - etaIndex2);
0535 long long deltaEtaLUT_02 = 0;
0536 unsigned int precDeltaEtaLUT_02 = 0;
0537 deltaEtaLUT_02 = m_gtScales->getLUT_DeltaEta(lutName_02, deltaEtaFW_02);
0538 precDeltaEtaLUT_02 = m_gtScales->getPrec_DeltaEta(lutName_02);
0539
0540 int deltaPhiFW_12 = abs(phiIndex1 - phiIndex2);
0541 if (deltaPhiFW_12 >= phiBound)
0542 deltaPhiFW_12 = 2 * phiBound - deltaPhiFW_12;
0543 std::string lutName_12 = lutObj1;
0544 lutName_12 += "-";
0545 lutName_12 += lutObj2;
0546 long long deltaPhiLUT_12 = m_gtScales->getLUT_DeltaPhi(lutName_12, deltaPhiFW_12);
0547 unsigned int precDeltaPhiLUT_12 = m_gtScales->getPrec_DeltaPhi(lutName_12);
0548
0549 int deltaEtaFW_12 = abs(etaIndex1 - etaIndex2);
0550 long long deltaEtaLUT_12 = 0;
0551 unsigned int precDeltaEtaLUT_12 = 0;
0552 deltaEtaLUT_12 = m_gtScales->getLUT_DeltaEta(lutName_12, deltaEtaFW_12);
0553 precDeltaEtaLUT_12 = m_gtScales->getPrec_DeltaEta(lutName_12);
0554
0555
0556 LogDebug("L1TGlobal") << "### Obj0 phiFW = " << phiIndex0 << " Obj1 phiFW = " << phiIndex1 << "\n"
0557 << " DeltaPhiFW = " << deltaPhiFW_01 << " LUT Name 01= " << lutName_01
0558 << " Prec = " << precDeltaPhiLUT_01 << "\n"
0559 << " LUT Name 02= " << lutName_02 << " Prec = " << precDeltaPhiLUT_02 << "\n"
0560 << " LUT Name 12= " << lutName_12 << " Prec = " << precDeltaPhiLUT_12 << "\n"
0561 << " DeltaPhiLUT_01 = " << deltaPhiLUT_01 << "\n"
0562 << " DeltaPhiLUT_02 = " << deltaPhiLUT_02 << "\n"
0563 << " DeltaPhiLUT_12 = " << deltaPhiLUT_12 << "\n"
0564 << "### Obj0 etaFW = " << etaIndex0 << " Obj1 etaFW = " << etaIndex1 << "\n"
0565 << " DeltaEtaFW = " << deltaEtaFW_01 << " LUT Name 01 = " << lutName_01
0566 << " Prec 01 = " << precDeltaEtaLUT_01 << "\n"
0567 << " LUT Name 02 = " << lutName_02 << " Prec 02 = " << precDeltaEtaLUT_02 << "\n"
0568 << " LUT Name 12 = " << lutName_12 << " Prec 12 = " << precDeltaEtaLUT_12 << "\n"
0569 << " DeltaEtaLUT_01 = " << deltaEtaLUT_01 << " DeltaEtaLUT_02 = " << deltaEtaLUT_02
0570 << " DeltaEtaLUT_12 = " << deltaEtaLUT_12 << std::endl;
0571
0572 if (corrPar.corrCutType & 0x9) {
0573
0574
0575
0576
0577 double cosDeltaPhiPhy_01 = cos(deltaPhiPhy_01);
0578 double coshDeltaEtaPhy_01 = cosh(deltaEtaPhy_01);
0579 double massSqPhy_01 = et0Phy * et1Phy * (coshDeltaEtaPhy_01 - cosDeltaPhiPhy_01);
0580
0581 long long cosDeltaPhiLUT_01 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_01, deltaPhiFW_01);
0582 unsigned int precCosLUT_01 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_01);
0583
0584 long long coshDeltaEtaLUT_01;
0585 coshDeltaEtaLUT_01 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_01, deltaEtaFW_01);
0586 unsigned int precCoshLUT_01 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_01);
0587 if (precCoshLUT_01 - precCosLUT_01 != 0)
0588 LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different precision" << std::endl;
0589
0590 double cosDeltaPhiPhy_02 = cos(deltaPhiPhy_02);
0591 double coshDeltaEtaPhy_02 = cosh(deltaEtaPhy_02);
0592 if (corrPar.corrCutType & 0x10)
0593 coshDeltaEtaPhy_02 = 1.;
0594 double massSqPhy_02 = et0Phy * et2Phy * (coshDeltaEtaPhy_02 - cosDeltaPhiPhy_02);
0595 long long cosDeltaPhiLUT_02 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_02, deltaPhiFW_02);
0596 unsigned int precCosLUT_02 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_02);
0597 long long coshDeltaEtaLUT_02;
0598 if (corrPar.corrCutType & 0x10) {
0599 coshDeltaEtaLUT_02 = 1 * pow(10, precCosLUT_02);
0600 } else {
0601 coshDeltaEtaLUT_02 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_02, deltaEtaFW_02);
0602 unsigned int precCoshLUT_02 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_02);
0603 if (precCoshLUT_02 - precCosLUT_02 != 0)
0604 LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different precision" << std::endl;
0605 }
0606
0607 double cosDeltaPhiPhy_12 = cos(deltaPhiPhy_12);
0608 double coshDeltaEtaPhy_12 = cosh(deltaEtaPhy_12);
0609 if (corrPar.corrCutType & 0x10)
0610 coshDeltaEtaPhy_12 = 1.;
0611 double massSqPhy_12 = et1Phy * et2Phy * (coshDeltaEtaPhy_12 - cosDeltaPhiPhy_12);
0612 long long cosDeltaPhiLUT_12 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_12, deltaPhiFW_12);
0613 unsigned int precCosLUT_12 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_12);
0614 long long coshDeltaEtaLUT_12;
0615 if (corrPar.corrCutType & 0x10) {
0616 coshDeltaEtaLUT_12 = 1 * pow(10, precCosLUT_12);
0617 } else {
0618 coshDeltaEtaLUT_12 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_12, deltaEtaFW_12);
0619 unsigned int precCoshLUT_12 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_12);
0620 if (precCoshLUT_12 - precCosLUT_12 != 0)
0621 LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different precision" << std::endl;
0622 }
0623
0624 std::string lutName = lutObj0;
0625 lutName += "-ET";
0626 long long ptObj0 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex0);
0627 unsigned int precPtLUTObj0 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0628
0629 lutName = lutObj1;
0630 lutName += "-ET";
0631 long long ptObj1 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex1);
0632 unsigned int precPtLUTObj1 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0633
0634 lutName = lutObj2;
0635 lutName += "-ET";
0636 long long ptObj2 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex2);
0637 unsigned int precPtLUTObj2 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0638
0639
0640 long long massSq_01 = ptObj0 * ptObj1 * (coshDeltaEtaLUT_01 - cosDeltaPhiLUT_01);
0641 long long massSq_02 = ptObj0 * ptObj2 * (coshDeltaEtaLUT_02 - cosDeltaPhiLUT_02);
0642 long long massSq_12 = ptObj1 * ptObj2 * (coshDeltaEtaLUT_12 - cosDeltaPhiLUT_12);
0643
0644
0645
0646 unsigned int preShift_01 = precPtLUTObj0 + precPtLUTObj1 + precCosLUT_01 - corrPar.precMassCut;
0647 unsigned int preShift_02 = precPtLUTObj0 + precPtLUTObj2 + precCosLUT_02 - corrPar.precMassCut;
0648 unsigned int preShift_12 = precPtLUTObj1 + precPtLUTObj2 + precCosLUT_12 - corrPar.precMassCut;
0649
0650 LogDebug("L1TGlobal") << "####################################\n";
0651 LogDebug("L1TGlobal") << " Testing the dimuon invariant mass between the FIRST PAIR 0-1 (" << lutObj0
0652 << "," << lutObj1 << ") \n"
0653 << " massSq/2 = " << massSq_01 << "\n"
0654 << " Precision Shift = " << preShift_01 << "\n"
0655 << " massSq (shift)= " << (massSq_01 / pow(10, preShift_01 + corrPar.precMassCut))
0656 << "\n"
0657 << " massSqPhy/2 = " << massSqPhy_01
0658 << " sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_01)) << std::endl;
0659
0660 LogDebug("L1TGlobal") << "####################################\n";
0661 LogDebug("L1TGlobal") << " Testing the dimuon invariant mass between the SECOND PAIR 0-2 (" << lutObj0
0662 << "," << lutObj2 << ") \n"
0663 << " massSq/2 = " << massSq_02 << "\n"
0664 << " Precision Shift = " << preShift_02 << "\n"
0665 << " massSq (shift)= " << (massSq_02 / pow(10, preShift_02 + corrPar.precMassCut))
0666 << "\n"
0667 << " massSqPhy/2 = " << massSqPhy_02
0668 << " sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_02)) << std::endl;
0669
0670 LogDebug("L1TGlobal") << "####################################\n";
0671 LogDebug("L1TGlobal") << " Testing the dimuon invariant mass between the THIRD PAIR 1-2 (" << lutObj1
0672 << "," << lutObj2 << ") \n"
0673 << " massSq/2 = " << massSq_12 << "\n"
0674 << " Precision Shift = " << preShift_12 << "\n"
0675 << " massSq (shift)= " << (massSq_12 / pow(10, preShift_12 + corrPar.precMassCut))
0676 << "\n"
0677 << " massSqPhy/2 = " << massSqPhy_12
0678 << " sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_12)) << std::endl;
0679
0680 LogDebug("L1TGlobal") << "\n ########### THREE-BODY INVARIANT MASS #########################\n";
0681 long long massSq = 0;
0682
0683 if (preShift_01 == preShift_02 && preShift_01 == preShift_12 && preShift_02 == preShift_12) {
0684 LogDebug("L1TGlobal") << "Check the preshift value: " << preShift_01 << " = " << preShift_02 << " = "
0685 << preShift_12 << std::endl;
0686 preShift = preShift_01;
0687 } else {
0688 LogDebug("L1TGlobal")
0689 << "Preshift values considered for the sum of the dimuon invariant masses are different!" << std::endl;
0690 }
0691
0692 if ((massSq_01 != massSq_02) && (massSq_01 != massSq_12) && (massSq_02 != massSq_12)) {
0693 massSq = massSq_01 + massSq_02 + massSq_12;
0694 LogDebug("L1TGlobal") << "massSq = " << massSq << std::endl;
0695 } else {
0696 LogDebug("L1TGlobal") << "Same pair of muons considered, three-body invariant mass do not computed"
0697 << std::endl;
0698 }
0699
0700 if (massSq >= 0 && massSq >= (long long)(corrPar.minMassCutValue * pow(10, preShift)) &&
0701 massSq <= (long long)(corrPar.maxMassCutValue * pow(10, preShift))) {
0702 LogDebug("L1TGlobal") << " Passed Invariant Mass Cut ["
0703 << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
0704 << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
0705 reqResult = true;
0706 } else {
0707 LogDebug("L1TGlobal") << " Failed Invariant Mass Cut ["
0708 << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
0709 << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
0710 reqResult = false;
0711 }
0712 }
0713
0714 if (reqResult && chrgCorrel) {
0715 condResult = true;
0716 (combinationsInCond()).push_back(objectsInComb);
0717 }
0718
0719 }
0720 }
0721 }
0722
0723 if (m_verbosity && condResult) {
0724 LogDebug("L1TGlobal") << " pass(es) the correlation condition.\n" << std::endl;
0725 }
0726 return condResult;
0727 }
0728
0729
0730
0731
0732
0733
0734
0735
0736
0737
0738 const bool l1t::CorrThreeBodyCondition::checkObjectParameter(const int iCondition, const l1t::L1Candidate& cand) const {
0739 return true;
0740 }
0741
0742 void l1t::CorrThreeBodyCondition::print(std::ostream& myCout) const {
0743 myCout << "Dummy Print for CorrThreeBodyCondition" << std::endl;
0744 m_gtCorrelationThreeBodyTemplate->print(myCout);
0745
0746 ConditionEvaluation::print(myCout);
0747 }