File indexing completed on 2024-06-06 04:26:58
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include <fstream>
0022 #include <iostream>
0023 #include <memory>
0024 #include <string>
0025 #include <vector>
0026 #include <stdlib.h>
0027 #include <cmath>
0028
0029
0030
0031 #include "FWCore/Framework/interface/Frameworkfwd.h"
0032 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0033
0034 #include "FWCore/Framework/interface/Event.h"
0035 #include "FWCore/Framework/interface/EventSetup.h"
0036 #include "FWCore/Framework/interface/MakerMacros.h"
0037
0038 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0041 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0042 #include "Geometry/HGCalCommonData/interface/HGCalParameters.h"
0043 #include "Geometry/HGCalCommonData/interface/HGCalCellUV.h"
0044 #include "Geometry/HGCalCommonData/interface/HGCalCell.h"
0045 #include "Geometry/HGCalCommonData/interface/HGCalCellOffset.h"
0046 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0047
0048 class HGCalCellOffsetTester : public edm::one::EDAnalyzer<> {
0049 public:
0050 explicit HGCalCellOffsetTester(const edm::ParameterSet&);
0051 ~HGCalCellOffsetTester() override = default;
0052 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0053
0054 void beginJob() override {}
0055 void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0056 void endJob() override {}
0057
0058 private:
0059 const double waferSize_;
0060 const int waferType_;
0061 const int placeIndex_;
0062 const int partial_;
0063 const double mouseBiteCut_;
0064 const double guardRingOffset_;
0065 const double sizeOffset_;
0066 std::ofstream outputFile;
0067 };
0068
0069 HGCalCellOffsetTester::HGCalCellOffsetTester(const edm::ParameterSet& iC)
0070 : waferSize_(iC.getParameter<double>("waferSize")),
0071 waferType_(iC.getParameter<int>("waferType")),
0072 placeIndex_(iC.getParameter<int>("cellPlacementIndex")),
0073 partial_(iC.getParameter<int>("cellType")),
0074 mouseBiteCut_(iC.getParameter<double>("mouseBiteCut")),
0075 guardRingOffset_(iC.getParameter<double>("guardRingOffset")),
0076 sizeOffset_(iC.getParameter<double>("sizeOffset")) {
0077 edm::LogVerbatim("HGCalGeom") << "Test positions for wafer of size " << waferSize_ << " Type " << waferType_
0078 << " Placement Index " << placeIndex_ << " GuardRing offset " << guardRingOffset_
0079 << " Mousebite cut " << mouseBiteCut_ << " SizeOffset " << sizeOffset_;
0080
0081 outputFile.open("nand.csv");
0082 if (!outputFile.is_open()) {
0083 edm::LogError("HGCalGeom") << "Could not open output file.";
0084 } else {
0085 outputFile << "x,y,u,v,\n";
0086 }
0087 }
0088
0089 void HGCalCellOffsetTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0090 edm::ParameterSetDescription desc;
0091 desc.add<double>("waferSize", 167.4408);
0092 desc.add<int>("waferType", 0);
0093 desc.add<int>("cellPlacementIndex", 11);
0094 desc.add<int>("cellType", 0);
0095 desc.add<double>("mouseBiteCut", 5.0);
0096 desc.add<double>("guardRingOffset", 0.9);
0097 desc.add<double>("sizeOffset", 0.435);
0098 descriptions.add("hgcalCellOffsetTester", desc);
0099 }
0100
0101
0102 void HGCalCellOffsetTester::analyze(const edm::Event&, const edm::EventSetup&) {
0103 const int nFine(12), nCoarse(8);
0104 int nCells = (waferType_ == 0) ? nFine : nCoarse;
0105 HGCalCellUV wafer(waferSize_, 0.0, nFine, nCoarse);
0106 HGCalCell wafer2(waferSize_, nFine, nCoarse);
0107 HGCalCellOffset offset(waferSize_, nFine, nCoarse, guardRingOffset_, mouseBiteCut_, sizeOffset_);
0108 edm::LogVerbatim("HGCalGeom") << "\nHGCalPartialCellTester:: nCells " << nCells << " and placement index "
0109 << placeIndex_ << "\n\n";
0110 for (int ui = 0; ui < 2 * nCells; ui++) {
0111 for (int vi = 0; vi < 2 * nCells; vi++) {
0112 if ((ui < 2 * nCells) && (vi < 2 * nCells) && ((vi - ui) < nCells) && ((ui - vi) <= nCells)) {
0113
0114 if (HGCalWaferMask::goodCell(ui, vi, partial_)) {
0115 std::pair<double, double> xy1 = wafer2.cellUV2XY2(ui, vi, placeIndex_, waferType_);
0116
0117 std::pair<int32_t, int32_t> uv1 =
0118 wafer.cellUVFromXY1(xy1.first, xy1.second, placeIndex_, waferType_, true, false);
0119
0120
0121
0122
0123
0124
0125
0126 std::pair<double, double> xyOffsetLD = offset.cellOffsetUV2XY1(ui, vi, placeIndex_, waferType_, partial_);
0127 auto area = offset.cellAreaUV(ui, vi, placeIndex_, waferType_, partial_, true);
0128
0129 outputFile << xyOffsetLD.first + xy1.first << "," << xyOffsetLD.second + xy1.second << "," << uv1.first << ","
0130 << uv1.second << "," << area << std::endl;
0131
0132 std::string comment = ((uv1.first != ui) || (uv1.second != vi))
0133 ? " ***** ERROR (u, v) from the methods dosent match *****"
0134 : "";
0135 edm::LogVerbatim("HGCalGeom") << "u = " << ui << " v = " << vi << " type = " << waferType_
0136 << " placement index " << placeIndex_ << " u " << uv1.first << " v "
0137 << uv1.second << " x " << xy1.first << " ,y " << xy1.second << " xoff "
0138 << xyOffsetLD.first << " ,yoff " << xyOffsetLD.second << " , area " << area
0139 << comment;
0140 }
0141 }
0142 }
0143 }
0144 }
0145
0146
0147 DEFINE_FWK_MODULE(HGCalCellOffsetTester);