File indexing completed on 2024-04-06 12:32:40
0001
0002 #include <iostream>
0003 #include <fstream>
0004 #include <vector>
0005 #include <map>
0006 #include <string>
0007
0008
0009 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0010 #include "DataFormats/HGCDigi/interface/HGCDigiCollections.h"
0011
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0019 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0020 #include "FWCore/ServiceRegistry/interface/Service.h"
0021
0022 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0023 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0024
0025 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0026 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0027 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0028
0029
0030 #include "TH1.h"
0031 #include "TH2.h"
0032
0033 class HGCalSiliconValidation : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0034 public:
0035 explicit HGCalSiliconValidation(const edm::ParameterSet& ps);
0036 ~HGCalSiliconValidation() override = default;
0037
0038 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0039
0040 protected:
0041 void beginJob() override {}
0042 void beginRun(edm::Run const&, edm::EventSetup const&) override;
0043 void endRun(edm::Run const&, edm::EventSetup const&) override {}
0044 void analyze(edm::Event const&, edm::EventSetup const&) override;
0045
0046 private:
0047 const std::string g4Label_, nameDetector_, hgcalHits_;
0048 const edm::InputTag hgcalDigis_;
0049 const int iSample_;
0050 const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> tok_hgcalgeom_;
0051 const edm::EDGetTokenT<edm::PCaloHitContainer> tok_hits_;
0052 const edm::EDGetTokenT<HGCalDigiCollection> tok_digi_;
0053
0054 TH1D *hsimE1_, *hsimE2_, *hsimTm_;
0055 TH1D *hsimLn_, *hdigEn_, *hdigLn_;
0056 TH2D *hsimOc_, *hsi2Oc_, *hdigOc_, *hdi2Oc_;
0057 };
0058
0059 HGCalSiliconValidation::HGCalSiliconValidation(const edm::ParameterSet& ps)
0060 : g4Label_(ps.getUntrackedParameter<std::string>("ModuleLabel", "g4SimHits")),
0061 nameDetector_(ps.getUntrackedParameter<std::string>("detectorName", "HGCalEESensitive")),
0062 hgcalHits_((ps.getUntrackedParameter<std::string>("HitCollection", "HGCHitsEE"))),
0063 hgcalDigis_(ps.getUntrackedParameter<edm::InputTag>("DigiCollection")),
0064 iSample_(ps.getUntrackedParameter<int>("Sample", 5)),
0065 tok_hgcalgeom_(esConsumes<HGCalGeometry, IdealGeometryRecord>(edm::ESInputTag{"", nameDetector_})),
0066 tok_hits_(consumes<edm::PCaloHitContainer>(edm::InputTag(g4Label_, hgcalHits_))),
0067 tok_digi_(consumes<HGCalDigiCollection>(hgcalDigis_)) {
0068 usesResource(TFileService::kSharedResource);
0069
0070 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation::Input for SimHit:"
0071 << edm::InputTag(g4Label_, hgcalHits_) << " Digits:" << hgcalDigis_
0072 << " Sample: " << iSample_;
0073 }
0074
0075 void HGCalSiliconValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0076 edm::ParameterSetDescription desc;
0077 desc.addUntracked<std::string>("ModuleLabel", "g4SimHits");
0078 desc.addUntracked<std::string>("detectorName", "HGCalEESensitive");
0079 desc.addUntracked<std::string>("HitCollection", "HGCHitsEE");
0080 desc.addUntracked<edm::InputTag>("DigiCollection", edm::InputTag("simHGCalUnsuppressedDigis", "EE"));
0081 desc.addUntracked<int>("Sample", 5);
0082 descriptions.add("hgcalSiliconAnalysisEE", desc);
0083 }
0084
0085 void HGCalSiliconValidation::beginRun(edm::Run const&, edm::EventSetup const& es) {
0086 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation::Booking the Histograms";
0087
0088
0089 edm::Service<TFileService> fs;
0090 hsimE1_ = fs->make<TH1D>("SimHitEn1", "Sim Hit Energy", 1000, 0.0, 1.0);
0091 hsimE2_ = fs->make<TH1D>("SimHitEn2", "Sim Hit Energy", 1000, 0.0, 1.0);
0092 hsimTm_ = fs->make<TH1D>("SimHitTime", "Sim Hit Time", 1000, 0.0, 500.0);
0093 hsimLn_ = fs->make<TH1D>("SimHitLong", "Sim Hit Long. Profile", 60, 0.0, 30.0);
0094 hsimOc_ = fs->make<TH2D>("SimHitOccup", "Sim Hit Occupnacy", 300, 0.0, 300.0, 60, 0.0, 30.0);
0095 hsi2Oc_ = fs->make<TH2D>("SimHitOccu2", "Sim Hit Occupnacy", 300, 300.0, 600.0, 300, 0.0, 300.0);
0096
0097 hdigEn_ = fs->make<TH1D>("DigiEnergy", "Digi ADC Sample", 1000, 0.0, 1000.0);
0098 hdigLn_ = fs->make<TH1D>("DigiLong", "Digi Long. Profile", 60, 0.0, 30.0);
0099 hdigOc_ = fs->make<TH2D>("DigiOccup", "Digi Occupnacy", 300, 0.0, 300.0, 60, 0.0, 30.0);
0100 hdi2Oc_ = fs->make<TH2D>("DigiOccu2", "Digi Occupnacy", 300, 300.0, 600.0, 300, 0.0, 300.0);
0101 }
0102
0103 void HGCalSiliconValidation::analyze(const edm::Event& e, const edm::EventSetup& iSetup) {
0104 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation:Run = " << e.id().run()
0105 << " Event = " << e.id().event();
0106
0107 const edm::ESHandle<HGCalGeometry>& geom = iSetup.getHandle(tok_hgcalgeom_);
0108 if (!geom.isValid()) {
0109 edm::LogWarning("HGCalValidation") << "Cannot get valid HGCalGeometry Object for " << nameDetector_;
0110 } else {
0111 const HGCalGeometry* geom0 = geom.product();
0112
0113
0114 const edm::Handle<edm::PCaloHitContainer>& hitsCalo = e.getHandle(tok_hits_);
0115 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation.: PCaloHitContainer obtained with flag "
0116 << hitsCalo.isValid();
0117 if (hitsCalo.isValid()) {
0118 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation: PCaloHit buffer " << hitsCalo->size();
0119 unsigned i(0);
0120 std::map<unsigned int, double> map_try;
0121 for (edm::PCaloHitContainer::const_iterator it = hitsCalo->begin(); it != hitsCalo->end(); ++it) {
0122 double energy = it->energy();
0123 double time = it->time();
0124 unsigned int id = it->id();
0125 GlobalPoint pos = geom0->getPosition(DetId(id));
0126 double r = pos.perp();
0127 double z = std::abs(pos.z());
0128 int lay = HGCSiliconDetId(id).layer();
0129 hsimE1_->Fill(energy);
0130 hsimTm_->Fill(time, energy);
0131 hsimOc_->Fill(r, lay, energy);
0132 hsi2Oc_->Fill(z, r, energy);
0133 hsimLn_->Fill(lay, energy);
0134 double ensum = (map_try.count(id) != 0) ? map_try[id] : 0;
0135 ensum += energy;
0136 map_try[id] = ensum;
0137 ++i;
0138 edm::LogVerbatim("HGCalValidation") << "HGCalHit[" << i << "] ID " << std::hex << " " << id << std::dec << " "
0139 << HGCSiliconDetId(id) << " E " << energy << " time " << time;
0140 }
0141 for (std::map<unsigned int, double>::iterator itr = map_try.begin(); itr != map_try.end(); ++itr) {
0142 hsimE2_->Fill((*itr).second);
0143 }
0144 }
0145
0146
0147 unsigned int kount(0);
0148 const edm::Handle<HGCalDigiCollection>& digicoll = e.getHandle(tok_digi_);
0149 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation.: HGCalDigiCollection obtained with flag "
0150 << digicoll.isValid();
0151 if (digicoll.isValid()) {
0152 edm::LogVerbatim("HGCalValidation") << "HGCalSiliconValidation: HGCalDigi buffer " << digicoll->size();
0153 for (HGCalDigiCollection::const_iterator it = digicoll->begin(); it != digicoll->end(); ++it) {
0154 HGCalDataFrame df(*it);
0155 double energy = df[iSample_].data();
0156 HGCSiliconDetId cell(df.id());
0157 GlobalPoint pos = geom0->getPosition(cell);
0158 double r = pos.perp();
0159 double z = std::abs(pos.z());
0160 int depth = cell.layer();
0161 hdigEn_->Fill(energy);
0162 hdigLn_->Fill(depth);
0163 hdigOc_->Fill(r, depth);
0164 hdi2Oc_->Fill(z, r);
0165 ++kount;
0166 edm::LogVerbatim("HGCalValidation") << "HGCalDigit[" << kount << "] ID " << cell << " E " << energy;
0167 }
0168 }
0169 }
0170 }
0171
0172 DEFINE_FWK_MODULE(HGCalSiliconValidation);
