File indexing completed on 2024-04-06 12:29:50
0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0003
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "FWCore/Framework/interface/MakerMacros.h"
0007
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013
0014 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0015 #include "DataFormats/HcalDetId/interface/HcalTestNumbering.h"
0016 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0017 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0018 #include "SimG4CMS/Calo/interface/CaloHitID.h"
0019 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0020 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0021
0022 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0023
0024 #include <memory>
0025 #include <iostream>
0026 #include <fstream>
0027 #include <vector>
0028 #include <map>
0029 #include <string>
0030
0031 class HcalSimHitDump : public edm::one::EDAnalyzer<> {
0032 public:
0033 HcalSimHitDump(const edm::ParameterSet& ps);
0034 ~HcalSimHitDump() override = default;
0035 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0036
0037 protected:
0038 void beginJob() override {}
0039 void endJob() override {}
0040 void analyze(const edm::Event& e, const edm::EventSetup& c) override;
0041
0042 void analyzeHits(std::vector<PCaloHit>&);
0043
0044 private:
0045 const std::string g4Label_, hitLab_;
0046 const int maxEvent_;
0047 const bool testNumber_;
0048 const edm::EDGetTokenT<edm::PCaloHitContainer> toks_calo_;
0049 int nevt_;
0050 };
0051
0052 HcalSimHitDump::HcalSimHitDump(const edm::ParameterSet& ps)
0053 : g4Label_(ps.getParameter<std::string>("ModuleLabel")),
0054 hitLab_(ps.getParameter<std::string>("HCCollection")),
0055 maxEvent_(ps.getParameter<int>("MaxEvent")),
0056 testNumber_(ps.getParameter<bool>("TestNumber")),
0057 toks_calo_(consumes<edm::PCaloHitContainer>(edm::InputTag(g4Label_, hitLab_))),
0058 nevt_(0) {
0059 edm::LogVerbatim("HitStudy") << "HcalSimHitDump::Module Label: " << g4Label_ << " Hits: " << hitLab_ << " MaxEvent "
0060 << maxEvent_ << " TestNumbering " << testNumber_;
0061 }
0062
0063 void HcalSimHitDump::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0064 edm::ParameterSetDescription desc;
0065 desc.add<std::string>("ModuleLabel", "g4SimHits");
0066 desc.add<std::string>("HCCollection", "HcalHits");
0067 desc.add<int>("MaxEvent", 10);
0068 desc.add<bool>("TestNumber", true);
0069 descriptions.add("hcalSimHitDump", desc);
0070 }
0071
0072 void HcalSimHitDump::analyze(const edm::Event& e, const edm::EventSetup&) {
0073 ++nevt_;
0074 edm::LogVerbatim("HitStudy") << "HcalSimHitDump::Serial # " << nevt_ << " Run # " << e.id().run() << " Event # "
0075 << e.id().event();
0076
0077 if (nevt_ <= maxEvent_) {
0078 std::vector<PCaloHit> hcHits;
0079 const edm::Handle<edm::PCaloHitContainer>& hitsCalo = e.getHandle(toks_calo_);
0080 if (hitsCalo.isValid()) {
0081 edm::LogVerbatim("HitStudy") << "HcalValidation: get valid hist for Hcal";
0082 std::vector<PCaloHit> caloHits;
0083 caloHits.insert(caloHits.end(), hitsCalo->begin(), hitsCalo->end());
0084 edm::LogVerbatim("HitStudy") << "HcalValidation: Hit buffer " << caloHits.size();
0085 analyzeHits(caloHits);
0086 }
0087 }
0088 }
0089
0090 void HcalSimHitDump::analyzeHits(std::vector<PCaloHit>& hits) {
0091
0092 for (unsigned int i = 0; i < hits.size(); i++) {
0093 double edep = hits[i].energy();
0094 double time = hits[i].time();
0095 unsigned int id_ = hits[i].id();
0096 if (testNumber_) {
0097 int det, z, depth, eta, phi, lay;
0098 HcalTestNumbering::unpackHcalIndex(id_, det, z, depth, eta, phi, lay);
0099 std::string sub("HX");
0100 if (det == 1)
0101 sub = "HB";
0102 else if (det == 2)
0103 sub = "HE";
0104 else if (det == 3)
0105 sub = "HO";
0106 else if (det == 4)
0107 sub = "HF";
0108 else if (det == 5)
0109 sub = "HT";
0110 int side = (z == 0) ? (-1) : (1);
0111 edm::LogVerbatim("HitStudy") << "[" << i << "] (" << sub << " " << side * eta << "," << phi << "," << depth << ","
0112 << lay << ") E " << edep << " T " << time;
0113 } else {
0114 edm::LogVerbatim("HitStudy") << "[" << i << "] " << HcalDetId(id_) << " E " << edep << " T " << time;
0115 }
0116 }
0117 }
0118
0119
0120 DEFINE_FWK_MODULE(HcalSimHitDump);