HcalCommonData/src/HcalSimParametersFromDD.cc

0001 #include "Geometry/HcalCommonData/interface/HcalSimParametersFromDD.h"
0002 #include "CondFormats/GeometryObjects/interface/HcalSimulationParameters.h"
0003 #include "DetectorDescription/Core/interface/DDFilter.h"
0004 #include "DetectorDescription/Core/interface/DDSplit.h"
0005 #include "DetectorDescription/Core/interface/DDValue.h"
0006 #include "DetectorDescription/Core/interface/DDutils.h"
0007 #include "DetectorDescription/DDCMS/interface/BenchmarkGrd.h"
0008 #include "DataFormats/Math/interface/GeantUnits.h"
0009
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include <iostream>
0012 #include <iomanip>
0013
0014 //#define EDM_ML_DEBUG
0015
0016 using namespace geant_units::operators;
0017
0018 bool HcalSimParametersFromDD::build(const DDCompactView* cpv, HcalSimulationParameters& php) {
0019   // Parameters for the fibers
0020   std::string attribute = "Volume";
0021   std::string value = "HF";
0022   DDSpecificsMatchesValueFilter filter1{DDValue(attribute, value, 0)};
0023   DDFilteredView fv1(*cpv, filter1);
0024
0025   // Names of sensitive volumes for HF
0026   php.hfNames_ = getNames(fv1);
0027   int nb(-1);
0028
0029   bool dodet = fv1.firstChild();
0030   if (dodet) {
0031     DDsvalues_type sv(fv1.mergedSpecifics());
0032
0033     // The level positions
0034     nb = -1;
0035     php.hfLevels_ = dbl_to_int(getDDDArray("Levels", sv, nb));
0036
0037     // Attenuation length
0038     nb = -1;
0039     php.attenuationLength_ = getDDDArray("attl", sv, nb);
0040
0041     // Limits on Lambda
0042     nb = 2;
0043     php.lambdaLimits_ = dbl_to_int(getDDDArray("lambLim", sv, nb));
0044
0045     // Fibre Lengths
0046     nb = 0;
0047     php.longFiberLength_ = getDDDArray("LongFL", sv, nb);
0048
0049     nb = 0;
0050     php.shortFiberLength_ = getDDDArray("ShortFL", sv, nb);
0051
0052   } else {
0053     throw cms::Exception("HcalSimParametersFromDD") << "Not found " << value << " for " << attribute << " but needed.";
0054   }
0055
0056   //Parameters for the PMT
0057   value = "HFPMT";
0058   DDSpecificsMatchesValueFilter filter2{DDValue(attribute, value, 0)};
0059   DDFilteredView fv2(*cpv, filter2);
0060   if (fv2.firstChild()) {
0061     DDsvalues_type sv(fv2.mergedSpecifics());
0062     int nb = -1;
0063     std::vector<double> neta = getDDDArray("indexPMTR", sv, nb);
0064     fillPMTs(neta, false, php);
0065     nb = -1;
0066     neta = getDDDArray("indexPMTL", sv, nb);
0067     fillPMTs(neta, true, php);
0068   } else {
0069     throw cms::Exception("HcalSimParametersFromDD") << "Not found " << value << " for " << attribute << " but needed.";
0070   }
0071
0072   //Names of special volumes (HFFibre, HFPMT, HFFibreBundles)
0073   fillNameVector(cpv, attribute, "HFFibre", php.hfFibreNames_);
0074   fillNameVector(cpv, attribute, "HFPMT", php.hfPMTNames_);
0075   fillNameVector(cpv, attribute, "HFFibreBundleStraight", php.hfFibreStraightNames_);
0076   fillNameVector(cpv, attribute, "HFFibreBundleConical", php.hfFibreConicalNames_);
0077
0078   // HCal materials
0079   attribute = "OnlyForHcalSimNumbering";
0080   DDSpecificsHasNamedValueFilter filter3{attribute};
0081   DDFilteredView fv3(*cpv, filter3);
0082   dodet = fv3.firstChild();
0083
0084   while (dodet) {
0085     const DDLogicalPart& log = fv3.logicalPart();
0086     if (!isItHF(log.name().name(), php)) {
0087       bool notIn = true;
0088       for (unsigned int i = 0; i < php.hcalMaterialNames_.size(); ++i) {
0089         if (!strcmp(php.hcalMaterialNames_[i].c_str(), log.material().name().name().c_str())) {
0090           notIn = false;
0091           break;
0092         }
0093       }
0094       if (notIn)
0095         php.hcalMaterialNames_.push_back(log.material().name().name());
0096     }
0097     dodet = fv2.next();
0098   }
0099
0100   return buildParameters(php);
0101 }
0102
0103 bool HcalSimParametersFromDD::build(const cms::DDCompactView& cpv, HcalSimulationParameters& php) {
0104 #ifdef EDM_ML_DEBUG
0105   edm::LogVerbatim("HCalGeom")
0106       << "Inside HcalSimParametersFromDD::build(const cms::DDCompactView*, HcalSimulationParameters&)";
0107 #endif
0108
0109   // HCal materials
0110   const cms::DDFilter filter("OnlyForHcalSimNumbering", "HCAL");
0111   cms::DDFilteredView fv(cpv, filter);
0112
0113   {
0114     BenchmarkGrd counter("HcalSimParametersFromDD get all vectors\n");
0115
0116     // The level positions
0117     php.hfLevels_ = fv.get<std::vector<int> >("hf", "Levels");
0118
0119     // Attenuation length
0120     php.attenuationLength_ = fv.get<std::vector<double> >("hf", "attl");
0121     std::for_each(
0122         php.attenuationLength_.begin(), php.attenuationLength_.end(), [](double& n) { n *= k_ScaleFromDD4hepInv; });
0123
0124     // Limits on Lambda
0125     php.lambdaLimits_ = fv.get<std::vector<int> >("hf", "lambLim");
0126
0127     // Fibre Lengths
0128     php.longFiberLength_ = fv.get<std::vector<double> >("hf", "LongFL");
0129     std::for_each(php.longFiberLength_.begin(), php.longFiberLength_.end(), [](double& n) { n *= k_ScaleFromDD4hep; });
0130     php.shortFiberLength_ = fv.get<std::vector<double> >("hf", "ShortFL");
0131     std::for_each(
0132         php.shortFiberLength_.begin(), php.shortFiberLength_.end(), [](double& n) { n *= k_ScaleFromDD4hep; });
0133
0134     //Parameters for the PMT
0135     std::vector<double> neta = fv.get<std::vector<double> >("hfpmt", "indexPMTR");
0136     fillPMTs(neta, false, php);
0137     neta = fv.get<std::vector<double> >("hfpmt", "indexPMTL");
0138     fillPMTs(neta, true, php);
0139
0140     // Parameters for the fibers
0141     fillNameVector(cpv, "HF", php.hfNames_);
0142
0143     //Names of special volumes (HFFibre, HFPMT, HFFibreBundles)
0144     fillNameVector(cpv, "HFFibre", php.hfFibreNames_);
0145     fillNameVector(cpv, "HFPMT", php.hfPMTNames_);
0146     fillNameVector(cpv, "HFFibreBundleStraight", php.hfFibreStraightNames_);
0147     fillNameVector(cpv, "HFFibreBundleConical", php.hfFibreConicalNames_);
0148   }
0149   {
0150     BenchmarkGrd counter("HcalSimParametersFromDD HCal materials OnlyForHcalSimNumbering, HCAL");
0151
0152     while (fv.firstChild()) {
0153       std::vector<int> copy = fv.copyNos();
0154       // idet = 3 for HB and 4 for HE (convention in the ddalgo code for HB/HE)
0155       int idet = (copy.size() > 1) ? (copy[1] / 1000) : 0;
0156       if ((idet == 3) || (idet == 4)) {
0157         std::string_view matName = dd4hep::dd::noNamespace(fv.materialName());
0158         if (std::find(std::begin(php.hcalMaterialNames_), std::end(php.hcalMaterialNames_), matName) ==
0159             std::end(php.hcalMaterialNames_)) {
0160           php.hcalMaterialNames_.emplace_back(matName);
0161         }
0162       }
0163     }
0164   }
0165
0166   return buildParameters(php);
0167 }
0168
0169 bool HcalSimParametersFromDD::buildParameters(const HcalSimulationParameters& php) {
0170 #ifdef EDM_ML_DEBUG
0171   std::stringstream ss0;
0172   for (unsigned int it = 0; it < php.hfNames_.size(); it++) {
0173     if (it / 10 * 10 == it)
0174       ss0 << "\n";
0175     ss0 << " [" << it << "] " << php.hfNames_[it];
0176   }
0177   edm::LogVerbatim("HCalGeom") << "HFNames: " << php.hfNames_.size() << ": " << ss0.str();
0178
0179   std::stringstream ss1;
0180   for (unsigned int it = 0; it < php.hfLevels_.size(); it++) {
0181     if (it / 10 * 10 == it)
0182       ss1 << "\n";
0183     ss1 << " [" << it << "] " << php.hfLevels_[it];
0184   }
0185   edm::LogVerbatim("HCalGeom") << "HF Volume Levels: " << php.hfLevels_.size() << " hfLevels: " << ss1.str();
0186
0187   std::stringstream ss2;
0188   for (unsigned int it = 0; it < php.attenuationLength_.size(); it++) {
0189     if (it / 10 * 10 == it)
0190       ss2 << "\n";
0191     ss2 << "  " << convertMmToCm(php.attenuationLength_[it]);
0192   }
0193   edm::LogVerbatim("HCalGeom") << "AttenuationLength: " << php.attenuationLength_.size()
0194                                << " attL(1/cm): " << ss2.str();
0195
0196   std::stringstream ss3;
0197   for (unsigned int it = 0; it < php.lambdaLimits_.size(); it++) {
0198     if (it / 10 * 10 == it)
0199       ss3 << "\n";
0200     ss3 << "  " << php.lambdaLimits_[it];
0201   }
0202   edm::LogVerbatim("HCalGeom") << php.lambdaLimits_.size() << " Limits on lambda " << ss3.str();
0203
0204   std::stringstream ss4;
0205   for (unsigned int it = 0; it < php.longFiberLength_.size(); it++) {
0206     if (it / 10 * 10 == it)
0207       ss4 << "\n";
0208     ss4 << "  " << convertMmToCm(php.longFiberLength_[it]);
0209   }
0210   edm::LogVerbatim("HCalGeom") << php.longFiberLength_.size() << " Long Fibre Length(cm):" << ss4.str();
0211
0212   std::stringstream ss5;
0213   for (unsigned int it = 0; it < php.shortFiberLength_.size(); it++) {
0214     if (it / 10 * 10 == it)
0215       ss5 << "\n";
0216     ss5 << "  " << convertMmToCm(php.shortFiberLength_[it]);
0217   }
0218   edm::LogVerbatim("HCalGeom") << php.shortFiberLength_.size() << " Short Fibre Length(cm):" << ss5.str();
0219
0220   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: gets the Index matches for " << php.pmtRight_.size() << " PMTs";
0221   for (unsigned int ii = 0; ii < php.pmtRight_.size(); ii++)
0222     edm::LogVerbatim("HCalGeom") << "rIndexR[" << ii << "] = " << php.pmtRight_[ii] << " fibreR[" << ii
0223                                  << "] = " << php.pmtFiberRight_[ii] << " rIndexL[" << ii << "] = " << php.pmtLeft_[ii]
0224                                  << " fibreL[" << ii << "] = " << php.pmtFiberLeft_[ii];
0225
0226   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: " << php.hfFibreNames_.size() << " names of HFFibre";
0227   for (unsigned int k = 0; k < php.hfFibreNames_.size(); ++k)
0228     edm::LogVerbatim("HCalGeom") << "[" << k << "] " << php.hfFibreNames_[k];
0229
0230   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: " << php.hfPMTNames_.size() << " names of HFPMT";
0231   for (unsigned int k = 0; k < php.hfPMTNames_.size(); ++k)
0232     edm::LogVerbatim("HCalGeom") << "[" << k << "] " << php.hfPMTNames_[k];
0233
0234   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: " << php.hfFibreStraightNames_.size()
0235                                << " names of HFFibreBundleStraight";
0236   for (unsigned int k = 0; k < php.hfFibreStraightNames_.size(); ++k)
0237     edm::LogVerbatim("HCalGeom") << "[" << k << "] " << php.hfFibreStraightNames_[k];
0238
0239   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: " << php.hfFibreConicalNames_.size()
0240                                << " names of FibreBundleConical";
0241   for (unsigned int k = 0; k < php.hfFibreConicalNames_.size(); ++k)
0242     edm::LogVerbatim("HCalGeom") << "[" << k << "] " << php.hfFibreConicalNames_[k];
0243
0244   edm::LogVerbatim("HCalGeom") << "HcalSimParameters: " << php.hcalMaterialNames_.size() << " names of HCAL materials";
0245   for (unsigned int k = 0; k < php.hcalMaterialNames_.size(); ++k)
0246     edm::LogVerbatim("HCalGeom") << "[" << k << "] " << php.hcalMaterialNames_[k];
0247 #endif
0248
0249   return true;
0250 }
0251
0252 void HcalSimParametersFromDD::fillNameVector(const DDCompactView* cpv,
0253                                              const std::string& attribute,
0254                                              const std::string& value,
0255                                              std::vector<std::string>& lvnames) {
0256   DDSpecificsMatchesValueFilter filter{DDValue(attribute, value, 0)};
0257   DDFilteredView fv(*cpv, filter);
0258   lvnames = getNames(fv);
0259 }
0260
0261 void HcalSimParametersFromDD::fillNameVector(const cms::DDCompactView& cpv,
0262                                              const std::string& value,
0263                                              std::vector<std::string>& lvnames) {
0264   {
0265     BenchmarkGrd counter("HcalSimParametersFromDD::fillNameVector");
0266     const cms::DDFilter filter("Volume", value);
0267     cms::DDFilteredView fv(cpv, filter);
0268     lvnames = getNames(fv);
0269   }
0270 }
0271
0272 void HcalSimParametersFromDD::fillPMTs(const std::vector<double>& neta, bool lOrR, HcalSimulationParameters& php) {
0273   {
0274     BenchmarkGrd counter("HcalSimParametersFromDD::fillPMTs");
0275     for (unsigned int ii = 0; ii < neta.size(); ii++) {
0276       int index = static_cast<int>(neta[ii]);
0277       int ir = -1, ifib = -1;
0278       if (index >= 0) {
0279         ir = index / 10;
0280         ifib = index % 10;
0281       }
0282       if (lOrR) {
0283         php.pmtLeft_.push_back(ir);
0284         php.pmtFiberLeft_.push_back(ifib);
0285       } else {
0286         php.pmtRight_.push_back(ir);
0287         php.pmtFiberRight_.push_back(ifib);
0288       }
0289     }
0290   }
0291 }
0292
0293 bool HcalSimParametersFromDD::isItHF(const std::string& name, const HcalSimulationParameters& php) {
0294   if (std::find(std::begin(php.hfNames_), std::end(php.hfNames_), name) != std::end(php.hfNames_))
0295     return true;
0296   if (std::find(std::begin(php.hfFibreNames_), std::end(php.hfFibreNames_), name) != std::end(php.hfFibreNames_))
0297     return true;
0298   if (std::find(std::begin(php.hfPMTNames_), std::end(php.hfPMTNames_), name) != std::end(php.hfPMTNames_))
0299     return true;
0300   if (std::find(std::begin(php.hfFibreStraightNames_), std::end(php.hfFibreStraightNames_), name) !=
0301       std::end(php.hfFibreStraightNames_))
0302     return true;
0303   if (std::find(std::begin(php.hfFibreConicalNames_), std::end(php.hfFibreConicalNames_), name) !=
0304       std::end(php.hfFibreConicalNames_))
0305     return true;
0306
0307   return false;
0308 }
0309
0310 std::vector<std::string> HcalSimParametersFromDD::getNames(DDFilteredView& fv) {
0311   std::vector<std::string> tmp;
0312   bool dodet = fv.firstChild();
0313   while (dodet) {
0314     const DDLogicalPart& log = fv.logicalPart();
0315     bool ok = true;
0316
0317     for (unsigned int i = 0; i < tmp.size(); ++i) {
0318       if (!strcmp(tmp[i].c_str(), log.name().name().c_str())) {
0319         ok = false;
0320         break;
0321       }
0322     }
0323     if (ok)
0324       tmp.push_back(log.name().name());
0325     dodet = fv.next();
0326   }
0327   return tmp;
0328 }
0329
0330 std::vector<std::string> HcalSimParametersFromDD::getNames(cms::DDFilteredView& fv) {
0331   std::vector<std::string> tmp;
0332   while (fv.firstChild()) {
0333     if (std::find(std::begin(tmp), std::end(tmp), fv.name()) == std::end(tmp))
0334       tmp.emplace_back(fv.name());
0335   }
0336   return tmp;
0337 }
0338
0339 std::vector<double> HcalSimParametersFromDD::getDDDArray(const std::string& str, const DDsvalues_type& sv, int& nmin) {
0340 #ifdef EDM_ML_DEBUG
0341   edm::LogVerbatim("HCalGeom") << "HcalSimParametersFromDD::getDDDArray called for " << str << " with nMin " << nmin;
0342 #endif
0343   DDValue value(str);
0344   if (DDfetch(&sv, value)) {
0345 #ifdef EDM_ML_DEBUG
0346     edm::LogVerbatim("HCalGeom") << value;
0347 #endif
0348     const std::vector<double>& fvec = value.doubles();
0349     int nval = fvec.size();
0350     if (nmin > 0) {
0351       if (nval < nmin) {
0352         edm::LogError("HCalGeom") << "# of " << str << " bins " << nval << " < " << nmin << " ==> illegal";
0353         throw cms::Exception("HcalSimParametersFromDD") << "nval < nmin for array " << str << "\n";
0354       }
0355     } else {
0356       if (nval < 1 && nmin != 0) {
0357         edm::LogError("HCalGeom") << "# of " << str << " bins " << nval << " < 1 ==> illegal (nmin=" << nmin << ")";
0358         throw cms::Exception("HcalSimParametersFromDD") << "nval < 1 for array " << str << "\n";
0359       }
0360     }
0361     nmin = nval;
0362     return fvec;
0363   } else {
0364     if (nmin != 0) {
0365       edm::LogError("HCalGeom") << "Cannot get array " << str;
0366       throw cms::Exception("HcalSimParametersFromDD") << "cannot get array " << str << "\n";
0367     } else {
0368       std::vector<double> fvec;
0369       return fvec;
0370     }
0371   }
0372 }