Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-01-25 02:56:25

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