File indexing completed on 2024-04-06 12:14:49
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019 #include <memory>
0020 #include <iostream>
0021 #include <fstream>
0022
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/EventSetup.h"
0029 #include "FWCore/Framework/interface/MakerMacros.h"
0030 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0033 #include "DataFormats/HcalDetId/interface/HcalTestNumbering.h"
0034 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0035 #include "Geometry/Records/interface/HcalRecNumberingRecord.h"
0036 #include "Geometry/HcalCommonData/interface/HcalDDDRecConstants.h"
0037 #include "Geometry/HcalCommonData/interface/HcalHitRelabeller.h"
0038
0039 class HcalHitRelabellerTester : public edm::one::EDAnalyzer<> {
0040 public:
0041 explicit HcalHitRelabellerTester(const edm::ParameterSet&);
0042 ~HcalHitRelabellerTester() override = default;
0043 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0044
0045 void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0046
0047 private:
0048 bool nd_;
0049 edm::ESGetToken<HcalDDDRecConstants, HcalRecNumberingRecord> token_;
0050 };
0051
0052 HcalHitRelabellerTester::HcalHitRelabellerTester(const edm::ParameterSet& ps)
0053 : nd_(ps.getParameter<bool>("neutralDensity")),
0054 token_{esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord>(edm::ESInputTag{})} {
0055 edm::LogVerbatim("HCalGeom") << "Construct HcalHitRelabellerTester with Neutral Density: " << nd_;
0056 }
0057
0058 void HcalHitRelabellerTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0059 edm::ParameterSetDescription desc;
0060 desc.add<bool>("neutralDensity", true);
0061 descriptions.add("hcalHitRelabellerTester", desc);
0062 }
0063
0064
0065 void HcalHitRelabellerTester::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0066 const HcalDDDRecConstants* theRecNumber = &iSetup.getData(token_);
0067 std::unique_ptr<HcalHitRelabeller> relabel = std::make_unique<HcalHitRelabeller>(nd_);
0068 relabel->setGeometry(theRecNumber);
0069 std::vector<int> etas = {-29, -26, -22, -19, -16, -13, -10, -7, -4, -1, 2, 5, 8, 11, 14, 17, 20, 23, 26, 29};
0070 std::vector<int> layers = {1, 2};
0071 const int iphi = 63;
0072 edm::LogVerbatim("HCalGeom") << "HcalHitRelabellerTester:: Testing " << etas.size() << " eta, " << layers.size()
0073 << " layer values and iphi = " << iphi;
0074 for (const auto& eta : etas) {
0075 int ieta = std::abs(eta);
0076 int det = (ieta <= 16) ? 1 : 2;
0077 int zside = (eta >= 0) ? 1 : -1;
0078 for (const auto& lay : layers) {
0079 if (ieta == 16)
0080 det = (lay <= 3) ? 1 : 2;
0081 int depth = theRecNumber->findDepth(det, ieta, iphi, zside, lay);
0082 if (depth > 0) {
0083 uint32_t id = HcalTestNumbering::packHcalIndex(det, zside, depth, ieta, iphi, lay);
0084 double wt = relabel->energyWt(id);
0085 edm::LogVerbatim("HCalGeom") << "Det " << det << " Z " << zside << " Eta" << ieta << " Phi " << iphi << " Lay "
0086 << lay << " Depth " << depth << " Layer0Wt " << wt;
0087 }
0088 }
0089 }
0090 }
0091
0092
0093 DEFINE_FWK_MODULE(HcalHitRelabellerTester);