File indexing completed on 2023-03-17 11:27:35
0001
0002
0003
0004
0005
0006
0007
0008 #include "Validation/GlobalHits/interface/GlobalHitsHistogrammer.h"
0009
0010
0011 GlobalHitsHistogrammer::GlobalHitsHistogrammer(const edm::ParameterSet &iPSet)
0012 : fName(""),
0013 verbosity(0),
0014 frequency(0),
0015 vtxunit(0),
0016 label(""),
0017 getAllProvenances(false),
0018 printProvenanceInfo(false),
0019 count(0) {
0020 std::string MsgLoggerCat = "GlobalHitsHistogrammer_GlobalHitsHistogrammer";
0021
0022
0023 fName = iPSet.getUntrackedParameter<std::string>("Name");
0024 verbosity = iPSet.getUntrackedParameter<int>("Verbosity");
0025 frequency = iPSet.getUntrackedParameter<int>("Frequency");
0026 vtxunit = iPSet.getUntrackedParameter<int>("VtxUnit");
0027 outputfile = iPSet.getParameter<std::string>("OutputFile");
0028 doOutput = iPSet.getParameter<bool>("DoOutput");
0029 edm::ParameterSet m_Prov = iPSet.getParameter<edm::ParameterSet>("ProvenanceLookup");
0030 getAllProvenances = m_Prov.getUntrackedParameter<bool>("GetAllProvenances");
0031 printProvenanceInfo = m_Prov.getUntrackedParameter<bool>("PrintProvenanceInfo");
0032
0033
0034 GlobalHitSrc_ = iPSet.getParameter<edm::InputTag>("GlobalHitSrc");
0035
0036 GlobalHitSrc_Token_ = consumes<PGlobalSimHit>(iPSet.getParameter<edm::InputTag>("GlobalHitSrc"));
0037
0038
0039
0040 verbosity %= 10;
0041
0042
0043 if (verbosity >= 0) {
0044 edm::LogInfo(MsgLoggerCat) << "\n===============================\n"
0045 << "Initialized as EDAnalyzer with parameter values:\n"
0046 << " Name = " << fName << "\n"
0047 << " Verbosity = " << verbosity << "\n"
0048 << " Frequency = " << frequency << "\n"
0049 << " VtxUnit = " << vtxunit << "\n"
0050 << " OutputFile = " << outputfile << "\n"
0051 << " DoOutput = " << doOutput << "\n"
0052 << " GetProv = " << getAllProvenances << "\n"
0053 << " PrintProv = " << printProvenanceInfo << "\n"
0054 << " GlobalHitSrc = " << GlobalHitSrc_.label() << ":" << GlobalHitSrc_.instance()
0055 << "\n"
0056 << "===============================\n";
0057 }
0058 }
0059
0060 void GlobalHitsHistogrammer::bookHistograms(DQMStore::IBooker &ibooker, edm::Run const &, edm::EventSetup const &) {
0061
0062 for (Int_t i = 0; i < 2; ++i) {
0063 meMCRGP[i] = nullptr;
0064 meMCG4Vtx[i] = nullptr;
0065 meGeantVtxX[i] = nullptr;
0066 meGeantVtxY[i] = nullptr;
0067 meGeantVtxZ[i] = nullptr;
0068 meMCG4Trk[i] = nullptr;
0069 meCaloEcal[i] = nullptr;
0070 meCaloEcalE[i] = nullptr;
0071 meCaloEcalToF[i] = nullptr;
0072 meCaloPreSh[i] = nullptr;
0073 meCaloPreShE[i] = nullptr;
0074 meCaloPreShToF[i] = nullptr;
0075 meCaloHcal[i] = nullptr;
0076 meCaloHcalE[i] = nullptr;
0077 meCaloHcalToF[i] = nullptr;
0078 meTrackerPx[i] = nullptr;
0079 meTrackerSi[i] = nullptr;
0080 meMuon[i] = nullptr;
0081 meMuonDtToF[i] = nullptr;
0082 meMuonCscToF[i] = nullptr;
0083 meMuonRpcFToF[i] = nullptr;
0084 meMuonRpcBToF[i] = nullptr;
0085 }
0086 meGeantTrkPt = nullptr;
0087 meGeantTrkE = nullptr;
0088 meCaloEcalPhi = nullptr;
0089 meCaloEcalEta = nullptr;
0090 meCaloPreShPhi = nullptr;
0091 meCaloPreShEta = nullptr;
0092 meCaloHcalPhi = nullptr;
0093 meCaloHcalEta = nullptr;
0094 meTrackerPxPhi = nullptr;
0095 meTrackerPxEta = nullptr;
0096 meTrackerPxBToF = nullptr;
0097 meTrackerPxBR = nullptr;
0098 meTrackerPxFToF = nullptr;
0099 meTrackerPxFZ = nullptr;
0100 meTrackerSiPhi = nullptr;
0101 meTrackerSiEta = nullptr;
0102 meTrackerSiBToF = nullptr;
0103 meTrackerSiBR = nullptr;
0104 meTrackerSiFToF = nullptr;
0105 meTrackerSiFZ = nullptr;
0106 meMuonPhi = nullptr;
0107 meMuonEta = nullptr;
0108 meMuonDtR = nullptr;
0109 meMuonCscZ = nullptr;
0110 meMuonRpcBR = nullptr;
0111 meMuonRpcFZ = nullptr;
0112
0113
0114 Char_t hname[200];
0115 Char_t htitle[200];
0116
0117
0118 ibooker.setCurrentFolder("GlobalHitsV/MCGeant");
0119 sprintf(hname, "hMCRGP1");
0120 sprintf(htitle, "RawGenParticles");
0121 meMCRGP[0] = ibooker.book1D(hname, htitle, 100, 0., 5000.);
0122 sprintf(hname, "hMCRGP2");
0123 meMCRGP[1] = ibooker.book1D(hname, htitle, 100, 0., 500.);
0124 for (Int_t i = 0; i < 2; ++i) {
0125 meMCRGP[i]->setAxisTitle("Number of Raw Generated Particles", 1);
0126 meMCRGP[i]->setAxisTitle("Count", 2);
0127 }
0128
0129 sprintf(hname, "hMCG4Vtx1");
0130 sprintf(htitle, "G4 Vertices");
0131 meMCG4Vtx[0] = ibooker.book1D(hname, htitle, 100, 0., 50000.);
0132 sprintf(hname, "hMCG4Vtx2");
0133 meMCG4Vtx[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0134 for (Int_t i = 0; i < 2; ++i) {
0135 meMCG4Vtx[i]->setAxisTitle("Number of Vertices", 1);
0136 meMCG4Vtx[i]->setAxisTitle("Count", 2);
0137 }
0138
0139 sprintf(hname, "hMCG4Trk1");
0140 sprintf(htitle, "G4 Tracks");
0141 meMCG4Trk[0] = ibooker.book1D(hname, htitle, 150, 0., 15000.);
0142 sprintf(hname, "hMCG4Trk2");
0143 meMCG4Trk[1] = ibooker.book1D(hname, htitle, 150, -0.5, 99.5);
0144 for (Int_t i = 0; i < 2; ++i) {
0145 meMCG4Trk[i]->setAxisTitle("Number of Tracks", 1);
0146 meMCG4Trk[i]->setAxisTitle("Count", 2);
0147 }
0148
0149 sprintf(hname, "hGeantVtxX1");
0150 sprintf(htitle, "Geant vertex x/micrometer");
0151 meGeantVtxX[0] = ibooker.book1D(hname, htitle, 100, -8000000., 8000000.);
0152 sprintf(hname, "hGeantVtxX2");
0153 meGeantVtxX[1] = ibooker.book1D(hname, htitle, 100, -50., 50.);
0154 for (Int_t i = 0; i < 2; ++i) {
0155 meGeantVtxX[i]->setAxisTitle("x of Vertex (um)", 1);
0156 meGeantVtxX[i]->setAxisTitle("Count", 2);
0157 }
0158
0159 sprintf(hname, "hGeantVtxY1");
0160 sprintf(htitle, "Geant vertex y/micrometer");
0161 meGeantVtxY[0] = ibooker.book1D(hname, htitle, 100, -8000000, 8000000.);
0162 sprintf(hname, "hGeantVtxY2");
0163 meGeantVtxY[1] = ibooker.book1D(hname, htitle, 100, -50., 50.);
0164 for (Int_t i = 0; i < 2; ++i) {
0165 meGeantVtxY[i]->setAxisTitle("y of Vertex (um)", 1);
0166 meGeantVtxY[i]->setAxisTitle("Count", 2);
0167 }
0168
0169 sprintf(hname, "hGeantVtxZ1");
0170 sprintf(htitle, "Geant vertex z/millimeter");
0171 meGeantVtxZ[0] = ibooker.book1D(hname, htitle, 100, -11000., 11000.);
0172 sprintf(hname, "hGeantVtxZ2");
0173 meGeantVtxZ[1] = ibooker.book1D(hname, htitle, 100, -250., 250.);
0174 for (Int_t i = 0; i < 2; ++i) {
0175 meGeantVtxZ[i]->setAxisTitle("z of Vertex (mm)", 1);
0176 meGeantVtxZ[i]->setAxisTitle("Count", 2);
0177 }
0178
0179 sprintf(hname, "hGeantTrkPt");
0180 sprintf(htitle, "Geant track pt/GeV");
0181 meGeantTrkPt = ibooker.book1D(hname, htitle, 100, 0., 200.);
0182 meGeantTrkPt->setAxisTitle("pT of Track (GeV)", 1);
0183 meGeantTrkPt->setAxisTitle("Count", 2);
0184
0185 sprintf(hname, "hGeantTrkE");
0186 sprintf(htitle, "Geant track E/GeV");
0187 meGeantTrkE = ibooker.book1D(hname, htitle, 100, 0., 5000.);
0188 meGeantTrkE->setAxisTitle("E of Track (GeV)", 1);
0189 meGeantTrkE->setAxisTitle("Count", 2);
0190
0191
0192 ibooker.setCurrentFolder("GlobalHitsV/ECals");
0193
0194 sprintf(hname, "hCaloEcal1");
0195 sprintf(htitle, "Ecal hits");
0196 meCaloEcal[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0197 sprintf(hname, "hCaloEcal2");
0198 meCaloEcal[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0199
0200 sprintf(hname, "hCaloEcalE1");
0201 sprintf(htitle, "Ecal hits, energy/GeV");
0202 meCaloEcalE[0] = ibooker.book1D(hname, htitle, 100, 0., 10.);
0203 sprintf(hname, "hCaloEcalE2");
0204 meCaloEcalE[1] = ibooker.book1D(hname, htitle, 100, 0., 0.1);
0205
0206 sprintf(hname, "hCaloEcalToF1");
0207 sprintf(htitle, "Ecal hits, ToF/ns");
0208 meCaloEcalToF[0] = ibooker.book1D(hname, htitle, 100, 0., 1000.);
0209 sprintf(hname, "hCaloEcalToF2");
0210 meCaloEcalToF[1] = ibooker.book1D(hname, htitle, 100, 0., 100.);
0211
0212 for (Int_t i = 0; i < 2; ++i) {
0213 meCaloEcal[i]->setAxisTitle("Number of Hits", 1);
0214 meCaloEcal[i]->setAxisTitle("Count", 2);
0215 meCaloEcalE[i]->setAxisTitle("Energy of Hits (GeV)", 1);
0216 meCaloEcalE[i]->setAxisTitle("Count", 2);
0217 meCaloEcalToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0218 meCaloEcalToF[i]->setAxisTitle("Count", 2);
0219 }
0220
0221 sprintf(hname, "hCaloEcalPhi");
0222 sprintf(htitle, "Ecal hits, phi/rad");
0223 meCaloEcalPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0224 meCaloEcalPhi->setAxisTitle("Phi of Hits (rad)", 1);
0225 meCaloEcalPhi->setAxisTitle("Count", 2);
0226
0227 sprintf(hname, "hCaloEcalEta");
0228 sprintf(htitle, "Ecal hits, eta");
0229 meCaloEcalEta = ibooker.book1D(hname, htitle, 100, -5.5, 5.5);
0230 meCaloEcalEta->setAxisTitle("Eta of Hits", 1);
0231 meCaloEcalEta->setAxisTitle("Count", 2);
0232
0233 sprintf(hname, "hCaloPreSh1");
0234 sprintf(htitle, "PreSh hits");
0235 meCaloPreSh[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0236 sprintf(hname, "hCaloPreSh2");
0237 meCaloPreSh[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0238
0239 sprintf(hname, "hCaloPreShE1");
0240 sprintf(htitle, "PreSh hits, energy/GeV");
0241 meCaloPreShE[0] = ibooker.book1D(hname, htitle, 100, 0., 10.);
0242 sprintf(hname, "hCaloPreShE2");
0243 meCaloPreShE[1] = ibooker.book1D(hname, htitle, 100, 0., 0.1);
0244
0245 sprintf(hname, "hCaloPreShToF1");
0246 sprintf(htitle, "PreSh hits, ToF/ns");
0247 meCaloPreShToF[0] = ibooker.book1D(hname, htitle, 100, 0., 1000.);
0248 sprintf(hname, "hCaloPreShToF2");
0249 meCaloPreShToF[1] = ibooker.book1D(hname, htitle, 100, 0., 100.);
0250
0251 for (Int_t i = 0; i < 2; ++i) {
0252 meCaloPreSh[i]->setAxisTitle("Number of Hits", 1);
0253 meCaloPreSh[i]->setAxisTitle("Count", 2);
0254 meCaloPreShE[i]->setAxisTitle("Energy of Hits (GeV)", 1);
0255 meCaloPreShE[i]->setAxisTitle("Count", 2);
0256 meCaloPreShToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0257 meCaloPreShToF[i]->setAxisTitle("Count", 2);
0258 }
0259
0260 sprintf(hname, "hCaloPreShPhi");
0261 sprintf(htitle, "PreSh hits, phi/rad");
0262 meCaloPreShPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0263 meCaloPreShPhi->setAxisTitle("Phi of Hits (rad)", 1);
0264 meCaloPreShPhi->setAxisTitle("Count", 2);
0265
0266 sprintf(hname, "hCaloPreShEta");
0267 sprintf(htitle, "PreSh hits, eta");
0268 meCaloPreShEta = ibooker.book1D(hname, htitle, 100, -5.5, 5.5);
0269 meCaloPreShEta->setAxisTitle("Eta of Hits", 1);
0270 meCaloPreShEta->setAxisTitle("Count", 2);
0271
0272
0273 ibooker.setCurrentFolder("GlobalHitsV/HCals");
0274 sprintf(hname, "hCaloHcal1");
0275 sprintf(htitle, "Hcal hits");
0276 meCaloHcal[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0277 sprintf(hname, "hCaloHcal2");
0278 meCaloHcal[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0279
0280 sprintf(hname, "hCaloHcalE1");
0281 sprintf(htitle, "Hcal hits, energy/GeV");
0282 meCaloHcalE[0] = ibooker.book1D(hname, htitle, 100, 0., 10.);
0283 sprintf(hname, "hCaloHcalE2");
0284 meCaloHcalE[1] = ibooker.book1D(hname, htitle, 100, 0., 0.1);
0285
0286 sprintf(hname, "hCaloHcalToF1");
0287 sprintf(htitle, "Hcal hits, ToF/ns");
0288 meCaloHcalToF[0] = ibooker.book1D(hname, htitle, 100, 0., 1000.);
0289 sprintf(hname, "hCaloHcalToF2");
0290 meCaloHcalToF[1] = ibooker.book1D(hname, htitle, 100, 0., 100.);
0291
0292 for (Int_t i = 0; i < 2; ++i) {
0293 meCaloHcal[i]->setAxisTitle("Number of Hits", 1);
0294 meCaloHcal[i]->setAxisTitle("Count", 2);
0295 meCaloHcalE[i]->setAxisTitle("Energy of Hits (GeV)", 1);
0296 meCaloHcalE[i]->setAxisTitle("Count", 2);
0297 meCaloHcalToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0298 meCaloHcalToF[i]->setAxisTitle("Count", 2);
0299 }
0300
0301 sprintf(hname, "hCaloHcalPhi");
0302 sprintf(htitle, "Hcal hits, phi/rad");
0303 meCaloHcalPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0304 meCaloHcalPhi->setAxisTitle("Phi of Hits (rad)", 1);
0305 meCaloHcalPhi->setAxisTitle("Count", 2);
0306
0307 sprintf(hname, "hCaloHcalEta");
0308 sprintf(htitle, "Hcal hits, eta");
0309 meCaloHcalEta = ibooker.book1D(hname, htitle, 100, -5.5, 5.5);
0310 meCaloHcalEta->setAxisTitle("Eta of Hits", 1);
0311 meCaloHcalEta->setAxisTitle("Count", 2);
0312
0313
0314 ibooker.setCurrentFolder("GlobalHitsV/SiPixels");
0315 sprintf(hname, "hTrackerPx1");
0316 sprintf(htitle, "Pixel hits");
0317 meTrackerPx[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0318 sprintf(hname, "hTrackerPx2");
0319 meTrackerPx[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0320 for (Int_t i = 0; i < 2; ++i) {
0321 meTrackerPx[i]->setAxisTitle("Number of Pixel Hits", 1);
0322 meTrackerPx[i]->setAxisTitle("Count", 2);
0323 }
0324
0325 sprintf(hname, "hTrackerPxPhi");
0326 sprintf(htitle, "Pixel hits phi/rad");
0327 meTrackerPxPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0328 meTrackerPxPhi->setAxisTitle("Phi of Hits (rad)", 1);
0329 meTrackerPxPhi->setAxisTitle("Count", 2);
0330
0331 sprintf(hname, "hTrackerPxEta");
0332 sprintf(htitle, "Pixel hits eta");
0333 meTrackerPxEta = ibooker.book1D(hname, htitle, 100, -3.5, 3.5);
0334 meTrackerPxEta->setAxisTitle("Eta of Hits", 1);
0335 meTrackerPxEta->setAxisTitle("Count", 2);
0336
0337 sprintf(hname, "hTrackerPxBToF");
0338 sprintf(htitle, "Pixel barrel hits, ToF/ns");
0339 meTrackerPxBToF = ibooker.book1D(hname, htitle, 100, 0., 40.);
0340 meTrackerPxBToF->setAxisTitle("Time of Flight of Hits (ns)", 1);
0341 meTrackerPxBToF->setAxisTitle("Count", 2);
0342
0343 sprintf(hname, "hTrackerPxBR");
0344 sprintf(htitle, "Pixel barrel hits, R/cm");
0345 meTrackerPxBR = ibooker.book1D(hname, htitle, 100, 0., 50.);
0346 meTrackerPxBR->setAxisTitle("R of Hits (cm)", 1);
0347 meTrackerPxBR->setAxisTitle("Count", 2);
0348
0349 sprintf(hname, "hTrackerPxFToF");
0350 sprintf(htitle, "Pixel forward hits, ToF/ns");
0351 meTrackerPxFToF = ibooker.book1D(hname, htitle, 100, 0., 50.);
0352 meTrackerPxFToF->setAxisTitle("Time of Flight of Hits (ns)", 1);
0353 meTrackerPxFToF->setAxisTitle("Count", 2);
0354
0355 sprintf(hname, "hTrackerPxFZ");
0356 sprintf(htitle, "Pixel forward hits, Z/cm");
0357 meTrackerPxFZ = ibooker.book1D(hname, htitle, 200, -100., 100.);
0358 meTrackerPxFZ->setAxisTitle("Z of Hits (cm)", 1);
0359 meTrackerPxFZ->setAxisTitle("Count", 2);
0360
0361
0362 ibooker.setCurrentFolder("GlobalHitsV/SiPixels");
0363 sprintf(hname, "hTrackerSi1");
0364 sprintf(htitle, "Silicon hits");
0365 meTrackerSi[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0366 sprintf(hname, "hTrackerSi2");
0367 meTrackerSi[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0368 for (Int_t i = 0; i < 2; ++i) {
0369 meTrackerSi[i]->setAxisTitle("Number of Silicon Hits", 1);
0370 meTrackerSi[i]->setAxisTitle("Count", 2);
0371 }
0372
0373 sprintf(hname, "hTrackerSiPhi");
0374 sprintf(htitle, "Silicon hits phi/rad");
0375 meTrackerSiPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0376 meTrackerSiPhi->setAxisTitle("Phi of Hits (rad)", 1);
0377 meTrackerSiPhi->setAxisTitle("Count", 2);
0378
0379 sprintf(hname, "hTrackerSiEta");
0380 sprintf(htitle, "Silicon hits eta");
0381 meTrackerSiEta = ibooker.book1D(hname, htitle, 100, -3.5, 3.5);
0382 meTrackerSiEta->setAxisTitle("Eta of Hits", 1);
0383 meTrackerSiEta->setAxisTitle("Count", 2);
0384
0385 sprintf(hname, "hTrackerSiBToF");
0386 sprintf(htitle, "Silicon barrel hits, ToF/ns");
0387 meTrackerSiBToF = ibooker.book1D(hname, htitle, 100, 0., 50.);
0388 meTrackerSiBToF->setAxisTitle("Time of Flight of Hits (ns)", 1);
0389 meTrackerSiBToF->setAxisTitle("Count", 2);
0390
0391 sprintf(hname, "hTrackerSiBR");
0392 sprintf(htitle, "Silicon barrel hits, R/cm");
0393 meTrackerSiBR = ibooker.book1D(hname, htitle, 100, 0., 200.);
0394 meTrackerSiBR->setAxisTitle("R of Hits (cm)", 1);
0395 meTrackerSiBR->setAxisTitle("Count", 2);
0396
0397 sprintf(hname, "hTrackerSiFToF");
0398 sprintf(htitle, "Silicon forward hits, ToF/ns");
0399 meTrackerSiFToF = ibooker.book1D(hname, htitle, 100, 0., 75.);
0400 meTrackerSiFToF->setAxisTitle("Time of Flight of Hits (ns)", 1);
0401 meTrackerSiFToF->setAxisTitle("Count", 2);
0402
0403 sprintf(hname, "hTrackerSiFZ");
0404 sprintf(htitle, "Silicon forward hits, Z/cm");
0405 meTrackerSiFZ = ibooker.book1D(hname, htitle, 200, -300., 300.);
0406 meTrackerSiFZ->setAxisTitle("Z of Hits (cm)", 1);
0407 meTrackerSiFZ->setAxisTitle("Count", 2);
0408
0409
0410 ibooker.setCurrentFolder("GlobalHitsV/Muons");
0411 sprintf(hname, "hMuon1");
0412 sprintf(htitle, "Muon hits");
0413 meMuon[0] = ibooker.book1D(hname, htitle, 100, 0., 10000.);
0414 sprintf(hname, "hMuon2");
0415 meMuon[1] = ibooker.book1D(hname, htitle, 100, -0.5, 99.5);
0416 for (Int_t i = 0; i < 2; ++i) {
0417 meMuon[i]->setAxisTitle("Number of Muon Hits", 1);
0418 meMuon[i]->setAxisTitle("Count", 2);
0419 }
0420
0421 sprintf(hname, "hMuonPhi");
0422 sprintf(htitle, "Muon hits phi/rad");
0423 meMuonPhi = ibooker.book1D(hname, htitle, 100, -3.2, 3.2);
0424 meMuonPhi->setAxisTitle("Phi of Hits (rad)", 1);
0425 meMuonPhi->setAxisTitle("Count", 2);
0426
0427 sprintf(hname, "hMuonEta");
0428 sprintf(htitle, "Muon hits eta");
0429 meMuonEta = ibooker.book1D(hname, htitle, 100, -3.5, 3.5);
0430 meMuonEta->setAxisTitle("Eta of Hits", 1);
0431 meMuonEta->setAxisTitle("Count", 2);
0432
0433 sprintf(hname, "hMuonCscToF1");
0434 sprintf(htitle, "Muon CSC hits, ToF/ns");
0435 meMuonCscToF[0] = ibooker.book1D(hname, htitle, 100, 0., 250.);
0436 sprintf(hname, "hMuonCscToF2");
0437 meMuonCscToF[1] = ibooker.book1D(hname, htitle, 100, 0., 50.);
0438 for (Int_t i = 0; i < 2; ++i) {
0439 meMuonCscToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0440 meMuonCscToF[i]->setAxisTitle("Count", 2);
0441 }
0442
0443 sprintf(hname, "hMuonCscZ");
0444 sprintf(htitle, "Muon CSC hits, Z/cm");
0445 meMuonCscZ = ibooker.book1D(hname, htitle, 200, -1500., 1500.);
0446 meMuonCscZ->setAxisTitle("Z of Hits (cm)", 1);
0447 meMuonCscZ->setAxisTitle("Count", 2);
0448
0449 sprintf(hname, "hMuonDtToF1");
0450 sprintf(htitle, "Muon DT hits, ToF/ns");
0451 meMuonDtToF[0] = ibooker.book1D(hname, htitle, 100, 0., 250.);
0452 sprintf(hname, "hMuonDtToF2");
0453 meMuonDtToF[1] = ibooker.book1D(hname, htitle, 100, 0., 50.);
0454 for (Int_t i = 0; i < 2; ++i) {
0455 meMuonDtToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0456 meMuonDtToF[i]->setAxisTitle("Count", 2);
0457 }
0458
0459 sprintf(hname, "hMuonDtR");
0460 sprintf(htitle, "Muon DT hits, R/cm");
0461 meMuonDtR = ibooker.book1D(hname, htitle, 100, 0., 1500.);
0462 meMuonDtR->setAxisTitle("R of Hits (cm)", 1);
0463 meMuonDtR->setAxisTitle("Count", 2);
0464
0465 sprintf(hname, "hMuonRpcFToF1");
0466 sprintf(htitle, "Muon RPC forward hits, ToF/ns");
0467 meMuonRpcFToF[0] = ibooker.book1D(hname, htitle, 100, 0., 250.);
0468 sprintf(hname, "hMuonRpcFToF2_4305");
0469 meMuonRpcFToF[1] = ibooker.book1D(hname, htitle, 100, 0., 50.);
0470 for (Int_t i = 0; i < 2; ++i) {
0471 meMuonRpcFToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0472 meMuonRpcFToF[i]->setAxisTitle("Count", 2);
0473 }
0474
0475 sprintf(hname, "hMuonRpcFZ");
0476 sprintf(htitle, "Muon RPC forward hits, Z/cm");
0477 meMuonRpcFZ = ibooker.book1D(hname, htitle, 201, -1500., 1500.);
0478 meMuonRpcFZ->setAxisTitle("Z of Hits (cm)", 1);
0479 meMuonRpcFZ->setAxisTitle("Count", 2);
0480
0481 sprintf(hname, "hMuonRpcBToF1");
0482 sprintf(htitle, "Muon RPC barrel hits, ToF/ns");
0483 meMuonRpcBToF[0] = ibooker.book1D(hname, htitle, 100, 0., 250.);
0484 sprintf(hname, "hMuonRpcBToF2");
0485 meMuonRpcBToF[1] = ibooker.book1D(hname, htitle, 100, 0., 50.);
0486 for (Int_t i = 0; i < 2; ++i) {
0487 meMuonRpcBToF[i]->setAxisTitle("Time of Flight of Hits (ns)", 1);
0488 meMuonRpcBToF[i]->setAxisTitle("Count", 2);
0489 }
0490
0491 sprintf(hname, "hMuonRpcBR");
0492 sprintf(htitle, "Muon RPC barrel hits, R/cm");
0493 meMuonRpcBR = ibooker.book1D(hname, htitle, 100, 0., 1500.);
0494 meMuonRpcBR->setAxisTitle("R of Hits (cm)", 1);
0495 meMuonRpcBR->setAxisTitle("Count", 2);
0496 }
0497
0498 GlobalHitsHistogrammer::~GlobalHitsHistogrammer() {}
0499
0500 void GlobalHitsHistogrammer::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0501 std::string MsgLoggerCat = "GlobalHitsHistogrammer_analyze";
0502
0503
0504 ++count;
0505
0506
0507 edm::RunNumber_t nrun = iEvent.id().run();
0508 edm::EventNumber_t nevt = iEvent.id().event();
0509
0510 if (verbosity > 0) {
0511 edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count << " events total)";
0512 } else if (verbosity == 0) {
0513 if (nevt % frequency == 0 || nevt == 1) {
0514 edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count
0515 << " events total)";
0516 }
0517 }
0518
0519
0520 if (getAllProvenances) {
0521 std::vector<const edm::StableProvenance *> AllProv;
0522 iEvent.getAllStableProvenance(AllProv);
0523
0524 if (verbosity >= 0)
0525 edm::LogInfo(MsgLoggerCat) << "Number of Provenances = " << AllProv.size();
0526
0527 if (printProvenanceInfo && (verbosity >= 0)) {
0528 TString eventout("\nProvenance info:\n");
0529
0530 for (unsigned int i = 0; i < AllProv.size(); ++i) {
0531 eventout += "\n ******************************";
0532 eventout += "\n Module : ";
0533 eventout += AllProv[i]->moduleLabel();
0534 eventout += "\n ProductID : ";
0535 eventout += AllProv[i]->productID().id();
0536 eventout += "\n ClassName : ";
0537 eventout += AllProv[i]->className();
0538 eventout += "\n InstanceName : ";
0539 eventout += AllProv[i]->productInstanceName();
0540 eventout += "\n BranchName : ";
0541 eventout += AllProv[i]->branchName();
0542 }
0543 eventout += "\n ******************************\n";
0544 edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0545 printProvenanceInfo = false;
0546 }
0547 getAllProvenances = false;
0548 }
0549
0550
0551 edm::Handle<PGlobalSimHit> srcGlobalHits;
0552 iEvent.getByToken(GlobalHitSrc_Token_, srcGlobalHits);
0553 if (!srcGlobalHits.isValid()) {
0554 edm::LogWarning(MsgLoggerCat) << "Unable to find PGlobalSimHit in event!";
0555 return;
0556 }
0557
0558 nPxlBrlHits = srcGlobalHits->getnPxlBrlHits();
0559 nPxlFwdHits = srcGlobalHits->getnPxlFwdHits();
0560 nPxlHits = nPxlBrlHits + nPxlFwdHits;
0561 nSiBrlHits = srcGlobalHits->getnSiBrlHits();
0562 nSiFwdHits = srcGlobalHits->getnSiFwdHits();
0563 nSiHits = nSiBrlHits + nSiFwdHits;
0564 nMuonDtHits = srcGlobalHits->getnMuonDtHits();
0565 nMuonCscHits = srcGlobalHits->getnMuonCscHits();
0566 nMuonRpcBrlHits = srcGlobalHits->getnMuonRpcBrlHits();
0567 nMuonRpcFwdHits = srcGlobalHits->getnMuonRpcFwdHits();
0568 nMuonHits = nMuonDtHits + nMuonCscHits + nMuonRpcBrlHits + nMuonRpcFwdHits;
0569
0570 for (Int_t i = 0; i < 2; ++i) {
0571 meMCRGP[i]->Fill((float)srcGlobalHits->getnRawGenPart());
0572 meMCG4Vtx[i]->Fill((float)srcGlobalHits->getnG4Vtx());
0573 meMCG4Trk[i]->Fill((float)srcGlobalHits->getnG4Trk());
0574 meCaloEcal[i]->Fill((float)srcGlobalHits->getnECalHits());
0575 meCaloPreSh[i]->Fill((float)srcGlobalHits->getnPreShHits());
0576 meCaloHcal[i]->Fill((float)srcGlobalHits->getnHCalHits());
0577 meTrackerPx[i]->Fill((float)nPxlHits);
0578 meTrackerSi[i]->Fill((float)nSiHits);
0579 meMuon[i]->Fill((float)nMuonHits);
0580 }
0581
0582
0583 std::vector<PGlobalSimHit::Vtx> G4Vtx = srcGlobalHits->getG4Vtx();
0584 for (unsigned int i = 0; i < G4Vtx.size(); ++i) {
0585 for (int j = 0; j < 2; ++j) {
0586 meGeantVtxX[j]->Fill(G4Vtx[i].x);
0587 meGeantVtxY[j]->Fill(G4Vtx[i].y);
0588 meGeantVtxZ[j]->Fill(G4Vtx[i].z);
0589 }
0590 }
0591
0592
0593 std::vector<PGlobalSimHit::Trk> G4Trk = srcGlobalHits->getG4Trk();
0594 for (unsigned int i = 0; i < G4Trk.size(); ++i) {
0595 meGeantTrkPt->Fill(G4Trk[i].pt);
0596 meGeantTrkE->Fill(G4Trk[i].e);
0597 }
0598
0599
0600 std::vector<PGlobalSimHit::CalHit> ECalHits = srcGlobalHits->getECalHits();
0601 for (unsigned int i = 0; i < ECalHits.size(); ++i) {
0602 for (Int_t j = 0; j < 2; ++j) {
0603 meCaloEcalE[j]->Fill(ECalHits[i].e);
0604 meCaloEcalToF[j]->Fill(ECalHits[i].tof);
0605 }
0606 meCaloEcalPhi->Fill(ECalHits[i].phi);
0607 meCaloEcalEta->Fill(ECalHits[i].eta);
0608 }
0609
0610
0611 std::vector<PGlobalSimHit::CalHit> PreShHits = srcGlobalHits->getPreShHits();
0612 for (unsigned int i = 0; i < PreShHits.size(); ++i) {
0613 for (Int_t j = 0; j < 2; ++j) {
0614 meCaloPreShE[j]->Fill(PreShHits[i].e);
0615 meCaloPreShToF[j]->Fill(PreShHits[i].tof);
0616 }
0617 meCaloPreShPhi->Fill(PreShHits[i].phi);
0618 meCaloPreShEta->Fill(PreShHits[i].eta);
0619 }
0620
0621
0622 std::vector<PGlobalSimHit::CalHit> HCalHits = srcGlobalHits->getHCalHits();
0623 for (unsigned int i = 0; i < HCalHits.size(); ++i) {
0624 for (Int_t j = 0; j < 2; ++j) {
0625 meCaloHcalE[j]->Fill(HCalHits[i].e);
0626 meCaloHcalToF[j]->Fill(HCalHits[i].tof);
0627 }
0628 meCaloHcalPhi->Fill(HCalHits[i].phi);
0629 meCaloHcalEta->Fill(HCalHits[i].eta);
0630 }
0631
0632
0633 std::vector<PGlobalSimHit::BrlHit> PxlBrlHits = srcGlobalHits->getPxlBrlHits();
0634 for (unsigned int i = 0; i < PxlBrlHits.size(); ++i) {
0635 meTrackerPxPhi->Fill(PxlBrlHits[i].phi);
0636 meTrackerPxEta->Fill(PxlBrlHits[i].eta);
0637 meTrackerPxBToF->Fill(PxlBrlHits[i].tof);
0638 meTrackerPxBR->Fill(PxlBrlHits[i].r);
0639 }
0640
0641
0642 std::vector<PGlobalSimHit::FwdHit> PxlFwdHits = srcGlobalHits->getPxlFwdHits();
0643 for (unsigned int i = 0; i < PxlFwdHits.size(); ++i) {
0644 meTrackerPxPhi->Fill(PxlFwdHits[i].phi);
0645 meTrackerPxEta->Fill(PxlFwdHits[i].eta);
0646 meTrackerPxFToF->Fill(PxlFwdHits[i].tof);
0647 meTrackerPxFZ->Fill(PxlFwdHits[i].z);
0648 }
0649
0650
0651 std::vector<PGlobalSimHit::BrlHit> SiBrlHits = srcGlobalHits->getSiBrlHits();
0652 for (unsigned int i = 0; i < SiBrlHits.size(); ++i) {
0653 meTrackerSiPhi->Fill(SiBrlHits[i].phi);
0654 meTrackerSiEta->Fill(SiBrlHits[i].eta);
0655 meTrackerSiBToF->Fill(SiBrlHits[i].tof);
0656 meTrackerSiBR->Fill(SiBrlHits[i].r);
0657 }
0658
0659
0660 std::vector<PGlobalSimHit::FwdHit> SiFwdHits = srcGlobalHits->getSiFwdHits();
0661 for (unsigned int i = 0; i < SiFwdHits.size(); ++i) {
0662 meTrackerSiPhi->Fill(SiFwdHits[i].phi);
0663 meTrackerSiEta->Fill(SiFwdHits[i].eta);
0664 meTrackerSiFToF->Fill(SiFwdHits[i].tof);
0665 meTrackerSiFZ->Fill(SiFwdHits[i].z);
0666 }
0667
0668
0669 std::vector<PGlobalSimHit::FwdHit> MuonCscHits = srcGlobalHits->getMuonCscHits();
0670 for (unsigned int i = 0; i < MuonCscHits.size(); ++i) {
0671 meMuonPhi->Fill(MuonCscHits[i].phi);
0672 meMuonEta->Fill(MuonCscHits[i].eta);
0673 for (Int_t j = 0; j < 2; ++j) {
0674 meMuonCscToF[j]->Fill(MuonCscHits[i].tof);
0675 }
0676 meMuonCscZ->Fill(MuonCscHits[i].z);
0677 }
0678
0679
0680 std::vector<PGlobalSimHit::BrlHit> MuonDtHits = srcGlobalHits->getMuonDtHits();
0681 for (unsigned int i = 0; i < MuonDtHits.size(); ++i) {
0682 meMuonPhi->Fill(MuonDtHits[i].phi);
0683 meMuonEta->Fill(MuonDtHits[i].eta);
0684 for (Int_t j = 0; j < 2; ++j) {
0685 meMuonDtToF[j]->Fill(MuonDtHits[i].tof);
0686 }
0687 meMuonDtR->Fill(MuonDtHits[i].r);
0688 }
0689
0690
0691 std::vector<PGlobalSimHit::FwdHit> MuonRpcFwdHits = srcGlobalHits->getMuonRpcFwdHits();
0692 for (unsigned int i = 0; i < MuonRpcFwdHits.size(); ++i) {
0693 meMuonPhi->Fill(MuonRpcFwdHits[i].phi);
0694 meMuonEta->Fill(MuonRpcFwdHits[i].eta);
0695 for (Int_t j = 0; j < 2; ++j) {
0696 meMuonRpcFToF[j]->Fill(MuonRpcFwdHits[i].tof);
0697 }
0698 meMuonRpcFZ->Fill(MuonRpcFwdHits[i].z);
0699 }
0700
0701
0702 std::vector<PGlobalSimHit::BrlHit> MuonRpcBrlHits = srcGlobalHits->getMuonRpcBrlHits();
0703 for (unsigned int i = 0; i < MuonRpcBrlHits.size(); ++i) {
0704 meMuonPhi->Fill(MuonRpcBrlHits[i].phi);
0705 meMuonEta->Fill(MuonRpcBrlHits[i].eta);
0706 for (Int_t j = 0; j < 2; ++j) {
0707 meMuonRpcBToF[j]->Fill(MuonRpcBrlHits[i].tof);
0708 }
0709 meMuonRpcBR->Fill(MuonRpcBrlHits[i].r);
0710 }
0711
0712 return;
0713 }