CMSSW Reference Manual

00001 // -*- C++ -*-
00002 //
00003 // Package:    DQMO/HcalMonitorClient/HcalDataCertification
00004 // Class:      HcalDataCertification
00005 // 
00013 //
00014 // Original Author:  "Igor Vodopiyanov"
00015 //         Created:  Nov-21 2008
00016 // $Id: HcalDataCertification.cc,v 1.1 2008/12/18 21:45:45 ivodop Exp $
00017 //
00018 //
00019 
00020 // system include files
00021 #include <memory>
00022 #include <stdio.h>
00023 #include <math.h>
00024 #include <sstream>
00025 #include <fstream>
00026 #include <exception>
00027 
00028 // user include files
00029 #include "FWCore/Framework/interface/Frameworkfwd.h"
00030 #include "FWCore/Framework/interface/EDAnalyzer.h"
00031 
00032 #include "FWCore/Framework/interface/Event.h"
00033 #include "FWCore/Framework/interface/MakerMacros.h"
00034 
00035 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00036 
00037 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00038 #include "DQMServices/Core/interface/DQMStore.h"
00039 #include "DQMServices/Core/interface/MonitorElement.h"
00040 #include "FWCore/ServiceRegistry/interface/Service.h"
00041 #include "PhysicsTools/UtilAlgos/interface/TFileService.h"
00042 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00043 
00044 //
00045 // class decleration
00046 //
00047 
00048 class HcalDataCertification : public edm::EDAnalyzer {
00049    public:
00050       explicit HcalDataCertification(const edm::ParameterSet&);
00051       ~HcalDataCertification();
00052 
00053    private:
00054       virtual void beginJob(const edm::EventSetup&) ;
00055       virtual void analyze(const edm::Event&, const edm::EventSetup&);
00056       virtual void endJob() ;
00057       virtual void beginRun(const edm::Run&, const edm::EventSetup&) ;
00058       virtual void endRun(const edm::Run&, const edm::EventSetup&) ;
00059 
00060    // ----------member data ---------------------------
00061 
00062    edm::ParameterSet conf_;
00063    DQMStore * dbe;
00064    edm::Service<TFileService> fs_;
00065 
00066 };
00067 
00068 //
00069 // constants, enums and typedefs
00070 //
00071 
00072 //
00073 // static data member definitions
00074 //
00075 
00076 //
00077 // constructors and destructor
00078 //
00079 
00080 HcalDataCertification::HcalDataCertification(const edm::ParameterSet& iConfig):conf_(iConfig)
00081 {
00082   // now do what ever initialization is needed
00083 }
00084 
00085 
00086 HcalDataCertification::~HcalDataCertification()
00087 { 
00088    // do anything here that needs to be done at desctruction time
00089    // (e.g. close files, deallocate resources etc.)
00090 }
00091 
00092 
00093 //
00094 // member functions
00095 //
00096 
00097 // ------------ method called to for each event  ------------
00098 void
00099 HcalDataCertification::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00100 {
00101   using namespace edm;
00102    
00103 #ifdef THIS_IS_AN_EVENT_EXAMPLE
00104   Handle<ExampleData> pIn;
00105   iEvent.getByLabel("example",pIn);
00106 #endif
00107    
00108 #ifdef THIS_IS_AN_EVENTSETUP_EXAMPLE
00109   ESHandle<SetupData> pSetup;
00110   iSetup.get<SetupRecord>().get(pSetup);
00111 #endif
00112 
00113 }
00114 
00115 // ------------ method called once each job just before starting event loop  ------------
00116 void 
00117 HcalDataCertification::beginJob(const edm::EventSetup&)
00118 {
00119   dbe = 0;
00120   dbe = edm::Service<DQMStore>().operator->();
00121   //std::cout<<"beginJob"<< std::endl;
00122 }
00123 
00124 // ------------ method called once each job just after ending the event loop  ------------
00125 void 
00126 HcalDataCertification::endJob() 
00127 {
00128   //std::cout << ">>> endJob " << std::endl;
00129 }
00130 
00131 // ------------ method called just before starting a new run  ------------
00132 void 
00133 HcalDataCertification::beginRun(const edm::Run& run, const edm::EventSetup& c)
00134 {
00135   //std::cout<<"beginRun"<<std::endl;
00136 }
00137 
00138 // ------------ method called right after a run ends ------------
00139 void 
00140 HcalDataCertification::endRun(const edm::Run& run, const edm::EventSetup& c)
00141 {
00142   
00143   dbe->setCurrentFolder("Hcal");
00144   std::string currDir = dbe->pwd();
00145   //std::cout << "--- Current Directory " << currDir << std::endl;
00146   std::vector<MonitorElement*> mes = dbe->getAllContents("");
00147   //std::cout << "found " << mes.size() << " monitoring elements:" << std::endl;
00148 
00149   dbe->setCurrentFolder("Hcal/EventInfo/Certification/");
00150   MonitorElement* Hcal_HB = dbe->bookFloat("Hcal_HB");
00151   MonitorElement* Hcal_HE = dbe->bookFloat("Hcal_HE");
00152   MonitorElement* Hcal_HF = dbe->bookFloat("Hcal_HF");
00153   MonitorElement* Hcal_HO = dbe->bookFloat("Hcal_HO");
00154 
00155   Hcal_HB->Fill(0);
00156   Hcal_HE->Fill(0);
00157   Hcal_HF->Fill(0);
00158   Hcal_HO->Fill(0);
00159 
00160   Int_t ndeadHB,ndeadHE,ndeadHF,ndeadHO,nhotHB,nhotHE,nhotHF,nhotHO;
00161   ndeadHB=ndeadHE=ndeadHF=ndeadHO=nhotHB=nhotHE=nhotHF=nhotHO=0;
00162   Int_t nenHB,nenHE,nenHF,nenHO;
00163   nenHB=nenHE=nenHF=nenHO=0;
00164 
00165   TH1F *hTP;
00166   TH2F *hhotHB,*hhotHF,*hhotHE,*hhotHO;
00167   TH2F *hdeadHB,*hdeadHF,*hdeadHE,*hdeadHO;
00168   TH2F *henergeticHB,*henergeticHE,*henergeticHF,*henergeticHO;
00169 
00170 
00171   // ************** Loop over Monitoring Elements and fill working histograms
00172   Int_t nWorkHist = 0;
00173   for(std::vector<MonitorElement*>::const_iterator ime = mes.begin(); ime!=mes.end(); ++ime) {
00174     std::string name = (*ime)->getName();
00175 
00176     // Hcal Nevents by Trig
00177 
00178     if (name == "# TP Digis") {
00179       hTP = (*ime)->getTH1F();
00180       //std::cout <<"TP Digis"<<std::endl;
00181       nWorkHist++;
00182     }
00183 
00184     // by SubDet
00185 
00186     if (name == "HB_OccupancyMap_HotCell_Threshold0") { // -----------Hot
00187       hhotHB = (*ime)->getTH2F();
00188       //std::cout <<"HB_OccupancyMap_HotCell_Threshold0"<<std::endl;
00189       nWorkHist++;
00190     }
00191     if (name == "HE_OccupancyMap_HotCell_Threshold0") { 
00192       hhotHE = (*ime)->getTH2F();
00193       //std::cout <<"HE_OccupancyMap_HotCell_Threshold0"<<std::endl;
00194       nWorkHist++;
00195     }
00196     if (name == "HF_OccupancyMap_HotCell_Threshold0") { 
00197       hhotHF = (*ime)->getTH2F();
00198       //std::cout <<"HF_OccupancyMap_HotCell_Threshold0"<<std::endl;
00199       nWorkHist++;
00200     }
00201     if (name == "HO_OccupancyMap_HotCell_Threshold0") { 
00202       hhotHO = (*ime)->getTH2F();
00203       //std::cout <<"HO_OccupancyMap_HotCell_Threshold0"<<std::endl;
00204       nWorkHist++;
00205     }
00206 
00207     if (name == "HB_HotCell_EnergyMap_Thresh0") { // -----------Energetic
00208       henergeticHB = (*ime)->getTH2F();
00209       //std::cout <<"HB_HotCell_EnergyMap_Thresh0"<<std::endl;
00210       nWorkHist++;
00211     }
00212     if (name == "HE_HotCell_EnergyMap_Thresh0") { 
00213       henergeticHE = (*ime)->getTH2F();
00214       //std::cout <<"HE_HotCell_EnergyMap_Thresh0"<<std::endl;
00215       nWorkHist++;
00216     }
00217     if (name == "HF_HotCell_EnergyMap_Thresh0") { 
00218       henergeticHF = (*ime)->getTH2F();
00219       //std::cout <<"HF_HotCell_EnergyMap_Thresh0"<<std::endl;
00220       nWorkHist++;
00221     }
00222     if (name == "HO_HotCell_EnergyMap_Thresh0") { 
00223       henergeticHO = (*ime)->getTH2F();
00224       //std::cout <<"HO_HotCell_EnergyMap_Thresh0"<<std::endl;
00225       nWorkHist++;
00226     }
00227 
00228     if (name == "HBProblemDeadCells") { // --------------------------Dead
00229       hdeadHB = (*ime)->getTH2F();
00230       //std::cout <<"HBProblemDeadCells"<<std::endl;
00231       nWorkHist++;
00232     }
00233     if (name == "HEProblemDeadCells") { 
00234       hdeadHE = (*ime)->getTH2F();
00235       //std::cout <<"HEProblemDeadCells"<<std::endl;
00236       nWorkHist++;
00237     }
00238     if (name == "HFProblemDeadCells") { 
00239       hdeadHF = (*ime)->getTH2F();
00240       //std::cout <<"HFProblemDeadCells"<<std::endl;
00241       nWorkHist++;
00242     }
00243     if (name == "HOProblemDeadCells") { 
00244       hdeadHO = (*ime)->getTH2F();
00245       //std::cout <<"HOProblemDeadCells"<<std::endl;
00246       nWorkHist++;
00247     }
00248 
00249   }  // ******** End loop over Monitoring Elements ------------------
00250 
00251 
00252   if (nWorkHist<13) {
00253     edm::LogPrint("HcalDataCertification")<<"N Hist Found ="<<nWorkHist<<" out of 13 => return"<<std::endl;
00254     return;
00255   }
00256 
00257   Int_t Nevents = (int) hTP->GetEntries();   // ----- Nevents
00258   if (Nevents<1) {
00259     edm::LogPrint("HcalDataCertification")<<"N events ="<<Nevents<<" => return"<<std::endl;
00260     return;
00261   }
00262 
00263   TH2F *hdeadratHB = (TH2F*) hdeadHB->Clone(); // ----- dead rate
00264   TH2F *hdeadratHE = (TH2F*) hdeadHE->Clone(); 
00265   TH2F *hdeadratHF = (TH2F*) hdeadHF->Clone(); 
00266   TH2F *hdeadratHO = (TH2F*) hdeadHO->Clone(); 
00267 
00268   TH2F *hhotratHB = (TH2F*) hhotHB->Clone(); // ----- hot rate
00269   TH2F *hhotratHE = (TH2F*) hhotHE->Clone(); 
00270   TH2F *hhotratHF = (TH2F*) hhotHF->Clone(); 
00271   TH2F *hhotratHO = (TH2F*) hhotHO->Clone(); 
00272 
00273   TH2F *henergeticcellHB = (TH2F*) henergeticHB->Clone(); // ----energetic rate
00274   TH2F *henergeticcellHE = (TH2F*) henergeticHE->Clone(); 
00275   TH2F *henergeticcellHF = (TH2F*) henergeticHF->Clone(); 
00276   TH2F *henergeticcellHO = (TH2F*) henergeticHO->Clone(); 
00277 
00278   hdeadratHB->Reset();hhotratHB->Reset();henergeticcellHB->Reset();
00279   hdeadratHE->Reset();hhotratHE->Reset();henergeticcellHE->Reset();
00280   hdeadratHF->Reset();hhotratHF->Reset();henergeticcellHF->Reset();
00281   hdeadratHO->Reset();hhotratHO->Reset();henergeticcellHO->Reset();
00282 
00283   // energetic cells
00284   for (int ii=1;ii<=henergeticHB->GetNbinsX();ii++) for (int jj=1;jj<=henergeticHB->GetNbinsY();jj++) { 
00285     if (henergeticHB->GetBinContent(ii,jj)>100) {
00286       double encell = henergeticHB->GetBinContent(ii,jj);
00287       double ensumcells = 0;
00288       int numcells = 0;
00289       int lmin=jj-2;
00290       int lmax=jj+2;
00291       for (int kk=TMath::Max(1,ii-2);kk<=TMath::Min(ii+2,henergeticHB->GetNbinsX());kk++) {
00292         for (int ll=TMath::Max(1,lmin);ll<=TMath::Min(lmax,henergeticHB->GetNbinsY());ll++) {
00293           if (henergeticHB->GetBinContent(kk,ll)>0 && (kk!=ii || ll!=jj)) {
00294             ensumcells += henergeticHB->GetBinContent(kk,ll);
00295             numcells++;
00296           }
00297         }
00298         if (numcells>0) if (1000*encell/(ensumcells/numcells+encell)>950) {
00299           henergeticcellHB->SetBinContent(ii,jj,1);
00300           //cout<<"energetic HB: "<<ii<<" / "<<jj<<"  E="<<encell<<"  "<<ensumcells<<"  "<<numcells<<std::endl;
00301         }
00302       }
00303     }
00304   }
00305   for (int ii=1;ii<=henergeticHE->GetNbinsX();ii++) for (int jj=1;jj<=henergeticHE->GetNbinsY();jj++) { // 
00306     if (henergeticHE->GetBinContent(ii,jj)>100) {
00307       double encell = henergeticHE->GetBinContent(ii,jj);
00308       double ensumcells = 0;
00309       int numcells = 0;
00310       int lmin=jj-2;
00311       int lmax=jj+2;
00312       if (ii<23 || ii>63) {lmin=jj-4; lmax=jj+4;}
00313       for (int kk=TMath::Max(1,ii-2);kk<=TMath::Min(ii+2,henergeticHE->GetNbinsX());kk++) {
00314         for (int ll=TMath::Max(1,lmin);ll<=TMath::Min(lmax,henergeticHE->GetNbinsY());ll++) {
00315           if (henergeticHE->GetBinContent(kk,ll)>0 && (kk!=ii || ll!=jj)) {
00316             ensumcells += henergeticHE->GetBinContent(kk,ll);
00317             numcells++;
00318           }
00319         }
00320       }
00321       if (numcells>0) if (1000*encell/(ensumcells/numcells+encell)>950) {
00322         henergeticcellHE->SetBinContent(ii,jj,1);
00323         //cout<<"energetic HE: "<<ii<<" / "<<jj<<"  E="<<encell<<"  "<<ensumcells<<"  "<<numcells<<std::endl;
00324       }
00325     }
00326   }
00327   for (int ii=1;ii<=henergeticHF->GetNbinsX();ii++) for (int jj=1;jj<=henergeticHF->GetNbinsY();jj++) { // 
00328     if (henergeticHF->GetBinContent(ii,jj)>100) {
00329       double encell = henergeticHF->GetBinContent(ii,jj);
00330       double ensumcells = 0;
00331       int numcells = 0;
00332       int lmin=jj-4;
00333       int lmax=jj+4;
00334       if (ii<3 || ii>83) {lmin=jj-8; lmax=jj+8;}
00335       for (int kk=TMath::Max(1,ii-2);kk<=TMath::Min(ii+2,henergeticHF->GetNbinsX());kk++) {
00336         for (int ll=TMath::Max(1,lmin);ll<=TMath::Min(lmax,henergeticHF->GetNbinsY());ll++) {
00337           if (henergeticHF->GetBinContent(kk,ll)>0 && (kk!=ii || ll!=jj)) {
00338             ensumcells += henergeticHF->GetBinContent(kk,ll);
00339             numcells++;
00340           }
00341         }
00342       }
00343       if (numcells>0) if (1000*encell/(ensumcells/numcells+encell)>950) {
00344         henergeticcellHF->SetBinContent(ii,jj,1);
00345         //cout<<"energetic HF: "<<ii<<" / "<<jj<<"  E="<<encell<<"  "<<ensumcells<<"  "<<numcells<<std::endl;
00346       }
00347     }
00348   }
00349   for (int ii=1;ii<=henergeticHO->GetNbinsX();ii++) for (int jj=1;jj<=henergeticHO->GetNbinsY();jj++) { // 
00350     if (henergeticHO->GetBinContent(ii,jj)>100) {
00351       double encell = henergeticHO->GetBinContent(ii,jj);
00352       double ensumcells = 0;
00353       int numcells = 0;
00354       int lmin=jj-2;
00355       int lmax=jj+2;
00356       for (int kk=TMath::Max(1,ii-2);kk<=TMath::Min(ii+2,henergeticHO->GetNbinsX());kk++) {
00357         for (int ll=TMath::Max(1,lmin);ll<=TMath::Min(lmax,henergeticHO->GetNbinsY());ll++) {
00358           if (henergeticHO->GetBinContent(kk,ll)>0 && (kk!=ii || ll!=jj)) {
00359             ensumcells += henergeticHO->GetBinContent(kk,ll);
00360             numcells++;
00361           }
00362         }
00363       }
00364       if (numcells>0) if (1000*encell/(ensumcells/numcells+encell)>950) {
00365         henergeticcellHO->SetBinContent(ii,jj,1);
00366         //cout<<"energetic HO: "<<ii<<" / "<<jj<<"  E="<<encell<<"  "<<ensumcells<<"  "<<numcells<<std::endl;
00367       }
00368     }
00369   }
00370   nenHB= (int) henergeticcellHB->Integral();
00371   nenHE= (int) henergeticcellHE->Integral();
00372   nenHF= (int) henergeticcellHF->Integral();
00373   nenHO= (int) henergeticcellHO->Integral();
00374   //std::cout<<"Energetic HB/HE/HF/HO= "<<nenHB<<" / "<<nenHE<<" / "<<nenHF<<" / "<<nenHO<<std::endl;
00375 
00376   // hot and dead cell rates
00377   for (int ii=1;ii<=hdeadHB->GetNbinsX();ii++) for (int jj=1;jj<=hdeadHB->GetNbinsY();jj++) { // --- ratio hist
00378     if (hdeadHB->GetBinContent(ii,jj)>0) {
00379       hdeadratHB->SetBinContent(ii,jj,hdeadHB->GetBinContent(ii,jj)/Nevents);
00380       hhotratHB->SetBinContent(ii,jj,hhotHB->GetBinContent(ii,jj)/Nevents);
00381     }
00382   }
00383   for (int ii=1;ii<=hdeadHE->GetNbinsX();ii++) for (int jj=1;jj<=hdeadHE->GetNbinsY();jj++) { 
00384     if (hdeadHE->GetBinContent(ii,jj)>0) {
00385       hdeadratHE->SetBinContent(ii,jj,hdeadHE->GetBinContent(ii,jj)/Nevents);
00386       hhotratHE->SetBinContent(ii,jj,hhotHE->GetBinContent(ii,jj)/Nevents);
00387     }
00388   }
00389   for (int ii=1;ii<=hdeadHF->GetNbinsX();ii++) for (int jj=1;jj<=hdeadHF->GetNbinsY();jj++) { 
00390     if (hdeadHF->GetBinContent(ii,jj)>0) {
00391       hdeadratHF->SetBinContent(ii,jj,hdeadHF->GetBinContent(ii,jj)/Nevents);
00392       hhotratHF->SetBinContent(ii,jj,hhotHF->GetBinContent(ii,jj)/Nevents);
00393     }
00394   }
00395   for (int ii=1;ii<=hdeadHO->GetNbinsX();ii++) for (int jj=1;jj<=hdeadHO->GetNbinsY();jj++) { 
00396     if (hdeadHO->GetBinContent(ii,jj)>0) {
00397       hdeadratHO->SetBinContent(ii,jj,hdeadHO->GetBinContent(ii,jj)/Nevents);
00398       hhotratHO->SetBinContent(ii,jj,hhotHO->GetBinContent(ii,jj)/Nevents);
00399     }
00400   }
00401 
00402   for (int ii=1;ii<=hhotHB->GetNbinsX();ii++) for (int jj=1;jj<=hhotHB->GetNbinsY();jj++) { // --- counters
00403     if (hdeadratHB->GetBinContent(ii,jj)>0.05) ndeadHB++;
00404     if (hhotratHB->GetBinContent(ii,jj)>0.05 || henergeticcellHB->GetBinContent(ii,jj)>0) nhotHB++;
00405   }
00406   for (int ii=1;ii<=hhotHE->GetNbinsX();ii++) for (int jj=1;jj<=hhotHE->GetNbinsY();jj++) { 
00407     if (hdeadratHE->GetBinContent(ii,jj)>0.05) ndeadHE++;
00408     if (hhotratHE->GetBinContent(ii,jj)>0.05 || henergeticcellHE->GetBinContent(ii,jj)>0) nhotHE++;
00409   }
00410   for (int ii=1;ii<=hhotHF->GetNbinsX();ii++) for (int jj=1;jj<=hhotHF->GetNbinsY();jj++) { 
00411     if (hdeadratHF->GetBinContent(ii,jj)>0.05) ndeadHF++;
00412     if (hhotratHF->GetBinContent(ii,jj)>0.05 || henergeticcellHF->GetBinContent(ii,jj)>0) nhotHF++;
00413   }
00414   for (int ii=1;ii<=hhotHO->GetNbinsX();ii++) for (int jj=1;jj<=hhotHO->GetNbinsY();jj++) { 
00415     if (hdeadratHO->GetBinContent(ii,jj)>0.05) ndeadHO++;
00416     if (hhotratHO->GetBinContent(ii,jj)>0.05 || henergeticcellHO->GetBinContent(ii,jj)>0) nhotHO++;
00417   }
00418 
00419   //std::cout<<"Dead HB/HE/HF/HO= "<<ndeadHB<<" / "<<ndeadHE<<" / "<<ndeadHF<<" / "<<ndeadHO<<std::endl;
00420   //std::cout<<"Hot  HB/HE/HF/HO= "<<nhotHB<<" / "<<nhotHE<<" / "<<nhotHF<<" / "<<nhotHO<<std::endl;
00421 
00422   Double_t valdeadHB= 1-ndeadHB/2592.0;
00423   Double_t valdeadHE= 1-ndeadHE/2592.0;
00424   Double_t valdeadHF= 1-ndeadHF/1728.0;
00425   Double_t valdeadHO= 1-ndeadHO/2160.0;
00426 
00427   Double_t valhotHB= 1-nhotHB/2592.0;
00428   Double_t valhotHE= 1-nhotHE/2592.0;
00429   Double_t valhotHF= 1-nhotHF/1728.0;
00430   Double_t valhotHO= 1-nhotHO/2160.0;
00431 
00432   Double_t valHB = valhotHB*valdeadHB;
00433   Double_t valHE = valhotHE*valdeadHE;
00434   Double_t valHF = valhotHF*valdeadHF;
00435   Double_t valHO = valhotHO*valdeadHO;
00436 
00437   Hcal_HB->Fill(valHB);
00438   Hcal_HE->Fill(valHE);
00439   Hcal_HF->Fill(valHF);
00440   Hcal_HO->Fill(valHO);
00441 
00442   //std::cout<<"Dead   HB/HE/HF/HO= "<<valdeadHB<<" / "<<valdeadHE<<" / "<<valdeadHF<<" / "<<valdeadHO<<std::endl;
00443   //std::cout<<"Hot    HB/HE/HF/HO= "<<valhotHB<<" / "<<valhotHE<<" / "<<valhotHF<<" / "<<valhotHO<<std::endl;
00444   //std::cout<<"ResVal HB/HE/HF/HO= "<<valHB<<" / "<<valHE<<" / "<<valHF<<" / "<<valHO<<std::endl;
00445 
00446 // ---------------------- end of certification
00447 
00448 }
00449 
00450 //define this as a plug-in
00451 DEFINE_FWK_MODULE(HcalDataCertification);