Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-10-30 00:11:40 

0001 #include "Geometry/HGCalCommonData/interface/HGCalGeomParameters.h"
0002 
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "FWCore/Utilities/interface/Exception.h"
0005 
0006 #include "DataFormats/Math/interface/GeantUnits.h"
0007 #include "DataFormats/Math/interface/Point3D.h"
0008 #include "DetectorDescription/Core/interface/DDConstant.h"
0009 #include "DetectorDescription/Core/interface/DDFilter.h"
0010 #include "DetectorDescription/Core/interface/DDFilteredView.h"
0011 #include "DetectorDescription/Core/interface/DDSolid.h"
0012 #include "DetectorDescription/Core/interface/DDValue.h"
0013 #include "DetectorDescription/Core/interface/DDutils.h"
0014 #include "DetectorDescription/RegressionTest/interface/DDErrorDetection.h"
0015 #include "Geometry/HGCalCommonData/interface/HGCalProperty.h"
0016 #include "Geometry/HGCalCommonData/interface/HGCalTileIndex.h"
0017 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0018 #include "Geometry/HGCalCommonData/interface/HGCalWaferIndex.h"
0019 #include "Geometry/HGCalCommonData/interface/HGCalWaferMask.h"
0020 #include "Geometry/HGCalCommonData/interface/HGCalWaferType.h"
0021 
0022 #include <algorithm>
0023 #include <sstream>
0024 #include <unordered_set>
0025 
0026 //#define EDM_ML_DEBUG
0027 using namespace geant_units::operators;
0028 
0029 const double tolerance = 0.001;
0030 const double tolmin = 1.e-20;
0031 
0032 HGCalGeomParameters::HGCalGeomParameters() : sqrt3_(std::sqrt(3.0)) {
0033 #ifdef EDM_ML_DEBUG
0034   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters::HGCalGeomParameters "
0035                                 << "constructor";
0036 #endif
0037 }
0038 
0039 void HGCalGeomParameters::loadGeometryHexagon(const DDFilteredView& _fv,
0040                                               HGCalParameters& php,
0041                                               const std::string& sdTag1,
0042                                               const DDCompactView* cpv,
0043                                               const std::string& sdTag2,
0044                                               const std::string& sdTag3,
0045                                               HGCalGeometryMode::WaferMode mode) {
0046   DDFilteredView fv = _fv;
0047   bool dodet(true);
0048   std::map<int, HGCalGeomParameters::layerParameters> layers;
0049   std::vector<HGCalParameters::hgtrform> trforms;
0050   std::vector<bool> trformUse;
0051 
0052   while (dodet) {
0053     const DDSolid& sol = fv.logicalPart().solid();
0054     // Layers first
0055     std::vector<int> copy = fv.copyNumbers();
0056     int nsiz = static_cast<int>(copy.size());
0057     int lay = (nsiz > 0) ? copy[nsiz - 1] : 0;
0058     int zp = (nsiz > 2) ? copy[nsiz - 3] : -1;
0059     if (zp != 1)
0060       zp = -1;
0061     if (lay == 0) {
0062       throw cms::Exception("DDException") << "Funny layer # " << lay << " zp " << zp << " in " << nsiz << " components";
0063     } else {
0064       if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
0065         php.layer_.emplace_back(lay);
0066       auto itr = layers.find(lay);
0067       if (itr == layers.end()) {
0068         double rin(0), rout(0);
0069         double zz = HGCalParameters::k_ScaleFromDDD * fv.translation().Z();
0070         if ((sol.shape() == DDSolidShape::ddpolyhedra_rz) || (sol.shape() == DDSolidShape::ddpolyhedra_rrz)) {
0071           const DDPolyhedra& polyhedra = static_cast<DDPolyhedra>(sol);
0072           const std::vector<double>& rmin = polyhedra.rMinVec();
0073           const std::vector<double>& rmax = polyhedra.rMaxVec();
0074           rin = 0.5 * HGCalParameters::k_ScaleFromDDD * (rmin[0] + rmin[1]);
0075           rout = 0.5 * HGCalParameters::k_ScaleFromDDD * (rmax[0] + rmax[1]);
0076         } else if (sol.shape() == DDSolidShape::ddtubs) {
0077           const DDTubs& tube = static_cast<DDTubs>(sol);
0078           rin = HGCalParameters::k_ScaleFromDDD * tube.rIn();
0079           rout = HGCalParameters::k_ScaleFromDDD * tube.rOut();
0080         }
0081         HGCalGeomParameters::layerParameters laypar(rin, rout, zz);
0082         layers[lay] = laypar;
0083       }
0084       DD3Vector x, y, z;
0085       fv.rotation().GetComponents(x, y, z);
0086       const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0087       const CLHEP::HepRotation hr(rotation);
0088       double xx = HGCalParameters::k_ScaleFromDDD * fv.translation().X();
0089       if (std::abs(xx) < tolerance)
0090         xx = 0;
0091       double yy = HGCalParameters::k_ScaleFromDDD * fv.translation().Y();
0092       if (std::abs(yy) < tolerance)
0093         yy = 0;
0094       double zz = HGCalParameters::k_ScaleFromDDD * fv.translation().Z();
0095       const CLHEP::Hep3Vector h3v(xx, yy, zz);
0096       HGCalParameters::hgtrform mytrf;
0097       mytrf.zp = zp;
0098       mytrf.lay = lay;
0099       mytrf.sec = 0;
0100       mytrf.subsec = 0;
0101       mytrf.h3v = h3v;
0102       mytrf.hr = hr;
0103       trforms.emplace_back(mytrf);
0104       trformUse.emplace_back(false);
0105     }
0106     dodet = fv.next();
0107   }
0108 
0109   // Then wafers
0110   // This assumes layers are build starting from 1 (which on 25 Jan 2016, they
0111   // were) to ensure that new copy numbers are always added to the end of the
0112   // list.
0113   std::unordered_map<int32_t, int32_t> copies;
0114   HGCalParameters::layer_map copiesInLayers(layers.size() + 1);
0115   std::vector<int32_t> wafer2copy;
0116   std::vector<HGCalGeomParameters::cellParameters> wafers;
0117   std::string attribute = "Volume";
0118   DDValue val1(attribute, sdTag2, 0.0);
0119   DDSpecificsMatchesValueFilter filter1{val1};
0120   DDFilteredView fv1(*cpv, filter1);
0121   bool ok = fv1.firstChild();
0122   if (!ok) {
0123     throw cms::Exception("DDException") << "Attribute " << val1 << " not found but needed.";
0124   } else {
0125     dodet = true;
0126     std::unordered_set<std::string> names;
0127     while (dodet) {
0128       const DDSolid& sol = fv1.logicalPart().solid();
0129       const std::string& name = fv1.logicalPart().name().name();
0130       std::vector<int> copy = fv1.copyNumbers();
0131       int nsiz = static_cast<int>(copy.size());
0132       int wafer = (nsiz > 0) ? copy[nsiz - 1] : 0;
0133       int layer = (nsiz > 1) ? copy[nsiz - 2] : 0;
0134       if (nsiz < 2) {
0135         throw cms::Exception("DDException") << "Funny wafer # " << wafer << " in " << nsiz << " components";
0136       } else if (layer > static_cast<int>(layers.size())) {
0137         edm::LogWarning("HGCalGeom") << "Funny wafer # " << wafer << " Layer " << layer << ":" << layers.size()
0138                                      << " among " << nsiz << " components";
0139       } else {
0140         auto itr = copies.find(wafer);
0141         auto cpy = copiesInLayers[layer].find(wafer);
0142         if (itr != copies.end() && cpy == copiesInLayers[layer].end()) {
0143           copiesInLayers[layer][wafer] = itr->second;
0144         }
0145         if (itr == copies.end()) {
0146           copies[wafer] = wafer2copy.size();
0147           copiesInLayers[layer][wafer] = wafer2copy.size();
0148           double xx = HGCalParameters::k_ScaleFromDDD * fv1.translation().X();
0149           if (std::abs(xx) < tolerance)
0150             xx = 0;
0151           double yy = HGCalParameters::k_ScaleFromDDD * fv1.translation().Y();
0152           if (std::abs(yy) < tolerance)
0153             yy = 0;
0154           wafer2copy.emplace_back(wafer);
0155           GlobalPoint p(xx, yy, HGCalParameters::k_ScaleFromDDD * fv1.translation().Z());
0156           HGCalGeomParameters::cellParameters cell(false, wafer, p);
0157           wafers.emplace_back(cell);
0158           if (names.count(name) == 0) {
0159             std::vector<double> zv, rv;
0160             if (mode == HGCalGeometryMode::Polyhedra) {
0161               const DDPolyhedra& polyhedra = static_cast<DDPolyhedra>(sol);
0162               zv = polyhedra.zVec();
0163               rv = polyhedra.rMaxVec();
0164             } else {
0165               const DDExtrudedPolygon& polygon = static_cast<DDExtrudedPolygon>(sol);
0166               zv = polygon.zVec();
0167               rv = polygon.xVec();
0168             }
0169             php.waferR_ = 2.0 * HGCalParameters::k_ScaleFromDDDToG4 * rv[0] * tan30deg_;
0170             php.waferSize_ = HGCalParameters::k_ScaleFromDDD * rv[0];
0171             double dz = 0.5 * HGCalParameters::k_ScaleFromDDDToG4 * (zv[1] - zv[0]);
0172 #ifdef EDM_ML_DEBUG
0173             edm::LogVerbatim("HGCalGeom")
0174                 << "Mode " << mode << " R " << php.waferSize_ << ":" << php.waferR_ << " z " << dz;
0175 #endif
0176             HGCalParameters::hgtrap mytr;
0177             mytr.lay = 1;
0178             mytr.bl = php.waferR_;
0179             mytr.tl = php.waferR_;
0180             mytr.h = php.waferR_;
0181             mytr.dz = dz;
0182             mytr.alpha = 0.0;
0183             mytr.cellSize = waferSize_;
0184             php.fillModule(mytr, false);
0185             names.insert(name);
0186           }
0187         }
0188       }
0189       dodet = fv1.next();
0190     }
0191   }
0192 
0193   // Finally the cells
0194   std::map<int, int> wafertype;
0195   std::map<int, HGCalGeomParameters::cellParameters> cellsf, cellsc;
0196   DDValue val2(attribute, sdTag3, 0.0);
0197   DDSpecificsMatchesValueFilter filter2{val2};
0198   DDFilteredView fv2(*cpv, filter2);
0199   ok = fv2.firstChild();
0200   if (!ok) {
0201     throw cms::Exception("DDException") << "Attribute " << val2 << " not found but needed.";
0202   } else {
0203     dodet = true;
0204     while (dodet) {
0205       const DDSolid& sol = fv2.logicalPart().solid();
0206       const std::string& name = sol.name().name();
0207       std::vector<int> copy = fv2.copyNumbers();
0208       int nsiz = static_cast<int>(copy.size());
0209       int cellx = (nsiz > 0) ? copy[nsiz - 1] : 0;
0210       int wafer = (nsiz > 1) ? copy[nsiz - 2] : 0;
0211       int cell = HGCalTypes::getUnpackedCell6(cellx);
0212       int type = HGCalTypes::getUnpackedCellType6(cellx);
0213       if (type != 1 && type != 2) {
0214         throw cms::Exception("DDException")
0215             << "Funny cell # " << cell << " type " << type << " in " << nsiz << " components";
0216       } else {
0217         auto ktr = wafertype.find(wafer);
0218         if (ktr == wafertype.end())
0219           wafertype[wafer] = type;
0220         bool newc(false);
0221         std::map<int, HGCalGeomParameters::cellParameters>::iterator itr;
0222         double cellsize = php.cellSize_[0];
0223         if (type == 1) {
0224           itr = cellsf.find(cell);
0225           newc = (itr == cellsf.end());
0226         } else {
0227           itr = cellsc.find(cell);
0228           newc = (itr == cellsc.end());
0229           cellsize = php.cellSize_[1];
0230         }
0231         if (newc) {
0232           bool half = (name.find("Half") != std::string::npos);
0233           double xx = HGCalParameters::k_ScaleFromDDD * fv2.translation().X();
0234           double yy = HGCalParameters::k_ScaleFromDDD * fv2.translation().Y();
0235           if (half) {
0236             math::XYZPointD p1(-2.0 * cellsize / 9.0, 0, 0);
0237             math::XYZPointD p2 = fv2.rotation()(p1);
0238             xx += (HGCalParameters::k_ScaleFromDDD * (p2.X()));
0239             yy += (HGCalParameters::k_ScaleFromDDD * (p2.Y()));
0240 #ifdef EDM_ML_DEBUG
0241             if (std::abs(p2.X()) < HGCalParameters::tol)
0242               p2.SetX(0.0);
0243             if (std::abs(p2.Z()) < HGCalParameters::tol)
0244               p2.SetZ(0.0);
0245             edm::LogVerbatim("HGCalGeom") << "Wafer " << wafer << " Type " << type << " Cell " << cellx << " local "
0246                                           << xx << ":" << yy << " new " << p1 << ":" << p2;
0247 #endif
0248           }
0249           HGCalGeomParameters::cellParameters cp(half, wafer, GlobalPoint(xx, yy, 0));
0250           if (type == 1) {
0251             cellsf[cell] = cp;
0252           } else {
0253             cellsc[cell] = cp;
0254           }
0255         }
0256       }
0257       dodet = fv2.next();
0258     }
0259   }
0260 
0261   loadGeometryHexagon(
0262       layers, trforms, trformUse, copies, copiesInLayers, wafer2copy, wafers, wafertype, cellsf, cellsc, php);
0263 }
0264 
0265 void HGCalGeomParameters::loadGeometryHexagon(const cms::DDCompactView* cpv,
0266                                               HGCalParameters& php,
0267                                               const std::string& sdTag1,
0268                                               const std::string& sdTag2,
0269                                               const std::string& sdTag3,
0270                                               HGCalGeometryMode::WaferMode mode) {
0271   const cms::DDFilter filter("Volume", sdTag1);
0272   cms::DDFilteredView fv((*cpv), filter);
0273   std::map<int, HGCalGeomParameters::layerParameters> layers;
0274   std::vector<HGCalParameters::hgtrform> trforms;
0275   std::vector<bool> trformUse;
0276   std::vector<std::pair<int, int> > trused;
0277 
0278   while (fv.firstChild()) {
0279     const std::vector<double>& pars = fv.parameters();
0280     // Layers first
0281     std::vector<int> copy = fv.copyNos();
0282     int nsiz = static_cast<int>(copy.size());
0283     int lay = (nsiz > 0) ? copy[0] : 0;
0284     int zp = (nsiz > 2) ? copy[2] : -1;
0285     if (zp != 1)
0286       zp = -1;
0287     if (lay == 0) {
0288       throw cms::Exception("DDException") << "Funny layer # " << lay << " zp " << zp << " in " << nsiz << " components";
0289     } else {
0290       if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
0291         php.layer_.emplace_back(lay);
0292       auto itr = layers.find(lay);
0293       double zz = HGCalParameters::k_ScaleFromDD4hep * fv.translation().Z();
0294       if (itr == layers.end()) {
0295         double rin(0), rout(0);
0296         if (dd4hep::isA<dd4hep::Polyhedra>(fv.solid())) {
0297           rin = 0.5 * HGCalParameters::k_ScaleFromDD4hep * (pars[5] + pars[8]);
0298           rout = 0.5 * HGCalParameters::k_ScaleFromDD4hep * (pars[6] + pars[9]);
0299         } else if (dd4hep::isA<dd4hep::Tube>(fv.solid())) {
0300           dd4hep::Tube tubeSeg(fv.solid());
0301           rin = HGCalParameters::k_ScaleFromDD4hep * tubeSeg.rMin();
0302           rout = HGCalParameters::k_ScaleFromDD4hep * tubeSeg.rMax();
0303         }
0304         HGCalGeomParameters::layerParameters laypar(rin, rout, zz);
0305         layers[lay] = laypar;
0306       }
0307       std::pair<int, int> layz(lay, zp);
0308       if (std::find(trused.begin(), trused.end(), layz) == trused.end()) {
0309         trused.emplace_back(layz);
0310         DD3Vector x, y, z;
0311         fv.rotation().GetComponents(x, y, z);
0312         const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0313         const CLHEP::HepRotation hr(rotation);
0314         double xx = HGCalParameters::k_ScaleFromDD4hep * fv.translation().X();
0315         if (std::abs(xx) < tolerance)
0316           xx = 0;
0317         double yy = HGCalParameters::k_ScaleFromDD4hep * fv.translation().Y();
0318         if (std::abs(yy) < tolerance)
0319           yy = 0;
0320         double zz = HGCalParameters::k_ScaleFromDD4hep * fv.translation().Z();
0321         const CLHEP::Hep3Vector h3v(xx, yy, zz);
0322         HGCalParameters::hgtrform mytrf;
0323         mytrf.zp = zp;
0324         mytrf.lay = lay;
0325         mytrf.sec = 0;
0326         mytrf.subsec = 0;
0327         mytrf.h3v = h3v;
0328         mytrf.hr = hr;
0329         trforms.emplace_back(mytrf);
0330         trformUse.emplace_back(false);
0331       }
0332     }
0333   }
0334 
0335   // Then wafers
0336   // This assumes layers are build starting from 1 (which on 25 Jan 2016, they
0337   // were) to ensure that new copy numbers are always added to the end of the
0338   // list.
0339   std::unordered_map<int32_t, int32_t> copies;
0340   HGCalParameters::layer_map copiesInLayers(layers.size() + 1);
0341   std::vector<int32_t> wafer2copy;
0342   std::vector<HGCalGeomParameters::cellParameters> wafers;
0343   const cms::DDFilter filter1("Volume", sdTag2);
0344   cms::DDFilteredView fv1((*cpv), filter1);
0345   bool ok = fv1.firstChild();
0346   if (!ok) {
0347     throw cms::Exception("DDException") << "Attribute " << sdTag2 << " not found but needed.";
0348   } else {
0349     bool dodet = true;
0350     std::unordered_set<std::string> names;
0351     while (dodet) {
0352       const std::string name = static_cast<std::string>(fv1.name());
0353       std::vector<int> copy = fv1.copyNos();
0354       int nsiz = static_cast<int>(copy.size());
0355       int wafer = (nsiz > 0) ? copy[0] : 0;
0356       int layer = (nsiz > 1) ? copy[1] : 0;
0357       if (nsiz < 2) {
0358         throw cms::Exception("DDException") << "Funny wafer # " << wafer << " in " << nsiz << " components";
0359       } else if (layer > static_cast<int>(layers.size())) {
0360         edm::LogWarning("HGCalGeom") << "Funny wafer # " << wafer << " Layer " << layer << ":" << layers.size()
0361                                      << " among " << nsiz << " components";
0362       } else {
0363         auto itr = copies.find(wafer);
0364         auto cpy = copiesInLayers[layer].find(wafer);
0365         if (itr != copies.end() && cpy == copiesInLayers[layer].end()) {
0366           copiesInLayers[layer][wafer] = itr->second;
0367         }
0368         if (itr == copies.end()) {
0369           copies[wafer] = wafer2copy.size();
0370           copiesInLayers[layer][wafer] = wafer2copy.size();
0371           double xx = HGCalParameters::k_ScaleFromDD4hep * fv1.translation().X();
0372           if (std::abs(xx) < tolerance)
0373             xx = 0;
0374           double yy = HGCalParameters::k_ScaleFromDD4hep * fv1.translation().Y();
0375           if (std::abs(yy) < tolerance)
0376             yy = 0;
0377           wafer2copy.emplace_back(wafer);
0378           GlobalPoint p(xx, yy, HGCalParameters::k_ScaleFromDD4hep * fv1.translation().Z());
0379           HGCalGeomParameters::cellParameters cell(false, wafer, p);
0380           wafers.emplace_back(cell);
0381           if (names.count(name) == 0) {
0382             double zv[2], rv;
0383             const std::vector<double>& pars = fv1.parameters();
0384             if (mode == HGCalGeometryMode::Polyhedra) {
0385               zv[0] = pars[4];
0386               zv[1] = pars[7];
0387               rv = pars[6];
0388             } else {
0389               zv[0] = pars[3];
0390               zv[1] = pars[9];
0391               rv = pars[4];
0392             }
0393             php.waferR_ = 2.0 * HGCalParameters::k_ScaleFromDD4hepToG4 * rv * tan30deg_;
0394             php.waferSize_ = HGCalParameters::k_ScaleFromDD4hep * rv;
0395             double dz = 0.5 * HGCalParameters::k_ScaleFromDD4hepToG4 * (zv[1] - zv[0]);
0396 #ifdef EDM_ML_DEBUG
0397             edm::LogVerbatim("HGCalGeom")
0398                 << "Mode " << mode << " R " << php.waferSize_ << ":" << php.waferR_ << " z " << dz;
0399 #endif
0400             HGCalParameters::hgtrap mytr;
0401             mytr.lay = 1;
0402             mytr.bl = php.waferR_;
0403             mytr.tl = php.waferR_;
0404             mytr.h = php.waferR_;
0405             mytr.dz = dz;
0406             mytr.alpha = 0.0;
0407             mytr.cellSize = waferSize_;
0408             php.fillModule(mytr, false);
0409             names.insert(name);
0410           }
0411         }
0412       }
0413       dodet = fv1.firstChild();
0414     }
0415   }
0416 
0417   // Finally the cells
0418   std::map<int, int> wafertype;
0419   std::map<int, HGCalGeomParameters::cellParameters> cellsf, cellsc;
0420   const cms::DDFilter filter2("Volume", sdTag3);
0421   cms::DDFilteredView fv2((*cpv), filter2);
0422   ok = fv2.firstChild();
0423   if (!ok) {
0424     throw cms::Exception("DDException") << "Attribute " << sdTag3 << " not found but needed.";
0425   } else {
0426     bool dodet = true;
0427     while (dodet) {
0428       const std::string name = static_cast<std::string>(fv2.name());
0429       std::vector<int> copy = fv2.copyNos();
0430       int nsiz = static_cast<int>(copy.size());
0431       int cellx = (nsiz > 0) ? copy[0] : 0;
0432       int wafer = (nsiz > 1) ? copy[1] : 0;
0433       int cell = HGCalTypes::getUnpackedCell6(cellx);
0434       int type = HGCalTypes::getUnpackedCellType6(cellx);
0435       if (type != 1 && type != 2) {
0436         throw cms::Exception("DDException")
0437             << "Funny cell # " << cell << " type " << type << " in " << nsiz << " components";
0438       } else {
0439         auto ktr = wafertype.find(wafer);
0440         if (ktr == wafertype.end())
0441           wafertype[wafer] = type;
0442         bool newc(false);
0443         std::map<int, HGCalGeomParameters::cellParameters>::iterator itr;
0444         double cellsize = php.cellSize_[0];
0445         if (type == 1) {
0446           itr = cellsf.find(cell);
0447           newc = (itr == cellsf.end());
0448         } else {
0449           itr = cellsc.find(cell);
0450           newc = (itr == cellsc.end());
0451           cellsize = php.cellSize_[1];
0452         }
0453         if (newc) {
0454           bool half = (name.find("Half") != std::string::npos);
0455           double xx = HGCalParameters::k_ScaleFromDD4hep * fv2.translation().X();
0456           double yy = HGCalParameters::k_ScaleFromDD4hep * fv2.translation().Y();
0457           if (half) {
0458             math::XYZPointD p1(-2.0 * cellsize / 9.0, 0, 0);
0459             math::XYZPointD p2 = fv2.rotation()(p1);
0460             xx += (HGCalParameters::k_ScaleFromDDD * (p2.X()));
0461             yy += (HGCalParameters::k_ScaleFromDDD * (p2.Y()));
0462 #ifdef EDM_ML_DEBUG
0463             if (std::abs(p2.X()) < HGCalParameters::tol)
0464               p2.SetX(0.0);
0465             if (std::abs(p2.Z()) < HGCalParameters::tol)
0466               p2.SetZ(0.0);
0467             edm::LogVerbatim("HGCalGeom") << "Wafer " << wafer << " Type " << type << " Cell " << cellx << " local "
0468                                           << xx << ":" << yy << " new " << p1 << ":" << p2;
0469 #endif
0470           }
0471           HGCalGeomParameters::cellParameters cp(half, wafer, GlobalPoint(xx, yy, 0));
0472           if (type == 1) {
0473             cellsf[cell] = cp;
0474           } else {
0475             cellsc[cell] = cp;
0476           }
0477         }
0478       }
0479       dodet = fv2.firstChild();
0480     }
0481   }
0482 
0483   loadGeometryHexagon(
0484       layers, trforms, trformUse, copies, copiesInLayers, wafer2copy, wafers, wafertype, cellsf, cellsc, php);
0485 }
0486 
0487 void HGCalGeomParameters::loadGeometryHexagon(const std::map<int, HGCalGeomParameters::layerParameters>& layers,
0488                                               std::vector<HGCalParameters::hgtrform>& trforms,
0489                                               std::vector<bool>& trformUse,
0490                                               const std::unordered_map<int32_t, int32_t>& copies,
0491                                               const HGCalParameters::layer_map& copiesInLayers,
0492                                               const std::vector<int32_t>& wafer2copy,
0493                                               const std::vector<HGCalGeomParameters::cellParameters>& wafers,
0494                                               const std::map<int, int>& wafertype,
0495                                               const std::map<int, HGCalGeomParameters::cellParameters>& cellsf,
0496                                               const std::map<int, HGCalGeomParameters::cellParameters>& cellsc,
0497                                               HGCalParameters& php) {
0498   if (((cellsf.size() + cellsc.size()) == 0) || (wafers.empty()) || (layers.empty())) {
0499     throw cms::Exception("DDException") << "HGCalGeomParameters: mismatch between geometry and specpar: cells "
0500                                         << cellsf.size() << ":" << cellsc.size() << " wafers " << wafers.size()
0501                                         << " layers " << layers.size();
0502   }
0503 
0504   for (unsigned int i = 0; i < layers.size(); ++i) {
0505     for (auto& layer : layers) {
0506       if (layer.first == static_cast<int>(i + php.firstLayer_)) {
0507         php.layerIndex_.emplace_back(i);
0508         php.rMinLayHex_.emplace_back(layer.second.rmin);
0509         php.rMaxLayHex_.emplace_back(layer.second.rmax);
0510         php.zLayerHex_.emplace_back(layer.second.zpos);
0511         break;
0512       }
0513     }
0514   }
0515 
0516   for (unsigned int i = 0; i < php.layer_.size(); ++i) {
0517     for (unsigned int i1 = 0; i1 < trforms.size(); ++i1) {
0518       if (!trformUse[i1] && php.layerGroup_[trforms[i1].lay - 1] == static_cast<int>(i + 1)) {
0519         trforms[i1].h3v *= static_cast<double>(HGCalParameters::k_ScaleFromDDD);
0520         trforms[i1].lay = (i + 1);
0521         trformUse[i1] = true;
0522         php.fillTrForm(trforms[i1]);
0523         int nz(1);
0524         for (unsigned int i2 = i1 + 1; i2 < trforms.size(); ++i2) {
0525           if (!trformUse[i2] && trforms[i2].zp == trforms[i1].zp &&
0526               php.layerGroup_[trforms[i2].lay - 1] == static_cast<int>(i + 1)) {
0527             php.addTrForm(trforms[i2].h3v);
0528             nz++;
0529             trformUse[i2] = true;
0530           }
0531         }
0532         if (nz > 0) {
0533           php.scaleTrForm(double(1.0 / nz));
0534         }
0535       }
0536     }
0537   }
0538 
0539   double rmin = HGCalParameters::k_ScaleFromDDD * php.waferR_;
0540   for (unsigned i = 0; i < wafer2copy.size(); ++i) {
0541     php.waferCopy_.emplace_back(wafer2copy[i]);
0542     php.waferPosX_.emplace_back(wafers[i].xyz.x());
0543     php.waferPosY_.emplace_back(wafers[i].xyz.y());
0544     auto ktr = wafertype.find(wafer2copy[i]);
0545     int typet = (ktr == wafertype.end()) ? 0 : (ktr->second);
0546     php.waferTypeT_.emplace_back(typet);
0547     double r = wafers[i].xyz.perp();
0548     int type(3);
0549     for (int k = 1; k < 4; ++k) {
0550       if ((r + rmin) <= php.boundR_[k]) {
0551         type = k;
0552         break;
0553       }
0554     }
0555     php.waferTypeL_.emplace_back(type);
0556   }
0557   php.copiesInLayers_ = copiesInLayers;
0558   php.nSectors_ = static_cast<int>(php.waferCopy_.size());
0559 
0560   std::vector<HGCalGeomParameters::cellParameters>::const_iterator itrf = wafers.end();
0561   for (unsigned int i = 0; i < cellsf.size(); ++i) {
0562     auto itr = cellsf.find(i);
0563     if (itr == cellsf.end()) {
0564       throw cms::Exception("DDException") << "HGCalGeomParameters: missing info for fine cell number " << i;
0565     } else {
0566       double xx = (itr->second).xyz.x();
0567       double yy = (itr->second).xyz.y();
0568       int waf = (itr->second).wafer;
0569       std::pair<double, double> xy = cellPosition(wafers, itrf, waf, xx, yy);
0570       php.cellFineX_.emplace_back(xy.first);
0571       php.cellFineY_.emplace_back(xy.second);
0572       php.cellFineHalf_.emplace_back((itr->second).half);
0573     }
0574   }
0575   itrf = wafers.end();
0576   for (unsigned int i = 0; i < cellsc.size(); ++i) {
0577     auto itr = cellsc.find(i);
0578     if (itr == cellsc.end()) {
0579       throw cms::Exception("DDException") << "HGCalGeomParameters: missing info for coarse cell number " << i;
0580     } else {
0581       double xx = (itr->second).xyz.x();
0582       double yy = (itr->second).xyz.y();
0583       int waf = (itr->second).wafer;
0584       std::pair<double, double> xy = cellPosition(wafers, itrf, waf, xx, yy);
0585       php.cellCoarseX_.emplace_back(xy.first);
0586       php.cellCoarseY_.emplace_back(xy.second);
0587       php.cellCoarseHalf_.emplace_back((itr->second).half);
0588     }
0589   }
0590   int depth(0);
0591   for (unsigned int i = 0; i < php.layerGroup_.size(); ++i) {
0592     bool first(true);
0593     for (unsigned int k = 0; k < php.layerGroup_.size(); ++k) {
0594       if (php.layerGroup_[k] == static_cast<int>(i + 1)) {
0595         if (first) {
0596           php.depth_.emplace_back(i + 1);
0597           php.depthIndex_.emplace_back(depth);
0598           php.depthLayerF_.emplace_back(k);
0599           ++depth;
0600           first = false;
0601         }
0602       }
0603     }
0604   }
0605   HGCalParameters::hgtrap mytr = php.getModule(0, false);
0606   mytr.bl *= HGCalParameters::k_ScaleFromDDD;
0607   mytr.tl *= HGCalParameters::k_ScaleFromDDD;
0608   mytr.h *= HGCalParameters::k_ScaleFromDDD;
0609   mytr.dz *= HGCalParameters::k_ScaleFromDDD;
0610   mytr.cellSize *= HGCalParameters::k_ScaleFromDDD;
0611   double dz = mytr.dz;
0612   php.fillModule(mytr, true);
0613   mytr.dz = 2 * dz;
0614   php.fillModule(mytr, true);
0615   mytr.dz = 3 * dz;
0616   php.fillModule(mytr, true);
0617 #ifdef EDM_ML_DEBUG
0618   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters finds " << php.zLayerHex_.size() << " layers";
0619   for (unsigned int i = 0; i < php.zLayerHex_.size(); ++i) {
0620     int k = php.layerIndex_[i];
0621     edm::LogVerbatim("HGCalGeom") << "Layer[" << i << ":" << k << ":" << php.layer_[k]
0622                                   << "] with r = " << php.rMinLayHex_[i] << ":" << php.rMaxLayHex_[i]
0623                                   << " at z = " << php.zLayerHex_[i];
0624   }
0625   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters has " << php.depthIndex_.size() << " depths";
0626   for (unsigned int i = 0; i < php.depthIndex_.size(); ++i) {
0627     int k = php.depthIndex_[i];
0628     edm::LogVerbatim("HGCalGeom") << "Reco Layer[" << i << ":" << k << "]  First Layer " << php.depthLayerF_[i]
0629                                   << " Depth " << php.depth_[k];
0630   }
0631   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters finds " << php.nSectors_ << " wafers";
0632   for (unsigned int i = 0; i < php.waferCopy_.size(); ++i)
0633     edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << ": " << php.waferCopy_[i] << "] type " << php.waferTypeL_[i]
0634                                   << ":" << php.waferTypeT_[i] << " at (" << php.waferPosX_[i] << ","
0635                                   << php.waferPosY_[i] << ",0)";
0636   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: wafer radius " << php.waferR_ << " and dimensions of the "
0637                                 << "wafers:";
0638   edm::LogVerbatim("HGCalGeom") << "Sim[0] " << php.moduleLayS_[0] << " dx " << php.moduleBlS_[0] << ":"
0639                                 << php.moduleTlS_[0] << " dy " << php.moduleHS_[0] << " dz " << php.moduleDzS_[0]
0640                                 << " alpha " << php.moduleAlphaS_[0];
0641   for (unsigned int k = 0; k < php.moduleLayR_.size(); ++k)
0642     edm::LogVerbatim("HGCalGeom") << "Rec[" << k << "] " << php.moduleLayR_[k] << " dx " << php.moduleBlR_[k] << ":"
0643                                   << php.moduleTlR_[k] << " dy " << php.moduleHR_[k] << " dz " << php.moduleDzR_[k]
0644                                   << " alpha " << php.moduleAlphaR_[k];
0645   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters finds " << php.cellFineX_.size() << " fine cells in a  wafer";
0646   for (unsigned int i = 0; i < php.cellFineX_.size(); ++i)
0647     edm::LogVerbatim("HGCalGeom") << "Fine Cell[" << i << "] at (" << php.cellFineX_[i] << "," << php.cellFineY_[i]
0648                                   << ",0)";
0649   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters finds " << php.cellCoarseX_.size()
0650                                 << " coarse cells in a wafer";
0651   for (unsigned int i = 0; i < php.cellCoarseX_.size(); ++i)
0652     edm::LogVerbatim("HGCalGeom") << "Coarse Cell[" << i << "] at (" << php.cellCoarseX_[i] << ","
0653                                   << php.cellCoarseY_[i] << ",0)";
0654   edm::LogVerbatim("HGCalGeom") << "Obtained " << php.trformIndex_.size() << " transformation matrices";
0655   for (unsigned int k = 0; k < php.trformIndex_.size(); ++k) {
0656     edm::LogVerbatim("HGCalGeom") << "Matrix[" << k << "] (" << std::hex << php.trformIndex_[k] << std::dec
0657                                   << ") Translation (" << php.trformTranX_[k] << ", " << php.trformTranY_[k] << ", "
0658                                   << php.trformTranZ_[k] << " Rotation (" << php.trformRotXX_[k] << ", "
0659                                   << php.trformRotYX_[k] << ", " << php.trformRotZX_[k] << ", " << php.trformRotXY_[k]
0660                                   << ", " << php.trformRotYY_[k] << ", " << php.trformRotZY_[k] << ", "
0661                                   << php.trformRotXZ_[k] << ", " << php.trformRotYZ_[k] << ", " << php.trformRotZZ_[k]
0662                                   << ")";
0663   }
0664   edm::LogVerbatim("HGCalGeom") << "Dump copiesInLayers for " << php.copiesInLayers_.size() << " layers";
0665   for (unsigned int k = 0; k < php.copiesInLayers_.size(); ++k) {
0666     const auto& theModules = php.copiesInLayers_[k];
0667     edm::LogVerbatim("HGCalGeom") << "Layer " << k << ":" << theModules.size();
0668     int k2(0);
0669     for (std::unordered_map<int, int>::const_iterator itr = theModules.begin(); itr != theModules.end(); ++itr, ++k2) {
0670       edm::LogVerbatim("HGCalGeom") << "[" << k2 << "] " << itr->first << ":" << itr->second;
0671     }
0672   }
0673 #endif
0674 }
0675 
0676 void HGCalGeomParameters::loadGeometryHexagon8(const DDFilteredView& _fv, HGCalParameters& php, int firstLayer) {
0677   DDFilteredView fv = _fv;
0678   bool dodet(true);
0679   std::map<int, HGCalGeomParameters::layerParameters> layers;
0680   std::map<std::pair<int, int>, HGCalParameters::hgtrform> trforms;
0681   int levelTop = 3 + std::max(php.levelT_[0], php.levelT_[1]);
0682 #ifdef EDM_ML_DEBUG
0683   int ntot(0);
0684 #endif
0685   while (dodet) {
0686 #ifdef EDM_ML_DEBUG
0687     ++ntot;
0688 #endif
0689     std::vector<int> copy = fv.copyNumbers();
0690     int nsiz = static_cast<int>(copy.size());
0691     if (nsiz < levelTop) {
0692       int lay = copy[nsiz - 1];
0693       int zside = (nsiz > php.levelZSide_) ? copy[php.levelZSide_] : -1;
0694       if (zside != 1)
0695         zside = -1;
0696       const DDSolid& sol = fv.logicalPart().solid();
0697 #ifdef EDM_ML_DEBUG
0698       edm::LogVerbatim("HGCalGeom") << sol.name() << " shape " << sol.shape() << " size " << nsiz << ":" << levelTop
0699                                     << " lay " << lay << " z " << zside;
0700 #endif
0701       if (lay == 0) {
0702         throw cms::Exception("DDException")
0703             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
0704       } else if (sol.shape() == DDSolidShape::ddtubs) {
0705         if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
0706           php.layer_.emplace_back(lay);
0707         const DDTubs& tube = static_cast<DDTubs>(sol);
0708         double rin = HGCalParameters::k_ScaleFromDDD * tube.rIn();
0709         double rout = HGCalParameters::k_ScaleFromDDD * tube.rOut();
0710         auto itr = layers.find(lay);
0711         if (itr == layers.end()) {
0712           double zp = HGCalParameters::k_ScaleFromDDD * fv.translation().Z();
0713           HGCalGeomParameters::layerParameters laypar(rin, rout, zp);
0714           layers[lay] = laypar;
0715         } else {
0716           (itr->second).rmin = std::min(rin, (itr->second).rmin);
0717           (itr->second).rmax = std::max(rout, (itr->second).rmax);
0718         }
0719         if (trforms.find(std::make_pair(lay, zside)) == trforms.end()) {
0720           DD3Vector x, y, z;
0721           fv.rotation().GetComponents(x, y, z);
0722           const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0723           const CLHEP::HepRotation hr(rotation);
0724           double xx =
0725               ((std::abs(fv.translation().X()) < tolerance) ? 0
0726                                                             : HGCalParameters::k_ScaleFromDDD * fv.translation().X());
0727           double yy =
0728               ((std::abs(fv.translation().Y()) < tolerance) ? 0
0729                                                             : HGCalParameters::k_ScaleFromDDD * fv.translation().Y());
0730           const CLHEP::Hep3Vector h3v(xx, yy, HGCalParameters::k_ScaleFromDDD * fv.translation().Z());
0731           HGCalParameters::hgtrform mytrf;
0732           mytrf.zp = zside;
0733           mytrf.lay = lay;
0734           mytrf.sec = 0;
0735           mytrf.subsec = 0;
0736           mytrf.h3v = h3v;
0737           mytrf.hr = hr;
0738           trforms[std::make_pair(lay, zside)] = mytrf;
0739         }
0740       }
0741     }
0742     dodet = fv.next();
0743   }
0744 #ifdef EDM_ML_DEBUG
0745   edm::LogVerbatim("HGCalGeom") << "Total # of views " << ntot;
0746 #endif
0747   loadGeometryHexagon8(layers, trforms, firstLayer, php);
0748 }
0749 
0750 void HGCalGeomParameters::loadGeometryHexagon8(const cms::DDCompactView* cpv,
0751                                                HGCalParameters& php,
0752                                                const std::string& sdTag1,
0753                                                int firstLayer) {
0754   const cms::DDFilter filter("Volume", sdTag1);
0755   cms::DDFilteredView fv((*cpv), filter);
0756   std::map<int, HGCalGeomParameters::layerParameters> layers;
0757   std::map<std::pair<int, int>, HGCalParameters::hgtrform> trforms;
0758   int levelTop = 3 + std::max(php.levelT_[0], php.levelT_[1]);
0759 #ifdef EDM_ML_DEBUG
0760   int ntot(0);
0761 #endif
0762   while (fv.firstChild()) {
0763 #ifdef EDM_ML_DEBUG
0764     ++ntot;
0765 #endif
0766     // Layers first
0767     int nsiz = static_cast<int>(fv.level());
0768     if (nsiz < levelTop) {
0769       std::vector<int> copy = fv.copyNos();
0770       int lay = copy[0];
0771       int zside = (nsiz > php.levelZSide_) ? copy[nsiz - php.levelZSide_ - 1] : -1;
0772       if (zside != 1)
0773         zside = -1;
0774 #ifdef EDM_ML_DEBUG
0775       edm::LogVerbatim("HGCalGeom") << fv.name() << " shape " << cms::dd::name(cms::DDSolidShapeMap, fv.shape())
0776                                     << " size " << nsiz << ":" << levelTop << " lay " << lay << " z " << zside << ":"
0777                                     << php.levelZSide_;
0778 #endif
0779       if (lay == 0) {
0780         throw cms::Exception("DDException")
0781             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
0782       } else if (fv.shape() == cms::DDSolidShape::ddtubs) {
0783         if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
0784           php.layer_.emplace_back(lay);
0785         const std::vector<double>& pars = fv.parameters();
0786         double rin = HGCalParameters::k_ScaleFromDD4hep * pars[0];
0787         double rout = HGCalParameters::k_ScaleFromDD4hep * pars[1];
0788         auto itr = layers.find(lay);
0789         if (itr == layers.end()) {
0790           double zp = HGCalParameters::k_ScaleFromDD4hep * fv.translation().Z();
0791           HGCalGeomParameters::layerParameters laypar(rin, rout, zp);
0792           layers[lay] = laypar;
0793         } else {
0794           (itr->second).rmin = std::min(rin, (itr->second).rmin);
0795           (itr->second).rmax = std::max(rout, (itr->second).rmax);
0796         }
0797         if (trforms.find(std::make_pair(lay, zside)) == trforms.end()) {
0798           DD3Vector x, y, z;
0799           fv.rotation().GetComponents(x, y, z);
0800           const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0801           const CLHEP::HepRotation hr(rotation);
0802           double xx = ((std::abs(fv.translation().X()) < tolerance)
0803                            ? 0
0804                            : HGCalParameters::k_ScaleFromDD4hep * fv.translation().X());
0805           double yy = ((std::abs(fv.translation().Y()) < tolerance)
0806                            ? 0
0807                            : HGCalParameters::k_ScaleFromDD4hep * fv.translation().Y());
0808           const CLHEP::Hep3Vector h3v(xx, yy, HGCalParameters::k_ScaleFromDD4hep * fv.translation().Z());
0809           HGCalParameters::hgtrform mytrf;
0810           mytrf.zp = zside;
0811           mytrf.lay = lay;
0812           mytrf.sec = 0;
0813           mytrf.subsec = 0;
0814           mytrf.h3v = h3v;
0815           mytrf.hr = hr;
0816           trforms[std::make_pair(lay, zside)] = mytrf;
0817         }
0818       }
0819     }
0820   }
0821 #ifdef EDM_ML_DEBUG
0822   edm::LogVerbatim("HGCalGeom") << "Total # of views " << ntot;
0823 #endif
0824   loadGeometryHexagon8(layers, trforms, firstLayer, php);
0825 }
0826 
0827 void HGCalGeomParameters::loadGeometryHexagonModule(const DDCompactView* cpv,
0828                                                     HGCalParameters& php,
0829                                                     const std::string& sdTag1,
0830                                                     const std::string& sdTag2,
0831                                                     int firstLayer) {
0832 #ifdef EDM_ML_DEBUG
0833   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters (DDD)::loadGeometryHexagonModule called with tags " << sdTag1
0834                                 << ":" << sdTag2 << " firstLayer " << firstLayer << ":" << php.firstMixedLayer_;
0835   int ntot1(0), ntot2(0);
0836 #endif
0837   std::map<int, HGCalGeomParameters::layerParameters> layers;
0838   std::map<std::pair<int, int>, double> zvals;
0839   std::map<std::pair<int, int>, HGCalParameters::hgtrform> trforms;
0840   int levelTop = php.levelT_[0];
0841 
0842   std::string attribute = "Volume";
0843   DDValue val1(attribute, sdTag2, 0.0);
0844   DDSpecificsMatchesValueFilter filter1{val1};
0845   DDFilteredView fv1(*cpv, filter1);
0846   bool dodet = fv1.firstChild();
0847   while (dodet) {
0848 #ifdef EDM_ML_DEBUG
0849     ++ntot1;
0850 #endif
0851     std::vector<int> copy = fv1.copyNumbers();
0852     int nsiz = static_cast<int>(copy.size());
0853     if (levelTop < nsiz) {
0854       int lay = copy[levelTop];
0855       int zside = (nsiz > php.levelZSide_) ? copy[php.levelZSide_] : -1;
0856       if (zside != 1)
0857         zside = -1;
0858       if (lay == 0) {
0859         throw cms::Exception("DDException")
0860             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
0861       } else {
0862         if (zvals.find(std::make_pair(lay, zside)) == zvals.end()) {
0863           zvals[std::make_pair(lay, zside)] = HGCalParameters::k_ScaleFromDDD * fv1.translation().Z();
0864 #ifdef EDM_ML_DEBUG
0865           std::ostringstream st1;
0866           st1 << "Name0 " << fv1.name() << " LTop " << levelTop << ":" << lay << " ZSide " << zside << " # of levels "
0867               << nsiz;
0868           for (const auto& c : copy)
0869             st1 << ":" << c;
0870           st1 << " Z " << zvals[std::make_pair(lay, zside)];
0871           edm::LogVerbatim("HGCalGeom") << st1.str();
0872 #endif
0873         }
0874       }
0875     }
0876     dodet = fv1.next();
0877   }
0878 
0879   DDValue val2(attribute, sdTag1, 0.0);
0880   DDSpecificsMatchesValueFilter filter2{val2};
0881   DDFilteredView fv2(*cpv, filter2);
0882   dodet = fv2.firstChild();
0883   while (dodet) {
0884     std::vector<int> copy = fv2.copyNumbers();
0885     int nsiz = static_cast<int>(copy.size());
0886 #ifdef EDM_ML_DEBUG
0887     ++ntot2;
0888     edm::LogVerbatim("HGCalGeom") << "loadGeometryHexagonModule:: nsiz " << nsiz << " Ltop " << levelTop;
0889 #endif
0890     if (levelTop < nsiz) {
0891       int lay = copy[levelTop];
0892       int zside = (nsiz > php.levelZSide_) ? copy[php.levelZSide_] : -1;
0893       if (zside != 1)
0894         zside = -1;
0895       const DDSolid& sol = fv2.logicalPart().solid();
0896 #ifdef EDM_ML_DEBUG
0897       std::ostringstream st2;
0898       st2 << "Name1 " << sol.name() << " shape " << sol.shape() << " LTop " << levelTop << ":" << lay << " ZSide "
0899           << zside << ":" << php.levelZSide_ << " # of levels " << nsiz;
0900       for (const auto& c : copy)
0901         st2 << ":" << c;
0902       edm::LogVerbatim("HGCalGeom") << st2.str();
0903 #endif
0904       if (lay == 0) {
0905         throw cms::Exception("DDException")
0906             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
0907       } else if ((sol.shape() == DDSolidShape::ddtubs) || (sol.shape() == DDSolidShape::ddbox)) {
0908         if (zvals.find(std::make_pair(lay, zside)) != zvals.end()) {
0909           if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
0910             php.layer_.emplace_back(lay);
0911           auto itr = layers.find(lay);
0912           if (itr == layers.end()) {
0913             double rin(0), rout(0);
0914             if (sol.shape() == DDSolidShape::ddtubs) {
0915               const DDTubs& tube = static_cast<DDTubs>(sol);
0916               rin = HGCalParameters::k_ScaleFromDDD * tube.rIn();
0917               rout = (php.firstMixedLayer_ > 0 && lay >= php.firstMixedLayer_)
0918                          ? php.radiusMixBoundary_[lay - php.firstMixedLayer_]
0919                          : HGCalParameters::k_ScaleFromDDD * tube.rOut();
0920             } else {
0921               const DDBox& box = static_cast<DDBox>(sol);
0922               rout = HGCalParameters::k_ScaleFromDDD * box.halfX();
0923             }
0924             double zp = zvals[std::make_pair(lay, 1)];
0925             HGCalGeomParameters::layerParameters laypar(rin, rout, zp);
0926             layers[lay] = laypar;
0927 #ifdef EDM_ML_DEBUG
0928             std::ostringstream st3;
0929             st3 << "Name1 " << fv2.name() << " LTop " << levelTop << ":" << lay << " ZSide " << zside << " # of levels "
0930                 << nsiz;
0931             for (const auto& c : copy)
0932               st3 << ":" << c;
0933             st3 << " R " << rin << ":" << rout;
0934             edm::LogVerbatim("HGCalGeom") << st3.str();
0935 #endif
0936           }
0937 
0938           if (trforms.find(std::make_pair(lay, zside)) == trforms.end()) {
0939             DD3Vector x, y, z;
0940             fv2.rotation().GetComponents(x, y, z);
0941             const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
0942             const CLHEP::HepRotation hr(rotation);
0943             double xx = ((std::abs(fv2.translation().X()) < tolerance)
0944                              ? 0
0945                              : HGCalParameters::k_ScaleFromDDD * fv2.translation().X());
0946             double yy = ((std::abs(fv2.translation().Y()) < tolerance)
0947                              ? 0
0948                              : HGCalParameters::k_ScaleFromDDD * fv2.translation().Y());
0949             const CLHEP::Hep3Vector h3v(xx, yy, zvals[std::make_pair(lay, zside)]);
0950             HGCalParameters::hgtrform mytrf;
0951             mytrf.zp = zside;
0952             mytrf.lay = lay;
0953             mytrf.sec = 0;
0954             mytrf.subsec = 0;
0955             mytrf.h3v = h3v;
0956             mytrf.hr = hr;
0957             trforms[std::make_pair(lay, zside)] = mytrf;
0958 #ifdef EDM_ML_DEBUG
0959             edm::LogVerbatim("HGCalGeom") << "Translation " << h3v;
0960 #endif
0961           }
0962         }
0963       }
0964     }
0965     dodet = fv2.next();
0966   }
0967 #ifdef EDM_ML_DEBUG
0968   edm::LogVerbatim("HGCalGeom") << "Total # of views " << ntot1 << ":" << ntot2;
0969 #endif
0970   loadGeometryHexagon8(layers, trforms, firstLayer, php);
0971 }
0972 
0973 void HGCalGeomParameters::loadGeometryHexagonModule(const cms::DDCompactView* cpv,
0974                                                     HGCalParameters& php,
0975                                                     const std::string& sdTag1,
0976                                                     const std::string& sdTag2,
0977                                                     int firstLayer) {
0978 #ifdef EDM_ML_DEBUG
0979   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters (DD4hep)::loadGeometryHexagonModule called with tags " << sdTag1
0980                                 << ":" << sdTag2 << " firstLayer " << firstLayer;
0981   int ntot1(0), ntot2(0);
0982 #endif
0983   std::map<int, HGCalGeomParameters::layerParameters> layers;
0984   std::map<std::pair<int, int>, HGCalParameters::hgtrform> trforms;
0985   std::map<std::pair<int, int>, double> zvals;
0986   int levelTop = php.levelT_[0];
0987 
0988   const cms::DDFilter filter1("Volume", sdTag2);
0989   cms::DDFilteredView fv1((*cpv), filter1);
0990   while (fv1.firstChild()) {
0991 #ifdef EDM_ML_DEBUG
0992     ++ntot1;
0993 #endif
0994     int nsiz = static_cast<int>(fv1.level());
0995     if (nsiz > levelTop) {
0996       std::vector<int> copy = fv1.copyNos();
0997       int lay = copy[nsiz - levelTop - 1];
0998       int zside = (nsiz > php.levelZSide_) ? copy[nsiz - php.levelZSide_ - 1] : -1;
0999       if (zside != 1)
1000         zside = -1;
1001       if (lay == 0) {
1002         throw cms::Exception("DDException")
1003             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
1004       } else {
1005         if (zvals.find(std::make_pair(lay, zside)) == zvals.end()) {
1006           zvals[std::make_pair(lay, zside)] = HGCalParameters::k_ScaleFromDD4hep * fv1.translation().Z();
1007 #ifdef EDM_ML_DEBUG
1008           std::ostringstream st1;
1009           st1 << "Name0 " << fv1.name() << " LTop " << levelTop << ":" << lay << " ZSide " << zside << " # of levels "
1010               << nsiz;
1011           for (const auto& c : copy)
1012             st1 << ":" << c;
1013           st1 << " Z " << zvals[std::make_pair(lay, zside)];
1014           edm::LogVerbatim("HGCalGeom") << st1.str();
1015 #endif
1016         }
1017       }
1018     }
1019   }
1020 
1021   const cms::DDFilter filter2("Volume", sdTag1);
1022   cms::DDFilteredView fv2((*cpv), filter2);
1023   while (fv2.firstChild()) {
1024     // Layers first
1025     int nsiz = static_cast<int>(fv2.level());
1026 #ifdef EDM_ML_DEBUG
1027     ++ntot2;
1028     edm::LogVerbatim("HGCalGeom") << "loadGeometryHexagonModule:: nsiz " << nsiz << " Ltop " << levelTop;
1029 #endif
1030     if (nsiz > levelTop) {
1031       std::vector<int> copy = fv2.copyNos();
1032       int lay = copy[nsiz - levelTop - 1];
1033       int zside = (nsiz > php.levelZSide_) ? copy[nsiz - php.levelZSide_ - 1] : -1;
1034       if (zside != 1)
1035         zside = -1;
1036 #ifdef EDM_ML_DEBUG
1037       std::ostringstream st2;
1038       st2 << "Name1 " << fv2.name() << "Shape " << cms::dd::name(cms::DDSolidShapeMap, fv2.shape()) << " LTop "
1039           << levelTop << ":" << lay << " ZSide " << zside << ":" << php.levelZSide_ << " # of levels " << nsiz;
1040       for (const auto& c : copy)
1041         st2 << ":" << c;
1042       edm::LogVerbatim("HGCalGeom") << st2.str();
1043 #endif
1044       if (lay == 0) {
1045         throw cms::Exception("DDException")
1046             << "Funny layer # " << lay << " zp " << zside << " in " << nsiz << " components";
1047       } else {
1048         if (zvals.find(std::make_pair(lay, zside)) != zvals.end()) {
1049           if (std::find(php.layer_.begin(), php.layer_.end(), lay) == php.layer_.end())
1050             php.layer_.emplace_back(lay);
1051           auto itr = layers.find(lay);
1052           if (itr == layers.end()) {
1053             const std::vector<double>& pars = fv2.parameters();
1054             double rin(0), rout(0);
1055             if (dd4hep::isA<dd4hep::Box>(fv2.solid())) {
1056               rout = HGCalParameters::k_ScaleFromDD4hep * pars[0];
1057             } else {
1058               rin = HGCalParameters::k_ScaleFromDD4hep * pars[0];
1059               rout = (php.firstMixedLayer_ > 0 && lay >= php.firstMixedLayer_)
1060                          ? php.radiusMixBoundary_[lay - php.firstMixedLayer_]
1061                          : HGCalParameters::k_ScaleFromDD4hep * pars[1];
1062             }
1063             double zp = zvals[std::make_pair(lay, 1)];
1064             HGCalGeomParameters::layerParameters laypar(rin, rout, zp);
1065             layers[lay] = laypar;
1066 #ifdef EDM_ML_DEBUG
1067             std::ostringstream st3;
1068             st3 << "Name2 " << fv2.name() << " LTop " << levelTop << ":" << lay << " ZSide " << zside << " # of levels "
1069                 << nsiz;
1070             for (const auto& c : copy)
1071               st3 << ":" << c;
1072             st3 << " R " << rin << ":" << rout;
1073             edm::LogVerbatim("HGCalGeom") << st3.str();
1074 #endif
1075           }
1076 
1077           if (trforms.find(std::make_pair(lay, zside)) == trforms.end()) {
1078             DD3Vector x, y, z;
1079             fv2.rotation().GetComponents(x, y, z);
1080             const CLHEP::HepRep3x3 rotation(x.X(), y.X(), z.X(), x.Y(), y.Y(), z.Y(), x.Z(), y.Z(), z.Z());
1081             const CLHEP::HepRotation hr(rotation);
1082             double xx = ((std::abs(fv2.translation().X()) < tolerance)
1083                              ? 0
1084                              : HGCalParameters::k_ScaleFromDD4hep * fv2.translation().X());
1085             double yy = ((std::abs(fv2.translation().Y()) < tolerance)
1086                              ? 0
1087                              : HGCalParameters::k_ScaleFromDD4hep * fv2.translation().Y());
1088             const CLHEP::Hep3Vector h3v(xx, yy, zvals[std::make_pair(lay, zside)]);
1089             HGCalParameters::hgtrform mytrf;
1090             mytrf.zp = zside;
1091             mytrf.lay = lay;
1092             mytrf.sec = 0;
1093             mytrf.subsec = 0;
1094             mytrf.h3v = h3v;
1095             mytrf.hr = hr;
1096             trforms[std::make_pair(lay, zside)] = mytrf;
1097 #ifdef EDM_ML_DEBUG
1098             edm::LogVerbatim("HGCalGeom") << "Translation " << h3v;
1099 #endif
1100           }
1101         }
1102       }
1103     }
1104   }
1105 #ifdef EDM_ML_DEBUG
1106   edm::LogVerbatim("HGCalGeom") << "Total # of views " << ntot1 << ":" << ntot2;
1107 #endif
1108   loadGeometryHexagon8(layers, trforms, firstLayer, php);
1109 }
1110 
1111 void HGCalGeomParameters::loadGeometryHexagon8(const std::map<int, HGCalGeomParameters::layerParameters>& layers,
1112                                                std::map<std::pair<int, int>, HGCalParameters::hgtrform>& trforms,
1113                                                const int& firstLayer,
1114                                                HGCalParameters& php) {
1115   double rmin(0), rmax(0);
1116   for (unsigned int i = 0; i < layers.size(); ++i) {
1117     for (auto& layer : layers) {
1118       if (layer.first == static_cast<int>(i + firstLayer)) {
1119         php.layerIndex_.emplace_back(i);
1120         php.rMinLayHex_.emplace_back(layer.second.rmin);
1121         php.rMaxLayHex_.emplace_back(layer.second.rmax);
1122         php.zLayerHex_.emplace_back(layer.second.zpos);
1123         if (i == 0) {
1124           rmin = layer.second.rmin;
1125           rmax = layer.second.rmax;
1126         } else {
1127           if (rmin > layer.second.rmin)
1128             rmin = layer.second.rmin;
1129           if (rmax < layer.second.rmax)
1130             rmax = layer.second.rmax;
1131         }
1132         break;
1133       }
1134     }
1135   }
1136   php.rLimit_.emplace_back(rmin);
1137   php.rLimit_.emplace_back(rmax);
1138   php.depth_ = php.layer_;
1139   php.depthIndex_ = php.layerIndex_;
1140   php.depthLayerF_ = php.layerIndex_;
1141 
1142   for (unsigned int i = 0; i < php.layer_.size(); ++i) {
1143     for (auto& trform : trforms) {
1144       if (trform.first.first == static_cast<int>(i + firstLayer)) {
1145         php.fillTrForm(trform.second);
1146       }
1147     }
1148   }
1149 
1150 #ifdef EDM_ML_DEBUG
1151   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Minimum/maximum R " << php.rLimit_[0] << ":" << php.rLimit_[1];
1152   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters finds " << php.zLayerHex_.size() << " layers";
1153   for (unsigned int i = 0; i < php.zLayerHex_.size(); ++i) {
1154     int k = php.layerIndex_[i];
1155     edm::LogVerbatim("HGCalGeom") << "Layer[" << i << ":" << k << ":" << php.layer_[k]
1156                                   << "] with r = " << php.rMinLayHex_[i] << ":" << php.rMaxLayHex_[i]
1157                                   << " at z = " << php.zLayerHex_[i];
1158   }
1159   edm::LogVerbatim("HGCalGeom") << "Obtained " << php.trformIndex_.size() << " transformation matrices";
1160   for (unsigned int k = 0; k < php.trformIndex_.size(); ++k) {
1161     edm::LogVerbatim("HGCalGeom") << "Matrix[" << k << "] (" << std::hex << php.trformIndex_[k] << std::dec
1162                                   << ") Translation (" << php.trformTranX_[k] << ", " << php.trformTranY_[k] << ", "
1163                                   << php.trformTranZ_[k] << " Rotation (" << php.trformRotXX_[k] << ", "
1164                                   << php.trformRotYX_[k] << ", " << php.trformRotZX_[k] << ", " << php.trformRotXY_[k]
1165                                   << ", " << php.trformRotYY_[k] << ", " << php.trformRotZY_[k] << ", "
1166                                   << php.trformRotXZ_[k] << ", " << php.trformRotYZ_[k] << ", " << php.trformRotZZ_[k]
1167                                   << ")";
1168   }
1169 #endif
1170 }
1171 
1172 void HGCalGeomParameters::loadSpecParsHexagon(const DDFilteredView& fv,
1173                                               HGCalParameters& php,
1174                                               const DDCompactView* cpv,
1175                                               const std::string& sdTag1,
1176                                               const std::string& sdTag2) {
1177   DDsvalues_type sv(fv.mergedSpecifics());
1178   php.boundR_ = getDDDArray("RadiusBound", sv, 4);
1179   rescale(php.boundR_, HGCalParameters::k_ScaleFromDDD);
1180   php.rLimit_ = getDDDArray("RadiusLimits", sv, 2);
1181   rescale(php.rLimit_, HGCalParameters::k_ScaleFromDDD);
1182   php.levelT_ = dbl_to_int(getDDDArray("LevelTop", sv, 0));
1183 
1184   // Grouping of layers
1185   php.layerGroup_ = dbl_to_int(getDDDArray("GroupingZFine", sv, 0));
1186   php.layerGroupM_ = dbl_to_int(getDDDArray("GroupingZMid", sv, 0));
1187   php.layerGroupO_ = dbl_to_int(getDDDArray("GroupingZOut", sv, 0));
1188   php.slopeMin_ = getDDDArray("Slope", sv, 1);
1189   const auto& dummy2 = getDDDArray("LayerOffset", sv, 0);
1190   if (!dummy2.empty())
1191     php.layerOffset_ = dummy2[0];
1192   else
1193     php.layerOffset_ = 0;
1194 
1195   // Wafer size
1196   std::string attribute = "Volume";
1197   DDSpecificsMatchesValueFilter filter1{DDValue(attribute, sdTag1, 0.0)};
1198   DDFilteredView fv1(*cpv, filter1);
1199   if (fv1.firstChild()) {
1200     DDsvalues_type sv(fv1.mergedSpecifics());
1201     const auto& dummy = getDDDArray("WaferSize", sv, 0);
1202     waferSize_ = dummy[0];
1203   }
1204 
1205   // Cell size
1206   DDSpecificsMatchesValueFilter filter2{DDValue(attribute, sdTag2, 0.0)};
1207   DDFilteredView fv2(*cpv, filter2);
1208   if (fv2.firstChild()) {
1209     DDsvalues_type sv(fv2.mergedSpecifics());
1210     php.cellSize_ = getDDDArray("CellSize", sv, 0);
1211   }
1212 
1213   loadSpecParsHexagon(php);
1214 }
1215 
1216 void HGCalGeomParameters::loadSpecParsHexagon(const cms::DDFilteredView& fv,
1217                                               HGCalParameters& php,
1218                                               const std::string& sdTag1,
1219                                               const std::string& sdTag2,
1220                                               const std::string& sdTag3,
1221                                               const std::string& sdTag4) {
1222   php.boundR_ = fv.get<std::vector<double> >(sdTag4, "RadiusBound");
1223   rescale(php.boundR_, HGCalParameters::k_ScaleFromDD4hep);
1224   php.rLimit_ = fv.get<std::vector<double> >(sdTag4, "RadiusLimits");
1225   rescale(php.rLimit_, HGCalParameters::k_ScaleFromDD4hep);
1226   php.levelT_ = dbl_to_int(fv.get<std::vector<double> >(sdTag4, "LevelTop"));
1227 
1228   // Grouping of layers
1229   php.layerGroup_ = dbl_to_int(fv.get<std::vector<double> >(sdTag1, "GroupingZFine"));
1230   php.layerGroupM_ = dbl_to_int(fv.get<std::vector<double> >(sdTag1, "GroupingZMid"));
1231   php.layerGroupO_ = dbl_to_int(fv.get<std::vector<double> >(sdTag1, "GroupingZOut"));
1232   php.slopeMin_ = fv.get<std::vector<double> >(sdTag4, "Slope");
1233   if (php.slopeMin_.empty())
1234     php.slopeMin_.emplace_back(0);
1235 
1236   // Wafer size
1237   const auto& dummy = fv.get<std::vector<double> >(sdTag2, "WaferSize");
1238   waferSize_ = dummy[0] * HGCalParameters::k_ScaleFromDD4hepToG4;
1239 
1240   // Cell size
1241   php.cellSize_ = fv.get<std::vector<double> >(sdTag3, "CellSize");
1242   rescale(php.cellSize_, HGCalParameters::k_ScaleFromDD4hepToG4);
1243 
1244   // Layer Offset
1245   const auto& dummy2 = fv.get<std::vector<double> >(sdTag1, "LayerOffset");
1246   if (!dummy2.empty()) {
1247     php.layerOffset_ = dummy2[0];
1248   } else {
1249     php.layerOffset_ = 0;
1250   }
1251 
1252   loadSpecParsHexagon(php);
1253 }
1254 
1255 void HGCalGeomParameters::loadSpecParsHexagon(const HGCalParameters& php) {
1256 #ifdef EDM_ML_DEBUG
1257   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: wafer radius ranges"
1258                                 << " for cell grouping " << php.boundR_[0] << ":" << php.boundR_[1] << ":"
1259                                 << php.boundR_[2] << ":" << php.boundR_[3];
1260   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Minimum/maximum R " << php.rLimit_[0] << ":" << php.rLimit_[1];
1261   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: LevelTop " << php.levelT_[0];
1262   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: minimum slope " << php.slopeMin_[0] << " and layer groupings "
1263                                 << "for the 3 ranges:";
1264   for (unsigned int k = 0; k < php.layerGroup_.size(); ++k)
1265     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.layerGroup_[k] << ":" << php.layerGroupM_[k] << ":"
1266                                   << php.layerGroupO_[k];
1267   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Wafer Size: " << waferSize_;
1268   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: " << php.cellSize_.size() << " cells of sizes:";
1269   for (unsigned int k = 0; k < php.cellSize_.size(); ++k)
1270     edm::LogVerbatim("HGCalGeom") << " [" << k << "] " << php.cellSize_[k];
1271   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: First Layer " << php.firstLayer_ << " and layer offset "
1272                                 << php.layerOffset_;
1273 #endif
1274 }
1275 
1276 void HGCalGeomParameters::loadSpecParsHexagon8(const DDFilteredView& fv, HGCalParameters& php) {
1277   DDsvalues_type sv(fv.mergedSpecifics());
1278   php.cellThickness_ = getDDDArray("CellThickness", sv, 3);
1279   rescale(php.cellThickness_, HGCalParameters::k_ScaleFromDDD);
1280   if ((php.mode_ == HGCalGeometryMode::Hexagon8Module) || (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
1281       (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell)) {
1282     php.waferThickness_ = getDDDArray("WaferThickness", sv, 3);
1283     rescale(php.waferThickness_, HGCalParameters::k_ScaleFromDDD);
1284   } else {
1285     for (unsigned int k = 0; k < php.cellThickness_.size(); ++k)
1286       php.waferThickness_.emplace_back(php.waferThick_);
1287   }
1288 
1289   php.radius100to200_ = getDDDArray("Radius100to200", sv, 5);
1290   php.radius200to300_ = getDDDArray("Radius200to300", sv, 5);
1291 
1292   const auto& dummy = getDDDArray("RadiusCuts", sv, 4);
1293   php.choiceType_ = static_cast<int>(dummy[0]);
1294   php.nCornerCut_ = static_cast<int>(dummy[1]);
1295   php.fracAreaMin_ = dummy[2];
1296   php.zMinForRad_ = HGCalParameters::k_ScaleFromDDD * dummy[3];
1297 
1298   php.radiusMixBoundary_ = fv.vector("RadiusMixBoundary");
1299   rescale(php.radiusMixBoundary_, HGCalParameters::k_ScaleFromDDD);
1300 
1301   php.slopeMin_ = getDDDArray("SlopeBottom", sv, 0);
1302   php.zFrontMin_ = getDDDArray("ZFrontBottom", sv, 0);
1303   rescale(php.zFrontMin_, HGCalParameters::k_ScaleFromDDD);
1304   php.rMinFront_ = getDDDArray("RMinFront", sv, 0);
1305   rescale(php.rMinFront_, HGCalParameters::k_ScaleFromDDD);
1306 
1307   php.slopeTop_ = getDDDArray("SlopeTop", sv, 0);
1308   php.zFrontTop_ = getDDDArray("ZFrontTop", sv, 0);
1309   rescale(php.zFrontTop_, HGCalParameters::k_ScaleFromDDD);
1310   php.rMaxFront_ = getDDDArray("RMaxFront", sv, 0);
1311   rescale(php.rMaxFront_, HGCalParameters::k_ScaleFromDDD);
1312 
1313   php.zRanges_ = fv.vector("ZRanges");
1314   rescale(php.zRanges_, HGCalParameters::k_ScaleFromDDD);
1315 
1316   const auto& dummy2 = getDDDArray("LayerOffset", sv, 1);
1317   php.layerOffset_ = dummy2[0];
1318   php.layerCenter_ = dbl_to_int(fv.vector("LayerCenter"));
1319 
1320   if (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell) {
1321     const auto& dummy3 = fv.vector("CalibCellRadius");
1322     php.calibCellRHD_ = HGCalParameters::k_ScaleFromDDD * dummy3[0];
1323     php.calibCellFullHD_ = dbl_to_int(fv.vector("CalibCellFullHD"));
1324     php.calibCellPartHD_ = dbl_to_int(fv.vector("CalibCellPartHD"));
1325     php.calibCellRLD_ = HGCalParameters::k_ScaleFromDDD * dummy3[1];
1326     php.calibCellFullLD_ = dbl_to_int(fv.vector("CalibCellFullLD"));
1327     php.calibCellPartLD_ = dbl_to_int(fv.vector("CalibCellPartLD"));
1328   } else {
1329     php.calibCellRHD_ = php.calibCellRLD_ = 0;
1330     php.calibCellFullHD_.clear();
1331     php.calibCellPartHD_.clear();
1332     php.calibCellFullLD_.clear();
1333     php.calibCellPartLD_.clear();
1334   }
1335 
1336   loadSpecParsHexagon8(php);
1337 
1338   // Read in parameters from Philip's file
1339   if (php.waferMaskMode_ > 1) {
1340     std::vector<int> layerType, waferIndex, waferProperties;
1341     std::vector<double> cassetteShift;
1342     if (php.waferMaskMode_ == siliconFileEE) {
1343       waferIndex = dbl_to_int(fv.vector("WaferIndexEE"));
1344       waferProperties = dbl_to_int(fv.vector("WaferPropertiesEE"));
1345     } else if (php.waferMaskMode_ == siliconCassetteEE) {
1346       waferIndex = dbl_to_int(fv.vector("WaferIndexEE"));
1347       waferProperties = dbl_to_int(fv.vector("WaferPropertiesEE"));
1348       cassetteShift = fv.vector("CassetteShiftEE");
1349     } else if (php.waferMaskMode_ == siliconFileHE) {
1350       waferIndex = dbl_to_int(fv.vector("WaferIndexHE"));
1351       waferProperties = dbl_to_int(fv.vector("WaferPropertiesHE"));
1352     } else if (php.waferMaskMode_ == siliconCassetteHE) {
1353       waferIndex = dbl_to_int(fv.vector("WaferIndexHE"));
1354       waferProperties = dbl_to_int(fv.vector("WaferPropertiesHE"));
1355       cassetteShift = fv.vector("CassetteShiftHE");
1356     }
1357     if ((php.mode_ == HGCalGeometryMode::Hexagon8Module) || (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
1358         (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell)) {
1359       if ((php.waferMaskMode_ == siliconFileEE) || (php.waferMaskMode_ == siliconCassetteEE)) {
1360         layerType = dbl_to_int(fv.vector("LayerTypesEE"));
1361       } else if ((php.waferMaskMode_ == siliconFileHE) || (php.waferMaskMode_ == siliconCassetteHE)) {
1362         layerType = dbl_to_int(fv.vector("LayerTypesHE"));
1363       }
1364     }
1365 
1366     php.cassetteShift_ = cassetteShift;
1367     rescale(php.cassetteShift_, HGCalParameters::k_ScaleFromDDD);
1368     loadSpecParsHexagon8(php, layerType, waferIndex, waferProperties);
1369   }
1370 }
1371 
1372 void HGCalGeomParameters::loadSpecParsHexagon8(const cms::DDFilteredView& fv,
1373                                                const cms::DDVectorsMap& vmap,
1374                                                HGCalParameters& php,
1375                                                const std::string& sdTag1) {
1376   php.cellThickness_ = fv.get<std::vector<double> >(sdTag1, "CellThickness");
1377   rescale(php.cellThickness_, HGCalParameters::k_ScaleFromDD4hep);
1378   if ((php.mode_ == HGCalGeometryMode::Hexagon8Module) || (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
1379       (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell)) {
1380     php.waferThickness_ = fv.get<std::vector<double> >(sdTag1, "WaferThickness");
1381     rescale(php.waferThickness_, HGCalParameters::k_ScaleFromDD4hep);
1382   } else {
1383     for (unsigned int k = 0; k < php.cellThickness_.size(); ++k)
1384       php.waferThickness_.emplace_back(php.waferThick_);
1385   }
1386 
1387   php.radius100to200_ = fv.get<std::vector<double> >(sdTag1, "Radius100to200");
1388   php.radius200to300_ = fv.get<std::vector<double> >(sdTag1, "Radius200to300");
1389 
1390   const auto& dummy = fv.get<std::vector<double> >(sdTag1, "RadiusCuts");
1391   if (dummy.size() > 3) {
1392     php.choiceType_ = static_cast<int>(dummy[0]);
1393     php.nCornerCut_ = static_cast<int>(dummy[1]);
1394     php.fracAreaMin_ = dummy[2];
1395     php.zMinForRad_ = HGCalParameters::k_ScaleFromDD4hep * dummy[3];
1396   } else {
1397     php.choiceType_ = php.nCornerCut_ = php.fracAreaMin_ = php.zMinForRad_ = 0;
1398   }
1399 
1400   php.slopeMin_ = fv.get<std::vector<double> >(sdTag1, "SlopeBottom");
1401   php.zFrontMin_ = fv.get<std::vector<double> >(sdTag1, "ZFrontBottom");
1402   rescale(php.zFrontMin_, HGCalParameters::k_ScaleFromDD4hep);
1403   php.rMinFront_ = fv.get<std::vector<double> >(sdTag1, "RMinFront");
1404   rescale(php.rMinFront_, HGCalParameters::k_ScaleFromDD4hep);
1405 
1406   php.slopeTop_ = fv.get<std::vector<double> >(sdTag1, "SlopeTop");
1407   php.zFrontTop_ = fv.get<std::vector<double> >(sdTag1, "ZFrontTop");
1408   rescale(php.zFrontTop_, HGCalParameters::k_ScaleFromDD4hep);
1409   php.rMaxFront_ = fv.get<std::vector<double> >(sdTag1, "RMaxFront");
1410   rescale(php.rMaxFront_, HGCalParameters::k_ScaleFromDD4hep);
1411   unsigned int kmax = (php.zFrontTop_.size() - php.slopeTop_.size());
1412   for (unsigned int k = 0; k < kmax; ++k)
1413     php.slopeTop_.emplace_back(0);
1414 
1415   const auto& dummy2 = fv.get<std::vector<double> >(sdTag1, "LayerOffset");
1416   if (!dummy2.empty()) {
1417     php.layerOffset_ = dummy2[0];
1418   } else {
1419     php.layerOffset_ = 0;
1420   }
1421 
1422   php.calibCellRHD_ = php.calibCellRLD_ = 0;
1423   php.calibCellFullHD_.clear();
1424   php.calibCellPartHD_.clear();
1425   php.calibCellFullLD_.clear();
1426   php.calibCellPartLD_.clear();
1427   for (auto const& it : vmap) {
1428     if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "RadiusMixBoundary")) {
1429       for (const auto& i : it.second)
1430         php.radiusMixBoundary_.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1431     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "ZRanges")) {
1432       for (const auto& i : it.second)
1433         php.zRanges_.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1434     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "LayerCenter")) {
1435       for (const auto& i : it.second)
1436         php.layerCenter_.emplace_back(std::round(i));
1437     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CalibCellRadius")) {
1438       php.calibCellRHD_ = HGCalParameters::k_ScaleFromDD4hep * it.second[0];
1439       php.calibCellRLD_ = HGCalParameters::k_ScaleFromDDD * it.second[1];
1440     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CalibCellFullHD")) {
1441       for (const auto& i : it.second)
1442         php.calibCellFullHD_.emplace_back(std::round(i));
1443     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CalibCellPartHD")) {
1444       for (const auto& i : it.second)
1445         php.calibCellPartHD_.emplace_back(std::round(i));
1446     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CalibCellFullLD")) {
1447       for (const auto& i : it.second)
1448         php.calibCellFullLD_.emplace_back(std::round(i));
1449     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CalibCellPartLD")) {
1450       for (const auto& i : it.second)
1451         php.calibCellPartLD_.emplace_back(std::round(i));
1452     }
1453   }
1454 
1455   loadSpecParsHexagon8(php);
1456 
1457   // Read in parameters from Philip's file
1458   if (php.waferMaskMode_ > 1) {
1459     std::vector<int> layerType, waferIndex, waferProperties;
1460     std::vector<double> cassetteShift;
1461     if ((php.waferMaskMode_ == siliconFileEE) || (php.waferMaskMode_ == siliconCassetteEE)) {
1462       for (auto const& it : vmap) {
1463         if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "WaferIndexEE")) {
1464           for (const auto& i : it.second)
1465             waferIndex.emplace_back(std::round(i));
1466         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "WaferPropertiesEE")) {
1467           for (const auto& i : it.second)
1468             waferProperties.emplace_back(std::round(i));
1469         }
1470       }
1471       if (php.waferMaskMode_ == siliconCassetteEE) {
1472         for (auto const& it : vmap) {
1473           if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CassetteShiftEE")) {
1474             for (const auto& i : it.second)
1475               cassetteShift.emplace_back(i);
1476           }
1477         }
1478       }
1479     } else if ((php.waferMaskMode_ == siliconFileHE) || (php.waferMaskMode_ == siliconCassetteHE)) {
1480       for (auto const& it : vmap) {
1481         if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "WaferIndexHE")) {
1482           for (const auto& i : it.second)
1483             waferIndex.emplace_back(std::round(i));
1484         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "WaferPropertiesHE")) {
1485           for (const auto& i : it.second)
1486             waferProperties.emplace_back(std::round(i));
1487         }
1488       }
1489       if (php.waferMaskMode_ == siliconCassetteHE) {
1490         for (auto const& it : vmap) {
1491           if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CassetteShiftHE")) {
1492             for (const auto& i : it.second)
1493               cassetteShift.emplace_back(i);
1494           }
1495         }
1496       }
1497     }
1498     if ((php.mode_ == HGCalGeometryMode::Hexagon8Module) || (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
1499         (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell)) {
1500       if ((php.waferMaskMode_ == siliconFileEE) || (php.waferMaskMode_ == siliconCassetteEE)) {
1501         for (auto const& it : vmap) {
1502           if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "LayerTypesEE")) {
1503             for (const auto& i : it.second)
1504               layerType.emplace_back(std::round(i));
1505           }
1506         }
1507       } else if ((php.waferMaskMode_ == siliconFileHE) || (php.waferMaskMode_ == siliconCassetteHE)) {
1508         for (auto const& it : vmap) {
1509           if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "LayerTypesHE")) {
1510             for (const auto& i : it.second)
1511               layerType.emplace_back(std::round(i));
1512           }
1513         }
1514       }
1515     }
1516 
1517     php.cassetteShift_ = cassetteShift;
1518     rescale(php.cassetteShift_, HGCalParameters::k_ScaleFromDD4hep);
1519     loadSpecParsHexagon8(php, layerType, waferIndex, waferProperties);
1520   }
1521 }
1522 
1523 void HGCalGeomParameters::loadSpecParsHexagon8(HGCalParameters& php) {
1524 #ifdef EDM_ML_DEBUG
1525   for (unsigned int k = 0; k < php.waferThickness_.size(); ++k)
1526     edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: wafer[" << k << "] Thickness " << php.waferThickness_[k];
1527   for (unsigned int k = 0; k < php.cellThickness_.size(); ++k)
1528     edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: cell[" << k << "] Thickness " << php.cellThickness_[k];
1529   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Polynomial "
1530                                 << "parameters for 120 to 200 micron "
1531                                 << "transition with " << php.radius100to200_.size() << " elements";
1532   for (unsigned int k = 0; k < php.radius100to200_.size(); ++k)
1533     edm::LogVerbatim("HGCalGeom") << "Element [" << k << "] " << php.radius100to200_[k];
1534   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Polynomial "
1535                                 << "parameters for 200 to 300 micron "
1536                                 << "transition with " << php.radius200to300_.size() << " elements";
1537   for (unsigned int k = 0; k < php.radius200to300_.size(); ++k)
1538     edm::LogVerbatim("HGCalGeom") << "Element [" << k << "] " << php.radius200to300_[k];
1539   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Parameters for the"
1540                                 << " transition " << php.choiceType_ << ":" << php.nCornerCut_ << ":"
1541                                 << php.fracAreaMin_ << ":" << php.zMinForRad_;
1542   for (unsigned int k = 0; k < php.radiusMixBoundary_.size(); ++k)
1543     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Mix[" << k << "] R = " << php.radiusMixBoundary_[k];
1544   for (unsigned int k = 0; k < php.zFrontMin_.size(); ++k)
1545     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Boundary[" << k << "] Bottom Z = " << php.zFrontMin_[k]
1546                                   << " Slope = " << php.slopeMin_[k] << " rMax = " << php.rMinFront_[k];
1547   for (unsigned int k = 0; k < php.zFrontTop_.size(); ++k)
1548     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Boundary[" << k << "] Top Z = " << php.zFrontTop_[k]
1549                                   << " Slope = " << php.slopeTop_[k] << " rMax = " << php.rMaxFront_[k];
1550   std::ostringstream st1;
1551   for (unsigned int k = 0; k < php.zRanges_.size(); ++k)
1552     st1 << ":" << php.zRanges_[k];
1553   edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Z-Boundary[" << php.zRanges_.size() << "] " << st1.str();
1554   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: LayerOffset " << php.layerOffset_ << " in array of size "
1555                                 << php.layerCenter_.size();
1556   for (unsigned int k = 0; k < php.layerCenter_.size(); ++k)
1557     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.layerCenter_[k];
1558   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Calibration cells in HD having radius " << php.calibCellRHD_
1559                                 << " in arrays of size " << php.calibCellFullHD_.size() << ":"
1560                                 << php.calibCellPartHD_.size();
1561   for (unsigned int k = 0; k < php.calibCellFullHD_.size(); ++k)
1562     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.calibCellFullHD_[k] << ":" << php.calibCellPartHD_[k];
1563   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Calibration cells in LD having radius " << php.calibCellRLD_
1564                                 << " in arrays of size " << php.calibCellFullLD_.size() << ":"
1565                                 << php.calibCellPartLD_.size();
1566   for (unsigned int k = 0; k < php.calibCellFullLD_.size(); ++k)
1567     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.calibCellFullLD_[k] << ":" << php.calibCellPartLD_[k];
1568 #endif
1569 }
1570 
1571 void HGCalGeomParameters::loadSpecParsHexagon8(HGCalParameters& php,
1572                                                const std::vector<int>& layerType,
1573                                                const std::vector<int>& waferIndex,
1574                                                const std::vector<int>& waferProperties) {
1575   // Store parameters from Philip's file
1576   for (unsigned int k = 0; k < layerType.size(); ++k) {
1577     php.layerType_.emplace_back(HGCalTypes::layerType(layerType[k]));
1578 #ifdef EDM_ML_DEBUG
1579     edm::LogVerbatim("HGCalGeom") << "Layer[" << k << "] Type " << layerType[k] << ":" << php.layerType_.back();
1580 #endif
1581   }
1582   for (unsigned int k = 0; k < php.layerType_.size(); ++k) {
1583     double cth = (php.layerType_[k] == HGCalTypes::WaferCenterR) ? cos(php.layerRotation_) : 1.0;
1584     double sth = (php.layerType_[k] == HGCalTypes::WaferCenterR) ? sin(php.layerRotation_) : 0.0;
1585     php.layerRotV_.emplace_back(std::make_pair(cth, sth));
1586 #ifdef EDM_ML_DEBUG
1587     edm::LogVerbatim("HGCalGeom") << "Layer[" << k << "] Type " << php.layerType_[k] << " cos|sin(Theta) "
1588                                   << php.layerRotV_.back().first << ":" << php.layerRotV_.back().second;
1589 #endif
1590   }
1591   for (unsigned int k = 0; k < waferIndex.size(); ++k) {
1592     int partial = HGCalProperty::waferPartial(waferProperties[k]);
1593     int orient =
1594         ((php.mode_ == HGCalGeometryMode::Hexagon8Cassette) || (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell))
1595             ? HGCalProperty::waferOrient(waferProperties[k])
1596             : HGCalWaferMask::getRotation(php.waferZSide_, partial, HGCalProperty::waferOrient(waferProperties[k]));
1597     php.waferInfoMap_[waferIndex[k]] = HGCalParameters::waferInfo(HGCalProperty::waferThick(waferProperties[k]),
1598                                                                   partial,
1599                                                                   orient,
1600                                                                   HGCalProperty::waferCassette(waferProperties[k]));
1601 #ifdef EDM_ML_DEBUG
1602     edm::LogVerbatim("HGCalGeom") << "[" << k << ":" << waferIndex[k] << ":"
1603                                   << HGCalWaferIndex::waferLayer(waferIndex[k]) << ":"
1604                                   << HGCalWaferIndex::waferU(waferIndex[k]) << ":"
1605                                   << HGCalWaferIndex::waferV(waferIndex[k]) << "]  Thickness type "
1606                                   << HGCalProperty::waferThick(waferProperties[k]) << " Partial type " << partial
1607                                   << " Orientation " << HGCalProperty::waferOrient(waferProperties[k]) << ":" << orient;
1608 #endif
1609   }
1610 }
1611 
1612 void HGCalGeomParameters::loadSpecParsTrapezoid(const DDFilteredView& fv, HGCalParameters& php) {
1613   DDsvalues_type sv(fv.mergedSpecifics());
1614   php.radiusMixBoundary_ = fv.vector("RadiusMixBoundary");
1615   rescale(php.radiusMixBoundary_, HGCalParameters::k_ScaleFromDDD);
1616 
1617   php.nPhiBinBH_ = dbl_to_int(getDDDArray("NPhiBinBH", sv, 0));
1618   php.layerFrontBH_ = dbl_to_int(getDDDArray("LayerFrontBH", sv, 0));
1619   php.rMinLayerBH_ = getDDDArray("RMinLayerBH", sv, 0);
1620   rescale(php.rMinLayerBH_, HGCalParameters::k_ScaleFromDDD);
1621   assert(php.nPhiBinBH_.size() > 1);
1622   php.nCellsFine_ = php.nPhiBinBH_[0];
1623   php.nCellsCoarse_ = php.nPhiBinBH_[1];
1624   assert(0 != php.nCellsFine_);
1625   assert(0 != php.nCellsCoarse_);
1626   php.cellSize_.emplace_back(2.0 * M_PI / php.nCellsFine_);
1627   php.cellSize_.emplace_back(2.0 * M_PI / php.nCellsCoarse_);
1628 
1629   php.slopeMin_ = getDDDArray("SlopeBottom", sv, 0);
1630   php.zFrontMin_ = getDDDArray("ZFrontBottom", sv, 0);
1631   rescale(php.zFrontMin_, HGCalParameters::k_ScaleFromDDD);
1632   php.rMinFront_ = getDDDArray("RMinFront", sv, 0);
1633   rescale(php.rMinFront_, HGCalParameters::k_ScaleFromDDD);
1634 
1635   php.slopeTop_ = getDDDArray("SlopeTop", sv, 0);
1636   php.zFrontTop_ = getDDDArray("ZFrontTop", sv, 0);
1637   rescale(php.zFrontTop_, HGCalParameters::k_ScaleFromDDD);
1638   php.rMaxFront_ = getDDDArray("RMaxFront", sv, 0);
1639   rescale(php.rMaxFront_, HGCalParameters::k_ScaleFromDDD);
1640 
1641   php.zRanges_ = fv.vector("ZRanges");
1642   rescale(php.zRanges_, HGCalParameters::k_ScaleFromDDD);
1643 
1644   // Offsets
1645   const auto& dummy2 = getDDDArray("LayerOffset", sv, 1);
1646   php.layerOffset_ = dummy2[0];
1647   php.layerCenter_ = dbl_to_int(fv.vector("LayerCenter"));
1648 
1649   loadSpecParsTrapezoid(php);
1650 #ifdef EDM_ML_DEBUG
1651   edm::LogVerbatim("HGCalGeom") << "WaferMaskMode " << php.waferMaskMode_ << " Compare " << scintillatorFile << ":"
1652                                 << scintillatorCassette << ":" << scintillatorFineCell;
1653 #endif
1654   // tile parameters from Katja's file
1655   if ((php.waferMaskMode_ == scintillatorFile) || (php.waferMaskMode_ == scintillatorCassette) ||
1656       (php.waferMaskMode_ == scintillatorFineCell)) {
1657     std::vector<int> tileIndx, tileProperty;
1658     std::vector<int> tileHEX1, tileHEX2, tileHEX3, tileHEX4, tileHEX5, tileHEX6;
1659     std::vector<double> tileRMin, tileRMax, tileRMinFine, tileRMaxFine;
1660     std::vector<int> tileRingMin, tileRingMax, tileRingMinFine, tileRingMaxFine;
1661     std::vector<double> cassetteShift;
1662     php.nPhiLayer_ = dbl_to_int(fv.vector("NPhiLayer"));
1663     tileIndx = dbl_to_int(fv.vector("TileIndex"));
1664     tileProperty = dbl_to_int(fv.vector("TileProperty"));
1665     tileHEX1 = dbl_to_int(fv.vector("TileHEX1"));
1666     tileHEX2 = dbl_to_int(fv.vector("TileHEX2"));
1667     tileHEX3 = dbl_to_int(fv.vector("TileHEX3"));
1668     tileHEX4 = dbl_to_int(fv.vector("TileHEX4"));
1669     tileRMin = fv.vector("TileRMin");
1670     tileRMax = fv.vector("TileRMax");
1671     rescale(tileRMin, HGCalParameters::k_ScaleFromDDD);
1672     rescale(tileRMax, HGCalParameters::k_ScaleFromDDD);
1673     tileRingMin = dbl_to_int(fv.vector("TileRingMin"));
1674     tileRingMax = dbl_to_int(fv.vector("TileRingMax"));
1675     if (php.waferMaskMode_ == scintillatorFineCell) {
1676       tileHEX5 = dbl_to_int(fv.vector("TileHEX5"));
1677       tileHEX6 = dbl_to_int(fv.vector("TileHEX6"));
1678       tileRMinFine = fv.vector("TileRMin6");
1679       tileRMaxFine = fv.vector("TileRMax6");
1680       rescale(tileRMinFine, HGCalParameters::k_ScaleFromDDD);
1681       rescale(tileRMaxFine, HGCalParameters::k_ScaleFromDDD);
1682       tileRingMinFine = dbl_to_int(fv.vector("TileRingMin6"));
1683       tileRingMaxFine = dbl_to_int(fv.vector("TileRingMax6"));
1684       php.nphiCassette_ = php.nCellsCoarse_ / php.cassettes_;
1685       php.nphiFineCassette_ = php.nCellsFine_ / php.cassettes_;
1686       std::vector<double> rectract = fv.vector("ScintRetract");
1687       rescale(rectract, HGCalParameters::k_ScaleFromDDD);
1688       double dphi = M_PI / php.cassettes_;
1689       for (int k = 0; k < php.cassettes_; ++k) {
1690         double phi = (2 * k + 1) * dphi;
1691         cassetteShift.emplace_back(rectract[k] * cos(phi));
1692         cassetteShift.emplace_back(rectract[k] * sin(phi));
1693       }
1694     } else if (php.waferMaskMode_ == scintillatorCassette) {
1695       if (php.cassettes_ > 0)
1696         php.nphiCassette_ = php.nCellsCoarse_ / php.cassettes_;
1697       cassetteShift = fv.vector("CassetteShiftHE");
1698     }
1699 
1700     rescale(cassetteShift, HGCalParameters::k_ScaleFromDDD);
1701     if (php.waferMaskMode_ == scintillatorFineCell)
1702       php.cassetteShiftTile_ = cassetteShift;
1703     else
1704       php.cassetteShift_ = cassetteShift;
1705     loadSpecParsTrapezoid(php,
1706                           tileIndx,
1707                           tileProperty,
1708                           tileHEX1,
1709                           tileHEX2,
1710                           tileHEX3,
1711                           tileHEX4,
1712                           tileHEX5,
1713                           tileHEX6,
1714                           tileRMin,
1715                           tileRMax,
1716                           tileRMinFine,
1717                           tileRMaxFine,
1718                           tileRingMin,
1719                           tileRingMax,
1720                           tileRingMinFine,
1721                           tileRingMaxFine);
1722   }
1723 }
1724 
1725 void HGCalGeomParameters::loadSpecParsTrapezoid(const cms::DDFilteredView& fv,
1726                                                 const cms::DDVectorsMap& vmap,
1727                                                 HGCalParameters& php,
1728                                                 const std::string& sdTag1) {
1729   for (auto const& it : vmap) {
1730     if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "RadiusMixBoundary")) {
1731       for (const auto& i : it.second)
1732         php.radiusMixBoundary_.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1733     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "ZRanges")) {
1734       for (const auto& i : it.second)
1735         php.zRanges_.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1736     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "LayerCenter")) {
1737       for (const auto& i : it.second)
1738         php.layerCenter_.emplace_back(std::round(i));
1739     }
1740   }
1741 
1742   php.nPhiBinBH_ = dbl_to_int(fv.get<std::vector<double> >(sdTag1, "NPhiBinBH"));
1743   php.layerFrontBH_ = dbl_to_int(fv.get<std::vector<double> >(sdTag1, "LayerFrontBH"));
1744   php.rMinLayerBH_ = fv.get<std::vector<double> >(sdTag1, "RMinLayerBH");
1745   rescale(php.rMinLayerBH_, HGCalParameters::k_ScaleFromDD4hep);
1746   assert(php.nPhiBinBH_.size() > 1);
1747   php.nCellsFine_ = php.nPhiBinBH_[0];
1748   php.nCellsCoarse_ = php.nPhiBinBH_[1];
1749   assert(0 != php.nCellsFine_);
1750   assert(0 != php.nCellsCoarse_);
1751   php.cellSize_.emplace_back(2.0 * M_PI / php.nCellsFine_);
1752   php.cellSize_.emplace_back(2.0 * M_PI / php.nCellsCoarse_);
1753 
1754   php.slopeMin_ = fv.get<std::vector<double> >(sdTag1, "SlopeBottom");
1755   php.zFrontMin_ = fv.get<std::vector<double> >(sdTag1, "ZFrontBottom");
1756   rescale(php.zFrontMin_, HGCalParameters::k_ScaleFromDD4hep);
1757   php.rMinFront_ = fv.get<std::vector<double> >(sdTag1, "RMinFront");
1758   rescale(php.rMinFront_, HGCalParameters::k_ScaleFromDD4hep);
1759 
1760   php.slopeTop_ = fv.get<std::vector<double> >(sdTag1, "SlopeTop");
1761   php.zFrontTop_ = fv.get<std::vector<double> >(sdTag1, "ZFrontTop");
1762   rescale(php.zFrontTop_, HGCalParameters::k_ScaleFromDD4hep);
1763   php.rMaxFront_ = fv.get<std::vector<double> >(sdTag1, "RMaxFront");
1764   rescale(php.rMaxFront_, HGCalParameters::k_ScaleFromDD4hep);
1765   unsigned int kmax = (php.zFrontTop_.size() - php.slopeTop_.size());
1766   for (unsigned int k = 0; k < kmax; ++k)
1767     php.slopeTop_.emplace_back(0);
1768 
1769   const auto& dummy2 = fv.get<std::vector<double> >(sdTag1, "LayerOffset");
1770   php.layerOffset_ = dummy2[0];
1771 
1772   loadSpecParsTrapezoid(php);
1773 
1774   // tile parameters from Katja's file
1775   if ((php.waferMaskMode_ == scintillatorFile) || (php.waferMaskMode_ == scintillatorCassette) ||
1776       (php.waferMaskMode_ == scintillatorFineCell)) {
1777     std::vector<int> tileIndx, tileProperty;
1778     std::vector<int> tileHEX1, tileHEX2, tileHEX3, tileHEX4, tileHEX5, tileHEX6;
1779     std::vector<double> tileRMin, tileRMax, tileRMinFine, tileRMaxFine;
1780     std::vector<int> tileRingMin, tileRingMax, tileRingMinFine, tileRingMaxFine;
1781     std::vector<double> cassetteShift;
1782     for (auto const& it : vmap) {
1783       if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileIndex")) {
1784         for (const auto& i : it.second)
1785           tileIndx.emplace_back(std::round(i));
1786       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileProperty")) {
1787         for (const auto& i : it.second)
1788           tileProperty.emplace_back(std::round(i));
1789       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "NPhiLayer")) {
1790         for (const auto& i : it.second)
1791           php.nPhiLayer_.emplace_back(std::round(i));
1792       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX1")) {
1793         for (const auto& i : it.second)
1794           tileHEX1.emplace_back(std::round(i));
1795       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX2")) {
1796         for (const auto& i : it.second)
1797           tileHEX2.emplace_back(std::round(i));
1798       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX3")) {
1799         for (const auto& i : it.second)
1800           tileHEX3.emplace_back(std::round(i));
1801       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX4")) {
1802         for (const auto& i : it.second)
1803           tileHEX4.emplace_back(std::round(i));
1804       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRMin")) {
1805         for (const auto& i : it.second)
1806           tileRMin.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1807       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRMax")) {
1808         for (const auto& i : it.second)
1809           tileRMax.emplace_back(HGCalParameters::k_ScaleFromDD4hep * i);
1810       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRingMin")) {
1811         for (const auto& i : it.second)
1812           tileRingMin.emplace_back(std::round(i));
1813       } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRingMax")) {
1814         for (const auto& i : it.second)
1815           tileRingMax.emplace_back(std::round(i));
1816       }
1817     }
1818     if (php.waferMaskMode_ == scintillatorFineCell) {
1819       php.nphiCassette_ = php.nCellsCoarse_ / php.cassettes_;
1820       php.nphiFineCassette_ = php.nCellsFine_ / php.cassettes_;
1821       std::vector<double> rectract;
1822       for (auto const& it : vmap) {
1823         if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX5")) {
1824           for (const auto& i : it.second)
1825             tileHEX5.emplace_back(std::round(i));
1826         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileHEX6")) {
1827           for (const auto& i : it.second)
1828             tileHEX6.emplace_back(std::round(i));
1829         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRMin6")) {
1830           for (const auto& i : it.second)
1831             tileRMinFine.emplace_back(i);
1832           rescale(tileRMinFine, HGCalParameters::k_ScaleFromDDD);
1833         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRMax6")) {
1834           for (const auto& i : it.second)
1835             tileRMaxFine.emplace_back(i);
1836           rescale(tileRMaxFine, HGCalParameters::k_ScaleFromDDD);
1837         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRingMin6")) {
1838           for (const auto& i : it.second)
1839             tileRingMinFine.emplace_back(std::round(i));
1840         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "TileRingMax6")) {
1841           for (const auto& i : it.second)
1842             tileRingMaxFine.emplace_back(std::round(i));
1843         } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "ScintRetract")) {
1844           for (const auto& i : it.second)
1845             rectract.emplace_back(i);
1846           double dphi = M_PI / php.cassettes_;
1847           for (int k = 0; k < php.cassettes_; ++k) {
1848             double phi = (2 * k + 1) * dphi;
1849             cassetteShift.emplace_back(rectract[k] * cos(phi));
1850             cassetteShift.emplace_back(rectract[k] * sin(phi));
1851           }
1852         }
1853       }
1854     } else if (php.waferMaskMode_ == scintillatorCassette) {
1855       for (auto const& it : vmap) {
1856         if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "CassetteShiftHE")) {
1857           for (const auto& i : it.second)
1858             cassetteShift.emplace_back(i);
1859         }
1860       }
1861     }
1862 
1863     rescale(cassetteShift, HGCalParameters::k_ScaleFromDD4hep);
1864     if (php.waferMaskMode_ == scintillatorFineCell)
1865       php.cassetteShiftTile_ = cassetteShift;
1866     else
1867       php.cassetteShift_ = cassetteShift;
1868     loadSpecParsTrapezoid(php,
1869                           tileIndx,
1870                           tileProperty,
1871                           tileHEX1,
1872                           tileHEX2,
1873                           tileHEX3,
1874                           tileHEX4,
1875                           tileHEX5,
1876                           tileHEX6,
1877                           tileRMin,
1878                           tileRMax,
1879                           tileRMinFine,
1880                           tileRMaxFine,
1881                           tileRingMin,
1882                           tileRingMax,
1883                           tileRingMinFine,
1884                           tileRingMaxFine);
1885   }
1886 }
1887 
1888 void HGCalGeomParameters::loadSpecParsTrapezoid(HGCalParameters& php) {
1889 #ifdef EDM_ML_DEBUG
1890   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters:nCells " << php.nCellsFine_ << ":" << php.nCellsCoarse_
1891                                 << " cellSize: " << php.cellSize_[0] << ":" << php.cellSize_[1];
1892   for (unsigned int k = 0; k < php.layerFrontBH_.size(); ++k)
1893     edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Type[" << k << "] Front Layer = " << php.layerFrontBH_[k]
1894                                   << " rMin = " << php.rMinLayerBH_[k];
1895   for (unsigned int k = 0; k < php.radiusMixBoundary_.size(); ++k) {
1896     edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Mix[" << k << "] R = " << php.radiusMixBoundary_[k]
1897                                   << " Nphi = " << php.scintCells(k + php.firstLayer_)
1898                                   << " dPhi = " << php.scintCellSize(k + php.firstLayer_);
1899   }
1900 
1901   for (unsigned int k = 0; k < php.zFrontMin_.size(); ++k)
1902     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Boundary[" << k << "] Bottom Z = " << php.zFrontMin_[k]
1903                                   << " Slope = " << php.slopeMin_[k] << " rMax = " << php.rMinFront_[k];
1904 
1905   for (unsigned int k = 0; k < php.zFrontTop_.size(); ++k)
1906     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Boundary[" << k << "] Top Z = " << php.zFrontTop_[k]
1907                                   << " Slope = " << php.slopeTop_[k] << " rMax = " << php.rMaxFront_[k];
1908 
1909   std::ostringstream st1;
1910   for (unsigned int k = 0; k < php.zRanges_.size(); ++k)
1911     st1 << ":" << php.zRanges_[k];
1912   edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Z-Boundary[" << php.zRanges_.size() << "] " << st1.str();
1913 
1914   edm::LogVerbatim("HGCalGeom") << "HGCalParameters: LayerOffset " << php.layerOffset_ << " in array of size "
1915                                 << php.layerCenter_.size();
1916   for (unsigned int k = 0; k < php.layerCenter_.size(); ++k)
1917     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.layerCenter_[k];
1918 #endif
1919 }
1920 
1921 void HGCalGeomParameters::loadSpecParsTrapezoid(HGCalParameters& php,
1922                                                 const std::vector<int>& tileIndx,
1923                                                 const std::vector<int>& tileProperty,
1924                                                 const std::vector<int>& tileHEX1,
1925                                                 const std::vector<int>& tileHEX2,
1926                                                 const std::vector<int>& tileHEX3,
1927                                                 const std::vector<int>& tileHEX4,
1928                                                 const std::vector<int>& tileHEX5,
1929                                                 const std::vector<int>& tileHEX6,
1930                                                 const std::vector<double>& tileRMin,
1931                                                 const std::vector<double>& tileRMax,
1932                                                 const std::vector<double>& tileRMinFine,
1933                                                 const std::vector<double>& tileRMaxFine,
1934                                                 const std::vector<int>& tileRingMin,
1935                                                 const std::vector<int>& tileRingMax,
1936                                                 const std::vector<int>& tileRingMinFine,
1937                                                 const std::vector<int>& tileRingMaxFine) {
1938   // tile parameters from Katja's file
1939   for (unsigned int k = 0; k < tileIndx.size(); ++k) {
1940     int hex5 = (k < tileHEX5.size()) ? tileHEX5[k] : 0;
1941     int hex6 = (k < tileHEX6.size()) ? tileHEX6[k] : 0;
1942     php.tileInfoMap_[tileIndx[k]] = HGCalParameters::tileInfo(HGCalProperty::tileType(tileProperty[k]),
1943                                                               HGCalProperty::tileSiPM(tileProperty[k]),
1944                                                               tileHEX1[k],
1945                                                               tileHEX2[k],
1946                                                               tileHEX3[k],
1947                                                               tileHEX4[k],
1948                                                               hex5,
1949                                                               hex6);
1950 #ifdef EDM_ML_DEBUG
1951     edm::LogVerbatim("HGCalGeom") << "Tile[" << k << ":" << tileIndx[k] << "] "
1952                                   << " Type " << HGCalProperty::tileType(tileProperty[k]) << " SiPM "
1953                                   << HGCalProperty::tileSiPM(tileProperty[k]) << " HEX " << std::hex << tileHEX1[k]
1954                                   << ":" << tileHEX2[k] << ":" << tileHEX3[k] << ":" << tileHEX4[k] << ":" << hex5
1955                                   << ":" << hex6 << std::dec;
1956 #endif
1957   }
1958 
1959   for (unsigned int k = 0; k < tileRMinFine.size(); ++k) {
1960     php.tileRingFineR_.emplace_back(tileRMinFine[k], tileRMaxFine[k]);
1961 #ifdef EDM_ML_DEBUG
1962     edm::LogVerbatim("HGCalGeom") << "TileRingFineR[" << k << "] " << tileRMinFine[k] << ":" << tileRMaxFine[k];
1963 #endif
1964   }
1965   for (unsigned int k = 0; k < tileRMin.size(); ++k) {
1966     php.tileRingR_.emplace_back(tileRMin[k], tileRMax[k]);
1967 #ifdef EDM_ML_DEBUG
1968     edm::LogVerbatim("HGCalGeom") << "TileRingR[" << k << "] " << tileRMin[k] << ":" << tileRMax[k];
1969 #endif
1970   }
1971 
1972   for (unsigned int k = 0; k < tileRingMinFine.size(); ++k) {
1973     php.tileRingFineRange_.emplace_back(tileRingMinFine[k], tileRingMaxFine[k]);
1974 #ifdef EDM_ML_DEBUG
1975     edm::LogVerbatim("HGCalGeom") << "TileRingFineRange[" << k << "] " << tileRingMinFine[k] << ":"
1976                                   << tileRingMaxFine[k];
1977 #endif
1978   }
1979   for (unsigned k = 0; k < tileRingMin.size(); ++k) {
1980     php.tileRingRange_.emplace_back(tileRingMin[k], tileRingMax[k]);
1981 #ifdef EDM_ML_DEBUG
1982     edm::LogVerbatim("HGCalGeom") << "TileRingRange[" << k << "] " << tileRingMin[k] << ":" << tileRingMax[k];
1983 #endif
1984   }
1985 #ifdef EDM_ML_DEBUG
1986   std::ostringstream st1, st2;
1987   st1 << "Scintillator CassetteShift for " << php.cassetteShiftTile_.size() << " tiles:";
1988   for (unsigned int k = 0; k < php.cassetteShiftTile_.size(); ++k)
1989     st1 << " " << php.cassetteShiftTile_[k];
1990   edm::LogVerbatim("HGCalGeom") << st1.str();
1991   st2 << "NPhi for scntillator " << php.nPhiLayer_.size() << " layers:";
1992   for (unsigned int k = 0; k < php.nPhiLayer_.size(); ++k)
1993     st2 << " " << php.nPhiLayer_[k];
1994   edm::LogVerbatim("HGCalGeom") << st2.str();
1995 #endif
1996 }
1997 
1998 void HGCalGeomParameters::loadWaferHexagon(HGCalParameters& php) {
1999   double waferW(HGCalParameters::k_ScaleFromDDD * waferSize_), rmin(HGCalParameters::k_ScaleFromDDD * php.waferR_);
2000   double rin(php.rLimit_[0]), rout(php.rLimit_[1]), rMaxFine(php.boundR_[1]);
2001 #ifdef EDM_ML_DEBUG
2002   edm::LogVerbatim("HGCalGeom") << "Input waferWidth " << waferW << ":" << rmin << " R Limits: " << rin << ":" << rout
2003                                 << " Fine " << rMaxFine;
2004 #endif
2005   // Clear the vectors
2006   php.waferCopy_.clear();
2007   php.waferTypeL_.clear();
2008   php.waferTypeT_.clear();
2009   php.waferPosX_.clear();
2010   php.waferPosY_.clear();
2011   double dx = 0.5 * waferW;
2012   double dy = 3.0 * dx * tan(30._deg);
2013   double rr = 2.0 * dx * tan(30._deg);
2014   int ncol = static_cast<int>(2.0 * rout / waferW) + 1;
2015   int nrow = static_cast<int>(rout / (waferW * tan(30._deg))) + 1;
2016   int ns2 = (2 * ncol + 1) * (2 * nrow + 1) * php.layer_.size();
2017   int incm(0), inrm(0);
2018   HGCalParameters::layer_map copiesInLayers(php.layer_.size() + 1);
2019   HGCalParameters::waferT_map waferTypes(ns2 + 1);
2020 #ifdef EDM_ML_DEBUG
2021   int kount(0), ntot(0);
2022   edm::LogVerbatim("HGCalGeom") << "Row " << nrow << " Column " << ncol;
2023 #endif
2024   for (int nr = -nrow; nr <= nrow; ++nr) {
2025     int inr = (nr >= 0) ? nr : -nr;
2026     for (int nc = -ncol; nc <= ncol; ++nc) {
2027       int inc = (nc >= 0) ? nc : -nc;
2028       if (inr % 2 == inc % 2) {
2029         double xpos = nc * dx;
2030         double ypos = nr * dy;
2031         std::pair<int, int> corner = HGCalGeomTools::waferCorner(xpos, ypos, dx, rr, rin, rout, true);
2032         double rpos = std::sqrt(xpos * xpos + ypos * ypos);
2033         int typet = (rpos < rMaxFine) ? 1 : 2;
2034         int typel(3);
2035         for (int k = 1; k < 4; ++k) {
2036           if ((rpos + rmin) <= php.boundR_[k]) {
2037             typel = k;
2038             break;
2039           }
2040         }
2041 #ifdef EDM_ML_DEBUG
2042         ++ntot;
2043 #endif
2044         if (corner.first > 0) {
2045           int copy = HGCalTypes::packTypeUV(typel, nc, nr);
2046           if (inc > incm)
2047             incm = inc;
2048           if (inr > inrm)
2049             inrm = inr;
2050 #ifdef EDM_ML_DEBUG
2051           kount++;
2052           edm::LogVerbatim("HGCalGeom") << kount << ":" << ntot << " Copy " << copy << " Type " << typel << ":" << typet
2053                                         << " Location " << corner.first << " Position " << xpos << ":" << ypos
2054                                         << " Layers " << php.layer_.size();
2055 #endif
2056           php.waferCopy_.emplace_back(copy);
2057           php.waferTypeL_.emplace_back(typel);
2058           php.waferTypeT_.emplace_back(typet);
2059           php.waferPosX_.emplace_back(xpos);
2060           php.waferPosY_.emplace_back(ypos);
2061           for (unsigned int il = 0; il < php.layer_.size(); ++il) {
2062             std::pair<int, int> corner =
2063                 HGCalGeomTools::waferCorner(xpos, ypos, dx, rr, php.rMinLayHex_[il], php.rMaxLayHex_[il], true);
2064             if (corner.first > 0) {
2065               auto cpy = copiesInLayers[php.layer_[il]].find(copy);
2066               if (cpy == copiesInLayers[php.layer_[il]].end())
2067                 copiesInLayers[php.layer_[il]][copy] =
2068                     ((corner.first == static_cast<int>(HGCalParameters::k_CornerSize)) ? php.waferCopy_.size() : -1);
2069             }
2070             if ((corner.first > 0) && (corner.first < static_cast<int>(HGCalParameters::k_CornerSize))) {
2071               int wl = HGCalWaferIndex::waferIndex(php.layer_[il], copy, 0, true);
2072               waferTypes[wl] = corner;
2073             }
2074           }
2075         }
2076       }
2077     }
2078   }
2079   php.copiesInLayers_ = copiesInLayers;
2080   php.waferTypes_ = waferTypes;
2081   php.nSectors_ = static_cast<int>(php.waferCopy_.size());
2082   php.waferUVMax_ = 0;
2083 #ifdef EDM_ML_DEBUG
2084   edm::LogVerbatim("HGCalGeom") << "HGCalWaferHexagon: # of columns " << incm << " # of rows " << inrm << " and "
2085                                 << kount << ":" << ntot << " wafers; R " << rin << ":" << rout;
2086   edm::LogVerbatim("HGCalGeom") << "Dump copiesInLayers for " << php.copiesInLayers_.size() << " layers";
2087   for (unsigned int k = 0; k < copiesInLayers.size(); ++k) {
2088     const auto& theModules = copiesInLayers[k];
2089     edm::LogVerbatim("HGCalGeom") << "Layer " << k << ":" << theModules.size();
2090     int k2(0);
2091     for (std::unordered_map<int, int>::const_iterator itr = theModules.begin(); itr != theModules.end(); ++itr, ++k2) {
2092       edm::LogVerbatim("HGCalGeom") << "[" << k2 << "] " << itr->first << ":" << itr->second;
2093     }
2094   }
2095 #endif
2096 }
2097 
2098 void HGCalGeomParameters::loadWaferHexagon8(HGCalParameters& php) {
2099   double waferW(php.waferSize_);
2100   double waferS(php.sensorSeparation_);
2101   auto wType = std::make_unique<HGCalWaferType>(php.radius100to200_,
2102                                                 php.radius200to300_,
2103                                                 HGCalParameters::k_ScaleToDDD * (waferW + waferS),
2104                                                 HGCalParameters::k_ScaleToDDD * php.zMinForRad_,
2105                                                 php.choiceType_,
2106                                                 php.nCornerCut_,
2107                                                 php.fracAreaMin_);
2108 
2109   double rout(php.rLimit_[1]);
2110 #ifdef EDM_ML_DEBUG
2111   edm::LogVerbatim("HGCalGeom") << "Input waferWidth " << waferW << ":" << waferS << " R Max: " << rout;
2112 #endif
2113   // Clear the vectors
2114   php.waferCopy_.clear();
2115   php.waferTypeL_.clear();
2116   php.waferTypeT_.clear();
2117   php.waferPosX_.clear();
2118   php.waferPosY_.clear();
2119   double r = 0.5 * (waferW + waferS);
2120   double R = 2.0 * r / sqrt3_;
2121   double dy = 0.75 * R;
2122   double r1 = 0.5 * waferW;
2123   double R1 = 2.0 * r1 / sqrt3_;
2124   int N = (r == 0) ? 3 : (static_cast<int>(0.5 * rout / r) + 3);
2125   int ns1 = (2 * N + 1) * (2 * N + 1);
2126   int ns2 = ns1 * php.zLayerHex_.size();
2127 #ifdef EDM_ML_DEBUG
2128   edm::LogVerbatim("HGCalGeom") << "wafer " << waferW << ":" << waferS << " r " << r << " dy " << dy << " N " << N
2129                                 << " sizes " << ns1 << ":" << ns2;
2130   std::vector<int> indtypes(ns1 + 1);
2131   indtypes.clear();
2132 #endif
2133   HGCalParameters::wafer_map wafersInLayers(ns1 + 1);
2134   HGCalParameters::wafer_map typesInLayers(ns2 + 1);
2135   HGCalParameters::waferT_map waferTypes(ns2 + 1);
2136   int ipos(0), lpos(0), uvmax(0), nwarn(0);
2137   std::vector<int> uvmx(php.zLayerHex_.size(), 0);
2138   for (int v = -N; v <= N; ++v) {
2139     for (int u = -N; u <= N; ++u) {
2140       int nr = 2 * v;
2141       int nc = -2 * u + v;
2142       double xpos = nc * r;
2143       double ypos = nr * dy;
2144       int indx = HGCalWaferIndex::waferIndex(0, u, v);
2145       php.waferCopy_.emplace_back(indx);
2146       php.waferPosX_.emplace_back(xpos);
2147       php.waferPosY_.emplace_back(ypos);
2148       wafersInLayers[indx] = ipos;
2149       ++ipos;
2150       std::pair<int, int> corner = HGCalGeomTools::waferCorner(xpos, ypos, r1, R1, 0, rout, false);
2151       if ((corner.first == static_cast<int>(HGCalParameters::k_CornerSize)) ||
2152           ((corner.first > 0) && php.defineFull_)) {
2153         uvmax = std::max(uvmax, std::max(std::abs(u), std::abs(v)));
2154       }
2155       for (unsigned int i = 0; i < php.zLayerHex_.size(); ++i) {
2156         int copy = i + php.layerOffset_;
2157         std::pair<double, double> xyoff = geomTools_.shiftXY(php.layerCenter_[copy], (waferW + waferS));
2158         int lay = php.layer_[php.layerIndex_[i]];
2159         double xpos0 = xpos + xyoff.first;
2160         double ypos0 = ypos + xyoff.second;
2161         double zpos = php.zLayerHex_[i];
2162         int kndx = HGCalWaferIndex::waferIndex(lay, u, v);
2163         int type(-1);
2164         if ((php.mode_ == HGCalGeometryMode::Hexagon8File) || (php.mode_ == HGCalGeometryMode::Hexagon8Module) ||
2165             (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) || (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell))
2166           type = wType->getType(kndx, php.waferInfoMap_);
2167         if (type < 0)
2168           type = wType->getType(HGCalParameters::k_ScaleToDDD * xpos0,
2169                                 HGCalParameters::k_ScaleToDDD * ypos0,
2170                                 HGCalParameters::k_ScaleToDDD * zpos);
2171         php.waferTypeL_.emplace_back(type);
2172         typesInLayers[kndx] = lpos;
2173         ++lpos;
2174 #ifdef EDM_ML_DEBUG
2175         indtypes.emplace_back(kndx);
2176 #endif
2177         std::pair<int, int> corner =
2178             HGCalGeomTools::waferCorner(xpos0, ypos0, r1, R1, php.rMinLayHex_[i], php.rMaxLayHex_[i], false);
2179 #ifdef EDM_ML_DEBUG
2180         if (((corner.first == 0) && std::abs(u) < 5 && std::abs(v) < 5) || (std::abs(u) < 2 && std::abs(v) < 2)) {
2181           edm::LogVerbatim("HGCalGeom") << "Layer " << lay << " R " << php.rMinLayHex_[i] << ":" << php.rMaxLayHex_[i]
2182                                         << " u " << u << " v " << v << " with " << corner.first << " corners";
2183         }
2184 #endif
2185         if ((corner.first == static_cast<int>(HGCalParameters::k_CornerSize)) ||
2186             ((corner.first > 0) && php.defineFull_)) {
2187           uvmx[i] = std::max(uvmx[i], std::max(std::abs(u), std::abs(v)));
2188         }
2189         if ((corner.first < static_cast<int>(HGCalParameters::k_CornerSize)) && (corner.first > 0)) {
2190 #ifdef EDM_ML_DEBUG
2191           edm::LogVerbatim("HGCalGeom") << "Layer " << lay << " u|v " << u << ":" << v << " with corner "
2192                                         << corner.first << ":" << corner.second;
2193 #endif
2194           int wl = HGCalWaferIndex::waferIndex(lay, u, v);
2195           if (php.waferMaskMode_ > 0) {
2196             bool v17OrLess = (php.mode_ < HGCalGeometryMode::Hexagon8CalibCell);
2197             std::pair<int, int> corner0 = HGCalWaferMask::getTypeMode(
2198                 xpos0, ypos0, r1, R1, php.rMinLayHex_[i], php.rMaxLayHex_[i], type, php.waferMaskMode_, v17OrLess);
2199             if ((php.mode_ == HGCalGeometryMode::Hexagon8File) || (php.mode_ == HGCalGeometryMode::Hexagon8Module) ||
2200                 (php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
2201                 (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell)) {
2202               auto itr = php.waferInfoMap_.find(wl);
2203               if (itr != php.waferInfoMap_.end()) {
2204                 int part = (itr->second).part;
2205                 int orient = (itr->second).orient;
2206                 bool ok = ((php.mode_ == HGCalGeometryMode::Hexagon8Cassette) ||
2207                            (php.mode_ == HGCalGeometryMode::Hexagon8CalibCell))
2208                               ? true
2209                               : HGCalWaferMask::goodTypeMode(
2210                                     xpos0, ypos0, r1, R1, php.rMinLayHex_[i], php.rMaxLayHex_[i], part, orient, false);
2211                 if (ok)
2212                   corner0 = std::make_pair(part, (HGCalTypes::k_OffsetRotation + orient));
2213 #ifdef EDM_ML_DEBUG
2214                 edm::LogVerbatim("HGCalGeom")
2215                     << "Layer:u:v " << i << ":" << lay << ":" << u << ":" << v << " Part " << corner0.first << ":"
2216                     << part << " Orient " << corner0.second << ":" << orient << " Position " << xpos0 << ":" << ypos0
2217                     << " delta " << r1 << ":" << R1 << " Limit " << php.rMinLayHex_[i] << ":" << php.rMaxLayHex_[i]
2218                     << " Compatibiliety Flag " << ok;
2219 #endif
2220                 if (!ok)
2221                   ++nwarn;
2222               }
2223             }
2224             waferTypes[wl] = corner0;
2225 #ifdef EDM_ML_DEBUG
2226             edm::LogVerbatim("HGCalGeom")
2227                 << "Layer " << lay << " u|v " << u << ":" << v << " Index " << std::hex << wl << std::dec << " pos "
2228                 << xpos0 << ":" << ypos0 << " R " << r1 << ":" << R1 << " Range " << php.rMinLayHex_[i] << ":"
2229                 << php.rMaxLayHex_[i] << type << ":" << php.waferMaskMode_ << " corner " << corner.first << ":"
2230                 << corner.second << " croner0 " << corner0.first << ":" << corner0.second;
2231 #endif
2232           } else {
2233             waferTypes[wl] = corner;
2234 #ifdef EDM_ML_DEBUG
2235             edm::LogVerbatim("HGCalGeom") << "Layer " << lay << " u|v " << u << ":" << v << " with corner "
2236                                           << corner.first << ":" << corner.second;
2237 #endif
2238           }
2239         }
2240       }
2241     }
2242   }
2243   if (nwarn > 0)
2244     edm::LogWarning("HGCalGeom") << "HGCalGeomParameters::loadWafer8: there are " << nwarn
2245                                  << " wafers with non-matching partial- orientation types";
2246   php.waferUVMax_ = uvmax;
2247   php.waferUVMaxLayer_ = uvmx;
2248   php.wafersInLayers_ = wafersInLayers;
2249   php.typesInLayers_ = typesInLayers;
2250   php.waferTypes_ = waferTypes;
2251   php.nSectors_ = static_cast<int>(php.waferCopy_.size());
2252   HGCalParameters::hgtrap mytr;
2253   mytr.lay = 1;
2254   mytr.bl = php.waferR_;
2255   mytr.tl = php.waferR_;
2256   mytr.h = php.waferR_;
2257   mytr.alpha = 0.0;
2258   mytr.cellSize = HGCalParameters::k_ScaleToDDD * php.waferSize_;
2259   for (auto const& dz : php.cellThickness_) {
2260     mytr.dz = 0.5 * HGCalParameters::k_ScaleToDDD * dz;
2261     php.fillModule(mytr, false);
2262   }
2263   for (unsigned k = 0; k < php.cellThickness_.size(); ++k) {
2264     HGCalParameters::hgtrap mytr = php.getModule(k, false);
2265     mytr.bl *= HGCalParameters::k_ScaleFromDDD;
2266     mytr.tl *= HGCalParameters::k_ScaleFromDDD;
2267     mytr.h *= HGCalParameters::k_ScaleFromDDD;
2268     mytr.dz *= HGCalParameters::k_ScaleFromDDD;
2269     mytr.cellSize *= HGCalParameters::k_ScaleFromDDD;
2270     php.fillModule(mytr, true);
2271   }
2272 #ifdef EDM_ML_DEBUG
2273   edm::LogVerbatim("HGCalGeom") << "HGCalGeomParameters: Total of " << php.waferCopy_.size() << " wafers";
2274   for (unsigned int k = 0; k < php.waferCopy_.size(); ++k) {
2275     int id = php.waferCopy_[k];
2276     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << std::hex << id << std::dec << ":"
2277                                   << HGCalWaferIndex::waferLayer(id) << ":" << HGCalWaferIndex::waferU(id) << ":"
2278                                   << HGCalWaferIndex::waferV(id) << " x " << php.waferPosX_[k] << " y "
2279                                   << php.waferPosY_[k] << " index " << php.wafersInLayers_[id];
2280   }
2281   edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Total of " << php.waferTypeL_.size() << " wafer types";
2282   for (unsigned int k = 0; k < php.waferTypeL_.size(); ++k) {
2283     int id = indtypes[k];
2284     edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << php.typesInLayers_[id] << ":" << php.waferTypeL_[k] << " ID "
2285                                   << std::hex << id << std::dec << ":" << HGCalWaferIndex::waferLayer(id) << ":"
2286                                   << HGCalWaferIndex::waferU(id) << ":" << HGCalWaferIndex::waferV(id);
2287   }
2288 #endif
2289 
2290   //Wafer offset
2291   php.xLayerHex_.clear();
2292   php.yLayerHex_.clear();
2293   double waferSize = php.waferSize_ + php.sensorSeparation_;
2294 #ifdef EDM_ML_DEBUG
2295   edm::LogVerbatim("HGCalGeom") << "WaferSize " << waferSize;
2296 #endif
2297   for (unsigned int k = 0; k < php.zLayerHex_.size(); ++k) {
2298     int copy = k + php.layerOffset_;
2299     std::pair<double, double> xyoff = geomTools_.shiftXY(php.layerCenter_[copy], waferSize);
2300     php.xLayerHex_.emplace_back(xyoff.first);
2301     php.yLayerHex_.emplace_back(xyoff.second);
2302 #ifdef EDM_ML_DEBUG
2303     edm::LogVerbatim("HGCalGeom") << "Layer[" << k << "] Off " << copy << ":" << php.layerCenter_[copy] << " Shift "
2304                                   << xyoff.first << ":" << xyoff.second;
2305 #endif
2306   }
2307 }
2308 
2309 void HGCalGeomParameters::loadCellParsHexagon(const DDCompactView* cpv, HGCalParameters& php) {
2310   // Special parameters for cell parameters
2311   std::string attribute = "OnlyForHGCalNumbering";
2312   DDSpecificsHasNamedValueFilter filter1{attribute};
2313   DDFilteredView fv1(*cpv, filter1);
2314   bool ok = fv1.firstChild();
2315 
2316   if (ok) {
2317     php.cellFine_ = dbl_to_int(cpv->vector("waferFine"));
2318     php.cellCoarse_ = dbl_to_int(cpv->vector("waferCoarse"));
2319   }
2320 
2321   loadCellParsHexagon(php);
2322 }
2323 
2324 void HGCalGeomParameters::loadCellParsHexagon(const cms::DDVectorsMap& vmap, HGCalParameters& php) {
2325   for (auto const& it : vmap) {
2326     if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "waferFine")) {
2327       for (const auto& i : it.second)
2328         php.cellFine_.emplace_back(std::round(i));
2329     } else if (dd4hep::dd::compareEqual(dd4hep::dd::noNamespace(it.first), "waferCoarse")) {
2330       for (const auto& i : it.second)
2331         php.cellCoarse_.emplace_back(std::round(i));
2332     }
2333   }
2334 
2335   loadCellParsHexagon(php);
2336 }
2337 
2338 void HGCalGeomParameters::loadCellParsHexagon(const HGCalParameters& php) {
2339 #ifdef EDM_ML_DEBUG
2340   edm::LogVerbatim("HGCalGeom") << "HGCalLoadCellPars: " << php.cellFine_.size() << " rows for fine cells";
2341   for (unsigned int k = 0; k < php.cellFine_.size(); ++k)
2342     edm::LogVerbatim("HGCalGeom") << "[" << k << "]: " << php.cellFine_[k];
2343   edm::LogVerbatim("HGCalGeom") << "HGCalLoadCellPars: " << php.cellCoarse_.size() << " rows for coarse cells";
2344   for (unsigned int k = 0; k < php.cellCoarse_.size(); ++k)
2345     edm::LogVerbatim("HGCalGeom") << "[" << k << "]: " << php.cellCoarse_[k];
2346 #endif
2347 }
2348 
2349 void HGCalGeomParameters::loadCellTrapezoid(HGCalParameters& php) {
2350   php.xLayerHex_.resize(php.zLayerHex_.size(), 0);
2351   php.yLayerHex_.resize(php.zLayerHex_.size(), 0);
2352 #ifdef EDM_ML_DEBUG
2353   edm::LogVerbatim("HGCalGeom") << "HGCalParameters: x|y|zLayerHex in array of size " << php.zLayerHex_.size();
2354   for (unsigned int k = 0; k < php.zLayerHex_.size(); ++k)
2355     edm::LogVerbatim("HGCalGeom") << "Layer[" << k << "] Shift " << php.xLayerHex_[k] << ":" << php.yLayerHex_[k] << ":"
2356                                   << php.zLayerHex_[k];
2357 #endif
2358   // Find the radius of each eta-partitions
2359   if ((php.mode_ == HGCalGeometryMode::TrapezoidFile) || (php.mode_ == HGCalGeometryMode::TrapezoidModule) ||
2360       (php.mode_ == HGCalGeometryMode::TrapezoidCassette) || (php.mode_ == HGCalGeometryMode::TrapezoidFineCell)) {
2361     //Ring radii for each partition
2362 #ifdef EDM_ML_DEBUG
2363     edm::LogVerbatim("HGCalGeom") << "HGCalParameters: Mode " << php.mode_ << ":"
2364                                   << HGCalGeometryMode::TrapezoidFineCell << " Sizes " << php.tileRingFineR_.size()
2365                                   << ":" << php.tileRingR_.size();
2366 #endif
2367     for (unsigned int k = 0; k < 2; ++k) {
2368       bool fine = ((k == 0) && (php.mode_ == HGCalGeometryMode::TrapezoidFineCell));
2369       unsigned int sizeR = (fine) ? php.tileRingFineR_.size() : php.tileRingR_.size();
2370       for (unsigned int kk = 0; kk < sizeR; ++kk) {
2371         if (fine)
2372           php.radiusLayer_[k].emplace_back(php.tileRingFineR_[kk].first);
2373         else
2374           php.radiusLayer_[k].emplace_back(php.tileRingR_[kk].first);
2375 #ifdef EDM_ML_DEBUG
2376         double zv = ((k == 0) ? (php.zLayerHex_[php.layerFrontBH_[1] - php.firstLayer_])
2377                               : (php.zLayerHex_[php.zLayerHex_.size() - 1]));
2378         double rv = php.radiusLayer_[k].back();
2379         double eta = -(std::log(std::tan(0.5 * std::atan(rv / zv))));
2380         edm::LogVerbatim("HGCalGeom") << "New [" << kk << "] new R = " << rv << " Eta = " << eta;
2381 #endif
2382       }
2383       if (fine)
2384         php.radiusLayer_[k].emplace_back(php.tileRingFineR_[php.tileRingFineR_.size() - 1].second);
2385       else
2386         php.radiusLayer_[k].emplace_back(php.tileRingR_[php.tileRingR_.size() - 1].second);
2387     }
2388     // Minimum and maximum radius index for each layer
2389 #ifdef EDM_ML_DEBUG
2390     edm::LogVerbatim("HGCalGeom") << "Till Ring Range size " << php.tileRingFineRange_.size() << ":"
2391                                   << php.tileRingRange_.size() << ":" << php.nPhiBinBH_.size() << ":"
2392                                   << php.zLayerHex_.size() << ":" << php.nPhiLayer_.size();
2393 #endif
2394     unsigned int k1(0), k2(0);
2395     for (unsigned int k = 0; k < php.zLayerHex_.size(); ++k) {
2396       if (php.nPhiLayer_[k] > 288) {
2397         php.iradMinBH_.emplace_back(1 + php.tileRingFineRange_[k1].first);
2398         php.iradMaxBH_.emplace_back(1 + php.tileRingFineRange_[k1].second);
2399         ++k1;
2400       } else {
2401         php.iradMinBH_.emplace_back(1 + php.tileRingRange_[k2].first);
2402         php.iradMaxBH_.emplace_back(1 + php.tileRingRange_[k2].second);
2403         ++k2;
2404       }
2405 #ifdef EDM_ML_DEBUG
2406       int kk = php.scintType(php.firstLayer_ + static_cast<int>(k));
2407       edm::LogVerbatim("HGCalGeom") << "New Layer " << k << " Type " << kk << " Low edge " << php.iradMinBH_.back()
2408                                     << " Top edge " << php.iradMaxBH_.back();
2409 #endif
2410     }
2411   } else {
2412     //Ring radii for each partition
2413     for (unsigned int k = 0; k < 2; ++k) {
2414       double rmax = ((k == 0) ? (php.rMaxLayHex_[php.layerFrontBH_[1] - php.firstLayer_] - 1)
2415                               : (php.rMaxLayHex_[php.rMaxLayHex_.size() - 1]));
2416       double rv = php.rMinLayerBH_[k];
2417       double zv = ((k == 0) ? (php.zLayerHex_[php.layerFrontBH_[1] - php.firstLayer_])
2418                             : (php.zLayerHex_[php.zLayerHex_.size() - 1]));
2419       php.radiusLayer_[k].emplace_back(rv);
2420 #ifdef EDM_ML_DEBUG
2421       double eta = -(std::log(std::tan(0.5 * std::atan(rv / zv))));
2422       edm::LogVerbatim("HGCalGeom") << "Old [" << k << "] rmax " << rmax << " Z = " << zv
2423                                     << " dEta = " << php.cellSize_[k] << "\nOld[0] new R = " << rv << " Eta = " << eta;
2424       int kount(1);
2425 #endif
2426       while (rv < rmax) {
2427         double eta = -(php.cellSize_[k] + std::log(std::tan(0.5 * std::atan(rv / zv))));
2428         rv = zv * std::tan(2.0 * std::atan(std::exp(-eta)));
2429         php.radiusLayer_[k].emplace_back(rv);
2430 #ifdef EDM_ML_DEBUG
2431         edm::LogVerbatim("HGCalGeom") << "Old [" << kount << "] new R = " << rv << " Eta = " << eta;
2432         ++kount;
2433 #endif
2434       }
2435     }
2436     // Find minimum and maximum radius index for each layer
2437     for (unsigned int k = 0; k < php.zLayerHex_.size(); ++k) {
2438       int kk = php.scintType(php.firstLayer_ + static_cast<int>(k));
2439       std::vector<double>::iterator low, high;
2440       low = std::lower_bound(php.radiusLayer_[kk].begin(), php.radiusLayer_[kk].end(), php.rMinLayHex_[k]);
2441 #ifdef EDM_ML_DEBUG
2442       edm::LogVerbatim("HGCalGeom") << "Old [" << k << "] RLow = " << php.rMinLayHex_[k] << " pos "
2443                                     << static_cast<int>(low - php.radiusLayer_[kk].begin());
2444 #endif
2445       if (low == php.radiusLayer_[kk].begin())
2446         ++low;
2447       int irlow = static_cast<int>(low - php.radiusLayer_[kk].begin());
2448       double drlow = php.radiusLayer_[kk][irlow] - php.rMinLayHex_[k];
2449 #ifdef EDM_ML_DEBUG
2450       edm::LogVerbatim("HGCalGeom") << "irlow " << irlow << " dr " << drlow << " min " << php.minTileSize_;
2451 #endif
2452       if (drlow < php.minTileSize_) {
2453         ++irlow;
2454 #ifdef EDM_ML_DEBUG
2455         drlow = php.radiusLayer_[kk][irlow] - php.rMinLayHex_[k];
2456         edm::LogVerbatim("HGCalGeom") << "Modified irlow " << irlow << " dr " << drlow;
2457 #endif
2458       }
2459       high = std::lower_bound(php.radiusLayer_[kk].begin(), php.radiusLayer_[kk].end(), php.rMaxLayHex_[k]);
2460 #ifdef EDM_ML_DEBUG
2461       edm::LogVerbatim("HGCalGeom") << "Old [" << k << "] RHigh = " << php.rMaxLayHex_[k] << " pos "
2462                                     << static_cast<int>(high - php.radiusLayer_[kk].begin());
2463 #endif
2464       if (high == php.radiusLayer_[kk].end())
2465         --high;
2466       int irhigh = static_cast<int>(high - php.radiusLayer_[kk].begin());
2467       double drhigh = php.rMaxLayHex_[k] - php.radiusLayer_[kk][irhigh - 1];
2468 #ifdef EDM_ML_DEBUG
2469       edm::LogVerbatim("HGCalGeom") << "irhigh " << irhigh << " dr " << drhigh << " min " << php.minTileSize_;
2470 #endif
2471       if (drhigh < php.minTileSize_) {
2472         --irhigh;
2473 #ifdef EDM_ML_DEBUG
2474         drhigh = php.rMaxLayHex_[k] - php.radiusLayer_[kk][irhigh - 1];
2475         edm::LogVerbatim("HGCalGeom") << "Modified irhigh " << irhigh << " dr " << drhigh;
2476 #endif
2477       }
2478       php.iradMinBH_.emplace_back(irlow);
2479       php.iradMaxBH_.emplace_back(irhigh);
2480 #ifdef EDM_ML_DEBUG
2481       edm::LogVerbatim("HGCalGeom") << "Old Layer " << k << " Type " << kk << " Low edge " << irlow << ":" << drlow
2482                                     << " Top edge " << irhigh << ":" << drhigh;
2483 #endif
2484     }
2485   }
2486 #ifdef EDM_ML_DEBUG
2487   for (unsigned int k = 0; k < 2; ++k) {
2488     edm::LogVerbatim("HGCalGeom") << "Type " << k << " with " << php.radiusLayer_[k].size() << " radii";
2489     for (unsigned int kk = 0; kk < php.radiusLayer_[k].size(); ++kk)
2490       edm::LogVerbatim("HGCalGeom") << "Ring[" << kk << "] " << php.radiusLayer_[k][kk];
2491   }
2492 #endif
2493 
2494   // Now define the volumes
2495   int im(0);
2496   php.waferUVMax_ = 0;
2497   HGCalParameters::hgtrap mytr;
2498   mytr.alpha = 0.0;
2499   for (unsigned int k = 0; k < php.zLayerHex_.size(); ++k) {
2500     if (php.iradMaxBH_[k] > php.waferUVMax_)
2501       php.waferUVMax_ = php.iradMaxBH_[k];
2502     int kk = ((php.firstLayer_ + static_cast<int>(k)) < php.layerFrontBH_[1]) ? 0 : 1;
2503     int irm = php.radiusLayer_[kk].size() - 1;
2504 #ifdef EDM_ML_DEBUG
2505     double rmin = php.radiusLayer_[kk][std::max((php.iradMinBH_[k] - 1), 0)];
2506     double rmax = php.radiusLayer_[kk][std::min(php.iradMaxBH_[k], irm)];
2507     edm::LogVerbatim("HGCalGeom") << "Layer " << php.firstLayer_ + k << ":" << kk << " Radius range "
2508                                   << php.iradMinBH_[k] << ":" << php.iradMaxBH_[k] << ":" << rmin << ":" << rmax;
2509 #endif
2510     mytr.lay = php.firstLayer_ + k;
2511     for (int irad = php.iradMinBH_[k]; irad <= php.iradMaxBH_[k]; ++irad) {
2512       double rmin = php.radiusLayer_[kk][std::max((irad - 1), 0)];
2513       double rmax = php.radiusLayer_[kk][std::min(irad, irm)];
2514       mytr.bl = 0.5 * rmin * php.scintCellSize(mytr.lay);
2515       mytr.tl = 0.5 * rmax * php.scintCellSize(mytr.lay);
2516       mytr.h = 0.5 * (rmax - rmin);
2517       mytr.dz = 0.5 * php.waferThick_;
2518       mytr.cellSize = 0.5 * (rmax + rmin) * php.scintCellSize(mytr.lay);
2519       php.fillModule(mytr, true);
2520       mytr.bl *= HGCalParameters::k_ScaleToDDD;
2521       mytr.tl *= HGCalParameters::k_ScaleToDDD;
2522       mytr.h *= HGCalParameters::k_ScaleToDDD;
2523       mytr.dz *= HGCalParameters::k_ScaleToDDD;
2524       mytr.cellSize *= HGCalParameters::k_ScaleFromDDD;
2525       php.fillModule(mytr, false);
2526       if (irad == php.iradMinBH_[k])
2527         php.firstModule_.emplace_back(im);
2528       ++im;
2529       if (irad == php.iradMaxBH_[k] - 1)
2530         php.lastModule_.emplace_back(im);
2531     }
2532   }
2533   php.nSectors_ = php.waferUVMax_;
2534 #ifdef EDM_ML_DEBUG
2535   edm::LogVerbatim("HGCalGeom") << "Maximum radius index " << php.waferUVMax_;
2536   for (unsigned int k = 0; k < php.firstModule_.size(); ++k)
2537     edm::LogVerbatim("HGCalGeom") << "Layer " << k + php.firstLayer_ << " Modules " << php.firstModule_[k] << ":"
2538                                   << php.lastModule_[k];
2539 #endif
2540 }
2541 
2542 std::vector<double> HGCalGeomParameters::getDDDArray(const std::string& str, const DDsvalues_type& sv, const int nmin) {
2543   DDValue value(str);
2544   if (DDfetch(&sv, value)) {
2545     const std::vector<double>& fvec = value.doubles();
2546     int nval = fvec.size();
2547     if (nmin > 0) {
2548       if (nval < nmin) {
2549         throw cms::Exception("DDException")
2550             << "HGCalGeomParameters:  # of " << str << " bins " << nval << " < " << nmin << " ==> illegal";
2551       }
2552     } else {
2553       if (nval < 1 && nmin == 0) {
2554         throw cms::Exception("DDException")
2555             << "HGCalGeomParameters: # of " << str << " bins " << nval << " < 1 ==> illegal"
2556             << " (nmin=" << nmin << ")";
2557       }
2558     }
2559     return fvec;
2560   } else {
2561     if (nmin >= 0) {
2562       throw cms::Exception("DDException") << "HGCalGeomParameters: cannot get array " << str;
2563     }
2564     std::vector<double> fvec;
2565     return fvec;
2566   }
2567 }
2568 
2569 std::pair<double, double> HGCalGeomParameters::cellPosition(
2570     const std::vector<HGCalGeomParameters::cellParameters>& wafers,
2571     std::vector<HGCalGeomParameters::cellParameters>::const_iterator& itrf,
2572     int wafer,
2573     double xx,
2574     double yy) {
2575   if (itrf == wafers.end()) {
2576     for (std::vector<HGCalGeomParameters::cellParameters>::const_iterator itr = wafers.begin(); itr != wafers.end();
2577          ++itr) {
2578       if (itr->wafer == wafer) {
2579         itrf = itr;
2580         break;
2581       }
2582     }
2583   }
2584   double dx(0), dy(0);
2585   if (itrf != wafers.end()) {
2586     dx = (xx - itrf->xyz.x());
2587     if (std::abs(dx) < tolerance)
2588       dx = 0;
2589     dy = (yy - itrf->xyz.y());
2590     if (std::abs(dy) < tolerance)
2591       dy = 0;
2592   }
2593   return std::make_pair(dx, dy);
2594 }
2595 
2596 void HGCalGeomParameters::rescale(std::vector<double>& v, const double s) {
2597   std::for_each(v.begin(), v.end(), [s](double& n) { n *= s; });
2598 }
2599 
2600 void HGCalGeomParameters::resetZero(std::vector<double>& v) {
2601   for (auto& n : v) {
2602     if (std::abs(n) < tolmin)
2603       n = 0;
2604   }
2605 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!