CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/DQM/HcalMonitorTasks/src/HcalCaloTowerMonitor.cc

Go to the documentation of this file.
00001 #include "DQM/HcalMonitorTasks/interface/HcalCaloTowerMonitor.h"
00002 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
00003 
00004 HcalCaloTowerMonitor::HcalCaloTowerMonitor(){}
00005 
00006 HcalCaloTowerMonitor::~HcalCaloTowerMonitor() {}
00007 
00008 void HcalCaloTowerMonitor::setup(const edm::ParameterSet& ps, DQMStore* dbe)
00009 {
00010   HcalBaseMonitor::setup(ps,dbe);
00011 
00012   baseFolder_ = rootFolder_+"CaloTowerMonitor";
00013   etaMax_ = ps.getUntrackedParameter<double>("MaxEta", 41.5);
00014   etaMin_ = ps.getUntrackedParameter<double>("MinEta", -41.5);
00015   etaBins_ = (int)(etaMax_ - etaMin_);
00016   std::cout << "CaloTower Monitor eta min/max set to " << etaMin_ << "/" << etaMax_ << std::endl;
00017 
00018   phiMax_ = ps.getUntrackedParameter<double>("MaxPhi", 73);
00019   phiMin_ = ps.getUntrackedParameter<double>("MinPhi", 0);
00020   phiBins_ = (int)(phiMax_ - phiMin_);
00021   std::cout << "CaloTower Monitor phi min/max set to " << phiMin_ << "/" << phiMax_ << std::endl;
00022 
00023   ievt_=0;
00024   // book histograms
00025   if (m_dbe)
00026     {
00027       m_dbe->setCurrentFolder(baseFolder_); 
00028       meEVT_ = m_dbe->bookInt("CaloTower Event Number");     
00029       meEVT_->Fill(ievt_); 
00030 
00031       // Calo Tower occupancy
00032       caloTowerOccMap=m_dbe->book2D("CaloTowerOccupancyMap", "Calo Tower Occupancy Map",
00033                                     etaBins_,etaMin_,etaMax_,
00034                                     phiBins_,phiMin_,phiMax_);
00035       // Calo Tower energy
00036       caloTowerEnergyMap=m_dbe->book2D("CaloTowerEnergyMap", "Calo Tower Energy Map",
00037                                        etaBins_,etaMin_,etaMax_,
00038                                        phiBins_,phiMin_,phiMax_);
00039 
00040       caloTowerTime = m_dbe->book1D("CaloTowerTime", "Calo Tower Time",
00041                                     10,0,10);
00042       caloTowerEnergy = m_dbe->book1D("CalotowerEnergy", "Calo Tower Energy",
00043                                       100,0,100);
00044       caloTowerMeanEnergyEta = m_dbe->bookProfile("CaloTowerMeanEnergyEta","Calo tower Mean Energy vs. Eta",
00045                                                   etaBins_,etaMin_,etaMax_,100,0,100);
00046 
00047       // Make plots of HCAL contributions to tower
00048       m_dbe->setCurrentFolder(baseFolder_+"/"+"HCAL");
00049       hcalOccMap=m_dbe->book2D("HcalOccupancyMap", "Calo Tower Occupancy Map",
00050                                etaBins_,etaMin_,etaMax_,
00051                                phiBins_,phiMin_,phiMax_);
00052       // Calo Tower energy
00053       hcalEnergyMap=m_dbe->book2D("HcalEnergyMap", "Calo Tower Energy Map",
00054                                   etaBins_,etaMin_,etaMax_,
00055                                   phiBins_,phiMin_,phiMax_);
00056 
00057       hcalTime = m_dbe->book1D("HcalTime", "Calo Tower Time",
00058                                     10,0,10);
00059       hcalEnergy = m_dbe->book1D("HcalEnergy", "Calo Tower Energy",
00060                                       100,0,100);
00061       hcalMeanEnergyEta = m_dbe->bookProfile("HcalMeanEnergyEta",
00062                                              "Calo tower Mean Energy vs. Eta",
00063                                              etaBins_,etaMin_,etaMax_,
00064                                              100,0,100);
00065       
00066 
00067       m_dbe->setCurrentFolder(baseFolder_+"/"+"ECAL");
00068       ecalOccMap=m_dbe->book2D("EcalOccupancyMap", "Calo Tower Occupancy Map",
00069                                etaBins_,etaMin_,etaMax_,
00070                                phiBins_,phiMin_,phiMax_);
00071       // Calo Tower energy
00072       ecalEnergyMap=m_dbe->book2D("EcalEnergyMap", "Calo Tower Energy Map",
00073                                   etaBins_,etaMin_,etaMax_,
00074                                   phiBins_,phiMin_,phiMax_);
00075 
00076       ecalTime = m_dbe->book1D("EcalTime", "Calo Tower Time",
00077                                10,0,10);
00078       ecalEnergy = m_dbe->book1D("EcalEnergy", "Calo Tower Energy",
00079                                  100,0,100);
00080       ecalMeanEnergyEta = m_dbe->bookProfile("EcalMeanEnergyEta",
00081                                              "Calo tower Mean Energy vs. Eta",
00082                                              etaBins_,etaMin_,etaMax_,
00083                                              100,0,100);
00084 
00085       // Comparison Plots
00086       m_dbe->setCurrentFolder(baseFolder_+"/"+"ComparisonPlots");
00087       time_HcalvsEcal=m_dbe->book2D("HcalvsEcalTime",
00088                                     "Hcal Time vs. Ecal time",
00089                                     10,0,10,10,0,10);
00090       time_CaloTowervsEcal=m_dbe->book2D("CaloTowervsEcalTime",
00091                                          "Calotower Time vs. Ecal time",
00092                                          10,0,10,10,0,10);
00093       time_CaloTowervsHcal=m_dbe->book2D("CaloTowervsHcalTime",
00094                                          "CaloTower Time vs. Hcal time",
00095                                          10,0,10,10,0,10);
00096       energy_HcalvsEcal=m_dbe->book2D("HcalvsEcalEnergy",
00097                                     "Hcal Energy vs. Ecal energy",
00098                                     100,0,100,100,0,100);
00099       energy_CaloTowervsEcal=m_dbe->book2D("CaloTowervsEcalEnergy",
00100                                          "Calotower Energy vs. Ecal energy",
00101                                          100,0,100,100,0,100);
00102       energy_CaloTowervsHcal=m_dbe->book2D("CaloTowervsHcalEnergy",
00103                                          "CaloTower Energy vs. Hcal energy",
00104                                          100,0,100,100,0,100);
00105     } // if (m_dbe)
00106   
00107 } //void HcalCaloTowerMonitor::setup(...)
00108 
00109 
00110 void HcalCaloTowerMonitor::processEvent(const CaloTowerCollection& calotower)
00111 {
00112   // fill histograms for each event
00113   if(!m_dbe)  
00114     { 
00115       if(fVerbosity) std::cout <<"<HcalCaloTowerMonitor::processEvent>    DQMStore not instantiated!!!\n"; 
00116       return; 
00117     }  // if(!m_dbe)
00118   
00119   if (fVerbosity) std::cout <<"Processing calotower"<<std::endl;
00120 
00121 
00122   // Store values of hcal, energy for forming averages at each eta
00123   float hcalenergy[90]={0.};
00124   float ecalenergy[90]={0.};
00125   int etacounts[90]={0};
00126 
00127   for (CaloTowerCollection::const_iterator it=calotower.begin();
00128        it!=calotower.end();++it)
00129     {
00130       //caloTowers give eta and phi; need to map to ieta,iphi
00131       // phi is mapped into 72 equally spaced segments (2pi/72=.087 radians each), starting at iphi=1
00132       // at ieta=21 (eta=1.74), iphi spaced in 10 degree increments 
00133       // at ieta=40 (eta=4.716), iphi spaced in 20 degree increments
00134       // However, segmentation values are similar
00135       // (iphi=1,2,3,4,5 for 5 degree segments,
00136       //  iphi 1,3,5,... for 10 degree segments, and iphi=3,7,... for 20 deg.)
00137       int iphi;
00138       double eta=it->eta();
00139       double phi=it->phi();
00140       // if phi<0; shift it by 2pi (so that phi runs from 0->2pi, rather than -pi->pi)
00141       if (phi<0)
00142         phi+=2*3.14159265359;
00143       iphi = int(phi/.087)+1;
00144       if (fabs(eta)>1.74)
00145         {
00146           if (fabs(eta)<=4.716)
00147             iphi=(iphi-1)/2*2+1;
00148           else
00149             iphi=(iphi-1)/4*4+3;
00150           
00151         }
00152       //ieta is more complicated -- it runs in segments of .087 for the first 20 segments, then the segmentation becomes non-uniform
00153       int ieta=getIeta(eta);
00154 
00155       // Add 42 to the indices.  Only index[1...] will have values within [HB...HF].  Index[0] and [84] store values outside the expected calorimeter range.
00156       hcalenergy[ieta+42]+=it->hadEnergy();
00157       ecalenergy[ieta+42]+=it->emEnergy();
00158       etacounts[ieta+42]++;
00159 
00160       /*
00161       std::cout <<"CALOTOWER eta = \t"<<eta<<"\tphi = "<<it->phi()<<std::endl;
00162       std::cout <<"\t ieta = \t"<<ieta<<"\tiphi = "<<iphi<<std::endl;
00163       */
00164 
00165       // Fill histograms
00166 
00167       // calotowers
00168       caloTowerOccMap->Fill(ieta,iphi);
00169       caloTowerEnergyMap->Fill(ieta,iphi,it->energy());
00170       caloTowerTime->Fill((it->ecalTime()+it->hcalTime())/2.);
00171       caloTowerEnergy->Fill(it->energy());
00172       //std::cout <<"calotower times: "<<it->ecalTime()<<"  "<<it->hcalTime()<<std::endlx;
00173 
00174       //hcal
00175       hcalOccMap->Fill(ieta,iphi);
00176       hcalEnergyMap->Fill(ieta,iphi,it->hadEnergy());
00177       hcalTime->Fill(it->hcalTime());
00178       hcalEnergy->Fill(it->hadEnergy());
00179 
00180       //ecal
00181       ecalOccMap->Fill(ieta,iphi);
00182       ecalEnergyMap->Fill(ieta,iphi,it->emEnergy());
00183       ecalTime->Fill(it->ecalTime());
00184       ecalEnergy->Fill(it->emEnergy());
00185 
00186       time_HcalvsEcal->Fill(it->ecalTime(),it->hcalTime());
00187       energy_HcalvsEcal->Fill(it->emEnergy(),it->hadEnergy());
00188     } // for CaloTowerCollection::const_iterator it...
00189 
00190   // Now form average for each eta bin;
00191   for (int i=0;i<90;++i)
00192     {
00193       if (etacounts[i]==0) continue; 
00194       caloTowerMeanEnergyEta->Fill(i-42, 1.*(hcalenergy[i]+ecalenergy[i])/etacounts[i]);
00195       hcalMeanEnergyEta->Fill(i-42, 1.*hcalenergy[i]/etacounts[i]);
00196       ecalMeanEnergyEta->Fill(i-42,1.*ecalenergy[i]/etacounts[i]);
00197     }
00198 }
00199 
00200 int HcalCaloTowerMonitor::getIeta(double eta)
00201 {
00202   // return Ieta index given physical eta value
00203 
00204   double myeta=fabs(eta);
00205   int neg=int(eta/myeta); // eta can be negative or positive
00206 
00207   if (myeta<=1.740) // first 20 bins are spaced evenly in increments of .087 radians, starting at ieta=1
00208     return (1+int(myeta/.087))*neg;
00209 
00210   if (fabs(eta)<=1.83)
00211     return int(21*neg);
00212   if (fabs(eta)<=1.93)
00213     return int(22*neg);
00214   if (fabs(eta)<=2.043)
00215     return int(23*neg);
00216   if (fabs(eta)<=2.172)
00217     return int(24*neg);
00218   if (fabs(eta)<=2.322)
00219     return int(25*neg);
00220              
00221   if (fabs(eta)<=2.5)
00222     return int(26*neg);
00223   if (fabs(eta)<=2.65)
00224     return int(27*neg);
00225 
00226   // 28 and 29 aren't quite right -- eta values depend on depth
00227   if (fabs(eta)<=2.853)
00228     return int(28*neg);
00229   if (fabs(eta)<=2.964)
00230     return int(29*neg);
00231   if (fabs(eta)<=3.139)
00232     return int(30*neg);
00233                  
00234   if (fabs(eta)<=3.314)
00235     return int(31*neg);
00236   if (fabs(eta)<=3.489)
00237     return int(32*neg);
00238   if (fabs(eta)<=3.664)
00239     return int(33*neg);
00240   if (fabs(eta)<=3.839)
00241     return int(34*neg);
00242   if (fabs(eta)<=4.013)
00243     return int(35*neg);
00244          
00245   if (fabs(eta)<=4.191)
00246     return int(36*neg);
00247   if (fabs(eta)<=4.363)
00248     return int(37*neg);
00249   if (fabs(eta)<=4.538)
00250     return int(38*neg);
00251   if (fabs(eta)<=4.716)
00252     return int(39*neg);
00253   if (fabs(eta)<=4.889)
00254     return int(40*neg);
00255   if (fabs(eta)<=5.191)
00256     return int(41*neg);
00257 
00258   // Anything with large eta is outside calorimeter; skip it?
00259   // Are there ZDC-based calotowers?
00260   return 42;
00261     
00262 } // int HcalCaloTowerMonitor::getIeta(double eta)
00263 
00264 
00265 void HcalCaloTowerMonitor::reset(){}