File indexing completed on 2023-01-13 23:39:43
0001 #include "DQMOffline/Hcal/interface/HcalRecHitsAnalyzer.h"
0002 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0003
0004 #include "FWCore/Utilities/interface/ESInputTag.h"
0005 #include "FWCore/Utilities/interface/Transition.h"
0006
0007 HcalRecHitsAnalyzer::HcalRecHitsAnalyzer(edm::ParameterSet const &conf)
0008 : topFolderName_(conf.getParameter<std::string>("TopFolderName")),
0009 hcalDDDRecConstantsToken_{esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord, edm::Transition::BeginRun>()},
0010 caloGeometryRunToken_{esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>()},
0011 caloGeometryEventToken_{esConsumes<CaloGeometry, CaloGeometryRecord>()},
0012 hcalTopologyToken_{esConsumes<HcalTopology, HcalRecNumberingRecord>()},
0013 hcalChannelQualityToken_{esConsumes<HcalChannelQuality, HcalChannelQualityRcd>(edm::ESInputTag("", "withTopo"))},
0014 hcalSeverityLevelComputerToken_{esConsumes<HcalSeverityLevelComputer, HcalSeverityLevelComputerRcd>()} {
0015
0016 outputFile_ = conf.getUntrackedParameter<std::string>("outputFile", "myfile.root");
0017
0018 if (!outputFile_.empty()) {
0019 edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will be saved to '" << outputFile_.c_str() << "'";
0020 } else {
0021 edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will NOT be saved";
0022 }
0023
0024 nevtot = 0;
0025
0026 hcalselector_ = conf.getUntrackedParameter<std::string>("hcalselector", "all");
0027 ecalselector_ = conf.getUntrackedParameter<std::string>("ecalselector", "yes");
0028 eventype_ = conf.getUntrackedParameter<std::string>("eventype", "single");
0029 sign_ = conf.getUntrackedParameter<std::string>("sign", "*");
0030
0031
0032
0033 hep17_ = conf.getUntrackedParameter<bool>("hep17");
0034
0035
0036 tok_hbhe_ = consumes<HBHERecHitCollection>(conf.getUntrackedParameter<edm::InputTag>("HBHERecHitCollectionLabel"));
0037 tok_hf_ = consumes<HFRecHitCollection>(conf.getUntrackedParameter<edm::InputTag>("HFRecHitCollectionLabel"));
0038 tok_ho_ = consumes<HORecHitCollection>(conf.getUntrackedParameter<edm::InputTag>("HORecHitCollectionLabel"));
0039 edm::InputTag EBRecHitCollectionLabel = conf.getParameter<edm::InputTag>("EBRecHitCollectionLabel");
0040 tok_EB_ = consumes<EBRecHitCollection>(EBRecHitCollectionLabel);
0041 edm::InputTag EERecHitCollectionLabel = conf.getParameter<edm::InputTag>("EERecHitCollectionLabel");
0042 tok_EE_ = consumes<EERecHitCollection>(EERecHitCollectionLabel);
0043
0044 subdet_ = 5;
0045 if (hcalselector_ == "noise")
0046 subdet_ = 0;
0047 if (hcalselector_ == "HB")
0048 subdet_ = 1;
0049 if (hcalselector_ == "HE")
0050 subdet_ = 2;
0051 if (hcalselector_ == "HO")
0052 subdet_ = 3;
0053 if (hcalselector_ == "HF")
0054 subdet_ = 4;
0055 if (hcalselector_ == "all")
0056 subdet_ = 5;
0057 if (hcalselector_ == "ZS")
0058 subdet_ = 6;
0059
0060 etype_ = 1;
0061 if (eventype_ == "multi")
0062 etype_ = 2;
0063
0064 iz = 1;
0065 if (sign_ == "-")
0066 iz = -1;
0067 if (sign_ == "*")
0068 iz = 0;
0069
0070 imc = 0;
0071 }
0072
0073 void HcalRecHitsAnalyzer::dqmBeginRun(const edm::Run &run, const edm::EventSetup &iSetup) {
0074 HcalDDDRecConstants const &hcons = iSetup.getData(hcalDDDRecConstantsToken_);
0075 maxDepthHB_ = hcons.getMaxDepth(0);
0076 maxDepthHE_ = hcons.getMaxDepth(1);
0077 maxDepthHF_ = std::max(hcons.getMaxDepth(2), 1);
0078 maxDepthHO_ = hcons.getMaxDepth(3);
0079
0080 CaloGeometry const &geo = iSetup.getData(caloGeometryRunToken_);
0081
0082 const HcalGeometry *gHB = static_cast<const HcalGeometry *>(geo.getSubdetectorGeometry(DetId::Hcal, HcalBarrel));
0083 const HcalGeometry *gHE = static_cast<const HcalGeometry *>(geo.getSubdetectorGeometry(DetId::Hcal, HcalEndcap));
0084 const HcalGeometry *gHO = static_cast<const HcalGeometry *>(geo.getSubdetectorGeometry(DetId::Hcal, HcalOuter));
0085 const HcalGeometry *gHF = static_cast<const HcalGeometry *>(geo.getSubdetectorGeometry(DetId::Hcal, HcalForward));
0086
0087 nChannels_[1] = gHB->getHxSize(1);
0088 nChannels_[2] = std::max(int(gHE->getHxSize(2)), 1);
0089 nChannels_[3] = gHO->getHxSize(3);
0090 nChannels_[4] = gHF->getHxSize(4);
0091
0092 nChannels_[0] = nChannels_[1] + nChannels_[2] + nChannels_[3] + nChannels_[4];
0093
0094
0095
0096
0097
0098
0099 maxDepthHF_ = (maxDepthHF_ > 2 ? 2 : maxDepthHF_);
0100
0101
0102 maxDepthAll_ = (maxDepthHB_ + maxDepthHO_ > maxDepthHE_ ? maxDepthHB_ + maxDepthHO_ : maxDepthHE_);
0103 maxDepthAll_ = (maxDepthAll_ > maxDepthHF_ ? maxDepthAll_ : maxDepthHF_);
0104
0105
0106
0107
0108 int NphiMax = hcons.getNPhi(0);
0109
0110 NphiMax = (hcons.getNPhi(1) > NphiMax ? hcons.getNPhi(1) : NphiMax);
0111 NphiMax = (hcons.getNPhi(2) > NphiMax ? hcons.getNPhi(2) : NphiMax);
0112 NphiMax = (hcons.getNPhi(3) > NphiMax ? hcons.getNPhi(3) : NphiMax);
0113
0114
0115 iphi_min_ = 0.5;
0116 iphi_max_ = NphiMax + 0.5;
0117 iphi_bins_ = (int)(iphi_max_ - iphi_min_);
0118
0119
0120
0121 int iEtaMax = (hcons.getEtaRange(0).second > hcons.getEtaRange(1).second ? hcons.getEtaRange(0).second
0122 : hcons.getEtaRange(1).second);
0123 iEtaMax = (iEtaMax > hcons.getEtaRange(2).second ? iEtaMax : hcons.getEtaRange(2).second);
0124 iEtaMax = (iEtaMax > hcons.getEtaRange(3).second ? iEtaMax : hcons.getEtaRange(3).second);
0125
0126
0127
0128 ieta_min_ = -iEtaMax - 1.5;
0129 ieta_max_ = iEtaMax + 1.5;
0130 ieta_bins_ = (int)(ieta_max_ - ieta_min_);
0131 }
0132
0133 void HcalRecHitsAnalyzer::bookHistograms(DQMStore::IBooker &ibooker,
0134 edm::Run const & ,
0135 edm::EventSetup const &)
0136
0137 {
0138 Char_t histo[200];
0139
0140 ibooker.setCurrentFolder(topFolderName_);
0141
0142
0143
0144
0145
0146
0147
0148
0149 for (int depth = 0; depth <= maxDepthHB_; depth++) {
0150 if (depth == 0) {
0151 sprintf(histo, "N_HB");
0152 } else {
0153 sprintf(histo, "N_HB_depth%d", depth);
0154 }
0155 int NBins = (int)(nChannels_[1] * 1.1);
0156 Nhb.push_back(ibooker.book1D(histo, histo, NBins, 0., (float)NBins));
0157 }
0158 for (int depth = 0; depth <= maxDepthHE_; depth++) {
0159 if (depth == 0) {
0160 sprintf(histo, "N_HE");
0161 } else {
0162 sprintf(histo, "N_HE_depth%d", depth);
0163 }
0164 int NBins = (int)(nChannels_[2] * 1.1);
0165 Nhe.push_back(ibooker.book1D(histo, histo, NBins, 0., (float)NBins));
0166 }
0167 for (int depth = 0; depth <= maxDepthHO_; depth++) {
0168 if (depth == 0) {
0169 sprintf(histo, "N_HO");
0170 } else {
0171 sprintf(histo, "N_HO_depth%d", depth);
0172 }
0173 int NBins = (int)(nChannels_[3] * 1.1);
0174 Nho.push_back(ibooker.book1D(histo, histo, NBins, 0., (float)NBins));
0175 }
0176 for (int depth = 0; depth <= maxDepthHF_; depth++) {
0177 if (depth == 0) {
0178 sprintf(histo, "N_HF");
0179 } else {
0180 sprintf(histo, "N_HF_depth%d", depth);
0181 }
0182 int NBins = (int)(nChannels_[4] * 1.1);
0183 Nhf.push_back(ibooker.book1D(histo, histo, NBins, 0., (float)NBins));
0184 }
0185
0186
0187 if (subdet_ == 6) {
0188 }
0189
0190
0191 else {
0192 for (int depth = 1; depth <= maxDepthAll_; depth++) {
0193 sprintf(histo, "emap_depth%d", depth);
0194 emap.push_back(ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_));
0195 }
0196 sprintf(histo, "emap_HO");
0197 emap_HO = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0198
0199
0200
0201 for (int depth = 1; depth <= maxDepthHB_; depth++) {
0202 sprintf(histo, "emean_vs_ieta_HB%d", depth);
0203 emean_vs_ieta_HB.push_back(ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0204
0205 sprintf(histo, "emean_vs_ieta_M0_HB%d", depth);
0206 emean_vs_ieta_HBM0.push_back(
0207 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0208
0209 sprintf(histo, "emean_vs_ieta_M3_HB%d", depth);
0210 emean_vs_ieta_HBM3.push_back(
0211 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0212 }
0213 for (int depth = 1; depth <= maxDepthHE_; depth++) {
0214 sprintf(histo, "emean_vs_ieta_HE%d", depth);
0215 emean_vs_ieta_HE.push_back(ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0216
0217 sprintf(histo, "emean_vs_ieta_M0_HE%d", depth);
0218 emean_vs_ieta_HEM0.push_back(
0219 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0220
0221 sprintf(histo, "emean_vs_ieta_M3_HE%d", depth);
0222 emean_vs_ieta_HEM3.push_back(
0223 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0224 }
0225
0226 if (hep17_) {
0227 for (int depth = 1; depth <= maxDepthHE_; depth++) {
0228 sprintf(histo, "emean_vs_ieta_HEP17_depth%d", depth);
0229 emean_vs_ieta_HEP17.push_back(
0230 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0231
0232 sprintf(histo, "emean_vs_ieta_M0_HEP17_depth%d", depth);
0233 emean_vs_ieta_HEP17M0.push_back(
0234 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0235
0236 sprintf(histo, "emean_vs_ieta_M3_HEP17_depth%d", depth);
0237 emean_vs_ieta_HEP17M3.push_back(
0238 ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0239 }
0240 }
0241
0242 for (int depth = 1; depth <= maxDepthHF_; depth++) {
0243 sprintf(histo, "emean_vs_ieta_HF%d", depth);
0244 emean_vs_ieta_HF.push_back(ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " "));
0245 }
0246 sprintf(histo, "emean_vs_ieta_HO");
0247 emean_vs_ieta_HO = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -10., 2000., " ");
0248
0249
0250
0251
0252 for (int depth = 1; depth <= maxDepthHB_; depth++) {
0253 sprintf(histo, "occupancy_map_HB%d", depth);
0254 occupancy_map_HB.push_back(
0255 ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_));
0256 }
0257
0258 for (int depth = 1; depth <= maxDepthHE_; depth++) {
0259 sprintf(histo, "occupancy_map_HE%d", depth);
0260 occupancy_map_HE.push_back(
0261 ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_));
0262 }
0263
0264 sprintf(histo, "occupancy_map_HO");
0265 occupancy_map_HO = ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_);
0266
0267 for (int depth = 1; depth <= maxDepthHF_; depth++) {
0268 sprintf(histo, "occupancy_map_HF%d", depth);
0269 occupancy_map_HF.push_back(
0270 ibooker.book2D(histo, histo, ieta_bins_, ieta_min_, ieta_max_, iphi_bins_, iphi_min_, iphi_max_));
0271 }
0272
0273
0274 for (int depth = 1; depth <= maxDepthHB_; depth++) {
0275 sprintf(histo, "occupancy_vs_ieta_HB%d", depth);
0276 occupancy_vs_ieta_HB.push_back(ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_));
0277 sprintf(histo, "nrechits_vs_iphi_HBP_d%d", depth);
0278 nrechits_vs_iphi_HBP.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0279 sprintf(histo, "nrechits_vs_iphi_HBM_d%d", depth);
0280 nrechits_vs_iphi_HBM.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0281 }
0282
0283 for (int depth = 1; depth <= maxDepthHE_; depth++) {
0284 sprintf(histo, "occupancy_vs_ieta_HE%d", depth);
0285 occupancy_vs_ieta_HE.push_back(ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_));
0286 sprintf(histo, "nrechits_vs_iphi_HEP_d%d", depth);
0287 nrechits_vs_iphi_HEP.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0288 sprintf(histo, "nrechits_vs_iphi_HEM_d%d", depth);
0289 nrechits_vs_iphi_HEM.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0290 }
0291
0292 sprintf(histo, "occupancy_vs_ieta_HO");
0293 occupancy_vs_ieta_HO = ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_);
0294 sprintf(histo, "nrechits_vs_iphi_HOP");
0295 nrechits_vs_iphi_HOP = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0296 sprintf(histo, "nrechits_vs_iphi_HOM");
0297 nrechits_vs_iphi_HOM = ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_);
0298
0299 for (int depth = 1; depth <= maxDepthHF_; depth++) {
0300 sprintf(histo, "occupancy_vs_ieta_HF%d", depth);
0301 occupancy_vs_ieta_HF.push_back(ibooker.book1D(histo, histo, ieta_bins_, ieta_min_, ieta_max_));
0302 sprintf(histo, "nrechits_vs_iphi_HFP_d%d", depth);
0303 nrechits_vs_iphi_HFP.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0304 sprintf(histo, "nrechits_vs_iphi_HFM_d%d", depth);
0305 nrechits_vs_iphi_HFM.push_back(ibooker.book1D(histo, histo, iphi_bins_, iphi_min_, iphi_max_));
0306 }
0307
0308
0309 sprintf(histo, "HcalRecHitTask_RecHit_StatusWord_HB");
0310 RecHit_StatusWord_HB = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0311
0312 sprintf(histo, "HcalRecHitTask_RecHit_StatusWord_HE");
0313 RecHit_StatusWord_HE = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0314
0315 sprintf(histo, "HcalRecHitTask_RecHit_StatusWord_HF");
0316 RecHit_StatusWord_HF = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0317
0318 sprintf(histo, "HcalRecHitTask_RecHit_StatusWord_HO");
0319 RecHit_StatusWord_HO = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0320
0321
0322 sprintf(histo, "HcalRecHitTask_RecHit_Aux_StatusWord_HB");
0323 RecHit_Aux_StatusWord_HB = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0324
0325 sprintf(histo, "HcalRecHitTask_RecHit_Aux_StatusWord_HE");
0326 RecHit_Aux_StatusWord_HE = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0327
0328 sprintf(histo, "HcalRecHitTask_RecHit_Aux_StatusWord_HF");
0329 RecHit_Aux_StatusWord_HF = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0330
0331 sprintf(histo, "HcalRecHitTask_RecHit_Aux_StatusWord_HO");
0332 RecHit_Aux_StatusWord_HO = ibooker.book1DD(histo, histo, 32, -0.5, 31.5);
0333
0334 }
0335
0336
0337 sprintf(histo, "HcalRecHitTask_RecHit_StatusWordCorr_HB");
0338 RecHit_StatusWordCorr_HB = ibooker.book2D(histo, histo, 2, -0.5, 1.5, 2, -0.5, 1.5);
0339
0340 sprintf(histo, "HcalRecHitTask_RecHit_StatusWordCorr_HE");
0341 RecHit_StatusWordCorr_HE = ibooker.book2D(histo, histo, 2, -0.5, 1.5, 2, -0.5, 1.5);
0342
0343
0344
0345
0346 if (subdet_ != 0 && imc != 0) {
0347 sprintf(histo, "HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths");
0348 meEnConeEtaProfile = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -100., 2000., " ");
0349
0350 sprintf(histo, "HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_E");
0351 meEnConeEtaProfile_E = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -100., 2000., " ");
0352
0353 sprintf(histo, "HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_EH");
0354 meEnConeEtaProfile_EH = ibooker.bookProfile(histo, histo, ieta_bins_, ieta_min_, ieta_max_, -100., 2000., " ");
0355 }
0356
0357
0358 if (subdet_ == 1 || subdet_ == 5) {
0359
0360
0361
0362 sprintf(histo, "HcalRecHitTask_severityLevel_HB");
0363 sevLvl_HB = ibooker.book1DD(histo, histo, 25, -0.5, 24.5);
0364
0365 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HB");
0366 meRecHitsEnergyHB = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0367
0368 sprintf(histo, "HcalRecHitTask_cleaned_energy_of_rechits_HB");
0369 meRecHitsCleanedEnergyHB = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0370
0371 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M0_HB");
0372 meRecHitsEnergyHBM0 = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0373
0374 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3_HB");
0375 meRecHitsEnergyHBM3 = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0376
0377 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M2vM0_HB");
0378 meRecHitsEnergyM2vM0HB = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0379
0380 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3vM0_HB");
0381 meRecHitsEnergyM3vM0HB = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0382
0383 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3vM2_HB");
0384 meRecHitsEnergyM3vM2HB = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0385
0386 sprintf(histo, "HcalRecHitTask_M2Log10Chi2_of_rechits_HB");
0387 meRecHitsM2Chi2HB = ibooker.book1D(histo, histo, 120, -2., 10.);
0388
0389 sprintf(histo, "HcalRecHitTask_timing_HB");
0390 meTimeHB = ibooker.book1DD(histo, histo, 70, -48., 92.);
0391
0392
0393 sprintf(histo, "HcalRecHitTask_timing_vs_energy_Low_HB");
0394 meTE_Low_HB = ibooker.book2D(histo, histo, 50, -5., 45., 70, -48., 92.);
0395
0396 sprintf(histo, "HcalRecHitTask_timing_vs_energy_HB");
0397 meTE_HB = ibooker.book2D(histo, histo, 150, -5., 295., 70, -48., 92.);
0398
0399 sprintf(histo, "HcalRecHitTask_timing_vs_energy_High_HB");
0400 meTE_High_HB = ibooker.book2D(histo, histo, 150, -5., 2995., 70, -48., 92.);
0401
0402 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_Low_HB");
0403 meTEprofileHB_Low = ibooker.bookProfile(histo, histo, 50, -5., 45., -48., 92., " ");
0404
0405 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_HB");
0406 meTEprofileHB = ibooker.bookProfile(histo, histo, 150, -5., 295., -48., 92., " ");
0407
0408 sprintf(histo, "HcalRecHitTask_Log10Chi2_vs_energy_profile_HB");
0409 meLog10Chi2profileHB = ibooker.bookProfile(histo, histo, 150, -5., 295., -2., 9.9, " ");
0410
0411 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_High_HB");
0412 meTEprofileHB_High = ibooker.bookProfile(histo, histo, 150, -5., 2995., -48., 92., " ");
0413 }
0414
0415
0416 if (subdet_ == 2 || subdet_ == 5) {
0417
0418
0419 sprintf(histo, "HcalRecHitTask_severityLevel_HE");
0420 sevLvl_HE = ibooker.book1DD(histo, histo, 25, -0.5, 24.5);
0421
0422 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HE");
0423 meRecHitsEnergyHE = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0424
0425 sprintf(histo, "HcalRecHitTask_cleaned_energy_of_rechits_HE");
0426 meRecHitsCleanedEnergyHE = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0427
0428 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M0_HE");
0429 meRecHitsEnergyHEM0 = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0430
0431 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3_HE");
0432 meRecHitsEnergyHEM3 = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0433
0434 if (hep17_) {
0435 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HEP17");
0436 meRecHitsEnergyHEP17.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0437
0438 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M0_HEP17");
0439 meRecHitsEnergyHEP17M0.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0440
0441 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3_HEP17");
0442 meRecHitsEnergyHEP17M3.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0443 for (int depth = 1; depth <= maxDepthHE_; depth++) {
0444 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HEP17_depth%d", depth);
0445 meRecHitsEnergyHEP17.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0446
0447 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M0_HEP17_depth%d", depth);
0448 meRecHitsEnergyHEP17M0.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0449
0450 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3_HEP17_depth%d", depth);
0451 meRecHitsEnergyHEP17M3.push_back(ibooker.book1D(histo, histo, 2010, -10., 2000.));
0452 }
0453 }
0454
0455 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M2vM0_HE");
0456 meRecHitsEnergyM2vM0HE = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0457
0458 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3vM0_HE");
0459 meRecHitsEnergyM3vM0HE = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0460
0461 sprintf(histo, "HcalRecHitTask_energy_of_rechits_M3vM2_HE");
0462 meRecHitsEnergyM3vM2HE = ibooker.book2D(histo, histo, 42, -10., 200., 42, -10., 200.);
0463
0464 sprintf(histo, "HcalRecHitTask_M2Log10Chi2_of_rechits_HE");
0465 meRecHitsM2Chi2HE = ibooker.book1D(histo, histo, 120, -2., 10.);
0466
0467 sprintf(histo, "HcalRecHitTask_timing_HE");
0468 meTimeHE = ibooker.book1DD(histo, histo, 70, -48., 92.);
0469
0470 sprintf(histo, "HcalRecHitTask_timing_vs_energy_Low_HE");
0471 meTE_Low_HE = ibooker.book2D(histo, histo, 80, -5., 75., 70, -48., 92.);
0472
0473 sprintf(histo, "HcalRecHitTask_timing_vs_energy_HE");
0474 meTE_HE = ibooker.book2D(histo, histo, 200, -5., 395., 70, -48., 92.);
0475
0476 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_Low_HE");
0477 meTEprofileHE_Low = ibooker.bookProfile(histo, histo, 80, -5., 75., -48., 92., " ");
0478
0479 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_HE");
0480 meTEprofileHE = ibooker.bookProfile(histo, histo, 200, -5., 395., -48., 92., " ");
0481
0482 sprintf(histo, "HcalRecHitTask_Log10Chi2_vs_energy_profile_HE");
0483 meLog10Chi2profileHE = ibooker.bookProfile(histo, histo, 200, -5., 395., -2., 9.9, " ");
0484 }
0485
0486
0487 if (subdet_ == 3 || subdet_ == 5) {
0488
0489
0490
0491 sprintf(histo, "HcalRecHitTask_severityLevel_HO");
0492 sevLvl_HO = ibooker.book1DD(histo, histo, 25, -0.5, 24.5);
0493
0494 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HO");
0495 meRecHitsEnergyHO = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0496
0497 sprintf(histo, "HcalRecHitTask_cleaned_energy_of_rechits_HO");
0498 meRecHitsCleanedEnergyHO = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0499
0500 sprintf(histo, "HcalRecHitTask_timing_HO");
0501 meTimeHO = ibooker.book1DD(histo, histo, 80, -80., 80.);
0502
0503 sprintf(histo, "HcalRecHitTask_timing_vs_energy_HO");
0504 meTE_HO = ibooker.book2D(histo, histo, 60, -5., 55., 80, -80., 80.);
0505
0506 sprintf(histo, "HcalRecHitTask_timing_vs_energy_High_HO");
0507 meTE_High_HO = ibooker.book2D(histo, histo, 100, -5., 995., 80, -80., 80.);
0508
0509 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_HO");
0510 meTEprofileHO = ibooker.bookProfile(histo, histo, 60, -5., 55., -80., 80., " ");
0511
0512 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_High_HO");
0513 meTEprofileHO_High = ibooker.bookProfile(histo, histo, 100, -5., 995., -80., 80., " ");
0514 }
0515
0516
0517 if (subdet_ == 4 || subdet_ == 5) {
0518
0519
0520
0521 sprintf(histo, "HcalRecHitTask_severityLevel_HF");
0522 sevLvl_HF = ibooker.book1DD(histo, histo, 25, -0.5, 24.5);
0523
0524 sprintf(histo, "HcalRecHitTask_energy_of_rechits_HF");
0525 meRecHitsEnergyHF = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0526
0527 sprintf(histo, "HcalRecHitTask_cleaned_energy_of_rechits_HF");
0528 meRecHitsCleanedEnergyHF = ibooker.book1DD(histo, histo, 2010, -10., 2000.);
0529
0530 sprintf(histo, "HcalRecHitTask_timing_HF");
0531 meTimeHF = ibooker.book1DD(histo, histo, 70, -48., 92.);
0532
0533 sprintf(histo, "HcalRecHitTask_timing_vs_energy_Low_HF");
0534 meTE_Low_HF = ibooker.book2D(histo, histo, 100, -5., 195., 70, -48., 92.);
0535
0536 sprintf(histo, "HcalRecHitTask_timing_vs_energy_HF");
0537 meTE_HF = ibooker.book2D(histo, histo, 200, -5., 995., 70, -48., 92.);
0538
0539 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_Low_HF");
0540 meTEprofileHF_Low = ibooker.bookProfile(histo, histo, 100, -5., 195., -48., 92., " ");
0541
0542 sprintf(histo, "HcalRecHitTask_timing_vs_energy_profile_HF");
0543 meTEprofileHF = ibooker.bookProfile(histo, histo, 200, -5., 995., -48., 92., " ");
0544 }
0545 }
0546
0547 void HcalRecHitsAnalyzer::analyze(edm::Event const &ev, edm::EventSetup const &iSetup) {
0548 using namespace edm;
0549
0550
0551
0552
0553
0554 int nrechits = 0;
0555 int nrechitsThresh = 0;
0556
0557
0558 double eEcalCone = 0.;
0559
0560
0561 double partR = 0.3;
0562
0563
0564
0565 double etaHot = 99999.;
0566 double phiHot = 99999.;
0567
0568
0569 geometry = &iSetup.getData(caloGeometryEventToken_);
0570
0571
0572 theHcalTopology = &iSetup.getData(hcalTopologyToken_);
0573
0574
0575 theHcalChStatus = &iSetup.getData(hcalChannelQualityToken_);
0576
0577
0578 theHcalSevLvlComputer = &iSetup.getData(hcalSeverityLevelComputerToken_);
0579
0580
0581 fillRecHitsTmp(subdet_, ev);
0582
0583
0584 if (subdet_ == 5 || subdet_ == 1) {
0585 for (unsigned int iv = 0; iv < hcalHBSevLvlVec.size(); iv++) {
0586 sevLvl_HB->Fill(hcalHBSevLvlVec[iv]);
0587 }
0588 }
0589
0590 if (subdet_ == 5 || subdet_ == 2) {
0591 for (unsigned int iv = 0; iv < hcalHESevLvlVec.size(); iv++) {
0592 sevLvl_HE->Fill(hcalHESevLvlVec[iv]);
0593 }
0594 }
0595
0596 if (subdet_ == 5 || subdet_ == 3) {
0597 for (unsigned int iv = 0; iv < hcalHOSevLvlVec.size(); iv++) {
0598 sevLvl_HO->Fill(hcalHOSevLvlVec[iv]);
0599 }
0600 }
0601
0602 if (subdet_ == 5 || subdet_ == 4) {
0603 for (unsigned int iv = 0; iv < hcalHFSevLvlVec.size(); iv++) {
0604 sevLvl_HF->Fill(hcalHFSevLvlVec[iv]);
0605 }
0606 }
0607
0608
0609
0610
0611
0612
0613 std::vector<int> nhb_v(maxDepthHB_ + 1, 0);
0614 std::vector<int> nhe_v(maxDepthHE_ + 1, 0);
0615 std::vector<int> nho_v(maxDepthHO_ + 1, 0);
0616 std::vector<int> nhf_v(maxDepthHF_ + 1, 0);
0617
0618 for (unsigned int i = 0; i < cen.size(); i++) {
0619 int sub = csub[i];
0620 int depth = cdepth[i];
0621 int ieta = cieta[i];
0622 int iphi = ciphi[i];
0623 double en = cen[i];
0624 double enM0 = cenM0[i];
0625 double enM3 = cenM3[i];
0626
0627
0628 uint32_t stwd = cstwd[i];
0629 uint32_t auxstwd = cauxstwd[i];
0630
0631
0632
0633 bool isHEP17 = (sub == 2) && (iphi >= 63) && (iphi <= 66) && (ieta > 0) && (hep17_);
0634
0635
0636
0637 if (depth < 1)
0638 continue;
0639 if (sub == 1 && depth > maxDepthHB_)
0640 continue;
0641 if (sub == 2 && depth > maxDepthHE_)
0642 continue;
0643 if (sub == 3 && depth > maxDepthHO_)
0644 continue;
0645 if (sub == 4 && depth > maxDepthHF_)
0646 continue;
0647
0648 if (sub == 1) {
0649 nhb_v[depth]++;
0650 nhb_v[0]++;
0651 }
0652 if (sub == 2) {
0653 nhe_v[depth]++;
0654 nhe_v[0]++;
0655 }
0656 if (sub == 3) {
0657 nho_v[depth]++;
0658 nho_v[0]++;
0659 }
0660 if (sub == 4) {
0661 nhf_v[depth]++;
0662 nhf_v[0]++;
0663 }
0664
0665 if (subdet_ == 6) {
0666 }
0667
0668 if (subdet_ != 6) {
0669 int ieta2 = ieta;
0670 int depth2 = depth;
0671 if (sub == 4) {
0672 if (ieta2 < 0)
0673 ieta2--;
0674 else
0675 ieta2++;
0676 }
0677 if (sub == 3)
0678 emap_HO->Fill(double(ieta2), double(iphi), en);
0679 else
0680 emap[depth2 - 1]->Fill(double(ieta2), double(iphi), en);
0681
0682
0683 if (depth == 1 || depth == 2) {
0684 int ieta1 = ieta;
0685 if (sub == 4) {
0686 if (ieta1 < 0)
0687 ieta1--;
0688 else
0689 ieta1++;
0690 }
0691 }
0692
0693 if (sub == 1) {
0694 emean_vs_ieta_HB[depth - 1]->Fill(double(ieta), en);
0695 emean_vs_ieta_HBM0[depth - 1]->Fill(double(ieta), enM0);
0696 emean_vs_ieta_HBM3[depth - 1]->Fill(double(ieta), enM3);
0697 occupancy_map_HB[depth - 1]->Fill(double(ieta), double(iphi));
0698 if (ieta > 0)
0699 nrechits_vs_iphi_HBP[depth - 1]->Fill(double(iphi));
0700 else
0701 nrechits_vs_iphi_HBM[depth - 1]->Fill(double(iphi));
0702 }
0703 if (sub == 2) {
0704 if (!isHEP17) {
0705 emean_vs_ieta_HE[depth - 1]->Fill(double(ieta), en);
0706 emean_vs_ieta_HEM0[depth - 1]->Fill(double(ieta), enM0);
0707 emean_vs_ieta_HEM3[depth - 1]->Fill(double(ieta), enM3);
0708 } else {
0709 emean_vs_ieta_HEP17[depth - 1]->Fill(double(ieta), en);
0710 emean_vs_ieta_HEP17M0[depth - 1]->Fill(double(ieta), enM0);
0711 emean_vs_ieta_HEP17M3[depth - 1]->Fill(double(ieta), enM3);
0712 }
0713 occupancy_map_HE[depth - 1]->Fill(double(ieta), double(iphi));
0714 if (ieta > 0)
0715 nrechits_vs_iphi_HEP[depth - 1]->Fill(double(iphi));
0716 else
0717 nrechits_vs_iphi_HEM[depth - 1]->Fill(double(iphi));
0718 }
0719 if (sub == 3) {
0720 emean_vs_ieta_HO->Fill(double(ieta), en);
0721 occupancy_map_HO->Fill(double(ieta), double(iphi));
0722 if (ieta > 0)
0723 nrechits_vs_iphi_HOP->Fill(double(iphi));
0724 else
0725 nrechits_vs_iphi_HOM->Fill(double(iphi));
0726 }
0727 if (sub == 4) {
0728 emean_vs_ieta_HF[depth - 1]->Fill(double(ieta), en);
0729 occupancy_map_HF[depth - 1]->Fill(double(ieta), double(iphi));
0730 if (ieta > 0)
0731 nrechits_vs_iphi_HFP[depth - 1]->Fill(double(iphi));
0732 else
0733 nrechits_vs_iphi_HFM[depth - 1]->Fill(double(iphi));
0734 }
0735 }
0736
0737
0738 uint32_t statadd;
0739
0740
0741 unsigned int sw27 = 27;
0742 unsigned int sw13 = 13;
0743
0744 uint32_t statadd27 = 0x1 << sw27;
0745 uint32_t statadd13 = 0x1 << sw13;
0746
0747 float status27 = 0;
0748 float status13 = 0;
0749
0750 if (stwd & statadd27)
0751 status27 = 1;
0752 if (stwd & statadd13)
0753 status13 = 1;
0754
0755 if (sub == 1) {
0756 RecHit_StatusWordCorr_HB->Fill(status13, status27);
0757 } else if (sub == 2) {
0758 RecHit_StatusWordCorr_HE->Fill(status13, status27);
0759 }
0760
0761 for (unsigned int isw = 0; isw < 32; isw++) {
0762 statadd = 0x1 << (isw);
0763 if (stwd & statadd) {
0764 if (sub == 1)
0765 RecHit_StatusWord_HB->Fill(isw);
0766 else if (sub == 2)
0767 RecHit_StatusWord_HE->Fill(isw);
0768 else if (sub == 3)
0769 RecHit_StatusWord_HO->Fill(isw);
0770 else if (sub == 4) {
0771 RecHit_StatusWord_HF->Fill(isw);
0772 }
0773 }
0774 }
0775
0776 for (unsigned int isw = 0; isw < 32; isw++) {
0777 statadd = 0x1 << (isw);
0778 if (auxstwd & statadd) {
0779 if (sub == 1)
0780 RecHit_Aux_StatusWord_HB->Fill(isw);
0781 else if (sub == 2)
0782 RecHit_Aux_StatusWord_HE->Fill(isw);
0783 else if (sub == 3)
0784 RecHit_Aux_StatusWord_HO->Fill(isw);
0785 else if (sub == 4)
0786 RecHit_Aux_StatusWord_HF->Fill(isw);
0787 }
0788 }
0789 }
0790
0791 for (int depth = 0; depth <= maxDepthHB_; depth++)
0792 Nhb[depth]->Fill(double(nhb_v[depth]));
0793 for (int depth = 0; depth <= maxDepthHE_; depth++)
0794 Nhe[depth]->Fill(double(nhe_v[depth]));
0795 for (int depth = 0; depth <= maxDepthHO_; depth++)
0796 Nho[depth]->Fill(double(nho_v[depth]));
0797 for (int depth = 0; depth <= maxDepthHF_; depth++)
0798 Nhf[depth]->Fill(double(nhf_v[depth]));
0799
0800
0801
0802
0803
0804 if ((subdet_ != 6) && (subdet_ != 0)) {
0805 double clusEta = 999.;
0806 double clusPhi = 999.;
0807 double clusEn = 0.;
0808
0809 double HcalCone = 0.;
0810
0811 int ietaMax = 9999;
0812
0813
0814
0815
0816
0817
0818
0819
0820
0821 for (unsigned int i = 0; i < cen.size(); i++) {
0822 int sub = csub[i];
0823 double eta = ceta[i];
0824 double phi = cphi[i];
0825 double ieta = cieta[i];
0826 double iphi = ciphi[i];
0827 double en = cen[i];
0828 double enM0 = cenM0[i];
0829 double enM3 = cenM3[i];
0830 double chi2 = cchi2[i];
0831 double chi2_log10 = 9.99;
0832 if (chi2 > 0.)
0833 chi2_log10 = log10(chi2);
0834 double t = ctime[i];
0835 double depth = cdepth[i];
0836 int sevlev = csevlev[i];
0837
0838 bool isHEP17 = (sub == 2) && (iphi >= 63) && (iphi <= 66) && (ieta > 0) && (hep17_);
0839
0840
0841
0842 double rhot = dR(etaHot, phiHot, eta, phi);
0843 if (rhot < partR && en > 1.) {
0844 clusEta = (clusEta * clusEn + eta * en) / (clusEn + en);
0845 clusPhi = phi12(clusPhi, clusEn, phi, en);
0846 clusEn += en;
0847 }
0848
0849 nrechits++;
0850
0851 if (en > 1.)
0852 nrechitsThresh++;
0853
0854
0855
0856 if (sub == 1 && (subdet_ == 1 || subdet_ == 5)) {
0857 meTimeHB->Fill(t);
0858 meRecHitsEnergyHB->Fill(en);
0859 if (sevlev <= 9)
0860 meRecHitsCleanedEnergyHB->Fill(en);
0861
0862 meRecHitsEnergyHBM0->Fill(enM0);
0863 meRecHitsEnergyHBM3->Fill(enM3);
0864
0865 meRecHitsEnergyM2vM0HB->Fill(enM0, en);
0866 meRecHitsEnergyM3vM0HB->Fill(enM0, enM3);
0867 meRecHitsEnergyM3vM2HB->Fill(en, enM3);
0868
0869 meRecHitsM2Chi2HB->Fill(chi2_log10);
0870 meLog10Chi2profileHB->Fill(en, chi2_log10);
0871
0872 meTE_Low_HB->Fill(en, t);
0873 meTE_HB->Fill(en, t);
0874 meTE_High_HB->Fill(en, t);
0875 meTEprofileHB_Low->Fill(en, t);
0876 meTEprofileHB->Fill(en, t);
0877 meTEprofileHB_High->Fill(en, t);
0878 }
0879 if (sub == 2 && (subdet_ == 2 || subdet_ == 5)) {
0880 meTimeHE->Fill(t);
0881 if (!isHEP17) {
0882 meRecHitsEnergyHE->Fill(en);
0883 if (sevlev <= 9)
0884 meRecHitsCleanedEnergyHE->Fill(en);
0885
0886 meRecHitsEnergyHEM0->Fill(enM0);
0887 meRecHitsEnergyHEM3->Fill(enM3);
0888 } else {
0889 meRecHitsEnergyHEP17[0]->Fill(en);
0890 meRecHitsEnergyHEP17M0[0]->Fill(enM0);
0891 meRecHitsEnergyHEP17M3[0]->Fill(enM3);
0892 meRecHitsEnergyHEP17[depth]->Fill(en);
0893 meRecHitsEnergyHEP17M0[depth]->Fill(enM0);
0894 meRecHitsEnergyHEP17M3[depth]->Fill(enM3);
0895 }
0896
0897 meRecHitsEnergyM2vM0HE->Fill(enM0, en);
0898 meRecHitsEnergyM3vM0HE->Fill(enM0, enM3);
0899 meRecHitsEnergyM3vM2HE->Fill(en, enM3);
0900
0901 meRecHitsM2Chi2HE->Fill(chi2_log10);
0902 meLog10Chi2profileHE->Fill(en, chi2_log10);
0903
0904 meTE_Low_HE->Fill(en, t);
0905 meTE_HE->Fill(en, t);
0906 meTEprofileHE_Low->Fill(en, t);
0907 meTEprofileHE->Fill(en, t);
0908 }
0909 if (sub == 4 && (subdet_ == 4 || subdet_ == 5)) {
0910 meTimeHF->Fill(t);
0911 meRecHitsEnergyHF->Fill(en);
0912 if (sevlev <= 9)
0913 meRecHitsCleanedEnergyHF->Fill(en);
0914
0915 meTE_Low_HF->Fill(en, t);
0916 meTE_HF->Fill(en, t);
0917 meTEprofileHF_Low->Fill(en, t);
0918 meTEprofileHF->Fill(en, t);
0919 }
0920 if (sub == 3 && (subdet_ == 3 || subdet_ == 5)) {
0921 meTimeHO->Fill(t);
0922 meRecHitsEnergyHO->Fill(en);
0923 if (sevlev <= 9)
0924 meRecHitsCleanedEnergyHO->Fill(en);
0925
0926 meTE_HO->Fill(en, t);
0927 meTE_High_HO->Fill(en, t);
0928 meTEprofileHO->Fill(en, t);
0929 meTEprofileHO_High->Fill(en, t);
0930 }
0931 }
0932
0933 if (imc != 0) {
0934
0935 meEnConeEtaProfile->Fill(double(ietaMax), HcalCone);
0936 meEnConeEtaProfile_E->Fill(double(ietaMax), eEcalCone);
0937 meEnConeEtaProfile_EH->Fill(double(ietaMax), HcalCone + eEcalCone);
0938 }
0939
0940
0941
0942
0943 }
0944
0945 nevtot++;
0946 }
0947
0948
0949 void HcalRecHitsAnalyzer::fillRecHitsTmp(int subdet_, edm::Event const &ev) {
0950 using namespace edm;
0951
0952
0953 csub.clear();
0954 cen.clear();
0955 cenM0.clear();
0956 cenM3.clear();
0957 cchi2.clear();
0958 ceta.clear();
0959 cphi.clear();
0960 ctime.clear();
0961 cieta.clear();
0962 ciphi.clear();
0963 cdepth.clear();
0964 cz.clear();
0965 cstwd.clear();
0966 cauxstwd.clear();
0967 csevlev.clear();
0968 hcalHBSevLvlVec.clear();
0969 hcalHESevLvlVec.clear();
0970 hcalHFSevLvlVec.clear();
0971 hcalHOSevLvlVec.clear();
0972
0973 if (subdet_ == 1 || subdet_ == 2 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
0974
0975 edm::Handle<HBHERecHitCollection> hbhecoll;
0976 if (ev.getByToken(tok_hbhe_, hbhecoll)) {
0977 for (HBHERecHitCollection::const_iterator j = hbhecoll->begin(); j != hbhecoll->end(); j++) {
0978 HcalDetId cell(j->id());
0979 const HcalGeometry *cellGeometry = dynamic_cast<const HcalGeometry *>(geometry->getSubdetectorGeometry(cell));
0980 double eta = cellGeometry->getPosition(cell).eta();
0981 double phi = cellGeometry->getPosition(cell).phi();
0982 double zc = cellGeometry->getPosition(cell).z();
0983 int sub = cell.subdet();
0984 int depth = cell.depth();
0985 int inteta = cell.ieta();
0986 int intphi = cell.iphi();
0987 double en = j->energy();
0988 double enM0 = j->eraw();
0989 double enM3 = j->eaux();
0990 double chi2 = j->chi2();
0991 double t = j->time();
0992 int stwd = j->flags();
0993 int auxstwd = j->aux();
0994
0995 int severityLevel = hcalSevLvl((CaloRecHit *)&*j);
0996 if (cell.subdet() == HcalBarrel) {
0997 hcalHBSevLvlVec.push_back(severityLevel);
0998 } else if (cell.subdet() == HcalEndcap) {
0999 hcalHESevLvlVec.push_back(severityLevel);
1000 }
1001
1002 if ((iz > 0 && eta > 0.) || (iz < 0 && eta < 0.) || iz == 0) {
1003 csub.push_back(sub);
1004 cen.push_back(en);
1005 cenM0.push_back(enM0);
1006 cenM3.push_back(enM3);
1007 cchi2.push_back(chi2);
1008 ceta.push_back(eta);
1009 cphi.push_back(phi);
1010 ctime.push_back(t);
1011 cieta.push_back(inteta);
1012 ciphi.push_back(intphi);
1013 cdepth.push_back(depth);
1014 cz.push_back(zc);
1015 cstwd.push_back(stwd);
1016 cauxstwd.push_back(auxstwd);
1017 csevlev.push_back(severityLevel);
1018 }
1019 }
1020 }
1021 }
1022
1023 if (subdet_ == 4 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
1024
1025 edm::Handle<HFRecHitCollection> hfcoll;
1026 if (ev.getByToken(tok_hf_, hfcoll)) {
1027 for (HFRecHitCollection::const_iterator j = hfcoll->begin(); j != hfcoll->end(); j++) {
1028 HcalDetId cell(j->id());
1029 auto cellGeometry = (geometry->getSubdetectorGeometry(cell))->getGeometry(cell);
1030 double eta = cellGeometry->getPosition().eta();
1031 double phi = cellGeometry->getPosition().phi();
1032 double zc = cellGeometry->getPosition().z();
1033 int sub = cell.subdet();
1034 int depth = cell.depth();
1035 int inteta = cell.ieta();
1036 int intphi = cell.iphi();
1037 double en = j->energy();
1038 double enM0 = 0.;
1039 double enM3 = 0.;
1040 double chi2 = 0.;
1041 double t = j->time();
1042 int stwd = j->flags();
1043 int auxstwd = j->aux();
1044
1045 int severityLevel = hcalSevLvl((CaloRecHit *)&*j);
1046 if (cell.subdet() == HcalForward) {
1047 hcalHFSevLvlVec.push_back(severityLevel);
1048 }
1049
1050 if ((iz > 0 && eta > 0.) || (iz < 0 && eta < 0.) || iz == 0) {
1051 csub.push_back(sub);
1052 cen.push_back(en);
1053 cenM0.push_back(enM0);
1054 cenM3.push_back(enM3);
1055 cchi2.push_back(chi2);
1056 ceta.push_back(eta);
1057 cphi.push_back(phi);
1058 ctime.push_back(t);
1059 cieta.push_back(inteta);
1060 ciphi.push_back(intphi);
1061 cdepth.push_back(depth);
1062 cz.push_back(zc);
1063 cstwd.push_back(stwd);
1064 cauxstwd.push_back(auxstwd);
1065 csevlev.push_back(severityLevel);
1066 }
1067 }
1068 }
1069 }
1070
1071
1072 if (subdet_ == 3 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
1073 edm::Handle<HORecHitCollection> hocoll;
1074 if (ev.getByToken(tok_ho_, hocoll)) {
1075 for (HORecHitCollection::const_iterator j = hocoll->begin(); j != hocoll->end(); j++) {
1076 HcalDetId cell(j->id());
1077 auto cellGeometry = (geometry->getSubdetectorGeometry(cell))->getGeometry(cell);
1078 double eta = cellGeometry->getPosition().eta();
1079 double phi = cellGeometry->getPosition().phi();
1080 double zc = cellGeometry->getPosition().z();
1081 int sub = cell.subdet();
1082 int depth = cell.depth();
1083 int inteta = cell.ieta();
1084 int intphi = cell.iphi();
1085 double t = j->time();
1086 double en = j->energy();
1087 double enM0 = 0.;
1088 double enM3 = 0.;
1089 double chi2 = 0.;
1090 int stwd = j->flags();
1091 int auxstwd = j->aux();
1092
1093 int severityLevel = hcalSevLvl((CaloRecHit *)&*j);
1094 if (cell.subdet() == HcalOuter) {
1095 hcalHOSevLvlVec.push_back(severityLevel);
1096 }
1097
1098 if ((iz > 0 && eta > 0.) || (iz < 0 && eta < 0.) || iz == 0) {
1099 csub.push_back(sub);
1100 cen.push_back(en);
1101 cenM0.push_back(enM0);
1102 cenM3.push_back(enM3);
1103 cchi2.push_back(chi2);
1104 ceta.push_back(eta);
1105 cphi.push_back(phi);
1106 ctime.push_back(t);
1107 cieta.push_back(inteta);
1108 ciphi.push_back(intphi);
1109 cdepth.push_back(depth);
1110 cz.push_back(zc);
1111 cstwd.push_back(stwd);
1112 cauxstwd.push_back(auxstwd);
1113 csevlev.push_back(severityLevel);
1114 }
1115 }
1116 }
1117 }
1118 }
1119
1120 double HcalRecHitsAnalyzer::dR(double eta1, double phi1, double eta2, double phi2) {
1121 double PI = 3.1415926535898;
1122 double deltaphi = phi1 - phi2;
1123 if (phi2 > phi1) {
1124 deltaphi = phi2 - phi1;
1125 }
1126 if (deltaphi > PI) {
1127 deltaphi = 2. * PI - deltaphi;
1128 }
1129 double deltaeta = eta2 - eta1;
1130 double tmp = sqrt(deltaeta * deltaeta + deltaphi * deltaphi);
1131 return tmp;
1132 }
1133
1134 double HcalRecHitsAnalyzer::phi12(double phi1, double en1, double phi2, double en2) {
1135
1136
1137 double tmp;
1138 double PI = 3.1415926535898;
1139 double a1 = phi1;
1140 double a2 = phi2;
1141
1142 if (a1 > 0.5 * PI && a2 < 0.)
1143 a2 += 2 * PI;
1144 if (a2 > 0.5 * PI && a1 < 0.)
1145 a1 += 2 * PI;
1146 tmp = (a1 * en1 + a2 * en2) / (en1 + en2);
1147 if (tmp > PI)
1148 tmp -= 2. * PI;
1149
1150 return tmp;
1151 }
1152
1153 double HcalRecHitsAnalyzer::dPhiWsign(double phi1, double phi2) {
1154
1155
1156
1157 double PI = 3.1415926535898;
1158 double a1 = phi1;
1159 double a2 = phi2;
1160 double tmp = a2 - a1;
1161 if (a1 * a2 < 0.) {
1162 if (a1 > 0.5 * PI)
1163 tmp += 2. * PI;
1164 if (a2 > 0.5 * PI)
1165 tmp -= 2. * PI;
1166 }
1167 return tmp;
1168 }
1169
1170 int HcalRecHitsAnalyzer::hcalSevLvl(const CaloRecHit *hit) {
1171 HcalDetId id = hit->detid();
1172 if (theHcalTopology->getMergePositionFlag() && id.subdet() == HcalEndcap) {
1173 id = theHcalTopology->idFront(id);
1174 }
1175
1176 const uint32_t recHitFlag = hit->flags();
1177 const uint32_t dbStatusFlag = theHcalChStatus->getValues(id)->getValue();
1178
1179 int severityLevel = theHcalSevLvlComputer->getSeverityLevel(id, recHitFlag, dbStatusFlag);
1180
1181 return severityLevel;
1182 }
1183
1184 DEFINE_FWK_MODULE(HcalRecHitsAnalyzer);