File indexing completed on 2024-04-06 12:14:59
0001
0002
0003
0004
0005
0006 #include <algorithm>
0007 #include <cmath>
0008 #include <map>
0009 #include <string>
0010 #include <unordered_set>
0011 #include <vector>
0012
0013 #include "DataFormats/Math/interface/angle_units.h"
0014 #include "DetectorDescription/Core/interface/DDAlgorithm.h"
0015 #include "DetectorDescription/Core/interface/DDAlgorithmFactory.h"
0016 #include "DetectorDescription/Core/interface/DDCurrentNamespace.h"
0017 #include "DetectorDescription/Core/interface/DDLogicalPart.h"
0018 #include "DetectorDescription/Core/interface/DDMaterial.h"
0019 #include "DetectorDescription/Core/interface/DDSolid.h"
0020 #include "DetectorDescription/Core/interface/DDSplit.h"
0021 #include "DetectorDescription/Core/interface/DDTypes.h"
0022 #include "DetectorDescription/Core/interface/DDutils.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 #include "FWCore/PluginManager/interface/PluginFactory.h"
0025 #include "Geometry/HGCalCommonData/interface/HGCalGeomTools.h"
0026 #include "Geometry/HGCalCommonData/interface/HGCalParameters.h"
0027 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0028
0029
0030 using namespace angle_units::operators;
0031
0032 class DDHGCalModuleAlgo : public DDAlgorithm {
0033 public:
0034
0035 DDHGCalModuleAlgo();
0036
0037 void initialize(const DDNumericArguments& nArgs,
0038 const DDVectorArguments& vArgs,
0039 const DDMapArguments& mArgs,
0040 const DDStringArguments& sArgs,
0041 const DDStringVectorArguments& vsArgs) override;
0042 void execute(DDCompactView& cpv) override;
0043
0044 protected:
0045 void constructLayers(const DDLogicalPart&, DDCompactView& cpv);
0046 double rMax(double z);
0047 void positionSensitive(DDLogicalPart& glog, double rin, double rout, DDCompactView& cpv);
0048
0049 private:
0050 static constexpr double tol_ = 0.00001;
0051
0052 std::vector<std::string> wafer_;
0053 std::vector<std::string> materials_;
0054 std::vector<std::string> names_;
0055 std::vector<double> thick_;
0056 std::vector<int> copyNumber_;
0057 std::vector<int> layers_;
0058 std::vector<double> layerThick_;
0059 std::vector<int> layerType_;
0060 std::vector<int> layerSense_;
0061 double zMinBlock_;
0062 double rMaxFine_;
0063 double waferW_;
0064 double waferGap_;
0065 int sectors_;
0066 std::vector<double> slopeB_;
0067 std::vector<double> slopeT_;
0068 std::vector<double> zFront_;
0069 std::vector<double> rMaxFront_;
0070 std::string idNameSpace_;
0071 std::unordered_set<int> copies_;
0072 };
0073
0074 DDHGCalModuleAlgo::DDHGCalModuleAlgo() {
0075 #ifdef EDM_ML_DEBUG
0076 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: Creating an instance";
0077 #endif
0078 }
0079
0080 void DDHGCalModuleAlgo::initialize(const DDNumericArguments& nArgs,
0081 const DDVectorArguments& vArgs,
0082 const DDMapArguments&,
0083 const DDStringArguments& sArgs,
0084 const DDStringVectorArguments& vsArgs) {
0085 wafer_ = vsArgs["WaferName"];
0086 #ifdef EDM_ML_DEBUG
0087 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << wafer_.size() << " wafers";
0088 for (unsigned int i = 0; i < wafer_.size(); ++i)
0089 edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << "] " << wafer_[i];
0090 #endif
0091 materials_ = vsArgs["MaterialNames"];
0092 names_ = vsArgs["VolumeNames"];
0093 thick_ = vArgs["Thickness"];
0094 copyNumber_.resize(materials_.size(), 1);
0095 #ifdef EDM_ML_DEBUG
0096 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << materials_.size() << " types of volumes";
0097 for (unsigned int i = 0; i < names_.size(); ++i)
0098 edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << names_[i] << " of thickness " << thick_[i]
0099 << " filled with " << materials_[i] << " first copy number " << copyNumber_[i];
0100 #endif
0101 layers_ = dbl_to_int(vArgs["Layers"]);
0102 layerThick_ = vArgs["LayerThick"];
0103 #ifdef EDM_ML_DEBUG
0104 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << layers_.size() << " blocks";
0105 for (unsigned int i = 0; i < layers_.size(); ++i)
0106 edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] of thickness " << layerThick_[i] << " with " << layers_[i]
0107 << " layers";
0108 #endif
0109 layerType_ = dbl_to_int(vArgs["LayerType"]);
0110 layerSense_ = dbl_to_int(vArgs["LayerSense"]);
0111 #ifdef EDM_ML_DEBUG
0112 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << layerType_.size() << " layers";
0113 for (unsigned int i = 0; i < layerType_.size(); ++i)
0114 edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerType_[i] << " sensitive class "
0115 << layerSense_[i];
0116 #endif
0117 zMinBlock_ = nArgs["zMinBlock"];
0118 rMaxFine_ = nArgs["rMaxFine"];
0119 waferW_ = nArgs["waferW"];
0120 waferGap_ = nArgs["waferGap"];
0121 sectors_ = (int)(nArgs["Sectors"]);
0122 #ifdef EDM_ML_DEBUG
0123 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: zStart " << zMinBlock_ << " rFineCoarse " << rMaxFine_
0124 << " wafer width " << waferW_ << " gap among wafers " << waferGap_ << " sectors "
0125 << sectors_;
0126 #endif
0127 slopeB_ = vArgs["SlopeBottom"];
0128 slopeT_ = vArgs["SlopeTop"];
0129 zFront_ = vArgs["ZFront"];
0130 rMaxFront_ = vArgs["RMaxFront"];
0131 #ifdef EDM_ML_DEBUG
0132 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: Bottom slopes " << slopeB_[0] << ":" << slopeB_[1] << " and "
0133 << slopeT_.size() << " slopes for top";
0134 for (unsigned int i = 0; i < slopeT_.size(); ++i)
0135 edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFront_[i] << " Rmax " << rMaxFront_[i] << " Slope "
0136 << slopeT_[i];
0137 #endif
0138 idNameSpace_ = DDCurrentNamespace::ns();
0139 #ifdef EDM_ML_DEBUG
0140 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: NameSpace " << idNameSpace_;
0141 #endif
0142 }
0143
0144
0145
0146
0147
0148 void DDHGCalModuleAlgo::execute(DDCompactView& cpv) {
0149 #ifdef EDM_ML_DEBUG
0150 edm::LogVerbatim("HGCalGeom") << "==>> Constructing DDHGCalModuleAlgo...";
0151 #endif
0152 copies_.clear();
0153 constructLayers(parent(), cpv);
0154 #ifdef EDM_ML_DEBUG
0155 edm::LogVerbatim("HGCalGeom") << copies_.size() << " different wafer copy numbers";
0156 #endif
0157 copies_.clear();
0158 #ifdef EDM_ML_DEBUG
0159 edm::LogVerbatim("HGCalGeom") << "<<== End of DDHGCalModuleAlgo construction";
0160 #endif
0161 }
0162
0163 void DDHGCalModuleAlgo::constructLayers(const DDLogicalPart& module, DDCompactView& cpv) {
0164 #ifdef EDM_ML_DEBUG
0165 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: \t\tInside Layers";
0166 #endif
0167 double zi(zMinBlock_);
0168 int laymin(0);
0169 const double tol(0.01);
0170 for (unsigned int i = 0; i < layers_.size(); i++) {
0171 double zo = zi + layerThick_[i];
0172 double routF = rMax(zi);
0173 int laymax = laymin + layers_[i];
0174 double zz = zi;
0175 double thickTot(0);
0176 for (int ly = laymin; ly < laymax; ++ly) {
0177 int ii = layerType_[ly];
0178 int copy = copyNumber_[ii];
0179 double rinB = (layerSense_[ly] == 0) ? (zo * slopeB_[0]) : (zo * slopeB_[1]);
0180 zz += (0.5 * thick_[ii]);
0181 thickTot += thick_[ii];
0182
0183 std::string name = "HGCal" + names_[ii] + std::to_string(copy);
0184 #ifdef EDM_ML_DEBUG
0185 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: Layer " << ly << ":" << ii << " Front " << zi << ", "
0186 << routF << " Back " << zo << ", " << rinB << " superlayer thickness "
0187 << layerThick_[i];
0188 #endif
0189 DDName matName(DDSplit(materials_[ii]).first, DDSplit(materials_[ii]).second);
0190 DDMaterial matter(matName);
0191 DDLogicalPart glog;
0192 if (layerSense_[ly] == 0) {
0193 double alpha = 1._pi / sectors_;
0194 double rmax = routF * cos(alpha) - tol;
0195 std::vector<double> pgonZ, pgonRin, pgonRout;
0196 pgonZ.emplace_back(-0.5 * thick_[ii]);
0197 pgonZ.emplace_back(0.5 * thick_[ii]);
0198 pgonRin.emplace_back(rinB);
0199 pgonRin.emplace_back(rinB);
0200 pgonRout.emplace_back(rmax);
0201 pgonRout.emplace_back(rmax);
0202 DDSolid solid =
0203 DDSolidFactory::polyhedra(DDName(name, idNameSpace_), sectors_, -alpha, 2._pi, pgonZ, pgonRin, pgonRout);
0204 glog = DDLogicalPart(solid.ddname(), matter, solid);
0205 #ifdef EDM_ML_DEBUG
0206 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << solid.name() << " polyhedra of " << sectors_
0207 << " sectors covering " << convertRadToDeg(-alpha) << ":"
0208 << (360.0 + convertRadToDeg(-alpha)) << " with " << pgonZ.size() << " sections";
0209 for (unsigned int k = 0; k < pgonZ.size(); ++k)
0210 edm::LogVerbatim("HGCalGeom") << "[" << k << "] z " << pgonZ[k] << " R " << pgonRin[k] << ":" << pgonRout[k];
0211 #endif
0212 } else {
0213 DDSolid solid = DDSolidFactory::tubs(DDName(name, idNameSpace_), 0.5 * thick_[ii], rinB, routF, 0.0, 2._pi);
0214 glog = DDLogicalPart(solid.ddname(), matter, solid);
0215 #ifdef EDM_ML_DEBUG
0216 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << solid.name() << " Tubs made of " << matName
0217 << " of dimensions " << rinB << ", " << routF << ", " << 0.5 * thick_[ii]
0218 << ", 0.0, 360.0";
0219 #endif
0220 positionSensitive(glog, rinB, routF, cpv);
0221 }
0222 DDTranslation r1(0, 0, zz);
0223 DDRotation rot;
0224 cpv.position(glog, module, copy, r1, rot);
0225 ++copyNumber_[ii];
0226 #ifdef EDM_ML_DEBUG
0227 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << glog.name() << " number " << copy << " positioned in "
0228 << module.name() << " at " << r1 << " with " << rot;
0229 #endif
0230 zz += (0.5 * thick_[ii]);
0231 }
0232 zi = zo;
0233 laymin = laymax;
0234 if (fabs(thickTot - layerThick_[i]) > tol_) {
0235 if (thickTot > layerThick_[i]) {
0236 edm::LogError("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " is smaller than thickness "
0237 << thickTot << " of all its components **** ERROR ****\n";
0238 } else {
0239 edm::LogWarning("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " does not match with "
0240 << thickTot << " of the components\n";
0241 }
0242 }
0243 }
0244 }
0245
0246 double DDHGCalModuleAlgo::rMax(double z) {
0247 double r(0);
0248 #ifdef EDM_ML_DEBUG
0249 unsigned int ik(0);
0250 #endif
0251 for (unsigned int k = 0; k < slopeT_.size(); ++k) {
0252 if (z < zFront_[k])
0253 break;
0254 r = rMaxFront_[k] + (z - zFront_[k]) * slopeT_[k];
0255 #ifdef EDM_ML_DEBUG
0256 ik = k;
0257 #endif
0258 }
0259 #ifdef EDM_ML_DEBUG
0260 edm::LogVerbatim("HGCalGeom") << "rMax : " << z << ":" << ik << ":" << r;
0261 #endif
0262 return r;
0263 }
0264
0265 void DDHGCalModuleAlgo::positionSensitive(DDLogicalPart& glog, double rin, double rout, DDCompactView& cpv) {
0266 double ww = (waferW_ + waferGap_);
0267 double dx = 0.5 * ww;
0268 double dy = 3.0 * dx * tan(30._deg);
0269 double rr = 2.0 * dx * tan(30._deg);
0270 int ncol = (int)(2.0 * rout / ww) + 1;
0271 int nrow = (int)(rout / (ww * tan(30._deg))) + 1;
0272 int incm(0), inrm(0);
0273 #ifdef EDM_ML_DEBUG
0274 int kount(0);
0275 edm::LogVerbatim("HGCalGeom") << glog.ddname() << " rout " << rout << " Row " << nrow << " Column " << ncol;
0276 #endif
0277 for (int nr = -nrow; nr <= nrow; ++nr) {
0278 int inr = (nr >= 0) ? nr : -nr;
0279 for (int nc = -ncol; nc <= ncol; ++nc) {
0280 int inc = (nc >= 0) ? nc : -nc;
0281 if (inr % 2 == inc % 2) {
0282 double xpos = nc * dx;
0283 double ypos = nr * dy;
0284 auto const& corner = HGCalGeomTools::waferCorner(xpos, ypos, dx, rr, rin, rout, true);
0285 if (corner.first == (int)(HGCalParameters::k_CornerSize)) {
0286 double rpos = std::sqrt(xpos * xpos + ypos * ypos);
0287 DDTranslation tran(xpos, ypos, 0.0);
0288 DDRotation rotation;
0289 int copy = HGCalTypes::packTypeUV(0, nc, nr);
0290 DDName name = (rpos < rMaxFine_) ? DDName(DDSplit(wafer_[0]).first, DDSplit(wafer_[0]).second)
0291 : DDName(DDSplit(wafer_[1]).first, DDSplit(wafer_[1]).second);
0292 cpv.position(name, glog.ddname(), copy, tran, rotation);
0293 if (inc > incm)
0294 incm = inc;
0295 if (inr > inrm)
0296 inrm = inr;
0297 #ifdef EDM_ML_DEBUG
0298 kount++;
0299 #endif
0300 if (copies_.count(copy) == 0)
0301 copies_.insert(copy);
0302 #ifdef EDM_ML_DEBUG
0303 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: " << name << " number " << copy << " positioned in "
0304 << glog.ddname() << " at " << tran << " with " << rotation;
0305 #endif
0306 }
0307 }
0308 }
0309 }
0310 #ifdef EDM_ML_DEBUG
0311 edm::LogVerbatim("HGCalGeom") << "DDHGCalModuleAlgo: # of columns " << incm << " # of rows " << inrm << " and "
0312 << kount << " wafers for " << glog.ddname();
0313 #endif
0314 }
0315
0316 DEFINE_EDM_PLUGIN(DDAlgorithmFactory, DDHGCalModuleAlgo, "hgcal:DDHGCalModuleAlgo");