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