SiStripCommissioningAnalysis/src/FastFedCablingAlgorithm.cc

0001 #include "DQM/SiStripCommissioningAnalysis/interface/FastFedCablingAlgorithm.h"
0002 #include "CondFormats/SiStripObjects/interface/FastFedCablingAnalysis.h"
0003 #include "DataFormats/SiStripCommon/interface/SiStripHistoTitle.h"
0004 #include "DataFormats/SiStripCommon/interface/SiStripEnumsAndStrings.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "TProfile.h"
0007 #include "TH1.h"
0008 #include <iostream>
0009 #include <sstream>
0010 #include <iomanip>
0011 #include <cmath>
0012
0013 using namespace sistrip;
0014
0015 // ----------------------------------------------------------------------------
0016 //
0017 FastFedCablingAlgorithm::FastFedCablingAlgorithm(const edm::ParameterSet& pset, FastFedCablingAnalysis* const anal)
0018     : CommissioningAlgorithm(anal), histo_(nullptr, "") {
0019   ;
0020 }
0021
0022 // ----------------------------------------------------------------------------
0023 //
0024 void FastFedCablingAlgorithm::extract(const std::vector<TH1*>& histos) {
0025   if (!anal()) {
0026     edm::LogWarning(mlCommissioning_) << "[FastFedCablingAlgorithm::" << __func__ << "]"
0027                                       << " NULL pointer to Analysis object!";
0028     return;
0029   }
0030
0031   // Check number of histograms
0032   if (histos.size() != 1) {
0033     anal()->addErrorCode(sistrip::numberOfHistos_);
0034   }
0035
0036   // Extract FED key from histo title
0037   if (!histos.empty()) {
0038     anal()->fedKey(extractFedKey(histos.front()));
0039   }
0040
0041   // Extract histograms
0042   std::vector<TH1*>::const_iterator ihis = histos.begin();
0043   for (; ihis != histos.end(); ihis++) {
0044     // Check for NULL pointer
0045     if (!(*ihis)) {
0046       continue;
0047     }
0048
0049     // Check name
0050     SiStripHistoTitle title((*ihis)->GetName());
0051     if (title.runType() != sistrip::FAST_CABLING) {
0052       anal()->addErrorCode(sistrip::unexpectedTask_);
0053       continue;
0054     }
0055
0056     // Extract cabling histo
0057     histo_.first = *ihis;
0058     histo_.second = (*ihis)->GetName();
0059   }
0060 }
0061
0062 // -----------------------------------------------------------------------------
0063 //
0064 void FastFedCablingAlgorithm::analyse() {
0065   if (!anal()) {
0066     edm::LogWarning(mlCommissioning_) << "[FastFedCablingAlgorithm::" << __func__ << "]"
0067                                       << " NULL pointer to base Analysis object!";
0068     return;
0069   }
0070
0071   CommissioningAnalysis* tmp = const_cast<CommissioningAnalysis*>(anal());
0072   FastFedCablingAnalysis* anal = dynamic_cast<FastFedCablingAnalysis*>(tmp);
0073   if (!anal) {
0074     edm::LogWarning(mlCommissioning_) << "[FastFedCablingAlgorithm::" << __func__ << "]"
0075                                       << " NULL pointer to derived Analysis object!";
0076     return;
0077   }
0078
0079   if (!histo_.first) {
0080     anal->addErrorCode(sistrip::nullPtr_);
0081     return;
0082   }
0083
0084   TProfile* histo = dynamic_cast<TProfile*>(histo_.first);
0085   if (!histo) {
0086     anal->addErrorCode(sistrip::nullPtr_);
0087     return;
0088   }
0089
0090   // Initialization
0091   uint16_t zero_entries = 0;
0092   uint16_t nbins = static_cast<uint16_t>(histo->GetNbinsX());
0093   std::vector<float> contents;
0094   std::vector<float> errors;
0095   std::vector<float> entries;
0096   contents.reserve(nbins);
0097   errors.reserve(nbins);
0098   entries.reserve(nbins);
0099
0100   // Copy histo contents to containers and find min/max
0101   anal->max_ = -1. * sistrip::invalid_;
0102   for (uint16_t ibin = 0; ibin < nbins; ibin++) {
0103     contents.push_back(histo->GetBinContent(ibin + 1));
0104     errors.push_back(histo->GetBinError(ibin + 1));
0105     entries.push_back(histo->GetBinEntries(ibin + 1));
0106     if (entries[ibin]) {
0107       if (contents[ibin] > anal->max_) {
0108         anal->max_ = contents[ibin];
0109       }
0110       if (contents[ibin] < anal->min_) {
0111         anal->min_ = contents[ibin];
0112       }
0113     } else {
0114       zero_entries++;
0115     }
0116   }
0117   if (anal->max_ < -1. * sistrip::valid_) {
0118     anal->max_ = sistrip::invalid_;
0119   }
0120
0121   // Check number of bins
0122   if (contents.size() != FastFedCablingAnalysis::nBitsForDcuId_ + FastFedCablingAnalysis::nBitsForLldCh_) {
0123     anal->addErrorCode(sistrip::numberOfBins_);
0124     return;
0125   }
0126
0127   // Check for bins with zero entries
0128   if (zero_entries) {
0129     anal->addErrorCode(sistrip::noEntries_);
0130     return;
0131   }
0132
0133   // Check min and max found
0134   if (anal->max_ > sistrip::valid_ || anal->min_ > sistrip::valid_) {
0135     return;
0136   }
0137
0138   // Calculate range and mid-range levels
0139   anal->range_ = anal->max_ - anal->min_;
0140   anal->midRange_ = anal->min_ + anal->range_ / 2.;
0141
0142   // Check if range is above threshold
0143   if (anal->range_ < FastFedCablingAnalysis::threshold_) {
0144     anal->addErrorCode(sistrip::smallDataRange_);
0145     return;
0146   }
0147
0148   // Identify samples to be either "low" or "high"
0149   std::vector<float> high;
0150   std::vector<float> low;
0151   for (uint16_t ibin = 0; ibin < nbins; ibin++) {
0152     if (entries[ibin]) {
0153       if (contents[ibin] < anal->midRange_) {
0154         low.push_back(contents[ibin]);
0155       } else {
0156         high.push_back(contents[ibin]);
0157       }
0158     }
0159   }
0160
0161   // Find median of high and low levels
0162   sort(high.begin(), high.end());
0163   sort(low.begin(), low.end());
0164   if (!high.empty()) {
0165     anal->highMedian_ = high[high.size() % 2 ? high.size() / 2 : high.size() / 2];
0166   }
0167   if (!low.empty()) {
0168     anal->lowMedian_ = low[low.size() % 2 ? low.size() / 2 : low.size() / 2];
0169   }
0170
0171   // Check if light levels above thresholds
0172   //if ( anal->highMedian_ < FastFedCablingAnalysis::dirtyThreshold_ ) { anal->addErrorCode(sistrip::invalidLightLevel_); }
0173   //if ( anal->lowMedian_ < FastFedCablingAnalysis::trimDacThreshold_ ) { anal->addErrorCode(sistrip::invalidTrimDacLevel_); }
0174
0175   // Find mean and rms in "low" samples
0176   anal->lowMean_ = 0.;
0177   anal->lowRms_ = 0.;
0178   for (uint16_t ibin = 0; ibin < low.size(); ibin++) {
0179     anal->lowMean_ += low[ibin];
0180     anal->lowRms_ += low[ibin] * low[ibin];
0181   }
0182   if (!low.empty()) {
0183     anal->lowMean_ = anal->lowMean_ / low.size();
0184     anal->lowRms_ = anal->lowRms_ / low.size();
0185   } else {
0186     anal->lowMean_ = 1. * sistrip::invalid_;
0187     anal->lowRms_ = 1. * sistrip::invalid_;
0188   }
0189   if (anal->lowMean_ < sistrip::valid_) {
0190     anal->lowRms_ = sqrt(fabs(anal->lowRms_ - anal->lowMean_ * anal->lowMean_));
0191   } else {
0192     anal->lowMean_ = 1. * sistrip::invalid_;
0193     anal->lowRms_ = 1. * sistrip::invalid_;
0194   }
0195
0196   // Find mean and rms in "high" samples
0197   anal->highMean_ = 0.;
0198   anal->highRms_ = 0.;
0199   for (uint16_t ibin = 0; ibin < high.size(); ibin++) {
0200     anal->highMean_ += high[ibin];
0201     anal->highRms_ += high[ibin] * high[ibin];
0202   }
0203   if (!high.empty()) {
0204     anal->highMean_ = anal->highMean_ / high.size();
0205     anal->highRms_ = anal->highRms_ / high.size();
0206   } else {
0207     anal->highMean_ = 1. * sistrip::invalid_;
0208     anal->highRms_ = 1. * sistrip::invalid_;
0209   }
0210   if (anal->highMean_ < sistrip::valid_) {
0211     anal->highRms_ = sqrt(fabs(anal->highRms_ - anal->highMean_ * anal->highMean_));
0212   } else {
0213     anal->highMean_ = 1. * sistrip::invalid_;
0214     anal->highRms_ = 1. * sistrip::invalid_;
0215   }
0216
0217   // Check if light levels above thresholds
0218   //if ( anal->highMean_ < FastFedCablingAnalysis::dirtyThreshold_ ) { anal->addErrorCode(sistrip::invalidLightLevel_); }
0219   //if ( anal->lowMean_ < FastFedCablingAnalysis::trimDacThreshold_ ) { anal->addErrorCode(sistrip::invalidTrimDacLevel_); }
0220
0221   // Recalculate range
0222   if (anal->highMean_ < 1. * sistrip::valid_ && anal->lowMean_ < 1. * sistrip::valid_) {
0223     anal->range_ = anal->highMean_ - anal->lowMean_;
0224     anal->midRange_ = anal->lowMean_ + anal->range_ / 2.;
0225   } else {
0226     anal->range_ = 1. * sistrip::invalid_;
0227     anal->midRange_ = 1. * sistrip::invalid_;
0228   }
0229
0230   // Check if updated range is valid and above threshold
0231   if (anal->range_ > 1. * sistrip::valid_ || anal->range_ < FastFedCablingAnalysis::threshold_) {
0232     anal->addErrorCode(sistrip::smallDataRange_);
0233     return;
0234   }
0235
0236   // Extract DCU id
0237   anal->dcuHardId_ = 0;
0238   for (uint16_t ibin = 0; ibin < FastFedCablingAnalysis::nBitsForDcuId_; ibin++) {
0239     if (entries[ibin]) {
0240       if (contents[ibin] > anal->midRange_) {
0241         anal->dcuHardId_ += 0xFFFFFFFF & (1 << ibin);
0242       }
0243     }
0244   }
0245   if (!anal->dcuHardId_) {
0246     anal->dcuHardId_ = sistrip::invalid32_;
0247   }
0248
0249   // Extract DCU id
0250   anal->lldCh_ = 0;
0251   for (uint16_t ibin = 0; ibin < FastFedCablingAnalysis::nBitsForLldCh_; ibin++) {
0252     if (entries[FastFedCablingAnalysis::nBitsForDcuId_ + ibin]) {
0253       if (contents[FastFedCablingAnalysis::nBitsForDcuId_ + ibin] > anal->midRange_) {
0254         anal->lldCh_ += (0x3 & (1 << ibin));
0255       }
0256     }
0257   }
0258   anal->lldCh_++;  // starts from 1
0259   if (!anal->lldCh_) {
0260     anal->lldCh_ = sistrip::invalid_;
0261   }
0262 }