File indexing completed on 2024-11-22 02:45:25
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include <cmath>
0022 #include <cstdlib>
0023 #include <fstream>
0024 #include <iostream>
0025 #include <memory>
0026 #include <string>
0027 #include <vector>
0028
0029
0030
0031 #include "FWCore/Framework/interface/Frameworkfwd.h"
0032 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0033
0034 #include <CLHEP/Vector/ThreeVector.h>
0035 #include "FWCore/Framework/interface/Event.h"
0036 #include "FWCore/Framework/interface/EventSetup.h"
0037 #include "FWCore/Framework/interface/MakerMacros.h"
0038 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0041 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0042 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0043 #include "Geometry/HGCalCommonData/interface/HGCalParameters.h"
0044 #include "Geometry/HGCalCommonData/interface/HGCalCellUV.h"
0045 #include "Geometry/HGCalCommonData/interface/HGCalCell.h"
0046 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0047 #include "Geometry/HGCalCommonData/interface/HGCalWaferType.h"
0048 #include "SimG4CMS/Calo/interface/HGCMouseBite.h"
0049 #include "SimG4CMS/Calo/interface/HGCGuardRing.h"
0050 #include "SimG4CMS/Calo/interface/HGCGuardRingPartial.h"
0051 #include "G4ThreeVector.hh"
0052
0053 class HGCalMouseBiteTester : public edm::one::EDAnalyzer<> {
0054 public:
0055 explicit HGCalMouseBiteTester(const edm::ParameterSet&);
0056 ~HGCalMouseBiteTester() override = default;
0057 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0058
0059 void beginJob() override {}
0060 void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0061 void endJob() override {}
0062
0063 private:
0064 const std::string nameSense_;
0065 const int waferU_;
0066 const int waferV_;
0067 const int nTrials_;
0068 const int layer_;
0069 const edm::ESGetToken<HGCalDDDConstants, IdealGeometryRecord> dddToken_;
0070 std::ofstream outputFile;
0071 };
0072
0073 HGCalMouseBiteTester::HGCalMouseBiteTester(const edm::ParameterSet& iC)
0074 : nameSense_(iC.getParameter<std::string>("nameSense")),
0075 waferU_(iC.getParameter<int>("waferU")),
0076 waferV_(iC.getParameter<int>("waferV")),
0077 nTrials_(iC.getParameter<int>("numbberOfTrials")),
0078 layer_(iC.getParameter<int>("layer")),
0079 dddToken_(esConsumes<HGCalDDDConstants, IdealGeometryRecord>(edm::ESInputTag{"", nameSense_})) {
0080 edm::LogVerbatim("HGCalGeom") << "Test Guard_Ring for wafer in layer" << layer_ << " U " << waferU_ << " V "
0081 << waferV_ << " with " << nTrials_ << " trials";
0082
0083 outputFile.open("full1.csv");
0084 if (!outputFile.is_open()) {
0085 edm::LogError("HGCalGeom") << "Could not open output file.";
0086 } else {
0087 outputFile << "x,y,u,v,\n";
0088 }
0089 }
0090
0091 void HGCalMouseBiteTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0092 edm::ParameterSetDescription desc;
0093 desc.add<std::string>("nameSense", "HGCalEESensitive");
0094 desc.add<int>("waferU", 1);
0095 desc.add<int>("waferV", 9);
0096 desc.add<int>("numbberOfTrials", 1000000);
0097 desc.add<int>("layer", 1);
0098 descriptions.add("hgcalMouseBiteTester", desc);
0099 }
0100
0101
0102 void HGCalMouseBiteTester::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0103 const HGCalDDDConstants& hgcons_ = iSetup.getData(dddToken_);
0104 double waferSize_(hgcons_.waferSize(false));
0105 int zside(1);
0106 int layertype = hgcons_.layerType(layer_);
0107 int frontBack = HGCalTypes::layerFrontBack(layertype);
0108 const std::vector<double> angle_{90.0, 30.0};
0109 int index = HGCalWaferIndex::waferIndex(layer_, waferU_, waferV_);
0110 int partialType_ = HGCalWaferType::getPartial(index, hgcons_.getParameter()->waferInfoMap_);
0111 int orient = HGCalWaferType::getOrient(index, hgcons_.getParameter()->waferInfoMap_);
0112 int placeIndex_ = HGCalCell::cellPlacementIndex(zside, frontBack, orient);
0113 int waferType_ = HGCalWaferType::getType(index, hgcons_.getParameter()->waferInfoMap_);
0114 double mouseBiteCut_ = hgcons_.mouseBite(false);
0115 bool v17OrLess = hgcons_.v17OrLess();
0116 HGCGuardRing guardRing_(hgcons_);
0117 HGCGuardRingPartial guardRingPartial_(hgcons_);
0118 HGCMouseBite mouseBite_(hgcons_, angle_, (waferSize_ * tan(30.0 * CLHEP::deg) - mouseBiteCut_), true);
0119 const int nFine(12), nCoarse(8);
0120 double r2 = 0.5 * waferSize_;
0121 double R2 = 2 * r2 / sqrt(3);
0122 int nCells = (waferType_ == 0) ? nFine : nCoarse;
0123 std::cout << "start" << std::endl;
0124 HGCalCellUV wafer(waferSize_, 0.0, nFine, nCoarse);
0125 HGCalCell wafer2(waferSize_, nFine, nCoarse);
0126 std::pair<double, double> xy = hgcons_.waferPosition(layer_, waferU_, waferV_, false, false);
0127 double x0 = (zside > 0) ? xy.first : -xy.first;
0128 double y0 = xy.second;
0129 std::ofstream guard_ring("Guard_ring.csv");
0130 std::ofstream guard_ring_partial("Guard_ring_partial.csv");
0131 std::ofstream mouse_bite("Mouse_bite.csv");
0132 std::ofstream selected("Selected.csv");
0133 edm::LogVerbatim("HGCalGeom") << "\nHGCalMouseBiteTester:: nCells " << nCells << " FrontBack " << frontBack
0134 << " Wafer Size " << waferSize_ << " and placement index " << placeIndex_
0135 << " WaferType " << waferType_ << " Partial " << partialType_ << " WaferX " << x0
0136 << " WaferY " << y0 << "\n\n";
0137 auto start_t = std::chrono::high_resolution_clock::now();
0138 std::cout << "v17 ? " << hgcons_.v17OrLess() << std::endl;
0139 for (int i = 0; i < nTrials_; i++) {
0140 double xi = (2 * r2 * static_cast<double>(rand()) / RAND_MAX) - r2;
0141 double yi = (2 * R2 * static_cast<double>(rand()) / RAND_MAX) - R2;
0142 bool goodPoint = true;
0143 int ug = 0;
0144 int vg = 0;
0145 if (partialType_ == 11 || partialType_ == 13 || partialType_ == 15 || partialType_ == 21 || partialType_ == 23 ||
0146 partialType_ == 25 || partialType_ == 0) {
0147 ug = 0;
0148 vg = 0;
0149 } else if (partialType_ == 12 || partialType_ == 14 || partialType_ == 16 || partialType_ == 22 ||
0150 partialType_ == 24) {
0151 ug = nCells + 1;
0152 vg = 2 * (nCells - 1);
0153 }
0154 std::pair<double, double> xyg = wafer2.cellUV2XY2(ug, vg, placeIndex_, waferType_);
0155 std::vector<std::pair<double, double> > wxy =
0156 HGCalWaferMask::waferXY(0, placeIndex_, waferSize_, 0.0, 0.0, 0.0, v17OrLess);
0157 for (unsigned int i = 0; i < (wxy.size() - 1); ++i) {
0158 double xp1 = wxy[i].first;
0159 double yp1 = wxy[i].second;
0160 double xp2 = wxy[i + 1].first;
0161 double yp2 = wxy[i + 1].second;
0162 if ((((xi - xp1) / (xp2 - xp1)) - ((yi - yp1) / (yp2 - yp1))) *
0163 (((xyg.first - xp1) / (xp2 - xp1)) - ((xyg.second - yp1) / (yp2 - yp1))) <=
0164 0) {
0165 goodPoint = false;
0166 }
0167 }
0168 if (goodPoint) {
0169 G4ThreeVector point(xi, yi, 0.0);
0170 std::pair<int32_t, int32_t> uv5;
0171 if (hgcons_.v17OrLess()) {
0172 uv5 = wafer.cellUVFromXY1(xi, yi, placeIndex_, waferType_, partialType_, true, false);
0173 } else {
0174 uv5 = wafer.cellUVFromXY2(xi, yi, placeIndex_, waferType_, partialType_, true, false);
0175 }
0176 if (guardRing_.exclude(point, zside, frontBack, layer_, waferU_, waferV_)) {
0177 guard_ring << xi << "," << yi << std::endl;
0178 }
0179
0180 if (guardRingPartial_.exclude(point, zside, frontBack, layer_, waferU_, waferV_)) {
0181 guard_ring_partial << xi << "," << yi << std::endl;
0182 } else if (mouseBite_.exclude(point, zside, layer_, waferU_, waferV_)) {
0183 mouse_bite << xi << "," << yi << std::endl;
0184 } else {
0185 selected << xi << "," << yi << std::endl;
0186 outputFile << xi << "," << yi << "," << uv5.first << "," << uv5.second << "," << std::endl;
0187 }
0188 }
0189 }
0190 guard_ring.close();
0191 guard_ring_partial.close();
0192 mouse_bite.close();
0193 selected.close();
0194 outputFile.close();
0195 auto end_t = std::chrono::high_resolution_clock::now();
0196 auto diff_t = end_t - start_t;
0197 edm::LogVerbatim("HGCalGeom") << "Execution time for " << nTrials_
0198 << " events = " << std::chrono::duration<double, std::milli>(diff_t).count() << " ms";
0199 }
0200
0201
0202 DEFINE_FWK_MODULE(HGCalMouseBiteTester);