Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-10-12 22:41:52

0001 // -*- C++ -*-
0002 //
0003 // Package:    HGCalWaferInFileTest
0004 // Class:      HGCalWaferInFileTest
0005 //
0006 /**\class HGCalWaferInFileTest HGCalWaferInFileTest.cc
0007  test/HGCalWaferInFileTest.cc
0008 
0009  Description: <one line class summary>
0010 
0011  Implementation:
0012      <Notes on implementation>
0013 */
0014 //
0015 // Original Author:  Sunanda Banerjee
0016 //         Created:  Mon 2020/06/24
0017 //
0018 //
0019 
0020 // system include files
0021 #include <sstream>
0022 #include <string>
0023 #include <vector>
0024 
0025 // user include files
0026 #include "FWCore/Framework/interface/Frameworkfwd.h"
0027 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/EventSetup.h"
0030 #include "FWCore/Framework/interface/MakerMacros.h"
0031 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0032 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0033 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0034 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0035 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0036 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0037 #include "Geometry/HGCalCommonData/interface/HGCalWaferIndex.h"
0038 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0039 
0040 class HGCalWaferInFileTest : public edm::one::EDAnalyzer<> {
0041 public:
0042   explicit HGCalWaferInFileTest(const edm::ParameterSet&);
0043 
0044   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0045 
0046   void beginJob() override {}
0047   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0048   void endJob() override {}
0049 
0050 private:
0051   std::vector<std::string> getPoints(
0052       double xpos, double ypos, double delX, double delY, double rin, double rout, int lay, int waferU, int waferV);
0053   std::vector<double> getCorners(double xpos, double ypos, double delX, double delY);
0054   std::pair<bool, std::string> getPoints(double xpos,
0055                                          double ypos,
0056                                          double delX,
0057                                          double delY,
0058                                          double rin,
0059                                          double rout,
0060                                          int part,
0061                                          int rotn,
0062                                          int layer,
0063                                          int waferU,
0064                                          int waferV);
0065 
0066   const std::string nameSense_, nameDetector_;
0067   const int verbosity_;
0068   const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> geomToken_;
0069 };
0070 
0071 HGCalWaferInFileTest::HGCalWaferInFileTest(const edm::ParameterSet& iC)
0072     : nameSense_(iC.getParameter<std::string>("NameSense")),
0073       nameDetector_(iC.getParameter<std::string>("NameDevice")),
0074       verbosity_(iC.getParameter<int>("Verbosity")),
0075       geomToken_(esConsumes<HGCalGeometry, IdealGeometryRecord>(edm::ESInputTag{"", nameSense_})) {
0076   edm::LogVerbatim("HGCalGeom") << "Test wafer characteristics for " << nameDetector_ << " using constants of "
0077                                 << nameSense_ << " with verbosity " << verbosity_;
0078 }
0079 
0080 void HGCalWaferInFileTest::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0081   edm::ParameterSetDescription desc;
0082   desc.add<std::string>("NameSense", "HGCalEESensitive");
0083   desc.add<std::string>("NameDevice", "HGCal EE");
0084   desc.add<int>("Verbosity", 0);
0085   descriptions.add("hgcalEEWaferInFileTest", desc);
0086 }
0087 
0088 // ------------ method called to produce the data  ------------
0089 void HGCalWaferInFileTest::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0090   const auto& geomR = iSetup.getData(geomToken_);
0091   const HGCalGeometry* geom = &geomR;
0092   const auto& hgdc = geom->topology().dddConstants();
0093   double delX = 0.5 * hgdc.getParameter()->waferSize_;
0094   double delY = 2.0 * delX / std::sqrt(3.0);
0095 
0096   edm::LogVerbatim("HGCalGeom") << nameDetector_ << "\nCheck Wafers in file are all valid for " << nameDetector_
0097                                 << std::endl;
0098   if (hgdc.waferHexagon8()) {
0099     DetId::Detector det = (nameSense_ == "HGCalHESiliconSensitive") ? DetId::HGCalHSi : DetId::HGCalEE;
0100     static std::vector<std::string> types = {"F", "b", "g", "gm", "a", "d", "dm", "c", "am", "bm", "X"};
0101     int layers = hgdc.layers(true);
0102     int layerf = hgdc.firstLayer();
0103     std::vector<int> miss(layers, 0);
0104     // See if all entries in the file are valid
0105     int bad1(0);
0106     for (unsigned int k = 0; k < hgdc.waferFileSize(); ++k) {
0107       int indx = hgdc.waferFileIndex(k);
0108       int layer = HGCalWaferIndex::waferLayer(indx);
0109       int waferU = HGCalWaferIndex::waferU(indx);
0110       int waferV = HGCalWaferIndex::waferV(indx);
0111       int type = std::get<0>(hgdc.waferFileInfo(k));
0112       HGCSiliconDetId id(det, 1, type, layer, waferU, waferV, 0, 0);
0113       if (!geom->topology().validModule(id, 3)) {
0114         int part = std::get<1>(hgdc.waferFileInfoFromIndex(indx));
0115         std::string typex = (part < static_cast<int>(types.size())) ? types[part] : "X";
0116         const auto& xy = hgdc.waferPosition(layer, waferU, waferV, true, true);
0117         const auto& rr = hgdc.rangeRLayer(layer, true);
0118         auto points = getPoints(xy.first, xy.second, delX, delY, rr.first, rr.second, layer, waferU, waferV);
0119         auto rpos = getCorners(xy.first, xy.second, delX, delY);
0120         std::ostringstream st1;
0121         st1 << "ID[" << k << "]: (" << (hgdc.getLayerOffset() + layer) << ", " << waferU << ", " << waferV << ", "
0122             << typex << ") at (" << std::setprecision(4) << xy.first << ", " << xy.second << ", "
0123             << hgdc.waferZ(layer, true) << ") not present with " << points.size() << " points:";
0124         for (auto point : points)
0125           st1 << " " << point;
0126         st1 << " in the region " << rr.first << ":" << rr.second << " Corners";
0127         for (auto point : rpos)
0128           st1 << " " << point;
0129         edm::LogVerbatim("HGCalGeom") << st1.str();
0130         ++bad1;
0131         if ((layer - layerf) < layers)
0132           ++miss[layer - layerf];
0133       }
0134     }
0135     edm::LogVerbatim("HGCalGeom") << "\n\nFinds " << bad1 << " invalid wafers among " << hgdc.waferFileSize()
0136                                   << " wafers in the list";
0137     for (unsigned int k = 0; k < miss.size(); ++k) {
0138       if (miss[k] > 0)
0139         edm::LogVerbatim("HGCalGeom") << "Layer[" << k << ":" << (layerf + k) << "] " << miss[k];
0140     }
0141     edm::LogVerbatim("HGCalGeom") << std::endl;
0142 
0143     // Now cross check the content (first type only)
0144     int allG(0), badT(0), badT1(0), badT2(0);
0145     for (unsigned int k = 0; k < hgdc.waferFileSize(); ++k) {
0146       int indx = hgdc.waferFileIndex(k);
0147       int type1 = std::get<0>(hgdc.waferFileInfo(k));
0148       int layer = HGCalWaferIndex::waferLayer(indx);
0149       int waferU = HGCalWaferIndex::waferU(indx);
0150       int waferV = HGCalWaferIndex::waferV(indx);
0151       int type2 = hgdc.waferType(layer, waferU, waferV, false);
0152       HGCSiliconDetId id(det, 1, type2, layer, waferU, waferV, 0, 0);
0153       if (geom->topology().validModule(id, 3)) {
0154         ++allG;
0155         bool typeOK = (type1 == type2);
0156         if (!typeOK) {
0157           ++badT;
0158           if (type1 == 0)
0159             ++badT1;
0160           else if (type2 == 0)
0161             ++badT2;
0162           const auto& xy = hgdc.waferPosition(layer, waferU, waferV, true, false);
0163           edm::LogVerbatim("HGCalGeom") << "ID[" << k << "]: (" << (hgdc.getLayerOffset() + layer) << ", " << waferU
0164                                         << ", " << waferV << ", " << type1 << ":" << type2 << ") at ("
0165                                         << std::setprecision(4) << xy.first << ", " << xy.second << ", "
0166                                         << hgdc.waferZ(layer, true) << ") failure flag " << typeOK;
0167         }
0168       }
0169     }
0170     edm::LogVerbatim("HGCalGeom") << "\n\nFinds " << badT << "[" << badT1 << ":" << badT2 << "] mismatch in type among "
0171                                   << allG << " wafers with the same indices\n";
0172 
0173     // Now cross check the content (partial and orientation)
0174     int allX(0), badG(0), badP(0), badP2(0), badR(0), badR2(0);
0175     std::vector<int> wrongP(layers, 0), wrongP2(layers, 0), wrongR(layers, 0), wrongR2(layers, 0);
0176     for (unsigned int k = 0; k < hgdc.waferFileSize(); ++k) {
0177       int indx = hgdc.waferFileIndex(k);
0178       int part1 = std::get<1>(hgdc.waferFileInfo(k));
0179       int rotn1 = std::get<2>(hgdc.waferFileInfo(k));
0180       int layer = HGCalWaferIndex::waferLayer(indx);
0181       int waferU = HGCalWaferIndex::waferU(indx);
0182       int waferV = HGCalWaferIndex::waferV(indx);
0183       int type2 = hgdc.waferType(layer, waferU, waferV, false);
0184       HGCSiliconDetId id(det, 1, type2, layer, waferU, waferV, 0, 0);
0185       if (geom->topology().validModule(id, 3)) {
0186         ++allX;
0187         int part2 = hgdc.waferTypeRotation(id.layer(), id.waferU(), id.waferV(), false, false).first;
0188         int rotn2 = hgdc.waferTypeRotation(id.layer(), id.waferU(), id.waferV(), false, false).second;
0189         bool partOK = ((part1 == part2) || ((part1 == HGCalTypes::WaferFull) && (part2 == HGCalTypes::WaferOut)));
0190         bool rotnOK = ((rotn1 == rotn2) || (part1 == HGCalTypes::WaferFull) || (part2 == HGCalTypes::WaferFull));
0191         bool partOK2 = (partOK) || (part2 < part1);
0192         bool rotnOK2(true);
0193         std::string miss;
0194         if (!partOK) {
0195           ++badP;
0196           if ((layer - layerf) < layers)
0197             ++wrongP[layer - layerf];
0198         }
0199         const auto& xy = hgdc.waferPosition(layer, waferU, waferV, true, false);
0200         const auto& rr = hgdc.rangeRLayer(layer, true);
0201 
0202         if (!partOK2) {
0203           ++badP2;
0204           if ((layer - layerf) < layers)
0205             ++wrongP2[layer - layerf];
0206           auto result =
0207               getPoints(xy.first, xy.second, delX, delY, rr.first, rr.second, part1, rotn1, layer, waferU, waferV);
0208           rotnOK2 = result.first;
0209           miss = result.second;
0210           if (!rotnOK2) {
0211             ++badR2;
0212             if ((layer - layerf) < layers)
0213               ++wrongR2[layer - layerf];
0214           }
0215         }
0216         if (!rotnOK) {
0217           ++badR;
0218           if ((layer - layerf) < layers)
0219             ++wrongR[layer - layerf];
0220         }
0221         if ((!partOK) || (!rotnOK)) {
0222           ++badG;
0223           if ((verbosity_ > 0) || ((!partOK2) || (!rotnOK2))) {
0224             std::string partx1 = (part1 < static_cast<int>(types.size())) ? types[part1] : "X";
0225             std::string partx2 = (part2 < static_cast<int>(types.size())) ? types[part2] : "X";
0226             auto points = getPoints(xy.first, xy.second, delX, delY, rr.first, rr.second, layer, waferU, waferV);
0227             auto rpos = getCorners(xy.first, xy.second, delX, delY);
0228             std::ostringstream st1;
0229             st1 << "ID[" << k << "]: (" << (hgdc.getLayerOffset() + layer) << ", " << waferU << ", " << waferV << ","
0230                 << type2 << ", " << partx1 << ":" << partx2 << ":" << part1 << ":" << part2 << ", " << rotn1 << ":"
0231                 << rotn2 << ") at (" << std::setprecision(4) << xy.first << ", " << xy.second << ", "
0232                 << hgdc.waferZ(layer, true) << ") failure flags (part = " << partOK << ":" << partOK2 << " rotn "
0233                 << rotnOK << ":" << rotnOK2 << " at " << miss << " with " << points.size() << " points:";
0234             for (auto point : points)
0235               st1 << " " << point;
0236             st1 << " in the region " << rr.first << ":" << rr.second << " Corners";
0237             for (auto point : rpos)
0238               st1 << " " << point;
0239             edm::LogVerbatim("HGCalGeom") << st1.str();
0240           }
0241         }
0242       }
0243     }
0244     edm::LogVerbatim("HGCalGeom") << "\n\nFinds " << badG << " (" << badP << ":" << badP2 << ":" << badR << ":" << badR2
0245                                   << ") mismatch in partial|orientation among " << allX
0246                                   << " wafers with the same indices";
0247     for (int k = 0; k < layers; ++k) {
0248       if ((wrongP[k] > 0) || (wrongR[k] > 0))
0249         edm::LogVerbatim("HGCalGeom") << "Layer[" << k << ":" << (layerf + k) << "] " << wrongP[k] << ":" << wrongP2[k]
0250                                       << ":" << wrongR[k] << ":" << wrongR2[k];
0251     }
0252     edm::LogVerbatim("HGCalGeom") << std::endl;
0253   }
0254 }
0255 
0256 std::vector<std::string> HGCalWaferInFileTest::getPoints(
0257     double xpos, double ypos, double delX, double delY, double rin, double rout, int layer, int waferU, int waferV) {
0258   std::vector<std::string> points;
0259   static const int corners = 6;
0260   static const int base = 10;
0261   static const std::string c0[corners] = {"A0", "B0", "C0", "D0", "E0", "F0"};
0262   double dx0[corners] = {
0263       0.0, HGCalTypes::c10 * delX, HGCalTypes::c10 * delX, 0.0, -HGCalTypes::c10 * delX, -HGCalTypes::c10 * delX};
0264   double dy0[corners] = {-HGCalTypes::c10 * delY,
0265                          -HGCalTypes::c50 * delY,
0266                          HGCalTypes::c50 * delY,
0267                          HGCalTypes::c10 * delY,
0268                          HGCalTypes::c50 * delY,
0269                          -HGCalTypes::c50 * delY};
0270   double xc[corners], yc[corners];
0271   int ncor(0), iok(0);
0272   for (int k = 0; k < corners; ++k) {
0273     xc[k] = xpos + dx0[k];
0274     yc[k] = ypos + dy0[k];
0275     double rpos = sqrt(xc[k] * xc[k] + yc[k] * yc[k]);
0276     if (rpos <= rout && rpos >= rin) {
0277       ++ncor;
0278       iok = iok * base + 1;
0279       points.emplace_back(c0[k]);
0280     } else {
0281       iok *= base;
0282     }
0283   }
0284   if (verbosity_ > 1)
0285     edm::LogVerbatim("HGCalGeom") << "I/p " << layer << ":" << waferU << ":" << waferV << ":" << xpos << ":" << ypos
0286                                   << ":" << delX << ":" << delY << ":" << rin << ":" << rout << " Corners " << ncor
0287                                   << " iok " << iok;
0288 
0289   static const int parts = 3;
0290   static const std::string c1[parts] = {"A1", "A2", "A3"};
0291   double dx1[parts] = {HGCalTypes::c22 * delX, HGCalTypes::c50 * delX, HGCalTypes::c77 * delX};
0292   double dy1[parts] = {-HGCalTypes::c88 * delY, -HGCalTypes::c75 * delY, -HGCalTypes::c61 * delY};
0293   if ((((iok / 10000) % 10) == 1) && (((iok / 100000) % 10) == 0)) {
0294     for (int k = 0; k < parts; ++k) {
0295       double xc1 = xpos + dx1[k];
0296       double yc1 = ypos + dy1[k];
0297       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0298       if (rpos <= rout && rpos >= rin) {
0299         points.emplace_back(c1[k]);
0300         break;
0301       }
0302     }
0303   }
0304   if ((((iok / 10000) % 10) == 0) && (((iok / 100000) % 10) == 1)) {
0305     for (int k1 = 0; k1 < parts; ++k1) {
0306       int k = parts - k1 - 1;
0307       double xc1 = xpos + dx1[k];
0308       double yc1 = ypos + dy1[k];
0309       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0310       if (rpos <= rout && rpos >= rin) {
0311         points.emplace_back(c1[k]);
0312         break;
0313       }
0314     }
0315   }
0316 
0317   static const std::string c2[parts] = {"B1", "B2", "B3"};
0318   double dx2[parts] = {HGCalTypes::c10 * delX, HGCalTypes::c10 * delX, HGCalTypes::c10 * delX};
0319   double dy2[parts] = {-HGCalTypes::c27 * delY, 0.0, HGCalTypes::c27 * delY};
0320   if ((((iok / 1000) % 10) == 1) && (((iok / 10000) % 10) == 0)) {
0321     for (int k = 0; k < parts; ++k) {
0322       double xc1 = xpos + dx2[k];
0323       double yc1 = ypos + dy2[k];
0324       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0325       if (rpos <= rout && rpos >= rin) {
0326         points.emplace_back(c2[k]);
0327         break;
0328       }
0329     }
0330   }
0331   if ((((iok / 1000) % 10) == 0) && (((iok / 10000) % 10) == 1)) {
0332     for (int k1 = 0; k1 < parts; ++k1) {
0333       int k = parts - k1 - 1;
0334       double xc1 = xpos + dx2[k];
0335       double yc1 = ypos + dy2[k];
0336       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0337       if (rpos <= rout && rpos >= rin) {
0338         points.emplace_back(c2[k]);
0339         break;
0340       }
0341     }
0342   }
0343 
0344   static const std::string c3[parts] = {"C1", "C2", "C3"};
0345   double dx3[parts] = {HGCalTypes::c77 * delX, HGCalTypes::c50 * delX, HGCalTypes::c22 * delX};
0346   double dy3[parts] = {HGCalTypes::c61 * delY, HGCalTypes::c75 * delY, HGCalTypes::c88 * delY};
0347   if ((((iok / 100) % 10) == 1) && (((iok / 1000) % 10) == 0)) {
0348     for (int k = 0; k < parts; ++k) {
0349       double xc1 = xpos + dx3[k];
0350       double yc1 = ypos + dy3[k];
0351       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0352       if (rpos <= rout && rpos >= rin) {
0353         points.emplace_back(c3[k]);
0354         break;
0355       }
0356     }
0357   }
0358   if ((((iok / 100) % 10) == 0) && (((iok / 1000) % 10) == 1)) {
0359     for (int k1 = 0; k1 < parts; ++k1) {
0360       int k = parts - k1 - 1;
0361       double xc1 = xpos + dx3[k];
0362       double yc1 = ypos + dy3[k];
0363       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0364       if (rpos <= rout && rpos >= rin) {
0365         points.emplace_back(c3[k]);
0366         break;
0367       }
0368     }
0369   }
0370 
0371   static const std::string c4[parts] = {"D1", "D2", "D3"};
0372   double dx4[parts] = {-HGCalTypes::c22 * delX, -HGCalTypes::c50 * delX, -HGCalTypes::c77 * delX};
0373   double dy4[parts] = {HGCalTypes::c88 * delY, HGCalTypes::c75 * delY, HGCalTypes::c61 * delY};
0374   if ((((iok / 10) % 10) == 1) && (((iok / 100) % 10) == 0)) {
0375     for (int k = 0; k < parts; ++k) {
0376       double xc1 = xpos + dx4[k];
0377       double yc1 = ypos + dy4[k];
0378       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0379       if (rpos <= rout && rpos >= rin) {
0380         points.emplace_back(c4[k]);
0381         break;
0382       }
0383     }
0384   }
0385   if ((((iok / 10) % 10) == 0) && (((iok / 100) % 10) == 1)) {
0386     for (int k1 = 0; k1 < parts; ++k1) {
0387       int k = parts - k1 - 1;
0388       double xc1 = xpos + dx4[k];
0389       double yc1 = ypos + dy4[k];
0390       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0391       if (rpos <= rout && rpos >= rin) {
0392         points.emplace_back(c4[k]);
0393         break;
0394       }
0395     }
0396   }
0397 
0398   static const std::string c5[parts] = {"E1", "E2", "E3"};
0399   double dx5[parts] = {-delX, -delX, -delX};
0400   double dy5[parts] = {HGCalTypes::c27 * delY, 0.0, -HGCalTypes::c27 * delY};
0401   if ((((iok / 1) % 10) == 1) && (((iok / 10) % 10) == 0)) {
0402     for (int k = 0; k < parts; ++k) {
0403       double xc1 = xpos + dx5[k];
0404       double yc1 = ypos + dy5[k];
0405       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0406       if (rpos <= rout && rpos >= rin) {
0407         points.emplace_back(c5[k]);
0408         break;
0409       }
0410     }
0411   }
0412   if ((((iok / 1) % 10) == 0) && (((iok / 10) % 10) == 1)) {
0413     for (int k1 = 0; k1 < parts; ++k1) {
0414       int k = parts - k1 - 1;
0415       double xc1 = xpos + dx5[k];
0416       double yc1 = ypos + dy5[k];
0417       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0418       if (rpos <= rout && rpos >= rin) {
0419         points.emplace_back(c5[k]);
0420         break;
0421       }
0422     }
0423   }
0424 
0425   static const std::string c6[parts] = {"F1", "F2", "F3"};
0426   double dx6[parts] = {-HGCalTypes::c77 * delX, -HGCalTypes::c50 * delX, -HGCalTypes::c22 * delX};
0427   double dy6[parts] = {-HGCalTypes::c61 * delY, -HGCalTypes::c75 * delY, -HGCalTypes::c88 * delY};
0428   if ((((iok / 100000) % 10) == 1) && (((iok / 1) % 10) == 0)) {
0429     for (int k = 0; k < parts; ++k) {
0430       double xc1 = xpos + dx6[k];
0431       double yc1 = ypos + dy6[k];
0432       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0433       if (rpos <= rout && rpos >= rin) {
0434         points.emplace_back(c6[k]);
0435         break;
0436       }
0437     }
0438   }
0439   if ((((iok / 100000) % 10) == 0) && (((iok / 1) % 10) == 1)) {
0440     for (int k1 = 0; k1 < parts; ++k1) {
0441       int k = parts - k1 - 1;
0442       double xc1 = xpos + dx6[k];
0443       double yc1 = ypos + dy6[k];
0444       double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0445       if (rpos <= rout && rpos >= rin) {
0446         points.emplace_back(c6[k]);
0447         break;
0448       }
0449     }
0450   }
0451 
0452   return points;
0453 }
0454 
0455 std::vector<double> HGCalWaferInFileTest::getCorners(double xpos, double ypos, double delX, double delY) {
0456   std::vector<double> points;
0457   static const int corners = 6;
0458   double dx0[corners] = {
0459       0.0, HGCalTypes::c10 * delX, HGCalTypes::c10 * delX, 0.0, -HGCalTypes::c10 * delX, -HGCalTypes::c10 * delX};
0460   double dy0[corners] = {-HGCalTypes::c10 * delY,
0461                          -HGCalTypes::c50 * delY,
0462                          HGCalTypes::c50 * delY,
0463                          HGCalTypes::c10 * delY,
0464                          HGCalTypes::c50 * delY,
0465                          -HGCalTypes::c50 * delY};
0466   for (int k = 0; k < corners; ++k) {
0467     double xc = xpos + dx0[k];
0468     double yc = ypos + dy0[k];
0469     double rpos = sqrt(xc * xc + yc * yc);
0470     points.emplace_back(rpos);
0471   }
0472   return points;
0473 }
0474 
0475 std::pair<bool, std::string> HGCalWaferInFileTest::getPoints(double xpos,
0476                                                              double ypos,
0477                                                              double delX,
0478                                                              double delY,
0479                                                              double rin,
0480                                                              double rout,
0481                                                              int part,
0482                                                              int rotn,
0483                                                              int layer,
0484                                                              int waferU,
0485                                                              int waferV) {
0486   std::string point;
0487   static const int corners = 6;
0488   static const int corner2 = 12;
0489   static const int base = 10;
0490   static const int base2 = 100;
0491   static const std::string c0[corners] = {"A0", "B0", "C0", "D0", "E0", "F0"};
0492   static const std::string c1[corners] = {"A2", "B2", "C2", "D2", "E2", "F2"};
0493   static const std::string c2[corner2] = {"A1", "A3", "B1", "B3", "C1", "C3", "D1", "D3", "E1", "E3", "F1", "F3"};
0494   double dx0[corners] = {
0495       0.0, HGCalTypes::c10 * delX, HGCalTypes::c10 * delX, 0.0, -HGCalTypes::c10 * delX, -HGCalTypes::c10 * delX};
0496   double dy0[corners] = {-HGCalTypes::c10 * delY,
0497                          -HGCalTypes::c50 * delY,
0498                          HGCalTypes::c50 * delY,
0499                          HGCalTypes::c10 * delY,
0500                          HGCalTypes::c50 * delY,
0501                          -HGCalTypes::c50 * delY};
0502   double dx1[corners] = {HGCalTypes::c50 * delX,
0503                          HGCalTypes::c10 * delX,
0504                          HGCalTypes::c50 * delX,
0505                          -HGCalTypes::c50 * delX,
0506                          -HGCalTypes::c10 * delX,
0507                          -HGCalTypes::c50 * delX};
0508   double dy1[corners] = {
0509       -HGCalTypes::c75 * delY, 0.0, HGCalTypes::c75 * delY, HGCalTypes::c75 * delY, 0.0, -HGCalTypes::c75 * delY};
0510   double dx2[corner2] = {HGCalTypes::c22 * delX,
0511                          HGCalTypes::c77 * delX,
0512                          HGCalTypes::c10 * delX,
0513                          HGCalTypes::c10 * delX,
0514                          HGCalTypes::c77 * delX,
0515                          HGCalTypes::c22 * delX,
0516                          -HGCalTypes::c22 * delX,
0517                          -HGCalTypes::c77 * delX,
0518                          -HGCalTypes::c10 * delX,
0519                          -HGCalTypes::c10 * delX,
0520                          -HGCalTypes::c77 * delX,
0521                          -HGCalTypes::c22 * delX};
0522   double dy2[corner2] = {-HGCalTypes::c88 * delY,
0523                          -HGCalTypes::c61 * delY,
0524                          -HGCalTypes::c27 * delY,
0525                          HGCalTypes::c27 * delY,
0526                          HGCalTypes::c61 * delY,
0527                          HGCalTypes::c88 * delY,
0528                          HGCalTypes::c88 * delY,
0529                          HGCalTypes::c61 * delY,
0530                          HGCalTypes::c27 * delY,
0531                          -HGCalTypes::c27 * delY,
0532                          -HGCalTypes::c61 * delY,
0533                          -HGCalTypes::c88 * delY};
0534   bool ok(true);
0535   switch (part) {
0536     case (HGCalTypes::WaferThree): {
0537       static const int nc0[corners] = {450, 150, 201, 312, 423, 534};
0538       int nc = nc0[rotn];
0539       for (int k1 = 0; k1 < 3; ++k1) {
0540         int k = nc % base;
0541         double xc1 = xpos + dx0[k];
0542         double yc1 = ypos + dy0[k];
0543         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0544         if (rpos <= rout && rpos >= rin) {
0545           point = c0[k];
0546           ok = false;
0547           break;
0548         }
0549         nc /= base;
0550       }
0551       break;
0552     }
0553     case (HGCalTypes::WaferSemi2): {
0554       static const int nc10[corners] = {450, 150, 201, 312, 423, 534};
0555       static const int nc11[corners] = {700, 902, 1104, 106, 308, 510};
0556       int nc = nc10[rotn];
0557       for (int k1 = 0; k1 < 3; ++k1) {
0558         int k = nc % base;
0559         double xc1 = xpos + dx0[k];
0560         double yc1 = ypos + dy0[k];
0561         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0562         if (rpos <= rout && rpos >= rin) {
0563           point = c0[k];
0564           ok = false;
0565           break;
0566         }
0567         nc /= base;
0568       }
0569       nc = nc11[rotn];
0570       for (int k1 = 0; k1 < 2; ++k1) {
0571         int k = nc % base2;
0572         double xc1 = xpos + dx2[k];
0573         double yc1 = ypos + dy2[k];
0574         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0575         if (rpos <= rout && rpos >= rin) {
0576           point = c2[k];
0577           ok = false;
0578           break;
0579         }
0580         nc /= base2;
0581       }
0582       break;
0583     }
0584     case (HGCalTypes::WaferSemi): {
0585       static const int nc20[corners] = {450, 150, 201, 312, 423, 534};
0586       static const int nc21[corners] = {30, 14, 25, 30, 41, 52};
0587       int nc = nc20[rotn];
0588       for (int k1 = 0; k1 < 3; ++k1) {
0589         int k = nc % base;
0590         double xc1 = xpos + dx0[k];
0591         double yc1 = ypos + dy0[k];
0592         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0593         if (rpos <= rout && rpos >= rin) {
0594           point = c0[k];
0595           ok = false;
0596           break;
0597         }
0598         nc /= base;
0599       }
0600       nc = nc21[rotn];
0601       for (int k1 = 0; k1 < 2; ++k1) {
0602         int k = nc % base;
0603         double xc1 = xpos + dx1[k];
0604         double yc1 = ypos + dy1[k];
0605         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0606         if (rpos <= rout && rpos >= rin) {
0607           point = c1[k];
0608           ok = false;
0609           break;
0610         }
0611         nc /= base;
0612       }
0613       break;
0614     }
0615     case (HGCalTypes::WaferHalf): {
0616       static const int nc3[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0617       int nc = nc3[rotn];
0618       for (int k1 = 0; k1 < 4; ++k1) {
0619         int k = nc % base;
0620         double xc1 = xpos + dx0[k];
0621         double yc1 = ypos + dy0[k];
0622         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0623         if (rpos <= rout && rpos >= rin) {
0624           point = c0[k];
0625           ok = false;
0626           break;
0627         }
0628         nc /= base;
0629       }
0630       break;
0631     }
0632     case (HGCalTypes::WaferChopTwoM): {
0633       static const int nc40[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0634       static const int nc41[corners] = {500, 702, 904, 1106, 108, 310};
0635       int nc = nc40[rotn];
0636       for (int k1 = 0; k1 < 4; ++k1) {
0637         int k = nc % base;
0638         double xc1 = xpos + dx0[k];
0639         double yc1 = ypos + dy0[k];
0640         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0641         if (rpos <= rout && rpos >= rin) {
0642           point = c0[k];
0643           ok = false;
0644           break;
0645         }
0646         nc /= base;
0647       }
0648       nc = nc41[rotn];
0649       for (int k1 = 0; k1 < 2; ++k1) {
0650         int k = nc % base2;
0651         double xc1 = xpos + dx2[k];
0652         double yc1 = ypos + dy2[k];
0653         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0654         if (rpos <= rout && rpos >= rin) {
0655           point = c2[k];
0656           ok = false;
0657           break;
0658         }
0659         nc /= base2;
0660       }
0661       break;
0662     }
0663     case (HGCalTypes::WaferChopTwo): {
0664       static const int nc50[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0665       static const int nc51[corners] = {20, 13, 24, 35, 40, 51};
0666       int nc = nc50[rotn];
0667       for (int k1 = 0; k1 < 4; ++k1) {
0668         int k = nc % base;
0669         double xc1 = xpos + dx0[k];
0670         double yc1 = ypos + dy0[k];
0671         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0672         if (rpos <= rout && rpos >= rin) {
0673           point = c0[k];
0674           ok = false;
0675           break;
0676         }
0677         nc /= base;
0678       }
0679       nc = nc51[rotn];
0680       for (int k1 = 0; k1 < 2; ++k1) {
0681         int k = nc % base;
0682         double xc1 = xpos + dx1[k];
0683         double yc1 = ypos + dy1[k];
0684         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0685         if (rpos <= rout && rpos >= rin) {
0686           point = c1[k];
0687           ok = false;
0688           break;
0689         }
0690         nc /= base;
0691       }
0692       break;
0693     }
0694     case (HGCalTypes::WaferFive): {
0695       static const int nc6[corners] = {23450, 13450, 24501, 35012, 40123, 51234};
0696       int nc = nc6[rotn];
0697       for (int k1 = 0; k1 < 5; ++k1) {
0698         int k = nc % base;
0699         double xc1 = xpos + dx0[k];
0700         double yc1 = ypos + dy0[k];
0701         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0702         if (rpos <= rout && rpos >= rin) {
0703           point = c0[k];
0704           ok = false;
0705           break;
0706         }
0707       }
0708       break;
0709     }
0710     default: {
0711       for (int k = 0; k < corners; ++k) {
0712         double xc1 = xpos + dx0[k];
0713         double yc1 = ypos + dy0[k];
0714         double rpos = sqrt(xc1 * xc1 + yc1 * yc1);
0715         if (rpos <= rout && rpos >= rin) {
0716           point = c0[k];
0717           ok = false;
0718           break;
0719         }
0720       }
0721       break;
0722     }
0723   }
0724   if (verbosity_ > 1)
0725     edm::LogVerbatim("HGCalGeom") << "I/p " << layer << ":" << waferU << ":" << waferV << ":" << xpos << ":" << ypos
0726                                   << ":" << delX << ":" << delY << ":" << rin << ":" << rout << " Results " << ok
0727                                   << " point " << point;
0728   return std::make_pair(ok, point);
0729 }
0730 
0731 // define this as a plug-in
0732 DEFINE_FWK_MODULE(HGCalWaferInFileTest);