Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-08 22:40:07

0001 /**
0002  * \class CorrThreeBodyCondition
0003  *
0004  * \orig author: Elisa Fontanesi - Boston University
0005  *               CorrCondition and CorrWithOverlapRemovalCondition classes used as a starting point
0006  *
0007  * Description: L1 Global Trigger three-body correlation conditions:                                                                                     
0008  *              evaluation of a three-body correlation condition (= three-muon invariant mass)
0009  *
0010  */
0011 
0012 // this class header
0013 #include "L1Trigger/L1TGlobal/interface/CorrCondition.h"
0014 #include "L1Trigger/L1TGlobal/interface/CorrThreeBodyCondition.h"
0015 
0016 // system include files
0017 #include <iostream>
0018 #include <iomanip>
0019 
0020 #include <string>
0021 #include <vector>
0022 #include <algorithm>
0023 
0024 // user include files
0025 //   base classes
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 // constructors
0040 //     default
0041 l1t::CorrThreeBodyCondition::CorrThreeBodyCondition() : ConditionEvaluation() {}
0042 
0043 //     from base template condition (from event setup usually)
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 // copy constructor
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 // destructor
0073 l1t::CorrThreeBodyCondition::~CorrThreeBodyCondition() {
0074   // empty
0075 }
0076 
0077 // equal operator
0078 l1t::CorrThreeBodyCondition& l1t::CorrThreeBodyCondition::operator=(const l1t::CorrThreeBodyCondition& cp) {
0079   copy(cp);
0080   return *this;
0081 }
0082 
0083 ///   set the pointer to uGT GlobalBoard
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 // try all object permutations and check spatial correlations, if required
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   bool reqObjResult = false;
0098 
0099   // number of objects in the condition (three) and their type
0100   int nObjInCond = 3;
0101   std::vector<GlobalObject> cndObjTypeVec(nObjInCond);
0102 
0103   // evaluate first the three subconditions (Type1s)
0104   const GtConditionCategory cond0Categ = m_gtCorrelationThreeBodyTemplate->cond0Category();
0105   const GtConditionCategory cond1Categ = m_gtCorrelationThreeBodyTemplate->cond1Category();
0106   const GtConditionCategory cond2Categ = m_gtCorrelationThreeBodyTemplate->cond2Category();
0107 
0108   const MuonTemplate* corrMuon = nullptr;
0109 
0110   CombinationsInCond cond0Comb;
0111   CombinationsInCond cond1Comb;
0112   CombinationsInCond cond2Comb;
0113 
0114   int cond0bx(0);
0115   int cond1bx(0);
0116   int cond2bx(0);
0117 
0118   // FIRST OBJECT
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     cond0bx = bxEval + (corrMuon->condRelativeBx());
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     // Interested only in three-muon correlations
0138     LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 0" << std::endl;
0139     return false;
0140   }
0141 
0142   // return if first subcondition is false
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   // SECOND OBJECT
0150   reqObjResult = false;
0151 
0152   if (cond1Categ == CondMuon) {
0153     LogDebug("L1TGlobal") << "\n --------------------- Second muon checks ---------------------" << std::endl;
0154     corrMuon = static_cast<const MuonTemplate*>(m_gtCond1);
0155     MuCondition muCondition(corrMuon, m_uGtB, 0, 0);
0156 
0157     muCondition.evaluateConditionStoreResult(bxEval);
0158     reqObjResult = muCondition.condLastResult();
0159 
0160     cond1Comb = (muCondition.getCombinationsInCond());
0161     cond1bx = bxEval + (corrMuon->condRelativeBx());
0162     cndObjTypeVec[1] = (corrMuon->objectType())[0];
0163 
0164     if (m_verbosity) {
0165       std::ostringstream myCout;
0166       muCondition.print(myCout);
0167       LogDebug("L1TGlobal") << myCout.str() << std::endl;
0168     }
0169   }
0170 
0171   else {
0172     // Interested only in three-muon correlations
0173     LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 1" << std::endl;
0174     return false;
0175   }
0176 
0177   // return if second subcondition is false
0178   if (!reqObjResult) {
0179     LogDebug("L1TGlobal") << "\n  Second subcondition false, third subcondition not evaluated and not printed."
0180                           << std::endl;
0181     return false;
0182   }
0183 
0184   // THIRD OBJECT
0185   reqObjResult = false;
0186 
0187   if (cond2Categ == CondMuon) {
0188     LogDebug("L1TGlobal") << "\n --------------------- Third muon checks ---------------------" << std::endl;
0189     corrMuon = static_cast<const MuonTemplate*>(m_gtCond2);
0190     MuCondition muCondition(corrMuon, m_uGtB, 0, 0);
0191 
0192     muCondition.evaluateConditionStoreResult(bxEval);
0193     reqObjResult = muCondition.condLastResult();
0194 
0195     cond2Comb = (muCondition.getCombinationsInCond());
0196     cond2bx = bxEval + (corrMuon->condRelativeBx());
0197     cndObjTypeVec[2] = (corrMuon->objectType())[0];
0198 
0199     if (m_verbosity) {
0200       std::ostringstream myCout;
0201       muCondition.print(myCout);
0202       LogDebug("L1TGlobal") << myCout.str() << std::endl;
0203     }
0204   }
0205 
0206   else {
0207     // Interested only in three-muon correlations
0208     LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 2" << std::endl;
0209     return false;
0210   }
0211 
0212   // return if third subcondition is false
0213   if (!reqObjResult) {
0214     return false;
0215   } else {
0216     LogDebug("L1TGlobal")
0217         << "\n"
0218         << "Found three objects satisfying subconditions: evaluate three-body correlation requirements.\n"
0219         << std::endl;
0220   }
0221 
0222   // since we have three good legs get the correlation parameters
0223   CorrelationThreeBodyTemplate::CorrelationThreeBodyParameter corrPar =
0224       *(m_gtCorrelationThreeBodyTemplate->correlationThreeBodyParameter());
0225 
0226   // vector to store the indices of the objects involved in the condition evaluation
0227   SingleCombInCond objectsInComb;
0228   objectsInComb.reserve(nObjInCond);
0229 
0230   // clear the m_combinationsInCond vector:
0231   // it will store the set of objects satisfying the condition evaluated as true
0232   (combinationsInCond()).clear();
0233 
0234   // pointers to objects
0235   const BXVector<const l1t::Muon*>* candMuVec = nullptr;
0236 
0237   // make the conversions of the indices, depending on the combination of objects involved
0238   int phiIndex0 = 0;
0239   double phi0Phy = 0.;
0240   int phiIndex1 = 0;
0241   double phi1Phy = 0.;
0242   int phiIndex2 = 0;
0243   double phi2Phy = 0.;
0244 
0245   int etaIndex0 = 0;
0246   double eta0Phy = 0.;
0247   int etaBin0 = 0;
0248   int etaIndex1 = 0;
0249   double eta1Phy = 0.;
0250   int etaBin1 = 0;
0251   int etaIndex2 = 0;
0252   double eta2Phy = 0.;
0253   int etaBin2 = 0;
0254 
0255   int etIndex0 = 0;
0256   int etBin0 = 0;
0257   double et0Phy = 0.;
0258   int etIndex1 = 0;
0259   int etBin1 = 0;
0260   double et1Phy = 0.;
0261   int etIndex2 = 0;
0262   int etBin2 = 0;
0263   double et2Phy = 0.;
0264 
0265   // Determine the number of phi bins to get cutoff at pi
0266   int phiBound = 0;
0267   if (cond0Categ == CondMuon || cond1Categ == CondMuon || cond2Categ == CondMuon) {
0268     const GlobalScales::ScaleParameters& par = m_gtScales->getMUScales();
0269     phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0270   } else {
0271     //Assumes all objects are on same phi scale
0272     const GlobalScales::ScaleParameters& par = m_gtScales->getEGScales();
0273     phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0274   }
0275   LogDebug("L1TGlobal") << "Phi Bound = " << phiBound << std::endl;
0276 
0277   // Keep track of objects for LUTS
0278   std::string lutObj0 = "NULL";
0279   std::string lutObj1 = "NULL";
0280   std::string lutObj2 = "NULL";
0281 
0282   LogTrace("L1TGlobal") << "  Number of objects satisfying the subcondition 0: " << (cond0Comb.size()) << std::endl;
0283   LogTrace("L1TGlobal") << "  Number of objects satisfying the subcondition 1: " << (cond1Comb.size()) << std::endl;
0284   LogTrace("L1TGlobal") << "  Number of objects satisfying the subcondition 2: " << (cond2Comb.size()) << std::endl;
0285 
0286   ////////////////////////////////
0287   // LOOP OVER ALL COMBINATIONS //
0288   ////////////////////////////////
0289   unsigned int preShift = 0;
0290 
0291   // *** Looking for a set of three objects
0292   for (std::vector<SingleCombInCond>::const_iterator it0Comb = cond0Comb.begin(); it0Comb != cond0Comb.end();
0293        it0Comb++) {
0294     // Type1s: there is 1 object only, no need for a loop, index 0 should be OK in (*it0Comb)[0]
0295     // ... but add protection to not crash
0296     LogDebug("L1TGlobal") << "Looking at first subcondition" << std::endl;
0297     int obj0Index = -1;
0298 
0299     if (!(*it0Comb).empty()) {
0300       obj0Index = (*it0Comb)[0];
0301     } else {
0302       LogTrace("L1TGlobal") << "\n  SingleCombInCond (*it0Comb).size() " << ((*it0Comb).size()) << std::endl;
0303       return false;
0304     }
0305 
0306     // FIRST OBJECT: Collect the information on the first leg of the correlation
0307     if (cond0Categ == CondMuon) {
0308       lutObj0 = "MU";
0309       candMuVec = m_uGtB->getCandL1Mu();
0310       phiIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwPhiAtVtx();  //(*candMuVec)[obj0Index]->phiIndex();
0311       etaIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwEtaAtVtx();
0312       etIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwPt();
0313       etaBin0 = etaIndex0;
0314       if (etaBin0 < 0)
0315         etaBin0 = m_gtScales->getMUScales().etaBins.size() + etaBin0;  //twos complement
0316       // LogDebug("L1TGlobal") << "Muon phi" << phiIndex0 << " eta " << etaIndex0 << " etaBin0 = " << etaBin0  << " et " << etIndex0 << std::endl;
0317 
0318       etBin0 = etIndex0;
0319       int ssize = m_gtScales->getMUScales().etBins.size();
0320       if (etBin0 >= ssize) {
0321         etBin0 = ssize - 1;
0322         LogTrace("L1TGlobal") << "MU0 hw et" << etBin0 << " out of scale range.  Setting to maximum.";
0323       }
0324 
0325       // Determine Floating Pt numbers for floating point caluclation
0326       std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex0);
0327       phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0328       binEdges = m_gtScales->getMUScales().etaBins.at(etaBin0);
0329       eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0330       binEdges = m_gtScales->getMUScales().etBins.at(etBin0);
0331       et0Phy = 0.5 * (binEdges.second + binEdges.first);
0332       LogDebug("L1TGlobal") << "Found all quantities for MU0" << std::endl;
0333     } else {
0334       // Interested only in three-muon correlations
0335       LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 0" << std::endl;
0336       return false;
0337     }
0338 
0339     // SECOND OBJECT: Now loop over the second leg to get its information
0340     for (std::vector<SingleCombInCond>::const_iterator it1Comb = cond1Comb.begin(); it1Comb != cond1Comb.end();
0341          it1Comb++) {
0342       LogDebug("L1TGlobal") << "Looking at second subdondition" << std::endl;
0343       int obj1Index = -1;
0344 
0345       if (!(*it1Comb).empty()) {
0346         obj1Index = (*it1Comb)[0];
0347       } else {
0348         LogTrace("L1TGlobal") << "\n  SingleCombInCond (*it1Comb).size() " << ((*it1Comb).size()) << std::endl;
0349         return false;
0350       }
0351 
0352       //If we are dealing with the same object type avoid the two legs
0353       // either being the same object
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         phiIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwPhiAtVtx();  //(*candMuVec)[obj0Index]->phiIndex();
0363         etaIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwEtaAtVtx();
0364         etIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwPt();
0365         etaBin1 = etaIndex1;
0366         if (etaBin1 < 0)
0367           etaBin1 = m_gtScales->getMUScales().etaBins.size() + etaBin1;
0368         // LogDebug("L1TGlobal") << "Muon phi" << phiIndex1 << " eta " << etaIndex1 << " etaBin1 = " << etaBin1  << " et " << etIndex1 << std::endl;
0369 
0370         etBin1 = etIndex1;
0371         int ssize = m_gtScales->getMUScales().etBins.size();
0372         if (etBin1 >= ssize) {
0373           LogTrace("L1TGlobal") << "MU1 hw et" << etBin1 << " out of scale range.  Setting to maximum.";
0374           etBin1 = ssize - 1;
0375         }
0376 
0377         // Determine Floating Pt numbers for floating point calculation
0378         std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex1);
0379         phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0380         binEdges = m_gtScales->getMUScales().etaBins.at(etaBin1);
0381         eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0382         binEdges = m_gtScales->getMUScales().etBins.at(etBin1);
0383         et1Phy = 0.5 * (binEdges.second + binEdges.first);
0384         LogDebug("L1TGlobal") << "Found all quantities for MU1" << std::endl;
0385       } else {
0386         // Interested only in three-muon correlations
0387         LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 1" << std::endl;
0388         return false;
0389       }
0390 
0391       // THIRD OBJECT: Finally loop over the third leg to get its information
0392       for (std::vector<SingleCombInCond>::const_iterator it2Comb = cond2Comb.begin(); it2Comb != cond2Comb.end();
0393            it2Comb++) {
0394         LogDebug("L1TGlobal") << "Looking at the third object for the three-body condition" << std::endl;
0395         int obj2Index = -1;
0396 
0397         if (!(*it2Comb).empty()) {
0398           obj2Index = (*it2Comb)[0];
0399         } else {
0400           LogTrace("L1TGlobal") << "\n  SingleCombInCond (*it2Comb).size() " << ((*it2Comb).size()) << std::endl;
0401           return false;
0402         }
0403 
0404         //If we are dealing with the same object type avoid the two legs
0405         // either being the same object
0406         if ((cndObjTypeVec[0] == cndObjTypeVec[2] && obj0Index == obj2Index && cond0bx == cond2bx) ||
0407             (cndObjTypeVec[1] == cndObjTypeVec[2] && obj1Index == obj2Index && cond1bx == cond2bx)) {
0408           LogDebug("L1TGlobal") << "Corr Condition looking at same leg...skip" << std::endl;
0409           continue;
0410         }
0411 
0412         if (cond2Categ == CondMuon) {
0413           lutObj2 = "MU";
0414           candMuVec = m_uGtB->getCandL1Mu();
0415           phiIndex2 = (candMuVec->at(cond2bx, obj2Index))->hwPhiAtVtx();  //(*candMuVec)[obj0Index]->phiIndex();
0416           etaIndex2 = (candMuVec->at(cond2bx, obj2Index))->hwEtaAtVtx();
0417           etIndex2 = (candMuVec->at(cond2bx, obj2Index))->hwPt();
0418           etaBin2 = etaIndex2;
0419           if (etaBin2 < 0)
0420             etaBin2 = m_gtScales->getMUScales().etaBins.size() + etaBin2;
0421 
0422           etBin2 = etIndex2;
0423           int ssize = m_gtScales->getMUScales().etBins.size();
0424           if (etBin2 >= ssize) {
0425             LogTrace("L1TGlobal") << "MU2 hw et" << etBin2 << " out of scale range.  Setting to maximum.";
0426             etBin2 = ssize - 1;
0427           }
0428 
0429           // Determine Floating Pt numbers for floating point calculation
0430           std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex2);
0431           phi2Phy = 0.5 * (binEdges.second + binEdges.first);
0432           binEdges = m_gtScales->getMUScales().etaBins.at(etaBin2);
0433           eta2Phy = 0.5 * (binEdges.second + binEdges.first);
0434           binEdges = m_gtScales->getMUScales().etBins.at(etBin2);
0435           et2Phy = 0.5 * (binEdges.second + binEdges.first);
0436           LogDebug("L1TGlobal") << "Found all quantities for MU2" << std::endl;
0437         }
0438 
0439         else {
0440           // Interested only in three-muon correlations
0441           LogDebug("L1TGlobal") << "CondMuon not satisfied for Leg 2" << std::endl;
0442           return false;
0443         };
0444 
0445         if (m_verbosity) {
0446           LogDebug("L1TGlobal") << "\n >>>>>> THREE-MUON EVENT!" << std::endl;
0447           LogDebug("L1TGlobal") << ">>>>>> Object involved in the three-body correlation condition are ["
0448                                 << l1TGtObjectEnumToString(cndObjTypeVec[0]) << ", "
0449                                 << l1TGtObjectEnumToString(cndObjTypeVec[1]) << ", "
0450                                 << l1TGtObjectEnumToString(cndObjTypeVec[2]) << "] with collection indices ["
0451                                 << obj0Index << ", " << obj1Index << obj2Index << "] "
0452                                 << " having: \n"
0453                                 << "     Et  values  = [" << etIndex0 << ", " << etIndex1 << ", " << etIndex2 << "]\n"
0454                                 << "     phi indices = [" << phiIndex0 << ", " << phiIndex1 << ", " << phiIndex2
0455                                 << "]\n"
0456                                 << "     eta indices = [" << etaIndex0 << ", " << etaIndex1 << ", " << etaIndex2
0457                                 << "]\n"
0458                                 << std::endl;
0459         }
0460 
0461         // Now perform the desired correlation on these three objects:
0462         //reqResult will be set true in case all checks were successful for a given combination of three muons
0463         bool reqResult = false;
0464 
0465         // Clear the vector containing indices of the objects of the combination involved in the condition evaluation
0466         objectsInComb.clear();
0467         objectsInComb.push_back(obj0Index);
0468         objectsInComb.push_back(obj1Index);
0469         objectsInComb.push_back(obj2Index);
0470 
0471         // Delta eta and phi calculations needed to evaluate the three-body invariant mass
0472         double deltaPhiPhy_01 = fabs(phi1Phy - phi0Phy);
0473         if (deltaPhiPhy_01 > M_PI)
0474           deltaPhiPhy_01 = 2. * M_PI - deltaPhiPhy_01;
0475         double deltaEtaPhy_01 = fabs(eta1Phy - eta0Phy);
0476 
0477         double deltaPhiPhy_02 = fabs(phi2Phy - phi0Phy);
0478         if (deltaPhiPhy_02 > M_PI)
0479           deltaPhiPhy_02 = 2. * M_PI - deltaPhiPhy_02;
0480         double deltaEtaPhy_02 = fabs(eta2Phy - eta0Phy);
0481 
0482         double deltaPhiPhy_12 = fabs(phi2Phy - phi1Phy);
0483         if (deltaPhiPhy_12 > M_PI)
0484           deltaPhiPhy_12 = 2. * M_PI - deltaPhiPhy_12;
0485         double deltaEtaPhy_12 = fabs(eta2Phy - eta1Phy);
0486 
0487         // Determine the integer based delta eta and delta phi
0488         int deltaPhiFW_01 = abs(phiIndex0 - phiIndex1);
0489         if (deltaPhiFW_01 >= phiBound)
0490           deltaPhiFW_01 = 2 * phiBound - deltaPhiFW_01;
0491         std::string lutName_01 = lutObj0;
0492         lutName_01 += "-";
0493         lutName_01 += lutObj1;
0494         long long deltaPhiLUT_01 = m_gtScales->getLUT_DeltaPhi(lutName_01, deltaPhiFW_01);
0495         unsigned int precDeltaPhiLUT_01 = m_gtScales->getPrec_DeltaPhi(lutName_01);
0496 
0497         int deltaEtaFW_01 = abs(etaIndex0 - etaIndex1);
0498         long long deltaEtaLUT_01 = 0;
0499         unsigned int precDeltaEtaLUT_01 = 0;
0500         deltaEtaLUT_01 = m_gtScales->getLUT_DeltaEta(lutName_01, deltaEtaFW_01);
0501         precDeltaEtaLUT_01 = m_gtScales->getPrec_DeltaEta(lutName_01);
0502         ///
0503         int deltaPhiFW_02 = abs(phiIndex0 - phiIndex2);
0504         if (deltaPhiFW_02 >= phiBound)
0505           deltaPhiFW_02 = 2 * phiBound - deltaPhiFW_02;
0506         std::string lutName_02 = lutObj0;
0507         lutName_02 += "-";
0508         lutName_02 += lutObj2;
0509         long long deltaPhiLUT_02 = m_gtScales->getLUT_DeltaPhi(lutName_02, deltaPhiFW_02);
0510         unsigned int precDeltaPhiLUT_02 = m_gtScales->getPrec_DeltaPhi(lutName_02);
0511 
0512         int deltaEtaFW_02 = abs(etaIndex0 - etaIndex2);
0513         long long deltaEtaLUT_02 = 0;
0514         unsigned int precDeltaEtaLUT_02 = 0;
0515         deltaEtaLUT_02 = m_gtScales->getLUT_DeltaEta(lutName_02, deltaEtaFW_02);
0516         precDeltaEtaLUT_02 = m_gtScales->getPrec_DeltaEta(lutName_02);
0517         ///
0518         int deltaPhiFW_12 = abs(phiIndex1 - phiIndex2);
0519         if (deltaPhiFW_12 >= phiBound)
0520           deltaPhiFW_12 = 2 * phiBound - deltaPhiFW_12;
0521         std::string lutName_12 = lutObj1;
0522         lutName_12 += "-";
0523         lutName_12 += lutObj2;
0524         long long deltaPhiLUT_12 = m_gtScales->getLUT_DeltaPhi(lutName_12, deltaPhiFW_12);
0525         unsigned int precDeltaPhiLUT_12 = m_gtScales->getPrec_DeltaPhi(lutName_12);
0526 
0527         int deltaEtaFW_12 = abs(etaIndex1 - etaIndex2);
0528         long long deltaEtaLUT_12 = 0;
0529         unsigned int precDeltaEtaLUT_12 = 0;
0530         deltaEtaLUT_12 = m_gtScales->getLUT_DeltaEta(lutName_12, deltaEtaFW_12);
0531         precDeltaEtaLUT_12 = m_gtScales->getPrec_DeltaEta(lutName_12);
0532         ///
0533 
0534         LogDebug("L1TGlobal") << "### Obj0 phiFW = " << phiIndex0 << " Obj1 phiFW = " << phiIndex1 << "\n"
0535                               << "    DeltaPhiFW = " << deltaPhiFW_01 << "    LUT Name 01= " << lutName_01
0536                               << " Prec = " << precDeltaPhiLUT_01 << "\n"
0537                               << "    LUT Name 02= " << lutName_02 << " Prec = " << precDeltaPhiLUT_02 << "\n"
0538                               << "    LUT Name 12= " << lutName_12 << " Prec = " << precDeltaPhiLUT_12 << "\n"
0539                               << "    DeltaPhiLUT_01 = " << deltaPhiLUT_01 << "\n"
0540                               << "    DeltaPhiLUT_02 = " << deltaPhiLUT_02 << "\n"
0541                               << "    DeltaPhiLUT_12 = " << deltaPhiLUT_12 << "\n"
0542                               << "### Obj0 etaFW = " << etaIndex0 << " Obj1 etaFW = " << etaIndex1 << "\n"
0543                               << "    DeltaEtaFW = " << deltaEtaFW_01 << "    LUT Name 01 = " << lutName_01
0544                               << " Prec 01 = " << precDeltaEtaLUT_01 << "\n"
0545                               << "    LUT Name 02 = " << lutName_02 << " Prec 02 = " << precDeltaEtaLUT_02 << "\n"
0546                               << "    LUT Name 12 = " << lutName_12 << " Prec 12 = " << precDeltaEtaLUT_12 << "\n"
0547                               << "    DeltaEtaLUT_01 = " << deltaEtaLUT_01 << "    DeltaEtaLUT_02 = " << deltaEtaLUT_02
0548                               << "    DeltaEtaLUT_12 = " << deltaEtaLUT_12 << std::endl;
0549 
0550         if (corrPar.corrCutType & 0x9) {
0551           //invariant mass calculation based for each pair on
0552           // M = sqrt(2*p1*p2(cosh(eta1-eta2) - cos(phi1 - phi2)))
0553           // but we calculate (1/2)M^2
0554           //
0555           double cosDeltaPhiPhy_01 = cos(deltaPhiPhy_01);
0556           double coshDeltaEtaPhy_01 = cosh(deltaEtaPhy_01);
0557           double massSqPhy_01 = et0Phy * et1Phy * (coshDeltaEtaPhy_01 - cosDeltaPhiPhy_01);
0558 
0559           long long cosDeltaPhiLUT_01 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_01, deltaPhiFW_01);
0560           unsigned int precCosLUT_01 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_01);
0561 
0562           long long coshDeltaEtaLUT_01;
0563           coshDeltaEtaLUT_01 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_01, deltaEtaFW_01);
0564           unsigned int precCoshLUT_01 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_01);
0565           if (precCoshLUT_01 - precCosLUT_01 != 0)
0566             LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different Precision" << std::endl;
0567 
0568           double cosDeltaPhiPhy_02 = cos(deltaPhiPhy_02);
0569           double coshDeltaEtaPhy_02 = cosh(deltaEtaPhy_02);
0570           if (corrPar.corrCutType & 0x10)
0571             coshDeltaEtaPhy_02 = 1.;
0572           double massSqPhy_02 = et0Phy * et2Phy * (coshDeltaEtaPhy_02 - cosDeltaPhiPhy_02);
0573           long long cosDeltaPhiLUT_02 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_02, deltaPhiFW_02);
0574           unsigned int precCosLUT_02 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_02);
0575           long long coshDeltaEtaLUT_02;
0576           if (corrPar.corrCutType & 0x10) {
0577             coshDeltaEtaLUT_02 = 1 * pow(10, precCosLUT_02);
0578           } else {
0579             coshDeltaEtaLUT_02 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_02, deltaEtaFW_02);
0580             unsigned int precCoshLUT_02 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_02);
0581             if (precCoshLUT_02 - precCosLUT_02 != 0)
0582               LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different Precision" << std::endl;
0583           }
0584 
0585           double cosDeltaPhiPhy_12 = cos(deltaPhiPhy_12);
0586           double coshDeltaEtaPhy_12 = cosh(deltaEtaPhy_12);
0587           if (corrPar.corrCutType & 0x10)
0588             coshDeltaEtaPhy_12 = 1.;
0589           double massSqPhy_12 = et1Phy * et2Phy * (coshDeltaEtaPhy_12 - cosDeltaPhiPhy_12);
0590           long long cosDeltaPhiLUT_12 = m_gtScales->getLUT_DeltaPhi_Cos(lutName_12, deltaPhiFW_12);
0591           unsigned int precCosLUT_12 = m_gtScales->getPrec_DeltaPhi_Cos(lutName_12);
0592           long long coshDeltaEtaLUT_12;
0593           if (corrPar.corrCutType & 0x10) {
0594             coshDeltaEtaLUT_12 = 1 * pow(10, precCosLUT_12);
0595           } else {
0596             coshDeltaEtaLUT_12 = m_gtScales->getLUT_DeltaEta_Cosh(lutName_12, deltaEtaFW_12);
0597             unsigned int precCoshLUT_12 = m_gtScales->getPrec_DeltaEta_Cosh(lutName_12);
0598             if (precCoshLUT_12 - precCosLUT_12 != 0)
0599               LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different Precision" << std::endl;
0600           }
0601 
0602           std::string lutName = lutObj0;
0603           lutName += "-ET";
0604           long long ptObj0 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex0);
0605           unsigned int precPtLUTObj0 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0606 
0607           lutName = lutObj1;
0608           lutName += "-ET";
0609           long long ptObj1 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex1);
0610           unsigned int precPtLUTObj1 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0611 
0612           lutName = lutObj2;
0613           lutName += "-ET";
0614           long long ptObj2 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex2);
0615           unsigned int precPtLUTObj2 = m_gtScales->getPrec_Pt("Mass_" + lutName);
0616 
0617           // Pt and Angles are at different precission.
0618           long long massSq_01 = ptObj0 * ptObj1 * (coshDeltaEtaLUT_01 - cosDeltaPhiLUT_01);
0619           long long massSq_02 = ptObj0 * ptObj2 * (coshDeltaEtaLUT_02 - cosDeltaPhiLUT_02);
0620           long long massSq_12 = ptObj1 * ptObj2 * (coshDeltaEtaLUT_12 - cosDeltaPhiLUT_12);
0621 
0622           //Note: There is an assumption here that Cos and Cosh have the same precission
0623           //unsigned int preShift_01 = precPtLUTObj0 + precPtLUTObj1 + precCosLUT - corrPar.precMassCut;
0624           unsigned int preShift_01 = precPtLUTObj0 + precPtLUTObj1 + precCosLUT_01 - corrPar.precMassCut;
0625           unsigned int preShift_02 = precPtLUTObj0 + precPtLUTObj2 + precCosLUT_02 - corrPar.precMassCut;
0626           unsigned int preShift_12 = precPtLUTObj1 + precPtLUTObj2 + precCosLUT_12 - corrPar.precMassCut;
0627 
0628           LogDebug("L1TGlobal") << "####################################\n";
0629           LogDebug("L1TGlobal")
0630               << "    Testing the dimuon invariant mass between the FIRST PAIR 0-1 (" << lutObj0 << "," << lutObj1
0631               << ") \n"
0632               //<< (long long)(corrPar.minMassCutValue * pow(10, preShift_01)) << ","
0633               //<< (long long)(corrPar.maxMassCutValue * pow(10, preShift_01))
0634               //<< "] with precision = " << corrPar.precMassCut << "\n"
0635               //<< "    deltaPhiLUT  = " << deltaPhiLUT_01 << "  cosLUT  = " << cosDeltaPhiLUT_01 << "\n"
0636               //<< "    deltaEtaLUT  = " << deltaEtaLUT_01 << "  coshLUT = " << coshDeltaEtaLUT_01 << "\n"
0637               //<< "    etIndex0     = " << etIndex0 << "    pt0LUT      = " << ptObj0
0638               //<< " PhyEt0 = " << et0Phy << "\n"
0639               //<< "    etIndex1     = " << etIndex1 << "    pt1LUT      = " << ptObj1
0640               //<< " PhyEt1 = " << et1Phy << "\n"
0641               << "    massSq/2     = " << massSq_01 << "\n"
0642               << "    Precision Shift = " << preShift_01 << "\n"
0643               << "    massSq   (shift)= " << (massSq_01 / pow(10, preShift_01 + corrPar.precMassCut))
0644               << "\n"
0645               //<< "    deltaPhiPhy  = " << deltaPhiPhy_01 << "  cos() = " << cosDeltaPhiPhy_01 << "\n"
0646               //<< "    deltaEtaPhy  = " << deltaEtaPhy_01 << "  cosh()= " << coshDeltaEtaPhy_01 << "\n"
0647               << "    massSqPhy/2  = " << massSqPhy_01 << "  sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_01))
0648               << std::endl;
0649 
0650           LogDebug("L1TGlobal") << "####################################\n";
0651           LogDebug("L1TGlobal")
0652               << "    Testing the dimuon invariant mass between the SECOND PAIR 0-2 (" << lutObj0 << "," << lutObj2
0653               << ") \n"
0654               //<< (long long)(corrPar.minMassCutValue * pow(10, preShift_02)) << ","
0655               //<< (long long)(corrPar.maxMassCutValue * pow(10, preShift_02))
0656               //<< "] with precision = " << corrPar.precMassCut << "\n"
0657               //<< "    deltaPhiLUT  = " << deltaPhiLUT_02 << "  cosLUT  = " << cosDeltaPhiLUT_02 << "\n"
0658               //<< "    deltaEtaLUT  = " << deltaEtaLUT_02 << "  coshLUT = " << coshDeltaEtaLUT_02 << "\n"
0659               //<< "    etIndex0     = " << etIndex0 << "    pt0LUT      = " << ptObj0
0660               //<< " PhyEt0 = " << et0Phy << "\n"
0661               //<< "    etIndex2     = " << etIndex2 << "    pt2LUT      = " << ptObj2
0662               //<< " PhyEt2 = " << et2Phy << "\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               //<< "    deltaPhiPhy  = " << deltaPhiPhy_02 << "  cos() = " << cosDeltaPhiPhy_02 << "\n"
0668               //<< "    deltaEtaPhy  = " << deltaEtaPhy_02 << "  cosh()= " << coshDeltaEtaPhy_02 << "\n"
0669               << "    massSqPhy/2  = " << massSqPhy_02 << "  sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_02))
0670               << std::endl;
0671 
0672           LogDebug("L1TGlobal") << "####################################\n";
0673           LogDebug("L1TGlobal")
0674               << "    Testing the dimuon invariant mass between the THIRD PAIR 1-2 (" << lutObj1 << "," << lutObj2
0675               << ") \n"
0676               //<< (long long)(corrPar.minMassCutValue * pow(10, preShift_12)) << ","
0677               //<< (long long)(corrPar.maxMassCutValue * pow(10, preShift_12))
0678               //<< "] with precision = " << corrPar.precMassCut << "\n"
0679               //<< "    deltaPhiLUT  = " << deltaPhiLUT_12 << "  cosLUT  = " << cosDeltaPhiLUT_12 << "\n"
0680               //<< "    deltaEtaLUT  = " << deltaEtaLUT_12 << "  coshLUT = " << coshDeltaEtaLUT_12 << "\n"
0681               //<< "    etIndex1     = " << etIndex1 << "    pt1LUT      = " << ptObj1
0682               //<< " PhyEt1 = " << et0Phy << "\n"
0683               //<< "    etIndex2     = " << etIndex2 << "    pt2LUT      = " << ptObj2
0684               //<< " PhyEt2 = " << et2Phy << "\n"
0685               << "    massSq/2     = " << massSq_12 << "\n"
0686               << "    Precision Shift = " << preShift_12 << "\n"
0687               << "    massSq   (shift)= " << (massSq_12 / pow(10, preShift_12 + corrPar.precMassCut))
0688               << "\n"
0689               //<< "    deltaPhiPhy  = " << deltaPhiPhy_12 << "  cos() = " << cosDeltaPhiPhy_12 << "\n"
0690               //<< "    deltaEtaPhy  = " << deltaEtaPhy_12 << "  cosh()= " << coshDeltaEtaPhy_12 << "\n"
0691               << "    massSqPhy/2  = " << massSqPhy_12 << "  sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy_12))
0692               << std::endl;
0693 
0694           LogDebug("L1TGlobal") << "\n ########### THREE-BODY INVARIANT MASS #########################\n";
0695           long long massSq = 0;
0696 
0697           if (preShift_01 == preShift_02 && preShift_01 == preShift_12 && preShift_02 == preShift_12) {
0698             LogDebug("L1TGlobal") << "Check the preshift value: " << preShift_01 << " = " << preShift_02 << " = "
0699                                   << preShift_12 << std::endl;
0700             preShift = preShift_01;
0701           } else {
0702             LogDebug("L1TGlobal")
0703                 << "Preshift values considered for the sum of the dimuon invariant masses are different!" << std::endl;
0704           }
0705 
0706           if ((massSq_01 != massSq_02) && (massSq_01 != massSq_12) && (massSq_02 != massSq_12)) {
0707             massSq = massSq_01 + massSq_02 + massSq_12;
0708             LogDebug("L1TGlobal") << "massSq = " << massSq << std::endl;
0709           } else {
0710             LogDebug("L1TGlobal") << "Same pair of muons considered, three-body invariant mass do not computed"
0711                                   << std::endl;
0712           }
0713 
0714           if (massSq >= 0 && massSq >= (long long)(corrPar.minMassCutValue * pow(10, preShift)) &&
0715               massSq <= (long long)(corrPar.maxMassCutValue * pow(10, preShift))) {
0716             LogDebug("L1TGlobal") << "    Passed Invariant Mass Cut ["
0717                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
0718                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
0719             reqResult = true;
0720           } else {
0721             LogDebug("L1TGlobal") << "    Failed Invariant Mass Cut ["
0722                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
0723                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
0724             reqResult = false;
0725           }
0726         }
0727 
0728         if (reqResult) {
0729           condResult = true;
0730           (combinationsInCond()).push_back(objectsInComb);
0731         }
0732 
0733       }  //end loop over third leg
0734     }    //end loop over second leg
0735   }      //end loop over first leg
0736 
0737   if (m_verbosity && condResult) {
0738     LogDebug("L1TGlobal") << " pass(es) the correlation condition.\n" << std::endl;
0739   }
0740   return condResult;
0741 }
0742 
0743 /**
0744  * checkObjectParameter - Compare a single particle with a numbered condition.
0745  *
0746  * @param iCondition The number of the condition.
0747  * @param cand The candidate to compare.
0748  *
0749  * @return The result of the comparison (false if a condition does not exist).
0750  */
0751 
0752 const bool l1t::CorrThreeBodyCondition::checkObjectParameter(const int iCondition, const l1t::L1Candidate& cand) const {
0753   return true;
0754 }
0755 
0756 void l1t::CorrThreeBodyCondition::print(std::ostream& myCout) const {
0757   myCout << "Dummy Print for CorrThreeBodyCondition" << std::endl;
0758   m_gtCorrelationThreeBodyTemplate->print(myCout);
0759 
0760   ConditionEvaluation::print(myCout);
0761 }