Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:02:26

0001 #include "CondFormats/PCLConfig/interface/AlignPCLThresholds.h"
0002 #include "CondFormats/PCLConfig/interface/AlignPCLThreshold.h"
0003 #include "FWCore/Utilities/interface/Exception.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include <iostream>
0006 #include <iomanip>  // std::setw
0007 
0008 //****************************************************************************//
0009 void AlignPCLThresholds::setAlignPCLThreshold(const std::string &AlignableId, const AlignPCLThreshold &Threshold) {
0010   m_thresholds[AlignableId] = Threshold;
0011 }
0012 
0013 //****************************************************************************//
0014 void AlignPCLThresholds::setAlignPCLThresholds(const int &Nrecords, const threshold_map &AlignPCLThresholds) {
0015   m_nrecords = Nrecords;
0016   m_thresholds = AlignPCLThresholds;
0017 }
0018 
0019 //****************************************************************************//
0020 void AlignPCLThresholds::setNRecords(const int &Nrecords) { m_nrecords = Nrecords; }
0021 
0022 //****************************************************************************//
0023 AlignPCLThreshold AlignPCLThresholds::getAlignPCLThreshold(const std::string &AlignableId) const {
0024   threshold_map::const_iterator it = m_thresholds.find(AlignableId);
0025 
0026   if (it != m_thresholds.end()) {
0027     return it->second;
0028   } else {
0029     throw cms::Exception("AlignPCLThresholds") << "No Thresholds defined for Alignable id " << AlignableId << "\n";
0030   }
0031 }
0032 
0033 //****************************************************************************//
0034 AlignPCLThreshold &AlignPCLThresholds::getAlignPCLThreshold(const std::string &AlignableId) {
0035   return m_thresholds[AlignableId];
0036 }
0037 
0038 float AlignPCLThresholds::getSigCut(const std::string &AlignableId, const coordType &type) const {
0039   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0040   switch (type) {
0041     case X:
0042       return a.getSigXcut();
0043     case Y:
0044       return a.getSigYcut();
0045     case Z:
0046       return a.getSigZcut();
0047     case theta_X:
0048       return a.getSigThetaXcut();
0049     case theta_Y:
0050       return a.getSigThetaYcut();
0051     case theta_Z:
0052       return a.getSigThetaZcut();
0053     default:
0054       throw cms::Exception("AlignPCLThresholds")
0055           << "Requested significance threshold for undefined coordinate" << type << "\n";
0056   }
0057 }
0058 
0059 //****************************************************************************//
0060 // overloaded method
0061 //****************************************************************************//
0062 std::array<float, 6> AlignPCLThresholds::getSigCut(const std::string &AlignableId) const {
0063   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0064   return {
0065       {a.getSigXcut(), a.getSigYcut(), a.getSigZcut(), a.getSigThetaXcut(), a.getSigThetaYcut(), a.getSigThetaZcut()}};
0066 }
0067 
0068 float AlignPCLThresholds::getCut(const std::string &AlignableId, const coordType &type) const {
0069   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0070   switch (type) {
0071     case X:
0072       return a.getXcut();
0073     case Y:
0074       return a.getYcut();
0075     case Z:
0076       return a.getZcut();
0077     case theta_X:
0078       return a.getThetaXcut();
0079     case theta_Y:
0080       return a.getThetaYcut();
0081     case theta_Z:
0082       return a.getThetaZcut();
0083     default:
0084       throw cms::Exception("AlignPCLThresholds")
0085           << "Requested significance threshold for undefined coordinate" << type << "\n";
0086   }
0087 }
0088 
0089 //****************************************************************************//
0090 // overloaded method
0091 //****************************************************************************//
0092 std::array<float, 6> AlignPCLThresholds::getCut(const std::string &AlignableId) const {
0093   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0094   return {{a.getXcut(), a.getYcut(), a.getZcut(), a.getThetaXcut(), a.getThetaYcut(), a.getThetaZcut()}};
0095 }
0096 
0097 //****************************************************************************//
0098 float AlignPCLThresholds::getMaxMoveCut(const std::string &AlignableId, const coordType &type) const {
0099   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0100   switch (type) {
0101     case X:
0102       return a.getMaxMoveXcut();
0103     case Y:
0104       return a.getMaxMoveYcut();
0105     case Z:
0106       return a.getMaxMoveZcut();
0107     case theta_X:
0108       return a.getMaxMoveThetaXcut();
0109     case theta_Y:
0110       return a.getMaxMoveThetaYcut();
0111     case theta_Z:
0112       return a.getMaxMoveThetaZcut();
0113     default:
0114       throw cms::Exception("AlignPCLThresholds")
0115           << "Requested significance threshold for undefined coordinate" << type << "\n";
0116   }
0117 }
0118 
0119 //****************************************************************************//
0120 // overloaded method
0121 //****************************************************************************//
0122 std::array<float, 6> AlignPCLThresholds::getMaxMoveCut(const std::string &AlignableId) const {
0123   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0124   return {{a.getMaxMoveXcut(),
0125            a.getMaxMoveYcut(),
0126            a.getMaxMoveZcut(),
0127            a.getMaxMoveThetaXcut(),
0128            a.getMaxMoveThetaYcut(),
0129            a.getMaxMoveThetaZcut()}};
0130 }
0131 
0132 //****************************************************************************//
0133 float AlignPCLThresholds::getMaxErrorCut(const std::string &AlignableId, const coordType &type) const {
0134   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0135   switch (type) {
0136     case X:
0137       return a.getErrorXcut();
0138     case Y:
0139       return a.getErrorYcut();
0140     case Z:
0141       return a.getErrorZcut();
0142     case theta_X:
0143       return a.getErrorThetaXcut();
0144     case theta_Y:
0145       return a.getErrorThetaYcut();
0146     case theta_Z:
0147       return a.getErrorThetaZcut();
0148     default:
0149       throw cms::Exception("AlignPCLThresholds")
0150           << "Requested significance threshold for undefined coordinate" << type << "\n";
0151   }
0152 }
0153 
0154 //****************************************************************************//
0155 // overloaded method
0156 //****************************************************************************//
0157 std::array<float, 6> AlignPCLThresholds::getMaxErrorCut(const std::string &AlignableId) const {
0158   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0159   return {{a.getErrorXcut(),
0160            a.getErrorYcut(),
0161            a.getErrorZcut(),
0162            a.getErrorThetaXcut(),
0163            a.getErrorThetaYcut(),
0164            a.getErrorThetaZcut()}};
0165 }
0166 
0167 //****************************************************************************//
0168 std::array<float, 4> AlignPCLThresholds::getExtraDOFCutsForAlignable(const std::string &AlignableId,
0169                                                                      const unsigned int i) const {
0170   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0171   return a.getExtraDOFCuts(i);
0172 }
0173 
0174 //****************************************************************************//
0175 std::string AlignPCLThresholds::getExtraDOFLabelForAlignable(const std::string &AlignableId,
0176                                                              const unsigned int i) const {
0177   AlignPCLThreshold a = getAlignPCLThreshold(AlignableId);
0178   return a.getExtraDOFLabel(i);
0179 }
0180 
0181 //****************************************************************************//
0182 void AlignPCLThresholds::printAll() const {
0183   edm::LogVerbatim("AlignPCLThresholds") << "AlignPCLThresholds::printAll()";
0184   edm::LogVerbatim("AlignPCLThresholds") << " ========================================================================="
0185                                             "==========================================";
0186   edm::LogVerbatim("AlignPCLThresholds") << "N records cut: " << this->getNrecords();
0187   for (auto it = m_thresholds.begin(); it != m_thresholds.end(); ++it) {
0188     edm::LogVerbatim("AlignPCLThresholds") << " ======================================================================="
0189                                               "============================================";
0190     edm::LogVerbatim("AlignPCLThresholds")
0191         << "key : " << it->first << " \n"
0192         << "- Xcut             : " << std::setw(4) << (it->second).getXcut() << std::setw(5) << "   um"
0193         << "| sigXcut          : " << std::setw(4) << (it->second).getSigXcut() << std::setw(1) << " "
0194         << "| maxMoveXcut      : " << std::setw(4) << (it->second).getMaxMoveXcut() << std::setw(5) << "   um"
0195         << "| ErrorXcut        : " << std::setw(4) << (it->second).getErrorXcut() << std::setw(5) << "   um\n"
0196 
0197         << "- thetaXcut        : " << std::setw(4) << (it->second).getThetaXcut() << std::setw(5) << " urad"
0198         << "| sigThetaXcut     : " << std::setw(4) << (it->second).getSigThetaXcut() << std::setw(1) << " "
0199         << "| maxMoveThetaXcut : " << std::setw(4) << (it->second).getMaxMoveThetaXcut() << std::setw(5) << " urad"
0200         << "| ErrorThetaXcut   : " << std::setw(4) << (it->second).getErrorThetaXcut() << std::setw(5) << " urad\n"
0201 
0202         << "- Ycut             : " << std::setw(4) << (it->second).getYcut() << std::setw(5) << "   um"
0203         << "| sigYcut          : " << std::setw(4) << (it->second).getSigXcut() << std::setw(1) << " "
0204         << "| maxMoveYcut      : " << std::setw(4) << (it->second).getMaxMoveYcut() << std::setw(5) << "   um"
0205         << "| ErrorYcut        : " << std::setw(4) << (it->second).getErrorYcut() << std::setw(5) << "   um\n"
0206 
0207         << "- thetaYcut        : " << std::setw(4) << (it->second).getThetaYcut() << std::setw(5) << " urad"
0208         << "| sigThetaYcut     : " << std::setw(4) << (it->second).getSigThetaYcut() << std::setw(1) << " "
0209         << "| maxMoveThetaYcut : " << std::setw(4) << (it->second).getMaxMoveThetaYcut() << std::setw(5) << " urad"
0210         << "| ErrorThetaYcut   : " << std::setw(4) << (it->second).getErrorThetaYcut() << std::setw(5) << " urad\n"
0211 
0212         << "- Zcut             : " << std::setw(4) << (it->second).getZcut() << std::setw(5) << "   um"
0213         << "| sigZcut          : " << std::setw(4) << (it->second).getSigZcut() << std::setw(1) << " "
0214         << "| maxMoveZcut      : " << std::setw(4) << (it->second).getMaxMoveZcut() << std::setw(5) << "   um"
0215         << "| ErrorZcut        : " << std::setw(4) << (it->second).getErrorZcut() << std::setw(5) << "   um\n"
0216 
0217         << "- thetaZcut        : " << std::setw(4) << (it->second).getThetaZcut() << std::setw(5) << " urad"
0218         << "| sigThetaZcut     : " << std::setw(4) << (it->second).getSigThetaZcut() << std::setw(1) << " "
0219         << "| maxMoveThetaZcut : " << std::setw(4) << (it->second).getMaxMoveThetaZcut() << std::setw(5) << " urad"
0220         << "| ErrorThetaZcut   : " << std::setw(4) << (it->second).getErrorThetaZcut() << std::setw(5) << " urad";
0221 
0222     if ((it->second).hasExtraDOF()) {
0223       for (unsigned int j = 0; j < (it->second).extraDOFSize(); j++) {
0224         std::array<float, 4> extraDOFCuts = getExtraDOFCutsForAlignable(it->first, j);
0225 
0226         edm::LogVerbatim("AlignPCLThresholds")
0227             << "Extra DOF " << j << " with label: " << getExtraDOFLabelForAlignable(it->first, j);
0228         edm::LogVerbatim("AlignPCLThresholds")
0229             << "- cut              : " << std::setw(4) << extraDOFCuts.at(0) << std::setw(5) << "    "
0230             << "| sigCut           : " << std::setw(4) << extraDOFCuts.at(1) << std::setw(1) << " "
0231             << "| maxMoveCut       : " << std::setw(4) << extraDOFCuts.at(2) << std::setw(5) << "    "
0232             << "| maxErrorCut      : " << std::setw(4) << extraDOFCuts.at(3) << std::setw(5) << "    ";
0233       }
0234     }
0235   }
0236 }
0237 
0238 //****************************************************************************//
0239 std::vector<std::string> AlignPCLThresholds::getAlignableList() const {
0240   std::vector<std::string> alignables_;
0241   alignables_.reserve(m_thresholds.size());
0242 
0243   for (auto it = m_thresholds.begin(); it != m_thresholds.end(); ++it) {
0244     alignables_.push_back(it->first);
0245   }
0246   return alignables_;
0247 }