File indexing completed on 2024-04-06 12:30:08
0001 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0002 #include "FWCore/Framework/interface/Frameworkfwd.h"
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009 #include "FWCore/ServiceRegistry/interface/Service.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012
0013 #include "SimDataFormats/ValidationFormats/interface/MaterialAccountingCalo.h"
0014
0015 #include <TH1F.h>
0016
0017 #include <iostream>
0018 #include <string>
0019 #include <vector>
0020
0021
0022
0023 class HGCalTBMBAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0024 public:
0025 HGCalTBMBAnalyzer(const edm::ParameterSet &);
0026 ~HGCalTBMBAnalyzer() override = default;
0027
0028 static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0029
0030 private:
0031 void analyze(edm::Event const &, edm::EventSetup const &) override;
0032
0033 const std::vector<std::string> listNames_;
0034 const edm::InputTag labelMBCalo_;
0035 const unsigned int nList_;
0036 const edm::EDGetTokenT<MaterialAccountingCaloCollection> tokMBCalo_;
0037 std::vector<TH1D *> me100_, me200_, me300_;
0038 };
0039
0040 HGCalTBMBAnalyzer::HGCalTBMBAnalyzer(const edm::ParameterSet &p)
0041 : listNames_(p.getParameter<std::vector<std::string>>("detectorNames")),
0042 labelMBCalo_(p.getParameter<edm::InputTag>("labelMBCalo")),
0043 nList_(listNames_.size()),
0044 tokMBCalo_(consumes<MaterialAccountingCaloCollection>(labelMBCalo_)) {
0045 edm::LogVerbatim("HGCSim") << "HGCalTBMBAnalyzer initialized for " << nList_ << " volumes";
0046 for (unsigned int k = 0; k < nList_; ++k)
0047 edm::LogVerbatim("HGCSim") << " [" << k << "] " << listNames_[k];
0048
0049 edm::Service<TFileService> tfile;
0050 if (!tfile.isAvailable())
0051 throw cms::Exception("BadConfig") << "TFileService unavailable: "
0052 << "please add it to config file";
0053 char name[20], title[80];
0054 TH1D *hist;
0055 for (unsigned int i = 0; i <= nList_; i++) {
0056 std::string named = (i == nList_) ? "Total" : listNames_[i];
0057 sprintf(name, "RadL%d", i);
0058 sprintf(title, "MB(X0) for (%s)", named.c_str());
0059 hist = tfile->make<TH1D>(name, title, 100000, 0.0, 100.0);
0060 hist->Sumw2(true);
0061 me100_.push_back(hist);
0062 sprintf(name, "IntL%d", i);
0063 sprintf(title, "MB(L0) for (%s)", named.c_str());
0064 hist = tfile->make<TH1D>(name, title, 100000, 0.0, 10.0);
0065 hist->Sumw2(true);
0066 me200_.push_back(hist);
0067 sprintf(name, "StepL%d", i);
0068 sprintf(title, "MB(Step) for (%s)", named.c_str());
0069 hist = tfile->make<TH1D>(name, title, 100000, 0.0, 50000.0);
0070 hist->Sumw2(true);
0071 me300_.push_back(hist);
0072 }
0073 edm::LogVerbatim("HGCSim") << "HGCalTBMBAnalyzer: Booking user histos done ===";
0074 }
0075
0076 void HGCalTBMBAnalyzer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0077 edm::ParameterSetDescription desc;
0078 std::vector<std::string> names = {"HGCalBeamWChamb",
0079 "HGCalBeamS1",
0080 "HGCalBeamS2",
0081 "HGCalBeamS3",
0082 "HGCalBeamS4",
0083 "HGCalBeamS5",
0084 "HGCalBeamS6",
0085 "HGCalBeamCK3",
0086 "HGCalBeamHaloCounter",
0087 "HGCalBeamMuonCounter",
0088 "HGCalEE",
0089 "HGCalHE",
0090 "HGCalAH"};
0091 desc.add<std::vector<std::string>>("detectorNames", names);
0092 desc.add<edm::InputTag>("labelMBCalo", edm::InputTag("g4SimHits", "HGCalTBMB"));
0093 descriptions.add("hgcalTBMBAnalyzer", desc);
0094 }
0095
0096 void HGCalTBMBAnalyzer::analyze(edm::Event const &iEvent, edm::EventSetup const &iSetup) {
0097 #ifdef EDM_ML_DEBUG
0098 edm::LogVerbatim("HGCSim") << "Run " << iEvent.id().run() << " Event " << iEvent.id().event() << " Luminosity "
0099 << iEvent.luminosityBlock() << " Bunch " << iEvent.bunchCrossing();
0100 #endif
0101
0102
0103 auto const &hgcalMBColl = iEvent.getHandle(tokMBCalo_);
0104 if (hgcalMBColl.isValid()) {
0105 auto hgcalMB = hgcalMBColl.product();
0106 #ifdef EDM_ML_DEBUG
0107 edm::LogVerbatim("HGCSim") << "Finds MaterialBudegetCollection with " << hgcalMB->size() << " entries";
0108 #endif
0109
0110 for (auto itr = hgcalMB->begin(); itr != hgcalMB->end(); ++itr) {
0111 for (uint32_t ii = 0; ii < itr->m_stepLen.size(); ++ii) {
0112 #ifdef EDM_ML_DEBUG
0113 edm::LogVerbatim("HGCSim") << "HGCalTBMBAnalyzer:index " << ii << " integrated step " << itr->m_stepLen[ii]
0114 << " X0 " << itr->m_radLen[ii] << " Lamda " << itr->m_intLen[ii];
0115 #endif
0116 if (ii < nList_) {
0117 me100_[ii]->Fill(itr->m_radLen[ii]);
0118 me200_[ii]->Fill(itr->m_intLen[ii]);
0119 me300_[ii]->Fill(itr->m_stepLen[ii]);
0120 }
0121 }
0122 }
0123 }
0124 }
0125
0126 DEFINE_FWK_MODULE(HGCalTBMBAnalyzer);