diff --git a/Geometry/HGCalCommonData/interface/HGCalGeometryMode.h b/Geometry/HGCalCommonData/interface/HGCalGeometryMode.h index a7e703f968867..d60efaa83a542 100644 --- a/Geometry/HGCalCommonData/interface/HGCalGeometryMode.h +++ b/Geometry/HGCalCommonData/interface/HGCalGeometryMode.h @@ -22,7 +22,15 @@ class HGCalStringToEnumParser { }; namespace HGCalGeometryMode { - enum GeometryMode { Square = 0, Hexagon = 1, HexagonFull = 2, Hexagon8 = 3, Hexagon8Full = 4, Trapezoid = 5 }; + enum GeometryMode { + Square = 0, + Hexagon = 1, + HexagonFull = 2, + Hexagon8 = 3, + Hexagon8Full = 4, + Trapezoid = 5, + HexagonFullPart = 6 + }; enum WaferMode { Polyhedra = 0, ExtrudedPolygon = 1 }; } // namespace HGCalGeometryMode diff --git a/Geometry/HGCalCommonData/interface/HGCalWaferType.h b/Geometry/HGCalCommonData/interface/HGCalWaferType.h index 671f940a5ab78..f24bfef885203 100644 --- a/Geometry/HGCalCommonData/interface/HGCalWaferType.h +++ b/Geometry/HGCalCommonData/interface/HGCalWaferType.h @@ -27,6 +27,7 @@ class HGCalWaferType { double cutFracArea); ~HGCalWaferType(); int getType(double xpos, double ypos, double zpos); + static int getType(int index, const std::vector& indices, const std::vector& types); std::pair rLimits(double zpos); private: diff --git a/Geometry/HGCalCommonData/plugins/DDHGCalEEFileAlgo.cc b/Geometry/HGCalCommonData/plugins/DDHGCalEEFileAlgo.cc new file mode 100644 index 0000000000000..878931820fbf3 --- /dev/null +++ b/Geometry/HGCalCommonData/plugins/DDHGCalEEFileAlgo.cc @@ -0,0 +1,422 @@ +/////////////////////////////////////////////////////////////////////////////// +// File: DDHGCalEEFileAlgo.cc +// Description: Geometry factory class for HGCal (EE and HESil) using +// information from the file +/////////////////////////////////////////////////////////////////////////////// + +#include "DataFormats/Math/interface/GeantUnits.h" +#include "DetectorDescription/Core/interface/DDAlgorithm.h" +#include "DetectorDescription/Core/interface/DDAlgorithmFactory.h" +#include "DetectorDescription/Core/interface/DDCurrentNamespace.h" +#include "DetectorDescription/Core/interface/DDLogicalPart.h" +#include "DetectorDescription/Core/interface/DDMaterial.h" +#include "DetectorDescription/Core/interface/DDSolid.h" +#include "DetectorDescription/Core/interface/DDSplit.h" +#include "DetectorDescription/Core/interface/DDTypes.h" +#include "DetectorDescription/Core/interface/DDutils.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/PluginManager/interface/PluginFactory.h" +#include "Geometry/HGCalCommonData/interface/HGCalGeomTools.h" +#include "Geometry/HGCalCommonData/interface/HGCalParameters.h" +#include "Geometry/HGCalCommonData/interface/HGCalWaferIndex.h" +#include "Geometry/HGCalCommonData/interface/HGCalWaferType.h" + +#include +#include +#include +#include +#include + +//#define EDM_ML_DEBUG +using namespace geant_units::operators; + +class DDHGCalEEFileAlgo : public DDAlgorithm { +public: + DDHGCalEEFileAlgo(); + + void initialize(const DDNumericArguments& nArgs, + const DDVectorArguments& vArgs, + const DDMapArguments& mArgs, + const DDStringArguments& sArgs, + const DDStringVectorArguments& vsArgs) override; + void execute(DDCompactView& cpv) override; + +protected: + void constructLayers(const DDLogicalPart&, DDCompactView& cpv); + void positionSensitive( + const DDLogicalPart& glog, double rin, double rout, double zpos, int layertype, int layer, DDCompactView& cpv); + +private: + HGCalGeomTools geomTools_; + + static constexpr double tol1_ = 0.01; + static constexpr double tol2_ = 0.00001; + + std::vector wafers_; // Wafers + std::vector materials_; // Materials + std::vector names_; // Names + std::vector thick_; // Thickness of the material + std::vector copyNumber_; // Initial copy numbers + std::vector layers_; // Number of layers in a section + std::vector layerThick_; // Thickness of each section + std::vector layerType_; // Type of the layer + std::vector layerSense_; // Content of a layer (sensitive?) + std::vector layerCenter_; // Centering of the wafers + int firstLayer_; // Copy # of the first sensitive layer + int absorbMode_; // Absorber mode + double zMinBlock_; // Starting z-value of the block + std::vector rad100to200_; // Parameters for 120-200mum trans. + std::vector rad200to300_; // Parameters for 200-300mum trans. + double zMinRadPar_; // Minimum z for radius parametriz. + std::vector waferIndex_; // Wafer index for the types + std::vector waferTypes_; // Wafer types + int choiceType_; // Type of parametrization to be used + int nCutRadPar_; // Cut off threshold for corners + double fracAreaMin_; // Minimum fractional conatined area + double waferSize_; // Width of the wafer + double waferSepar_; // Sensor separation + int sectors_; // Sectors + std::vector slopeB_; // Slope at the lower R + std::vector zFrontB_; // Starting Z values for the slopes + std::vector rMinFront_; // Corresponding rMin's + std::vector slopeT_; // Slopes at the larger R + std::vector zFrontT_; // Starting Z values for the slopes + std::vector rMaxFront_; // Corresponding rMax's + std::string nameSpace_; // Namespace of this and ALL sub-parts + std::unordered_set copies_; // List of copy #'s + double alpha_, cosAlpha_; +}; + +DDHGCalEEFileAlgo::DDHGCalEEFileAlgo() { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: Creating an instance"; +#endif +} + +void DDHGCalEEFileAlgo::initialize(const DDNumericArguments& nArgs, + const DDVectorArguments& vArgs, + const DDMapArguments&, + const DDStringArguments& sArgs, + const DDStringVectorArguments& vsArgs) { + wafers_ = vsArgs["WaferNames"]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << wafers_.size() << " wafers"; + for (unsigned int i = 0; i < wafers_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << "] " << wafers_[i]; +#endif + materials_ = vsArgs["MaterialNames"]; + names_ = vsArgs["VolumeNames"]; + thick_ = vArgs["Thickness"]; + copyNumber_.resize(materials_.size(), 1); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << materials_.size() << " types of volumes"; + for (unsigned int i = 0; i < names_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << names_[i] << " of thickness " << thick_[i] + << " filled with " << materials_[i] << " first copy number " << copyNumber_[i]; +#endif + layers_ = dbl_to_int(vArgs["Layers"]); + layerThick_ = vArgs["LayerThick"]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "There are " << layers_.size() << " blocks"; + for (unsigned int i = 0; i < layers_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] of thickness " << layerThick_[i] << " with " << layers_[i] + << " layers"; +#endif + layerType_ = dbl_to_int(vArgs["LayerType"]); + layerSense_ = dbl_to_int(vArgs["LayerSense"]); + firstLayer_ = (int)(nArgs["FirstLayer"]); + absorbMode_ = (int)(nArgs["AbsorberMode"]); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "First Layer " << firstLayer_ << " and " + << "Absober mode " << absorbMode_; +#endif + layerCenter_ = dbl_to_int(vArgs["LayerCenter"]); +#ifdef EDM_ML_DEBUG + for (unsigned int i = 0; i < layerCenter_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "LayerCenter [" << i << "] " << layerCenter_[i]; +#endif + if (firstLayer_ > 0) { + for (unsigned int i = 0; i < layerType_.size(); ++i) { + if (layerSense_[i] > 0) { + int ii = layerType_[i]; + copyNumber_[ii] = firstLayer_; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "First copy number for layer type " << i << ":" << ii << " with " + << materials_[ii] << " changed to " << copyNumber_[ii]; +#endif + break; + } + } + } +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "There are " << layerType_.size() << " layers"; + for (unsigned int i = 0; i < layerType_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerType_[i] << " sensitive class " + << layerSense_[i]; +#endif + zMinBlock_ = nArgs["zMinBlock"]; + rad100to200_ = vArgs["rad100to200"]; + rad200to300_ = vArgs["rad200to300"]; + zMinRadPar_ = nArgs["zMinForRadPar"]; + choiceType_ = (int)(nArgs["choiceType"]); + nCutRadPar_ = (int)(nArgs["nCornerCut"]); + fracAreaMin_ = nArgs["fracAreaMin"]; + waferSize_ = nArgs["waferSize"]; + waferSepar_ = nArgs["SensorSeparation"]; + sectors_ = (int)(nArgs["Sectors"]); + alpha_ = (1._pi) / sectors_; + cosAlpha_ = cos(alpha_); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "zStart " << zMinBlock_ << " radius for wafer type separation uses " + << rad100to200_.size() << " parameters; zmin " << zMinRadPar_ << " cutoff " + << choiceType_ << ":" << nCutRadPar_ << ":" << fracAreaMin_ << " wafer width " + << waferSize_ << " separations " << waferSepar_ << " sectors " << sectors_ << ":" + << convertRadToDeg(alpha_) << ":" << cosAlpha_; + for (unsigned int k = 0; k < rad100to200_.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] 100-200 " << rad100to200_[k] << " 200-300 " << rad200to300_[k]; +#endif + waferIndex_ = dbl_to_int(vArgs["WaferIndex"]); + waferTypes_ = dbl_to_int(vArgs["WaferTypes"]); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "waferTypes with " << waferTypes_.size() << " entries"; + for (unsigned int k = 0; k < waferTypes_.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << waferIndex_[k] << " (" + << HGCalWaferIndex::waferLayer(waferIndex_[k]) << ", " + << HGCalWaferIndex::waferU(waferIndex_[k]) << ", " + << HGCalWaferIndex::waferV(waferIndex_[k]) << ") : " << waferTypes_[k]; +#endif + slopeB_ = vArgs["SlopeBottom"]; + zFrontB_ = vArgs["ZFrontBottom"]; + rMinFront_ = vArgs["RMinFront"]; + slopeT_ = vArgs["SlopeTop"]; + zFrontT_ = vArgs["ZFrontTop"]; + rMaxFront_ = vArgs["RMaxFront"]; +#ifdef EDM_ML_DEBUG + for (unsigned int i = 0; i < slopeB_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontB_[i] << " Rmin " << rMinFront_[i] + << " Slope " << slopeB_[i]; + for (unsigned int i = 0; i < slopeT_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontT_[i] << " Rmax " << rMaxFront_[i] + << " Slope " << slopeT_[i]; +#endif + nameSpace_ = DDCurrentNamespace::ns(); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: NameSpace " << nameSpace_; +#endif +} + +//////////////////////////////////////////////////////////////////// +// DDHGCalEEFileAlgo methods... +//////////////////////////////////////////////////////////////////// + +void DDHGCalEEFileAlgo::execute(DDCompactView& cpv) { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "==>> Constructing DDHGCalEEFileAlgo..."; + copies_.clear(); +#endif + constructLayers(parent(), cpv); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << copies_.size() << " different wafer copy numbers"; + int k(0); + for (std::unordered_set::const_iterator itr = copies_.begin(); itr != copies_.end(); ++itr, ++k) { + edm::LogVerbatim("HGCalGeom") << "Copy [" << k << "] : " << (*itr); + } + copies_.clear(); + edm::LogVerbatim("HGCalGeom") << "<<== End of DDHGCalEEFileAlgo construction..."; +#endif +} + +void DDHGCalEEFileAlgo::constructLayers(const DDLogicalPart& module, DDCompactView& cpv) { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: \t\tInside Layers"; +#endif + double zi(zMinBlock_); + int laymin(0); + for (unsigned int i = 0; i < layers_.size(); i++) { + double zo = zi + layerThick_[i]; + double routF = HGCalGeomTools::radius(zi, zFrontT_, rMaxFront_, slopeT_); + int laymax = laymin + layers_[i]; + double zz = zi; + double thickTot(0); + for (int ly = laymin; ly < laymax; ++ly) { + int ii = layerType_[ly]; + int copy = copyNumber_[ii]; + double hthick = 0.5 * thick_[ii]; + double rinB = HGCalGeomTools::radius(zo, zFrontB_, rMinFront_, slopeB_); + zz += hthick; + thickTot += thick_[ii]; + + std::string name = names_[ii] + std::to_string(copy); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: Layer " << ly << ":" << ii << " Front " << zi << ", " + << routF << " Back " << zo << ", " << rinB << " superlayer thickness " + << layerThick_[i]; +#endif + DDName matName(DDSplit(materials_[ii]).first, DDSplit(materials_[ii]).second); + DDMaterial matter(matName); + DDLogicalPart glog; + if (layerSense_[ly] < 1) { + std::vector pgonZ, pgonRin, pgonRout; + if (layerSense_[ly] == 0 || absorbMode_ == 0) { + double rmax = routF * cosAlpha_ - tol1_; + pgonZ.emplace_back(-hthick); + pgonZ.emplace_back(hthick); + pgonRin.emplace_back(rinB); + pgonRin.emplace_back(rinB); + pgonRout.emplace_back(rmax); + pgonRout.emplace_back(rmax); + } else { + HGCalGeomTools::radius(zz - hthick, + zz + hthick, + zFrontB_, + rMinFront_, + slopeB_, + zFrontT_, + rMaxFront_, + slopeT_, + -layerSense_[ly], + pgonZ, + pgonRin, + pgonRout); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: z " << (zz - hthick) << ":" << (zz + hthick) << " with " + << pgonZ.size() << " palnes"; + for (unsigned int isec = 0; isec < pgonZ.size(); ++isec) + edm::LogVerbatim("HGCalGeom") + << "[" << isec << "] z " << pgonZ[isec] << " R " << pgonRin[isec] << ":" << pgonRout[isec]; +#endif + for (unsigned int isec = 0; isec < pgonZ.size(); ++isec) { + pgonZ[isec] -= zz; + pgonRout[isec] = pgonRout[isec] * cosAlpha_ - tol1_; + } + } + DDSolid solid = + DDSolidFactory::polyhedra(DDName(name, nameSpace_), sectors_, -alpha_, 2._pi, pgonZ, pgonRin, pgonRout); + glog = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << solid.name() << " polyhedra of " << sectors_ + << " sectors covering " << convertRadToDeg(-alpha_) << ":" + << convertRadToDeg(-alpha_ + 2._pi) << " with " << pgonZ.size() + << " sections and filled with " << matName << ":" << &matter; + for (unsigned int k = 0; k < pgonZ.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] z " << pgonZ[k] << " R " << pgonRin[k] << ":" << pgonRout[k]; +#endif + } else { + DDSolid solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rinB, routF, 0.0, 2._pi); + glog = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << solid.name() << " Tubs made of " << matName << ":" + << &matter << " of dimensions " << rinB << ", " << routF << ", " << hthick + << ", 0.0, 360.0 and position " << glog.name() << " number " << copy << ":" + << layerCenter_[copy - 1]; +#endif + positionSensitive(glog, rinB, routF, zz, layerSense_[ly], (copy - 1), cpv); + } + DDTranslation r1(0, 0, zz); + DDRotation rot; + cpv.position(glog, module, copy, r1, rot); + ++copyNumber_[ii]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << glog.name() << " number " << copy << " positioned in " + << module.name() << " at " << r1 << " with " << rot; +#endif + zz += hthick; + } // End of loop over layers in a block + zi = zo; + laymin = laymax; + // Make consistency check of all the partitions of the block + if (std::abs(thickTot - layerThick_[i]) >= tol2_) { + if (thickTot > layerThick_[i]) { + edm::LogError("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " is smaller than " << thickTot + << ": thickness of all its components **** ERROR ****"; + } else { + edm::LogWarning("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " does not match with " + << thickTot << " of the components"; + } + } + } // End of loop over blocks +} + +void DDHGCalEEFileAlgo::positionSensitive( + const DDLogicalPart& glog, double rin, double rout, double zpos, int layertype, int layer, DDCompactView& cpv) { + static const double sqrt3 = std::sqrt(3.0); + int layercenter = layerCenter_[layer]; + double r = 0.5 * (waferSize_ + waferSepar_); + double R = 2.0 * r / sqrt3; + double dy = 0.75 * R; + int N = (int)(0.5 * rout / r) + 2; + const auto& xyoff = geomTools_.shiftXY(layercenter, (waferSize_ + waferSepar_)); +#ifdef EDM_ML_DEBUG + int ium(0), ivm(0), iumAll(0), ivmAll(0), kount(0), ntot(0), nin(0); + std::vector ntype(6, 0); + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << glog.ddname() << " rin:rout " << rin << ":" << rout + << " zpos " << zpos << " N " << N << " for maximum u, v; r " << r << " R " << R + << " dy " << dy << " Shift " << xyoff.first << ":" << xyoff.second << " WaferSize " + << (waferSize_ + waferSepar_); +#endif + for (int u = -N; u <= N; ++u) { + int iu = std::abs(u); + for (int v = -N; v <= N; ++v) { + int iv = std::abs(v); + int nr = 2 * v; + int nc = -2 * u + v; + double xpos = xyoff.first + nc * r; + double ypos = xyoff.second + nr * dy; + const auto& corner = HGCalGeomTools::waferCorner(xpos, ypos, r, R, rin, rout, false); +#ifdef EDM_ML_DEBUG + ++ntot; + if (((corner.first <= 0) && std::abs(u) < 5 && std::abs(v) < 5) || (std::abs(u) < 2 && std::abs(v) < 2)) { + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << glog.ddname() << " R " << rin << ":" << rout + << "\n Z " << zpos << " LayerType " << layertype << " u " << u << " v " << v + << " with " << corner.first << " corners"; + } +#endif + int type = HGCalWaferType::getType(HGCalWaferIndex::waferIndex(layer, u, v, false), waferIndex_, waferTypes_); + if (corner.first > 0 && type >= 0) { + int copy = type * 1000000 + iv * 100 + iu; + if (u < 0) + copy += 10000; + if (v < 0) + copy += 100000; +#ifdef EDM_ML_DEBUG + if (iu > ium) + ium = iu; + if (iv > ivm) + ivm = iv; + kount++; + if (copies_.count(copy) == 0) + copies_.insert(copy); +#endif + if (corner.first == (int)(HGCalParameters::k_CornerSize)) { +#ifdef EDM_ML_DEBUG + if (iu > iumAll) + iumAll = iu; + if (iv > ivmAll) + ivmAll = iv; + ++nin; +#endif + DDTranslation tran(xpos, ypos, 0.0); + DDRotation rotation; + if (layertype > 1) + type += 3; + DDName name = DDName(DDSplit(wafers_[type]).first, DDSplit(wafers_[type]).second); + cpv.position(name, glog.ddname(), copy, tran, rotation); +#ifdef EDM_ML_DEBUG + ++ntype[type]; + edm::LogVerbatim("HGCalGeom") << " DDHGCalEEFileAlgo: " << name << " number " << copy << " type " << layertype + << ":" << type << " positioned in " << glog.ddname() << " at " << tran + << " with " << rotation; +#endif + } + } + } + } +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: Maximum # of u " << ium << ":" << iumAll << " # of v " << ivm + << ":" << ivmAll << " and " << nin << ":" << kount << ":" << ntot << " wafers (" + << ntype[0] << ":" << ntype[1] << ":" << ntype[2] << ":" << ntype[3] << ":" << ntype[4] + << ":" << ntype[5] << ") for " << glog.ddname() << " R " << rin << ":" << rout; +#endif +} + +DEFINE_EDM_PLUGIN(DDAlgorithmFactory, DDHGCalEEFileAlgo, "hgcal:DDHGCalEEFileAlgo"); diff --git a/Geometry/HGCalCommonData/plugins/DDHGCalHEFileAlgo.cc b/Geometry/HGCalCommonData/plugins/DDHGCalHEFileAlgo.cc new file mode 100644 index 0000000000000..b4eb97ec50536 --- /dev/null +++ b/Geometry/HGCalCommonData/plugins/DDHGCalHEFileAlgo.cc @@ -0,0 +1,606 @@ +/////////////////////////////////////////////////////////////////////////////// +// File: DDHGCalHEFileAlgo.cc +// Description: Geometry factory class for HGCal (Mix) +/////////////////////////////////////////////////////////////////////////////// + +#include "DataFormats/Math/interface/CMSUnits.h" +#include "DetectorDescription/Core/interface/DDAlgorithm.h" +#include "DetectorDescription/Core/interface/DDAlgorithmFactory.h" +#include "DetectorDescription/Core/interface/DDCurrentNamespace.h" +#include "DetectorDescription/Core/interface/DDLogicalPart.h" +#include "DetectorDescription/Core/interface/DDMaterial.h" +#include "DetectorDescription/Core/interface/DDSolid.h" +#include "DetectorDescription/Core/interface/DDSplit.h" +#include "DetectorDescription/Core/interface/DDTypes.h" +#include "DetectorDescription/Core/interface/DDutils.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/PluginManager/interface/PluginFactory.h" +#include "Geometry/HGCalCommonData/interface/HGCalGeomTools.h" +#include "Geometry/HGCalCommonData/interface/HGCalParameters.h" +#include "Geometry/HGCalCommonData/interface/HGCalWaferIndex.h" +#include "Geometry/HGCalCommonData/interface/HGCalWaferType.h" + +#include +#include +#include +#include +#include + +//#define EDM_ML_DEBUG +using namespace cms_units::operators; + +class DDHGCalHEFileAlgo : public DDAlgorithm { +public: + DDHGCalHEFileAlgo(); + + void initialize(const DDNumericArguments& nArgs, + const DDVectorArguments& vArgs, + const DDMapArguments& mArgs, + const DDStringArguments& sArgs, + const DDStringVectorArguments& vsArgs) override; + void execute(DDCompactView& cpv) override; + +protected: + void constructLayers(const DDLogicalPart&, DDCompactView& cpv); + void positionMix(const DDLogicalPart& glog, + const std::string& name, + int copy, + double thick, + const DDMaterial& matter, + double rin, + double rmid, + double routF, + double zz, + DDCompactView& cpv); + void positionSensitive( + const DDLogicalPart& glog, double rin, double rout, double zpos, int layertype, int layer, DDCompactView& cpv); + +private: + HGCalGeomTools geomTools_; + + static constexpr double tol1_ = 0.01; + static constexpr double tol2_ = 0.00001; + + std::vector wafers_; // Wafers + std::vector materials_; // Materials + std::vector names_; // Names + std::vector thick_; // Thickness of the material + std::vector copyNumber_; // Initial copy numbers + std::vector layers_; // Number of layers in a section + std::vector layerThick_; // Thickness of each section + std::vector rMixLayer_; // Partition between Si/Sci part + std::vector layerType_; // Type of the layer + std::vector layerSense_; // Content of a layer (sensitive?) + int firstLayer_; // Copy # of the first sensitive layer + int absorbMode_; // Absorber mode + std::vector materialsTop_; // Materials of top layers + std::vector namesTop_; // Names of top layers + std::vector layerThickTop_; // Thickness of the top sections + std::vector layerTypeTop_; // Type of the Top layer + std::vector copyNumberTop_; // Initial copy numbers (top section) + std::vector materialsBot_; // Materials of bottom layers + std::vector namesBot_; // Names of bottom layers + std::vector layerThickBot_; // Thickness of the bottom sections + std::vector layerTypeBot_; // Type of the bottom layers + std::vector copyNumberBot_; // Initial copy numbers (bot section) + std::vector layerSenseBot_; // Content of bottom layer (sensitive?) + std::vector layerCenter_; // Centering of the wafers + + double zMinBlock_; // Starting z-value of the block + std::vector rad100to200_; // Parameters for 120-200mum trans. + std::vector rad200to300_; // Parameters for 200-300mum trans. + double zMinRadPar_; // Minimum z for radius parametriz. + std::vector waferIndex_; // Wafer index for the types + std::vector waferTypes_; // Wafer types + int choiceType_; // Type of parametrization to be used + int nCutRadPar_; // Cut off threshold for corners + double fracAreaMin_; // Minimum fractional conatined area + double waferSize_; // Width of the wafer + double waferSepar_; // Sensor separation + int sectors_; // Sectors + std::vector slopeB_; // Slope at the lower R + std::vector zFrontB_; // Starting Z values for the slopes + std::vector rMinFront_; // Corresponding rMin's + std::vector slopeT_; // Slopes at the larger R + std::vector zFrontT_; // Starting Z values for the slopes + std::vector rMaxFront_; // Corresponding rMax's + std::string nameSpace_; // Namespace of this and ALL sub-parts + std::unordered_set copies_; // List of copy #'s + double alpha_, cosAlpha_; +}; + +DDHGCalHEFileAlgo::DDHGCalHEFileAlgo() { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Creating an instance"; +#endif +} + +void DDHGCalHEFileAlgo::initialize(const DDNumericArguments& nArgs, + const DDVectorArguments& vArgs, + const DDMapArguments&, + const DDStringArguments& sArgs, + const DDStringVectorArguments& vsArgs) { + wafers_ = vsArgs["WaferNames"]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << wafers_.size() << " wafers"; + for (unsigned int i = 0; i < wafers_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << "] " << wafers_[i]; +#endif + materials_ = vsArgs["MaterialNames"]; + names_ = vsArgs["VolumeNames"]; + thick_ = vArgs["Thickness"]; + copyNumber_.resize(materials_.size(), 1); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << materials_.size() << " types of volumes"; + for (unsigned int i = 0; i < names_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << names_[i] << " of thickness " << thick_[i] + << " filled with " << materials_[i] << " first copy number " << copyNumber_[i]; +#endif + layers_ = dbl_to_int(vArgs["Layers"]); + layerThick_ = vArgs["LayerThick"]; + rMixLayer_ = vArgs["LayerRmix"]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "There are " << layers_.size() << " blocks"; + for (unsigned int i = 0; i < layers_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] of thickness " << layerThick_[i] << " Rmid " << rMixLayer_[i] + << " with " << layers_[i] << " layers"; +#endif + layerType_ = dbl_to_int(vArgs["LayerType"]); + layerSense_ = dbl_to_int(vArgs["LayerSense"]); + firstLayer_ = (int)(nArgs["FirstLayer"]); + absorbMode_ = (int)(nArgs["AbsorberMode"]); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "First Layer " << firstLayer_ << " and " + << "Absober mode " << absorbMode_; +#endif + layerCenter_ = dbl_to_int(vArgs["LayerCenter"]); +#ifdef EDM_ML_DEBUG + for (unsigned int i = 0; i < layerCenter_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "LayerCenter [" << i << "] " << layerCenter_[i]; +#endif + if (firstLayer_ > 0) { + for (unsigned int i = 0; i < layerType_.size(); ++i) { + if (layerSense_[i] > 0) { + int ii = layerType_[i]; + copyNumber_[ii] = firstLayer_; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "First copy number for layer type " << i << ":" << ii << " with " + << materials_[ii] << " changed to " << copyNumber_[ii]; +#endif + break; + } + } + } +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "There are " << layerType_.size() << " layers"; + for (unsigned int i = 0; i < layerType_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerType_[i] << " sensitive class " + << layerSense_[i]; +#endif + materialsTop_ = vsArgs["TopMaterialNames"]; + namesTop_ = vsArgs["TopVolumeNames"]; + layerThickTop_ = vArgs["TopLayerThickness"]; + layerTypeTop_ = dbl_to_int(vArgs["TopLayerType"]); + copyNumberTop_.resize(materialsTop_.size(), 1); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << materialsTop_.size() << " types of volumes in the top part"; + for (unsigned int i = 0; i < materialsTop_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << namesTop_[i] << " of thickness " << layerThickTop_[i] + << " filled with " << materialsTop_[i] << " first copy number " << copyNumberTop_[i]; + edm::LogVerbatim("HGCalGeom") << "There are " << layerTypeTop_.size() << " layers in the top part"; + for (unsigned int i = 0; i < layerTypeTop_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerTypeTop_[i]; +#endif + materialsBot_ = vsArgs["BottomMaterialNames"]; + namesBot_ = vsArgs["BottomVolumeNames"]; + layerTypeBot_ = dbl_to_int(vArgs["BottomLayerType"]); + layerSenseBot_ = dbl_to_int(vArgs["BottomLayerSense"]); + layerThickBot_ = vArgs["BottomLayerThickness"]; + copyNumberBot_.resize(materialsBot_.size(), 1); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << materialsBot_.size() + << " types of volumes in the bottom part"; + for (unsigned int i = 0; i < materialsBot_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << namesBot_[i] << " of thickness " << layerThickBot_[i] + << " filled with " << materialsBot_[i] << " first copy number " << copyNumberBot_[i]; + edm::LogVerbatim("HGCalGeom") << "There are " << layerTypeBot_.size() << " layers in the bottom part"; + for (unsigned int i = 0; i < layerTypeBot_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerTypeBot_[i] + << " sensitive class " << layerSenseBot_[i]; +#endif + zMinBlock_ = nArgs["zMinBlock"]; + rad100to200_ = vArgs["rad100to200"]; + rad200to300_ = vArgs["rad200to300"]; + zMinRadPar_ = nArgs["zMinForRadPar"]; + choiceType_ = (int)(nArgs["choiceType"]); + nCutRadPar_ = (int)(nArgs["nCornerCut"]); + fracAreaMin_ = nArgs["fracAreaMin"]; + waferSize_ = nArgs["waferSize"]; + waferSepar_ = nArgs["SensorSeparation"]; + sectors_ = (int)(nArgs["Sectors"]); + alpha_ = (1._pi) / sectors_; + cosAlpha_ = cos(alpha_); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: zStart " << zMinBlock_ + << " radius for wafer type separation uses " << rad100to200_.size() + << " parameters; zmin " << zMinRadPar_ << " cutoff " << choiceType_ << ":" + << nCutRadPar_ << ":" << fracAreaMin_ << " wafer width " << waferSize_ + << " separations " << waferSepar_ << " sectors " << sectors_ << ":" + << convertRadToDeg(alpha_) << ":" << cosAlpha_; + for (unsigned int k = 0; k < rad100to200_.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] 100-200 " << rad100to200_[k] << " 200-300 " << rad200to300_[k]; +#endif + waferIndex_ = dbl_to_int(vArgs["WaferIndex"]); + waferTypes_ = dbl_to_int(vArgs["WaferTypes"]); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "waferTypes with " << waferTypes_.size() << " entries"; + for (unsigned int k = 0; k < waferTypes_.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << waferIndex_[k] << " (" + << HGCalWaferIndex::waferLayer(waferIndex_[k]) << ", " + << HGCalWaferIndex::waferU(waferIndex_[k]) << ", " + << HGCalWaferIndex::waferV(waferIndex_[k]) << ") : " << waferTypes_[k]; +#endif + slopeB_ = vArgs["SlopeBottom"]; + zFrontB_ = vArgs["ZFrontBottom"]; + rMinFront_ = vArgs["RMinFront"]; + slopeT_ = vArgs["SlopeTop"]; + zFrontT_ = vArgs["ZFrontTop"]; + rMaxFront_ = vArgs["RMaxFront"]; +#ifdef EDM_ML_DEBUG + for (unsigned int i = 0; i < slopeB_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontB_[i] << " Rmin " << rMinFront_[i] + << " Slope " << slopeB_[i]; + for (unsigned int i = 0; i < slopeT_.size(); ++i) + edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontT_[i] << " Rmax " << rMaxFront_[i] + << " Slope " << slopeT_[i]; +#endif + nameSpace_ = DDCurrentNamespace::ns(); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: NameSpace " << nameSpace_; +#endif +} + +//////////////////////////////////////////////////////////////////// +// DDHGCalHEFileAlgo methods... +//////////////////////////////////////////////////////////////////// + +void DDHGCalHEFileAlgo::execute(DDCompactView& cpv) { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "==>> Constructing DDHGCalHEFileAlgo..."; + copies_.clear(); +#endif + constructLayers(parent(), cpv); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << copies_.size() << " different wafer copy numbers"; + int k(0); + for (std::unordered_set::const_iterator itr = copies_.begin(); itr != copies_.end(); ++itr, ++k) { + edm::LogVerbatim("HGCalGeom") << "Copy [" << k << "] : " << (*itr); + } + copies_.clear(); + edm::LogVerbatim("HGCalGeom") << "<<== End of DDHGCalHEFileAlgo construction..."; +#endif +} + +void DDHGCalHEFileAlgo::constructLayers(const DDLogicalPart& module, DDCompactView& cpv) { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: \t\tInside Layers"; +#endif + double zi(zMinBlock_); + int laymin(0); + for (unsigned int i = 0; i < layers_.size(); i++) { + double zo = zi + layerThick_[i]; + double routF = HGCalGeomTools::radius(zi, zFrontT_, rMaxFront_, slopeT_); + int laymax = laymin + layers_[i]; + double zz = zi; + double thickTot(0); + for (int ly = laymin; ly < laymax; ++ly) { + int ii = layerType_[ly]; + int copy = copyNumber_[ii]; + double hthick = 0.5 * thick_[ii]; + double rinB = HGCalGeomTools::radius(zo, zFrontB_, rMinFront_, slopeB_); + zz += hthick; + thickTot += thick_[ii]; + + std::string name = names_[ii] + std::to_string(copy); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Layer " << ly << ":" << ii << " Front " << zi << ", " + << routF << " Back " << zo << ", " << rinB << " superlayer thickness " + << layerThick_[i]; +#endif + DDName matName(DDSplit(materials_[ii]).first, DDSplit(materials_[ii]).second); + DDMaterial matter(matName); + DDLogicalPart glog; + if (layerSense_[ly] < 1) { + std::vector pgonZ, pgonRin, pgonRout; + if (layerSense_[ly] == 0 || absorbMode_ == 0) { + double rmax = + (std::min(routF, HGCalGeomTools::radius(zz + hthick, zFrontT_, rMaxFront_, slopeT_)) * cosAlpha_) - tol1_; + pgonZ.emplace_back(-hthick); + pgonZ.emplace_back(hthick); + pgonRin.emplace_back(rinB); + pgonRin.emplace_back(rinB); + pgonRout.emplace_back(rmax); + pgonRout.emplace_back(rmax); + } else { + HGCalGeomTools::radius(zz - hthick, + zz + hthick, + zFrontB_, + rMinFront_, + slopeB_, + zFrontT_, + rMaxFront_, + slopeT_, + -layerSense_[ly], + pgonZ, + pgonRin, + pgonRout); + for (unsigned int isec = 0; isec < pgonZ.size(); ++isec) { + pgonZ[isec] -= zz; + pgonRout[isec] = pgonRout[isec] * cosAlpha_ - tol1_; + } + } + DDSolid solid = + DDSolidFactory::polyhedra(DDName(name, nameSpace_), sectors_, -alpha_, 2._pi, pgonZ, pgonRin, pgonRout); + glog = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " polyhedra of " << sectors_ + << " sectors covering " << convertRadToDeg(-alpha_) << ":" + << convertRadToDeg(-alpha_ + 2._pi) << " with " << pgonZ.size() << " sections"; + for (unsigned int k = 0; k < pgonZ.size(); ++k) + edm::LogVerbatim("HGCalGeom") << "[" << k << "] z " << pgonZ[k] << " R " << pgonRin[k] << ":" << pgonRout[k]; +#endif + } else { + DDSolid solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rinB, routF, 0.0, 2._pi); + glog = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " Tubs made of " << matName + << " of dimensions " << rinB << ", " << routF << ", " << hthick + << ", 0.0, 360.0 and positioned in: " << glog.name() << " number " << copy; +#endif + positionMix(glog, name, copy, thick_[ii], matter, rinB, rMixLayer_[i], routF, zz, cpv); + } + DDTranslation r1(0, 0, zz); + DDRotation rot; + cpv.position(glog, module, copy, r1, rot); + ++copyNumber_[ii]; +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << glog.name() << " number " << copy << " positioned in " + << module.name() << " at " << r1 << " with " << rot; +#endif + zz += hthick; + } // End of loop over layers in a block + zi = zo; + laymin = laymax; + // Make consistency check of all the partitions of the block + if (std::abs(thickTot - layerThick_[i]) >= tol2_) { + if (thickTot > layerThick_[i]) { + edm::LogError("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " is smaller than " << thickTot + << ": thickness of all its components **** ERROR ****"; + } else { + edm::LogWarning("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " does not match with " + << thickTot << " of the components"; + } + } + } // End of loop over blocks +} + +void DDHGCalHEFileAlgo::positionMix(const DDLogicalPart& glog, + const std::string& nameM, + int copyM, + double thick, + const DDMaterial& matter, + double rin, + double rmid, + double rout, + double zz, + DDCompactView& cpv) { + DDLogicalPart glog1; + DDTranslation tran; + DDRotation rot; + for (unsigned int ly = 0; ly < layerTypeTop_.size(); ++ly) { + int ii = layerTypeTop_[ly]; + copyNumberTop_[ii] = copyM; + } + for (unsigned int ly = 0; ly < layerTypeBot_.size(); ++ly) { + int ii = layerTypeBot_[ly]; + copyNumberBot_[ii] = copyM; + } + double hthick = 0.5 * thick; + // Make the top part first + std::string name = nameM + "Top"; + DDSolid solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rmid, rout, 0.0, 2._pi); + glog1 = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " Tubs made of " << matter.name() + << " of dimensions " << rmid << ", " << rout << ", " << hthick << ", 0.0, 360.0"; +#endif + cpv.position(glog1, glog, 1, tran, rot); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << glog1.name() << " number 1 positioned in " << glog.name() + << " at " << tran << " with " << rot; +#endif + double thickTot(0), zpos(-hthick); + for (unsigned int ly = 0; ly < layerTypeTop_.size(); ++ly) { + int ii = layerTypeTop_[ly]; + int copy = copyNumberTop_[ii]; + double hthickl = 0.5 * layerThickTop_[ii]; + thickTot += layerThickTop_[ii]; + name = namesTop_[ii] + std::to_string(copy); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Layer " << ly << ":" << ii << " R " << rmid << ":" << rout + << " Thick " << layerThickTop_[ii]; +#endif + DDName matName(DDSplit(materialsTop_[ii]).first, DDSplit(materialsTop_[ii]).second); + DDMaterial matter1(matName); + solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthickl, rmid, rout, 0.0, 2._pi); + DDLogicalPart glog2 = DDLogicalPart(solid.ddname(), matter1, solid); +#ifdef EDM_ML_DEBUG + double eta1 = -log(tan(0.5 * atan(rmid / zz))); + double eta2 = -log(tan(0.5 * atan(rout / zz))); + edm::LogVerbatim("HGCalGeom") << name << " z|rin|rout " << zz << ":" << rmid << ":" << rout << " eta " << eta1 + << ":" << eta2; + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " Tubs made of " << matName + << " of dimensions " << rmid << ", " << rout << ", " << hthickl << ", 0.0, 360.0"; +#endif + zpos += hthickl; + DDTranslation r1(0, 0, zpos); + cpv.position(glog2, glog1, copy, r1, rot); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Position " << glog2.name() << " number " << copy << " in " + << glog1.name() << " at " << r1 << " with " << rot; +#endif + ++copyNumberTop_[ii]; + zpos += hthickl; + } + if (std::abs(thickTot - thick) >= tol2_) { + if (thickTot > thick) { + edm::LogError("HGCalGeom") << "Thickness of the partition " << thick << " is smaller than " << thickTot + << ": thickness of all its components in the top part **** ERROR ****"; + } else { + edm::LogWarning("HGCalGeom") << "Thickness of the partition " << thick << " does not match with " << thickTot + << " of the components in top part"; + } + } + + // Make the bottom part next + name = nameM + "Bottom"; + solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rin, rmid, 0.0, 2._pi); + glog1 = DDLogicalPart(solid.ddname(), matter, solid); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " Tubs made of " << matter.name() + << " of dimensions " << rin << ", " << rmid << ", " << hthick << ", 0.0, 360.0"; +#endif + cpv.position(glog1, glog, 1, tran, rot); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << glog1.name() << " number 1 positioned in " << glog.name() + << " at " << tran << " with " << rot; +#endif + thickTot = 0; + zpos = -hthick; + for (unsigned int ly = 0; ly < layerTypeBot_.size(); ++ly) { + int ii = layerTypeBot_[ly]; + int copy = copyNumberBot_[ii]; + double hthickl = 0.5 * layerThickBot_[ii]; + thickTot += layerThickBot_[ii]; + name = namesBot_[ii] + std::to_string(copy); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Layer " << ly << ":" << ii << " R " << rin << ":" << rmid + << " Thick " << layerThickBot_[ii]; +#endif + DDName matName(DDSplit(materialsBot_[ii]).first, DDSplit(materialsBot_[ii]).second); + DDMaterial matter1(matName); + solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthickl, rin, rmid, 0.0, 2._pi); + DDLogicalPart glog2 = DDLogicalPart(solid.ddname(), matter1, solid); +#ifdef EDM_ML_DEBUG + double eta1 = -log(tan(0.5 * atan(rin / zz))); + double eta2 = -log(tan(0.5 * atan(rmid / zz))); + edm::LogVerbatim("HGCalGeom") << name << " z|rin|rout " << zz << ":" << rin << ":" << rmid << " eta " << eta1 << ":" + << eta2; + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << solid.name() << " Tubs made of " << matName + << " of dimensions " << rin << ", " << rmid << ", " << hthickl << ", 0.0, 360.0"; +#endif + zpos += hthickl; + DDTranslation r1(0, 0, zpos); + cpv.position(glog2, glog1, copy, r1, rot); +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Position " << glog2.name() << " number " << copy << " in " + << glog1.name() << " at " << r1 << " with " << rot; +#endif + if (layerSenseBot_[ly] != 0) { +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: z " << (zz + zpos) << " Center " << copy << ":" + << (copy - firstLayer_) << ":" << layerCenter_[copy - firstLayer_]; +#endif + positionSensitive(glog2, rin, rmid, zz + zpos, layerSenseBot_[ly], (copy - firstLayer_), cpv); + } + zpos += hthickl; + ++copyNumberBot_[ii]; + } + if (std::abs(thickTot - thick) >= tol2_) { + if (thickTot > thick) { + edm::LogError("HGCalGeom") << "Thickness of the partition " << thick << " is smaller than " << thickTot + << ": thickness of all its components in the top part **** ERROR ****"; + } else { + edm::LogWarning("HGCalGeom") << "Thickness of the partition " << thick << " does not match with " << thickTot + << " of the components in top part"; + } + } +} + +void DDHGCalHEFileAlgo::positionSensitive( + const DDLogicalPart& glog, double rin, double rout, double zpos, int layertype, int layer, DDCompactView& cpv) { + static const double sqrt3 = std::sqrt(3.0); + int layercenter = layerCenter_[layer]; + double r = 0.5 * (waferSize_ + waferSepar_); + double R = 2.0 * r / sqrt3; + double dy = 0.75 * R; + int N = (int)(0.5 * rout / r) + 2; + const auto& xyoff = geomTools_.shiftXY(layercenter, (waferSize_ + waferSepar_)); +#ifdef EDM_ML_DEBUG + int ium(0), ivm(0), iumAll(0), ivmAll(0), kount(0), ntot(0), nin(0); + std::vector ntype(6, 0); + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: " << glog.ddname() << " rin:rout " << rin << ":" << rout + << " zpos " << zpos << " N " << N << " for maximum u, v Offset; Shift " << xyoff.first + << ":" << xyoff.second << " WaferSize " << (waferSize_ + waferSepar_); +#endif + for (int u = -N; u <= N; ++u) { + int iu = std::abs(u); + for (int v = -N; v <= N; ++v) { + int iv = std::abs(v); + int nr = 2 * v; + int nc = -2 * u + v; + double xpos = xyoff.first + nc * r; + double ypos = xyoff.second + nr * dy; + const auto& corner = HGCalGeomTools::waferCorner(xpos, ypos, r, R, rin, rout, false); +#ifdef EDM_ML_DEBUG + ++ntot; +#endif + int type = HGCalWaferType::getType(HGCalWaferIndex::waferIndex(layer, u, v, false), waferIndex_, waferTypes_); + if (corner.first > 0 && type >= 0) { + int copy = type * 1000000 + iv * 100 + iu; + if (u < 0) + copy += 10000; + if (v < 0) + copy += 100000; +#ifdef EDM_ML_DEBUG + if (iu > ium) + ium = iu; + if (iv > ivm) + ivm = iv; + kount++; + if (copies_.count(copy) == 0) + copies_.insert(copy); +#endif + if (corner.first == (int)(HGCalParameters::k_CornerSize)) { +#ifdef EDM_ML_DEBUG + if (iu > iumAll) + iumAll = iu; + if (iv > ivmAll) + ivmAll = iv; + ++nin; +#endif + DDTranslation tran(xpos, ypos, 0.0); + DDRotation rotation; + if (layertype > 1) + type += 3; + DDName name = DDName(DDSplit(wafers_[type]).first, DDSplit(wafers_[type]).second); + cpv.position(name, glog.ddname(), copy, tran, rotation); +#ifdef EDM_ML_DEBUG + ++ntype[type]; + edm::LogVerbatim("HGCalGeom") << " DDHGCalHEFileAlgo: " << name << " number " << copy << " type " << layertype + << ":" << type << " positioned in " << glog.ddname() << " at " << tran + << " with " << rotation; +#endif + } + } + } + } +#ifdef EDM_ML_DEBUG + edm::LogVerbatim("HGCalGeom") << "DDHGCalHEFileAlgo: Maximum # of u " << ium << ":" << iumAll << " # of v " << ivm + << ":" << ivmAll << " and " << nin << ":" << kount << ":" << ntot << " wafers (" + << ntype[0] << ":" << ntype[1] << ":" << ntype[2] << ":" << ntype[3] << ":" << ntype[4] + << ":" << ntype[5] << ") for " << glog.ddname() << " R " << rin << ":" << rout; +#endif +} + +DEFINE_EDM_PLUGIN(DDAlgorithmFactory, DDHGCalHEFileAlgo, "hgcal:DDHGCalHEFileAlgo"); diff --git a/Geometry/HGCalCommonData/python/testHGCalV10XML_cfi.py b/Geometry/HGCalCommonData/python/testHGCalV10XML_cfi.py index 963397d6e7113..26287e36d6ed2 100644 --- a/Geometry/HGCalCommonData/python/testHGCalV10XML_cfi.py +++ b/Geometry/HGCalCommonData/python/testHGCalV10XML_cfi.py @@ -15,7 +15,7 @@ 'Geometry/CMSCommonData/data/cmsCalo.xml', 'Geometry/CMSCommonData/data/beampipe/2026/v1/beampipe.xml', 'Geometry/CMSCommonData/data/cmsBeam/2026/v1/cmsBeam.xml', - 'Geometry/CMSCommonData/data/cavernData/2017/v1/cavernData.xml', + 'Geometry/CMSCommonData/data/cavernData/2021/v1/cavernData.xml', 'Geometry/EcalCommonData/data/eregalgo/2026/v2/eregalgo.xml', 'Geometry/EcalCommonData/data/ectkcable/2026/v1/ectkcable.xml', 'Geometry/EcalCommonData/data/ectkcablemat/2026/v1/ectkcablemat.xml', diff --git a/Geometry/HGCalCommonData/python/testHGCalV11XML_cfi.py b/Geometry/HGCalCommonData/python/testHGCalV11XML_cfi.py index 8d966ef5f005e..1684f9bcda875 100644 --- a/Geometry/HGCalCommonData/python/testHGCalV11XML_cfi.py +++ b/Geometry/HGCalCommonData/python/testHGCalV11XML_cfi.py @@ -15,7 +15,7 @@ 'Geometry/CMSCommonData/data/cmsCalo.xml', 'Geometry/CMSCommonData/data/beampipe/2026/v1/beampipe.xml', 'Geometry/CMSCommonData/data/cmsBeam/2026/v1/cmsBeam.xml', - 'Geometry/CMSCommonData/data/cavernData/2017/v1/cavernData.xml', + 'Geometry/CMSCommonData/data/cavernData/2021/v1/cavernData.xml', 'Geometry/CMSCommonData/data/muonBase/2026/v2/muonBase.xml', 'Geometry/CMSCommonData/data/cmsMuon.xml', 'Geometry/CMSCommonData/data/mgnt.xml', diff --git a/Geometry/HGCalCommonData/python/testHGCalV12XML_cfi.py b/Geometry/HGCalCommonData/python/testHGCalV12XML_cfi.py index e2569dacae092..b3bfe1e321169 100644 --- a/Geometry/HGCalCommonData/python/testHGCalV12XML_cfi.py +++ b/Geometry/HGCalCommonData/python/testHGCalV12XML_cfi.py @@ -15,7 +15,7 @@ 'Geometry/CMSCommonData/data/cmsCalo.xml', 'Geometry/CMSCommonData/data/beampipe/2026/v1/beampipe.xml', 'Geometry/CMSCommonData/data/cmsBeam/2026/v1/cmsBeam.xml', - 'Geometry/CMSCommonData/data/cavernData/2017/v1/cavernData.xml', + 'Geometry/CMSCommonData/data/cavernData/2021/v1/cavernData.xml', 'Geometry/CMSCommonData/data/muonBase/2026/v4/muonBase.xml', 'Geometry/CMSCommonData/data/cmsMuon.xml', 'Geometry/CMSCommonData/data/mgnt.xml', diff --git a/Geometry/HGCalCommonData/python/testHGCalV14XML_cfi.py b/Geometry/HGCalCommonData/python/testHGCalV14XML_cfi.py new file mode 100644 index 0000000000000..3cf71c81486c8 --- /dev/null +++ b/Geometry/HGCalCommonData/python/testHGCalV14XML_cfi.py @@ -0,0 +1,89 @@ +import FWCore.ParameterSet.Config as cms + +# This config was generated automatically using generate2026Geometry.py +# If you notice a mistake, please update the generating script, not just this config + +XMLIdealGeometryESSource = cms.ESSource("XMLIdealGeometryESSource", + geomXMLFiles = cms.vstring( + 'Geometry/CMSCommonData/data/materials.xml', + 'Geometry/CMSCommonData/data/rotations.xml', + 'Geometry/CMSCommonData/data/extend/cmsextent.xml', + 'Geometry/CMSCommonData/data/cms/2026/v4/cms.xml', + 'Geometry/CMSCommonData/data/eta3/etaMax.xml', + 'Geometry/CMSCommonData/data/cmsMother.xml', + 'Geometry/CMSCommonData/data/caloBase/2026/v4/caloBase.xml', + 'Geometry/CMSCommonData/data/cmsCalo.xml', + 'Geometry/CMSCommonData/data/beampipe/2026/v1/beampipe.xml', + 'Geometry/CMSCommonData/data/cmsBeam/2026/v1/cmsBeam.xml', + 'Geometry/CMSCommonData/data/cavernData/2021/v1/cavernData.xml', + 'Geometry/CMSCommonData/data/muonBase/2026/v5/muonBase.xml', + 'Geometry/CMSCommonData/data/cmsMuon.xml', + 'Geometry/CMSCommonData/data/mgnt.xml', + 'Geometry/CMSCommonData/data/muonMB.xml', + 'Geometry/CMSCommonData/data/muonMagnet.xml', + 'Geometry/EcalCommonData/data/eregalgo/2026/v2/eregalgo.xml', + 'Geometry/EcalCommonData/data/ectkcable/2026/v1/ectkcable.xml', + 'Geometry/EcalCommonData/data/ectkcablemat/2026/v1/ectkcablemat.xml', + 'Geometry/EcalCommonData/data/ebalgo.xml', + 'Geometry/EcalCommonData/data/ebcon.xml', + 'Geometry/EcalCommonData/data/ebrot.xml', + 'Geometry/HcalCommonData/data/hcalrotations.xml', + 'Geometry/HcalCommonData/data/hcal/v2/hcalalgo.xml', + 'Geometry/HcalCommonData/data/hcalbarrelalgo.xml', + 'Geometry/HcalCommonData/data/hcalcablealgo/v2/hcalcablealgo.xml', + 'Geometry/HcalCommonData/data/hcalouteralgo.xml', + 'Geometry/HcalCommonData/data/hcalforwardalgo.xml', + 'Geometry/HcalCommonData/data/hcalSimNumbering/NoHE/hcalSimNumbering.xml', + 'Geometry/HcalCommonData/data/hcalRecNumbering/NoHE/hcalRecNumbering.xml', + 'Geometry/HcalCommonData/data/average/hcalforwardmaterial.xml', + 'Geometry/HGCalCommonData/data/hgcalMaterial/v1/hgcalMaterial.xml', + 'Geometry/HGCalCommonData/data/hgcal/v13/hgcal.xml', + 'Geometry/HGCalCommonData/data/hgcalcell/v9/hgcalcell.xml', + 'Geometry/HGCalCommonData/data/hgcalwafer/v9/hgcalwafer.xml', + 'Geometry/HGCalCommonData/data/hgcalEE/v14/hgcalEE.xml', + 'Geometry/HGCalCommonData/data/hgcalHEsil/v14/hgcalHEsil.xml', + 'Geometry/HGCalCommonData/data/hgcalHEmix/v14/hgcalHEmix.xml', + 'Geometry/HGCalCommonData/data/hgcalCons/v13/hgcalCons.xml', + 'Geometry/HGCalCommonData/data/hgcalConsData/v13/hgcalConsData.xml', + 'Geometry/ForwardCommonData/data/forwardshield/2017/v1/forwardshield.xml', + 'Geometry/MuonCommonData/data/mbCommon/2021/v1/mbCommon.xml', + 'Geometry/MuonCommonData/data/mb1/2015/v2/mb1.xml', + 'Geometry/MuonCommonData/data/mb2/2015/v2/mb2.xml', + 'Geometry/MuonCommonData/data/mb3/2015/v2/mb3.xml', + 'Geometry/MuonCommonData/data/mb4/2015/v2/mb4.xml', + 'Geometry/MuonCommonData/data/mb4Shield/2021/v1/mb4Shield.xml', + 'Geometry/MuonCommonData/data/muonYoke/2021/v3/muonYoke.xml', + 'Geometry/MuonCommonData/data/csc/2021/v2/csc.xml', + 'Geometry/MuonCommonData/data/mfshield/2017/v2/mfshield.xml', + 'Geometry/MuonCommonData/data/mf/2026/v6/mf.xml', + 'Geometry/MuonCommonData/data/rpcf/2026/v3/rpcf.xml', + 'Geometry/MuonCommonData/data/gemf/TDR_BaseLine/gemf.xml', + 'Geometry/MuonCommonData/data/gem11/TDR_BaseLine/gem11.xml', + 'Geometry/MuonCommonData/data/gem21/TDR_Dev/gem21.xml', + 'Geometry/MuonCommonData/data/mfshield/2026/v3/mfshield.xml', + 'Geometry/MuonCommonData/data/me0/TDR_Dev/v3/me0.xml', + )+ + cms.vstring( + 'Geometry/MuonCommonData/data/muonNumbering/TDR_DeV/v1/muonNumbering.xml', + 'Geometry/EcalSimData/data/PhaseII/ecalsens.xml', + 'Geometry/HcalCommonData/data/hcalsens/NoHE/hcalsenspmf.xml', + 'Geometry/HcalSimData/data/hf.xml', + 'Geometry/HcalSimData/data/hfpmt.xml', + 'Geometry/HcalSimData/data/hffibrebundle.xml', + 'Geometry/HcalSimData/data/CaloUtil.xml', + 'Geometry/HGCalSimData/data/hgcsensv9.xml', + 'Geometry/MuonSimData/data/PhaseII/ME0EtaPart/muonSens.xml', + 'Geometry/DTGeometryBuilder/data/dtSpecsFilter.xml', + 'Geometry/CSCGeometryBuilder/data/cscSpecsFilter.xml', + 'Geometry/CSCGeometryBuilder/data/cscSpecs.xml', + 'Geometry/RPCGeometryBuilder/data/2026/v1/RPCSpecs.xml', + 'Geometry/GEMGeometryBuilder/data/v7/GEMSpecsFilter.xml', + 'Geometry/GEMGeometryBuilder/data/v7/GEMSpecs.xml', + 'Geometry/HcalSimData/data/HcalProdCuts.xml', + 'Geometry/EcalSimData/data/EcalProdCuts.xml', + 'Geometry/HGCalSimData/data/hgcProdCutsv9.xml', + 'Geometry/CMSCommonData/data/FieldParameters.xml', + 'Geometry/MuonSimData/data/PhaseII/muonProdCuts.xml', + ), + rootNodeName = cms.string('cms:OCMS') +) diff --git a/Geometry/HGCalCommonData/src/HGCalGeometryMode.cc b/Geometry/HGCalCommonData/src/HGCalGeometryMode.cc index ab7df7cc0f9c6..bfd8eb9a8a22a 100644 --- a/Geometry/HGCalCommonData/src/HGCalGeometryMode.cc +++ b/Geometry/HGCalCommonData/src/HGCalGeometryMode.cc @@ -8,6 +8,7 @@ HGCalStringToEnumParser::HGCalStringToEnumParse enumMap["HGCalGeometryMode::Hexagon8"] = HGCalGeometryMode::Hexagon8; enumMap["HGCalGeometryMode::Hexagon8Full"] = HGCalGeometryMode::Hexagon8Full; enumMap["HGCalGeometryMode::Trapezoid"] = HGCalGeometryMode::Trapezoid; + enumMap["HGCalGeometryMode::HexagonFullPart"] = HGCalGeometryMode::HexagonFullPart; } template <> diff --git a/Geometry/HGCalCommonData/src/HGCalWaferType.cc b/Geometry/HGCalCommonData/src/HGCalWaferType.cc index 402b119db2a2f..c5e905170e8c1 100644 --- a/Geometry/HGCalCommonData/src/HGCalWaferType.cc +++ b/Geometry/HGCalCommonData/src/HGCalWaferType.cc @@ -103,6 +103,12 @@ int HGCalWaferType::getType(double xpos, double ypos, double zpos) { return type; } +int HGCalWaferType::getType(int index, const std::vector& indices, const std::vector& types) { + auto itr = static_cast(std::find(std::begin(indices), std::end(indices), index) - std::begin(indices)); + int type = (itr < indices.size()) ? types[itr] : -1; + return type; +} + std::pair HGCalWaferType::rLimits(double zpos) { double zz = std::abs(zpos); if (zz < zMin_)