GlobalCaloTrigger/src/L1GctHtMissLut.cc

0001 #include "L1Trigger/GlobalCaloTrigger/interface/L1GctHtMissLut.h"
0002
0003 #include "CondFormats/L1TObjects/interface/L1CaloEtScale.h"
0004
0005 #include <cmath>
0006
0007 //DEFINE STATICS
0008 const int L1GctHtMissLut::NAddress = 2 * L1GctHtMissLut::kHxOrHyMissComponentNBits;
0009 const int L1GctHtMissLut::NData = L1GctHtMissLut::kHtMissMagnitudeNBits + L1GctHtMissLut::kHtMissAngleNBits;
0010
0011 L1GctHtMissLut::L1GctHtMissLut(const L1CaloEtScale* const scale, const double lsb)
0012     : L1GctLut<NAddress, NData>(), m_etScale(scale), m_componentLsb(lsb) {
0013   if (scale != nullptr)
0014     m_setupOk = true;
0015 }
0016
0017 L1GctHtMissLut::L1GctHtMissLut() : L1GctLut<NAddress, NData>(), m_etScale(nullptr), m_componentLsb(1.0) {}
0018
0019 L1GctHtMissLut::L1GctHtMissLut(const L1GctHtMissLut& lut)
0020     : L1GctLut<NAddress, NData>(), m_etScale(lut.etScale()), m_componentLsb(lut.componentLsb()) {
0021   if (m_etScale != nullptr)
0022     m_setupOk = true;
0023 }
0024
0025 L1GctHtMissLut::~L1GctHtMissLut() {}
0026
0027 uint16_t L1GctHtMissLut::value(const uint16_t lutAddress) const {
0028   uint16_t result = 0;
0029
0030   if (lutAddress != 0) {
0031     static const int maxComponent = 1 << kHxOrHyMissComponentNBits;
0032     static const int componentMask = maxComponent - 1;
0033
0034     static const int magnitudeMask = (1 << kHtMissMagnitudeNBits) - 1;
0035     static const int angleMask = (1 << kHtMissAngleNBits) - 1;
0036
0037     // Extract the bits corresponding to hx and hy components
0038     int hxCompGct = static_cast<int>(lutAddress >> kHxOrHyMissComponentNBits) & componentMask;
0039     int hyCompGct = static_cast<int>(lutAddress) & componentMask;
0040
0041     // These are twos-complement integers - if the MSB is set, the value is negative
0042     if (hxCompGct >= maxComponent / 2)
0043       hxCompGct -= maxComponent;
0044     if (hyCompGct >= maxComponent / 2)
0045       hyCompGct -= maxComponent;
0046
0047     // Convert to GeV. Add 0.5 to each component to compensate for truncation errors.
0048     double hxCompGeV = m_componentLsb * (static_cast<double>(hxCompGct) + 0.5);
0049     double hyCompGeV = m_componentLsb * (static_cast<double>(hyCompGct) + 0.5);
0050
0051     // Convert to magnitude and angle
0052     double htMissMag = sqrt(hxCompGeV * hxCompGeV + hyCompGeV * hyCompGeV);
0053     double htMissAng = atan2(hyCompGeV, hxCompGeV);
0054     if (htMissAng < 0.0)
0055       htMissAng += 2.0 * M_PI;
0056
0057     // Convert back to integer
0058     int htMissMagBits = static_cast<int>(m_etScale->rank(htMissMag)) & magnitudeMask;
0059     int htMissAngBits = static_cast<int>(htMissAng * 9.0 / M_PI) & angleMask;
0060
0061     // Form the lut output
0062     result = (htMissMagBits << kHtMissAngleNBits) | htMissAngBits;
0063   }
0064
0065   return result;
0066 }
0067
0068 std::vector<double> L1GctHtMissLut::getThresholdsGeV() const { return m_etScale->getThresholds(); }
0069
0070 std::vector<unsigned> L1GctHtMissLut::getThresholdsGct() const {
0071   std::vector<unsigned> result;
0072   std::vector<double> thresholdsGeV = m_etScale->getThresholds();
0073   for (std::vector<double>::const_iterator thr = thresholdsGeV.begin(); thr != thresholdsGeV.end(); thr++) {
0074     result.push_back(static_cast<unsigned>((*thr) / (m_componentLsb)));
0075   }
0076   return result;
0077 }
0078
0079 L1GctHtMissLut L1GctHtMissLut::operator=(const L1GctHtMissLut& lut) {
0080   const L1GctHtMissLut& temp(lut);
0081   return temp;
0082 }
0083
0084 std::ostream& operator<<(std::ostream& os, const L1GctHtMissLut& lut) {
0085   os << "===L1GctHtMissLut===" << std::endl;
0086   std::vector<double> thresholds = lut.m_etScale->getThresholds();
0087   std::vector<double>::const_iterator thr = thresholds.begin();
0088   os << "Thresholds are: " << *(thr++);
0089   for (; thr != thresholds.end(); thr++) {
0090     os << ", " << *thr;
0091   }
0092   os << std::endl;
0093   os << "Max values for input to et scale " << lut.m_etScale->linScaleMax() << " and for output "
0094      << lut.m_etScale->rankScaleMax() << std::endl;
0095   os << "LSB used for conversion is " << lut.m_componentLsb << " GeV" << std::endl;
0096   os << "\n===Lookup table contents===\n" << std::endl;
0097   const L1GctLut<L1GctHtMissLut::NAddress, L1GctHtMissLut::NData>* temp = &lut;
0098   os << *temp;
0099   return os;
0100 }
0101
0102 template class L1GctLut<L1GctHtMissLut::NAddress, L1GctHtMissLut::NData>;