ESTrendTask.cc

Go to the documentation of this file.

#include <iostream>
#include <fstream>
#include <vector>
#include <cmath>
 
#include "FWCore/Framework/interface/Run.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "DataFormats/EcalDetId/interface/ESDetId.h"
#include "DQMServices/Core/interface/DQMStore.h"
#include "DQM/EcalPreshowerMonitorModule/interface/ESTrendTask.h"
#include "DataFormats/EcalRawData/interface/ESDCCHeaderBlock.h"
#include "DataFormats/EcalRawData/interface/ESKCHIPBlock.h"
 
using namespace cms;
using namespace edm;
using namespace std;
 
ESTrendTask::ESTrendTask(const ParameterSet& ps) {
  prefixME_ = ps.getUntrackedParameter<string>("prefixME", "EcalPreshower");
  rechittoken_ = consumes<ESRecHitCollection>(ps.getParameter<InputTag>("RecHitLabel"));
  dccCollections_ = consumes<ESRawDataCollection>(ps.getParameter<InputTag>("ESDCCCollections"));
 
  for (int i = 0; i < 2; ++i)
    for (int j = 0; j < 2; ++j) {
      hESRecHitTrend_[i][j] = nullptr;
      hESRecHitTrendHr_[i][j] = nullptr;
    }
 
  hESSLinkErrTrend_ = nullptr;
  hESFiberErrTrend_ = nullptr;
  hESSLinkErrTrendHr_ = nullptr;
  hESFiberErrTrendHr_ = nullptr;
 
  start_time_ = 0;
  current_time_ = 0;
  last_time_ = 0;
 
  ievt_ = 0;
}
 
void ESTrendTask::dqmBeginRun(const Run& r, const EventSetup& c) {
  start_time_ = (r.beginTime()).unixTime();
 
  //std::cout << "start time : " << start_time_ << std::endl;
}
 
void ESTrendTask::bookHistograms(DQMStore::IBooker& iBooker, Run const&, EventSetup const&) {
  char histo[200];
 
  iBooker.setCurrentFolder(prefixME_ + "/ESTrendTask");
 
  for (int i = 0; i < 2; ++i)
    for (int j = 0; j < 2; ++j) {
      int iz = (i == 0) ? 1 : -1;
      sprintf(histo, "ES Trending RH Occ per 5 mins Z %d P %d", iz, j + 1);
      hESRecHitTrend_[i][j] = iBooker.bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
      hESRecHitTrend_[i][j]->setAxisTitle("Elapse time (Minutes)", 1);
      hESRecHitTrend_[i][j]->setAxisTitle("ES RecHit Occupancy / 5 minutes", 2);
 
      sprintf(histo, "ES Trending RH Occ per hour Z %d P %d", iz, j + 1);
      hESRecHitTrendHr_[i][j] = iBooker.bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
      hESRecHitTrendHr_[i][j]->setAxisTitle("Elapse time (Hours)", 1);
      hESRecHitTrendHr_[i][j]->setAxisTitle("ES RecHit Occupancy / hour", 2);
    }
 
  sprintf(histo, "ES Trending SLink CRC Error per 5 mins");
  hESSLinkErrTrend_ = iBooker.bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
  hESSLinkErrTrend_->setAxisTitle("Elapse time (Minutes)", 1);
  hESSLinkErrTrend_->setAxisTitle("ES SLink CRC Err / 5 minutes", 2);
 
  sprintf(histo, "ES Trending Fiber Error per 5 mins");
  hESFiberErrTrend_ = iBooker.bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
  hESFiberErrTrend_->setAxisTitle("Elapse time (Minutes)", 1);
  hESFiberErrTrend_->setAxisTitle("ES Fiber Err / 5 minutes", 2);
 
  sprintf(histo, "ES Trending SLink CRC Error per hour");
  hESSLinkErrTrendHr_ = iBooker.bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
  hESSLinkErrTrendHr_->setAxisTitle("Elapse time (Hours)", 1);
  hESSLinkErrTrendHr_->setAxisTitle("ES SLink CRC Err / hour", 2);
 
  sprintf(histo, "ES Trending Fiber Error per hour");
  hESFiberErrTrendHr_ = iBooker.bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
  hESFiberErrTrendHr_->setAxisTitle("Elapse time (Hours)", 1);
  hESFiberErrTrendHr_->setAxisTitle("ES Fiber Err / hour", 2);
}
 
void ESTrendTask::analyze(const Event& e, const EventSetup& c) {
  ievt_++;
 
  // Collect time information
  updateTime(e);
 
  long int diff_current_start = current_time_ - start_time_;
  long int diff_last_start = last_time_ - start_time_;
  //LogInfo("ESTrendTask") << "time difference is negative in " << ievt_ << " events\n"
  //<< "\tcurrent - start time = " << diff_current_start
  //<< ", \tlast - start time = " << diff_last_start << endl;
 
  //  std::cout << "current_time : " << current_time_ << ", diff : " << diff_current_start << std::endl;
 
  // Calculate time interval and bin width
  //  int minuteBinWidth = int(nBasicClusterMinutely_->getTProfile()->GetXaxis()->GetBinWidth(1));
  int minuteBinWidth = 5;
  long int minuteBinDiff = diff_current_start / 60 / minuteBinWidth - diff_last_start / 60 / minuteBinWidth;
  long int minuteDiff = (current_time_ - last_time_) / 60;
 
  //  int hourBinWidth = int(nBasicClusterHourly_->getTProfile()->GetXaxis()->GetBinWidth(1));
  int hourBinWidth = 1;
  long int hourBinDiff = diff_current_start / 3600 / hourBinWidth - diff_last_start / 3600 / hourBinWidth;
  long int hourDiff = (current_time_ - last_time_) / 3600;
 
  if (minuteDiff >= minuteBinWidth) {
    while (minuteDiff >= minuteBinWidth)
      minuteDiff -= minuteBinWidth;
  }
  if (hourDiff >= hourBinWidth) {
    while (hourDiff >= hourBinWidth)
      hourDiff -= hourBinWidth;
  }
 
  // ES DCC
  Int_t slinkCRCErr = 0;
  Int_t fiberErr = 0;
  vector<int> fiberStatus;
  Handle<ESRawDataCollection> dccs;
  if (e.getByToken(dccCollections_, dccs)) {
    for (ESRawDataCollection::const_iterator dccItr = dccs->begin(); dccItr != dccs->end(); ++dccItr) {
      ESDCCHeaderBlock dcc = (*dccItr);
 
      if (dcc.getDCCErrors() == 101)
        slinkCRCErr++;
 
      fiberStatus = dcc.getFEChannelStatus();
      for (unsigned int i = 0; i < fiberStatus.size(); ++i) {
        if (fiberStatus[i] == 4 || fiberStatus[i] == 8 || fiberStatus[i] == 10 || fiberStatus[i] == 11 ||
            fiberStatus[i] == 12)
          fiberErr++;
      }
    }
  }
 
  shift2Right(hESSLinkErrTrend_->getTProfile(), minuteBinDiff);
  hESSLinkErrTrend_->Fill(minuteDiff, slinkCRCErr);
 
  shift2Right(hESFiberErrTrend_->getTProfile(), minuteBinDiff);
  hESFiberErrTrend_->Fill(minuteDiff, fiberErr);
 
  shift2Right(hESSLinkErrTrendHr_->getTProfile(), hourBinDiff);
  hESSLinkErrTrendHr_->Fill(hourDiff, slinkCRCErr);
 
  shift2Right(hESFiberErrTrendHr_->getTProfile(), hourBinDiff);
  hESFiberErrTrendHr_->Fill(hourDiff, fiberErr);
 
  // ES RecHits
  int zside, plane;
  int nrh[2][2];
  for (int i = 0; i < 2; i++)
    for (int j = 0; j < 2; j++) {
      nrh[i][j] = 0;
    }
 
  Handle<ESRecHitCollection> ESRecHit;
  if (e.getByToken(rechittoken_, ESRecHit)) {
    for (ESRecHitCollection::const_iterator hitItr = ESRecHit->begin(); hitItr != ESRecHit->end(); ++hitItr) {
      ESDetId id = ESDetId(hitItr->id());
 
      zside = id.zside();
      plane = id.plane();
 
      int i = (zside == 1) ? 0 : 1;
      int j = plane - 1;
 
      nrh[i][j]++;
    }
  } else {
    LogWarning("ESTrendTask") << "RecHitCollection not available";
  }
 
  for (int i = 0; i < 2; ++i)
    for (int j = 0; j < 2; ++j) {
      shift2Right(hESRecHitTrend_[i][j]->getTProfile(), minuteBinDiff);
      hESRecHitTrend_[i][j]->Fill(minuteDiff, nrh[i][j] / (1072 * 32.));
 
      shift2Right(hESRecHitTrendHr_[i][j]->getTProfile(), hourBinDiff);
      hESRecHitTrendHr_[i][j]->Fill(hourDiff, nrh[i][j] / (1072 * 32.));
    }
}
 
void ESTrendTask::updateTime(const edm::Event& e) {
  last_time_ = current_time_;
  current_time_ = e.time().unixTime();
}
 
void ESTrendTask::shift2Right(TProfile* p, int bins) {
  if (bins <= 0)
    return;
 
  if (!p->GetSumw2())
    p->Sumw2();
  int nBins = p->GetXaxis()->GetNbins();
 
  // by shifting n bin to the right, the number of entries are
  // reduced by the number in n bins including the overflow bin.
  double nentries = p->GetEntries();
  for (int i = 0; i < bins; i++)
    nentries -= p->GetBinEntries(nBins + 1 - bins);
  p->SetEntries(nentries);
 
  // the last bin goes to overflow
  // each bin moves to the right
 
  TArrayD* sumw2 = p->GetSumw2();
 
  for (int i = nBins + 1; i > bins; i--) {
    // GetBinContent return binContent/binEntries
    p->SetBinContent(i, p->GetBinContent(i - bins) * p->GetBinEntries(i - bins));
    p->SetBinEntries(i, p->GetBinEntries(i - bins));
    sumw2->SetAt(sumw2->GetAt(i - bins), i);
  }
}
 
void ESTrendTask::shift2Left(TProfile* p, int bins) {
  if (bins <= 0)
    return;
 
  if (!p->GetSumw2())
    p->Sumw2();
  int nBins = p->GetXaxis()->GetNbins();
 
  // by shifting n bin to the left, the number of entries are
  // reduced by the number in n bins including the underflow bin.
  double nentries = p->GetEntries();
  for (int i = 0; i < bins; i++)
    nentries -= p->GetBinEntries(i);
  p->SetEntries(nentries);
 
  // the first bin goes to underflow
  // each bin moves to the right
 
  TArrayD* sumw2 = p->GetSumw2();
 
  for (int i = 0; i <= nBins + 1 - bins; i++) {
    // GetBinContent return binContent/binEntries
    p->SetBinContent(i, p->GetBinContent(i + bins) * p->GetBinEntries(i + bins));
    p->SetBinEntries(i, p->GetBinEntries(i + bins));
    sumw2->SetAt(sumw2->GetAt(i + bins), i);
  }
}
 
DEFINE_FWK_MODULE(ESTrendTask);