Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-08 22:57:06

0001 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0002 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0003 #include "Geometry/HGCalCommonData/interface/HGCalGeomTools.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 
0006 #include <algorithm>
0007 #include <sstream>
0008 
0009 //#define EDM_ML_DEBUG
0010 
0011 bool HGCalWaferMask::maskCell(int u, int v, int n, int ncor, int fcor, int corners) {
0012   /*
0013 Masks each cell (or not) according to its wafer and cell position (detId) and to the user needs (corners).
0014 Each wafer has k_CornerSize corners which are defined in anti-clockwise order starting from the corner at the top, which is always #0. 'ncor' denotes the number of corners inside the physical region. 'fcor' is the defined to be the first corner that appears inside the detector's physical volume in anti-clockwise order. 
0015 The argument 'corners' controls the types of wafers the user wants: for instance, corners=3 masks all wafers that have at least 3 corners inside the physical region. 
0016  */
0017   bool mask(false);
0018   if (ncor < corners) {
0019     mask = true;
0020   } else {
0021     if (ncor == HGCalGeomTools::k_fourCorners) {
0022       switch (fcor) {
0023         case (0): {
0024           mask = (v >= n);
0025           break;
0026         }
0027         case (1): {
0028           mask = (u >= n);
0029           break;
0030         }
0031         case (2): {
0032           mask = (u > v);
0033           break;
0034         }
0035         case (3): {
0036           mask = (v < n);
0037           break;
0038         }
0039         case (4): {
0040           mask = (u < n);
0041           break;
0042         }
0043         default: {
0044           mask = (u <= v);
0045           break;
0046         }
0047       }
0048     } else {
0049       switch (fcor) {
0050         case (0): {
0051           if (ncor == HGCalGeomTools::k_threeCorners) {
0052             mask = !((u > 2 * v) && (v < n));
0053           } else {
0054             mask = ((u >= n) && (v >= n) && ((u + v) > (3 * n - 2)));
0055           }
0056           break;
0057         }
0058         case (1): {
0059           if (ncor == HGCalGeomTools::k_threeCorners) {
0060             mask = !((u + v) < n);
0061           } else {
0062             mask = ((u >= n) && (u > v) && ((2 * u - v) > 2 * n));
0063           }
0064           break;
0065         }
0066         case (2): {
0067           if (ncor == HGCalGeomTools::k_threeCorners) {
0068             mask = !((u < n) && (v > u) && (v > (2 * u - 1)));
0069           } else {
0070             mask = ((u > 2 * v) && (v < n));
0071           }
0072           break;
0073         }
0074         case (3): {
0075           if (ncor == HGCalGeomTools::k_threeCorners) {
0076             mask = !((v >= u) && ((2 * v - u) > (2 * n - 2)));
0077           } else {
0078             mask = ((u + v) < n);
0079           }
0080           break;
0081         }
0082         case (4): {
0083           if (ncor == HGCalGeomTools::k_threeCorners) {
0084             mask = !((u >= n) && (v >= n) && ((u + v) > (3 * n - 2)));
0085           } else {
0086             mask = ((u < n) && (v > u) && (v > (2 * u - 1)));
0087           }
0088           break;
0089         }
0090         default: {
0091           if (ncor == HGCalGeomTools::k_threeCorners) {
0092             mask = !((u >= n) && (u > v) && ((2 * u - v) > 2 * n));
0093           } else {
0094             mask = ((v >= u) && ((2 * v - u) > (2 * n - 2)));
0095           }
0096           break;
0097         }
0098       }
0099     }
0100   }
0101 #ifdef EDM_ML_DEBUG
0102   edm::LogVerbatim("HGCalGeom") << "Corners: " << ncor << ":" << fcor << " N " << n << " u " << u << " v " << v
0103                                 << " Mask " << mask;
0104 #endif
0105   return mask;
0106 }
0107 
0108 bool HGCalWaferMask::goodCell(int u, int v, int n, int type, int rotn) {
0109   // for V15 and V16
0110   bool good(false);
0111   int n2 = n / 2;
0112   int n4 = n / 4;
0113   int n3 = (n + 1) / 3;
0114   switch (type) {
0115     case (HGCalTypes::WaferFull): {  //WaferFull
0116       good = true;
0117       break;
0118     }
0119     case (HGCalTypes::WaferFive): {  //WaferFive
0120       switch (rotn) {
0121         case (HGCalTypes::WaferCorner0): {
0122           int u2 = u / 2;
0123           good = ((v - u2) < n);
0124           break;
0125         }
0126         case (HGCalTypes::WaferCorner1): {
0127           good = ((v + u) < (3 * n - 1));
0128           break;
0129         }
0130         case (HGCalTypes::WaferCorner2): {
0131           int v2 = (v + 1) / 2;
0132           good = ((u - v2) < n);
0133           break;
0134         }
0135         case (HGCalTypes::WaferCorner3): {
0136           int u2 = (u + 1) / 2;
0137           good = (u2 <= v);
0138           break;
0139         }
0140         case (HGCalTypes::WaferCorner4): {
0141           good = ((v + u) >= n);
0142           break;
0143         }
0144         default: {
0145           int v2 = v / 2;
0146           good = (u > v2);
0147           break;
0148         }
0149       }
0150       break;
0151     }
0152     case (HGCalTypes::WaferChopTwo): {  //WaferChopTwo
0153       switch (rotn) {
0154         case (HGCalTypes::WaferCorner0): {
0155           good = (v < (3 * n2));
0156           break;
0157         }
0158         case (HGCalTypes::WaferCorner1): {
0159           good = (u < (3 * n2));
0160           break;
0161         }
0162         case (HGCalTypes::WaferCorner2): {
0163           good = ((u - v) <= n2);
0164           break;
0165         }
0166         case (HGCalTypes::WaferCorner3): {
0167           good = (v >= n2);
0168           break;
0169         }
0170         case (HGCalTypes::WaferCorner4): {
0171           good = (u >= n2);
0172           break;
0173         }
0174         default: {
0175           good = ((v - u) < n2);
0176           break;
0177         }
0178       }
0179       break;
0180     }
0181     case (HGCalTypes::WaferChopTwoM): {  //WaferChopTwoM
0182       switch (rotn) {
0183         case (HGCalTypes::WaferCorner0): {
0184           good = (v < (5 * n4));
0185           break;
0186         }
0187         case (HGCalTypes::WaferCorner1): {
0188           good = (u < (5 * n4));
0189           break;
0190         }
0191         case (HGCalTypes::WaferCorner2): {
0192           good = ((u - v) <= n4);
0193           break;
0194         }
0195         case (HGCalTypes::WaferCorner3): {
0196           good = (v >= (3 * n4));
0197           break;
0198         }
0199         case (HGCalTypes::WaferCorner4): {
0200           good = (u >= (3 * n4));
0201           break;
0202         }
0203         default: {
0204           good = ((v - u) < n4);
0205           break;
0206         }
0207       }
0208       break;
0209     }
0210     case (HGCalTypes::WaferHalf): {  //WaferHalf
0211       switch (rotn) {
0212         case (HGCalTypes::WaferCorner0): {
0213           good = (v < n);
0214           break;
0215         }
0216         case (HGCalTypes::WaferCorner1): {
0217           good = (u < n);
0218           break;
0219         }
0220         case (HGCalTypes::WaferCorner2): {
0221           good = (v >= u);
0222           break;
0223         }
0224         case (HGCalTypes::WaferCorner3): {
0225           good = (v >= n);
0226           break;
0227         }
0228         case (HGCalTypes::WaferCorner4): {
0229           good = (u >= n);
0230           break;
0231         }
0232         default: {
0233           good = (u > v);
0234           break;
0235         }
0236       }
0237       break;
0238     }
0239     case (HGCalTypes::WaferSemi): {  //WaferSemi
0240       switch (rotn) {
0241         case (HGCalTypes::WaferCorner0): {
0242           good = ((u + v) < (2 * n));
0243           break;
0244         }
0245         case (HGCalTypes::WaferCorner1): {
0246           good = ((2 * u - v) < n);
0247           break;
0248         }
0249         case (HGCalTypes::WaferCorner2): {
0250           good = ((2 * v - u) >= n);
0251           break;
0252         }
0253         case (HGCalTypes::WaferCorner3): {
0254           good = ((u + v) >= (2 * n));
0255           break;
0256         }
0257         case (HGCalTypes::WaferCorner4): {
0258           good = ((2 * u - v) > n);
0259           break;
0260         }
0261         default: {
0262           good = ((2 * v - u) < n);
0263           break;
0264         }
0265       }
0266       break;
0267     }
0268     case (HGCalTypes::WaferThree): {  //WaferThree
0269       switch (rotn) {
0270         case (HGCalTypes::WaferCorner0): {
0271           good = ((v + u) < n);
0272           break;
0273         }
0274         case (HGCalTypes::WaferCorner1): {
0275           int v2 = v / 2;
0276           good = (u <= v2);
0277           break;
0278         }
0279         case (HGCalTypes::WaferCorner2): {
0280           int u2 = (u / 2);
0281           good = ((v - u2) >= n);
0282           break;
0283         }
0284         case (HGCalTypes::WaferCorner3): {
0285           good = ((v + u) >= (3 * n - 1));
0286           break;
0287         }
0288         case (HGCalTypes::WaferCorner4): {
0289           int v2 = ((v + 1) / 2);
0290           good = ((u - v2) >= n);
0291           break;
0292         }
0293         default: {
0294           int u2 = ((u + 1) / 2);
0295           good = (v < u2);
0296           break;
0297         }
0298       }
0299       break;
0300     }
0301     case (HGCalTypes::WaferSemi2): {  //WaferSemi2
0302       switch (rotn) {
0303         case (HGCalTypes::WaferCorner0): {
0304           good = ((u + v) < (4 * n3));
0305           break;
0306         }
0307         case (HGCalTypes::WaferCorner1): {
0308           good = ((2 * u - v) < n2);
0309           break;
0310         }
0311         case (HGCalTypes::WaferCorner2): {
0312           int u2 = ((u + 1) / 2);
0313           good = ((v - u2) >= (3 * n4));
0314           break;
0315         }
0316         case (HGCalTypes::WaferCorner3): {
0317           good = ((u + v) >= (5 * n2));
0318           break;
0319         }
0320         case (HGCalTypes::WaferCorner4): {
0321           good = ((2 * u - v) > (3 * n2));
0322           break;
0323         }
0324         default: {
0325           int u2 = ((n == 8) ? ((u + 1) / 2) : (u / 2));
0326           good = ((v - u2) < n4);
0327           break;
0328         }
0329       }
0330       break;
0331     }
0332     case (HGCalTypes::WaferFive2): {  //WaferFive2
0333       switch (rotn) {
0334         case (HGCalTypes::WaferCorner0): {
0335           good = ((2 * v - u) <= (3 * n2));
0336           break;
0337         }
0338         case (HGCalTypes::WaferCorner1): {
0339           good = ((u + v) < (5 * n2));
0340           break;
0341         }
0342         case (HGCalTypes::WaferCorner2): {
0343           good = ((2 * u - v) >= (3 * n2));
0344           break;
0345         }
0346         case (HGCalTypes::WaferCorner3): {
0347           good = ((2 * v - u) >= n3);
0348           break;
0349         }
0350         case (HGCalTypes::WaferCorner4): {
0351           good = ((u + v) > (4 * n3));
0352           break;
0353         }
0354         default: {
0355           good = ((2 * u - v) >= n2);
0356           break;
0357         }
0358       }
0359       break;
0360     }
0361     case (HGCalTypes::WaferHalf2): {  //WaferHalf2
0362       switch (rotn) {
0363         case (HGCalTypes::WaferCorner0): {
0364           good = (v <= (3 * n4));
0365           break;
0366         }
0367         case (HGCalTypes::WaferCorner1): {
0368           good = (u < (3 * n4));
0369           break;
0370         }
0371         case (HGCalTypes::WaferCorner2): {
0372           good = ((v - u) >= n4);
0373           break;
0374         }
0375         case (HGCalTypes::WaferCorner3): {
0376           good = (v >= (5 * n4));
0377           break;
0378         }
0379         case (HGCalTypes::WaferCorner4): {
0380           good = (u >= (5 * n4));
0381           break;
0382         }
0383         default: {
0384           good = ((u - v) >= n4);
0385           break;
0386         }
0387       }
0388       break;
0389     }
0390   }
0391 #ifdef EDM_ML_DEBUG
0392   edm::LogVerbatim("HGCalGeom") << "u|v " << u << ":" << v << " N " << n << " type " << type << " rot " << rotn
0393                                 << " good " << good;
0394 #endif
0395   return good;
0396 }
0397 
0398 bool HGCalWaferMask::goodCell(int u, int v, int waferType) {
0399   // for V17
0400   bool good(false);
0401   switch (waferType) {
0402     case (HGCalTypes::WaferFull): {  //WaferFull
0403       good = true;
0404       break;
0405     }
0406     case (HGCalTypes::WaferLDTop): {
0407       good = (u * HGCalTypes::edgeWaferLDTop[0] + v * HGCalTypes::edgeWaferLDTop[1] <= HGCalTypes::edgeWaferLDTop[2]);
0408       break;
0409     }
0410     case (HGCalTypes::WaferLDBottom): {
0411       good = (u * HGCalTypes::edgeWaferLDBottom[0] + v * HGCalTypes::edgeWaferLDBottom[1] <=
0412               HGCalTypes::edgeWaferLDBottom[2]);
0413       break;
0414     }
0415     case (HGCalTypes::WaferLDLeft): {
0416       good =
0417           (u * HGCalTypes::edgeWaferLDLeft[0] + v * HGCalTypes::edgeWaferLDLeft[1] <= HGCalTypes::edgeWaferLDLeft[2]);
0418       break;
0419     }
0420     case (HGCalTypes::WaferLDRight): {
0421       good = (u * HGCalTypes::edgeWaferLDRight[0] + v * HGCalTypes::edgeWaferLDRight[1] <=
0422               HGCalTypes::edgeWaferLDRight[2]);
0423       break;
0424     }
0425     case (HGCalTypes::WaferLDFive): {
0426       good =
0427           (u * HGCalTypes::edgeWaferLDFive[0] + v * HGCalTypes::edgeWaferLDFive[1] <= HGCalTypes::edgeWaferLDFive[2]);
0428       break;
0429     }
0430     case (HGCalTypes::WaferLDThree): {
0431       good = (u * HGCalTypes::edgeWaferLDThree[0] + v * HGCalTypes::edgeWaferLDThree[1] <=
0432               HGCalTypes::edgeWaferLDThree[2]);
0433       break;
0434     }
0435     case (HGCalTypes::WaferHDTop): {
0436       good = (u * HGCalTypes::edgeWaferHDTop[0] + v * HGCalTypes::edgeWaferHDTop[1] <= HGCalTypes::edgeWaferHDTop[2]);
0437       break;
0438     }
0439     case (HGCalTypes::WaferHDBottom): {
0440       good = (u * HGCalTypes::edgeWaferHDBottom[0] + v * HGCalTypes::edgeWaferHDBottom[1] <=
0441               HGCalTypes::edgeWaferHDBottom[2]);
0442       break;
0443     }
0444     case (HGCalTypes::WaferHDLeft): {
0445       good =
0446           (u * HGCalTypes::edgeWaferHDLeft[0] + v * HGCalTypes::edgeWaferHDLeft[1] <= HGCalTypes::edgeWaferHDLeft[2]);
0447       break;
0448     }
0449     case (HGCalTypes::WaferHDRight): {
0450       good = (u * HGCalTypes::edgeWaferHDRight[0] + v * HGCalTypes::edgeWaferHDRight[1] <=
0451               HGCalTypes::edgeWaferHDRight[2]);
0452       break;
0453     }
0454     case (HGCalTypes::WaferHDFive): {
0455       good =
0456           (u * HGCalTypes::edgeWaferHDFive[0] + v * HGCalTypes::edgeWaferHDFive[1] <= HGCalTypes::edgeWaferHDFive[2]);
0457       break;
0458     }
0459   }
0460 #ifdef EDM_ML_DEBUG
0461   edm::LogVerbatim("HGCalGeom") << "u|v " << u << ":" << v << " WaferType " << waferType << " good " << good;
0462 #endif
0463   return good;
0464 }
0465 
0466 int HGCalWaferMask::getRotation(int zside, int type, int rotn) {
0467   // Needs extension for V17
0468   if (rotn >= HGCalTypes::WaferCornerMax)
0469     rotn = HGCalTypes::WaferCorner0;
0470   int newrotn(rotn);
0471   if ((zside < 0) && (type != HGCalTypes::WaferFull)) {
0472     if ((type == HGCalTypes::WaferFive) || (type == HGCalTypes::WaferFive2)) {  //WaferFive/WaferFive2
0473       static constexpr int rot1[HGCalTypes::WaferCornerMax] = {HGCalTypes::WaferCorner4,
0474                                                                HGCalTypes::WaferCorner3,
0475                                                                HGCalTypes::WaferCorner2,
0476                                                                HGCalTypes::WaferCorner1,
0477                                                                HGCalTypes::WaferCorner0,
0478                                                                HGCalTypes::WaferCorner5};
0479       newrotn = rot1[rotn];
0480     } else if ((type == HGCalTypes::WaferThree) || (type == HGCalTypes::WaferSemi) ||
0481                (type == HGCalTypes::WaferSemi2)) {  //WaferThree/WaferSemi/WaferSemi2
0482       static constexpr int rot2[HGCalTypes::WaferCornerMax] = {HGCalTypes::WaferCorner2,
0483                                                                HGCalTypes::WaferCorner1,
0484                                                                HGCalTypes::WaferCorner0,
0485                                                                HGCalTypes::WaferCorner5,
0486                                                                HGCalTypes::WaferCorner4,
0487                                                                HGCalTypes::WaferCorner3};
0488       newrotn = rot2[rotn];
0489     } else {  //WaferHalf/WaferChopTwo/WaferChopTwoM/WaferHalf2
0490       static constexpr int rot3[HGCalTypes::WaferCornerMax] = {HGCalTypes::WaferCorner3,
0491                                                                HGCalTypes::WaferCorner2,
0492                                                                HGCalTypes::WaferCorner1,
0493                                                                HGCalTypes::WaferCorner0,
0494                                                                HGCalTypes::WaferCorner5,
0495                                                                HGCalTypes::WaferCorner4};
0496       newrotn = rot3[rotn];
0497     }
0498   }
0499 #ifdef EDM_ML_DEBUG
0500   edm::LogVerbatim("HGCalGeom") << "zside " << zside << " type " << type << " rotn " << rotn << ":" << newrotn;
0501 #endif
0502   return newrotn;
0503 }
0504 
0505 std::pair<int, int> HGCalWaferMask::getTypeMode(const double& xpos,
0506                                                 const double& ypos,
0507                                                 const double& delX,
0508                                                 const double& delY,
0509                                                 const double& rin,
0510                                                 const double& rout,
0511                                                 const int& wType,
0512                                                 const int& mode,
0513                                                 bool debug) {
0514   // No need to extend this for V17 -- use flat file information only
0515   int ncor(0), iok(0);
0516   int type(HGCalTypes::WaferFull), rotn(HGCalTypes::WaferCorner0);
0517 
0518   static constexpr int corners = 6;
0519   static constexpr int base = 10;
0520   double rin2 = rin * rin;
0521   double rout2 = rout * rout;
0522   double dx0[corners] = {HGCalTypes::c00 * delX,
0523                          HGCalTypes::c10 * delX,
0524                          HGCalTypes::c10 * delX,
0525                          HGCalTypes::c00 * delX,
0526                          -HGCalTypes::c10 * delX,
0527                          -HGCalTypes::c10 * delX};
0528   double dy0[corners] = {-HGCalTypes::c10 * delY,
0529                          -HGCalTypes::c50 * delY,
0530                          HGCalTypes::c50 * delY,
0531                          HGCalTypes::c10 * delY,
0532                          HGCalTypes::c50 * delY,
0533                          -HGCalTypes::c50 * delY};
0534   double xc[corners], yc[corners];
0535   for (int k = 0; k < corners; ++k) {
0536     xc[k] = xpos + dx0[k];
0537     yc[k] = ypos + dy0[k];
0538     double rpos2 = (xc[k] * xc[k] + yc[k] * yc[k]);
0539     if (rpos2 <= rout2 && rpos2 >= rin2) {
0540       ++ncor;
0541       iok = iok * base + 1;
0542     } else {
0543       iok *= base;
0544     }
0545   }
0546   if (debug)
0547     edm::LogVerbatim("HGCalGeom") << "I/p " << xpos << ":" << ypos << ":" << delX << ":" << delY << ":" << rin << ":"
0548                                   << rout << ":" << wType << ":" << mode << " Corners " << ncor << " iok " << iok;
0549 
0550   static constexpr int ipat5[corners] = {101111, 110111, 111011, 111101, 111110, 11111};
0551   static constexpr int ipat4[corners] = {100111, 110011, 111001, 111100, 11110, 1111};
0552   static constexpr int ipat3[corners] = {100011, 110001, 111000, 11100, 1110, 111};
0553   static constexpr int ipat2[corners] = {11, 100001, 110000, 11000, 1100, 110};
0554   double dx1[corners] = {HGCalTypes::c50 * delX,
0555                          HGCalTypes::c10 * delX,
0556                          HGCalTypes::c50 * delX,
0557                          -HGCalTypes::c50 * delX,
0558                          -HGCalTypes::c10 * delX,
0559                          -HGCalTypes::c50 * delX};
0560   double dy1[corners] = {-HGCalTypes::c75 * delY,
0561                          HGCalTypes::c00 * delY,
0562                          HGCalTypes::c75 * delY,
0563                          HGCalTypes::c75 * delY,
0564                          HGCalTypes::c00 * delY,
0565                          -HGCalTypes::c75 * delY};
0566   double dx2[corners] = {HGCalTypes::c50 * delX,
0567                          -HGCalTypes::c50 * delX,
0568                          -HGCalTypes::c10 * delX,
0569                          -HGCalTypes::c50 * delX,
0570                          HGCalTypes::c50 * delX,
0571                          HGCalTypes::c10 * delX};
0572   double dy2[corners] = {HGCalTypes::c75 * delY,
0573                          HGCalTypes::c75 * delY,
0574                          HGCalTypes::c00 * delY,
0575                          -HGCalTypes::c75 * delY,
0576                          -HGCalTypes::c75 * delY,
0577                          HGCalTypes::c00 * delY};
0578   double dx3[corners] = {HGCalTypes::c22 * delX,
0579                          HGCalTypes::c10 * delX,
0580                          HGCalTypes::c77 * delX,
0581                          -HGCalTypes::c22 * delX,
0582                          -HGCalTypes::c10 * delX,
0583                          -HGCalTypes::c77 * delX};
0584   double dy3[corners] = {-HGCalTypes::c88 * delY,
0585                          -HGCalTypes::c27 * delY,
0586                          HGCalTypes::c61 * delY,
0587                          HGCalTypes::c88 * delY,
0588                          HGCalTypes::c27 * delY,
0589                          -HGCalTypes::c61 * delY};
0590   double dx4[corners] = {HGCalTypes::c22 * delX,
0591                          -HGCalTypes::c77 * delX,
0592                          -HGCalTypes::c10 * delX,
0593                          -HGCalTypes::c22 * delX,
0594                          HGCalTypes::c77 * delX,
0595                          HGCalTypes::c10 * delX};
0596   double dy4[corners] = {HGCalTypes::c88 * delY,
0597                          HGCalTypes::c61 * delY,
0598                          -HGCalTypes::c27 * delY,
0599                          -HGCalTypes::c88 * delY,
0600                          -HGCalTypes::c61 * delY,
0601                          HGCalTypes::c27 * delY};
0602   double dx5[corners] = {-HGCalTypes::c50 * delX,
0603                          -HGCalTypes::c10 * delX,
0604                          -HGCalTypes::c50 * delX,
0605                          HGCalTypes::c50 * delX,
0606                          HGCalTypes::c10 * delX,
0607                          HGCalTypes::c50 * delX};
0608   double dy5[corners] = {HGCalTypes::c75 * delY,
0609                          HGCalTypes::c00 * delY,
0610                          -HGCalTypes::c75 * delY,
0611                          -HGCalTypes::c75 * delY,
0612                          HGCalTypes::c00 * delY,
0613                          HGCalTypes::c75 * delY};
0614   double dx6[corners] = {-HGCalTypes::c77 * delX,
0615                          -HGCalTypes::c10 * delX,
0616                          -HGCalTypes::c22 * delX,
0617                          HGCalTypes::c77 * delX,
0618                          HGCalTypes::c10 * delX,
0619                          HGCalTypes::c22 * delX};
0620   double dy6[corners] = {HGCalTypes::c61 * delY,
0621                          -HGCalTypes::c27 * delY,
0622                          -HGCalTypes::c88 * delY,
0623                          -HGCalTypes::c61 * delY,
0624                          HGCalTypes::c27 * delY,
0625                          HGCalTypes::c88 * delY};
0626   double dx7[corners] = {-HGCalTypes::c22 * delX,
0627                          -HGCalTypes::c10 * delX,
0628                          -HGCalTypes::c77 * delX,
0629                          HGCalTypes::c22 * delX,
0630                          HGCalTypes::c10 * delX,
0631                          HGCalTypes::c77 * delX};
0632   double dy7[corners] = {HGCalTypes::c88 * delY,
0633                          HGCalTypes::c27 * delY,
0634                          -HGCalTypes::c61 * delY,
0635                          -HGCalTypes::c88 * delY,
0636                          -HGCalTypes::c27 * delY,
0637                          HGCalTypes::c61 * delY};
0638   double dx8[corners] = {HGCalTypes::c77 * delX,
0639                          HGCalTypes::c10 * delX,
0640                          HGCalTypes::c22 * delX,
0641                          -HGCalTypes::c77 * delX,
0642                          -HGCalTypes::c10 * delX,
0643                          -HGCalTypes::c22 * delX};
0644   double dy8[corners] = {-HGCalTypes::c61 * delY,
0645                          HGCalTypes::c27 * delY,
0646                          HGCalTypes::c88 * delY,
0647                          HGCalTypes::c61 * delY,
0648                          -HGCalTypes::c27 * delY,
0649                          -HGCalTypes::c88 * delY};
0650   double dx9[corners] = {-HGCalTypes::c22 * delX,
0651                          HGCalTypes::c77 * delX,
0652                          HGCalTypes::c10 * delX,
0653                          HGCalTypes::c22 * delX,
0654                          -HGCalTypes::c77 * delX,
0655                          -HGCalTypes::c10 * delX};
0656   double dy9[corners] = {-HGCalTypes::c88 * delY,
0657                          -HGCalTypes::c61 * delY,
0658                          HGCalTypes::c27 * delY,
0659                          HGCalTypes::c88 * delY,
0660                          HGCalTypes::c61 * delY,
0661                          -HGCalTypes::c27 * delY};
0662 
0663   if (ncor == HGCalGeomTools::k_allCorners) {
0664   } else if (ncor == HGCalGeomTools::k_fiveCorners) {
0665     rotn = static_cast<int>(std::find(ipat5, ipat5 + 6, iok) - ipat5);
0666     type = HGCalTypes::WaferFive;
0667   } else if (ncor == HGCalGeomTools::k_fourCorners) {
0668     rotn = static_cast<int>(std::find(ipat4, ipat4 + 6, iok) - ipat4);
0669     type = HGCalTypes::WaferHalf;
0670     double rpos12 = ((xpos + dx1[rotn]) * (xpos + dx1[rotn]) + (ypos + dy1[rotn]) * (ypos + dy1[rotn]));
0671     double rpos22(0);
0672     if (rpos12 <= rout2 && rpos12 >= rin2) {
0673       rpos22 = ((xpos + dx2[rotn]) * (xpos + dx2[rotn]) + (ypos + dy2[rotn]) * (ypos + dy2[rotn]));
0674       if (rpos22 <= rout2 && rpos22 >= rin2)
0675         type = HGCalTypes::WaferChopTwo;
0676     }
0677     if (debug)
0678       edm::LogVerbatim("HGCalGeom") << "Test for Chop2 " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0679                                     << type;
0680     if ((type == HGCalTypes::WaferHalf) && (wType == 0)) {
0681       rpos12 = ((xpos + dx3[rotn]) * (xpos + dx3[rotn]) + (ypos + dy3[rotn]) * (ypos + dy3[rotn]));
0682       if (rpos12 <= rout2 && rpos12 >= rin2) {
0683         rpos22 = ((xpos + dx4[rotn]) * (xpos + dx4[rotn]) + (ypos + dy4[rotn]) * (ypos + dy4[rotn]));
0684         if (rpos22 <= rout2 && rpos22 >= rin2)
0685           type = HGCalTypes::WaferChopTwoM;
0686       }
0687       if (debug)
0688         edm::LogVerbatim("HGCalGeom") << "Test for Chop2M " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0689                                       << type;
0690     }
0691   } else if (ncor == HGCalGeomTools::k_threeCorners) {
0692     rotn = static_cast<int>(std::find(ipat3, ipat3 + 6, iok) - ipat3);
0693     type = HGCalTypes::WaferThree;
0694     double rpos12 = ((xpos + dx7[rotn]) * (xpos + dx7[rotn]) + (ypos + dy7[rotn]) * (ypos + dy7[rotn]));
0695     double rpos22(0);
0696     if (rpos12 <= rout2 && rpos12 >= rin2) {
0697       rpos22 = ((xpos + dx8[rotn]) * (xpos + dx8[rotn]) + (ypos + dy8[rotn]) * (ypos + dy8[rotn]));
0698       if (rpos22 <= rout2 && rpos22 >= rin2)
0699         type = HGCalTypes::WaferFive2;
0700     }
0701     if (debug)
0702       edm::LogVerbatim("HGCalGeom") << "Test for Five2 " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0703                                     << type;
0704     if ((type == HGCalTypes::WaferThree) && (wType == 0)) {
0705       rpos12 = ((xpos + dx1[rotn]) * (xpos + dx1[rotn]) + (ypos + dy1[rotn]) * (ypos + dy1[rotn]));
0706       if (rpos12 <= rout2 && rpos12 >= rin2) {
0707         rpos22 = ((xpos + dx5[rotn]) * (xpos + dx5[rotn]) + (ypos + dy5[rotn]) * (ypos + dy5[rotn]));
0708         if (rpos22 <= rout2 && rpos22 >= rin2)
0709           type = HGCalTypes::WaferSemi;
0710       }
0711       if (debug)
0712         edm::LogVerbatim("HGCalGeom") << "Test for Semi " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0713                                       << type;
0714     }
0715     if ((type == HGCalTypes::WaferThree) && (wType == 0)) {
0716       rpos12 = ((xpos + dx3[rotn]) * (xpos + dx3[rotn]) + (ypos + dy3[rotn]) * (ypos + dy3[rotn]));
0717       if (rpos12 <= rout2 && rpos12 >= rin2) {
0718         rpos22 = ((xpos + dx6[rotn]) * (xpos + dx6[rotn]) + (ypos + dy6[rotn]) * (ypos + dy6[rotn]));
0719         if (rpos22 <= rout2 && rpos22 >= rin2)
0720           type = HGCalTypes::WaferSemi2;
0721       }
0722       if (debug)
0723         edm::LogVerbatim("HGCalGeom") << "Test for SemiM " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0724                                       << type;
0725     }
0726   } else if (ncor == HGCalGeomTools::k_twoCorners) {
0727     rotn = static_cast<int>(std::find(ipat2, ipat2 + 6, iok) - ipat2);
0728     type = HGCalTypes::WaferOut;
0729     double rpos12 = ((xpos + dx7[rotn]) * (xpos + dx7[rotn]) + (ypos + dy7[rotn]) * (ypos + dy7[rotn]));
0730     double rpos22(0);
0731     if (rpos12 <= rout2 && rpos12 >= rin2) {
0732       rpos22 = ((xpos + dx9[rotn]) * (xpos + dx9[rotn]) + (ypos + dy9[rotn]) * (ypos + dy9[rotn]));
0733       if (rpos22 <= rout2 && rpos22 >= rin2)
0734         type = HGCalTypes::WaferHalf2;
0735       else
0736         rotn = HGCalTypes::WaferCorner0;
0737     }
0738     if (debug)
0739       edm::LogVerbatim("HGCalGeom") << "Test for Half2 " << std::sqrt(rpos12) << ":" << std::sqrt(rpos22) << " Type "
0740                                     << type;
0741   } else {
0742     type = HGCalTypes::WaferOut;
0743   }
0744 
0745   if (debug)
0746     edm::LogVerbatim("HGCalGeom") << "I/p " << xpos << ":" << ypos << ":" << delX << ":" << delY << ":" << rin << ":"
0747                                   << rout << ":" << wType << ":" << mode << " o/p " << iok << ":" << ncor << ":" << type
0748                                   << ":" << rotn;
0749   return ((mode == 0) ? std::make_pair(ncor, rotn) : std::make_pair(type, (rotn + HGCalTypes::k_OffsetRotation)));
0750 }
0751 
0752 bool HGCalWaferMask::goodTypeMode(
0753     double xpos, double ypos, double delX, double delY, double rin, double rout, int part, int rotn, bool debug) {
0754   // Needs extension for V17
0755   if (part < 0 || part > HGCalTypes::WaferSizeMax)
0756     return false;
0757   if (rotn < 0 || rotn > HGCalTypes::WaferCornerMax)
0758     return false;
0759   double rin2 = rin * rin;
0760   double rout2 = rout * rout;
0761   double rpos2(0);
0762   static constexpr int corners = HGCalTypes::WaferCornerMax;
0763   static constexpr int corner2 = 2 * HGCalTypes::WaferCornerMax;
0764   static constexpr int base = 10;
0765   static constexpr int base2 = 100;
0766   double dx0[corners] = {HGCalTypes::c00 * delX,
0767                          HGCalTypes::c10 * delX,
0768                          HGCalTypes::c10 * delX,
0769                          HGCalTypes::c00 * delX,
0770                          -HGCalTypes::c10 * delX,
0771                          -HGCalTypes::c10 * delX};
0772   double dy0[corners] = {-HGCalTypes::c10 * delY,
0773                          -HGCalTypes::c50 * delY,
0774                          HGCalTypes::c50 * delY,
0775                          HGCalTypes::c10 * delY,
0776                          HGCalTypes::c50 * delY,
0777                          -HGCalTypes::c50 * delY};
0778   double dx1[corners] = {HGCalTypes::c50 * delX,
0779                          HGCalTypes::c10 * delX,
0780                          HGCalTypes::c50 * delX,
0781                          -HGCalTypes::c50 * delX,
0782                          -HGCalTypes::c10 * delX,
0783                          -HGCalTypes::c50 * delX};
0784   double dy1[corners] = {-HGCalTypes::c75 * delY,
0785                          HGCalTypes::c00 * delY,
0786                          HGCalTypes::c75 * delY,
0787                          HGCalTypes::c75 * delY,
0788                          HGCalTypes::c00 * delY,
0789                          -HGCalTypes::c75 * delY};
0790   double dx2[corner2] = {HGCalTypes::c22 * delX,
0791                          HGCalTypes::c77 * delX,
0792                          HGCalTypes::c10 * delX,
0793                          HGCalTypes::c10 * delX,
0794                          HGCalTypes::c77 * delX,
0795                          HGCalTypes::c22 * delX,
0796                          -HGCalTypes::c22 * delX,
0797                          -HGCalTypes::c77 * delX,
0798                          -HGCalTypes::c10 * delX,
0799                          -HGCalTypes::c10 * delX,
0800                          -HGCalTypes::c77 * delX,
0801                          -HGCalTypes::c22 * delX};
0802   double dy2[corner2] = {-HGCalTypes::c88 * delY,
0803                          -HGCalTypes::c61 * delY,
0804                          -HGCalTypes::c27 * delY,
0805                          HGCalTypes::c27 * delY,
0806                          HGCalTypes::c61 * delY,
0807                          HGCalTypes::c88 * delY,
0808                          HGCalTypes::c88 * delY,
0809                          HGCalTypes::c61 * delY,
0810                          HGCalTypes::c27 * delY,
0811                          -HGCalTypes::c27 * delY,
0812                          -HGCalTypes::c61 * delY,
0813                          -HGCalTypes::c88 * delY};
0814   bool ok(true);
0815   int ncf(-1);
0816   switch (part) {
0817     case (HGCalTypes::WaferThree): {
0818       static constexpr int nc0[corners] = {450, 150, 201, 312, 423, 534};
0819       int nc = nc0[rotn];
0820       for (int k1 = 0; k1 < 3; ++k1) {
0821         int k = nc % base;
0822         double xc1 = xpos + dx0[k];
0823         double yc1 = ypos + dy0[k];
0824         rpos2 = (xc1 * xc1 + yc1 * yc1);
0825         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0826           ok = false;
0827           ncf = k;
0828           break;
0829         }
0830         nc /= base;
0831       }
0832       break;
0833     }
0834     case (HGCalTypes::WaferSemi2): {
0835       static constexpr int nc10[corners] = {450, 150, 201, 312, 423, 534};
0836       static constexpr int nc11[corners] = {700, 902, 1104, 106, 308, 510};
0837       int nc = nc10[rotn];
0838       for (int k1 = 0; k1 < 3; ++k1) {
0839         int k = nc % base;
0840         double xc1 = xpos + dx0[k];
0841         double yc1 = ypos + dy0[k];
0842         rpos2 = (xc1 * xc1 + yc1 * yc1);
0843         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0844           ok = false;
0845           ncf = k;
0846           break;
0847         }
0848         nc /= base;
0849       }
0850       nc = nc11[rotn];
0851       for (int k1 = 0; k1 < 2; ++k1) {
0852         int k = nc % base2;
0853         double xc1 = xpos + dx2[k];
0854         double yc1 = ypos + dy2[k];
0855         rpos2 = (xc1 * xc1 + yc1 * yc1);
0856         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0857           ok = false;
0858           ncf = k + base2;
0859           break;
0860         }
0861         nc /= base2;
0862       }
0863       break;
0864     }
0865     case (HGCalTypes::WaferSemi): {
0866       static constexpr int nc20[corners] = {450, 150, 201, 312, 423, 534};
0867       static constexpr int nc21[corners] = {30, 14, 25, 30, 41, 52};
0868       int nc = nc20[rotn];
0869       for (int k1 = 0; k1 < 3; ++k1) {
0870         int k = nc % base;
0871         double xc1 = xpos + dx0[k];
0872         double yc1 = ypos + dy0[k];
0873         rpos2 = (xc1 * xc1 + yc1 * yc1);
0874         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0875           ok = false;
0876           ncf = k;
0877           break;
0878         }
0879         nc /= base;
0880       }
0881       nc = nc21[rotn];
0882       for (int k1 = 0; k1 < 2; ++k1) {
0883         int k = nc % base;
0884         double xc1 = xpos + dx1[k];
0885         double yc1 = ypos + dy1[k];
0886         rpos2 = (xc1 * xc1 + yc1 * yc1);
0887         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0888           ok = false;
0889           ncf = k + base2;
0890           break;
0891         }
0892         nc /= base;
0893       }
0894       break;
0895     }
0896     case (HGCalTypes::WaferHalf): {
0897       static constexpr int nc3[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0898       int nc = nc3[rotn];
0899       for (int k1 = 0; k1 < 4; ++k1) {
0900         int k = nc % base;
0901         double xc1 = xpos + dx0[k];
0902         double yc1 = ypos + dy0[k];
0903         rpos2 = (xc1 * xc1 + yc1 * yc1);
0904         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0905           ok = false;
0906           ncf = k;
0907           break;
0908         }
0909         nc /= base;
0910       }
0911       break;
0912     }
0913     case (HGCalTypes::WaferChopTwoM): {
0914       static constexpr int nc40[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0915       static constexpr int nc41[corners] = {500, 702, 904, 1106, 108, 310};
0916       int nc = nc40[rotn];
0917       for (int k1 = 0; k1 < 4; ++k1) {
0918         int k = nc % base;
0919         double xc1 = xpos + dx0[k];
0920         double yc1 = ypos + dy0[k];
0921         rpos2 = (xc1 * xc1 + yc1 * yc1);
0922         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0923           ok = false;
0924           ncf = k;
0925           break;
0926         }
0927         nc /= base;
0928       }
0929       nc = nc41[rotn];
0930       for (int k1 = 0; k1 < 2; ++k1) {
0931         int k = nc % base2;
0932         double xc1 = xpos + dx2[k];
0933         double yc1 = ypos + dy2[k];
0934         rpos2 = (xc1 * xc1 + yc1 * yc1);
0935         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0936           ok = false;
0937           ncf = k + base2;
0938           break;
0939         }
0940         nc /= base2;
0941       }
0942       break;
0943     }
0944     case (HGCalTypes::WaferChopTwo): {
0945       static constexpr int nc50[corners] = {3450, 1450, 2501, 3012, 4123, 5234};
0946       static constexpr int nc51[corners] = {20, 13, 24, 35, 40, 51};
0947       int nc = nc50[rotn];
0948       for (int k1 = 0; k1 < 4; ++k1) {
0949         int k = nc % base;
0950         double xc1 = xpos + dx0[k];
0951         double yc1 = ypos + dy0[k];
0952         rpos2 = (xc1 * xc1 + yc1 * yc1);
0953         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0954           ok = false;
0955           ncf = k;
0956           break;
0957         }
0958         nc /= base;
0959       }
0960       nc = nc51[rotn];
0961       for (int k1 = 0; k1 < 2; ++k1) {
0962         int k = nc % base;
0963         double xc1 = xpos + dx1[k];
0964         double yc1 = ypos + dy1[k];
0965         rpos2 = (xc1 * xc1 + yc1 * yc1);
0966         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0967           ok = false;
0968           ncf = k + base2;
0969           break;
0970         }
0971         nc /= base;
0972       }
0973       break;
0974     }
0975     case (HGCalTypes::WaferFive): {
0976       static constexpr int nc6[corners] = {23450, 13450, 24501, 35012, 40123, 51234};
0977       int nc = nc6[rotn];
0978       for (int k1 = 0; k1 < 5; ++k1) {
0979         int k = nc % base;
0980         double xc1 = xpos + dx0[k];
0981         double yc1 = ypos + dy0[k];
0982         rpos2 = (xc1 * xc1 + yc1 * yc1);
0983         if ((rpos2 > rout2) || (rpos2 < rin2)) {
0984           ok = false;
0985           ncf = k;
0986           break;
0987         }
0988       }
0989       break;
0990     }
0991     case (HGCalTypes::WaferFive2): {
0992       static constexpr int nc60[corners] = {450, 150, 201, 312, 423, 534};
0993       static constexpr int nc61[corners] = {601, 803, 1005, 7, 209, 411};
0994       int nc = nc60[rotn];
0995       for (int k1 = 0; k1 < 3; ++k1) {
0996         int k = nc % base;
0997         double xc1 = xpos + dx0[k];
0998         double yc1 = ypos + dy0[k];
0999         rpos2 = (xc1 * xc1 + yc1 * yc1);
1000         if ((rpos2 > rout2) || (rpos2 < rin2)) {
1001           ok = false;
1002           ncf = k;
1003           break;
1004         }
1005         nc /= base;
1006       }
1007       nc = nc61[rotn];
1008       for (int k1 = 0; k1 < 2; ++k1) {
1009         int k = nc % base2;
1010         double xc1 = xpos + dx2[k];
1011         double yc1 = ypos + dy2[k];
1012         rpos2 = (xc1 * xc1 + yc1 * yc1);
1013         if ((rpos2 > rout2) || (rpos2 < rin2)) {
1014           ok = false;
1015           ncf = k + base2;
1016           break;
1017         }
1018         nc /= base2;
1019       }
1020       break;
1021     }
1022     case (HGCalTypes::WaferHalf2): {
1023       static constexpr int nc70[corners] = {45, 50, 1, 12, 23, 34};
1024       static constexpr int nc71[corners] = {611, 801, 1003, 5, 207, 409};
1025       int nc = nc70[rotn];
1026       for (int k1 = 0; k1 < 2; ++k1) {
1027         int k = nc % base;
1028         double xc1 = xpos + dx0[k];
1029         double yc1 = ypos + dy0[k];
1030         rpos2 = (xc1 * xc1 + yc1 * yc1);
1031         if ((rpos2 > rout2) || (rpos2 < rin2)) {
1032           ok = false;
1033           ncf = k;
1034           break;
1035         }
1036         nc /= base;
1037       }
1038       nc = nc71[rotn];
1039       for (int k1 = 0; k1 < 2; ++k1) {
1040         int k = nc % base2;
1041         double xc1 = xpos + dx2[k];
1042         double yc1 = ypos + dy2[k];
1043         rpos2 = (xc1 * xc1 + yc1 * yc1);
1044         if ((rpos2 > rout2) || (rpos2 < rin2)) {
1045           ok = false;
1046           ncf = k + base2;
1047           break;
1048         }
1049         nc /= base2;
1050       }
1051       break;
1052     }
1053     default: {
1054       for (int k = 0; k < corners; ++k) {
1055         double xc1 = xpos + dx0[k];
1056         double yc1 = ypos + dy0[k];
1057         rpos2 = (xc1 * xc1 + yc1 * yc1);
1058         if ((rpos2 > rout2) || (rpos2 < rin2)) {
1059           ok = false;
1060           ncf = k;
1061           break;
1062         }
1063       }
1064       break;
1065     }
1066   }
1067   if (debug || (!ok))
1068     edm::LogVerbatim("HGCalGeom") << "I/p "
1069                                   << ":" << xpos << ":" << ypos << ":" << delX << ":" << delY << ":" << rin << ":"
1070                                   << rout << ":" << part << ":" << rotn << " Results " << ok << ":" << ncf << " R "
1071                                   << rin2 << ":" << rout2 << ":" << rpos2;
1072   return ok;
1073 }
1074 
1075 std::vector<std::pair<double, double> > HGCalWaferMask::waferXY(
1076     int part, int ori, int zside, double delX, double delY, double xpos, double ypos) {
1077   // Good for V15 and V16 versions
1078   std::vector<std::pair<double, double> > xy;
1079   int orient = getRotation(-zside, part, ori);
1080 #ifdef EDM_ML_DEBUG
1081   edm::LogVerbatim("HGCalGeom") << "Part " << part << " zSide " << zside << " Orient " << ori << ":" << orient;
1082 #endif
1083   /*
1084     The exact contour of partial wafers are obtained by joining points on
1085     the circumference of a full wafer.
1086     Numbering the points along the edges of a hexagonal wafer, starting from
1087     the bottom corner:
1088 
1089                                    3
1090                                15     18
1091                              9           8
1092                           19               14
1093                         4                     2 
1094                        16                    23
1095                        10                     7
1096                        20                    13
1097                         5                     1
1098                           17               22
1099                             11           6
1100                                21     12
1101                                    0
1102 
1103     Depending on the wafer type and orientation index, the corners
1104     are chosen in the variable *np*
1105   */
1106   double dx[24] = {HGCalTypes::c00 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c00 * delX,
1107                    -HGCalTypes::c10 * delX, -HGCalTypes::c10 * delX, HGCalTypes::c50 * delX,  HGCalTypes::c10 * delX,
1108                    HGCalTypes::c50 * delX,  -HGCalTypes::c50 * delX, -HGCalTypes::c10 * delX, -HGCalTypes::c50 * delX,
1109                    HGCalTypes::c22 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c77 * delX,  -HGCalTypes::c22 * delX,
1110                    -HGCalTypes::c10 * delX, -HGCalTypes::c77 * delX, HGCalTypes::c22 * delX,  -HGCalTypes::c77 * delX,
1111                    -HGCalTypes::c10 * delX, -HGCalTypes::c22 * delX, HGCalTypes::c77 * delX,  HGCalTypes::c10 * delX};
1112   double dy[24] = {-HGCalTypes::c10 * delY, -HGCalTypes::c50 * delY, HGCalTypes::c50 * delY,  HGCalTypes::c10 * delY,
1113                    HGCalTypes::c50 * delY,  -HGCalTypes::c50 * delY, -HGCalTypes::c75 * delY, HGCalTypes::c00 * delY,
1114                    HGCalTypes::c75 * delY,  HGCalTypes::c75 * delY,  HGCalTypes::c00 * delY,  -HGCalTypes::c75 * delY,
1115                    -HGCalTypes::c88 * delY, -HGCalTypes::c27 * delY, HGCalTypes::c61 * delY,  HGCalTypes::c88 * delY,
1116                    HGCalTypes::c27 * delY,  -HGCalTypes::c61 * delY, HGCalTypes::c88 * delY,  HGCalTypes::c61 * delY,
1117                    -HGCalTypes::c27 * delY, -HGCalTypes::c88 * delY, -HGCalTypes::c61 * delY, HGCalTypes::c27 * delY};
1118   if (part == HGCalTypes::WaferFull) {
1119     int np[7] = {0, 1, 2, 3, 4, 5, 0};
1120     for (int k = 0; k < 7; ++k)
1121       xy.push_back(std::make_pair((xpos + dx[np[k]]), (ypos + dy[np[k]])));
1122   } else if (part == HGCalTypes::WaferFive) {
1123     int np[6][6] = {{0, 2, 3, 4, 5, 0},
1124                     {1, 3, 4, 5, 0, 1},
1125                     {2, 4, 5, 0, 1, 2},
1126                     {3, 5, 0, 1, 2, 3},
1127                     {4, 0, 1, 2, 3, 4},
1128                     {5, 1, 2, 3, 4, 5}};
1129     for (int k = 0; k < 6; ++k) {
1130       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1131 #ifdef EDM_ML_DEBUG
1132       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1133                                     << dy[np[orient][k]];
1134 #endif
1135     }
1136   } else if (part == HGCalTypes::WaferHalf) {
1137     int np[6][5] = {
1138         {0, 3, 4, 5, 0}, {1, 4, 5, 0, 1}, {2, 5, 0, 1, 2}, {3, 0, 1, 2, 3}, {4, 1, 2, 3, 4}, {5, 2, 3, 4, 5}};
1139     for (int k = 0; k < 5; ++k) {
1140       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1141 #ifdef EDM_ML_DEBUG
1142       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1143                                     << dy[np[orient][k]];
1144 #endif
1145     }
1146   } else if (part == HGCalTypes::WaferThree) {
1147     int np[6][4] = {{0, 4, 5, 0}, {1, 5, 0, 1}, {2, 0, 1, 2}, {3, 1, 2, 3}, {4, 2, 3, 4}, {5, 3, 4, 5}};
1148     for (int k = 0; k < 4; ++k) {
1149       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1150 #ifdef EDM_ML_DEBUG
1151       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1152                                     << dy[np[orient][k]];
1153 #endif
1154     }
1155   } else if (part == HGCalTypes::WaferChopTwo) {
1156     int np[6][7] = {{6, 8, 3, 4, 5, 0, 6},
1157                     {7, 9, 4, 5, 0, 1, 7},
1158                     {8, 10, 5, 0, 1, 2, 8},
1159                     {9, 11, 0, 1, 2, 3, 9},
1160                     {10, 6, 1, 2, 3, 4, 10},
1161                     {11, 7, 2, 3, 4, 5, 11}};
1162     for (int k = 0; k < 7; ++k) {
1163       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1164 #ifdef EDM_ML_DEBUG
1165       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1166                                     << dy[np[orient][k]];
1167 #endif
1168     }
1169   } else if (part == HGCalTypes::WaferSemi) {
1170     int np[6][6] = {{6, 9, 4, 5, 0, 6},
1171                     {7, 10, 5, 0, 1, 7},
1172                     {8, 11, 0, 1, 2, 8},
1173                     {9, 6, 1, 2, 3, 9},
1174                     {10, 7, 2, 3, 4, 10},
1175                     {11, 8, 3, 4, 5, 11}};
1176     for (int k = 0; k < 6; ++k) {
1177       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1178 #ifdef EDM_ML_DEBUG
1179       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1180                                     << dy[np[orient][k]];
1181 #endif
1182     }
1183   } else if (part == HGCalTypes::WaferChopTwoM) {
1184     int np[6][7] = {{12, 18, 3, 4, 5, 0, 12},
1185                     {13, 19, 4, 5, 0, 1, 13},
1186                     {14, 20, 5, 0, 1, 2, 14},
1187                     {15, 21, 0, 1, 2, 3, 15},
1188                     {16, 22, 1, 2, 3, 4, 16},
1189                     {17, 23, 2, 3, 4, 5, 17}};
1190     for (int k = 0; k < 7; ++k) {
1191       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1192 #ifdef EDM_ML_DEBUG
1193       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1194                                     << dy[np[orient][k]];
1195 #endif
1196     }
1197   } else if (part == HGCalTypes::WaferSemi2) {
1198     int np[6][6] = {{12, 19, 4, 5, 0, 12},
1199                     {13, 20, 5, 0, 1, 13},
1200                     {14, 21, 0, 1, 2, 14},
1201                     {15, 22, 1, 2, 3, 15},
1202                     {16, 23, 2, 3, 4, 16},
1203                     {17, 18, 3, 4, 5, 17}};
1204     for (int k = 0; k < 6; ++k) {
1205       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1206 #ifdef EDM_ML_DEBUG
1207       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1208                                     << dy[np[orient][k]];
1209 #endif
1210     }
1211   } else if (part == HGCalTypes::WaferFive2) {
1212     int np[6][6] = {{22, 15, 4, 5, 0, 22},
1213                     {23, 16, 5, 0, 1, 23},
1214                     {18, 17, 0, 1, 2, 18},
1215                     {19, 12, 1, 2, 3, 19},
1216                     {20, 13, 2, 3, 4, 20},
1217                     {21, 14, 3, 4, 5, 21}};
1218     for (int k = 0; k < 6; ++k) {
1219       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1220 #ifdef EDM_ML_DEBUG
1221       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1222                                     << dy[np[orient][k]];
1223 #endif
1224     }
1225   } else if (part == HGCalTypes::WaferHalf2) {
1226     int np[6][5] = {{21, 15, 4, 5, 21},
1227                     {22, 16, 5, 0, 22},
1228                     {23, 17, 0, 1, 23},
1229                     {18, 12, 1, 2, 18},
1230                     {19, 13, 2, 3, 19},
1231                     {20, 14, 3, 4, 20}};
1232     for (int k = 0; k < 5; ++k) {
1233       xy.push_back(std::make_pair((xpos + dx[np[orient][k]]), (ypos + dy[np[orient][k]])));
1234 #ifdef EDM_ML_DEBUG
1235       edm::LogVerbatim("HGCalGeom") << k << ":" << np[orient][k] << ":" << dx[np[orient][k]] << ":"
1236                                     << dy[np[orient][k]];
1237 #endif
1238     }
1239   }
1240 #ifdef EDM_ML_DEBUG
1241   edm::LogVerbatim("HGCalGeom") << "I/p: " << part << ":" << ori << ":" << zside << ":" << delX << ":" << delY << ":"
1242                                 << xpos << ":" << ypos << " O/p having " << xy.size() << " points:";
1243   std::ostringstream st1;
1244   for (unsigned int i = 0; i < xy.size(); ++i)
1245     st1 << " [" << i << "] " << xy[i].first << ":" << xy[i].second;
1246   edm::LogVerbatim("HGCalGeom") << st1.str();
1247 #endif
1248   return xy;
1249 }
1250 
1251 std::vector<std::pair<double, double> > HGCalWaferMask::waferXY(
1252     int part, int place, double delX, double delY, double xpos, double ypos) {
1253   std::vector<std::pair<double, double> > xy;
1254   // Good for V17 version and uses partial wafer type & placement index
1255 #ifdef EDM_ML_DEBUG
1256   edm::LogVerbatim("HGCalGeom") << "Part " << part << " Placement Index " << place;
1257 #endif
1258   /*
1259     The exact contour of partial wafers are obtained by joining points on
1260     the circumference of a full wafer.
1261     Numbering the points along the edges of a hexagonal wafer, starting from
1262     the bottom corner:
1263 
1264                                    3
1265                                15     18
1266                              9           8
1267                           19               14
1268                         4                     2 
1269                        16                    23
1270                        10                     7
1271                        20                    13
1272                         5                     1
1273                           17               22
1274                             11           6
1275                                21     12
1276                                    0
1277 
1278     Depending on the wafer type and placement index, the corners
1279     are chosen in the variable *np*
1280   */
1281   double dx[24] = {HGCalTypes::c00 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c00 * delX,
1282                    -HGCalTypes::c10 * delX, -HGCalTypes::c10 * delX, HGCalTypes::c50 * delX,  HGCalTypes::c10 * delX,
1283                    HGCalTypes::c50 * delX,  -HGCalTypes::c50 * delX, -HGCalTypes::c10 * delX, -HGCalTypes::c50 * delX,
1284                    HGCalTypes::c22 * delX,  HGCalTypes::c10 * delX,  HGCalTypes::c77 * delX,  -HGCalTypes::c22 * delX,
1285                    -HGCalTypes::c10 * delX, -HGCalTypes::c77 * delX, HGCalTypes::c22 * delX,  -HGCalTypes::c77 * delX,
1286                    -HGCalTypes::c10 * delX, -HGCalTypes::c22 * delX, HGCalTypes::c77 * delX,  HGCalTypes::c10 * delX};
1287   double dy[24] = {-HGCalTypes::c10 * delY, -HGCalTypes::c50 * delY, HGCalTypes::c50 * delY,  HGCalTypes::c10 * delY,
1288                    HGCalTypes::c50 * delY,  -HGCalTypes::c50 * delY, -HGCalTypes::c75 * delY, HGCalTypes::c00 * delY,
1289                    HGCalTypes::c75 * delY,  HGCalTypes::c75 * delY,  HGCalTypes::c00 * delY,  -HGCalTypes::c75 * delY,
1290                    -HGCalTypes::c88 * delY, -HGCalTypes::c27 * delY, HGCalTypes::c61 * delY,  HGCalTypes::c88 * delY,
1291                    HGCalTypes::c27 * delY,  -HGCalTypes::c61 * delY, HGCalTypes::c88 * delY,  HGCalTypes::c61 * delY,
1292                    -HGCalTypes::c27 * delY, -HGCalTypes::c88 * delY, -HGCalTypes::c61 * delY, HGCalTypes::c27 * delY};
1293   if (part == HGCalTypes::WaferFull) {
1294     int np[7] = {0, 1, 2, 3, 4, 5, 0};
1295     for (int k = 0; k < 7; ++k)
1296       xy.push_back(std::make_pair((xpos + dx[np[k]]), (ypos + dy[np[k]])));
1297   } else if (part == HGCalTypes::WaferLDTop) {
1298     int np[12][5] = {{0, 1, 4, 5, 0},
1299                      {1, 2, 5, 0, 1},
1300                      {2, 3, 0, 1, 2},
1301                      {3, 4, 1, 2, 3},
1302                      {4, 5, 2, 3, 4},
1303                      {5, 0, 3, 4, 5},
1304                      {0, 1, 2, 5, 0},
1305                      {5, 0, 1, 4, 5},
1306                      {4, 5, 0, 3, 4},
1307                      {3, 4, 5, 2, 3},
1308                      {2, 3, 4, 1, 2},
1309                      {1, 2, 3, 0, 1}};
1310     for (int k = 0; k < 5; ++k) {
1311       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1312 #ifdef EDM_ML_DEBUG
1313       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1314 #endif
1315     }
1316   } else if (part == HGCalTypes::WaferLDBottom) {
1317     int np[12][5] = {{1, 2, 3, 4, 1},
1318                      {2, 3, 4, 5, 2},
1319                      {3, 4, 5, 0, 3},
1320                      {4, 5, 0, 1, 4},
1321                      {5, 0, 1, 2, 5},
1322                      {0, 1, 2, 3, 0},
1323                      {5, 2, 3, 4, 5},
1324                      {4, 1, 2, 3, 4},
1325                      {3, 0, 1, 2, 3},
1326                      {2, 5, 0, 1, 2},
1327                      {1, 4, 5, 0, 1},
1328                      {0, 3, 4, 5, 0}};
1329     for (int k = 0; k < 5; ++k) {
1330       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1331 #ifdef EDM_ML_DEBUG
1332       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1333 #endif
1334     }
1335   } else if (part == HGCalTypes::WaferLDLeft) {
1336     int np[12][6] = {{0, 1, 2, 8, 11, 0},
1337                      {1, 2, 3, 9, 6, 1},
1338                      {2, 3, 4, 10, 7, 2},
1339                      {3, 4, 5, 11, 8, 3},
1340                      {4, 5, 0, 6, 9, 4},
1341                      {5, 0, 1, 7, 10, 5},
1342                      {0, 6, 9, 4, 5, 0},
1343                      {5, 11, 8, 3, 4, 5},
1344                      {4, 10, 7, 2, 3, 4},
1345                      {3, 9, 6, 1, 2, 3},
1346                      {2, 8, 11, 0, 1, 2},
1347                      {1, 7, 10, 5, 0, 1}};
1348     for (int k = 0; k < 6; ++k) {
1349       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1350 #ifdef EDM_ML_DEBUG
1351       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1352 #endif
1353     }
1354   } else if (part == HGCalTypes::WaferLDRight) {
1355     int np[12][6] = {{5, 11, 8, 3, 4, 5},
1356                      {0, 6, 9, 4, 5, 0},
1357                      {1, 7, 10, 5, 0, 1},
1358                      {2, 8, 11, 0, 1, 2},
1359                      {3, 9, 6, 1, 2, 3},
1360                      {4, 10, 7, 2, 3, 4},
1361                      {1, 2, 3, 9, 6, 1},
1362                      {0, 1, 2, 8, 11, 0},
1363                      {5, 0, 1, 7, 10, 5},
1364                      {4, 5, 0, 6, 9, 4},
1365                      {3, 4, 5, 11, 8, 3},
1366                      {2, 3, 4, 10, 7, 2}};
1367     for (int k = 0; k < 6; ++k) {
1368       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1369 #ifdef EDM_ML_DEBUG
1370       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1371 #endif
1372     }
1373   } else if (part == HGCalTypes::WaferLDFive) {
1374     int np[12][6] = {{0, 1, 2, 3, 5, 0},
1375                      {1, 2, 3, 4, 0, 1},
1376                      {2, 3, 4, 5, 1, 2},
1377                      {3, 4, 5, 0, 2, 3},
1378                      {4, 5, 0, 1, 3, 4},
1379                      {5, 0, 1, 2, 4, 5},
1380                      {0, 1, 3, 4, 5, 0},
1381                      {5, 0, 2, 3, 4, 5},
1382                      {4, 5, 1, 2, 3, 4},
1383                      {3, 4, 0, 1, 2, 3},
1384                      {2, 3, 5, 0, 1, 2},
1385                      {1, 2, 4, 5, 0, 1}};
1386     for (int k = 0; k < 6; ++k) {
1387       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1388 #ifdef EDM_ML_DEBUG
1389       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1390 #endif
1391     }
1392   } else if (part == HGCalTypes::WaferLDThree) {
1393     int np[12][4] = {{5, 3, 4, 5},
1394                      {0, 4, 5, 0},
1395                      {1, 5, 0, 1},
1396                      {2, 0, 1, 2},
1397                      {3, 1, 2, 3},
1398                      {4, 2, 3, 4},
1399                      {1, 2, 3, 1},
1400                      {0, 1, 2, 0},
1401                      {5, 0, 1, 5},
1402                      {4, 5, 0, 4},
1403                      {3, 4, 5, 3},
1404                      {2, 3, 4, 2}};
1405     for (int k = 0; k < 4; ++k) {
1406       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1407 #ifdef EDM_ML_DEBUG
1408       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1409 #endif
1410     }
1411   } else if (part == HGCalTypes::WaferHDTop) {
1412     int np[12][5] = {{0, 22, 16, 5, 0},
1413                      {1, 23, 17, 0, 1},
1414                      {2, 18, 12, 1, 2},
1415                      {3, 19, 13, 2, 3},
1416                      {4, 20, 14, 3, 4},
1417                      {5, 21, 15, 4, 5},
1418                      {0, 1, 23, 17, 0},
1419                      {5, 0, 22, 16, 5},
1420                      {4, 5, 21, 15, 4},
1421                      {3, 4, 20, 14, 3},
1422                      {2, 3, 19, 13, 2},
1423                      {1, 2, 18, 12, 1}};
1424     for (int k = 0; k < 5; ++k) {
1425       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1426 #ifdef EDM_ML_DEBUG
1427       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1428 #endif
1429     }
1430   } else if (part == HGCalTypes::WaferHDBottom) {
1431     int np[12][7] = {{1, 2, 3, 4, 16, 22, 1},
1432                      {2, 3, 4, 5, 17, 23, 2},
1433                      {3, 4, 5, 0, 12, 18, 3},
1434                      {4, 5, 0, 1, 13, 19, 4},
1435                      {5, 0, 1, 2, 14, 20, 5},
1436                      {0, 1, 2, 3, 15, 21, 0},
1437                      {5, 17, 23, 2, 3, 4, 5},
1438                      {4, 16, 22, 1, 2, 3, 4},
1439                      {3, 15, 21, 0, 1, 2, 3},
1440                      {2, 14, 20, 5, 0, 1, 2},
1441                      {1, 13, 19, 4, 5, 0, 1},
1442                      {0, 12, 18, 3, 4, 5, 0}};
1443     for (int k = 0; k < 7; ++k) {
1444       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1445 #ifdef EDM_ML_DEBUG
1446       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1447 #endif
1448     }
1449   } else if (part == HGCalTypes::WaferHDLeft) {
1450     int np[12][6] = {{0, 1, 2, 14, 21, 0},
1451                      {1, 2, 3, 15, 22, 1},
1452                      {2, 3, 4, 16, 23, 2},
1453                      {3, 4, 5, 17, 18, 3},
1454                      {4, 5, 0, 12, 19, 4},
1455                      {5, 0, 1, 13, 20, 5},
1456                      {0, 12, 19, 4, 5, 0},
1457                      {5, 17, 18, 3, 4, 5},
1458                      {4, 16, 23, 2, 3, 4},
1459                      {3, 15, 22, 1, 2, 3},
1460                      {2, 14, 21, 0, 1, 2},
1461                      {1, 13, 20, 5, 0, 1}};
1462     for (int k = 0; k < 6; ++k) {
1463       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1464 #ifdef EDM_ML_DEBUG
1465       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1466 #endif
1467     }
1468   } else if (part == HGCalTypes::WaferHDRight) {
1469     int np[12][6] = {{5, 17, 18, 3, 4, 5},
1470                      {0, 12, 19, 4, 5, 0},
1471                      {1, 13, 20, 5, 0, 1},
1472                      {2, 14, 21, 0, 1, 2},
1473                      {3, 15, 22, 1, 2, 3},
1474                      {4, 16, 23, 2, 3, 4},
1475                      {1, 2, 3, 15, 22, 1},
1476                      {0, 1, 2, 14, 21, 0},
1477                      {5, 0, 1, 13, 20, 5},
1478                      {4, 5, 0, 12, 19, 4},
1479                      {3, 4, 5, 17, 18, 3},
1480                      {2, 3, 4, 16, 23, 2}};
1481     for (int k = 0; k < 6; ++k) {
1482       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1483 #ifdef EDM_ML_DEBUG
1484       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1485 #endif
1486     }
1487   } else if (part == HGCalTypes::WaferHDFive) {
1488     int np[12][6] = {{0, 1, 2, 18, 17, 0},
1489                      {1, 2, 3, 19, 12, 1},
1490                      {2, 3, 4, 20, 13, 2},
1491                      {3, 4, 5, 21, 14, 3},
1492                      {4, 5, 0, 22, 15, 4},
1493                      {5, 0, 1, 23, 16, 5},
1494                      {0, 22, 15, 4, 5, 0},
1495                      {5, 21, 14, 3, 4, 5},
1496                      {4, 20, 13, 2, 3, 4},
1497                      {3, 19, 12, 1, 2, 3},
1498                      {2, 18, 17, 0, 1, 2},
1499                      {1, 23, 16, 5, 0, 1}};
1500     for (int k = 0; k < 6; ++k) {
1501       xy.push_back(std::make_pair((xpos + dx[np[place][k]]), (ypos + dy[np[place][k]])));
1502 #ifdef EDM_ML_DEBUG
1503       edm::LogVerbatim("HGCalGeom") << k << ":" << np[place][k] << ":" << dx[np[place][k]] << ":" << dy[np[place][k]];
1504 #endif
1505     }
1506   }
1507 #ifdef EDM_ML_DEBUG
1508   edm::LogVerbatim("HGCalGeom") << "I/p: " << part << ":" << place << ":" << delX << ":" << delY << ":" << xpos << ":"
1509                                 << ypos << " O/p having " << xy.size() << " points:";
1510   std::ostringstream st1;
1511   for (unsigned int i = 0; i < xy.size(); ++i)
1512     st1 << " [" << i << "] " << xy[i].first << ":" << xy[i].second;
1513   edm::LogVerbatim("HGCalGeom") << st1.str();
1514 #endif
1515   return xy;
1516 }