Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-11-06 04:13:40

0001 #include "Geometry/HcalCommonData/interface/HcalGeomParameters.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 "DetectorDescription/Core/interface/DDutils.h"
0008 #include "DetectorDescription/Core/interface/DDValue.h"
0009 #include "DetectorDescription/Core/interface/DDFilter.h"
0010 #include "DetectorDescription/Core/interface/DDSolid.h"
0011 #include "DetectorDescription/Core/interface/DDConstant.h"
0012 #include "DetectorDescription/RegressionTest/interface/DDErrorDetection.h"
0013 #include "CondFormats/GeometryObjects/interface/HcalParameters.h"
0014 
0015 //#define EDM_ML_DEBUG
0016 using namespace geant_units::operators;
0017 
0018 static const double tan10deg = std::tan(10._deg);
0019 
0020 void HcalGeomParameters::getConstRHO(std::vector<double>& rHO) const {
0021   rHO.emplace_back(rminHO_);
0022   for (double i : etaHO_)
0023     rHO.emplace_back(i);
0024 }
0025 
0026 std::vector<int> HcalGeomParameters::getModHalfHBHE(const int type) const {
0027   std::vector<int> modHalf;
0028   if (type == 0) {
0029     modHalf.emplace_back(nmodHB_);
0030     modHalf.emplace_back(nzHB_);
0031   } else {
0032     modHalf.emplace_back(nmodHE_);
0033     modHalf.emplace_back(nzHE_);
0034   }
0035   return modHalf;
0036 }
0037 
0038 void HcalGeomParameters::loadGeometry(const DDFilteredView& _fv, HcalParameters& php) {
0039   DDFilteredView fv = _fv;
0040   bool dodet = true;
0041   bool hf = false;
0042   clear(php);
0043 
0044   while (dodet) {
0045     DDTranslation t = fv.translation();
0046     std::vector<int> copy = fv.copyNumbers();
0047     const DDSolid& sol = fv.logicalPart().solid();
0048     int idet = 0, lay = -1;
0049     int nsiz = static_cast<int>(copy.size());
0050     if (nsiz > 0)
0051       lay = copy[nsiz - 1] / 10;
0052     if (nsiz > 1)
0053       idet = copy[nsiz - 2] / 1000;
0054 #ifdef EDM_ML_DEBUG
0055     edm::LogVerbatim("HCalGeom") << "Name " << fv.logicalPart().solid().name() << " Copy " << copy.size();
0056 #endif
0057     double dx = 0, dy = 0, dz = 0, dx1 = 0, dx2 = 0;
0058     double alp(0);
0059     if (sol.shape() == DDSolidShape::ddbox) {
0060       const DDBox& box = static_cast<DDBox>(sol);
0061       dx = HcalGeomParameters::k_ScaleFromDDDToG4 * box.halfX();
0062       dy = HcalGeomParameters::k_ScaleFromDDDToG4 * box.halfY();
0063       dz = HcalGeomParameters::k_ScaleFromDDDToG4 * box.halfZ();
0064     } else if (sol.shape() == DDSolidShape::ddtrap) {
0065       const DDTrap& trp = static_cast<DDTrap>(sol);
0066       dx1 = HcalGeomParameters::k_ScaleFromDDDToG4 * trp.x1();
0067       dx2 = HcalGeomParameters::k_ScaleFromDDDToG4 * trp.x2();
0068       dx = 0.25 * HcalGeomParameters::k_ScaleFromDDDToG4 * (trp.x1() + trp.x2() + trp.x3() + trp.x4());
0069       dy = 0.5 * HcalGeomParameters::k_ScaleFromDDDToG4 * (trp.y1() + trp.y2());
0070       dz = HcalGeomParameters::k_ScaleFromDDDToG4 * trp.halfZ();
0071       alp = 0.5 * (trp.alpha1() + trp.alpha2());
0072     } else if (sol.shape() == DDSolidShape::ddtubs) {
0073       const DDTubs& tub = static_cast<DDTubs>(sol);
0074       dx = HcalGeomParameters::k_ScaleFromDDDToG4 * tub.rIn();
0075       dy = HcalGeomParameters::k_ScaleFromDDDToG4 * tub.rOut();
0076       dz = HcalGeomParameters::k_ScaleFromDDDToG4 * tub.zhalf();
0077     }
0078     if (idet == 3) {
0079       // HB
0080 #ifdef EDM_ML_DEBUG
0081       edm::LogVerbatim("HCalGeom") << "HB " << sol.name() << " Shape " << sol.shape() << " Layer " << lay << " R "
0082                                    << t.Rho();
0083 #endif
0084       if (lay >= 0 && lay < maxLayer_) {
0085         ib_[lay]++;
0086         rb_[lay] += (HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho());
0087         if (thkb_[lay] <= 0) {
0088           if (lay < 17)
0089             thkb_[lay] = dx;
0090           else
0091             thkb_[lay] = std::min(dx, dy);
0092         }
0093         if (lay < 17) {
0094           bool found = false;
0095           for (double k : rxb_) {
0096             if (std::abs(k - (HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho())) < 0.01) {
0097               found = true;
0098               break;
0099             }
0100           }
0101           if (!found) {
0102             rxb_.emplace_back(HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho());
0103             php.rhoxHB.emplace_back(HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho() * std::cos(t.phi()));
0104             php.zxHB.emplace_back(HcalGeomParameters::k_ScaleFromDDDToG4 * std::abs(t.z()));
0105             php.dyHB.emplace_back(2. * dy);
0106             php.dxHB.emplace_back(2. * dz);
0107             php.layHB.emplace_back(lay);
0108           }
0109         }
0110       }
0111       if (lay == 2) {
0112         int iz = copy[nsiz - 5];
0113         int fi = copy[nsiz - 4];
0114         unsigned int it1 = find(iz, izb_);
0115         if (it1 == izb_.size())
0116           izb_.emplace_back(iz);
0117         unsigned int it2 = find(fi, phib_);
0118         if (it2 == phib_.size())
0119           phib_.emplace_back(fi);
0120       }
0121       if (lay == 18) {
0122         int ifi = -1, ich = -1;
0123         if (nsiz > 2)
0124           ifi = copy[nsiz - 3];
0125         if (nsiz > 3)
0126           ich = copy[nsiz - 4];
0127         double z1 = std::abs((HcalGeomParameters::k_ScaleFromDDDToG4 * t.z()) + dz);
0128         double z2 = std::abs((HcalGeomParameters::k_ScaleFromDDDToG4 * t.z()) - dz);
0129         if (std::abs(z1 - z2) < 0.01)
0130           z1 = 0;
0131         if (ifi == 1 && ich == 4) {
0132           if (z1 > z2) {
0133             double tmp = z1;
0134             z1 = z2;
0135             z2 = tmp;
0136           }
0137           bool sok = true;
0138           for (unsigned int kk = 0; kk < php.zHO.size(); kk++) {
0139             if (std::abs(z2 - php.zHO[kk]) < 0.01) {
0140               sok = false;
0141               break;
0142             } else if (z2 < php.zHO[kk]) {
0143               php.zHO.resize(php.zHO.size() + 2);
0144               for (unsigned int kz = php.zHO.size() - 1; kz > kk + 1; kz = kz - 2) {
0145                 php.zHO[kz] = php.zHO[kz - 2];
0146                 php.zHO[kz - 1] = php.zHO[kz - 3];
0147               }
0148               php.zHO[kk + 1] = z2;
0149               php.zHO[kk] = z1;
0150               sok = false;
0151               break;
0152             }
0153           }
0154           if (sok) {
0155             php.zHO.emplace_back(z1);
0156             php.zHO.emplace_back(z2);
0157           }
0158 #ifdef EDM_ML_DEBUG
0159           edm::LogVerbatim("HCalGeom") << "Detector " << idet << " Lay " << lay << " fi " << ifi << " " << ich << " z "
0160                                        << z1 << " " << z2;
0161 #endif
0162         }
0163       }
0164     } else if (idet == 4) {
0165       // HE
0166 #ifdef EDM_ML_DEBUG
0167       edm::LogVerbatim("HCalGeom") << "HE " << sol.name() << " Shape " << sol.shape() << " Layer " << lay << " Z "
0168                                    << t.z();
0169 #endif
0170       if (lay >= 0 && lay < maxLayer_) {
0171         ie_[lay]++;
0172         ze_[lay] += std::abs(HcalGeomParameters::k_ScaleFromDDDToG4 * t.z());
0173         if (thke_[lay] <= 0)
0174           thke_[lay] = dz;
0175         double rinHE = HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho() * cos(alp) - dy;
0176         double routHE = HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho() * cos(alp) + dy;
0177         rminHE_[lay] += rinHE;
0178         rmaxHE_[lay] += routHE;
0179         bool found = false;
0180         for (double k : php.zxHE) {
0181           if (std::abs(k - std::abs(HcalGeomParameters::k_ScaleFromDDDToG4 * t.z())) < 0.01) {
0182             found = true;
0183             break;
0184           }
0185         }
0186         if (!found) {
0187           php.zxHE.emplace_back(HcalGeomParameters::k_ScaleFromDDDToG4 * std::abs(t.z()));
0188           php.rhoxHE.emplace_back(HcalGeomParameters::k_ScaleFromDDDToG4 * t.Rho() * std::cos(t.phi()));
0189           php.dyHE.emplace_back(dy * std::cos(t.phi()));
0190           dx1 -= 0.5 * (HcalGeomParameters::k_ScaleFromDDDToG4 * t.rho() - dy) * std::cos(t.phi()) * tan10deg;
0191           dx2 -= 0.5 * (HcalGeomParameters::k_ScaleFromDDDToG4 * t.rho() + dy) * std::cos(t.phi()) * tan10deg;
0192           php.dx1HE.emplace_back(-dx1);
0193           php.dx2HE.emplace_back(-dx2);
0194           php.layHE.emplace_back(lay);
0195         }
0196       }
0197       if (copy[nsiz - 1] == kHELayer1_ || copy[nsiz - 1] == kHELayer2_) {
0198         int iz = copy[nsiz - 7];
0199         int fi = copy[nsiz - 5];
0200         unsigned int it1 = find(iz, ize_);
0201         if (it1 == ize_.size())
0202           ize_.emplace_back(iz);
0203         unsigned int it2 = find(fi, phie_);
0204         if (it2 == phie_.size())
0205           phie_.emplace_back(fi);
0206       }
0207     } else if (idet == 5) {
0208       // HF
0209       if (!hf) {
0210         const std::vector<double>& paras = sol.parameters();
0211 #ifdef EDM_ML_DEBUG
0212         edm::LogVerbatim("HCalGeom") << "HF " << sol.name() << " Shape " << sol.shape() << " Z " << t.z() << " with "
0213                                      << paras.size() << " Parameters";
0214         for (unsigned j = 0; j < paras.size(); j++)
0215           edm::LogVerbatim("HCalGeom") << "HF Parameter[" << j << "] = " << paras[j];
0216 #endif
0217         if (sol.shape() == DDSolidShape::ddpolycone_rrz) {
0218           int nz = (int)(paras.size()) - 3;
0219           dzVcal_ = 0.5 * HcalGeomParameters::k_ScaleFromDDDToG4 * (paras[nz] - paras[3]);
0220           hf = true;
0221         } else if (sol.shape() == DDSolidShape::ddtubs || sol.shape() == DDSolidShape::ddcons) {
0222           dzVcal_ = HcalGeomParameters::k_ScaleFromDDDToG4 * paras[0];
0223           hf = true;
0224         }
0225       }
0226 #ifdef EDM_ML_DEBUG
0227     } else {
0228       edm::LogVerbatim("HCalGeom") << "Unknown Detector " << idet << " for " << sol.name() << " Shape " << sol.shape()
0229                                    << " R " << t.Rho() << " Z " << t.z();
0230 #endif
0231     }
0232     dodet = fv.next();
0233   }
0234 
0235   loadfinal(php);
0236 }
0237 
0238 void HcalGeomParameters::loadGeometry(const cms::DDCompactView& cpv, HcalParameters& php) {
0239   const cms::DDFilter filter("OnlyForHcalSimNumbering", "HCAL");
0240   cms::DDFilteredView fv(cpv, filter);
0241   clear(php);
0242   bool hf(false);
0243   while (fv.firstChild()) {
0244     auto t = fv.translation();
0245     std::vector<double> paras = fv.parameters();
0246     std::vector<int> copy = fv.copyNos();
0247     int idet = 0, lay = -1;
0248     int nsiz = static_cast<int>(copy.size());
0249     if (nsiz > 0)
0250       lay = copy[0] / 10;
0251     if (nsiz > 1)
0252       idet = copy[1] / 1000;
0253 #ifdef EDM_ML_DEBUG
0254     edm::LogVerbatim("HCalGeom") << "Name " << fv.name() << " Copy " << copy.size();
0255     for (unsigned int n = 0; n < copy.size(); ++n)
0256       edm::LogVerbatim("HCalGeom") << "[" << n << "] " << copy[n];
0257     edm::LogVerbatim("HCalGeom") << "Detector " << idet << " Layer " << lay << " parameters: " << paras.size();
0258     for (unsigned int n = 0; n < paras.size(); ++n)
0259       edm::LogVerbatim("HCalGeom") << "[" << n << "] " << paras[n];
0260 #endif
0261     double dx = 0, dy = 0, dz = 0, dx1 = 0, dx2 = 0;
0262     double alp(0);
0263     if (dd4hep::isA<dd4hep::Box>(fv.solid())) {
0264       dx = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[0];
0265       dy = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[1];
0266       dz = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[2];
0267     } else if (dd4hep::isA<dd4hep::Trap>(fv.solid())) {
0268       dx1 = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[4];
0269       dx2 = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[5];
0270       dx = 0.25 * HcalGeomParameters::k_ScaleFromDD4hepToG4 * (paras[4] + paras[5] + paras[8] + paras[9]);
0271       dy = 0.5 * HcalGeomParameters::k_ScaleFromDD4hepToG4 * (paras[3] + paras[7]);
0272       dz = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[0];
0273       alp = 0.5 * (paras[6] + paras[10]);
0274     } else if (dd4hep::isA<dd4hep::Tube>(fv.solid())) {
0275       dx = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[0];
0276       dy = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[1];
0277       dz = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[2];
0278     }
0279     if (idet == 3) {
0280       // HB
0281 #ifdef EDM_ML_DEBUG
0282       edm::LogVerbatim("HCalGeom") << "HB " << fv.name() << " Shape " << cms::dd::name(cms::DDSolidShapeMap, fv.shape())
0283                                    << " Layer " << lay << " R " << t.Rho();
0284 #endif
0285       if (lay >= 0 && lay < maxLayer_) {
0286         ib_[lay]++;
0287         rb_[lay] += (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho());
0288         if (thkb_[lay] <= 0) {
0289           if (lay < 17)
0290             thkb_[lay] = dx;
0291           else
0292             thkb_[lay] = std::min(dx, dy);
0293         }
0294         if (lay < 17) {
0295           bool found = false;
0296           for (double k : rxb_) {
0297             if (std::abs(k - (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho())) < 0.01) {
0298               found = true;
0299               break;
0300             }
0301           }
0302           if (!found) {
0303             rxb_.emplace_back(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho());
0304             php.rhoxHB.emplace_back(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() * std::cos(t.phi()));
0305             php.zxHB.emplace_back(std::abs(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z()));
0306             php.dyHB.emplace_back(2. * dy);
0307             php.dxHB.emplace_back(2. * dz);
0308             php.layHB.emplace_back(lay);
0309           }
0310         }
0311       }
0312       if (lay == 2) {
0313         int iz = copy[4];
0314         int fi = copy[3];
0315         unsigned int it1 = find(iz, izb_);
0316         if (it1 == izb_.size())
0317           izb_.emplace_back(iz);
0318         unsigned int it2 = find(fi, phib_);
0319         if (it2 == phib_.size())
0320           phib_.emplace_back(fi);
0321       }
0322       if (lay == 18) {
0323         int ifi = -1, ich = -1;
0324         if (nsiz > 2)
0325           ifi = copy[2];
0326         if (nsiz > 3)
0327           ich = copy[3];
0328         double z1 = std::abs((HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z()) + dz);
0329         double z2 = std::abs((HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z()) - dz);
0330         if (std::abs(z1 - z2) < 0.01)
0331           z1 = 0;
0332         if (ifi == 1 && ich == 4) {
0333           if (z1 > z2) {
0334             double tmp = z1;
0335             z1 = z2;
0336             z2 = tmp;
0337           }
0338           bool sok = true;
0339           for (unsigned int kk = 0; kk < php.zHO.size(); kk++) {
0340             if (std::abs(z2 - php.zHO[kk]) < 0.01) {
0341               sok = false;
0342               break;
0343             } else if (z2 < php.zHO[kk]) {
0344               php.zHO.resize(php.zHO.size() + 2);
0345               for (unsigned int kz = php.zHO.size() - 1; kz > kk + 1; kz = kz - 2) {
0346                 php.zHO[kz] = php.zHO[kz - 2];
0347                 php.zHO[kz - 1] = php.zHO[kz - 3];
0348               }
0349               php.zHO[kk + 1] = z2;
0350               php.zHO[kk] = z1;
0351               sok = false;
0352               break;
0353             }
0354           }
0355           if (sok) {
0356             php.zHO.emplace_back(z1);
0357             php.zHO.emplace_back(z2);
0358           }
0359 #ifdef EDM_ML_DEBUG
0360           edm::LogVerbatim("HCalGeom") << "Detector " << idet << " Lay " << lay << " fi " << ifi << " " << ich << " z "
0361                                        << z1 << " " << z2;
0362 #endif
0363         }
0364       }
0365     } else if (idet == 4) {
0366       // HE
0367 #ifdef EDM_ML_DEBUG
0368       edm::LogVerbatim("HCalGeom") << "HE " << fv.name() << " Shape " << cms::dd::name(cms::DDSolidShapeMap, fv.shape())
0369                                    << " Layer " << lay << " Z " << t.z();
0370 #endif
0371       if (lay >= 0 && lay < maxLayer_) {
0372         ie_[lay]++;
0373         ze_[lay] += std::abs(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z());
0374         if (thke_[lay] <= 0)
0375           thke_[lay] = dz;
0376         double rinHE = HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() * cos(alp) - dy;
0377         double routHE = HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() * cos(alp) + dy;
0378         rminHE_[lay] += rinHE;
0379         rmaxHE_[lay] += routHE;
0380         bool found = false;
0381         for (double k : php.zxHE) {
0382           if (std::abs(k - std::abs(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z())) < 0.01) {
0383             found = true;
0384             break;
0385           }
0386         }
0387         if (!found) {
0388           php.zxHE.emplace_back(std::abs(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z()));
0389           php.rhoxHE.emplace_back(HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() * std::cos(t.phi()));
0390           php.dyHE.emplace_back(dy * std::cos(t.phi()));
0391           dx1 -= 0.5 * (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() - dy) * std::cos(t.phi()) * tan10deg;
0392           dx2 -= 0.5 * (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho() + dy) * std::cos(t.phi()) * tan10deg;
0393           php.dx1HE.emplace_back(-dx1);
0394           php.dx2HE.emplace_back(-dx2);
0395           php.layHE.emplace_back(lay);
0396         }
0397       }
0398       if (copy[0] == kHELayer1_ || copy[0] == kHELayer2_) {
0399         int iz = copy[6];
0400         int fi = copy[4];
0401         unsigned int it1 = find(iz, ize_);
0402         if (it1 == ize_.size())
0403           ize_.emplace_back(iz);
0404         unsigned int it2 = find(fi, phie_);
0405         if (it2 == phie_.size())
0406           phie_.emplace_back(fi);
0407       }
0408     } else if (idet == 5) {
0409       // HF
0410       if (!hf) {
0411 #ifdef EDM_ML_DEBUG
0412         edm::LogVerbatim("HCalGeom") << "HF " << fv.name() << " Shape "
0413                                      << cms::dd::name(cms::DDSolidShapeMap, fv.shape()) << " Z "
0414                                      << HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z() << " with " << paras.size()
0415                                      << " Parameters";
0416         for (unsigned j = 0; j < paras.size(); j++)
0417           edm::LogVerbatim("HCalGeom") << "HF Parameter[" << j << "] = " << paras[j];
0418 #endif
0419         if (dd4hep::isA<dd4hep::Polycone>(fv.solid())) {
0420           int nz = (int)(paras.size()) - 3;
0421           dzVcal_ = 0.5 * HcalGeomParameters::k_ScaleFromDD4hepToG4 * (paras[nz] - paras[3]);
0422           hf = true;
0423         } else if (dd4hep::isA<dd4hep::Tube>(fv.solid()) || dd4hep::isA<dd4hep::ConeSegment>(fv.solid())) {
0424           dzVcal_ = HcalGeomParameters::k_ScaleFromDD4hepToG4 * paras[2];
0425           hf = true;
0426         }
0427       }
0428 #ifdef EDM_ML_DEBUG
0429     } else {
0430       edm::LogVerbatim("HCalGeom") << "Unknown Detector " << idet << " for " << fv.name() << " Shape "
0431                                    << cms::dd::name(cms::DDSolidShapeMap, fv.shape()) << " R "
0432                                    << (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.Rho()) << " Z "
0433                                    << (HcalGeomParameters::k_ScaleFromDD4hepToG4 * t.z());
0434 #endif
0435     }
0436   }
0437   loadfinal(php);
0438 }
0439 
0440 unsigned int HcalGeomParameters::find(int element, std::vector<int>& array) const {
0441   unsigned int id = array.size();
0442   for (unsigned int i = 0; i < array.size(); i++) {
0443     if (element == array[i]) {
0444       id = i;
0445       break;
0446     }
0447   }
0448   return id;
0449 }
0450 
0451 double HcalGeomParameters::getEta(double r, double z) const {
0452   double tmp = 0;
0453   if (z != 0)
0454     tmp = -log(tan(0.5 * atan(r / z)));
0455 #ifdef EDM_ML_DEBUG
0456   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters::getEta " << r << " " << z << " ==> " << tmp;
0457 #endif
0458   return tmp;
0459 }
0460 
0461 void HcalGeomParameters::clear(HcalParameters& php) {
0462   php.rhoxHB.clear();
0463   php.zxHB.clear();
0464   php.dyHB.clear();
0465   php.dxHB.clear();
0466   php.layHB.clear();
0467   php.layHE.clear();
0468   php.zxHE.clear();
0469   php.rhoxHE.clear();
0470   php.dyHE.clear();
0471   php.dx1HE.clear();
0472   php.dx2HE.clear();
0473 
0474   // Initialize all variables
0475   nzHB_ = nmodHB_ = 0;
0476   nzHE_ = nmodHE_ = 0;
0477   for (int i = 0; i < 4; ++i)
0478     etaHO_[i] = 0;
0479   zVcal_ = dzVcal_ = dlShort_ = 0;
0480   rminHO_ = dzVcal_ = -1.;
0481   for (int i = 0; i < maxLayer_; ++i) {
0482     rb_.emplace_back(0.0);
0483     ze_.emplace_back(0.0);
0484     thkb_.emplace_back(-1.0);
0485     thke_.emplace_back(-1.0);
0486     ib_.emplace_back(0);
0487     ie_.emplace_back(0);
0488     rminHE_.emplace_back(0.0);
0489     rmaxHE_.emplace_back(0.0);
0490   }
0491 }
0492 
0493 void HcalGeomParameters::loadfinal(HcalParameters& php) {
0494   int ibmx = 0, iemx = 0;
0495   for (int i = 0; i < maxLayer_; i++) {
0496     if (ib_[i] > 0) {
0497       rb_[i] /= static_cast<double>(ib_[i]);
0498       ibmx = i + 1;
0499     }
0500     if (ie_[i] > 0) {
0501       ze_[i] /= static_cast<double>(ie_[i]);
0502       rminHE_[i] /= static_cast<double>(ie_[i]);
0503       rmaxHE_[i] /= static_cast<double>(ie_[i]);
0504       iemx = i + 1;
0505     }
0506 #ifdef EDM_ML_DEBUG
0507     edm::LogVerbatim("HCalGeom") << "Index " << i << " Barrel " << ib_[i] << " " << rb_[i] << " Endcap " << ie_[i]
0508                                  << " " << ze_[i] << ":" << rminHE_[i] << ":" << rmaxHE_[i];
0509 #endif
0510   }
0511   for (int i = 4; i >= 0; i--) {
0512     if (ib_[i] == 0) {
0513       rb_[i] = rb_[i + 1];
0514       thkb_[i] = thkb_[i + 1];
0515     }
0516     if (ie_[i] == 0) {
0517       ze_[i] = ze_[i + 1];
0518       thke_[i] = thke_[i + 1];
0519     }
0520 #ifdef EDM_ML_DEBUG
0521     if (ib_[i] == 0 || ie_[i] == 0)
0522       edm::LogVerbatim("HCalGeom") << "Index " << i << " Barrel " << ib_[i] << " " << rb_[i] << " Endcap " << ie_[i]
0523                                    << " " << ze_[i];
0524 #endif
0525   }
0526 
0527 #ifdef EDM_ML_DEBUG
0528   for (unsigned int k = 0; k < php.layHB.size(); ++k)
0529     edm::LogVerbatim("HCalGeom") << "HB: " << php.layHB[k] << " R " << rxb_[k] << " " << php.rhoxHB[k] << " Z "
0530                                  << php.zxHB[k] << " DY " << php.dyHB[k] << " DZ " << php.dxHB[k];
0531   for (unsigned int k = 0; k < php.layHE.size(); ++k)
0532     edm::LogVerbatim("HCalGeom") << "HE: " << php.layHE[k] << " R " << php.rhoxHE[k] << " Z " << php.zxHE[k]
0533                                  << " X1|X2 " << php.dx1HE[k] << "|" << php.dx2HE[k] << " DY " << php.dyHE[k];
0534   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters: Maximum Layer for HB " << ibmx << " for HE " << iemx
0535                                << " extent " << dzVcal_;
0536 #endif
0537 
0538   if (ibmx > 0) {
0539     php.rHB.resize(ibmx);
0540     php.drHB.resize(ibmx);
0541     for (int i = 0; i < ibmx; i++) {
0542       php.rHB[i] = rb_[i];
0543       php.drHB[i] = thkb_[i];
0544 #ifdef EDM_ML_DEBUG
0545       edm::LogVerbatim("HCalGeom") << "HcalGeomParameters: php.rHB[" << i << "] = " << php.rHB[i] << " php.drHB[" << i
0546                                    << "] = " << php.drHB[i];
0547 #endif
0548     }
0549   }
0550   if (iemx > 0) {
0551     php.zHE.resize(iemx);
0552     php.dzHE.resize(iemx);
0553     for (int i = 0; i < iemx; i++) {
0554       php.zHE[i] = ze_[i];
0555       php.dzHE[i] = thke_[i];
0556 #ifdef EDM_ML_DEBUG
0557       edm::LogVerbatim("HCalGeom") << "HcalGeomParameters: php.zHE[" << i << "] = " << php.zHE[i] << " php.dzHE[" << i
0558                                    << "] = " << php.dzHE[i];
0559 #endif
0560     }
0561   }
0562 
0563   nzHB_ = static_cast<int>(izb_.size());
0564   nmodHB_ = static_cast<int>(phib_.size());
0565 #ifdef EDM_ML_DEBUG
0566   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters::loadGeometry: " << nzHB_ << " barrel half-sectors";
0567   for (int i = 0; i < nzHB_; i++)
0568     edm::LogVerbatim("HCalGeom") << "Section " << i << " Copy number " << izb_[i];
0569   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters::loadGeometry: " << nmodHB_ << " barrel modules";
0570   for (int i = 0; i < nmodHB_; i++)
0571     edm::LogVerbatim("HCalGeom") << "Module " << i << " Copy number " << phib_[i];
0572 #endif
0573 
0574   nzHE_ = static_cast<int>(ize_.size());
0575   nmodHE_ = static_cast<int>(phie_.size());
0576 #ifdef EDM_ML_DEBUG
0577   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters::loadGeometry: " << nzHE_ << " endcap half-sectors";
0578   for (int i = 0; i < nzHE_; i++)
0579     edm::LogVerbatim("HCalGeom") << "Section " << i << " Copy number " << ize_[i];
0580   edm::LogVerbatim("HCalGeom") << "HcalGeomParameters::loadGeometry: " << nmodHE_ << " endcap modules";
0581   for (int i = 0; i < nmodHE_; i++)
0582     edm::LogVerbatim("HCalGeom") << "Module " << i << " Copy number " << phie_[i];
0583 #endif
0584 
0585 #ifdef EDM_ML_DEBUG
0586   edm::LogVerbatim("HCalGeom") << "HO has Z of size " << php.zHO.size();
0587   for (unsigned int kk = 0; kk < php.zHO.size(); kk++)
0588     edm::LogVerbatim("HCalGeom") << "ZHO[" << kk << "] = " << php.zHO[kk];
0589 #endif
0590   if (ibmx > 17 && php.zHO.size() > 4) {
0591     rminHO_ = php.rHB[17] - 100.0;
0592     etaHO_[0] = getEta(0.5 * (php.rHB[17] + php.rHB[18]), php.zHO[1]);
0593     etaHO_[1] = getEta(php.rHB[18] + php.drHB[18], php.zHO[2]);
0594     etaHO_[2] = getEta(php.rHB[18] - php.drHB[18], php.zHO[3]);
0595     etaHO_[3] = getEta(php.rHB[18] + php.drHB[18], php.zHO[4]);
0596   }
0597 #ifdef EDM_ML_DEBUG
0598   edm::LogVerbatim("HCalGeom") << "HO Eta boundaries " << etaHO_[0] << " " << etaHO_[1] << " " << etaHO_[2] << " "
0599                                << etaHO_[3];
0600   edm::LogVerbatim("HCalGeom") << "HO Parameters " << rminHO_ << " " << php.zHO.size();
0601   for (unsigned int i = 0; i < php.zHO.size(); ++i)
0602     edm::LogVerbatim("HCalGeom") << " zho[" << i << "] = " << php.zHO[i];
0603 #endif
0604 }