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, ix, iy, strip;
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       ix    = id.six();
00265       iy    = id.siy();
00266       strip = id.strip();
00267       
00268       int i = (zside==1)? 0:1;
00269       int j = plane-1;
00270       
00271       nrh[i][j]++;
00272     }
00273   } else {
00274     LogWarning("ESTrendTask") << rechitlabel_ << " not available";
00275   }
00276   
00277   for (int i=0; i<2; ++i) 
00278     for (int j=0; j<2; ++j) {
00279       shift2Right(hESRecHitTrend_[i][j]->getTProfile(), minuteBinDiff);
00280       hESRecHitTrend_[i][j]->Fill(minuteDiff, nrh[i][j]/(1072*32.));
00281 
00282       shift2Right(hESRecHitTrendHr_[i][j]->getTProfile(), hourBinDiff);
00283       hESRecHitTrendHr_[i][j]->Fill(hourDiff, nrh[i][j]/(1072*32.));
00284     }
00285 
00286 }
00287 
00288 void ESTrendTask::updateTime(const edm::Event& e) {
00289 
00290   last_time_ = current_time_;
00291   current_time_ = e.time().unixTime();
00292 
00293 }
00294 
00295 void ESTrendTask::shift2Right(TProfile* p, int bins) {
00296 
00297   if(bins <= 0) return;
00298 
00299   if(!p->GetSumw2()) p->Sumw2();
00300   int nBins = p->GetXaxis()->GetNbins();
00301 
00302   // by shifting n bin to the right, the number of entries are
00303   // reduced by the number in n bins including the overflow bin.
00304   double nentries = p->GetEntries();
00305   for(int i=0; i<bins; i++) nentries -= p->GetBinEntries(nBins+1-bins);
00306   p->SetEntries(nentries);
00307   
00308   // the last bin goes to overflow
00309   // each bin moves to the right
00310 
00311   TArrayD* sumw2 = p->GetSumw2();
00312 
00313   for(int i=nBins+1; i>bins; i--) {
00314     // GetBinContent return binContent/binEntries
00315     p->SetBinContent(i, p->GetBinContent(i-bins)*p->GetBinEntries(i-bins));
00316     p->SetBinEntries(i,p->GetBinEntries(i-bins));
00317     sumw2->SetAt(sumw2->GetAt(i-bins),i);
00318   }
00319 
00320 }
00321 
00322 void ESTrendTask::shift2Left(TProfile* p, int bins) {
00323 
00324   if(bins <= 0) return;
00325 
00326   if(!p->GetSumw2()) p->Sumw2();
00327   int nBins = p->GetXaxis()->GetNbins();
00328 
00329   // by shifting n bin to the left, the number of entries are
00330   // reduced by the number in n bins including the underflow bin.
00331   double nentries = p->GetEntries();
00332   for(int i=0; i<bins; i++) nentries -= p->GetBinEntries(i);
00333   p->SetEntries(nentries);
00334   
00335   // the first bin goes to underflow
00336   // each bin moves to the right
00337 
00338   TArrayD* sumw2 = p->GetSumw2();
00339 
00340   for(int i=0; i<=nBins+1-bins; i++) {
00341     // GetBinContent return binContent/binEntries
00342     p->SetBinContent(i, p->GetBinContent(i+bins)*p->GetBinEntries(i+bins));
00343     p->SetBinEntries(i,p->GetBinEntries(i+bins));
00344     sumw2->SetAt(sumw2->GetAt(i+bins),i);
00345   }
00346 
00347 }
00348 
00349 DEFINE_FWK_MODULE(ESTrendTask);
CMSSW_4_4_3_patch1/src/DQM/EcalPreshowerMonitorModule/src/ESTrendTask.cc
