Reference Manual

Go to the documentation of this file.
00001 #include <iostream>
00002 #include <fstream>
00003 #include <vector>
00004 #include <math.h>
00005 
00006 #include "FWCore/Framework/interface/Run.h"
00007 #include "FWCore/ServiceRegistry/interface/Service.h"
00008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00009 #include "FWCore/Framework/interface/ESHandle.h"
00010 #include "FWCore/Framework/interface/MakerMacros.h"
00011 #include "FWCore/Framework/interface/Frameworkfwd.h"
00012 #include "DataFormats/EcalDetId/interface/ESDetId.h"
00013 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00014 #include "DQMServices/Core/interface/MonitorElement.h"
00015 #include "DQMServices/Core/interface/DQMStore.h"
00016 #include "DQM/EcalPreshowerMonitorModule/interface/ESTrendTask.h"
00017 #include "DataFormats/EcalRawData/interface/ESDCCHeaderBlock.h"
00018 #include "DataFormats/EcalRawData/interface/ESKCHIPBlock.h"
00019 #include "DataFormats/EcalRawData/interface/EcalRawDataCollections.h"
00020 
00021 using namespace cms;
00022 using namespace edm;
00023 using namespace std;
00024 
00025 ESTrendTask::ESTrendTask(const ParameterSet& ps) {
00026 
00027   init_ = false;
00028 
00029   dqmStore_ = Service<DQMStore>().operator->();
00030 
00031   prefixME_       = ps.getUntrackedParameter<string>("prefixME", "EcalPreshower"); 
00032   enableCleanup_  = ps.getUntrackedParameter<bool>("enableCleanup", false);
00033   mergeRuns_      = ps.getUntrackedParameter<bool>("mergeRuns", false);
00034   rechitlabel_    = ps.getParameter<InputTag>("RecHitLabel");
00035   dccCollections_ = ps.getParameter<InputTag>("ESDCCCollections");
00036 
00037   for (int i=0; i<2; ++i)
00038     for (int j=0; j<2; ++j) {
00039       hESRecHitTrend_[i][j] = 0;
00040       hESRecHitTrendHr_[i][j] = 0;
00041     }
00042 
00043   hESSLinkErrTrend_ = 0;
00044   hESFiberErrTrend_ = 0;
00045   hESSLinkErrTrendHr_ = 0;
00046   hESFiberErrTrendHr_ = 0;
00047 
00048   start_time_ = 0;
00049   current_time_ = 0;
00050   last_time_ = 0;
00051 
00052 }
00053 
00054 ESTrendTask::~ESTrendTask() {
00055 }
00056 
00057 void ESTrendTask::beginJob(void) {
00058 
00059   ievt_ = 0;
00060 
00061   if ( dqmStore_ ) {
00062     dqmStore_->setCurrentFolder(prefixME_ + "/ESTrendTask");
00063     dqmStore_->rmdir(prefixME_ + "/ESTrendTask");
00064   }
00065 
00066 }
00067 
00068 void ESTrendTask::beginRun(const Run& r, const EventSetup& c) {
00069 
00070   if ( ! mergeRuns_ ) this->reset();
00071 
00072   start_time_ = (r.beginTime()).unixTime();
00073 
00074   //std::cout << "start time : " << start_time_ << std::endl;
00075 
00076 }
00077 
00078 void ESTrendTask::endRun(const Run& r, const EventSetup& c) {
00079 
00080 }
00081 
00082 void ESTrendTask::reset(void) {
00083 
00084   for (int i=0 ; i<2; ++i)  
00085     for (int j=0 ; j<2; ++j) {
00086       if (hESRecHitTrend_[i][j]) hESRecHitTrend_[i][j]->Reset();
00087       if (hESRecHitTrendHr_[i][j]) hESRecHitTrendHr_[i][j]->Reset();
00088     }
00089 
00090   if (hESSLinkErrTrend_) hESSLinkErrTrend_->Reset();
00091   if (hESFiberErrTrend_) hESFiberErrTrend_->Reset();
00092   if (hESSLinkErrTrendHr_) hESSLinkErrTrendHr_->Reset();
00093   if (hESFiberErrTrendHr_) hESFiberErrTrendHr_->Reset();
00094 
00095 }
00096 
00097 void ESTrendTask::setup(void) {
00098 
00099   init_ = true;
00100 
00101   char histo[200];
00102 
00103   if ( dqmStore_ ) {
00104     dqmStore_->setCurrentFolder(prefixME_ + "/ESTrendTask");
00105 
00106     for (int i=0 ; i<2; ++i)  
00107       for (int j=0 ; j<2; ++j) {
00108         int iz = (i==0)? 1:-1;
00109         sprintf(histo, "ES Trending RH Occ per 5 mins Z %d P %d", iz, j+1);
00110         hESRecHitTrend_[i][j] = dqmStore_->bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
00111         hESRecHitTrend_[i][j]->setAxisTitle("Elapse time (Minutes)", 1);
00112         hESRecHitTrend_[i][j]->setAxisTitle("ES RecHit Occupancy / 5 minutes", 2);
00113 
00114         sprintf(histo, "ES Trending RH Occ per hour Z %d P %d", iz, j+1);
00115         hESRecHitTrendHr_[i][j] = dqmStore_->bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
00116         hESRecHitTrendHr_[i][j]->setAxisTitle("Elapse time (Hours)", 1);
00117         hESRecHitTrendHr_[i][j]->setAxisTitle("ES RecHit Occupancy / hour", 2);
00118       }
00119 
00120     sprintf(histo, "ES Trending SLink CRC Error per 5 mins");
00121     hESSLinkErrTrend_ = dqmStore_->bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
00122     hESSLinkErrTrend_->setAxisTitle("Elapse time (Minutes)", 1);
00123     hESSLinkErrTrend_->setAxisTitle("ES SLink CRC Err / 5 minutes", 2);
00124 
00125     sprintf(histo, "ES Trending Fiber Error per 5 mins");
00126     hESFiberErrTrend_ = dqmStore_->bookProfile(histo, histo, 36, 0.0, 180.0, 100, 0.0, 1.0e6, "s");
00127     hESFiberErrTrend_->setAxisTitle("Elapse time (Minutes)", 1);
00128     hESFiberErrTrend_->setAxisTitle("ES Fiber Err / 5 minutes", 2);
00129 
00130     sprintf(histo, "ES Trending SLink CRC Error per hour");
00131     hESSLinkErrTrendHr_ = dqmStore_->bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
00132     hESSLinkErrTrendHr_->setAxisTitle("Elapse time (Hours)", 1);
00133     hESSLinkErrTrendHr_->setAxisTitle("ES SLink CRC Err / hour", 2);
00134 
00135     sprintf(histo, "ES Trending Fiber Error per hour");
00136     hESFiberErrTrendHr_ = dqmStore_->bookProfile(histo, histo, 24, 0.0, 24.0, 100, 0.0, 1.0e6, "s");
00137     hESFiberErrTrendHr_->setAxisTitle("Elapse time (Hours)", 1);
00138     hESFiberErrTrendHr_->setAxisTitle("ES Fiber Err / hour", 2);
00139   }
00140 
00141 }
00142 
00143 void ESTrendTask::cleanup(void) {
00144 
00145   if ( ! init_ ) return;
00146 
00147   if ( dqmStore_ ) {
00148     dqmStore_->setCurrentFolder(prefixME_ + "/ESTrendTask");
00149 
00150     for (int i=0 ; i<2; ++i)  
00151       for (int j=0 ; j<2; ++j) {
00152         if (hESRecHitTrend_[i][j]) dqmStore_->removeElement(hESRecHitTrend_[i][j]->getName());
00153         hESRecHitTrend_[i][j] = 0;
00154         if (hESRecHitTrendHr_[i][j]) dqmStore_->removeElement(hESRecHitTrendHr_[i][j]->getName());
00155         hESRecHitTrendHr_[i][j] = 0;
00156       }
00157 
00158     if (hESSLinkErrTrend_) dqmStore_->removeElement(hESSLinkErrTrend_->getName());
00159     hESSLinkErrTrend_ = 0;
00160     if (hESFiberErrTrend_) dqmStore_->removeElement(hESFiberErrTrend_->getName());
00161     hESFiberErrTrend_ = 0;
00162     if (hESSLinkErrTrendHr_) dqmStore_->removeElement(hESSLinkErrTrendHr_->getName());
00163     hESSLinkErrTrendHr_ = 0;
00164     if (hESFiberErrTrendHr_) dqmStore_->removeElement(hESFiberErrTrendHr_->getName());
00165     hESFiberErrTrendHr_ = 0;
00166   }
00167 
00168   init_ = false;
00169 
00170 }
00171 
00172 void ESTrendTask::endJob(void) {
00173 
00174   LogInfo("ESTrendTask") << "analyzed " << ievt_ << " events";
00175 
00176   if ( enableCleanup_ ) this->cleanup();
00177 
00178 }
00179 
00180 void ESTrendTask::analyze(const Event& e, const EventSetup& c) {
00181 
00182   if ( ! init_ ) this->setup();
00183 
00184   ievt_++;
00185 
00186   // Collect time information
00187   updateTime(e);
00188 
00189   long int diff_current_start = current_time_ - start_time_;
00190   long int diff_last_start    = last_time_ - start_time_;
00191   //LogInfo("ESTrendTask") << "time difference is negative in " << ievt_ << " events\n"
00192   //<< "\tcurrent - start time = " << diff_current_start
00193   //<< ", \tlast - start time = " << diff_last_start << endl;
00194 
00195   //  std::cout << "current_time : " << current_time_ << ", diff : " << diff_current_start << std::endl;
00196 
00197   // Calculate time interval and bin width
00198   //  int minuteBinWidth = int(nBasicClusterMinutely_->getTProfile()->GetXaxis()->GetBinWidth(1));
00199   int minuteBinWidth = 5;
00200   long int minuteBinDiff = diff_current_start/60/minuteBinWidth - diff_last_start/60/minuteBinWidth;
00201   long int minuteDiff = (current_time_ - last_time_)/60;
00202 
00203   //  int hourBinWidth = int(nBasicClusterHourly_->getTProfile()->GetXaxis()->GetBinWidth(1));
00204   int hourBinWidth = 1;
00205   long int hourBinDiff = diff_current_start/3600/hourBinWidth - diff_last_start/3600/hourBinWidth;
00206   long int hourDiff = (current_time_ - last_time_)/3600;
00207 
00208   if (minuteDiff >= minuteBinWidth) {
00209     while (minuteDiff >= minuteBinWidth) minuteDiff -= minuteBinWidth;
00210   }
00211   if (hourDiff >= hourBinWidth) {
00212     while (hourDiff >= hourBinWidth) hourDiff -= hourBinWidth;
00213   }
00214 
00215   // ES DCC
00216   Int_t slinkCRCErr = 0;
00217   Int_t fiberErr = 0;
00218   vector<int> fiberStatus;
00219   Handle<ESRawDataCollection> dccs;
00220   if ( e.getByLabel(dccCollections_, dccs) ) {
00221     for (ESRawDataCollection::const_iterator dccItr = dccs->begin(); dccItr != dccs->end(); ++dccItr) {
00222       ESDCCHeaderBlock dcc = (*dccItr);
00223 
00224       if (dcc.getDCCErrors() == 101) slinkCRCErr++;
00225 
00226       fiberStatus = dcc.getFEChannelStatus();
00227       for (unsigned int i=0; i<fiberStatus.size(); ++i) {
00228         if (fiberStatus[i]==4 || fiberStatus[i]==8 || fiberStatus[i]==10 || fiberStatus[i]==11 || fiberStatus[i]==12)
00229           fiberErr++;
00230       }
00231       
00232     }
00233   }
00234 
00235   shift2Right(hESSLinkErrTrend_->getTProfile(), minuteBinDiff);
00236   hESSLinkErrTrend_->Fill(minuteDiff, slinkCRCErr);
00237 
00238   shift2Right(hESFiberErrTrend_->getTProfile(), minuteBinDiff);
00239   hESFiberErrTrend_->Fill(minuteDiff, fiberErr);
00240 
00241   shift2Right(hESSLinkErrTrendHr_->getTProfile(), hourBinDiff);
00242   hESSLinkErrTrendHr_->Fill(hourDiff, slinkCRCErr);
00243 
00244   shift2Right(hESFiberErrTrendHr_->getTProfile(), hourBinDiff);
00245   hESFiberErrTrendHr_->Fill(hourDiff, fiberErr);
00246 
00247   // ES RecHits
00248   int zside, plane;
00249   int nrh[2][2];
00250   for (int i = 0; i < 2; i++ ) 
00251     for( int j = 0; j < 2; j++) {
00252       nrh[i][j] = 0;
00253     }
00254 
00255   Handle<EcalRecHitCollection> ESRecHit;
00256   if ( e.getByLabel(rechitlabel_, ESRecHit) ) {
00257     
00258     for (ESRecHitCollection::const_iterator hitItr = ESRecHit->begin(); hitItr != ESRecHit->end(); ++hitItr) {
00259       
00260       ESDetId id = ESDetId(hitItr->id());
00261       
00262       zside = id.zside();
00263       plane = id.plane();
00264       
00265       int i = (zside==1)? 0:1;
00266       int j = plane-1;
00267       
00268       nrh[i][j]++;
00269     }
00270   } else {
00271     LogWarning("ESTrendTask") << rechitlabel_ << " not available";
00272   }
00273   
00274   for (int i=0; i<2; ++i) 
00275     for (int j=0; j<2; ++j) {
00276       shift2Right(hESRecHitTrend_[i][j]->getTProfile(), minuteBinDiff);
00277       hESRecHitTrend_[i][j]->Fill(minuteDiff, nrh[i][j]/(1072*32.));
00278 
00279       shift2Right(hESRecHitTrendHr_[i][j]->getTProfile(), hourBinDiff);
00280       hESRecHitTrendHr_[i][j]->Fill(hourDiff, nrh[i][j]/(1072*32.));
00281     }
00282 
00283 }
00284 
00285 void ESTrendTask::updateTime(const edm::Event& e) {
00286 
00287   last_time_ = current_time_;
00288   current_time_ = e.time().unixTime();
00289 
00290 }
00291 
00292 void ESTrendTask::shift2Right(TProfile* p, int bins) {
00293 
00294   if(bins <= 0) return;
00295 
00296   if(!p->GetSumw2()) p->Sumw2();
00297   int nBins = p->GetXaxis()->GetNbins();
00298 
00299   // by shifting n bin to the right, the number of entries are
00300   // reduced by the number in n bins including the overflow bin.
00301   double nentries = p->GetEntries();
00302   for(int i=0; i<bins; i++) nentries -= p->GetBinEntries(nBins+1-bins);
00303   p->SetEntries(nentries);
00304   
00305   // the last bin goes to overflow
00306   // each bin moves to the right
00307 
00308   TArrayD* sumw2 = p->GetSumw2();
00309 
00310   for(int i=nBins+1; i>bins; i--) {
00311     // GetBinContent return binContent/binEntries
00312     p->SetBinContent(i, p->GetBinContent(i-bins)*p->GetBinEntries(i-bins));
00313     p->SetBinEntries(i,p->GetBinEntries(i-bins));
00314     sumw2->SetAt(sumw2->GetAt(i-bins),i);
00315   }
00316 
00317 }
00318 
00319 void ESTrendTask::shift2Left(TProfile* p, int bins) {
00320 
00321   if(bins <= 0) return;
00322 
00323   if(!p->GetSumw2()) p->Sumw2();
00324   int nBins = p->GetXaxis()->GetNbins();
00325 
00326   // by shifting n bin to the left, the number of entries are
00327   // reduced by the number in n bins including the underflow bin.
00328   double nentries = p->GetEntries();
00329   for(int i=0; i<bins; i++) nentries -= p->GetBinEntries(i);
00330   p->SetEntries(nentries);
00331   
00332   // the first bin goes to underflow
00333   // each bin moves to the right
00334 
00335   TArrayD* sumw2 = p->GetSumw2();
00336 
00337   for(int i=0; i<=nBins+1-bins; i++) {
00338     // GetBinContent return binContent/binEntries
00339     p->SetBinContent(i, p->GetBinContent(i+bins)*p->GetBinEntries(i+bins));
00340     p->SetBinEntries(i,p->GetBinEntries(i+bins));
00341     sumw2->SetAt(sumw2->GetAt(i+bins),i);
00342   }
00343 
00344 }
00345 
00346 DEFINE_FWK_MODULE(ESTrendTask);
/data/refman/pasoursint/CMSSW_5_3_10_patch2/src/DQM/EcalPreshowerMonitorModule/src/ESTrendTask.cc
