Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:55:07

0001 /**
0002  * \class CorrCondition
0003  *
0004  *
0005  * Description: evaluation of a correlation condition.
0006  *
0007  * Implementation:
0008  *    <TODO: enter implementation details>
0009  *
0010  *\new features: Dragana Pilipovic
0011  *                - updated for invariant mass over delta R condition*
0012  */
0013 
0014 // this class header
0015 #include "L1Trigger/L1TGlobal/interface/CorrCondition.h"
0016 
0017 // system include files
0018 #include <iostream>
0019 #include <iomanip>
0020 
0021 #include <string>
0022 #include <vector>
0023 #include <algorithm>
0024 
0025 // user include files
0026 //   base classes
0027 #include "L1Trigger/L1TGlobal/interface/CorrelationTemplate.h"
0028 #include "L1Trigger/L1TGlobal/interface/ConditionEvaluation.h"
0029 
0030 #include "L1Trigger/L1TGlobal/interface/MuCondition.h"
0031 #include "L1Trigger/L1TGlobal/interface/CaloCondition.h"
0032 #include "L1Trigger/L1TGlobal/interface/EnergySumCondition.h"
0033 #include "L1Trigger/L1TGlobal/interface/MuonTemplate.h"
0034 #include "L1Trigger/L1TGlobal/interface/CaloTemplate.h"
0035 #include "L1Trigger/L1TGlobal/interface/EnergySumTemplate.h"
0036 #include "L1Trigger/L1TGlobal/interface/GlobalScales.h"
0037 
0038 #include "DataFormats/L1Trigger/interface/L1Candidate.h"
0039 
0040 #include "L1Trigger/L1TGlobal/interface/GlobalBoard.h"
0041 
0042 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0043 #include "FWCore/MessageLogger/interface/MessageDrop.h"
0044 
0045 // constructors
0046 //     default
0047 l1t::CorrCondition::CorrCondition() : ConditionEvaluation() {}
0048 
0049 //     from base template condition (from event setup usually)
0050 l1t::CorrCondition::CorrCondition(const GlobalCondition* corrTemplate,
0051                                   const GlobalCondition* cond0Condition,
0052                                   const GlobalCondition* cond1Condition,
0053                                   const GlobalBoard* ptrGTB)
0054     : ConditionEvaluation(),
0055       m_gtCorrelationTemplate(static_cast<const CorrelationTemplate*>(corrTemplate)),
0056       m_gtCond0(cond0Condition),
0057       m_gtCond1(cond1Condition),
0058       m_uGtB(ptrGTB) {}
0059 
0060 // copy constructor
0061 void l1t::CorrCondition::copy(const l1t::CorrCondition& cp) {
0062   m_gtCorrelationTemplate = cp.gtCorrelationTemplate();
0063   m_uGtB = cp.getuGtB();
0064 
0065   m_condMaxNumberObjects = cp.condMaxNumberObjects();
0066   m_condLastResult = cp.condLastResult();
0067   m_combinationsInCond = cp.getCombinationsInCond();
0068 
0069   m_verbosity = cp.m_verbosity;
0070 }
0071 
0072 l1t::CorrCondition::CorrCondition(const l1t::CorrCondition& cp) : ConditionEvaluation() { copy(cp); }
0073 
0074 // destructor
0075 l1t::CorrCondition::~CorrCondition() {
0076   // empty
0077 }
0078 
0079 // equal operator
0080 l1t::CorrCondition& l1t::CorrCondition::operator=(const l1t::CorrCondition& cp) {
0081   copy(cp);
0082   return *this;
0083 }
0084 
0085 // methods
0086 void l1t::CorrCondition::setGtCorrelationTemplate(const CorrelationTemplate* caloTempl) {
0087   m_gtCorrelationTemplate = caloTempl;
0088 }
0089 
0090 ///   set the pointer to uGT GlobalBoard
0091 void l1t::CorrCondition::setuGtB(const GlobalBoard* ptrGTB) { m_uGtB = ptrGTB; }
0092 
0093 void l1t::CorrCondition::setScales(const GlobalScales* sc) { m_gtScales = sc; }
0094 
0095 // try all object permutations and check spatial correlations, if required
0096 const bool l1t::CorrCondition::evaluateCondition(const int bxEval) const {
0097   // std::cout << "m_isDebugEnabled = " << m_isDebugEnabled << std::endl;
0098   // std::cout << "m_verbosity = " << m_verbosity << std::endl;
0099 
0100   //std::ostringstream myCout;
0101   //m_gtCorrelationTemplate->print(myCout);
0102   //LogDebug("L1TGlobal")
0103   //   << "Correlation Condition Evaluation \n" << myCout.str() << std::endl;
0104 
0105   bool condResult = false;
0106   bool reqObjResult = false;
0107 
0108   // number of objects in condition (it is 2, no need to retrieve from
0109   // condition template) and their type
0110   int nObjInCond = 2;
0111   std::vector<GlobalObject> cndObjTypeVec(nObjInCond);
0112 
0113   // evaluate first the two sub-conditions (Type1s)
0114 
0115   const GtConditionCategory cond0Categ = m_gtCorrelationTemplate->cond0Category();
0116   const GtConditionCategory cond1Categ = m_gtCorrelationTemplate->cond1Category();
0117 
0118   //Decide if we have a mixed (muon + cal) condition
0119   bool convertCaloScales = false;
0120   if ((cond0Categ == CondMuon && (cond1Categ == CondCalo || cond1Categ == CondEnergySum)) ||
0121       (cond1Categ == CondMuon && (cond0Categ == CondCalo || cond0Categ == CondEnergySum)))
0122     convertCaloScales = true;
0123 
0124   const MuonTemplate* corrMuon = nullptr;
0125   const CaloTemplate* corrCalo = nullptr;
0126   const EnergySumTemplate* corrEnergySum = nullptr;
0127 
0128   // FIXME copying is slow...
0129   CombinationsInCond cond0Comb;
0130   CombinationsInCond cond1Comb;
0131 
0132   int cond0bx(0);
0133   int cond1bx(0);
0134 
0135   switch (cond0Categ) {
0136     case CondMuon: {
0137       corrMuon = static_cast<const MuonTemplate*>(m_gtCond0);
0138       MuCondition muCondition(
0139           corrMuon, m_uGtB, 0, 0);  //BLW these are counts that don't seem to be used...perhaps remove
0140 
0141       muCondition.evaluateConditionStoreResult(bxEval);
0142       reqObjResult = muCondition.condLastResult();
0143 
0144       cond0Comb = (muCondition.getCombinationsInCond());
0145       cond0bx = bxEval + (corrMuon->condRelativeBx());
0146       cndObjTypeVec[0] = (corrMuon->objectType())[0];
0147 
0148       if (m_verbosity) {
0149         std::ostringstream myCout;
0150         muCondition.print(myCout);
0151 
0152         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0153       }
0154     } break;
0155     case CondCalo: {
0156       corrCalo = static_cast<const CaloTemplate*>(m_gtCond0);
0157 
0158       CaloCondition caloCondition(
0159           corrCalo, m_uGtB, 0, 0, 0, 0);  //BLW these are counters that don't seem to be used...perhaps remove.
0160 
0161       caloCondition.evaluateConditionStoreResult(bxEval);
0162       reqObjResult = caloCondition.condLastResult();
0163 
0164       cond0Comb = (caloCondition.getCombinationsInCond());
0165       cond0bx = bxEval + (corrCalo->condRelativeBx());
0166       cndObjTypeVec[0] = (corrCalo->objectType())[0];
0167 
0168       if (m_verbosity) {
0169         std::ostringstream myCout;
0170         caloCondition.print(myCout);
0171 
0172         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0173       }
0174     } break;
0175     case CondEnergySum: {
0176       corrEnergySum = static_cast<const EnergySumTemplate*>(m_gtCond0);
0177       EnergySumCondition eSumCondition(corrEnergySum, m_uGtB);
0178 
0179       eSumCondition.evaluateConditionStoreResult(bxEval);
0180       reqObjResult = eSumCondition.condLastResult();
0181 
0182       cond0Comb = (eSumCondition.getCombinationsInCond());
0183       cond0bx = bxEval + (corrEnergySum->condRelativeBx());
0184       cndObjTypeVec[0] = (corrEnergySum->objectType())[0];
0185 
0186       if (m_verbosity) {
0187         std::ostringstream myCout;
0188         eSumCondition.print(myCout);
0189 
0190         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0191       }
0192     } break;
0193     default: {
0194       // should not arrive here, there are no correlation conditions defined for this object
0195       return false;
0196     } break;
0197   }
0198 
0199   // return if first subcondition is false
0200   if (!reqObjResult) {
0201     LogDebug("L1TGlobal") << "\n  First sub-condition false, second sub-condition not evaluated and not printed."
0202                           << std::endl;
0203     return false;
0204   }
0205 
0206   // second object
0207   reqObjResult = false;
0208 
0209   switch (cond1Categ) {
0210     case CondMuon: {
0211       corrMuon = static_cast<const MuonTemplate*>(m_gtCond1);
0212       MuCondition muCondition(
0213           corrMuon, m_uGtB, 0, 0);  //BLW these are counts that don't seem to be used...perhaps remove
0214 
0215       muCondition.evaluateConditionStoreResult(bxEval);
0216       reqObjResult = muCondition.condLastResult();
0217 
0218       cond1Comb = (muCondition.getCombinationsInCond());
0219       cond1bx = bxEval + (corrMuon->condRelativeBx());
0220 
0221       cndObjTypeVec[1] = (corrMuon->objectType())[0];
0222 
0223       if (m_verbosity) {
0224         std::ostringstream myCout;
0225         muCondition.print(myCout);
0226 
0227         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0228       }
0229     } break;
0230     case CondCalo: {
0231       corrCalo = static_cast<const CaloTemplate*>(m_gtCond1);
0232       CaloCondition caloCondition(
0233           corrCalo, m_uGtB, 0, 0, 0, 0);  //BLW these are counters that don't seem to be used...perhaps remove.
0234 
0235       caloCondition.evaluateConditionStoreResult(bxEval);
0236       reqObjResult = caloCondition.condLastResult();
0237 
0238       cond1Comb = (caloCondition.getCombinationsInCond());
0239       cond1bx = bxEval + (corrCalo->condRelativeBx());
0240       cndObjTypeVec[1] = (corrCalo->objectType())[0];
0241 
0242       if (m_verbosity) {
0243         std::ostringstream myCout;
0244         caloCondition.print(myCout);
0245 
0246         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0247       }
0248 
0249     } break;
0250     case CondEnergySum: {
0251       corrEnergySum = static_cast<const EnergySumTemplate*>(m_gtCond1);
0252 
0253       EnergySumCondition eSumCondition(corrEnergySum, m_uGtB);
0254 
0255       eSumCondition.evaluateConditionStoreResult(bxEval);
0256       reqObjResult = eSumCondition.condLastResult();
0257 
0258       cond1Comb = (eSumCondition.getCombinationsInCond());
0259       cond1bx = bxEval + (corrEnergySum->condRelativeBx());
0260       cndObjTypeVec[1] = (corrEnergySum->objectType())[0];
0261 
0262       if (m_verbosity) {
0263         std::ostringstream myCout;
0264         eSumCondition.print(myCout);
0265 
0266         LogDebug("L1TGlobal") << myCout.str() << std::endl;
0267       }
0268     } break;
0269     default: {
0270       // should not arrive here, there are no correlation conditions defined for this object
0271       return false;
0272     } break;
0273   }
0274 
0275   // return if second sub-condition is false
0276   if (!reqObjResult) {
0277     return false;
0278   } else {
0279     LogDebug("L1TGlobal") << "\n"
0280                           << "    Both sub-conditions true for object requirements."
0281                           << "    Evaluate correlation requirements.\n"
0282                           << std::endl;
0283   }
0284 
0285   // since we have two good legs get the correlation parameters
0286   CorrelationTemplate::CorrelationParameter corrPar = *(m_gtCorrelationTemplate->correlationParameter());
0287 
0288   // vector to store the indices of the calorimeter objects
0289   // from the combination evaluated in the condition
0290   SingleCombInCond objectsInComb;
0291   objectsInComb.reserve(nObjInCond);
0292 
0293   // clear the m_combinationsInCond vector
0294   (combinationsInCond()).clear();
0295 
0296   // pointers to objects
0297   const BXVector<const l1t::Muon*>* candMuVec = nullptr;
0298   const BXVector<const l1t::L1Candidate*>* candCaloVec = nullptr;
0299   const BXVector<const l1t::EtSum*>* candEtSumVec = nullptr;
0300 
0301   bool etSumCond = false;
0302 
0303   // make the conversions of the indices, depending on the combination of objects involved
0304   // (via pair index)
0305 
0306   int phiIndex0 = 0;
0307   double phi0Phy = 0.;
0308   int phiIndex1 = 0;
0309   double phi1Phy = 0.;
0310 
0311   int etaIndex0 = 0;
0312   double eta0Phy = 0.;
0313   int etaIndex1 = 0;
0314   double eta1Phy = 0.;
0315   int etaBin0 = 0;
0316   int etaBin1 = 0;
0317 
0318   int etIndex0 = 0;
0319   int etBin0 = 0;
0320   double et0Phy = 0.;
0321   int etIndex1 = 0;
0322   int etBin1 = 0;
0323   double et1Phy = 0.;
0324 
0325   int uptIndex0 = 0;     // Added displaced muons
0326   int uptBin0 = 0;       // Added displaced muons
0327   double upt0Phy = 0.0;  // Added displaced muons
0328   int uptIndex1 = 0;     // Added displaced muons
0329   int uptBin1 = 0;       // Added displaced muons
0330   double upt1Phy = 0.0;  // Added displaced muons
0331 
0332   int chrg0 = -1;
0333   int chrg1 = -1;
0334 
0335   // Determine the number of phi bins to get cutoff at pi
0336   int phiBound = 0;
0337   if (cond0Categ == CondMuon || cond1Categ == CondMuon) {
0338     const GlobalScales::ScaleParameters& par = m_gtScales->getMUScales();
0339     //phiBound = par.phiBins.size()/2;
0340     phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0341   } else {
0342     //Assumes all calorimeter objects are on same phi scale
0343     const GlobalScales::ScaleParameters& par = m_gtScales->getEGScales();
0344     //phiBound = par.phiBins.size()/2;
0345     phiBound = (int)((par.phiMax - par.phiMin) / par.phiStep) / 2;
0346   }
0347   LogDebug("L1TGlobal") << "Phi Bound = " << phiBound << std::endl;
0348 
0349   //          BUILD THE 1/dR2 LUT BASED ON CONDITION
0350   int iRSq = 0;
0351   int jRSq = 0;
0352   int iRSqmax = 227;
0353   int jRSqmax = 145;
0354   double tempRsq = 0.0;
0355   double tempdiffeta = 0.0435;
0356   double tempdiffphi = 0.02181661564992912;
0357   double resolution = 1.0;
0358   unsigned int precInvRsq = 0x5;
0359   if (cond0Categ == CondMuon || cond1Categ == CondMuon) {
0360     LogTrace("L1TGlobal") << "  Building 1/dR2 for Muon " << std::endl;
0361     tempdiffeta = tempdiffeta / 2.0;
0362     resolution = 0.5;
0363   } else {
0364     LogTrace("L1TGlobal") << "  Building 1/dR2 for Calo " << std::endl;
0365     iRSqmax = 231;
0366     jRSqmax = 73;
0367     tempdiffphi = 0.04363323129985824;
0368   }
0369   long long InvDeltaRSqLUT[iRSqmax][jRSqmax];
0370   double temp_InvDeltaRSq[iRSqmax][jRSqmax];
0371   if (corrPar.corrCutType & 0x80) {  //Only build the 1/dR2 LUT if necessary
0372     for (iRSq = 0; iRSq < iRSqmax; iRSq = iRSq + 1) {
0373       for (jRSq = 0; jRSq < jRSqmax; jRSq = jRSq + 1) {
0374         tempRsq = (tempdiffeta * iRSq) * (tempdiffeta * iRSq) + (tempdiffphi * jRSq) * (tempdiffphi * jRSq);
0375         if (tempRsq > 0.0) {
0376           temp_InvDeltaRSq[iRSq][jRSq] = 1.0 / tempRsq;
0377           InvDeltaRSqLUT[iRSq][jRSq] = (long long)round(pow(10, precInvRsq) * temp_InvDeltaRSq[iRSq][jRSq]);
0378         } else {
0379           temp_InvDeltaRSq[iRSq][jRSq] = 0.0;
0380           InvDeltaRSqLUT[iRSq][jRSq] = (long long)0;
0381         }
0382       }
0383     }
0384   }
0385 
0386   // Keep track of objects for LUTS
0387   std::string lutObj0 = "NULL";
0388   std::string lutObj1 = "NULL";
0389 
0390   LogTrace("L1TGlobal") << "  Sub-condition 0: std::vector<SingleCombInCond> size: " << (cond0Comb.size()) << std::endl;
0391   LogTrace("L1TGlobal") << "  Sub-condition 1: std::vector<SingleCombInCond> size: " << (cond1Comb.size()) << std::endl;
0392 
0393   // loop over all combinations which produced individually "true" as Type1s
0394   //
0395   // BLW: Optimization issue: potentially making the same comparison twice
0396   //                          if both legs are the same object type.
0397   for (std::vector<SingleCombInCond>::const_iterator it0Comb = cond0Comb.begin(); it0Comb != cond0Comb.end();
0398        it0Comb++) {
0399     // Type1s: there is 1 object only, no need for a loop, index 0 should be OK in (*it0Comb)[0]
0400     // ... but add protection to not crash
0401     int obj0Index = -1;
0402 
0403     if (!(*it0Comb).empty()) {
0404       obj0Index = (*it0Comb)[0];
0405     } else {
0406       LogTrace("L1TGlobal") << "\n  SingleCombInCond (*it0Comb).size() " << ((*it0Comb).size()) << std::endl;
0407       return false;
0408     }
0409 
0410     // Collect the information on the first leg of the correlation
0411     switch (cond0Categ) {
0412       case CondMuon: {
0413         lutObj0 = "MU";
0414         candMuVec = m_uGtB->getCandL1Mu();
0415         phiIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwPhiAtVtx();  //(*candMuVec)[obj0Index]->phiIndex();
0416         etaIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwEtaAtVtx();
0417         etIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwPt();
0418         uptIndex0 = (candMuVec->at(cond0bx, obj0Index))->hwPtUnconstrained();  // Added for displaced muons
0419         chrg0 = (candMuVec->at(cond0bx, obj0Index))->hwCharge();
0420         int etaBin0 = etaIndex0;
0421         if (etaBin0 < 0)
0422           etaBin0 = m_gtScales->getMUScales().etaBins.size() + etaBin0;  //twos complement
0423 
0424         etBin0 = etIndex0;
0425         uptBin0 = uptIndex0;  // Added for displaced muons
0426         int ssize = m_gtScales->getMUScales().etBins.size();
0427         if (etBin0 >= ssize) {
0428           etBin0 = ssize - 1;
0429           LogTrace("L1TGlobal") << "muon0 hw et" << etBin0 << " out of scale range.  Setting to maximum.";
0430         }
0431 
0432         // Determine Floating Pt numbers for floating point caluclation
0433         std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex0);
0434         phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0435         binEdges = m_gtScales->getMUScales().etaBins.at(etaBin0);
0436         eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0437         binEdges = m_gtScales->getMUScales().etBins.at(etBin0);
0438         et0Phy = 0.5 * (binEdges.second + binEdges.first);
0439         if (corrPar.corrCutType & 0x40)  // Added for displaced muons
0440         {
0441           binEdges = m_gtScales->getMUScales().uptBins.at(uptBin0);
0442           upt0Phy = 0.5 * (binEdges.second + binEdges.first);
0443         }
0444         LogDebug("L1TGlobal") << "Found all quantities for the muon 0" << std::endl;
0445       } break;
0446 
0447         // Calorimeter Objects (EG, Jet, Tau)
0448       case CondCalo: {
0449         switch (cndObjTypeVec[0]) {
0450           case gtEG: {
0451             lutObj0 = "EG";
0452             candCaloVec = m_uGtB->getCandL1EG();
0453             phiIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPhi();
0454             etaIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwEta();
0455             etIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPt();
0456             etaBin0 = etaIndex0;
0457             if (etaBin0 < 0)
0458               etaBin0 = m_gtScales->getEGScales().etaBins.size() + etaBin0;
0459             //          LogDebug("L1TGlobal") << "EG0 phi" << phiIndex0 << " eta " << etaIndex0 << " etaBin0 = " << etaBin0 << " et " << etIndex0 << std::endl;
0460 
0461             etBin0 = etIndex0;
0462             int ssize = m_gtScales->getEGScales().etBins.size();
0463             if (etBin0 >= ssize) {
0464               etBin0 = ssize - 1;
0465               LogTrace("L1TGlobal") << "EG0 hw et" << etBin0 << " out of scale range.  Setting to maximum.";
0466             }
0467 
0468             // Determine Floating Pt numbers for floating point caluclation
0469             std::pair<double, double> binEdges = m_gtScales->getEGScales().phiBins.at(phiIndex0);
0470             phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0471             binEdges = m_gtScales->getEGScales().etaBins.at(etaBin0);
0472             eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0473             binEdges = m_gtScales->getEGScales().etBins[etBin0];
0474             et0Phy = 0.5 * (binEdges.second + binEdges.first);
0475 
0476           } break;
0477           case gtJet: {
0478             lutObj0 = "JET";
0479             candCaloVec = m_uGtB->getCandL1Jet();
0480             phiIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPhi();
0481             etaIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwEta();
0482             etIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPt();
0483             etaBin0 = etaIndex0;
0484             if (etaBin0 < 0)
0485               etaBin0 = m_gtScales->getJETScales().etaBins.size() + etaBin0;
0486 
0487             etBin0 = etIndex0;
0488             int ssize = m_gtScales->getJETScales().etBins.size();
0489             if (etBin0 >= ssize) {
0490               //edm::LogWarning("L1TGlobal")
0491               //<< "jet0 hw et" << etBin0 << " out of scale range.  Setting to maximum.";
0492               etBin0 = ssize - 1;
0493             }
0494 
0495             // Determine Floating Pt numbers for floating point caluclation
0496             std::pair<double, double> binEdges = m_gtScales->getJETScales().phiBins.at(phiIndex0);
0497             phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0498             binEdges = m_gtScales->getJETScales().etaBins.at(etaBin0);
0499             eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0500             binEdges = m_gtScales->getJETScales().etBins.at(etBin0);
0501             et0Phy = 0.5 * (binEdges.second + binEdges.first);
0502 
0503           } break;
0504           case gtTau: {
0505             candCaloVec = m_uGtB->getCandL1Tau();
0506             phiIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPhi();
0507             etaIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwEta();
0508             etIndex0 = (candCaloVec->at(cond0bx, obj0Index))->hwPt();
0509             etaBin0 = etaIndex0;
0510             if (etaBin0 < 0)
0511               etaBin0 = m_gtScales->getTAUScales().etaBins.size() + etaBin0;
0512 
0513             etBin0 = etIndex0;
0514             int ssize = m_gtScales->getTAUScales().etBins.size();
0515             if (etBin0 >= ssize) {
0516               etBin0 = ssize - 1;
0517               LogTrace("L1TGlobal") << "tau0 hw et" << etBin0 << " out of scale range.  Setting to maximum.";
0518             }
0519 
0520             // Determine Floating Pt numbers for floating point caluclation
0521             std::pair<double, double> binEdges = m_gtScales->getTAUScales().phiBins.at(phiIndex0);
0522             phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0523             binEdges = m_gtScales->getTAUScales().etaBins.at(etaBin0);
0524             eta0Phy = 0.5 * (binEdges.second + binEdges.first);
0525             binEdges = m_gtScales->getTAUScales().etBins.at(etBin0);
0526             et0Phy = 0.5 * (binEdges.second + binEdges.first);
0527             lutObj0 = "TAU";
0528           } break;
0529           default: {
0530           } break;
0531         }  //end switch on calo type.
0532 
0533         //If needed convert calo scales to muon scales for comparison
0534         if (convertCaloScales) {
0535           std::string lutName = lutObj0;
0536           lutName += "-MU";
0537           long long tst = m_gtScales->getLUT_CalMuEta(lutName, etaBin0);
0538           LogDebug("L1TGlobal") << lutName << "  EtaCal = " << etaIndex0 << " etaBin0 = " << etaBin0
0539                                 << " EtaMu = " << tst << std::endl;
0540           etaIndex0 = tst;
0541 
0542           tst = m_gtScales->getLUT_CalMuPhi(lutName, phiIndex0);
0543           LogDebug("L1TGlobal") << lutName << "  PhiCal = " << phiIndex0 << " PhiMu = " << tst << std::endl;
0544           phiIndex0 = tst;
0545         }
0546 
0547       } break;
0548 
0549         // Energy Sums
0550       case CondEnergySum: {
0551         etSumCond = true;
0552         //Stupid mapping between enum types for energy sums.
0553         l1t::EtSum::EtSumType type;
0554         switch (cndObjTypeVec[0]) {
0555           case gtETM:
0556             type = l1t::EtSum::EtSumType::kMissingEt;
0557             lutObj0 = "ETM";
0558             break;
0559           case gtETT:
0560             type = l1t::EtSum::EtSumType::kTotalEt;
0561             lutObj0 = "ETT";
0562             break;
0563           case gtETTem:
0564             type = l1t::EtSum::EtSumType::kTotalEtEm;
0565             lutObj0 =
0566                 "ETTem";  //should this be just ETT (share LUTs?) Can't be used for CorrCond anyway since now directional information
0567             break;
0568           case gtHTM:
0569             type = l1t::EtSum::EtSumType::kMissingHt;
0570             lutObj0 = "HTM";
0571             break;
0572           case gtHTT:
0573             type = l1t::EtSum::EtSumType::kTotalHt;
0574             lutObj0 = "HTT";
0575             break;
0576           case gtETMHF:
0577             type = l1t::EtSum::EtSumType::kMissingEtHF;
0578             lutObj0 = "ETMHF";
0579             break;
0580           case gtMinBiasHFP0:
0581           case gtMinBiasHFM0:
0582           case gtMinBiasHFP1:
0583           case gtMinBiasHFM1:
0584             type = l1t::EtSum::EtSumType::kMinBiasHFP0;
0585             lutObj0 =
0586                 "MinBias";  //??Fix?? Not a valid LUT type Can't be used for CorrCond anyway since now directional information
0587             break;
0588           default:
0589             edm::LogError("L1TGlobal") << "\n  Error: "
0590                                        << "Unmatched object type from template to EtSumType, cndObjTypeVec[0] = "
0591                                        << cndObjTypeVec[0] << std::endl;
0592             type = l1t::EtSum::EtSumType::kTotalEt;
0593             break;
0594         }
0595 
0596         candEtSumVec = m_uGtB->getCandL1EtSum();
0597 
0598         for (int iEtSum = 0; iEtSum < (int)candEtSumVec->size(cond0bx); iEtSum++) {
0599           if ((candEtSumVec->at(cond0bx, iEtSum))->getType() == type) {
0600             phiIndex0 = (candEtSumVec->at(cond0bx, iEtSum))->hwPhi();
0601             etaIndex0 = (candEtSumVec->at(cond0bx, iEtSum))->hwEta();
0602             etIndex0 = (candEtSumVec->at(cond0bx, iEtSum))->hwPt();
0603 
0604             //  Get the floating point numbers
0605             if (cndObjTypeVec[0] == gtETM) {
0606               std::pair<double, double> binEdges = m_gtScales->getETMScales().phiBins.at(phiIndex0);
0607               phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0608               eta0Phy = 0.;  //No Eta for Energy Sums
0609 
0610               etBin0 = etIndex0;
0611               int ssize = m_gtScales->getETMScales().etBins.size();
0612               assert(ssize > 0);
0613               if (etBin0 >= ssize) {
0614                 etBin0 = ssize - 1;
0615               }
0616 
0617               binEdges = m_gtScales->getETMScales().etBins.at(etBin0);
0618               et0Phy = 0.5 * (binEdges.second + binEdges.first);
0619             } else if (cndObjTypeVec[0] == gtHTM) {
0620               std::pair<double, double> binEdges = m_gtScales->getHTMScales().phiBins.at(phiIndex0);
0621               phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0622               eta0Phy = 0.;  //No Eta for Energy Sums
0623 
0624               etBin0 = etIndex0;
0625               int ssize = m_gtScales->getHTMScales().etBins.size();
0626               assert(ssize > 0);
0627               if (etBin0 >= ssize) {
0628                 etBin0 = ssize - 1;
0629               }
0630 
0631               binEdges = m_gtScales->getHTMScales().etBins.at(etBin0);
0632               et0Phy = 0.5 * (binEdges.second + binEdges.first);
0633             } else if (cndObjTypeVec[0] == gtETMHF) {
0634               std::pair<double, double> binEdges = m_gtScales->getETMHFScales().phiBins.at(phiIndex0);
0635               phi0Phy = 0.5 * (binEdges.second + binEdges.first);
0636               eta0Phy = 0.;  //No Eta for Energy Sums
0637 
0638               etBin0 = etIndex0;
0639               int ssize = m_gtScales->getETMHFScales().etBins.size();
0640               assert(ssize > 0);
0641               if (etBin0 >= ssize) {
0642                 etBin0 = ssize - 1;
0643               }
0644 
0645               binEdges = m_gtScales->getETMHFScales().etBins.at(etBin0);
0646               et0Phy = 0.5 * (binEdges.second + binEdges.first);
0647             }
0648 
0649             //If needed convert calo scales to muon scales for comparison (only phi for energy sums)
0650             if (convertCaloScales) {
0651               std::string lutName = lutObj0;
0652               lutName += "-MU";
0653               long long tst = m_gtScales->getLUT_CalMuPhi(lutName, phiIndex0);
0654               LogDebug("L1TGlobal") << lutName << "  PhiCal = " << phiIndex0 << " PhiMu = " << tst << std::endl;
0655               phiIndex0 = tst;
0656             }
0657 
0658           }  //check it is the EtSum we want
0659         }    // loop over Etsums
0660 
0661       } break;
0662 
0663       default: {
0664         // should not arrive here, there are no correlation conditions defined for this object
0665         LogDebug("L1TGlobal") << "Error could not find the Cond Category for Leg 0" << std::endl;
0666         return false;
0667       } break;
0668     }  //end switch on first leg type
0669 
0670     // Now loop over the second leg to get its information
0671     for (std::vector<SingleCombInCond>::const_iterator it1Comb = cond1Comb.begin(); it1Comb != cond1Comb.end();
0672          it1Comb++) {
0673       LogDebug("L1TGlobal") << "Looking at second Condition" << std::endl;
0674       // Type1s: there is 1 object only, no need for a loop (*it1Comb)[0]
0675       // ... but add protection to not crash
0676       int obj1Index = -1;
0677 
0678       if (!(*it1Comb).empty()) {
0679         obj1Index = (*it1Comb)[0];
0680       } else {
0681         LogTrace("L1TGlobal") << "\n  SingleCombInCond (*it1Comb).size() " << ((*it1Comb).size()) << std::endl;
0682         return false;
0683       }
0684 
0685       //If we are dealing with the same object type avoid the two legs
0686       // either being the same object
0687       if (cndObjTypeVec[0] == cndObjTypeVec[1] && obj0Index == obj1Index && cond0bx == cond1bx) {
0688         LogDebug("L1TGlobal") << "Corr Condition looking at same leg...skip" << std::endl;
0689         continue;
0690       }
0691 
0692       switch (cond1Categ) {
0693         case CondMuon: {
0694           lutObj1 = "MU";
0695           candMuVec = m_uGtB->getCandL1Mu();
0696           phiIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwPhiAtVtx();  //(*candMuVec)[obj0Index]->phiIndex();
0697           etaIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwEtaAtVtx();
0698           etIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwPt();
0699           uptIndex1 = (candMuVec->at(cond1bx, obj1Index))->hwPtUnconstrained();  // Added for displaced muons
0700           chrg1 = (candMuVec->at(cond1bx, obj1Index))->hwCharge();
0701           etaBin1 = etaIndex1;
0702           if (etaBin1 < 0)
0703             etaBin1 = m_gtScales->getMUScales().etaBins.size() + etaBin1;
0704           //           LogDebug("L1TGlobal") << "Muon phi" << phiIndex1 << " eta " << etaIndex1 << " etaBin1 = " << etaBin1  << " et " << etIndex1 << std::endl;
0705 
0706           etBin1 = etIndex1;
0707           uptBin1 = uptIndex1;  // Added for displaced muons
0708           int ssize = m_gtScales->getMUScales().etBins.size();
0709           if (etBin1 >= ssize) {
0710             LogTrace("L1TGlobal") << "muon2 hw et" << etBin1 << " out of scale range.  Setting to maximum.";
0711             etBin1 = ssize - 1;
0712           }
0713 
0714           // Determine Floating Pt numbers for floating point caluclation
0715           std::pair<double, double> binEdges = m_gtScales->getMUScales().phiBins.at(phiIndex1);
0716           phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0717           binEdges = m_gtScales->getMUScales().etaBins.at(etaBin1);
0718           eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0719           binEdges = m_gtScales->getMUScales().etBins.at(etBin1);
0720           et1Phy = 0.5 * (binEdges.second + binEdges.first);
0721           if (corrPar.corrCutType & 0x40)  // Added for displaced muons
0722           {
0723             binEdges = m_gtScales->getMUScales().uptBins.at(uptBin1);
0724             upt1Phy = 0.5 * (binEdges.second + binEdges.first);
0725           }
0726           LogDebug("L1TGlobal") << "Found all quantities for the muon 1" << std::endl;
0727         } break;
0728         case CondCalo: {
0729           switch (cndObjTypeVec[1]) {
0730             case gtEG: {
0731               candCaloVec = m_uGtB->getCandL1EG();
0732               phiIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPhi();
0733               etaIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwEta();
0734               etIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPt();
0735               etaBin1 = etaIndex1;
0736               if (etaBin1 < 0)
0737                 etaBin1 = m_gtScales->getEGScales().etaBins.size() + etaBin1;
0738 
0739               etBin1 = etIndex1;
0740               int ssize = m_gtScales->getEGScales().etBins.size();
0741               if (etBin1 >= ssize) {
0742                 LogTrace("L1TGlobal") << "EG1 hw et" << etBin1 << " out of scale range.  Setting to maximum.";
0743                 etBin1 = ssize - 1;
0744               }
0745 
0746               // Determine Floating Pt numbers for floating point caluclation
0747               std::pair<double, double> binEdges = m_gtScales->getEGScales().phiBins.at(phiIndex1);
0748               phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0749               binEdges = m_gtScales->getEGScales().etaBins.at(etaBin1);
0750               eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0751               binEdges = m_gtScales->getEGScales().etBins.at(etBin1);
0752               et1Phy = 0.5 * (binEdges.second + binEdges.first);
0753               lutObj1 = "EG";
0754             } break;
0755             case gtJet: {
0756               candCaloVec = m_uGtB->getCandL1Jet();
0757               phiIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPhi();
0758               etaIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwEta();
0759               etIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPt();
0760               etaBin1 = etaIndex1;
0761               if (etaBin1 < 0)
0762                 etaBin1 = m_gtScales->getJETScales().etaBins.size() + etaBin1;
0763 
0764               etBin1 = etIndex1;
0765               int ssize = m_gtScales->getJETScales().etBins.size();
0766               assert(ssize);
0767               if (etBin1 >= ssize) {
0768                 //edm::LogWarning("L1TGlobal")
0769                 //<< "jet2 hw et" << etBin1 << " out of scale range.  Setting to maximum.";
0770                 etBin1 = ssize - 1;
0771               }
0772 
0773               // Determine Floating Pt numbers for floating point caluclation
0774               std::pair<double, double> binEdges = m_gtScales->getJETScales().phiBins.at(phiIndex1);
0775               phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0776               binEdges = m_gtScales->getJETScales().etaBins.at(etaBin1);
0777               eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0778               //CRASHES HERE:
0779               binEdges = m_gtScales->getJETScales().etBins.at(etBin1);
0780               et1Phy = 0.5 * (binEdges.second + binEdges.first);
0781               lutObj1 = "JET";
0782             } break;
0783             case gtTau: {
0784               candCaloVec = m_uGtB->getCandL1Tau();
0785               phiIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPhi();
0786               etaIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwEta();
0787               etIndex1 = (candCaloVec->at(cond1bx, obj1Index))->hwPt();
0788               etaBin1 = etaIndex1;
0789               if (etaBin1 < 0)
0790                 etaBin1 = m_gtScales->getTAUScales().etaBins.size() + etaBin1;
0791 
0792               etBin1 = etIndex1;
0793               int ssize = m_gtScales->getTAUScales().etBins.size();
0794               if (etBin1 >= ssize) {
0795                 LogTrace("L1TGlobal") << "tau2 hw et" << etBin1 << " out of scale range.  Setting to maximum.";
0796                 etBin1 = ssize - 1;
0797               }
0798 
0799               // Determine Floating Pt numbers for floating point caluclation
0800               std::pair<double, double> binEdges = m_gtScales->getTAUScales().phiBins.at(phiIndex1);
0801               phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0802               binEdges = m_gtScales->getTAUScales().etaBins.at(etaBin1);
0803               eta1Phy = 0.5 * (binEdges.second + binEdges.first);
0804               binEdges = m_gtScales->getTAUScales().etBins.at(etBin1);
0805               et1Phy = 0.5 * (binEdges.second + binEdges.first);
0806               lutObj1 = "TAU";
0807             } break;
0808             default: {
0809             } break;
0810           }  //end switch on calo type.
0811 
0812           //If needed convert calo scales to muon scales for comparison
0813           if (convertCaloScales) {
0814             std::string lutName = lutObj1;
0815             lutName += "-MU";
0816             long long tst = m_gtScales->getLUT_CalMuEta(lutName, etaBin1);
0817             LogDebug("L1TGlobal") << lutName << "  EtaCal = " << etaIndex1 << " EtaMu = " << tst << std::endl;
0818             etaIndex1 = tst;
0819 
0820             tst = m_gtScales->getLUT_CalMuPhi(lutName, phiIndex1);
0821             LogDebug("L1TGlobal") << lutName << "  PhiCal = " << phiIndex1 << " PhiMu = " << tst << std::endl;
0822             phiIndex1 = tst;
0823           }
0824 
0825         } break;
0826         case CondEnergySum: {
0827           LogDebug("L1TGlobal") << "Looking at second Condition as Energy Sum: " << cndObjTypeVec[1] << std::endl;
0828           etSumCond = true;
0829 
0830           //Stupid mapping between enum types for energy sums.
0831           l1t::EtSum::EtSumType type;
0832           switch (cndObjTypeVec[1]) {
0833             case gtETM:
0834               type = l1t::EtSum::EtSumType::kMissingEt;
0835               lutObj1 = "ETM";
0836               break;
0837             case gtETT:
0838               type = l1t::EtSum::EtSumType::kTotalEt;
0839               lutObj1 = "ETT";
0840               break;
0841             case gtETTem:
0842               type = l1t::EtSum::EtSumType::kTotalEtEm;
0843               lutObj1 = "ETTem";
0844               break;
0845             case gtHTM:
0846               type = l1t::EtSum::EtSumType::kMissingHt;
0847               lutObj1 = "HTM";
0848               break;
0849             case gtHTT:
0850               type = l1t::EtSum::EtSumType::kTotalHt;
0851               lutObj1 = "HTT";
0852               break;
0853             case gtETMHF:
0854               type = l1t::EtSum::EtSumType::kMissingEtHF;
0855               lutObj1 = "ETMHF";
0856               break;
0857             case gtMinBiasHFP0:
0858             case gtMinBiasHFM0:
0859             case gtMinBiasHFP1:
0860             case gtMinBiasHFM1:
0861               type = l1t::EtSum::EtSumType::kMinBiasHFP0;
0862               lutObj1 = "MinBias";
0863               break;
0864             default:
0865               edm::LogError("L1TGlobal") << "\n  Error: "
0866                                          << "Unmatched object type from template to EtSumType, cndObjTypeVec[1] = "
0867                                          << cndObjTypeVec[1] << std::endl;
0868               type = l1t::EtSum::EtSumType::kTotalEt;
0869               break;
0870           }
0871 
0872           candEtSumVec = m_uGtB->getCandL1EtSum();
0873 
0874           LogDebug("L1TGlobal") << "obj " << lutObj1 << " Vector Size " << candEtSumVec->size(cond1bx) << std::endl;
0875           for (int iEtSum = 0; iEtSum < (int)candEtSumVec->size(cond1bx); iEtSum++) {
0876             if ((candEtSumVec->at(cond1bx, iEtSum))->getType() == type) {
0877               phiIndex1 = (candEtSumVec->at(cond1bx, iEtSum))->hwPhi();
0878               etaIndex1 = (candEtSumVec->at(cond1bx, iEtSum))->hwEta();
0879               etIndex1 = (candEtSumVec->at(cond1bx, iEtSum))->hwPt();
0880 
0881               // Determine Floating Pt numbers for floating point caluclation
0882 
0883               if (cndObjTypeVec[1] == gtETM) {
0884                 std::pair<double, double> binEdges = m_gtScales->getETMScales().phiBins.at(phiIndex1);
0885                 phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0886                 eta1Phy = 0.;  //No Eta for Energy Sums
0887 
0888                 etBin1 = etIndex1;
0889                 int ssize = m_gtScales->getETMScales().etBins.size();
0890                 assert(ssize > 0);
0891                 if (etBin1 >= ssize) {
0892                   etBin1 = ssize - 1;
0893                 }
0894 
0895                 binEdges = m_gtScales->getETMScales().etBins.at(etBin1);
0896                 et1Phy = 0.5 * (binEdges.second + binEdges.first);
0897               } else if (cndObjTypeVec[1] == gtHTM) {
0898                 std::pair<double, double> binEdges = m_gtScales->getHTMScales().phiBins.at(phiIndex1);
0899                 phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0900                 eta1Phy = 0.;  //No Eta for Energy Sums
0901 
0902                 etBin1 = etIndex1;
0903                 int ssize = m_gtScales->getHTMScales().etBins.size();
0904                 assert(ssize > 0);
0905                 if (etBin1 >= ssize) {
0906                   etBin1 = ssize - 1;
0907                 }
0908 
0909                 binEdges = m_gtScales->getHTMScales().etBins.at(etBin1);
0910                 et1Phy = 0.5 * (binEdges.second + binEdges.first);
0911               } else if (cndObjTypeVec[1] == gtETMHF) {
0912                 std::pair<double, double> binEdges = m_gtScales->getETMHFScales().phiBins.at(phiIndex1);
0913                 phi1Phy = 0.5 * (binEdges.second + binEdges.first);
0914                 eta1Phy = 0.;  //No Eta for Energy Sums
0915 
0916                 etBin1 = etIndex1;
0917                 int ssize = m_gtScales->getETMHFScales().etBins.size();
0918                 assert(ssize > 0);
0919                 if (etBin1 >= ssize) {
0920                   etBin1 = ssize - 1;
0921                 }
0922 
0923                 binEdges = m_gtScales->getETMHFScales().etBins.at(etBin1);
0924                 et1Phy = 0.5 * (binEdges.second + binEdges.first);
0925               }
0926 
0927               //If needed convert calo scales to muon scales for comparison (only phi for energy sums)
0928               if (convertCaloScales) {
0929                 std::string lutName = lutObj1;
0930                 lutName += "-MU";
0931                 long long tst = m_gtScales->getLUT_CalMuPhi(lutName, phiIndex1);
0932                 LogDebug("L1TGlobal") << lutName << "  PhiCal = " << phiIndex1 << " PhiMu = " << tst << std::endl;
0933                 phiIndex1 = tst;
0934               }
0935 
0936             }  //check it is the EtSum we want
0937           }    // loop over Etsums
0938 
0939         } break;
0940         default: {
0941           // should not arrive here, there are no correlation conditions defined for this object
0942           LogDebug("L1TGlobal") << "Error could not find the Cond Category for Leg 0" << std::endl;
0943           return false;
0944         } break;
0945       }  //end switch on second leg
0946 
0947       if (m_verbosity) {
0948         LogDebug("L1TGlobal") << "    Correlation pair [" << l1t::GlobalObjectEnumToString(cndObjTypeVec[0]) << ", "
0949                               << l1t::GlobalObjectEnumToString(cndObjTypeVec[1]) << "] with collection indices ["
0950                               << obj0Index << ", " << obj1Index << "] "
0951                               << " has: \n"
0952                               << "     Et  value   = [" << etIndex0 << ", " << etIndex1 << "]\n"
0953                               << "     phi indices = [" << phiIndex0 << ", " << phiIndex1 << "]\n"
0954                               << "     eta indices = [" << etaIndex0 << ", " << etaIndex1 << "]\n"
0955                               << "     chrg        = [" << chrg0 << ", " << chrg1 << "]\n";
0956       }
0957 
0958       // Now perform the desired correlation on these two objects. Assume true until we find a contradition
0959       bool reqResult = true;
0960 
0961       // clear the indices in the combination
0962       objectsInComb.clear();
0963 
0964       objectsInComb.push_back(obj0Index);
0965       objectsInComb.push_back(obj1Index);
0966 
0967       // if we get here all checks were successful for this combination
0968       // set the general result for evaluateCondition to "true"
0969 
0970       // These all require some delta eta and phi calculations.  Do them first...for now real calculation but need to
0971       // revise this to line up with firmware calculations.
0972       double deltaPhiPhy = fabs(phi1Phy - phi0Phy);
0973       if (deltaPhiPhy > M_PI)
0974         deltaPhiPhy = 2. * M_PI - deltaPhiPhy;
0975       double deltaEtaPhy = fabs(eta1Phy - eta0Phy);
0976 
0977       // Determine the integer based delta eta and delta phi
0978       int deltaPhiFW = abs(phiIndex0 - phiIndex1);
0979       if (deltaPhiFW >= phiBound)
0980         deltaPhiFW = 2 * phiBound - deltaPhiFW;
0981       std::string lutName = lutObj0;
0982       lutName += "-";
0983       lutName += lutObj1;
0984       long long deltaPhiLUT = m_gtScales->getLUT_DeltaPhi(lutName, deltaPhiFW);
0985       unsigned int precDeltaPhiLUT = m_gtScales->getPrec_DeltaPhi(lutName);
0986 
0987       int deltaEtaFW = abs(etaIndex0 - etaIndex1);
0988       long long deltaEtaLUT = 0;
0989       unsigned int precDeltaEtaLUT = 0;
0990       if (!etSumCond) {
0991         deltaEtaLUT = m_gtScales->getLUT_DeltaEta(lutName, deltaEtaFW);
0992         precDeltaEtaLUT = m_gtScales->getPrec_DeltaEta(lutName);
0993       }
0994 
0995       //
0996       LogDebug("L1TGlobal") << "Obj0 phiFW = " << phiIndex0 << " Obj1 phiFW = " << phiIndex1 << "\n"
0997                             << "    DeltaPhiFW = " << deltaPhiFW << "\n"
0998                             << "    LUT Name = " << lutName << " Prec = " << precDeltaPhiLUT
0999                             << "  DeltaPhiLUT = " << deltaPhiLUT << "\n"
1000                             << "Obj0 etaFW = " << etaIndex0 << " Obj1 etaFW = " << etaIndex1 << "\n"
1001                             << "    DeltaEtaFW = " << deltaEtaFW << "\n"
1002                             << "    LUT Name = " << lutName << " Prec = " << precDeltaEtaLUT
1003                             << "  DeltaEtaLUT = " << deltaEtaLUT << std::endl;
1004 
1005       // If there is a delta eta, check it.
1006       if (corrPar.corrCutType & 0x1) {
1007         unsigned int preShift = precDeltaEtaLUT - corrPar.precEtaCut;
1008         LogDebug("L1TGlobal") << "    Testing Delta Eta Cut (" << lutObj0 << "," << lutObj1 << ") ["
1009                               << (long long)(corrPar.minEtaCutValue * pow(10, preShift)) << ","
1010                               << (long long)(corrPar.maxEtaCutValue * pow(10, preShift))
1011                               << "] with precision = " << corrPar.precEtaCut << "\n"
1012                               << "    deltaEtaLUT = " << deltaEtaLUT << "\n"
1013                               << "    Precision Shift = " << preShift << "\n"
1014                               << "    deltaEta (shift)= " << (deltaEtaLUT / pow(10, preShift + corrPar.precEtaCut))
1015                               << "\n"
1016                               << "    deltaEtaPhy = " << deltaEtaPhy << std::endl;
1017 
1018         //if(preShift>0) deltaEtaLUT /= pow(10,preShift);
1019         if (deltaEtaLUT >= (long long)(corrPar.minEtaCutValue * pow(10, preShift)) &&
1020             deltaEtaLUT <= (long long)(corrPar.maxEtaCutValue * pow(10, preShift))) {
1021           LogDebug("L1TGlobal") << "    Passed Delta Eta Cut ["
1022                                 << (long long)(corrPar.minEtaCutValue * pow(10, preShift)) << ","
1023                                 << (long long)(corrPar.maxEtaCutValue * pow(10, preShift)) << "]" << std::endl;
1024 
1025         } else {
1026           LogDebug("L1TGlobal") << "    Failed Delta Eta Cut ["
1027                                 << (long long)(corrPar.minEtaCutValue * pow(10, preShift)) << ","
1028                                 << (long long)(corrPar.maxEtaCutValue * pow(10, preShift)) << "]" << std::endl;
1029           reqResult = false;
1030         }
1031       }
1032 
1033       //if there is a delta phi check it.
1034       if (corrPar.corrCutType & 0x2) {
1035         unsigned int preShift = precDeltaPhiLUT - corrPar.precPhiCut;
1036         LogDebug("L1TGlobal") << "    Testing Delta Phi Cut (" << lutObj0 << "," << lutObj1 << ") ["
1037                               << (long long)(corrPar.minPhiCutValue * pow(10, preShift)) << ","
1038                               << (long long)(corrPar.maxPhiCutValue * pow(10, preShift))
1039                               << "] with precision = " << corrPar.precPhiCut << "\n"
1040                               << "    deltaPhiLUT = " << deltaPhiLUT << "\n"
1041                               << "    Precision Shift = " << preShift << "\n"
1042                               << "    deltaPhi (shift)= " << (deltaPhiLUT / pow(10, preShift + corrPar.precPhiCut))
1043                               << "\n"
1044                               << "    deltaPhiPhy = " << deltaPhiPhy << std::endl;
1045 
1046         //if(preShift>0) deltaPhiLUT /= pow(10,preShift);
1047         if (deltaPhiLUT >= (long long)(corrPar.minPhiCutValue * pow(10, preShift)) &&
1048             deltaPhiLUT <= (long long)(corrPar.maxPhiCutValue * pow(10, preShift))) {
1049           LogDebug("L1TGlobal") << "    Passed Delta Phi Cut ["
1050                                 << (long long)(corrPar.minPhiCutValue * pow(10, preShift)) << ","
1051                                 << (long long)(corrPar.maxPhiCutValue * pow(10, preShift)) << "]" << std::endl;
1052 
1053         } else {
1054           LogDebug("L1TGlobal") << "    Failed Delta Phi Cut ["
1055                                 << (long long)(corrPar.minPhiCutValue * pow(10, preShift)) << ","
1056                                 << (long long)(corrPar.maxPhiCutValue * pow(10, preShift)) << "]" << std::endl;
1057           reqResult = false;
1058         }
1059       }
1060 
1061       if (corrPar.corrCutType & 0x4) {
1062         //Assumes Delta Eta and Delta Phi LUTs have the same precision
1063         unsigned int preShift = 2 * precDeltaPhiLUT - corrPar.precDRCut;
1064         double deltaRSqPhy = deltaPhiPhy * deltaPhiPhy + deltaEtaPhy * deltaEtaPhy;
1065         long long deltaRSq = deltaEtaLUT * deltaEtaLUT + deltaPhiLUT * deltaPhiLUT;
1066 
1067         LogDebug("L1TGlobal") << "    Testing Delta R Cut (" << lutObj0 << "," << lutObj1 << ") ["
1068                               << (long long)(corrPar.minDRCutValue * pow(10, preShift)) << ","
1069                               << (long long)(corrPar.maxDRCutValue * pow(10, preShift))
1070                               << "] with precision = " << corrPar.precDRCut << "\n"
1071                               << "    deltaPhiLUT = " << deltaPhiLUT << "\n"
1072                               << "    deltaEtaLUT = " << deltaEtaLUT << "\n"
1073                               << "    deltaRSqLUT = " << deltaRSq << "\n"
1074                               << "    Precision Shift = " << preShift << "\n"
1075                               << "    deltaRSqLUT (shift)= " << (deltaRSq / pow(10, preShift + corrPar.precDRCut))
1076                               << "\n"
1077                               << "    deltaRSqPhy = " << deltaRSqPhy << std::endl;
1078 
1079         //if(preShift>0) deltaRSq /= pow(10,preShift);
1080         if (deltaRSq >= (long long)(corrPar.minDRCutValue * pow(10, preShift)) &&
1081             deltaRSq <= (long long)(corrPar.maxDRCutValue * pow(10, preShift))) {
1082           LogDebug("L1TGlobal") << "    Passed Delta R Cut [" << (long long)(corrPar.minDRCutValue * pow(10, preShift))
1083                                 << "," << (long long)(corrPar.maxDRCutValue * pow(10, preShift)) << "]" << deltaRSq
1084                                 << std::endl;
1085 
1086         } else {
1087           LogDebug("L1TGlobal") << "    Failed Delta R Cut [" << (int)(corrPar.minDRCutValue * pow(10, preShift)) << ","
1088                                 << (long long)(corrPar.maxDRCutValue * pow(10, preShift)) << "]" << deltaRSq
1089                                 << std::endl;
1090           reqResult = false;
1091         }
1092       }
1093 
1094       if (corrPar.corrCutType & 0x20) {
1095         // Two body pt: pt^2 = pt1^2+pt2^2+2*pt1*pt2*(cos(phi1)*cos(phi2)+sin(phi1)*sin(phi2)).
1096 
1097         LogDebug("L1TGlobal") << " corrPar.corrCutType: " << corrPar.corrCutType << "\n";
1098 
1099         //calculate math sins and cosines for debugging
1100         double cosPhi1Phy = cos(phi0Phy);
1101         double sinPhi1Phy = sin(phi0Phy);
1102         double cosPhi2Phy = cos(phi1Phy);
1103         double sinPhi2Phy = sin(phi1Phy);
1104 
1105         double tbptSqPhy = et0Phy * et0Phy + et1Phy * et1Phy +
1106                            2 * et0Phy * et1Phy * (cosPhi1Phy * cosPhi2Phy + sinPhi1Phy * sinPhi2Phy);
1107         // get values from LUT's
1108 
1109         const std::string& lutName0 = lutObj0;
1110         unsigned int precCosLUT0 = m_gtScales->getPrec_Cos(lutName0);
1111         unsigned int precSinLUT0 = m_gtScales->getPrec_Sin(lutName0);
1112 
1113         const std::string& lutName1 = lutObj1;
1114         unsigned int precCosLUT1 = m_gtScales->getPrec_Cos(lutName1);
1115         unsigned int precSinLUT1 = m_gtScales->getPrec_Sin(lutName1);
1116 
1117         if (precCosLUT0 - precCosLUT1 != 0)
1118           LogDebug("L1TGlobal") << "Warning: Cos LUTs for TwoBodyPt on different Precision" << std::endl;
1119         if (precSinLUT0 - precSinLUT1 != 0)
1120           LogDebug("L1TGlobal") << "Warning: Sin LUTs for TwoBodyPt on different Precision" << std::endl;
1121         if (precSinLUT0 - precCosLUT1 != 0)
1122           LogDebug("L1TGlobal") << "Warning: Sin and Cos LUTs for TwoBodyPt on different Precision" << std::endl;
1123         if (precSinLUT1 - precCosLUT0 != 0)
1124           LogDebug("L1TGlobal") << "Warning: Sin and Cos LUTs for TwoBodyPt on different Precision" << std::endl;
1125 
1126         long long cosPhi1LUT = m_gtScales->getLUT_Cos(lutName0, phiIndex0);
1127         long long sinPhi1LUT = m_gtScales->getLUT_Sin(lutName0, phiIndex0);
1128 
1129         long long cosPhi2LUT = m_gtScales->getLUT_Cos(lutName1, phiIndex1);
1130         long long sinPhi2LUT = m_gtScales->getLUT_Sin(lutName1, phiIndex1);
1131 
1132         // now get pt LUTs
1133         std::string lutName = lutObj0;
1134         lutName += "-ET";
1135         long long ptObj0 = m_gtScales->getLUT_Pt("TwoBody_" + lutName, etIndex0);
1136         unsigned int precPtLUTObj0 = m_gtScales->getPrec_Pt("TwoBody_" + lutName);
1137 
1138         lutName = lutObj1;
1139         lutName += "-ET";
1140         long long ptObj1 = m_gtScales->getLUT_Pt("TwoBody_" + lutName, etIndex1);
1141         unsigned int precPtLUTObj1 = m_gtScales->getPrec_Pt("TwoBody_" + lutName);
1142 
1143         LogTrace("L1TGlobal") << " TBPT Trig precisions:\t " << precCosLUT0 << "\t" << precCosLUT1 << "\t"
1144                               << precSinLUT0 << "\t" << precSinLUT1;
1145         LogTrace("L1TGlobal") << " TBPT Pt precisions:\t " << precPtLUTObj0 << "\t" << precPtLUTObj1;
1146         LogTrace("L1TGlobal") << " TBPT Pt cut:\t " << corrPar.minTBPTCutValue << "\tPrecTBPTCut\t"
1147                               << corrPar.precTBPTCut;
1148         LogTrace("L1TGlobal") << " TBPT Pt1*Pt1 -- Phys:\t " << et0Phy * et0Phy << "\tHW:\t"
1149                               << ptObj0 * ptObj0 * (pow(10, 6));
1150         LogTrace("L1TGlobal") << " TBPT Pt2*Pt2 -- Phys:\t " << et1Phy * et1Phy << "\tHW:\t"
1151                               << ptObj1 * ptObj1 * (pow(10, 6));
1152         LogTrace("L1TGlobal") << " TBPT 2Pt1*Pt2 -- Phys:\t " << 2 * et0Phy * et1Phy << "\tHW:\t"
1153                               << 2 * (ptObj0 * pow(10, 0)) * (ptObj1 * pow(10, 0));
1154         LogTrace("L1TGlobal") << " TBPT Trig -- Phys:\t " << cosPhi1Phy * cosPhi2Phy + sinPhi1Phy * sinPhi2Phy
1155                               << "\tHW:\t" << cosPhi1LUT * cosPhi2LUT + sinPhi1LUT * sinPhi2LUT;
1156 
1157         //double tbptSqPhy =   et0Phy*et0Phy             + et1Phy*et1Phy + 2*et0Phy*et1Phy*(cosPhi1Phy*cosPhi2Phy + sinPhi1Phy*sinPhi2Phy);
1158         long long tbptSqHW = ptObj0 * ptObj0 * (pow(10, 2 * precCosLUT0)) +
1159                              ptObj1 * ptObj1 * (pow(10, 2 * precCosLUT0)) +
1160                              2 * ptObj0 * ptObj1 * (cosPhi1LUT * cosPhi2LUT + sinPhi1LUT * sinPhi2LUT);
1161 
1162         unsigned int preShift = precPtLUTObj0 + precPtLUTObj1 + 2 * precCosLUT0;
1163 
1164         LogTrace("L1TGlobal") << "TBPT Result -- Phys: " << tbptSqPhy << "\tHW: " << tbptSqHW << "\tShifted\t"
1165                               << tbptSqHW / pow(10, preShift) << std::endl;
1166 
1167         preShift = preShift - corrPar.precTBPTCut;
1168 
1169         LogDebug("L1TGlobal")
1170             << "    Testing Two Body Pt Cut (" << lutObj0 << "," << lutObj1 << ") ["
1171             << (long long)(corrPar.minTBPTCutValue * pow(10, preShift)) << ","
1172             << (long long)(corrPar.maxTBPTCutValue * pow(10, preShift)) << "] with precision = " << corrPar.precTBPTCut
1173             << "\n"
1174             << "    etIndex0     = " << etIndex0 << "    pt0LUT      = " << ptObj0 << " PhyEt0 = " << et0Phy << "\n"
1175             << "    etIndex1     = " << etIndex1 << "    pt1LUT      = " << ptObj1 << " PhyEt1 = " << et1Phy << "\n"
1176             << "    Precision Shift = " << preShift << "\n"
1177             << "    Sin(phi1): LUT/Phys\t " << sinPhi1LUT << " / " << sinPhi1Phy << "\n"
1178             << "    Sin(phi2): LUT/Phys\t " << sinPhi2LUT << " / " << sinPhi2Phy << "\n"
1179             << "    Cos(phi1): LUT/Phys\t " << cosPhi1LUT << " / " << cosPhi1Phy << "\n"
1180             << "    Cos(phi2): LUT/Phys\t " << cosPhi2LUT << " / " << cosPhi2Phy
1181             << "\n"
1182 
1183             //    << "    deltaPhiLUT = " << deltaPhiLUT << "\n"
1184             //    << "    deltaEtaLUT = " << deltaEtaLUT << "\n"
1185             //    << "    deltaRSqLUT = " << deltaRSq <<  "\n"
1186             //    << "    Precision Shift = " << preShift << "\n"
1187             //    << "    deltaRSqLUT (shift)= " << (deltaRSq/pow(10,preShift+corrPar.precDRCut))   << "\n"
1188             //    << "    deltaRSqPhy = " << deltaRSqPhy
1189             << std::endl;
1190 
1191         if (tbptSqHW > 0. && tbptSqHW >= (long long)(corrPar.minTBPTCutValue * pow(10, preShift))) {
1192           LogDebug("L1TGlobal") << "    Passed Two Body pT Cut ["
1193                                 << (long long)(corrPar.minTBPTCutValue * pow(10, preShift)) << "]"
1194                                 << "\twith value: " << tbptSqHW << "\n"
1195                                 << "\tPhysics Cut[" << corrPar.minTBPTCutValue / pow(10, corrPar.precTBPTCut)
1196                                 << "]\tPhysics Value: " << tbptSqPhy << std::endl;
1197 
1198         } else {
1199           LogDebug("L1TGlobal") << "    Failed Two Body pT Cut ["
1200                                 << (long long)(corrPar.minTBPTCutValue * pow(10, preShift)) << "]"
1201                                 << "\t with value: " << tbptSqHW << "\n"
1202                                 << "\tPhysics Cut[" << corrPar.minTBPTCutValue / pow(10, corrPar.precTBPTCut)
1203                                 << "]\tPhysics Value: " << tbptSqPhy << std::endl;
1204           reqResult = false;
1205         }
1206       }
1207 
1208       if (corrPar.corrCutType & 0x8 || corrPar.corrCutType & 0x10 || corrPar.corrCutType & 0x40 ||
1209           corrPar.corrCutType & 0x80) {  // added 0x40 for massUpt; added 0x80 for mass/dR
1210         //invariant mass calculation based on
1211         // M = sqrt(2*p1*p2(cosh(eta1-eta2) - cos(phi1 - phi2)))
1212         // but we calculate (1/2)M^2
1213         //
1214         double cosDeltaPhiPhy = cos(deltaPhiPhy);
1215         double coshDeltaEtaPhy = cosh(deltaEtaPhy);
1216         if (corrPar.corrCutType & 0x10)
1217           coshDeltaEtaPhy = 1.;
1218         double massSqPhy = et0Phy * et1Phy * (coshDeltaEtaPhy - cosDeltaPhiPhy);
1219         if (corrPar.corrCutType & 0x40)  // Added for displaced muons
1220           massSqPhy = upt0Phy * upt1Phy * (coshDeltaEtaPhy - cosDeltaPhiPhy);
1221 
1222         long long cosDeltaPhiLUT = m_gtScales->getLUT_DeltaPhi_Cos(lutName, deltaPhiFW);
1223         unsigned int precCosLUT = m_gtScales->getPrec_DeltaPhi_Cos(lutName);
1224 
1225         long long coshDeltaEtaLUT;
1226         if (corrPar.corrCutType & 0x10) {
1227           coshDeltaEtaLUT = 1 * pow(10, precCosLUT);
1228         } else {
1229           coshDeltaEtaLUT = m_gtScales->getLUT_DeltaEta_Cosh(lutName, deltaEtaFW);
1230           unsigned int precCoshLUT = m_gtScales->getPrec_DeltaEta_Cosh(lutName);
1231           if (precCoshLUT - precCosLUT != 0)
1232             LogDebug("L1TGlobal") << "Warning: Cos and Cosh LUTs on different Precision" << std::endl;
1233         }
1234         // if (corrPar.corrCutType & 0x10) coshDeltaEtaLUT=1*pow(10,precCosLUT);
1235 
1236         std::string lutName = lutObj0;
1237         lutName += "-ET";
1238         long long ptObj0 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex0);
1239         unsigned int precPtLUTObj0 = m_gtScales->getPrec_Pt("Mass_" + lutName);
1240 
1241         lutName = lutObj1;
1242         lutName += "-ET";
1243         long long ptObj1 = m_gtScales->getLUT_Pt("Mass_" + lutName, etIndex1);
1244         unsigned int precPtLUTObj1 = m_gtScales->getPrec_Pt("Mass_" + lutName);
1245 
1246         if (corrPar.corrCutType & 0x40)  // Added for displaced muons
1247         {
1248           lutName = lutObj0;
1249           lutName += "-UPT";
1250           ptObj0 = m_gtScales->getLUT_Upt("Mass_" + lutName, uptIndex0);
1251           precPtLUTObj0 = m_gtScales->getPrec_Upt("Mass_" + lutName);
1252 
1253           lutName = lutObj1;
1254           lutName += "-UPT";
1255           ptObj1 = m_gtScales->getLUT_Upt("Mass_" + lutName, uptIndex1);
1256           precPtLUTObj1 = m_gtScales->getPrec_Upt("Mass_" + lutName);
1257         }
1258 
1259         // Pt and Angles are at different precission.
1260         long long massSq = ptObj0 * ptObj1 * (coshDeltaEtaLUT - cosDeltaPhiLUT);
1261 
1262         //Note: There is an assumption here that Cos and Cosh have the same precission
1263         unsigned int preShift = precPtLUTObj0 + precPtLUTObj1 + precCosLUT - corrPar.precMassCut;
1264 
1265         LogDebug("L1TGlobal") << "    Testing Invariant Mass (" << lutObj0 << "," << lutObj1 << ") ["
1266                               << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
1267                               << (long long)(corrPar.maxMassCutValue * pow(10, preShift))
1268                               << "] with precision = " << corrPar.precMassCut << "\n"
1269                               << "    deltaPhiLUT  = " << deltaPhiLUT << "  cosLUT  = " << cosDeltaPhiLUT << "\n"
1270                               << "    deltaEtaLUT  = " << deltaEtaLUT << "  coshLUT = " << coshDeltaEtaLUT << "\n"
1271                               << "    etIndex0     = " << etIndex0 << "    pt0LUT      = " << ptObj0
1272                               << " PhyEt0 = " << et0Phy << "\n"
1273                               << "    etIndex1     = " << etIndex1 << "    pt1LUT      = " << ptObj1
1274                               << " PhyEt1 = " << et1Phy << "\n"
1275                               << "    massSq/2     = " << massSq << "\n"
1276                               << "    Precision Shift = " << preShift << "\n"
1277                               << "    massSq   (shift)= " << (massSq / pow(10, preShift + corrPar.precMassCut)) << "\n"
1278                               << "    deltaPhiPhy  = " << deltaPhiPhy << "  cos() = " << cosDeltaPhiPhy << "\n"
1279                               << "    deltaEtaPhy  = " << deltaEtaPhy << "  cosh()= " << coshDeltaEtaPhy << "\n"
1280                               << "    massSqPhy/2  = " << massSqPhy
1281                               << "  sqrt(|massSq|) = " << sqrt(fabs(2. * massSqPhy)) << std::endl;
1282 
1283         //if(preShift>0) massSq /= pow(10,preShift);
1284         if (corrPar.corrCutType & 0x80) {  //deal with the Invariant Mass Over Delta R cut
1285           // The following is the most precise indexing for 1/dr2 LUT iEta and jPhi - it is based on the physical VALUES for deta and dphi:
1286           //             unsigned int inverseDeltaRIndexEta = int(round(abs(deltaEtaPhy / tempdiffeta)));
1287           //             unsigned int inverseDeltaRIndexPhi = int(round(abs(deltaPhiPhy/ tempdiffphi)));
1288           // To match the FW we instead must use these INDIVIDUAL indecies for eta and phi:
1289           // (NOTE: The following treatment of delta eta and delta phi indecies matches the most precise case above):
1290           int iEta0 = int(trunc(etaIndex0 * resolution));
1291           int iEta1 = int(trunc(etaIndex1 * resolution));
1292           unsigned int inverseDeltaRIndexEta = abs(iEta0 - iEta1);
1293           int jPhi0 = int(trunc(phiIndex0 * resolution));
1294           int jPhi1 = int(trunc(phiIndex1 * resolution));
1295           unsigned int inverseDeltaRIndexPhi = abs(jPhi0 - jPhi1);
1296           if (abs(phiIndex0 - phiIndex1) >= phiBound)
1297             inverseDeltaRIndexPhi = int(trunc(2 * phiBound * resolution)) - inverseDeltaRIndexPhi;
1298           unsigned int precInvDRSqLUT = 0x5;
1299           long long LInverseDeltaRSqLUT = 0x0;
1300           if (inverseDeltaRIndexEta + inverseDeltaRIndexPhi > 0)
1301             LInverseDeltaRSqLUT = InvDeltaRSqLUT[inverseDeltaRIndexEta][inverseDeltaRIndexPhi];
1302           long long massSqOverDeltaRSq =
1303               (long long)massSq * LInverseDeltaRSqLUT;  // keep full precision, do not '/pow(10,precInvDRSqLUT);'
1304           if ((inverseDeltaRIndexEta + inverseDeltaRIndexPhi == 0) ||
1305               ((massSqOverDeltaRSq >= 0) &&
1306                (massSqOverDeltaRSq >=
1307                 (long long)(corrPar.minMassCutValue * pow(10, preShift) *
1308                             pow(10, precInvDRSqLUT))))) {  // only check for the minimum threshold in case of invMass/dR
1309             LogDebug("L1TGlobal") << "    Passed Invariant Mass Over Delta R  Cut ["
1310                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift) * pow(10, precInvDRSqLUT))
1311                                   << ","
1312                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift) * pow(10, precInvDRSqLUT))
1313                                   << "]     massSqOverDeltaRSq " << massSqOverDeltaRSq << std::endl;
1314           } else {
1315             LogDebug("L1TGlobal") << "    Failed Invariant Mass OverDeltaR Cut ["
1316                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift) * pow(10, precInvDRSqLUT))
1317                                   << ","
1318                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift) * pow(10, precInvDRSqLUT))
1319                                   << "]     massSqOverDeltaRSq " << massSqOverDeltaRSq << std::endl;
1320             reqResult = false;
1321           }
1322           //Done with Invariant Mass Over Delta R vs Mass Cut choice
1323         } else {  //do the InvMassCut choice
1324           if (massSq >= 0 && massSq >= (long long)(corrPar.minMassCutValue * pow(10, preShift)) &&
1325               massSq <= (long long)(corrPar.maxMassCutValue * pow(10, preShift))) {
1326             LogDebug("L1TGlobal") << "    Passed Invariant Mass Cut ["
1327                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
1328                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
1329           } else {
1330             LogDebug("L1TGlobal") << "    Failed Invariant Mass Cut ["
1331                                   << (long long)(corrPar.minMassCutValue * pow(10, preShift)) << ","
1332                                   << (long long)(corrPar.maxMassCutValue * pow(10, preShift)) << "]" << std::endl;
1333             reqResult = false;
1334           }  //Done with Invariant Mass Cut
1335         }    //Done with choice of Invariant Mass Cut vs InvMass/dR
1336       }      //Done with any type of Mass Cut
1337 
1338       // For Muon-Muon Correlation Check the Charge Correlation if requested
1339       bool chrgCorrel = true;
1340       if (cond0Categ == CondMuon && cond1Categ == CondMuon) {
1341         // Check for like-sign
1342         if (corrPar.chargeCorrelation == 2 && chrg0 != chrg1)
1343           chrgCorrel = false;
1344         // Check for opp-sign
1345         if (corrPar.chargeCorrelation == 4 && chrg0 == chrg1)
1346           chrgCorrel = false;
1347       }
1348 
1349       if (reqResult & chrgCorrel) {
1350         condResult = true;
1351         (combinationsInCond()).push_back(objectsInComb);
1352       }
1353 
1354     }  //end loop over second leg
1355 
1356   }  //end loop over first leg
1357 
1358   if (m_verbosity && condResult) {
1359     LogDebug("L1TGlobal") << " pass(es) the correlation condition.\n" << std::endl;
1360   }
1361 
1362   return condResult;
1363 }
1364 
1365 // load calo candidates
1366 const l1t::L1Candidate* l1t::CorrCondition::getCandidate(const int bx, const int indexCand) const {
1367   // objectType() gives the type for nrObjects() only,
1368   // but in a CondCalo all objects have the same type
1369   // take type from the type of the first object
1370   switch ((m_gtCorrelationTemplate->objectType())[0]) {
1371     case gtEG:
1372       return (m_uGtB->getCandL1EG())->at(bx, indexCand);
1373       break;
1374 
1375     case gtJet:
1376       return (m_uGtB->getCandL1Jet())->at(bx, indexCand);
1377       break;
1378 
1379     case gtTau:
1380       return (m_uGtB->getCandL1Tau())->at(bx, indexCand);
1381       break;
1382     default:
1383       return nullptr;
1384       break;
1385   }
1386 
1387   return nullptr;
1388 }
1389 
1390 /**
1391  * checkObjectParameter - Compare a single particle with a numbered condition.
1392  *
1393  * @param iCondition The number of the condition.
1394  * @param cand The candidate to compare.
1395  *
1396  * @return The result of the comparison (false if a condition does not exist).
1397  */
1398 
1399 const bool l1t::CorrCondition::checkObjectParameter(const int iCondition, const l1t::L1Candidate& cand) const {
1400   return true;
1401 }
1402 
1403 void l1t::CorrCondition::print(std::ostream& myCout) const {
1404   myCout << "Dummy Print for CorrCondition" << std::endl;
1405   m_gtCorrelationTemplate->print(myCout);
1406 
1407   ConditionEvaluation::print(myCout);
1408 }