File indexing completed on 2023-03-17 13:03:41
0001 #include <cmath>
0002 #include <iostream>
0003 #include <string>
0004 #include <vector>
0005
0006 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0007
0008 #include "FWCore/Framework/interface/Frameworkfwd.h"
0009 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015
0016 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0017 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0018 #include "DataFormats/Math/interface/angle_units.h"
0019 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0020
0021 using namespace angle_units::operators;
0022
0023 class HGCalGeometryRotTest : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
0024 public:
0025 explicit HGCalGeometryRotTest(const edm::ParameterSet&);
0026 ~HGCalGeometryRotTest() override = default;
0027 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0028
0029 void beginRun(edm::Run const&, edm::EventSetup const&) override;
0030 void analyze(edm::Event const& iEvent, edm::EventSetup const&) override {}
0031 void endRun(edm::Run const&, edm::EventSetup const&) override {}
0032
0033 private:
0034 const std::string nameDetector_;
0035 const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> geomToken_;
0036 const std::vector<int> layers_;
0037 const std::vector<int> waferU_, waferV_, cellU_, cellV_, types_;
0038 };
0039
0040 HGCalGeometryRotTest::HGCalGeometryRotTest(const edm::ParameterSet& iC)
0041 : nameDetector_(iC.getParameter<std::string>("detectorName")),
0042 geomToken_(esConsumes<HGCalGeometry, IdealGeometryRecord, edm::Transition::BeginRun>(
0043 edm::ESInputTag{"", nameDetector_})),
0044 layers_(iC.getParameter<std::vector<int>>("layers")),
0045 waferU_(iC.getParameter<std::vector<int>>("waferUs")),
0046 waferV_(iC.getParameter<std::vector<int>>("waferVs")),
0047 cellU_(iC.getParameter<std::vector<int>>("cellUs")),
0048 cellV_(iC.getParameter<std::vector<int>>("cellVs")),
0049 types_(iC.getParameter<std::vector<int>>("types")) {}
0050
0051 void HGCalGeometryRotTest::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0052 edm::ParameterSetDescription desc;
0053 std::vector<int> layer = {27, 28, 29, 30, 31, 32};
0054 std::vector<int> waferU = {-2, -3, 1, -8, 2, 8};
0055 std::vector<int> waferV = {0, -2, 3, 0, 9, 0};
0056 std::vector<int> cellU = {16, 4, 8, 11, 11, 5};
0057 std::vector<int> cellV = {20, 10, 17, 13, 9, 2};
0058 std::vector<int> type = {0, 0, 0, 2, 2, 2};
0059 desc.add<std::string>("detectorName", "HGCalHESiliconSensitive");
0060 desc.add<std::vector<int>>("layers", layer);
0061 desc.add<std::vector<int>>("waferUs", waferU);
0062 desc.add<std::vector<int>>("waferVs", waferV);
0063 desc.add<std::vector<int>>("cellUs", cellU);
0064 desc.add<std::vector<int>>("cellVs", cellV);
0065 desc.add<std::vector<int>>("types", type);
0066 descriptions.add("hgcalGeometryRotTest", desc);
0067 }
0068
0069 void HGCalGeometryRotTest::beginRun(const edm::Run&, const edm::EventSetup& iSetup) {
0070 const auto& geomR = iSetup.getData(geomToken_);
0071 const HGCalGeometry* geom = &geomR;
0072 DetId::Detector det = (nameDetector_ == "HGCalEESensitive") ? DetId::HGCalEE : DetId::HGCalHSi;
0073 int layerF = *(layers_.begin());
0074 int layerL = *(--layers_.end());
0075 int layerOff = geom->topology().dddConstants().getLayerOffset();
0076 edm::LogVerbatim("HGCalGeom") << nameDetector_ << " with layers in the range " << layerF << ":" << layerL
0077 << " Offset " << layerOff << " and for " << waferU_.size() << " wafers and cells";
0078
0079 for (unsigned int k = 0; k < waferU_.size(); ++k) {
0080 for (auto lay : layers_) {
0081 HGCSiliconDetId detId(det, 1, types_[k], lay - layerOff, waferU_[k], waferV_[k], cellU_[k], cellV_[k]);
0082 GlobalPoint global = geom->getPosition(DetId(detId));
0083 double phi2 = global.phi();
0084 auto xy = geom->topology().dddConstants().waferPosition(lay - layerOff, waferU_[k], waferV_[k], true, false);
0085 double phi1 = std::atan2(xy.second, xy.first);
0086 edm::LogVerbatim("HGCalGeom") << "Layer: " << lay << " U " << waferU_[k] << " V " << waferV_[k] << " Position ("
0087 << xy.first << ", " << xy.second << ") phi " << convertRadToDeg(phi1);
0088 edm::LogVerbatim("HGCalGeom") << detId << " Position " << global << " phi " << convertRadToDeg(phi2);
0089 }
0090 }
0091 }
0092
0093
0094 #include "FWCore/Framework/interface/MakerMacros.h"
0095 DEFINE_FWK_MODULE(HGCalGeometryRotTest);