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