Hcal/src/CaloTowersAnalyzer.cc

0001 #include "DQMOffline/Hcal/interface/CaloTowersAnalyzer.h"
0002 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0003 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0004 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0005 #include "FWCore/Utilities/interface/Transition.h"
0006
0007 CaloTowersAnalyzer::CaloTowersAnalyzer(edm::ParameterSet const &conf)
0008     : hcalDDDRecConstantsToken_{esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord, edm::Transition::BeginRun>()} {
0009   tok_towers_ = consumes<CaloTowerCollection>(conf.getUntrackedParameter<edm::InputTag>("CaloTowerCollectionLabel"));
0010   // DQM ROOT output
0011
0012   outputFile_ = conf.getUntrackedParameter<std::string>("outputFile", "myfile.root");
0013
0014   hcalselector_ = conf.getUntrackedParameter<std::string>("hcalselector", "all");
0015
0016   useAllHistos_ = conf.getUntrackedParameter<bool>("useAllHistos", false);
0017 }
0018
0019 void CaloTowersAnalyzer::dqmBeginRun(const edm::Run &run, const edm::EventSetup &es) {
0020   HcalDDDRecConstants const &hcons = es.getData(hcalDDDRecConstantsToken_);
0021
0022   // Get Phi segmentation from geometry, use the max phi number so that all iphi
0023   // values are included.
0024
0025   int NphiMax = hcons.getNPhi(0);
0026
0027   NphiMax = (hcons.getNPhi(1) > NphiMax ? hcons.getNPhi(1) : NphiMax);
0028   NphiMax = (hcons.getNPhi(2) > NphiMax ? hcons.getNPhi(2) : NphiMax);
0029   NphiMax = (hcons.getNPhi(3) > NphiMax ? hcons.getNPhi(3) : NphiMax);
0030
0031   // Center the iphi bins on the integers
0032   iphi_min_ = 0.5;
0033   iphi_max_ = NphiMax + 0.5;
0034   iphi_bins_ = (int)(iphi_max_ - iphi_min_);
0035
0036   // Retain classic behavior, all plots have same ieta range.
0037
0038   int iEtaMax = (hcons.getEtaRange(0).second > hcons.getEtaRange(1).second ? hcons.getEtaRange(0).second
0039                                                                            : hcons.getEtaRange(1).second);
0040   iEtaMax = (iEtaMax > hcons.getEtaRange(2).second ? iEtaMax : hcons.getEtaRange(2).second);
0041   iEtaMax = (iEtaMax > hcons.getEtaRange(3).second ? iEtaMax : hcons.getEtaRange(3).second);
0042
0043   // Give an empty bin around the subdet ieta range to make it clear that all
0044   // ieta rings have been included
0045   ieta_min_ = -iEtaMax - 1.5;
0046   ieta_max_ = iEtaMax + 1.5;
0047   ieta_bins_ = (int)(ieta_max_ - ieta_min_);
0048 }
0049
0050 void CaloTowersAnalyzer::bookHistograms(DQMStore::IBooker &ibooker,
0051                                         edm::Run const & /* iRun*/,
0052                                         edm::EventSetup const & /* iSetup */) {
0053   etaMin[0] = 0.;
0054   etaMax[0] = 1.4;
0055   etaMin[1] = 1.4;
0056   etaMax[1] = 2.9;
0057   etaMin[2] = 2.9;
0058   etaMax[2] = 5.2;
0059
0060   isub = 0;
0061   if (hcalselector_ == "HB")
0062     isub = 1;
0063   if (hcalselector_ == "HE")
0064     isub = 2;
0065   if (hcalselector_ == "HF")
0066     isub = 3;
0067
0068   if (!outputFile_.empty()) {
0069     edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will be saved to '" << outputFile_.c_str() << "'";
0070   } else {
0071     edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will NOT be saved";
0072   }
0073
0074   Char_t histo[100];
0075
0076   ibooker.setCurrentFolder("CaloTowersD/CaloTowersTask");
0077
0078   // These two histos are not drawn by our macros, but they are used
0079   // in the EndJob for norms and such so I am leaving them alone for now
0080   //-------------------------------------------------------------------------------------------
0081   sprintf(histo, "Ntowers_per_event_vs_ieta");
0082   Ntowers_vs_ieta = ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_);
0083
0084   sprintf(histo, "CaloTowersTask_map_Nentries");
0085   mapEnergy_N = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0086   //-------------------------------------------------------------------------------------------
0087
0088   // These the single pion scan histos
0089   //-------------------------------------------------------------------------------------------
0090   // The first three are not used
0091   if (useAllHistos_) {
0092     sprintf(histo, "emean_vs_ieta_E");
0093     emean_vs_ieta_E = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0094     sprintf(histo, "emean_vs_ieta_H");
0095     emean_vs_ieta_H = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0096     sprintf(histo, "emean_vs_ieta_EH");
0097     emean_vs_ieta_EH = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0098   }
0099   // These are drawn
0100   sprintf(histo, "emean_vs_ieta_E1");
0101   emean_vs_ieta_E1 = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0102   sprintf(histo, "emean_vs_ieta_H1");
0103   emean_vs_ieta_H1 = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0104   sprintf(histo, "emean_vs_ieta_EH1");
0105   emean_vs_ieta_EH1 = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, 2100, -100., 2000.);
0106   //-------------------------------------------------------------------------------------------
0107
0108   // Map energy histos are not drawn
0109   if (useAllHistos_) {
0110     sprintf(histo, "CaloTowersTask_map_energy_E");
0111     mapEnergy_E = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0112     sprintf(histo, "CaloTowersTask_map_energy_H");
0113     mapEnergy_H = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0114     sprintf(histo, "CaloTowersTask_map_energy_EH");
0115     mapEnergy_EH = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0116   }
0117
0118   // All ECAL cell histos are used
0119   // XXX: ECAL 0-25 [0-26, 26 bins]   HCAL 0-4 [0-5, 5 bins]
0120   sprintf(histo, "number_of_bad_cells_Ecal_EB");
0121   numBadCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0122   sprintf(histo, "number_of_bad_cells_Ecal_EE");
0123   numBadCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0124   sprintf(histo, "number_of_recovered_cells_Ecal_EB");
0125   numRcvCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0126   sprintf(histo, "number_of_recovered_cells_Ecal_EE");
0127   numRcvCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0128   sprintf(histo, "number_of_problematic_cells_Ecal_EB");
0129   numPrbCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0130   sprintf(histo, "number_of_problematic_cells_Ecal_EE");
0131   numPrbCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0132
0133   // Occupancy vs. ieta is drawn, occupancy map is needed to draw it
0134   sprintf(histo, "CaloTowersTask_map_occupancy");
0135   occupancy_map = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0136
0137   sprintf(histo, "CaloTowersTask_occupancy_vs_ieta");
0138   occupancy_vs_ieta = ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_);
0139
0140   if (isub == 1 || isub == 0) {
0141     // All cell histos are used
0142     sprintf(histo, "number_of_bad_cells_Hcal_HB");
0143     numBadCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0144     sprintf(histo, "number_of_recovered_cells_Hcal_HB");
0145     numRcvCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0146     sprintf(histo, "number_of_problematic_cells_Hcal_HB");
0147     numPrbCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0148
0149     // These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL,
0150     // MET and SET
0151     //-------------------------------------------------------------------------------------------
0152     sprintf(histo, "CaloTowersTask_energy_HCAL_HB");
0153     meEnergyHcal_HB = ibooker.book1D(histo, histo, 4100, -200, 8000);
0154
0155     sprintf(histo, "CaloTowersTask_energy_ECAL_HB");
0156     meEnergyEcal_HB = ibooker.book1D(histo, histo, 3100, -200, 6000);
0157
0158     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HB");
0159     meNumFiredTowers_HB = ibooker.book1D(histo, histo, 1000, 0, 2000);
0160
0161     sprintf(histo, "CaloTowersTask_MET_HB");
0162     MET_HB = ibooker.book1D(histo, histo, 3000, 0., 3000.);
0163
0164     sprintf(histo, "CaloTowersTask_SET_HB");
0165     SET_HB = ibooker.book1D(histo, histo, 8000, 0., 8000.);
0166     //-------------------------------------------------------------------------------------------
0167
0168     // Timing histos and profiles -- all six are necessary
0169     //-------------------------------------------------------------------------------------------
0170     sprintf(histo, "CaloTowersTask_EM_Timing_HB");
0171     emTiming_HB = ibooker.book1D(histo, histo, 110, -120., 100.);
0172
0173     sprintf(histo, "CaloTowersTask_HAD_Timing_HB");
0174     hadTiming_HB = ibooker.book1D(histo, histo, 70, -48., 92.);
0175
0176     // Energy-Timing histos are divided into low, medium and high to reduce
0177     // memory usage EM
0178     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_Low_HB");
0179     emEnergyTiming_Low_HB = ibooker.book2D(histo, histo, 40, 0., 40., 110, -120., 100.);
0180
0181     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HB");
0182     emEnergyTiming_HB = ibooker.book2D(histo, histo, 200, 0., 400., 110, -120., 100.);
0183
0184     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_High_HB");
0185     emEnergyTiming_High_HB = ibooker.book2D(histo, histo, 200, 0., 3000., 110, -120., 100.);
0186
0187     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_Low_HB");
0188     emEnergyTiming_profile_Low_HB = ibooker.bookProfile(histo, histo, 40, 0., 40., 110, -120., 100.);
0189
0190     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HB");
0191     emEnergyTiming_profile_HB = ibooker.bookProfile(histo, histo, 200, 0., 400., 110, -120., 100.);
0192
0193     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_High_HB");
0194     emEnergyTiming_profile_High_HB = ibooker.bookProfile(histo, histo, 200, 0., 3000., 110, -120., 100.);
0195
0196     // HAD
0197     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HB");
0198     hadEnergyTiming_Low_HB = ibooker.book2D(histo, histo, 40, 0., 40., 70, -48., 92.);
0199
0200     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HB");
0201     hadEnergyTiming_HB = ibooker.book2D(histo, histo, 100, 0., 200., 70, -48., 92.);
0202
0203     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_High_HB");
0204     hadEnergyTiming_High_HB = ibooker.book2D(histo, histo, 300, 0., 3000., 70, -48., 92.);
0205
0206     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HB");
0207     hadEnergyTiming_profile_Low_HB = ibooker.bookProfile(histo, histo, 40, 0., 40., 70, -48., 92.);
0208
0209     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HB");
0210     hadEnergyTiming_profile_HB = ibooker.bookProfile(histo, histo, 100, 0., 200., 70, -48., 92.);
0211
0212     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_High_HB");
0213     hadEnergyTiming_profile_High_HB = ibooker.bookProfile(histo, histo, 300, 0., 3000., 70, -48., 92.);
0214     //-------------------------------------------------------------------------------------------
0215
0216     sprintf(histo, "CaloTowersTask_Iphi_HCAL_component_of_tower_HBP");
0217     meIphiHcalTower_HBP = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0218     sprintf(histo, "CaloTowersTask_Iphi_HCAL_component_of_tower_HBM");
0219     meIphiHcalTower_HBM = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0220
0221     //-------------------------------------------------------------------------------------------
0222
0223     // Everything else is not drawn
0224     if (useAllHistos_) {
0225       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HB");
0226       meEnergyHcalvsEcal_HB = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0227
0228       sprintf(histo, "CaloTowersTask_energy_OUTER_HB");
0229       meEnergyHO_HB = ibooker.book1D(histo, histo, 1640, -200, 8000);
0230
0231       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HB");
0232       meEnergyEcalTower_HB = ibooker.book1D(histo, histo, 440, -200, 2000);
0233
0234       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HB");
0235       meEnergyHcalTower_HB = ibooker.book1D(histo, histo, 440, -200, 2000);
0236
0237       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HB");
0238       meTotEnergy_HB = ibooker.book1D(histo, histo, 400, 0., 2000.);
0239
0240       sprintf(histo, "CaloTowersTask_map_energy_HB");
0241       mapEnergy_HB = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0242       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HB");
0243       mapEnergyHcal_HB =
0244           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0245       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HB");
0246       mapEnergyEcal_HB =
0247           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0248
0249       sprintf(histo, "CaloTowersTask_phi_MET_HB");
0250       phiMET_HB = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0251     }
0252   }
0253
0254   if (isub == 2 || isub == 0) {
0255     // All cell histos are used
0256     sprintf(histo, "number_of_bad_cells_Hcal_HE");
0257     numBadCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0258     sprintf(histo, "number_of_recovered_cells_Hcal_HE");
0259     numRcvCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0260     sprintf(histo, "number_of_problematic_cells_Hcal_HE");
0261     numPrbCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0262
0263     // These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL,
0264     // MET and SET
0265     //-------------------------------------------------------------------------------------------
0266     sprintf(histo, "CaloTowersTask_energy_HCAL_HE");
0267     meEnergyHcal_HE = ibooker.book1D(histo, histo, 1240, -200, 6000);
0268
0269     sprintf(histo, "CaloTowersTask_energy_ECAL_HE");
0270     meEnergyEcal_HE = ibooker.book1D(histo, histo, 1240, -200, 6000);
0271
0272     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HE");
0273     meNumFiredTowers_HE = ibooker.book1D(histo, histo, 1000, 0, 2000);
0274
0275     sprintf(histo, "CaloTowersTask_MET_HE");
0276     MET_HE = ibooker.book1D(histo, histo, 1000, 0., 1000.);
0277
0278     sprintf(histo, "CaloTowersTask_SET_HE");
0279     SET_HE = ibooker.book1D(histo, histo, 2000, 0., 2000.);
0280     //-------------------------------------------------------------------------------------------
0281
0282     // Timing histos and profiles -- all six are necessary
0283     //-------------------------------------------------------------------------------------------
0284     sprintf(histo, "CaloTowersTask_EM_Timing_HE");
0285     emTiming_HE = ibooker.book1D(histo, histo, 110, -120., 100.);
0286
0287     sprintf(histo, "CaloTowersTask_HAD_Timing_HE");
0288     hadTiming_HE = ibooker.book1D(histo, histo, 70, -48., 92.);
0289
0290     // Energy-Timing histos are divided into low and normal to reduce memory
0291     // usage EM
0292     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_Low_HE");
0293     emEnergyTiming_Low_HE = ibooker.book2D(histo, histo, 160, 0., 160., 110, -120., 100.);
0294
0295     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HE");
0296     emEnergyTiming_HE = ibooker.book2D(histo, histo, 200, 0., 800., 110, -120., 100.);
0297
0298     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_Low_HE");
0299     emEnergyTiming_profile_Low_HE = ibooker.bookProfile(histo, histo, 160, 0., 160., 110, -120., 100.);
0300
0301     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HE");
0302     emEnergyTiming_profile_HE = ibooker.bookProfile(histo, histo, 200, 0., 800., 110, -120., 100.);
0303
0304     // HAD
0305     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HE");
0306     hadEnergyTiming_Low_HE = ibooker.book2D(histo, histo, 160, 0., 160., 70, -48., 92.);
0307
0308     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HE");
0309     hadEnergyTiming_HE = ibooker.book2D(histo, histo, 200, 0., 800., 70, -48., 92.);
0310
0311     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HE");
0312     hadEnergyTiming_profile_Low_HE = ibooker.bookProfile(histo, histo, 160, 0., 160., 70, -48., 92.);
0313
0314     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HE");
0315     hadEnergyTiming_profile_HE = ibooker.bookProfile(histo, histo, 200, 0., 800., 70, -48., 92.);
0316     //-------------------------------------------------------------------------------------------
0317
0318     sprintf(histo, "CaloTowersTask_Iphi_HCAL_component_of_tower_HEP");
0319     meIphiHcalTower_HEP = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0320     sprintf(histo, "CaloTowersTask_Iphi_HCAL_component_of_tower_HEM");
0321     meIphiHcalTower_HEM = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0322
0323     //-------------------------------------------------------------------------------------------
0324
0325     // Everything else is not drawn
0326     if (useAllHistos_) {
0327       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HE");
0328       meEnergyHcalvsEcal_HE = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0329
0330       sprintf(histo, "CaloTowersTask_energy_OUTER_HE");
0331       meEnergyHO_HE = ibooker.book1D(histo, histo, 440, -200, 2000);
0332
0333       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HE");
0334       meEnergyEcalTower_HE = ibooker.book1D(histo, histo, 1100, -200, 2000);
0335
0336       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HE");
0337       meEnergyHcalTower_HE = ibooker.book1D(histo, histo, 1100, -200, 2000);
0338
0339       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HE");
0340       meTotEnergy_HE = ibooker.book1D(histo, histo, 400, 0., 2000.);
0341
0342       sprintf(histo, "CaloTowersTask_map_energy_HE");
0343       mapEnergy_HE = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0344       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HE");
0345       mapEnergyHcal_HE =
0346           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0347       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HE");
0348       mapEnergyEcal_HE =
0349           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0350
0351       sprintf(histo, "CaloTowersTask_phi_MET_HE");
0352       phiMET_HE = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0353     }
0354   }
0355
0356   if (isub == 3 || isub == 0) {
0357     // All cell histos are used
0358     sprintf(histo, "number_of_bad_cells_Hcal_HF");
0359     numBadCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0360     sprintf(histo, "number_of_recovered_cells_Hcal_HF");
0361     numRcvCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0362     sprintf(histo, "number_of_problematic_cells_Hcal_HF");
0363     numPrbCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0364
0365     // These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL,
0366     // MET and SET
0367     //-------------------------------------------------------------------------------------------
0368     sprintf(histo, "CaloTowersTask_energy_HCAL_HF");
0369     meEnergyHcal_HF = ibooker.book1D(histo, histo, 5001, -20., 1000000.);
0370
0371     sprintf(histo, "CaloTowersTask_energy_ECAL_HF");
0372     meEnergyEcal_HF = ibooker.book1D(histo, histo, 3501, -20., 70000.);
0373
0374     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HF");
0375     meNumFiredTowers_HF = ibooker.book1D(histo, histo, 1000, 0., 2000.);
0376
0377     sprintf(histo, "CaloTowersTask_MET_HF");
0378     MET_HF = ibooker.book1D(histo, histo, 500, 0., 500.);
0379
0380     sprintf(histo, "CaloTowersTask_SET_HF");
0381     SET_HF = ibooker.book1D(histo, histo, 500, 0., 5000.);
0382     //-------------------------------------------------------------------------------------------
0383
0384     // Timing histos and profiles -- all six are necessary
0385     //-------------------------------------------------------------------------------------------
0386     sprintf(histo, "CaloTowersTask_EM_Timing_HF");
0387     emTiming_HF = ibooker.book1D(histo, histo, 110, -120., 100.);
0388
0389     sprintf(histo, "CaloTowersTask_HAD_Timing_HF");
0390     hadTiming_HF = ibooker.book1D(histo, histo, 70, -48., 92.);
0391
0392     // EM
0393     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HF");
0394     emEnergyTiming_HF = ibooker.book2D(histo, histo, 150, 0., 300., 110, -120., 100.);
0395
0396     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HF");
0397     emEnergyTiming_profile_HF = ibooker.bookProfile(histo, histo, 150, 0., 300., 110, -120., 100.);
0398
0399     // HAD (requires two different sets of histograms to lower RAM usage)
0400     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HF");
0401     hadEnergyTiming_Low_HF = ibooker.book2D(histo, histo, 40, 0., 40., 70, -48., 92.);
0402
0403     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HF");
0404     hadEnergyTiming_HF = ibooker.book2D(histo, histo, 200, 0., 600., 70, -48., 92.);
0405
0406     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HF");
0407     hadEnergyTiming_profile_Low_HF = ibooker.bookProfile(histo, histo, 40, 0., 40., 70, -48., 92.);
0408
0409     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HF");
0410     hadEnergyTiming_profile_HF = ibooker.bookProfile(histo, histo, 200, 0., 600., 70, -48., 92.);
0411     //-------------------------------------------------------------------------------------------
0412
0413     sprintf(histo, "CaloTowersTask_Iphi_tower_HFP");
0414     meIphiCaloTower_HFP = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0415     sprintf(histo, "CaloTowersTask_Iphi_tower_HFM");
0416     meIphiCaloTower_HFM = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0417
0418     //-------------------------------------------------------------------------------------------
0419
0420     // Everything else is not drawn
0421     if (useAllHistos_) {
0422       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HF");
0423       meEnergyHcalvsEcal_HF = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0424
0425       sprintf(histo, "CaloTowersTask_energy_OUTER_HF");
0426       meEnergyHO_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0427
0428       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HF");
0429       meEnergyEcalTower_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0430
0431       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HF");
0432       meEnergyHcalTower_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0433
0434       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HF");
0435       meTotEnergy_HF = ibooker.book1D(histo, histo, 400, 0., 2000.);
0436
0437       sprintf(histo, "CaloTowersTask_map_energy_HF");
0438       mapEnergy_HF = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0439       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HF");
0440       mapEnergyHcal_HF =
0441           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0442       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HF");
0443       mapEnergyEcal_HF =
0444           ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0445
0446       sprintf(histo, "CaloTowersTask_phi_MET_HF");
0447       phiMET_HF = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0448     }
0449   }
0450 }
0451
0452 CaloTowersAnalyzer::~CaloTowersAnalyzer() {}
0453
0454 void CaloTowersAnalyzer::analyze(edm::Event const &event, edm::EventSetup const &) {
0455   nevent++;
0456
0457   edm::Handle<CaloTowerCollection> towers;
0458   event.getByToken(tok_towers_, towers);
0459   CaloTowerCollection::const_iterator cal;
0460
0461   double met;
0462   double phimet;
0463
0464   // ieta scan
0465   // double partR  = 0.3;
0466   // double Rmin   = 9999.;
0467   // double Econe  = 0.;
0468   // double Hcone  = 0.;
0469   // double Ee1    = 0.;
0470   // double Eh1    = 0.;
0471   double ieta_MC = 9999;
0472   double iphi_MC = 9999;
0473
0474   // HB
0475   double sumEnergyHcal_HB = 0.;
0476   double sumEnergyEcal_HB = 0.;
0477   double sumEnergyHO_HB = 0.;
0478   Int_t numFiredTowers_HB = 0;
0479   double metx_HB = 0.;
0480   double mety_HB = 0.;
0481   double metz_HB = 0.;
0482   double sEt_HB = 0.;
0483   // HE
0484   double sumEnergyHcal_HE = 0.;
0485   double sumEnergyEcal_HE = 0.;
0486   double sumEnergyHO_HE = 0.;
0487   Int_t numFiredTowers_HE = 0;
0488   double metx_HE = 0.;
0489   double mety_HE = 0.;
0490   double metz_HE = 0.;
0491   double sEt_HE = 0.;
0492   // HF
0493   double sumEnergyHcal_HF = 0.;
0494   double sumEnergyEcal_HF = 0.;
0495   double sumEnergyHO_HF = 0.;
0496   Int_t numFiredTowers_HF = 0;
0497   double metx_HF = 0.;
0498   double mety_HF = 0.;
0499   double metz_HF = 0.;
0500   double sEt_HF = 0.;
0501
0502   for (cal = towers->begin(); cal != towers->end(); ++cal) {
0503     double eE = cal->emEnergy();
0504     double eH = cal->hadEnergy();
0505     double eHO = cal->outerEnergy();
0506     double etaT = cal->eta();
0507     //    double phiT   = cal->phi();
0508     double en = cal->energy();
0509     double etT = cal->et();
0510     double had_tm = cal->hcalTime();
0511     double em_tm = cal->ecalTime();
0512
0513     int numBadEcalCells = cal->numBadEcalCells();
0514     int numRcvEcalCells = cal->numRecoveredEcalCells();
0515     int numPrbEcalCells = cal->numProblematicEcalCells();
0516
0517     int numBadHcalCells = cal->numBadHcalCells();
0518     int numRcvHcalCells = cal->numRecoveredHcalCells();
0519     int numPrbHcalCells = cal->numProblematicHcalCells();
0520
0521     math::RhoEtaPhiVector mom(cal->et(), cal->eta(), cal->phi());
0522     //  Vector mom  = cal->momentum();
0523
0524     // cell properties
0525     CaloTowerDetId idT = cal->id();
0526     int ieta = idT.ieta();
0527     int iphi = idT.iphi();
0528
0529     // ecal:  0 EcalBarrel  1 EcalEndcap
0530     // hcal:  0 hcalBarrel  1 HcalEndcap  2 HcalForward
0531     std::vector<int> inEcals(2), inHcals(3);
0532     unsigned int constitSize = cal->constituentsSize();
0533     for (unsigned int ic = 0; ic < constitSize; ic++) {
0534       DetId detId = cal->constituent(ic);
0535       if (detId.det() == DetId::Ecal) {
0536         if (detId.subdetId() == EcalBarrel)
0537           inEcals[0] = 1;
0538         else if (detId.subdetId() == EcalEndcap)
0539           inEcals[1] = 1;
0540       }
0541       if (detId.det() == DetId::Hcal) {
0542         if (HcalDetId(detId).subdet() == HcalBarrel)
0543           inHcals[0] = 1;
0544         else if (HcalDetId(detId).subdet() == HcalEndcap)
0545           inHcals[1] = 1;
0546         else if (HcalDetId(detId).subdet() == HcalForward)
0547           inHcals[2] = 1;
0548       }
0549     }
0550     // All cell histos are used
0551     if (inEcals[0]) {
0552       numBadCellsEcal_EB->Fill(numBadEcalCells);
0553       numRcvCellsEcal_EB->Fill(numRcvEcalCells);
0554       numPrbCellsEcal_EB->Fill(numPrbEcalCells);
0555     }
0556     if (inEcals[1]) {
0557       numBadCellsEcal_EE->Fill(numBadEcalCells);
0558       numRcvCellsEcal_EE->Fill(numRcvEcalCells);
0559       numPrbCellsEcal_EE->Fill(numPrbEcalCells);
0560     }
0561
0562     // Ntowers is used in EndJob, occupancy_map is used for occupancy vs ieta
0563     Ntowers_vs_ieta->Fill(double(ieta), 1.);
0564     occupancy_map->Fill(double(ieta), double(iphi));
0565
0566     if ((isub == 0 || isub == 1) && (fabs(etaT) < etaMax[0] && fabs(etaT) >= etaMin[0])) {
0567       // All cell histos are used
0568       numBadCellsHcal_HB->Fill(numBadHcalCells);
0569       numRcvCellsHcal_HB->Fill(numRcvHcalCells);
0570       numPrbCellsHcal_HB->Fill(numPrbHcalCells);
0571
0572       // Map energy histos are not used
0573       if (useAllHistos_) {
0574         mapEnergy_HB->Fill(double(ieta), double(iphi), en);
0575         mapEnergyHcal_HB->Fill(double(ieta), double(iphi), eH);
0576         mapEnergyEcal_HB->Fill(double(ieta), double(iphi), eE);
0577       }
0578       //      std::cout << " e_ecal = " << eE << std::endl;
0579
0580       if (eH > 0.) {  // iphi of HCAL component of calotower
0581         if (ieta > 0)
0582           meIphiHcalTower_HBP->Fill(double(iphi));
0583         else
0584           meIphiHcalTower_HBM->Fill(double(iphi));
0585       }
0586
0587       //  simple sums
0588       sumEnergyHcal_HB += eH;
0589       sumEnergyEcal_HB += eE;
0590       sumEnergyHO_HB += eHO;
0591
0592       numFiredTowers_HB++;
0593
0594       // Not used
0595       if (useAllHistos_) {
0596         meEnergyEcalTower_HB->Fill(eE);
0597         meEnergyHcalTower_HB->Fill(eH);
0598       }
0599
0600       // MET, SET & phimet
0601       //  double  etT = cal->et();
0602       metx_HB += mom.x();
0603       mety_HB += mom.y();  // etT * sin(phiT);
0604       sEt_HB += etT;
0605
0606       // Timing (all histos are used)
0607       emTiming_HB->Fill(em_tm);
0608       hadTiming_HB->Fill(had_tm);
0609
0610       emEnergyTiming_Low_HB->Fill(eE, em_tm);
0611       emEnergyTiming_HB->Fill(eE, em_tm);
0612       emEnergyTiming_High_HB->Fill(eE, em_tm);
0613       emEnergyTiming_profile_Low_HB->Fill(eE, em_tm);
0614       emEnergyTiming_profile_HB->Fill(eE, em_tm);
0615       emEnergyTiming_profile_High_HB->Fill(eE, em_tm);
0616
0617       hadEnergyTiming_Low_HB->Fill(eH, had_tm);
0618       hadEnergyTiming_HB->Fill(eH, had_tm);
0619       hadEnergyTiming_High_HB->Fill(eH, had_tm);
0620       hadEnergyTiming_profile_Low_HB->Fill(eH, had_tm);
0621       hadEnergyTiming_profile_HB->Fill(eH, had_tm);
0622       hadEnergyTiming_profile_High_HB->Fill(eH, had_tm);
0623     }
0624
0625     if ((isub == 0 || isub == 2) && (fabs(etaT) < etaMax[1] && fabs(etaT) >= etaMin[1])) {
0626       // All cell histos are used
0627       numBadCellsHcal_HE->Fill(numBadHcalCells);
0628       numRcvCellsHcal_HE->Fill(numRcvHcalCells);
0629       numPrbCellsHcal_HE->Fill(numPrbHcalCells);
0630
0631       // Map energy histos are not used
0632       if (useAllHistos_) {
0633         mapEnergy_HE->Fill(double(ieta), double(iphi), en);
0634         mapEnergyHcal_HE->Fill(double(ieta), double(iphi), eH);
0635         mapEnergyEcal_HE->Fill(double(ieta), double(iphi), eE);
0636       }
0637       //      std::cout << " e_ecal = " << eE << std::endl;
0638
0639       if (eH > 0.) {  // iphi of HCAL component of calotower
0640         if (ieta > 0)
0641           meIphiHcalTower_HEP->Fill(double(iphi));
0642         else
0643           meIphiHcalTower_HEM->Fill(double(iphi));
0644       }
0645
0646       //  simple sums
0647       sumEnergyHcal_HE += eH;
0648       sumEnergyEcal_HE += eE;
0649       sumEnergyHO_HE += eHO;
0650
0651       numFiredTowers_HE++;
0652
0653       // Not used
0654       if (useAllHistos_) {
0655         meEnergyEcalTower_HE->Fill(eE);
0656         meEnergyHcalTower_HE->Fill(eH);
0657       }
0658       // MET, SET & phimet
0659       //  double  etT = cal->et();
0660       metx_HE += mom.x();
0661       mety_HE += mom.y();  // etT * sin(phiT);
0662       sEt_HE += etT;
0663
0664       // Timing (all histos are used)
0665       emTiming_HE->Fill(em_tm);
0666       hadTiming_HE->Fill(had_tm);
0667
0668       emEnergyTiming_Low_HE->Fill(eE, em_tm);
0669       emEnergyTiming_HE->Fill(eE, em_tm);
0670       emEnergyTiming_profile_Low_HE->Fill(eE, em_tm);
0671       emEnergyTiming_profile_HE->Fill(eE, em_tm);
0672
0673       hadEnergyTiming_Low_HE->Fill(eH, had_tm);
0674       hadEnergyTiming_HE->Fill(eH, had_tm);
0675       hadEnergyTiming_profile_Low_HE->Fill(eH, had_tm);
0676       hadEnergyTiming_profile_HE->Fill(eH, had_tm);
0677     }
0678
0679     if ((isub == 0 || isub == 3) && (fabs(etaT) < etaMax[2] && fabs(etaT) >= etaMin[2])) {
0680       // All cell histos are used
0681       numBadCellsHcal_HF->Fill(numBadHcalCells);
0682       numRcvCellsHcal_HF->Fill(numRcvHcalCells);
0683       numPrbCellsHcal_HF->Fill(numPrbHcalCells);
0684
0685       // Map energy histos are not used
0686       if (useAllHistos_) {
0687         mapEnergy_HF->Fill(double(ieta), double(iphi), en);
0688         mapEnergyHcal_HF->Fill(double(ieta), double(iphi), eH);
0689         mapEnergyEcal_HF->Fill(double(ieta), double(iphi), eE);
0690       }
0691       //      std::cout << " e_ecal = " << eE << std::endl;
0692
0693       if (ieta > 0)
0694         meIphiCaloTower_HFP->Fill(double(iphi));
0695       else
0696         meIphiCaloTower_HFM->Fill(double(iphi));
0697
0698       //  simple sums
0699       sumEnergyHcal_HF += eH;
0700       sumEnergyEcal_HF += eE;
0701       sumEnergyHO_HF += eHO;
0702
0703       numFiredTowers_HF++;
0704
0705       // Not used
0706       if (useAllHistos_) {
0707         meEnergyEcalTower_HF->Fill(eE);
0708         meEnergyHcalTower_HF->Fill(eH);
0709       }
0710       // MET, SET & phimet
0711       //  double  etT = cal->et();
0712       metx_HF += mom.x();
0713       mety_HF += mom.y();  // etT * sin(phiT);
0714       sEt_HF += etT;
0715
0716       // Timing (all histos are used)
0717       emTiming_HF->Fill(em_tm);
0718       hadTiming_HF->Fill(had_tm);
0719       emEnergyTiming_HF->Fill(eE, em_tm);
0720       emEnergyTiming_profile_HF->Fill(eE, em_tm);
0721
0722       hadEnergyTiming_Low_HF->Fill(eH, had_tm);
0723       hadEnergyTiming_HF->Fill(eH, had_tm);
0724       hadEnergyTiming_profile_Low_HF->Fill(eH, had_tm);
0725       hadEnergyTiming_profile_HF->Fill(eH, had_tm);
0726     }
0727
0728   }  // end of Towers cycle
0729
0730   // These are the six single pion histos; only the second set is used
0731
0732   mapEnergy_N->Fill(double(ieta_MC), double(iphi_MC), 1.);
0733
0734   if (isub == 0 || isub == 1) {
0735     met = sqrt(metx_HB * metx_HB + mety_HB * mety_HB);
0736     Vector metv(metx_HB, mety_HB, metz_HB);
0737     phimet = metv.phi();
0738
0739     // Five oldest drawn histos first; the rest are not used
0740     meEnergyHcal_HB->Fill(sumEnergyHcal_HB);
0741     meEnergyEcal_HB->Fill(sumEnergyEcal_HB);
0742     meNumFiredTowers_HB->Fill(numFiredTowers_HB);
0743     MET_HB->Fill(met);
0744     SET_HB->Fill(sEt_HB);
0745
0746     if (useAllHistos_) {
0747       meEnergyHcalvsEcal_HB->Fill(sumEnergyEcal_HB, sumEnergyHcal_HB);
0748       meEnergyHO_HB->Fill(sumEnergyHO_HB);
0749       meTotEnergy_HB->Fill(sumEnergyHcal_HB + sumEnergyEcal_HB + sumEnergyHO_HB);
0750       phiMET_HB->Fill(phimet);
0751     }
0752   }
0753
0754   if (isub == 0 || isub == 2) {
0755     met = sqrt(metx_HE * metx_HE + mety_HE * mety_HE);
0756     Vector metv(metx_HE, mety_HE, metz_HE);
0757     phimet = metv.phi();
0758
0759     // Five oldest drawn histos first; the rest are not used
0760     meEnergyHcal_HE->Fill(sumEnergyHcal_HE);
0761     meEnergyEcal_HE->Fill(sumEnergyEcal_HE);
0762     meNumFiredTowers_HE->Fill(numFiredTowers_HE);
0763     MET_HE->Fill(met);
0764     SET_HE->Fill(sEt_HE);
0765
0766     if (useAllHistos_) {
0767       meEnergyHcalvsEcal_HE->Fill(sumEnergyEcal_HE, sumEnergyHcal_HE);
0768       meEnergyHO_HE->Fill(sumEnergyHO_HE);
0769       meTotEnergy_HE->Fill(sumEnergyHcal_HE + sumEnergyEcal_HE + sumEnergyHO_HE);
0770       phiMET_HE->Fill(phimet);
0771     }
0772   }
0773
0774   if (isub == 0 || isub == 3) {
0775     met = sqrt(metx_HF * metx_HF + mety_HF * mety_HF);
0776     Vector metv(metx_HF, mety_HF, metz_HF);
0777     phimet = metv.phi();
0778
0779     // Five oldest drawn histos first; the rest are not used
0780     meEnergyHcal_HF->Fill(sumEnergyHcal_HF);
0781     meEnergyEcal_HF->Fill(sumEnergyEcal_HF);
0782     meNumFiredTowers_HF->Fill(numFiredTowers_HF);
0783     MET_HF->Fill(met);
0784     SET_HF->Fill(sEt_HF);
0785
0786     if (useAllHistos_) {
0787       meEnergyHcalvsEcal_HF->Fill(sumEnergyEcal_HF, sumEnergyHcal_HF);
0788       meEnergyHO_HF->Fill(sumEnergyHO_HF);
0789       meTotEnergy_HF->Fill(sumEnergyHcal_HF + sumEnergyEcal_HF + sumEnergyHO_HF);
0790       phiMET_HF->Fill(phimet);
0791     }
0792   }
0793 }
0794
0795 double CaloTowersAnalyzer::dR(double eta1, double phi1, double eta2, double phi2) {
0796   double PI = 3.1415926535898;
0797   double deltaphi = phi1 - phi2;
0798   if (phi2 > phi1) {
0799     deltaphi = phi2 - phi1;
0800   }
0801   if (deltaphi > PI) {
0802     deltaphi = 2. * PI - deltaphi;
0803   }
0804   double deltaeta = eta2 - eta1;
0805   double tmp = sqrt(deltaeta * deltaeta + deltaphi * deltaphi);
0806   return tmp;
0807 }
0808
0809 DEFINE_FWK_MODULE(CaloTowersAnalyzer);