CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch12/src/DQM/HcalMonitorTasks/src/HcalDetDiagLEDMonitor.cc

Go to the documentation of this file.
00001 #include "DQMServices/Core/interface/MonitorElement.h"
00002 // this is to retrieve HCAL digi's
00003 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
00004 // to retrive trigger information (local runs only)
00005 #include "TBDataFormats/HcalTBObjects/interface/HcalTBTriggerData.h"
00006 // to retrive GMT information, for cosmic runs muon triggers can be used as pedestal (global runs only)
00007 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTReadoutCollection.h"
00008 // to retrive trigger desition words, to select pedestal (from hcal point of view) triggers (global runs only)
00009 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
00010 
00011 #include "CalibCalorimetry/HcalAlgos/interface/HcalLogicalMapGenerator.h"
00012 #include "CondFormats/HcalObjects/interface/HcalLogicalMap.h"
00013 #include "DQM/HcalMonitorTasks/interface/HcalBaseDQMonitor.h"
00014 #include "DQM/HcalMonitorTasks/interface/HcalEtaPhiHists.h"
00015 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
00016 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
00017 #include "CondFormats/HcalObjects/interface/HcalElectronicsMap.h"
00018 
00019 #include "TFile.h"
00020 #include "TTree.h"
00021 #include "TSystem.h"
00022 #include <math.h>
00023 
00024 class HcalDetDiagLEDData{
00025 public: 
00026    HcalDetDiagLEDData(){ 
00027              IsRefetence=false;
00028              status=0;
00029              reset();
00030           }
00031    void   reset(){
00032              Xe=XXe=Xt=XXt=n=0;
00033              overflow=0;
00034              undeflow=0;
00035           }
00036    void   add_statistics(double *data,int nTS){
00037              double e=GetEnergy(data,nTS);
00038              double t=GetTime(data,nTS);
00039              if(e<20) undeflow++; else if(e>10000) overflow++; else{
00040                 n++; Xe+=e; XXe+=e*e; Xt+=t; XXt+=t*t;
00041              }     
00042           }
00043    void   set_reference(float val,float rms){
00044              ref_led=val; ref_rms=rms;
00045              IsRefetence=true;
00046           }       
00047    void   change_status(int val){
00048              status|=val;
00049           }       
00050    int    get_status(){
00051              return status;
00052           }       
00053    bool   get_reference(double *val,double *rms){
00054              *val=ref_led; *rms=ref_rms;
00055              return IsRefetence;
00056           }       
00057    bool   get_average_led(double *ave,double *rms){
00058              if(n>0){ *ave=Xe/n; *rms=sqrt(XXe/n-(Xe*Xe)/(n*n));} else return false;
00059              return true; 
00060           }
00061    bool   get_average_time(double *ave,double *rms){
00062              if(n>0){ *ave=Xt/n; *rms=sqrt(XXt/n-(Xt*Xt)/(n*n));} else return false;
00063              return true; 
00064           }
00065    int    get_statistics(){
00066              return (int)n;
00067           } 
00068    int    get_overflow(){
00069              return overflow;
00070           }   
00071    int    get_undeflow(){
00072              return undeflow;
00073           }   
00074 private:   
00075    double GetEnergy(double *data,int n){
00076              int MaxI=0; double Energy,MaxE=0;
00077              for(int j=0;j<n;++j) if(MaxE<data[j]){ MaxE=data[j]; MaxI=j; }
00078              Energy=data[MaxI];
00079              if(MaxI>0) Energy+=data[MaxI-1];
00080              if(MaxI>1) Energy+=data[MaxI-2];
00081              if(MaxI<(n-1)) Energy+=data[MaxI+1];
00082              if(MaxI<(n-2)) Energy+=data[MaxI+2];
00083              return Energy;
00084           }
00085    double GetTime(double *data,int n=10){
00086              int MaxI=0; double Time=-9999,SumT=0,MaxT=-10;
00087              for(int j=0;j<n;++j) if(MaxT<data[j]){ MaxT=data[j]; MaxI=j; }
00088              Time=MaxI*data[MaxI];
00089              SumT=data[MaxI];
00090              if(MaxI>0){ Time+=(MaxI-1)*data[MaxI-1]; SumT+=data[MaxI-1]; }
00091              if(MaxI<(n-1)){ Time+=(MaxI+1)*data[MaxI+1]; SumT+=data[MaxI+1]; }
00092              Time=Time/SumT;
00093              return Time;
00094    }      
00095    int   overflow;
00096    int   undeflow;
00097    double Xe,XXe,Xt,XXt,n;
00098    bool  IsRefetence;
00099    float ref_led;
00100    float ref_rms;
00101    int   status;
00102 };
00103 
00104 class HcalDetDiagLEDMonitor:public HcalBaseDQMonitor {
00105 public:
00106   HcalDetDiagLEDMonitor(const edm::ParameterSet& ps); 
00107   ~HcalDetDiagLEDMonitor(); 
00108 
00109   void beginRun(const edm::Run& run, const edm::EventSetup& c);
00110   void setup();
00111   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup);// const HcalDbService& cond)
00112   void endRun(const edm::Run& run, const edm::EventSetup& c);
00113   void reset();
00114   void cleanup(); 
00115   void fillHistos();
00116   int  GetStatistics(){ return ievt_; }
00117 private:
00118   HcalLogicalMapGenerator *gen;
00119   HcalElectronicsMap      emap;
00120   HcalLogicalMap          *lmap;
00121   void SaveReference();
00122   void LoadReference();
00123   void CheckStatus();
00124   
00125   HcalDetDiagLEDData* GetCalib(std::string sd,int eta,int phi){
00126     int SD=0,ETA=0,PHI=0;
00127     if(sd.compare("HB")==0) SD=1; 
00128     if(sd.compare("HE")==0) SD=2; 
00129     if(sd.compare("HO")==0) SD=3; 
00130     if(sd.compare("HF")==0) SD=4; 
00131     if(SD==1 || SD==2){
00132       if(eta>0) ETA=1; else ETA=-1;
00133       if(phi==71 ||phi==72 || phi==1 || phi==2) PHI=71; else PHI=((phi-3)/4)*4+3;
00134     }else if(SD==3){
00135       if(abs(eta)<=4){
00136         ETA=0;
00137         if(phi==71 ||phi==72 || phi==1 || phi==2 || phi==3 || phi==4) PHI=71; else PHI=((phi-5)/6)*6+5;
00138       }else{
00139         if(abs(eta)>4  && abs(eta)<=10)  ETA=1;
00140         if(abs(eta)>10 && abs(eta)<=15)  ETA=2;
00141         if(eta<0) ETA=-ETA;
00142         if(phi==71 ||phi==72 || (phi>=1 && phi<=10)) PHI=71; else PHI=((phi-11)/12)*12+11;
00143       }
00144     }else if(SD==4){
00145       if(eta>0) ETA=1; else ETA=-1;
00146       if(phi>=1  && phi<=18) PHI=1;
00147       if(phi>=19 && phi<=36) PHI=19;
00148       if(phi>=37 && phi<=54) PHI=37;
00149       if(phi>=55 && phi<=72) PHI=55;
00150     }
00151     return &calib_data[SD][ETA+2][PHI-1];
00152   };
00153   int         ievt_;
00154   int         run_number;
00155   int         dataset_seq_number;
00156   bool        IsReference;
00157   
00158   double      LEDMeanTreshold;
00159   double      LEDRmsTreshold;
00160    
00161   std::string ReferenceData;
00162   std::string ReferenceRun;
00163   std::string OutputFilePath;
00164 
00165   MonitorElement *meEVT_,*meRUN_;
00166   MonitorElement *RefRun_;
00167   MonitorElement *Energy;
00168   MonitorElement *Time;
00169   MonitorElement *EnergyHF;
00170   MonitorElement *TimeHF;
00171   MonitorElement *Time2Dhbhehf;
00172   MonitorElement *Time2Dho;
00173   MonitorElement *Energy2Dhbhehf;
00174   MonitorElement *Energy2Dho;
00175   MonitorElement *EnergyRMS;
00176   MonitorElement *TimeRMS;
00177   MonitorElement *EnergyRMSHF;
00178   MonitorElement *TimeRMSHF;
00179   MonitorElement *EnergyCorr;
00180   MonitorElement *HBPphi;
00181   MonitorElement *HBMphi;
00182   MonitorElement *HEPphi;
00183   MonitorElement *HEMphi;
00184   MonitorElement *HFPphi;
00185   MonitorElement *HFMphi;
00186   MonitorElement *HO0phi;
00187   MonitorElement *HO1Pphi;
00188   MonitorElement *HO2Pphi;
00189   MonitorElement *HO1Mphi;
00190   MonitorElement *HO2Mphi;
00191 
00192   HcalDetDiagLEDData hb_data[85][72][4];
00193   HcalDetDiagLEDData he_data[85][72][4];
00194   HcalDetDiagLEDData ho_data[85][72][4];
00195   HcalDetDiagLEDData hf_data[85][72][4];
00196   HcalDetDiagLEDData calib_data[5][5][72];
00197   
00198   EtaPhiHists *ChannelsLEDEnergy;
00199   EtaPhiHists *ChannelsLEDEnergyRef;
00200   EtaPhiHists *ChannelStatusMissingChannels;
00201   EtaPhiHists *ChannelStatusUnstableChannels;
00202   EtaPhiHists *ChannelStatusUnstableLEDsignal;
00203   EtaPhiHists *ChannelStatusLEDMean;
00204   EtaPhiHists *ChannelStatusLEDRMS;
00205   EtaPhiHists *ChannelStatusTimeMean;
00206   EtaPhiHists *ChannelStatusTimeRMS;
00207 
00208   edm::InputTag digiLabel_;
00209   edm::InputTag calibDigiLabel_;
00210 
00211   void fill_channel_status(std::string subdet,int eta,int phi,int depth,int type,double status);
00212   void   fill_energy(std::string subdet,int eta,int phi,int depth,double e,int type);
00213   double get_energy(std::string subdet,int eta,int phi,int depth,int type);
00214 };
00215 
00217 static const float adc2fC[128]={-0.5,0.5,1.5,2.5,3.5,4.5,5.5,6.5,7.5,8.5,9.5, 10.5,11.5,12.5,
00218                    13.5,15.,17.,19.,21.,23.,25.,27.,29.5,32.5,35.5,38.5,42.,46.,50.,54.5,59.5,
00219                    64.5,59.5,64.5,69.5,74.5,79.5,84.5,89.5,94.5,99.5,104.5,109.5,114.5,119.5,
00220                    124.5,129.5,137.,147.,157.,167.,177.,187.,197.,209.5,224.5,239.5,254.5,272.,
00221                    292.,312.,334.5,359.5,384.5,359.5,384.5,409.5,434.5,459.5,484.5,509.5,534.5,
00222                    559.5,584.5,609.5,634.5,659.5,684.5,709.5,747.,797.,847.,897.,947.,997.,
00223                    1047.,1109.5,1184.5,1259.5,1334.5,1422.,1522.,1622.,1734.5,1859.5,1984.5,
00224                    1859.5,1984.5,2109.5,2234.5,2359.5,2484.5,2609.5,2734.5,2859.5,2984.5,
00225                    3109.5,3234.5,3359.5,3484.5,3609.5,3797.,4047.,4297.,4547.,4797.,5047.,
00226                    5297.,5609.5,5984.5,6359.5,6734.5,7172.,7672.,8172.,8734.5,9359.5,9984.5};
00228 
00229 
00230 
00231 
00232 HcalDetDiagLEDMonitor::HcalDetDiagLEDMonitor(const edm::ParameterSet& ps) {
00233   ievt_=0;
00234   dataset_seq_number=1;
00235   run_number=-1;
00236   IsReference=false;
00237 
00238   Online_                = ps.getUntrackedParameter<bool>("online",false);
00239   mergeRuns_             = ps.getUntrackedParameter<bool>("mergeRuns",false);
00240   enableCleanup_         = ps.getUntrackedParameter<bool>("enableCleanup",false);
00241   debug_                 = ps.getUntrackedParameter<int>("debug",0);
00242   prefixME_              = ps.getUntrackedParameter<std::string>("subSystemFolder","Hcal/");
00243   if (prefixME_.substr(prefixME_.size()-1,prefixME_.size())!="/")
00244     prefixME_.append("/");
00245   subdir_                = ps.getUntrackedParameter<std::string>("TaskFolder","DetDiagLEDMonitor_Hcal");
00246   if (subdir_.size()>0 && subdir_.substr(subdir_.size()-1,subdir_.size())!="/")
00247     subdir_.append("/");
00248   subdir_=prefixME_+subdir_;
00249   AllowedCalibTypes_     = ps.getUntrackedParameter<std::vector<int> > ("AllowedCalibTypes");
00250   skipOutOfOrderLS_      = ps.getUntrackedParameter<bool>("skipOutOfOrderLS","false");
00251   NLumiBlocks_           = ps.getUntrackedParameter<int>("NLumiBlocks",4000);
00252   makeDiagnostics_       = ps.getUntrackedParameter<bool>("makeDiagnostics",false);
00253 
00254   LEDMeanTreshold  = ps.getUntrackedParameter<double>("LEDMeanTreshold" , 0.1);
00255   LEDRmsTreshold   = ps.getUntrackedParameter<double>("LEDRmsTreshold"  , 0.1);
00256   
00257   ReferenceData    = ps.getUntrackedParameter<std::string>("LEDReferenceData" ,"");
00258   OutputFilePath   = ps.getUntrackedParameter<std::string>("OutputFilePath", "");
00259 
00260   digiLabel_       = ps.getUntrackedParameter<edm::InputTag>("digiLabel", edm::InputTag("hcalDigis"));
00261   calibDigiLabel_  = ps.getUntrackedParameter<edm::InputTag>("calibDigiLabel",edm::InputTag("hcalDigis"));
00262 }
00263 
00264 HcalDetDiagLEDMonitor::~HcalDetDiagLEDMonitor(){}
00265 
00266 void HcalDetDiagLEDMonitor::cleanup(){
00267   if(dbe_){
00268     dbe_->setCurrentFolder(subdir_);
00269     dbe_->removeContents();
00270     dbe_ = 0;
00271   }
00272 } 
00273 void HcalDetDiagLEDMonitor::reset(){}
00274 
00275 void HcalDetDiagLEDMonitor::beginRun(const edm::Run& run, const edm::EventSetup& c)
00276 {
00277   if (debug_>1) std::cout <<"HcalDetDiagLEDMonitor::beginRun"<<std::endl;
00278   HcalBaseDQMonitor::beginRun(run,c);
00279 
00280   if (tevt_==0) this->setup(); // set up histograms if they have not been created before
00281   if (mergeRuns_==false)
00282     this->reset();
00283 
00284   return;
00285 } // void HcalNDetDiagLEDMonitor::beginRun(...)
00286 
00287 void HcalDetDiagLEDMonitor::setup(){
00288      // Call base class setup
00289      HcalBaseDQMonitor::setup();
00290      if (!dbe_) return;
00291 
00292      std::string name;
00293      dbe_->setCurrentFolder(subdir_);   
00294      meEVT_ = dbe_->bookInt("HcalDetDiagLEDMonitor Event Number");
00295      meRUN_ = dbe_->bookInt("HcalDetDiagLEDMonitor Run Number");
00296      ReferenceRun="UNKNOWN";
00297      LoadReference();
00298      dbe_->setCurrentFolder(subdir_);
00299      RefRun_= dbe_->bookString("HcalDetDiagLEDMonitor Reference Run",ReferenceRun);
00300      dbe_->setCurrentFolder(subdir_+"Summary Plots");
00301      
00302      name="HBHEHO LED Energy Distribution";               Energy         = dbe_->book1D(name,name,200,0,3000);
00303      name="HBHEHO LED Timing Distribution";               Time           = dbe_->book1D(name,name,200,0,10);
00304      name="HBHEHO LED Energy RMS_div_Energy Distribution";EnergyRMS      = dbe_->book1D(name,name,200,0,0.2);
00305      name="HBHEHO LED Timing RMS Distribution";           TimeRMS        = dbe_->book1D(name,name,200,0,0.4);
00306      name="HF LED Energy Distribution";                   EnergyHF       = dbe_->book1D(name,name,200,0,3000);
00307      name="HF LED Timing Distribution";                   TimeHF         = dbe_->book1D(name,name,200,0,10);
00308      name="HF LED Energy RMS_div_Energy Distribution";    EnergyRMSHF    = dbe_->book1D(name,name,200,0,0.5);
00309      name="HF LED Timing RMS Distribution";               TimeRMSHF      = dbe_->book1D(name,name,200,0,0.4);
00310      name="LED Energy Corr(PinDiod) Distribution";        EnergyCorr     = dbe_->book1D(name,name,200,0,10);
00311      name="LED Timing HBHEHF";                            Time2Dhbhehf   = dbe_->book2D(name,name,87,-43,43,74,0,73);
00312      name="LED Timing HO";                                Time2Dho       = dbe_->book2D(name,name,33,-16,16,74,0,73);
00313      name="LED Energy HBHEHF";                            Energy2Dhbhehf = dbe_->book2D(name,name,87,-43,43,74,0,73);
00314      name="LED Energy HO";                                Energy2Dho     = dbe_->book2D(name,name,33,-16,16,74,0,73);
00315 
00316      name="HBP Average over HPD LED Ref";          HBPphi = dbe_->book2D(name,name,180,1,73,400,0,2);
00317      name="HBM Average over HPD LED Ref";          HBMphi = dbe_->book2D(name,name,180,1,73,400,0,2);
00318      name="HEP Average over HPD LED Ref";          HEPphi = dbe_->book2D(name,name,180,1,73,400,0,2);
00319      name="HEM Average over HPD LED Ref";          HEMphi = dbe_->book2D(name,name,180,1,73,400,0,2);
00320      name="HFP Average over RM LED Ref";           HFPphi = dbe_->book2D(name,name,180,1,37,400,0,2);
00321      name="HFM Average over RM LED Ref";           HFMphi = dbe_->book2D(name,name,180,1,37,400,0,2);
00322      name="HO0 Average over HPD LED Ref";          HO0phi = dbe_->book2D(name,name,180,1,49,400,0,2);
00323      name="HO1P Average over HPD LED Ref";         HO1Pphi= dbe_->book2D(name,name,180,1,49,400,0,2);
00324      name="HO2P Average over HPD LED Ref";         HO2Pphi= dbe_->book2D(name,name,180,1,49,400,0,2);
00325      name="HO1M Average over HPD LED Ref";         HO1Mphi= dbe_->book2D(name,name,180,1,49,400,0,2);
00326      name="HO2M Average over HPD LED Ref";         HO2Mphi= dbe_->book2D(name,name,180,1,49,400,0,2);
00327 
00328      ChannelsLEDEnergy = new EtaPhiHists();
00329      ChannelsLEDEnergy->setup(dbe_," Channel LED Energy");
00330      ChannelsLEDEnergyRef = new EtaPhiHists();
00331      ChannelsLEDEnergyRef->setup(dbe_," Channel LED Energy Reference");
00332      
00333      dbe_->setCurrentFolder(subdir_+"channel status");
00334      ChannelStatusMissingChannels = new EtaPhiHists();
00335      ChannelStatusMissingChannels->setup(dbe_," Missing Channels");
00336      ChannelStatusUnstableChannels = new EtaPhiHists();
00337      ChannelStatusUnstableChannels->setup(dbe_," Unstable Channels");
00338      ChannelStatusUnstableLEDsignal = new EtaPhiHists();
00339      ChannelStatusUnstableLEDsignal->setup(dbe_," Unstable LED");
00340      ChannelStatusLEDMean = new EtaPhiHists();
00341      ChannelStatusLEDMean->setup(dbe_," LED Mean");
00342      ChannelStatusLEDRMS = new EtaPhiHists();
00343      ChannelStatusLEDRMS->setup(dbe_," LED RMS");
00344      ChannelStatusTimeMean = new EtaPhiHists();
00345      ChannelStatusTimeMean->setup(dbe_," Time Mean");
00346      ChannelStatusTimeRMS = new EtaPhiHists();
00347      ChannelStatusTimeRMS->setup(dbe_," Time RMS");
00348 
00349      gen=new HcalLogicalMapGenerator();
00350      lmap =new HcalLogicalMap(gen->createMap());
00351      emap=lmap->generateHcalElectronicsMap();
00352      return;
00353 } 
00354 
00355 void HcalDetDiagLEDMonitor::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup){
00356 int  eta,phi,depth,nTS;
00357    if(!dbe_) return; 
00358    bool LEDEvent=false;
00359    bool LocalRun=false;
00360    // for local runs 
00361 
00362    edm::Handle<HcalTBTriggerData> trigger_data;
00363    iEvent.getByType(trigger_data);
00364    if(trigger_data.isValid()){
00365       if(trigger_data->triggerWord()==6){ LEDEvent=true;LocalRun=true;}
00366    } 
00367    if(!LocalRun) return;  
00368    if(!LEDEvent) return; 
00369    
00370    HcalBaseDQMonitor::analyze(iEvent, iSetup);
00371    meEVT_->Fill(++ievt_);
00372    run_number=iEvent.id().run();
00373    meRUN_->Fill(iEvent.id().run());
00374 
00375    double data[20];
00376 
00377    edm::Handle<HBHEDigiCollection> hbhe; 
00378    iEvent.getByLabel(digiLabel_, hbhe);
00379    if(hbhe.isValid()) for(HBHEDigiCollection::const_iterator digi=hbhe->begin();digi!=hbhe->end();digi++){
00380      eta=digi->id().ieta(); phi=digi->id().iphi(); depth=digi->id().depth(); nTS=digi->size();
00381      if(digi->id().subdet()==HcalBarrel){
00382        for(int i=0;i<nTS;i++) data[i]=adc2fC[digi->sample(i).adc()&0xff]-2.5;
00383        hb_data[eta+42][phi-1][depth-1].add_statistics(data,nTS);
00384      }   
00385      if(digi->id().subdet()==HcalEndcap){
00386        for(int i=0;i<nTS;i++) data[i]=adc2fC[digi->sample(i).adc()&0xff]-2.5;
00387        he_data[eta+42][phi-1][depth-1].add_statistics(data,nTS);
00388      }
00389    }   
00390 
00391    edm::Handle<HODigiCollection> ho; 
00392    iEvent.getByLabel(digiLabel_,ho);
00393    if(ho.isValid()) for(HODigiCollection::const_iterator digi=ho->begin();digi!=ho->end();digi++){
00394      eta=digi->id().ieta(); phi=digi->id().iphi(); depth=digi->id().depth(); nTS=digi->size();
00395      for(int i=0;i<nTS;i++) data[i]=adc2fC[digi->sample(i).adc()&0xff]-2.5;
00396      ho_data[eta+42][phi-1][depth-1].add_statistics(data,nTS);
00397    }   
00398 
00399    edm::Handle<HFDigiCollection> hf;
00400    iEvent.getByLabel(digiLabel_,hf);
00401    if(hf.isValid()) for(HFDigiCollection::const_iterator digi=hf->begin();digi!=hf->end();digi++){
00402      eta=digi->id().ieta(); phi=digi->id().iphi(); depth=digi->id().depth(); nTS=digi->size();
00403      for(int i=0;i<nTS;i++) data[i]=adc2fC[digi->sample(i).adc()&0xff]-2.5;
00404      hf_data[eta+42][phi-1][depth-1].add_statistics(data,nTS);
00405    }   
00406  
00407    edm::Handle<HcalCalibDigiCollection> calib;
00408    iEvent.getByLabel(calibDigiLabel_, calib);
00409    if(calib.isValid())for(HcalCalibDigiCollection::const_iterator digi=calib->begin();digi!=calib->end();digi++){
00410      if(digi->id().cboxChannel()!=0 || digi->id().hcalSubdet()==0) continue; 
00411      nTS=digi->size();
00412      double e=0; 
00413      for(int i=0;i<nTS;i++){ data[i]=adc2fC[digi->sample(i).adc()&0xff]; e+=data[i];}
00414      if(e<15000) calib_data[digi->id().hcalSubdet()][digi->id().ieta()+2][digi->id().iphi()-1].add_statistics(data,nTS);
00415    }   
00416   
00417    if(((ievt_)%500)==0){
00418        fillHistos();
00419        CheckStatus(); 
00420    }
00421    return;
00422 }
00423 
00424 
00425 void HcalDetDiagLEDMonitor::fillHistos(){
00426   std::string subdet[4]={"HB","HE","HO","HF"};
00427     Energy->Reset();
00428    Time->Reset();
00429    EnergyRMS->Reset();
00430    TimeRMS->Reset();
00431    EnergyHF->Reset();
00432    TimeHF->Reset();
00433    EnergyRMSHF->Reset();
00434    TimeRMSHF->Reset();
00435    EnergyCorr->Reset();
00436    Time2Dhbhehf->Reset();
00437    Time2Dho->Reset();
00438    Energy2Dhbhehf->Reset();
00439    Energy2Dho->Reset();
00440    HBPphi->Reset();
00441    HBMphi->Reset();
00442    HEPphi->Reset();
00443    HEMphi->Reset();
00444    HFPphi->Reset();
00445    HFMphi->Reset();
00446    HO0phi->Reset();
00447    HO1Pphi->Reset();
00448    HO2Pphi->Reset();
00449    HO1Mphi->Reset();
00450    HO2Mphi->Reset();
00451    
00452    // HB histograms
00453    for(int eta=-16;eta<=16;eta++) for(int phi=1;phi<=72;phi++){ 
00454       double T=0,nT=0,E=0,nE=0;
00455       for(int depth=1;depth<=2;depth++){
00456          if(hb_data[eta+42][phi-1][depth-1].get_statistics()>100){
00457             double ave=0;
00458             double rms=0;
00459             double time=0;
00460             double time_rms=0;
00461             hb_data[eta+42][phi-1][depth-1].get_average_led(&ave,&rms);
00462             hb_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00463             Energy->Fill(ave);
00464             if(ave>0)EnergyRMS->Fill(rms/ave);
00465             Time->Fill(time);
00466             TimeRMS->Fill(time_rms);
00467             T+=time; nT++; E+=ave; nE++;
00468             if(GetCalib("HB",eta,phi)->get_statistics()>100){
00469               double ave_calib=0;
00470               double rms_calib=0;
00471               GetCalib("HB",eta,phi)->get_average_led(&ave_calib,&rms_calib);
00472               fill_energy("HB",eta,phi,depth,ave/ave_calib,1);
00473               EnergyCorr->Fill(ave_calib/ave);
00474             }
00475          }
00476       } 
00477       if(nT>0){Time2Dhbhehf->Fill(eta,phi,T/nT);Energy2Dhbhehf->Fill(eta,phi,E/nE); }
00478    } 
00479    // HE histograms
00480    for(int eta=-29;eta<=29;eta++) for(int phi=1;phi<=72;phi++){
00481       double T=0,nT=0,E=0,nE=0;
00482       for(int depth=1;depth<=3;depth++){
00483          if(he_data[eta+42][phi-1][depth-1].get_statistics()>100){
00484             double ave=0;
00485             double rms=0;
00486             double time=0;
00487             double time_rms=0;
00488             he_data[eta+42][phi-1][depth-1].get_average_led(&ave,&rms);
00489             he_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00490             Energy->Fill(ave);
00491             if(ave>0)EnergyRMS->Fill(rms/ave);
00492             Time->Fill(time);
00493             T+=time; nT++; E+=ave; nE++;
00494             TimeRMS->Fill(time_rms);
00495             if(GetCalib("HE",eta,phi)->get_statistics()>100){
00496               double ave_calib=0;
00497               double rms_calib=0;
00498               GetCalib("HE",eta,phi)->get_average_led(&ave_calib,&rms_calib);
00499               fill_energy("HE",eta,phi,depth,ave/ave_calib,1);
00500               EnergyCorr->Fill(ave_calib/ave);
00501             }
00502          }
00503       }
00504       if(nT>0 && abs(eta)>16 ){Time2Dhbhehf->Fill(eta,phi,T/nT);   Energy2Dhbhehf->Fill(eta,phi,E/nE); }         
00505       if(nT>0 && abs(eta)>20 ){Time2Dhbhehf->Fill(eta,phi+1,T/nT); Energy2Dhbhehf->Fill(eta,phi+1,E/nE);}        
00506    } 
00507    // HF histograms
00508    for(int eta=-42;eta<=42;eta++) for(int phi=1;phi<=72;phi++){
00509       double T=0,nT=0,E=0,nE=0;
00510       for(int depth=1;depth<=2;depth++){
00511          if(hf_data[eta+42][phi-1][depth-1].get_statistics()>100){
00512            double ave=0;
00513            double rms=0;
00514            double time=0;
00515            double time_rms=0;
00516            hf_data[eta+42][phi-1][depth-1].get_average_led(&ave,&rms);
00517            hf_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00518            EnergyHF->Fill(ave);
00519            if(ave>0)EnergyRMSHF->Fill(rms/ave);
00520            TimeHF->Fill(time);
00521            T+=time; nT++; E+=ave; nE++;
00522            TimeRMSHF->Fill(time_rms);
00523            if(GetCalib("HF",eta,phi)->get_statistics()>100){
00524              double ave_calib=0;
00525              double rms_calib=0;
00526              GetCalib("HF",eta,phi)->get_average_led(&ave_calib,&rms_calib);
00527              fill_energy("HF",eta,phi,depth,ave/ave_calib,1);
00528              EnergyCorr->Fill(ave_calib/ave);
00529            }
00530          }
00531       } 
00532       if(nT>0 && abs(eta)>29 ){ Time2Dhbhehf->Fill(eta,phi,T/nT); Time2Dhbhehf->Fill(eta,phi+1,T/nT);}   
00533       if(nT>0 && abs(eta)>29 ){ Energy2Dhbhehf->Fill(eta,phi,E/nE); Energy2Dhbhehf->Fill(eta,phi+1,E/nE);}       
00534    } 
00535    // HO histograms
00536    for(int eta=-15;eta<=15;eta++) for(int phi=1;phi<=72;phi++){
00537       double T=0,nT=0,E=0,nE=0;
00538       for(int depth=4;depth<=4;depth++){
00539          if(ho_data[eta+42][phi-1][depth-1].get_statistics()>100){
00540             double ave=0;
00541             double rms=0;
00542             double time=0;
00543             double time_rms=0;
00544             ho_data[eta+42][phi-1][depth-1].get_average_led(&ave,&rms);
00545             ho_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00546             Energy->Fill(ave);
00547             if(ave>0)EnergyRMS->Fill(rms/ave);
00548             Time->Fill(time);
00549             T+=time; nT++; E+=ave; nE++;
00550             TimeRMS->Fill(time_rms);
00551             if(GetCalib("HO",eta,phi)->get_statistics()>100){
00552               double ave_calib=0;
00553               double rms_calib=0;
00554               GetCalib("HO",eta,phi)->get_average_led(&ave_calib,&rms_calib);
00555               fill_energy("HO",eta,phi,depth,ave/ave_calib,1);
00556               EnergyCorr->Fill(ave_calib/ave);
00557             }
00558          }
00559       }
00560       if(nT>0){ Time2Dho->Fill(eta,phi,T/nT); Energy2Dho->Fill(eta,phi+1,E/nE) ;}
00561    } 
00562 
00563    double ave = 0.,rms = 0.,ave_calib = 0.,rms_calib = 0.;
00564    // HB Ref histograms
00565    for(int eta=-16;eta<=16;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=2;depth++){
00566       if(hb_data[eta+42][phi-1][depth-1].get_reference(&ave,&rms) && GetCalib("HB",eta,phi)->get_reference(&ave_calib,&rms_calib)){
00567             fill_energy("HB",eta,phi,depth,ave/ave_calib,2);
00568       }
00569    } 
00570    // HE Ref histograms
00571    for(int eta=-29;eta<=29;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=3;depth++){
00572       if(he_data[eta+42][phi-1][depth-1].get_reference(&ave,&rms) && GetCalib("HE",eta,phi)->get_reference(&ave_calib,&rms_calib)){
00573             fill_energy("HE",eta,phi,depth,ave/ave_calib,2);
00574       }
00575    } 
00576    // HO Ref histograms
00577    for(int eta=-15;eta<=15;eta++) for(int phi=1;phi<=72;phi++) for(int depth=4;depth<=4;depth++){
00578       if(ho_data[eta+42][phi-1][depth-1].get_reference(&ave,&rms) && GetCalib("HO",eta,phi)->get_reference(&ave_calib,&rms_calib)){
00579             fill_energy("HO",eta,phi,depth,ave/ave_calib,2);
00580       }
00581    } 
00582    // HF Ref histograms
00583    for(int eta=-42;eta<=42;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=2;depth++){
00584       if(hf_data[eta+42][phi-1][depth-1].get_reference(&ave,&rms) && GetCalib("HF",eta,phi)->get_reference(&ave_calib,&rms_calib)){
00585             fill_energy("HF",eta,phi,depth,ave/ave_calib,2);
00586       }
00587    } 
00588 
00589   //fill RM histograms: this part is incomplete, will be modefied later 
00590   double hbp[18][4],nhbp[18][4],hbm[18][4],nhbm[18][4];
00591   double hep[18][4],nhep[18][4],hem[18][4],nhem[18][4];
00592   double hfp[18][4],nhfp[18][4],hfm[18][4],nhfm[18][4];
00593   double ho0[18][4],nho0[18][4];
00594   double ho1p[18][4],nho1p[18][4];
00595   double ho2p[18][4],nho2p[18][4];
00596   double ho1m[18][4],nho1m[18][4];
00597   double ho2m[18][4],nho2m[18][4];
00598   for(int i=0;i<18;i++) for(int j=0;j<4;j++)
00599    hbp[i][j]=nhbp[i][j]=hbm[i][j]=nhbm[i][j]=hep[i][j]=nhep[i][j]=hem[i][j]=nhem[i][j]=hfp[i][j]=nhfp[i][j]=hfm[i][j]=nhfm[i][j]=0;
00600   for(int i=0;i<18;i++) for(int j=0;j<4;j++)
00601    ho0[i][j]=nho0[i][j]=ho1p[i][j]=nho1p[i][j]=ho2p[i][j]=nho2p[i][j]=ho1m[i][j]=nho1m[i][j]=ho2m[i][j]=nho2m[i][j]=0;
00602 
00603    std::vector <HcalElectronicsId> AllElIds = emap.allElectronicsIdPrecision();
00604    for(std::vector <HcalElectronicsId>::iterator eid = AllElIds.begin(); eid != AllElIds.end(); eid++){
00605       DetId detid=emap.lookup(*eid);
00606       if(detid.det()!=DetId::Hcal) continue;
00607       HcalGenericDetId gid(emap.lookup(*eid));
00608       if(!(!(gid.null()) && 
00609             (gid.genericSubdet()==HcalGenericDetId::HcalGenBarrel ||
00610              gid.genericSubdet()==HcalGenericDetId::HcalGenEndcap  ||
00611              gid.genericSubdet()==HcalGenericDetId::HcalGenForward ||
00612              gid.genericSubdet()==HcalGenericDetId::HcalGenOuter))) continue;
00613       int sd=0,eta,phi,depth; 
00614       if(gid.genericSubdet()==HcalGenericDetId::HcalGenBarrel)      sd=0;
00615       else if(gid.genericSubdet()==HcalGenericDetId::HcalGenEndcap) sd=1;
00616       else if(gid.genericSubdet()==HcalGenericDetId::HcalGenOuter)  sd=2;
00617       else if(gid.genericSubdet()==HcalGenericDetId::HcalGenForward)sd=3;
00618       HcalDetId hid(detid);
00619       eta=hid.ieta();
00620       phi=hid.iphi();
00621       depth=hid.depth(); 
00622 
00623       double ave =get_energy(subdet[sd],eta,phi,depth,1);
00624       double ref =get_energy(subdet[sd],eta,phi,depth,2);
00625 
00626       HcalFrontEndId  lmap_entry=lmap->getHcalFrontEndId(hid);
00627       int rbx; 
00628       if(sd==0 || sd==1 || sd==3){
00629            sscanf(&(lmap_entry.rbx().c_str())[3],"%d",&rbx);
00630       }else{
00631            if(abs(eta)<5) sscanf(&(lmap_entry.rbx().c_str())[3],"%d",&rbx);
00632            if(abs(eta)>=5) sscanf(&(lmap_entry.rbx().c_str())[4],"%d",&rbx);           
00633       }
00634       if(ave>0 && ref>0){
00635            if(sd==0 && eta>0){ hbp[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhbp[rbx-1][lmap_entry.rm()-1]++; }
00636            if(sd==0 && eta<0){ hbm[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhbm[rbx-1][lmap_entry.rm()-1]++; }
00637            if(sd==1 && eta>0){ hep[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhep[rbx-1][lmap_entry.rm()-1]++; }
00638            if(sd==1 && eta<0){ hem[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhem[rbx-1][lmap_entry.rm()-1]++; }
00639            if(sd==3 && eta>0){ hfp[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhfp[rbx-1][lmap_entry.rm()-1]++; }
00640            if(sd==3 && eta<0){ hfm[rbx-1][lmap_entry.rm()-1]+=ave/ref; nhfm[rbx-1][lmap_entry.rm()-1]++; }
00641            if(sd==2 && abs(eta)<5){ ho0[rbx-1][lmap_entry.rm()-1]+=ave/ref; nho0[rbx-1][lmap_entry.rm()-1]++; }
00642            if(sd==2 && eta>=5 && eta<=10){ ho1p[rbx-1][lmap_entry.rm()-1]+=ave/ref; nho1p[rbx-1][lmap_entry.rm()-1]++; }
00643            if(sd==2 && eta>=11 && eta<=15){ ho2p[rbx-1][lmap_entry.rm()-1]+=ave/ref; nho2p[rbx-1][lmap_entry.rm()-1]++; }
00644            if(sd==2 && eta>=-10 && eta<=-5){ ho1m[rbx-1][lmap_entry.rm()-1]+=ave/ref; nho1m[rbx-1][lmap_entry.rm()-1]++; }
00645            if(sd==2 && eta>=-15 && eta<=-11){ ho2m[rbx-1][lmap_entry.rm()-1]+=ave/ref; nho2m[rbx-1][lmap_entry.rm()-1]++; }
00646       }
00647   }  
00648   for(int i=0;i<18;i++)for(int j=0;j<4;j++){
00649      int phi=i*4+j+1; 
00650      if(nhbp[i][j]>1) HBPphi->Fill(phi+0.5,hbp[i][j]/nhbp[i][j]);
00651      if(nhbm[i][j]>1) HBMphi->Fill(phi+0.5,hbm[i][j]/nhbm[i][j]);
00652      if(nhep[i][j]>1) HEPphi->Fill(phi+0.5,hep[i][j]/nhep[i][j]);
00653      if(nhem[i][j]>1) HEMphi->Fill(phi+0.5,hem[i][j]/nhem[i][j]);
00654   }   
00655   for(int i=0;i<12;i++)for(int j=0;j<3;j++){
00656      int phi=i*3+j+1; 
00657      if(nhfp[i][j]>1) HFPphi->Fill(phi+0.5,hfp[i][j]/nhfp[i][j]);
00658      if(nhfm[i][j]>1) HFMphi->Fill(phi+0.5,hfm[i][j]/nhfm[i][j]);
00659   } 
00660   for(int i=0;i<12;i++)for(int j=0;j<4;j++){
00661      int phi=i*4+j+1; 
00662      if(nho0[i][j]>1) HO0phi->Fill(phi+0.5,ho0[i][j]/nho0[i][j]);
00663      if(nho1p[i][j]>1) HO1Pphi->Fill(phi+0.5,ho1p[i][j]/nho1p[i][j]);
00664      if(nho2p[i][j]>1) HO2Pphi->Fill(phi+0.5,ho2p[i][j]/nho2p[i][j]);
00665      if(nho1m[i][j]>1) HO1Mphi->Fill(phi+0.5,ho1m[i][j]/nho1m[i][j]);
00666      if(nho2m[i][j]>1) HO2Mphi->Fill(phi+0.5,ho2m[i][j]/nho2m[i][j]);
00667   } 
00668 } 
00669 
00670 void HcalDetDiagLEDMonitor::SaveReference(){
00671 double led,rms,time,time_rms;
00672 int    Eta,Phi,Depth,Statistic,Status=0;
00673 char   Subdet[10],str[500];
00674        sprintf(str,"%sHcalDetDiagLEDData_run%06i_%i.root",OutputFilePath.c_str(),run_number,dataset_seq_number);
00675        TFile *theFile = new TFile(str, "RECREATE");
00676        if(!theFile->IsOpen()) return;
00677        theFile->cd();
00678        sprintf(str,"%d",run_number); TObjString run(str);    run.Write("run number");
00679        sprintf(str,"%d",ievt_);      TObjString events(str); events.Write("Total events processed");
00680        
00681        TTree *tree   =new TTree("HCAL LED data","HCAL LED data");
00682        if(tree==0)   return;
00683        tree->Branch("Subdet",   &Subdet,         "Subdet/C");
00684        tree->Branch("eta",      &Eta,            "Eta/I");
00685        tree->Branch("phi",      &Phi,            "Phi/I");
00686        tree->Branch("depth",    &Depth,          "Depth/I");
00687        tree->Branch("statistic",&Statistic,      "Statistic/I");
00688        tree->Branch("status",   &Status,         "Status/I");
00689        tree->Branch("led",      &led,            "led/D");
00690        tree->Branch("rms",      &rms,            "rms/D");
00691        tree->Branch("time",     &time,           "time/D");
00692        tree->Branch("time_rms", &time_rms,       "time_rms/D");
00693        sprintf(Subdet,"HB");
00694        for(int eta=-16;eta<=16;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=2;depth++){
00695           if((Statistic=hb_data[eta+42][phi-1][depth-1].get_statistics())>100){
00696              Eta=eta; Phi=phi; Depth=depth;
00697              Status=hb_data[eta+42][phi-1][depth-1].get_status();
00698              hb_data[eta+42][phi-1][depth-1].get_average_led(&led,&rms);
00699              hb_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00700              tree->Fill();
00701           }
00702        } 
00703        sprintf(Subdet,"HE");
00704        for(int eta=-29;eta<=29;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=3;depth++){
00705          if((Statistic=he_data[eta+42][phi-1][depth-1].get_statistics())>100){
00706             Eta=eta; Phi=phi; Depth=depth;
00707             Status=he_data[eta+42][phi-1][depth-1].get_status();
00708             he_data[eta+42][phi-1][depth-1].get_average_led(&led,&rms);
00709             he_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00710             tree->Fill();
00711          }
00712        } 
00713        sprintf(Subdet,"HO");
00714        for(int eta=-15;eta<=15;eta++) for(int phi=1;phi<=72;phi++) for(int depth=4;depth<=4;depth++){
00715          if((Statistic=ho_data[eta+42][phi-1][depth-1].get_statistics())>100){
00716              Eta=eta; Phi=phi; Depth=depth;
00717              Status=ho_data[eta+42][phi-1][depth-1].get_status();
00718              ho_data[eta+42][phi-1][depth-1].get_average_led(&led,&rms);
00719              ho_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00720              tree->Fill();
00721          }
00722        } 
00723        sprintf(Subdet,"HF");
00724        for(int eta=-42;eta<=42;eta++) for(int phi=1;phi<=72;phi++) for(int depth=1;depth<=2;depth++){
00725          if((Statistic=hf_data[eta+42][phi-1][depth-1].get_statistics())>100){
00726              Eta=eta; Phi=phi; Depth=depth;
00727              Status=hf_data[eta+42][phi-1][depth-1].get_status();
00728              hf_data[eta+42][phi-1][depth-1].get_average_led(&led,&rms);
00729              hf_data[eta+42][phi-1][depth-1].get_average_time(&time,&time_rms);
00730              tree->Fill();
00731          }
00732        }
00733        sprintf(Subdet,"CALIB_HB");
00734        for(int eta=-1;eta<=1;eta++) for(int phi=1;phi<=72;phi++){
00735           if((calib_data[1][eta+2][phi-1].get_statistics())>100){
00736              Eta=eta; Phi=phi; Depth=0;
00737              Status=calib_data[1][eta+2][phi-1].get_status();
00738              calib_data[1][eta+2][phi-1].get_average_led(&led,&rms);
00739              calib_data[1][eta+2][phi-1].get_average_time(&time,&time_rms);
00740              tree->Fill();
00741           }
00742        } 
00743        sprintf(Subdet,"CALIB_HE");
00744        for(int eta=-1;eta<=1;eta++) for(int phi=1;phi<=72;phi++){
00745           if((calib_data[2][eta+2][phi-1].get_statistics())>100){
00746              Eta=eta; Phi=phi; Depth=0;
00747              Status=calib_data[2][eta+2][phi-1].get_status();
00748              calib_data[2][eta+2][phi-1].get_average_led(&led,&rms);
00749              calib_data[2][eta+2][phi-1].get_average_time(&time,&time_rms);
00750              tree->Fill();
00751           }
00752        } 
00753        sprintf(Subdet,"CALIB_HO");
00754        for(int eta=-2;eta<=2;eta++) for(int phi=1;phi<=72;phi++){
00755           if((calib_data[3][eta+2][phi-1].get_statistics())>100){
00756              Eta=eta; Phi=phi; Depth=0;
00757              Status=calib_data[3][eta+2][phi-1].get_status();
00758              calib_data[3][eta+2][phi-1].get_average_led(&led,&rms);
00759              calib_data[3][eta+2][phi-1].get_average_time(&time,&time_rms);
00760              tree->Fill();
00761           }
00762        } 
00763        sprintf(Subdet,"CALIB_HF");
00764        for(int eta=-2;eta<=2;eta++) for(int phi=1;phi<=72;phi++){
00765           if((calib_data[4][eta+2][phi-1].get_statistics())>100){
00766              Eta=eta; Phi=phi; Depth=0;
00767              Status=calib_data[4][eta+2][phi-1].get_status();
00768              calib_data[4][eta+2][phi-1].get_average_led(&led,&rms);
00769              calib_data[4][eta+2][phi-1].get_average_time(&time,&time_rms);
00770              tree->Fill();
00771           }
00772        } 
00773        theFile->Write();
00774        theFile->Close();
00775        dataset_seq_number++;
00776 }
00777 
00778 void HcalDetDiagLEDMonitor::LoadReference(){
00779 double led,rms;
00780 int Eta,Phi,Depth;
00781 char subdet[10];
00782 TFile *f;
00783    if(gSystem->AccessPathName(ReferenceData.c_str())) return;
00784    f = new TFile(ReferenceData.c_str(),"READ");
00785    if(!f->IsOpen()) return ;
00786    TObjString *STR=(TObjString *)f->Get("run number");
00787    if(STR){ std::string Ref(STR->String()); ReferenceRun=Ref;}
00788    TTree*  t=(TTree*)f->Get("HCAL LED data");
00789    if(!t) return;
00790    t->SetBranchAddress("Subdet",   subdet);
00791    t->SetBranchAddress("eta",      &Eta);
00792    t->SetBranchAddress("phi",      &Phi);
00793    t->SetBranchAddress("depth",    &Depth);
00794    t->SetBranchAddress("led",      &led);
00795    t->SetBranchAddress("rms",      &rms);
00796    for(int ievt=0;ievt<t->GetEntries();ievt++){
00797 
00798      t->GetEntry(ievt);
00799      if(strcmp(subdet,"HB")==0) hb_data[Eta+42][Phi-1][Depth-1].set_reference(led,rms);
00800      if(strcmp(subdet,"HE")==0) he_data[Eta+42][Phi-1][Depth-1].set_reference(led,rms);
00801      if(strcmp(subdet,"HO")==0) ho_data[Eta+42][Phi-1][Depth-1].set_reference(led,rms);
00802      if(strcmp(subdet,"HF")==0) hf_data[Eta+42][Phi-1][Depth-1].set_reference(led,rms);
00803      if(strcmp(subdet,"CALIB_HB")==0) calib_data[1][Eta+2][Phi-1].set_reference(led,rms);
00804      if(strcmp(subdet,"CALIB_HE")==0) calib_data[2][Eta+2][Phi-1].set_reference(led,rms);
00805      if(strcmp(subdet,"CALIB_HO")==0) calib_data[3][Eta+2][Phi-1].set_reference(led,rms);
00806      if(strcmp(subdet,"CALIB_HF")==0) calib_data[4][Eta+2][Phi-1].set_reference(led,rms);
00807    }
00808    f->Close();
00809    IsReference=true;
00810 } 
00811 void HcalDetDiagLEDMonitor::CheckStatus(){
00812    for(int i=0;i<4;i++){
00813       ChannelStatusMissingChannels->depth[i]->Reset();
00814       ChannelStatusUnstableChannels->depth[i]->Reset();
00815       ChannelStatusUnstableLEDsignal->depth[i]->Reset();
00816       ChannelStatusLEDMean->depth[i]->Reset();
00817       ChannelStatusLEDRMS->depth[i]->Reset();
00818       ChannelStatusTimeMean->depth[i]->Reset();
00819       ChannelStatusTimeRMS->depth[i]->Reset();
00820    }
00821   
00822    std::vector <HcalElectronicsId> AllElIds = emap.allElectronicsIdPrecision();
00823    for (std::vector <HcalElectronicsId>::iterator eid = AllElIds.begin(); eid != AllElIds.end(); eid++) {
00824       DetId detid=emap.lookup(*eid);
00825       if (detid.det()!=DetId::Hcal) continue;
00826       HcalGenericDetId gid(emap.lookup(*eid));
00827       if(!(!(gid.null()) && 
00828             (gid.genericSubdet()==HcalGenericDetId::HcalGenBarrel ||
00829              gid.genericSubdet()==HcalGenericDetId::HcalGenEndcap  ||
00830              gid.genericSubdet()==HcalGenericDetId::HcalGenForward ||
00831              gid.genericSubdet()==HcalGenericDetId::HcalGenOuter))) continue;
00832       int eta=0,phi=0,depth=0;
00833       
00834       HcalDetId hid(detid);
00835       eta=hid.ieta();
00836       phi=hid.iphi();
00837       depth=hid.depth(); 
00838       
00839       double AVE_TIME=Time->getMean();
00840       if(detid.subdetId()==HcalBarrel){
00841          int stat=hb_data[eta+42][phi-1][depth-1].get_statistics()+
00842                    hb_data[eta+42][phi-1][depth-1].get_overflow()+hb_data[eta+42][phi-1][depth-1].get_undeflow();
00843          if(stat==0){ 
00844              fill_channel_status("HB",eta,phi,depth,1,1); 
00845              hb_data[eta+42][phi-1][depth-1].change_status(1); 
00846          }
00847          if(stat>0 && stat!=(ievt_)){ 
00848              fill_channel_status("HB",eta,phi,depth,2,(double)stat/(double)(ievt_)); 
00849              hb_data[eta+42][phi-1][depth-1].change_status(2); 
00850          }
00851          if(hb_data[eta+42][phi-1][depth-1].get_statistics()>100){ 
00852              double ave=0;
00853              double rms=0;
00854              hb_data[eta+42][phi-1][depth-1].get_average_time(&ave,&rms);
00855              if((AVE_TIME-ave)>0.75 || (AVE_TIME-ave)<-0.75){
00856                 fill_channel_status("HB",eta,phi,depth,6,AVE_TIME-ave); 
00857                 hb_data[eta+42][phi-1][depth-1].change_status(8); 
00858              }
00859          }  
00860          stat=hb_data[eta+42][phi-1][depth-1].get_undeflow();     
00861          if(stat>0){ 
00862              fill_channel_status("HB",eta,phi,depth,3,(double)stat/(double)(ievt_)); 
00863              hb_data[eta+42][phi-1][depth-1].change_status(4); 
00864          }    
00865       } 
00866       if(detid.subdetId()==HcalEndcap){
00867          int stat=he_data[eta+42][phi-1][depth-1].get_statistics()+
00868                    he_data[eta+42][phi-1][depth-1].get_overflow()+he_data[eta+42][phi-1][depth-1].get_undeflow();
00869          if(stat==0){ 
00870              fill_channel_status("HE",eta,phi,depth,1,1); 
00871              he_data[eta+42][phi-1][depth-1].change_status(1); 
00872          }
00873          if(stat>0 && stat!=(ievt_)){ 
00874              fill_channel_status("HE",eta,phi,depth,2,(double)stat/(double)(ievt_)); 
00875              he_data[eta+42][phi-1][depth-1].change_status(2); 
00876          }
00877          if(he_data[eta+42][phi-1][depth-1].get_statistics()>100){ 
00878              double ave=0;
00879              double rms=0;
00880              he_data[eta+42][phi-1][depth-1].get_average_time(&ave,&rms);
00881              if((AVE_TIME-ave)>0.75 || (AVE_TIME-ave)<-0.75){ 
00882                 fill_channel_status("HE",eta,phi,depth,6,AVE_TIME-ave); 
00883                 he_data[eta+42][phi-1][depth-1].change_status(8); 
00884              }  
00885          }  
00886          stat=he_data[eta+42][phi-1][depth-1].get_undeflow();     
00887          if(stat>0){ 
00888              fill_channel_status("HE",eta,phi,depth,3,(double)stat/(double)(ievt_)); 
00889              he_data[eta+42][phi-1][depth-1].change_status(4); 
00890          }  
00891       } 
00892       if(detid.subdetId()==HcalOuter){
00893          int stat=ho_data[eta+42][phi-1][depth-1].get_statistics()+
00894                    ho_data[eta+42][phi-1][depth-1].get_overflow()+ho_data[eta+42][phi-1][depth-1].get_undeflow();
00895          if(stat==0){ 
00896              fill_channel_status("HO",eta,phi,depth,1,1); 
00897              ho_data[eta+42][phi-1][depth-1].change_status(1); 
00898          }
00899          if(stat>0 && stat!=(ievt_)){ 
00900              fill_channel_status("HO",eta,phi,depth,2,(double)stat/(double)(ievt_)); 
00901              ho_data[eta+42][phi-1][depth-1].change_status(2); 
00902          }
00903          if(ho_data[eta+42][phi-1][depth-1].get_statistics()>100){ 
00904              double ave=0;
00905              double rms=0;
00906              ho_data[eta+42][phi-1][depth-1].get_average_time(&ave,&rms);
00907              if((AVE_TIME-ave)>0.75 || (AVE_TIME-ave)<-0.75){
00908                 fill_channel_status("HO",eta,phi,depth,6,AVE_TIME-ave); 
00909                 ho_data[eta+42][phi-1][depth-1].change_status(8);
00910              } 
00911          }  
00912          stat=ho_data[eta+42][phi-1][depth-1].get_undeflow();     
00913          if(stat>0){ 
00914              fill_channel_status("HO",eta,phi,depth,3,(double)stat/(double)(ievt_)); 
00915              ho_data[eta+42][phi-1][depth-1].change_status(4); 
00916          }  
00917       } 
00918       if(detid.subdetId()==HcalForward){
00919          AVE_TIME=TimeHF->getMean();
00920          int stat=hf_data[eta+42][phi-1][depth-1].get_statistics()+
00921                    hf_data[eta+42][phi-1][depth-1].get_overflow()+hf_data[eta+42][phi-1][depth-1].get_undeflow();
00922          if(stat==0){ 
00923              fill_channel_status("HF",eta,phi,depth,1,1); 
00924              hf_data[eta+42][phi-1][depth-1].change_status(1); 
00925          }
00926          if(stat>0 && stat!=(ievt_)){ 
00927              fill_channel_status("HF",eta,phi,depth,2,(double)stat/(double)(ievt_)); 
00928              hf_data[eta+42][phi-1][depth-1].change_status(2); 
00929          }
00930          if(hf_data[eta+42][phi-1][depth-1].get_statistics()>100){ 
00931              double ave=0;
00932              double rms=0;
00933              hf_data[eta+42][phi-1][depth-1].get_average_time(&ave,&rms);
00934              if((AVE_TIME-ave)>0.75 || (AVE_TIME-ave)<-0.75){
00935                 fill_channel_status("HF",eta,phi,depth,6,AVE_TIME-ave); 
00936                 hf_data[eta+42][phi-1][depth-1].change_status(8);
00937              } 
00938          }  
00939          stat=hf_data[eta+42][phi-1][depth-1].get_undeflow();     
00940          if(stat>0){ 
00941              fill_channel_status("HF",eta,phi,depth,3,(double)stat/(double)(ievt_)); 
00942              hf_data[eta+42][phi-1][depth-1].change_status(4); 
00943          }  
00944       } 
00945    }
00946 }
00947 void HcalDetDiagLEDMonitor::fill_energy(std::string subdet,int eta,int phi,int depth,double e,int type){ 
00948   int subdetval=-1;
00949   if (subdet.compare("HB")==0) subdetval=(int)HcalBarrel;
00950   else if (subdet.compare("HE")==0) subdetval=(int)HcalEndcap;
00951   else if (subdet.compare("HO")==0) subdetval=(int)HcalOuter;
00952   else if (subdet.compare("HF")==0) subdetval=(int)HcalForward;
00953   else return;
00954 
00955   int ietabin=CalcEtaBin(subdetval, eta, depth)+1;
00956   if(type==1) ChannelsLEDEnergy->depth[depth-1]   ->setBinContent(ietabin,phi,e);
00957   else if(type==2) ChannelsLEDEnergyRef->depth[depth-1]->setBinContent(ietabin,phi,e);
00958 }
00959 
00960 double HcalDetDiagLEDMonitor::get_energy(std::string subdet,int eta,int phi,int depth,int type){
00961   int subdetval=-1;
00962   if (subdet.compare("HB")==0) subdetval=(int)HcalBarrel;
00963   else if (subdet.compare("HE")==0) subdetval=(int)HcalEndcap;
00964   else if (subdet.compare("HO")==0) subdetval=(int)HcalOuter;
00965   else if (subdet.compare("HF")==0) subdetval=(int)HcalForward;
00966   else return -1.0;
00967 
00968   int ietabin=CalcEtaBin(subdetval, eta, depth)+1;
00969   if(type==1) return ChannelsLEDEnergy->depth[depth-1]  ->getBinContent(ietabin, phi);
00970   else if(type==2) return ChannelsLEDEnergyRef->depth[depth-1] ->getBinContent(ietabin,phi);
00971   return -1.0;
00972 }
00973 
00974 void HcalDetDiagLEDMonitor::fill_channel_status(std::string subdet,int eta,int phi,int depth,int type,double status){
00975   int subdetval=-1;
00976   if (subdet.compare("HB")==0) subdetval=(int)HcalBarrel;
00977   else if (subdet.compare("HE")==0) subdetval=(int)HcalEndcap;
00978   else if (subdet.compare("HO")==0) subdetval=(int)HcalOuter;
00979   else if (subdet.compare("HF")==0) subdetval=(int)HcalForward;
00980   else return;
00981   int ietabin=CalcEtaBin(subdetval, eta, depth)+1;
00982 
00983    if(type==1) ChannelStatusMissingChannels->depth[depth-1]  ->setBinContent(ietabin,phi,status);
00984    if(type==2) ChannelStatusUnstableChannels->depth[depth-1] ->setBinContent(ietabin,phi,status);
00985    if(type==3) ChannelStatusUnstableLEDsignal->depth[depth-1]->setBinContent(ietabin,phi,status);
00986    if(type==4) ChannelStatusLEDMean->depth[depth-1]          ->setBinContent(ietabin,phi,status);
00987    if(type==5) ChannelStatusLEDRMS->depth[depth-1]           ->setBinContent(ietabin,phi,status);
00988    if(type==6) ChannelStatusTimeMean->depth[depth-1]         ->setBinContent(ietabin,phi,status);
00989    if(type==7) ChannelStatusTimeRMS->depth[depth-1]          ->setBinContent(ietabin,phi,status);
00990 }
00991 void HcalDetDiagLEDMonitor::endRun(const edm::Run& run, const edm::EventSetup& c){   
00992    if(ievt_>=100){
00993       fillHistos();
00994       CheckStatus();
00995       SaveReference(); 
00996    }   
00997 } 
00998 DEFINE_FWK_MODULE (HcalDetDiagLEDMonitor);
00999