File indexing completed on 2022-04-04 22:01:45
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021 #include <vector>
0022 #include <string>
0023
0024 #include "FWCore/Framework/interface/Frameworkfwd.h"
0025 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0026
0027 #include "FWCore/Framework/interface/Event.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029
0030 #include "FWCore/Utilities/interface/InputTag.h"
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032
0033 #include "SimDataFormats/CaloHit/interface/HFShowerPhoton.h"
0034 #include "SimDataFormats/CaloHit/interface/HFShowerLibraryEventInfo.h"
0035
0036 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0037 #include "FWCore/ServiceRegistry/interface/Service.h"
0038 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0039 #include "TFile.h"
0040 #include "TTree.h"
0041 #include "TBranch.h"
0042 #include "TH1F.h"
0043 #include "TH1I.h"
0044
0045 class AnalyzeTuples : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0046 public:
0047 explicit AnalyzeTuples(const edm::ParameterSet&);
0048 ~AnalyzeTuples() override = default;
0049
0050 private:
0051 void beginJob() override;
0052 void analyze(const edm::Event&, const edm::EventSetup&) override;
0053 void endJob() override;
0054 void loadEventInfo(TBranch*);
0055 void getRecord(int type, int record);
0056 TFile* hf;
0057 TBranch *emBranch, *hadBranch;
0058
0059 int nMomBin, totEvents, evtPerBin;
0060 float libVers, listVersion;
0061 std::vector<double> pmom;
0062 std::vector<HFShowerPhoton> photon;
0063
0064 TH1I* hNPELongElec[12];
0065 TH1I* hNPEShortElec[12];
0066 TH1I* hNPELongPion[12];
0067 TH1I* hNPEShortPion[12];
0068 };
0069
0070 AnalyzeTuples::AnalyzeTuples(const edm::ParameterSet& iConfig) {
0071 usesResource(TFileService::kSharedResource);
0072
0073 edm::LogVerbatim("HcalSim") << "analyzetuples a buraya girdi";
0074 edm::ParameterSet m_HS = iConfig.getParameter<edm::ParameterSet>("HFShowerLibrary");
0075 edm::FileInPath fp = m_HS.getParameter<edm::FileInPath>("FileName");
0076 std::string pTreeName = fp.fullPath();
0077 std::string emName = m_HS.getParameter<std::string>("TreeEMID");
0078 std::string hadName = m_HS.getParameter<std::string>("TreeHadID");
0079 std::string branchEvInfo = m_HS.getUntrackedParameter<std::string>(
0080 "BranchEvt", "HFShowerLibraryEventInfos_hfshowerlib_HFShowerLibraryEventInfo");
0081 std::string branchPre = m_HS.getUntrackedParameter<std::string>("BranchPre", "HFShowerPhotons_hfshowerlib_");
0082 std::string branchPost = m_HS.getUntrackedParameter<std::string>("BranchPost", "_R.obj");
0083
0084 if (pTreeName.find(".") == 0)
0085 pTreeName.erase(0, 2);
0086 const char* nTree = pTreeName.c_str();
0087 hf = TFile::Open(nTree);
0088
0089 if (!hf->IsOpen()) {
0090 edm::LogError("HFShower") << "HFShowerLibrary: opening " << nTree << " failed";
0091 throw cms::Exception("Unknown", "HFShowerLibrary") << "Opening of " << pTreeName << " fails\n";
0092 } else {
0093 edm::LogInfo("HFShower") << "HFShowerLibrary: opening " << nTree << " successfully";
0094 }
0095
0096 TTree* event = (TTree*)hf->Get("Events");
0097 if (event) {
0098 std::string info = branchEvInfo + branchPost;
0099 TBranch* evtInfo = event->GetBranch(info.c_str());
0100 if (evtInfo) {
0101 loadEventInfo(evtInfo);
0102 } else {
0103 edm::LogError("HFShower") << "HFShowerLibrary: HFShowerLibrayEventInfo"
0104 << " Branch does not exist in Event";
0105 throw cms::Exception("Unknown", "HFShowerLibrary") << "Event information absent\n";
0106 }
0107 } else {
0108 edm::LogError("HFShower") << "HFShowerLibrary: Events Tree does not "
0109 << "exist";
0110 throw cms::Exception("Unknown", "HFShowerLibrary") << "Events tree absent\n";
0111 }
0112
0113 edm::LogInfo("HFShower") << "HFShowerLibrary: Library " << libVers << " ListVersion " << listVersion
0114 << " Events Total " << totEvents << " and " << evtPerBin << " per bin";
0115 edm::LogInfo("HFShower") << "HFShowerLibrary: Energies (GeV) with " << nMomBin << " bins";
0116 for (int i = 0; i < nMomBin; i++)
0117 edm::LogInfo("HFShower") << "HFShowerLibrary: pmom[" << i << "] = " << pmom[i] << " GeV";
0118
0119 std::string nameBr = branchPre + emName + branchPost;
0120 emBranch = event->GetBranch(nameBr.c_str());
0121 nameBr = branchPre + hadName + branchPost;
0122 hadBranch = event->GetBranch(nameBr.c_str());
0123 edm::LogInfo("HFShower") << "HFShowerLibrary:Branch " << emName << " has " << emBranch->GetEntries()
0124 << " entries and Branch " << hadName << " has " << hadBranch->GetEntries() << " entries";
0125 edm::LogInfo("HFShower") << "HFShowerLibrary::No packing information -"
0126 << " Assume x, y, z are not in packed form";
0127 }
0128
0129 void AnalyzeTuples::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0130 for (int ibin = 0; ibin < 12; ++ibin) {
0131 int min = evtPerBin * (ibin);
0132 int max = evtPerBin * (ibin + 1);
0133 for (int i = min; i < max; ++i) {
0134 getRecord(0, i);
0135 int npe_long = 0;
0136 int npe_short = 0;
0137 edm::LogVerbatim("HcalSim") << "phptons size" << photon.size();
0138 for (int j = 0; j < int(photon.size()); ++j) {
0139
0140 if (photon[j].z() < 0) {
0141
0142 ++npe_short;
0143 } else {
0144
0145 ++npe_long;
0146 edm::LogVerbatim("HcalSim") << photon[j].z();
0147 }
0148 }
0149 hNPELongElec[ibin]->Fill(npe_long);
0150 edm::LogVerbatim("HcalSim") << ibin << npe_long;
0151 hNPEShortElec[ibin]->Fill(npe_short);
0152 }
0153 }
0154 for (int ibin = 0; ibin < 12; ++ibin) {
0155 int min = evtPerBin * (ibin);
0156 int max = evtPerBin * (ibin + 1);
0157 for (int i = min; i < max; ++i) {
0158 getRecord(1, i);
0159 int npe_long = 0;
0160 int npe_short = 0;
0161 for (int j = 0; j < int(photon.size()); ++j) {
0162
0163 if (photon[j].z() < 0) {
0164
0165 ++npe_short;
0166 } else {
0167
0168 ++npe_long;
0169 }
0170 }
0171 hNPELongPion[ibin]->Fill(npe_long);
0172 hNPEShortPion[ibin]->Fill(npe_short);
0173 }
0174 }
0175 }
0176
0177 void AnalyzeTuples::beginJob() {
0178 edm::Service<TFileService> fs;
0179 TFileDirectory HFDir = fs->mkdir("HF");
0180 char title[128];
0181 for (int i = 0; i < int(pmom.size()); ++i) {
0182 sprintf(title, "NPELongElec_Mom_%i", int(pmom[i]));
0183 int maxBin = int(pmom[i] + 50);
0184 hNPELongElec[i] = HFDir.make<TH1I>(title, "NPE Long", 140, 0, 140);
0185 sprintf(title, "NPEShortElec_Mom_%i", int(pmom[i]));
0186 hNPEShortElec[i] = HFDir.make<TH1I>(title, "NPE Short", maxBin, 0, maxBin);
0187 sprintf(title, "NPELongPion_Mom_%i", int(pmom[i]));
0188 hNPELongPion[i] = HFDir.make<TH1I>(title, "NPE Long", maxBin, 0, maxBin);
0189 sprintf(title, "NPEShortPion_Mom_%i", int(pmom[i]));
0190 hNPEShortPion[i] = HFDir.make<TH1I>(title, "NPE Short", maxBin, 0, maxBin);
0191 }
0192 }
0193
0194 void AnalyzeTuples::endJob() {}
0195
0196 void AnalyzeTuples::loadEventInfo(TBranch* branch) {
0197 std::vector<HFShowerLibraryEventInfo> eventInfoCollection;
0198 branch->SetAddress(&eventInfoCollection);
0199 branch->GetEntry(0);
0200 edm::LogInfo("HFShower") << "HFShowerLibrary::loadEventInfo loads "
0201 << " EventInfo Collection of size " << eventInfoCollection.size() << " records";
0202 totEvents = eventInfoCollection[0].totalEvents();
0203 nMomBin = eventInfoCollection[0].numberOfBins();
0204 evtPerBin = eventInfoCollection[0].eventsPerBin();
0205 libVers = eventInfoCollection[0].showerLibraryVersion();
0206 listVersion = eventInfoCollection[0].physListVersion();
0207 pmom = eventInfoCollection[0].energyBins();
0208 }
0209
0210 void AnalyzeTuples::getRecord(int type, int record) {
0211 int nrc = record - 1;
0212 photon.clear();
0213 if (type > 0) {
0214 hadBranch->SetAddress(&photon);
0215 hadBranch->GetEntry(nrc);
0216 } else {
0217 emBranch->SetAddress(&photon);
0218 emBranch->GetEntry(nrc);
0219 }
0220 #ifdef EDM_ML_DEBUG
0221 int nPhoton = photon.size();
0222 edm::LogVerbatim("HFShower") << "HFShowerLibrary::getRecord: Record " << record << " of type " << type << " with "
0223 << nPhoton << " photons";
0224 for (int j = 0; j < nPhoton; j++)
0225 edm::LogVerbatim("HFShower") << "Photon " << j << " " << photon[j];
0226 #endif
0227 }
0228
0229
0230 DEFINE_FWK_MODULE(AnalyzeTuples);