File indexing completed on 2024-11-28 23:10:56
0001 #include "CondFormats/GeometryObjects/interface/EcalSimulationParameters.h"
0002 #include "DetectorDescription/Core/interface/DDFilteredView.h"
0003 #include "DetectorDescription/Core/interface/DDFilter.h"
0004 #include "DetectorDescription/Core/interface/DDValue.h"
0005 #include "DetectorDescription/Core/interface/DDutils.h"
0006 #include "DetectorDescription/DDCMS/interface/DDFilteredView.h"
0007 #include "Geometry/EcalCommonData/interface/EcalSimParametersFromDD.h"
0008 #include <iostream>
0009 #include <iomanip>
0010
0011
0012
0013 template <typename T>
0014 void myPrint(std::string value, const std::vector<T>& vec) {
0015 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD: " << vec.size() << " entries for " << value << ":";
0016 unsigned int i(0);
0017 for (const auto& e : vec) {
0018 edm::LogVerbatim("EcalGeom") << " (" << i << ") " << e;
0019 ++i;
0020 }
0021 }
0022
0023 bool EcalSimParametersFromDD::build(const DDCompactView* cpv, const std::string& name, EcalSimulationParameters& php) {
0024 #ifdef EDM_ML_DEBUG
0025 edm::LogVerbatim("EcalGeom")
0026 << "Inside EcalSimParametersFromDD::build(const DDCompactView*, const std::string&, EcalSimulationParameters&)";
0027 #endif
0028
0029 std::string attribute = "ReadOutName";
0030 DDSpecificsMatchesValueFilter filter{DDValue(attribute, name, 0)};
0031 DDFilteredView fv(*cpv, filter);
0032 bool dodet = fv.firstChild();
0033 DDsvalues_type sv(fv.mergedSpecifics());
0034
0035
0036 php.useWeight_ = true;
0037 std::vector<double> tempD = getDDDArray("EnergyWeight", sv);
0038 #ifdef EDM_ML_DEBUG
0039 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " EnergyWeight " << tempD.empty();
0040 #endif
0041 if (!tempD.empty()) {
0042 if (tempD[0] < 0.1)
0043 php.useWeight_ = false;
0044 }
0045 tempD = getDDDArray("nxtalEta", sv);
0046 #ifdef EDM_ML_DEBUG
0047 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " nxtalEta " << tempD.empty();
0048 #endif
0049 if (tempD.empty())
0050 php.nxtalEta_ = 0;
0051 else
0052 php.nxtalEta_ = static_cast<int>(tempD[0]);
0053 tempD = getDDDArray("nxtalPhi", sv);
0054 #ifdef EDM_ML_DEBUG
0055 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " nxtalPhi " << tempD.empty();
0056 #endif
0057 if (tempD.empty())
0058 php.nxtalPhi_ = 0;
0059 else
0060 php.nxtalPhi_ = static_cast<int>(tempD[0]);
0061 tempD = getDDDArray("PhiBaskets", sv);
0062 #ifdef EDM_ML_DEBUG
0063 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " PhiBaskets " << tempD.empty();
0064 #endif
0065 if (tempD.empty())
0066 php.phiBaskets_ = 0;
0067 else
0068 php.phiBaskets_ = static_cast<int>(tempD[0]);
0069 php.etaBaskets_ = dbl_to_int(getDDDArray("EtaBaskets", sv));
0070 #ifdef EDM_ML_DEBUG
0071 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " EtaBaskets " << tempD.empty();
0072 #endif
0073 tempD = getDDDArray("ncrys", sv);
0074 if (tempD.empty())
0075 php.ncrys_ = 0;
0076 else
0077 php.ncrys_ = static_cast<int>(tempD[0]);
0078 tempD = getDDDArray("nmods", sv);
0079 #ifdef EDM_ML_DEBUG
0080 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " ncrys " << tempD.empty();
0081 #endif
0082 if (tempD.empty())
0083 php.nmods_ = 0;
0084 else
0085 php.nmods_ = static_cast<int>(tempD[0]);
0086
0087 std::vector<std::string> tempS = getStringArray("Depth1Name", sv);
0088 #ifdef EDM_ML_DEBUG
0089 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " Depth1Name " << tempS.empty();
0090 #endif
0091 if (!tempS.empty())
0092 php.depth1Name_ = tempS[0];
0093 else
0094 php.depth1Name_ = " ";
0095 tempS = getStringArray("Depth2Name", sv);
0096 #ifdef EDM_ML_DEBUG
0097 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:" << name << " Depth2Name " << tempS.empty();
0098 #endif
0099 if (!tempS.empty())
0100 php.depth2Name_ = tempS[0];
0101 else
0102 php.depth2Name_ = " ";
0103
0104
0105 while (dodet) {
0106 if (std::find(php.lvNames_.begin(), php.lvNames_.end(), fv.logicalPart().name().name()) == php.lvNames_.end()) {
0107 php.matNames_.emplace_back(fv.logicalPart().material().name().name());
0108 php.lvNames_.emplace_back(fv.logicalPart().name().name());
0109 const DDSolid& sol = fv.logicalPart().solid();
0110 const std::vector<double>& paras = sol.parameters();
0111 double dz = (sol.shape() == DDSolidShape::ddtrap) ? (2 * paras[0]) : 0.0;
0112 php.dzs_.emplace_back(dz);
0113 }
0114 dodet = fv.next();
0115 }
0116 return this->buildParameters(php);
0117 }
0118
0119 bool EcalSimParametersFromDD::build(const cms::DDCompactView* cpv,
0120 const std::string& name,
0121 EcalSimulationParameters& php) {
0122 #ifdef EDM_ML_DEBUG
0123 edm::LogVerbatim("EcalGeom") << "Inside EcalSimParametersFromDD::build(const cms::DDCompactView*, const std::string, "
0124 "EcalSimulationParameters&)";
0125 #endif
0126
0127 std::string attribute = "ReadOutName";
0128 cms::DDFilteredView fv(cpv->detector(), cpv->detector()->worldVolume());
0129
0130
0131 std::string specName = ((name == "EcalHitsEE") ? "ecal_ee" : ((name == "EcalHitsES") ? "ecal_sf" : "ecal_eb"));
0132
0133 php.useWeight_ = true;
0134 std::vector<double> tempD = fv.get<std::vector<double> >(specName, "EnergyWeight");
0135 if (!tempD.empty()) {
0136 if (tempD[0] < 0.1)
0137 php.useWeight_ = false;
0138 }
0139 tempD = fv.get<std::vector<double> >(specName, "nxtalEta");
0140 if (tempD.empty())
0141 php.nxtalEta_ = 0;
0142 else
0143 php.nxtalEta_ = static_cast<int>(tempD[0]);
0144 tempD = fv.get<std::vector<double> >(specName, "nxtalPhi");
0145 if (tempD.empty())
0146 php.nxtalPhi_ = 0;
0147 else
0148 php.nxtalPhi_ = static_cast<int>(tempD[0]);
0149 tempD = fv.get<std::vector<double> >(specName, "PhiBaskets");
0150 if (tempD.empty())
0151 php.phiBaskets_ = 0;
0152 else
0153 php.phiBaskets_ = static_cast<int>(tempD[0]);
0154 php.etaBaskets_ = dbl_to_int(fv.get<std::vector<double> >(specName, "EtaBaskets"));
0155 tempD = fv.get<std::vector<double> >(specName, "ncrys");
0156 if (tempD.empty())
0157 php.ncrys_ = 0;
0158 else
0159 php.ncrys_ = static_cast<int>(tempD[0]);
0160 tempD = fv.get<std::vector<double> >(specName, "nmods");
0161 if (tempD.empty())
0162 php.nmods_ = 0;
0163 else
0164 php.nmods_ = static_cast<int>(tempD[0]);
0165
0166 std::vector<std::string> tempS = fv.get<std::vector<std::string> >(specName, "Depth1Name");
0167 if (!tempS.empty())
0168 php.depth1Name_ = tempS[0];
0169 else
0170 php.depth1Name_ = " ";
0171 tempS = fv.get<std::vector<std::string> >(specName, "Depth2Name");
0172 if (!tempS.empty())
0173 php.depth2Name_ = tempS[0];
0174 else
0175 php.depth2Name_ = " ";
0176
0177
0178 cms::DDSpecParRefs refs;
0179 const cms::DDSpecParRegistry& mypar = cpv->specpars();
0180 mypar.filter(refs, attribute, name);
0181 fv.mergedSpecifics(refs);
0182 while (fv.firstChild()) {
0183 const std::string name{dd4hep::dd::noNamespace(fv.name()).data(), dd4hep::dd::noNamespace(fv.name()).size()};
0184 const std::string matName{dd4hep::dd::noNamespace(fv.materialName()).data(),
0185 dd4hep::dd::noNamespace(fv.materialName()).size()};
0186 if (std::find(php.lvNames_.begin(), php.lvNames_.end(), name) == php.lvNames_.end()) {
0187 php.matNames_.emplace_back(matName);
0188 php.lvNames_.emplace_back(name);
0189 const std::vector<double>& paras = fv.parameters();
0190 double dz = (dd4hep::isA<dd4hep::Trap>(fv.solid())) ? ((2.0 * paras[0]) / dd4hep::mm) : 0.0;
0191 php.dzs_.emplace_back(dz);
0192 }
0193 };
0194
0195 return this->buildParameters(php);
0196 }
0197
0198 bool EcalSimParametersFromDD::buildParameters(const EcalSimulationParameters& php) {
0199 #ifdef EDM_ML_DEBUG
0200 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:: nxtalEta:" << php.nxtalEta_
0201 << " nxtalPhi:" << php.nxtalPhi_ << " phiBaskets:" << php.phiBaskets_
0202 << " ncrys:" << php.ncrys_ << " nmods: " << php.nmods_ << " useWeight:" << php.useWeight_
0203 << " DeothNames:" << php.depth1Name_ << ":" << php.depth2Name_;
0204 myPrint("etaBaskets", php.etaBaskets_);
0205 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:: " << php.lvNames_.size() << " lvNames, "
0206 << php.matNames_.size() << " matNames and " << php.dzs_.size() << "dzs";
0207 #endif
0208
0209 return true;
0210 }
0211
0212 std::vector<double> EcalSimParametersFromDD::getDDDArray(const std::string& str, const DDsvalues_type& sv) {
0213 #ifdef EDM_ML_DEBUG
0214 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:getDDDArray called for " << str;
0215 #endif
0216 DDValue value(str);
0217 if (DDfetch(&sv, value)) {
0218 #ifdef EDM_ML_DEBUG
0219 edm::LogVerbatim("EcalGeom") << value;
0220 #endif
0221 const std::vector<double>& fvec = value.doubles();
0222 return fvec;
0223 } else {
0224 std::vector<double> fvec;
0225 return fvec;
0226 }
0227 }
0228
0229 std::vector<std::string> EcalSimParametersFromDD::getStringArray(const std::string& str, const DDsvalues_type& sv) {
0230 #ifdef EDM_ML_DEBUG
0231 edm::LogVerbatim("EcalGeom") << "EcalSimParametersFromDD:getStringArray called for " << str;
0232 #endif
0233 DDValue value(str);
0234 if (DDfetch(&sv, value)) {
0235 #ifdef EDM_ML_DEBUG
0236 edm::LogVerbatim("EcalGeom") << value;
0237 #endif
0238 const std::vector<std::string>& fvec = value.strings();
0239 return fvec;
0240 } else {
0241 std::vector<std::string> fvec;
0242 return fvec;
0243 }
0244 }