Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:14:58

0001 ///////////////////////////////////////////////////////////////////////////////
0002 // File: DDHGCalHEAlgo.cc
0003 // Description: Geometry factory class for HGCal (Mix)
0004 ///////////////////////////////////////////////////////////////////////////////
0005 
0006 #include "DataFormats/Math/interface/angle_units.h"
0007 #include "DetectorDescription/Core/interface/DDAlgorithm.h"
0008 #include "DetectorDescription/Core/interface/DDAlgorithmFactory.h"
0009 #include "DetectorDescription/Core/interface/DDCurrentNamespace.h"
0010 #include "DetectorDescription/Core/interface/DDLogicalPart.h"
0011 #include "DetectorDescription/Core/interface/DDMaterial.h"
0012 #include "DetectorDescription/Core/interface/DDSolid.h"
0013 #include "DetectorDescription/Core/interface/DDSplit.h"
0014 #include "DetectorDescription/Core/interface/DDTypes.h"
0015 #include "DetectorDescription/Core/interface/DDutils.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/PluginManager/interface/PluginFactory.h"
0018 #include "Geometry/HGCalCommonData/interface/HGCalGeomTools.h"
0019 #include "Geometry/HGCalCommonData/interface/HGCalParameters.h"
0020 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0021 #include "Geometry/HGCalCommonData/interface/HGCalWaferType.h"
0022 
0023 #include <cmath>
0024 #include <memory>
0025 #include <string>
0026 #include <unordered_set>
0027 #include <vector>
0028 
0029 //#define EDM_ML_DEBUG
0030 using namespace angle_units::operators;
0031 
0032 class DDHGCalHEAlgo : public DDAlgorithm {
0033 public:
0034   // Constructor and Destructor
0035   DDHGCalHEAlgo();  // const std::string & name);
0036   ~DDHGCalHEAlgo() override;
0037 
0038   void initialize(const DDNumericArguments& nArgs,
0039                   const DDVectorArguments& vArgs,
0040                   const DDMapArguments& mArgs,
0041                   const DDStringArguments& sArgs,
0042                   const DDStringVectorArguments& vsArgs) override;
0043   void execute(DDCompactView& cpv) override;
0044 
0045 protected:
0046   void constructLayers(const DDLogicalPart&, DDCompactView& cpv);
0047   void positionMix(const DDLogicalPart& glog,
0048                    const std::string& name,
0049                    int copy,
0050                    double thick,
0051                    const DDMaterial& matter,
0052                    double rin,
0053                    double rmid,
0054                    double routF,
0055                    double zz,
0056                    DDCompactView& cpv);
0057   void positionSensitive(const DDLogicalPart& glog,
0058                          double rin,
0059                          double rout,
0060                          double zpos,
0061                          int layertype,
0062                          int layercenter,
0063                          DDCompactView& cpv);
0064 
0065 private:
0066   HGCalGeomTools geomTools_;
0067   std::unique_ptr<HGCalWaferType> waferType_;
0068 
0069   static constexpr double tol1_ = 0.01;
0070   static constexpr double tol2_ = 0.00001;
0071 
0072   std::vector<std::string> wafers_;        // Wafers
0073   std::vector<std::string> materials_;     // Materials
0074   std::vector<std::string> names_;         // Names
0075   std::vector<double> thick_;              // Thickness of the material
0076   std::vector<int> copyNumber_;            // Initial copy numbers
0077   std::vector<int> layers_;                // Number of layers in a section
0078   std::vector<double> layerThick_;         // Thickness of each section
0079   std::vector<double> rMixLayer_;          // Partition between Si/Sci part
0080   std::vector<int> layerType_;             // Type of the layer
0081   std::vector<int> layerSense_;            // Content of a layer (sensitive?)
0082   int firstLayer_;                         // Copy # of the first sensitive layer
0083   int absorbMode_;                         // Absorber mode
0084   int sensitiveMode_;                      // Sensitive mode
0085   std::vector<std::string> materialsTop_;  // Materials of top layers
0086   std::vector<std::string> namesTop_;      // Names of top layers
0087   std::vector<double> layerThickTop_;      // Thickness of the top sections
0088   std::vector<int> layerTypeTop_;          // Type of the Top layer
0089   std::vector<int> copyNumberTop_;         // Initial copy numbers (top section)
0090   std::vector<std::string> materialsBot_;  // Materials of bottom layers
0091   std::vector<std::string> namesBot_;      // Names of bottom layers
0092   std::vector<double> layerThickBot_;      // Thickness of the bottom sections
0093   std::vector<int> layerTypeBot_;          // Type of the bottom layers
0094   std::vector<int> copyNumberBot_;         // Initial copy numbers (bot section)
0095   std::vector<int> layerSenseBot_;         // Content of bottom layer (sensitive?)
0096   std::vector<int> layerCenter_;           // Centering of the wafers
0097 
0098   double zMinBlock_;                 // Starting z-value of the block
0099   std::vector<double> rad100to200_;  // Parameters for 120-200mum trans.
0100   std::vector<double> rad200to300_;  // Parameters for 200-300mum trans.
0101   double zMinRadPar_;                // Minimum z for radius parametriz.
0102   int choiceType_;                   // Type of parametrization to be used
0103   int nCutRadPar_;                   // Cut off threshold for corners
0104   double fracAreaMin_;               // Minimum fractional conatined area
0105   double waferSize_;                 // Width of the wafer
0106   double waferSepar_;                // Sensor separation
0107   int sectors_;                      // Sectors
0108   std::vector<double> slopeB_;       // Slope at the lower R
0109   std::vector<double> zFrontB_;      // Starting Z values for the slopes
0110   std::vector<double> rMinFront_;    // Corresponding rMin's
0111   std::vector<double> slopeT_;       // Slopes at the larger R
0112   std::vector<double> zFrontT_;      // Starting Z values for the slopes
0113   std::vector<double> rMaxFront_;    // Corresponding rMax's
0114   std::string nameSpace_;            // Namespace of this and ALL sub-parts
0115   std::unordered_set<int> copies_;   // List of copy #'s
0116   double alpha_, cosAlpha_;
0117 };
0118 
0119 DDHGCalHEAlgo::DDHGCalHEAlgo() {
0120 #ifdef EDM_ML_DEBUG
0121   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Creating an instance";
0122 #endif
0123 }
0124 
0125 DDHGCalHEAlgo::~DDHGCalHEAlgo() {}
0126 
0127 void DDHGCalHEAlgo::initialize(const DDNumericArguments& nArgs,
0128                                const DDVectorArguments& vArgs,
0129                                const DDMapArguments&,
0130                                const DDStringArguments& sArgs,
0131                                const DDStringVectorArguments& vsArgs) {
0132   wafers_ = vsArgs["WaferNames"];
0133 #ifdef EDM_ML_DEBUG
0134   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << wafers_.size() << " wafers";
0135   for (unsigned int i = 0; i < wafers_.size(); ++i)
0136     edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << "] " << wafers_[i];
0137 #endif
0138   materials_ = vsArgs["MaterialNames"];
0139   names_ = vsArgs["VolumeNames"];
0140   thick_ = vArgs["Thickness"];
0141   copyNumber_.resize(materials_.size(), 1);
0142 #ifdef EDM_ML_DEBUG
0143   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << materials_.size() << " types of volumes";
0144   for (unsigned int i = 0; i < names_.size(); ++i)
0145     edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << names_[i] << " of thickness " << thick_[i]
0146                                   << " filled with " << materials_[i] << " first copy number " << copyNumber_[i];
0147 #endif
0148   layers_ = dbl_to_int(vArgs["Layers"]);
0149   layerThick_ = vArgs["LayerThick"];
0150   rMixLayer_ = vArgs["LayerRmix"];
0151 #ifdef EDM_ML_DEBUG
0152   edm::LogVerbatim("HGCalGeom") << "There are " << layers_.size() << " blocks";
0153   for (unsigned int i = 0; i < layers_.size(); ++i)
0154     edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] of thickness " << layerThick_[i] << " Rmid " << rMixLayer_[i]
0155                                   << " with " << layers_[i] << " layers";
0156 #endif
0157   layerType_ = dbl_to_int(vArgs["LayerType"]);
0158   layerSense_ = dbl_to_int(vArgs["LayerSense"]);
0159   firstLayer_ = (int)(nArgs["FirstLayer"]);
0160   absorbMode_ = (int)(nArgs["AbsorberMode"]);
0161   sensitiveMode_ = (int)(nArgs["SensitiveMode"]);
0162 #ifdef EDM_ML_DEBUG
0163   edm::LogVerbatim("HGCalGeom") << "First Layer " << firstLayer_ << " and "
0164                                 << "Absober:Sensitive mode " << absorbMode_ << ":" << sensitiveMode_;
0165 #endif
0166   layerCenter_ = dbl_to_int(vArgs["LayerCenter"]);
0167 #ifdef EDM_ML_DEBUG
0168   for (unsigned int i = 0; i < layerCenter_.size(); ++i)
0169     edm::LogVerbatim("HGCalGeom") << "LayerCenter [" << i << "] " << layerCenter_[i];
0170 #endif
0171   if (firstLayer_ > 0) {
0172     for (unsigned int i = 0; i < layerType_.size(); ++i) {
0173       if (layerSense_[i] > 0) {
0174         int ii = layerType_[i];
0175         copyNumber_[ii] = firstLayer_;
0176 #ifdef EDM_ML_DEBUG
0177         edm::LogVerbatim("HGCalGeom") << "First copy number for layer type " << i << ":" << ii << " with "
0178                                       << materials_[ii] << " changed to " << copyNumber_[ii];
0179 #endif
0180         break;
0181       }
0182     }
0183   }
0184 #ifdef EDM_ML_DEBUG
0185   edm::LogVerbatim("HGCalGeom") << "There are " << layerType_.size() << " layers";
0186   for (unsigned int i = 0; i < layerType_.size(); ++i)
0187     edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerType_[i] << " sensitive class "
0188                                   << layerSense_[i];
0189 #endif
0190   materialsTop_ = vsArgs["TopMaterialNames"];
0191   namesTop_ = vsArgs["TopVolumeNames"];
0192   layerThickTop_ = vArgs["TopLayerThickness"];
0193   layerTypeTop_ = dbl_to_int(vArgs["TopLayerType"]);
0194   copyNumberTop_.resize(materialsTop_.size(), 1);
0195 #ifdef EDM_ML_DEBUG
0196   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << materialsTop_.size() << " types of volumes in the top part";
0197   for (unsigned int i = 0; i < materialsTop_.size(); ++i)
0198     edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << namesTop_[i] << " of thickness " << layerThickTop_[i]
0199                                   << " filled with " << materialsTop_[i] << " first copy number " << copyNumberTop_[i];
0200   edm::LogVerbatim("HGCalGeom") << "There are " << layerTypeTop_.size() << " layers in the top part";
0201   for (unsigned int i = 0; i < layerTypeTop_.size(); ++i)
0202     edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerTypeTop_[i];
0203 #endif
0204   materialsBot_ = vsArgs["BottomMaterialNames"];
0205   namesBot_ = vsArgs["BottomVolumeNames"];
0206   layerTypeBot_ = dbl_to_int(vArgs["BottomLayerType"]);
0207   layerSenseBot_ = dbl_to_int(vArgs["BottomLayerSense"]);
0208   layerThickBot_ = vArgs["BottomLayerThickness"];
0209   copyNumberBot_.resize(materialsBot_.size(), 1);
0210 #ifdef EDM_ML_DEBUG
0211   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << materialsBot_.size() << " types of volumes in the bottom part";
0212   for (unsigned int i = 0; i < materialsBot_.size(); ++i)
0213     edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << namesBot_[i] << " of thickness " << layerThickBot_[i]
0214                                   << " filled with " << materialsBot_[i] << " first copy number " << copyNumberBot_[i];
0215   edm::LogVerbatim("HGCalGeom") << "There are " << layerTypeBot_.size() << " layers in the bottom part";
0216   for (unsigned int i = 0; i < layerTypeBot_.size(); ++i)
0217     edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerTypeBot_[i]
0218                                   << " sensitive class " << layerSenseBot_[i];
0219 #endif
0220   zMinBlock_ = nArgs["zMinBlock"];
0221   rad100to200_ = vArgs["rad100to200"];
0222   rad200to300_ = vArgs["rad200to300"];
0223   zMinRadPar_ = nArgs["zMinForRadPar"];
0224   choiceType_ = (int)(nArgs["choiceType"]);
0225   nCutRadPar_ = (int)(nArgs["nCornerCut"]);
0226   fracAreaMin_ = nArgs["fracAreaMin"];
0227   waferSize_ = nArgs["waferSize"];
0228   waferSepar_ = nArgs["SensorSeparation"];
0229   sectors_ = (int)(nArgs["Sectors"]);
0230   alpha_ = (1._pi) / sectors_;
0231   cosAlpha_ = cos(alpha_);
0232 #ifdef EDM_ML_DEBUG
0233   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: zStart " << zMinBlock_ << " radius for wafer type separation uses "
0234                                 << rad100to200_.size() << " parameters; zmin " << zMinRadPar_ << " cutoff "
0235                                 << choiceType_ << ":" << nCutRadPar_ << ":" << fracAreaMin_ << " wafer width "
0236                                 << waferSize_ << " separations " << waferSepar_ << " sectors " << sectors_ << ":"
0237                                 << convertRadToDeg(alpha_) << ":" << cosAlpha_;
0238   for (unsigned int k = 0; k < rad100to200_.size(); ++k)
0239     edm::LogVerbatim("HGCalGeom") << "[" << k << "] 100-200 " << rad100to200_[k] << " 200-300 " << rad200to300_[k];
0240 #endif
0241   slopeB_ = vArgs["SlopeBottom"];
0242   zFrontB_ = vArgs["ZFrontBottom"];
0243   rMinFront_ = vArgs["RMinFront"];
0244   slopeT_ = vArgs["SlopeTop"];
0245   zFrontT_ = vArgs["ZFrontTop"];
0246   rMaxFront_ = vArgs["RMaxFront"];
0247 #ifdef EDM_ML_DEBUG
0248   for (unsigned int i = 0; i < slopeB_.size(); ++i)
0249     edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontB_[i] << " Rmin " << rMinFront_[i]
0250                                   << " Slope " << slopeB_[i];
0251   for (unsigned int i = 0; i < slopeT_.size(); ++i)
0252     edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << zFrontT_[i] << " Rmax " << rMaxFront_[i]
0253                                   << " Slope " << slopeT_[i];
0254 #endif
0255   nameSpace_ = DDCurrentNamespace::ns();
0256 #ifdef EDM_ML_DEBUG
0257   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: NameSpace " << nameSpace_;
0258 #endif
0259 
0260   waferType_ = std::make_unique<HGCalWaferType>(
0261       rad100to200_, rad200to300_, (waferSize_ + waferSepar_), zMinRadPar_, choiceType_, nCutRadPar_, fracAreaMin_);
0262 }
0263 
0264 ////////////////////////////////////////////////////////////////////
0265 // DDHGCalHEAlgo methods...
0266 ////////////////////////////////////////////////////////////////////
0267 
0268 void DDHGCalHEAlgo::execute(DDCompactView& cpv) {
0269 #ifdef EDM_ML_DEBUG
0270   edm::LogVerbatim("HGCalGeom") << "==>> Constructing DDHGCalHEAlgo...";
0271   copies_.clear();
0272 #endif
0273   constructLayers(parent(), cpv);
0274 #ifdef EDM_ML_DEBUG
0275   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << copies_.size() << " different wafer copy numbers";
0276   int k(0);
0277   for (std::unordered_set<int>::const_iterator itr = copies_.begin(); itr != copies_.end(); ++itr, ++k) {
0278     edm::LogVerbatim("HGCalGeom") << "Copy [" << k << "] : " << (*itr);
0279   }
0280   copies_.clear();
0281   edm::LogVerbatim("HGCalGeom") << "<<== End of DDHGCalHEAlgo construction...";
0282 #endif
0283 }
0284 
0285 void DDHGCalHEAlgo::constructLayers(const DDLogicalPart& module, DDCompactView& cpv) {
0286 #ifdef EDM_ML_DEBUG
0287   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: \t\tInside Layers";
0288 #endif
0289   double zi(zMinBlock_);
0290   int laymin(0);
0291   for (unsigned int i = 0; i < layers_.size(); i++) {
0292     double zo = zi + layerThick_[i];
0293     double routF = HGCalGeomTools::radius(zi, zFrontT_, rMaxFront_, slopeT_);
0294     int laymax = laymin + layers_[i];
0295     double zz = zi;
0296     double thickTot(0);
0297     for (int ly = laymin; ly < laymax; ++ly) {
0298       int ii = layerType_[ly];
0299       int copy = copyNumber_[ii];
0300       double hthick = 0.5 * thick_[ii];
0301       double rinB = HGCalGeomTools::radius(zo, zFrontB_, rMinFront_, slopeB_);
0302       zz += hthick;
0303       thickTot += thick_[ii];
0304 
0305       std::string name = names_[ii] + std::to_string(copy);
0306 #ifdef EDM_ML_DEBUG
0307       edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Layer " << ly << ":" << ii << " Front " << zi << ", " << routF
0308                                     << " Back " << zo << ", " << rinB << " superlayer thickness " << layerThick_[i];
0309 #endif
0310       DDName matName(DDSplit(materials_[ii]).first, DDSplit(materials_[ii]).second);
0311       DDMaterial matter(matName);
0312       DDLogicalPart glog;
0313       if (layerSense_[ly] < 1) {
0314         std::vector<double> pgonZ, pgonRin, pgonRout;
0315         if (layerSense_[ly] == 0 || absorbMode_ == 0) {
0316           double rmax =
0317               (std::min(routF, HGCalGeomTools::radius(zz + hthick, zFrontT_, rMaxFront_, slopeT_)) * cosAlpha_) - tol1_;
0318           pgonZ.emplace_back(-hthick);
0319           pgonZ.emplace_back(hthick);
0320           pgonRin.emplace_back(rinB);
0321           pgonRin.emplace_back(rinB);
0322           pgonRout.emplace_back(rmax);
0323           pgonRout.emplace_back(rmax);
0324         } else {
0325           HGCalGeomTools::radius(zz - hthick,
0326                                  zz + hthick,
0327                                  zFrontB_,
0328                                  rMinFront_,
0329                                  slopeB_,
0330                                  zFrontT_,
0331                                  rMaxFront_,
0332                                  slopeT_,
0333                                  -layerSense_[ly],
0334                                  pgonZ,
0335                                  pgonRin,
0336                                  pgonRout);
0337           for (unsigned int isec = 0; isec < pgonZ.size(); ++isec) {
0338             pgonZ[isec] -= zz;
0339             pgonRout[isec] = pgonRout[isec] * cosAlpha_ - tol1_;
0340           }
0341         }
0342         DDSolid solid =
0343             DDSolidFactory::polyhedra(DDName(name, nameSpace_), sectors_, -alpha_, 2._pi, pgonZ, pgonRin, pgonRout);
0344         glog = DDLogicalPart(solid.ddname(), matter, solid);
0345 #ifdef EDM_ML_DEBUG
0346         edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " polyhedra of " << sectors_
0347                                       << " sectors covering " << convertRadToDeg(-alpha_) << ":"
0348                                       << convertRadToDeg(-alpha_ + 2._pi) << " with " << pgonZ.size() << " sections";
0349         for (unsigned int k = 0; k < pgonZ.size(); ++k)
0350           edm::LogVerbatim("HGCalGeom") << "[" << k << "] z " << pgonZ[k] << " R " << pgonRin[k] << ":" << pgonRout[k];
0351 #endif
0352       } else {
0353         double rins = (sensitiveMode_ < 1) ? rinB : HGCalGeomTools::radius(zz + hthick, zFrontB_, rMinFront_, slopeB_);
0354         double routs =
0355             (sensitiveMode_ < 1) ? routF : HGCalGeomTools::radius(zz - hthick, zFrontT_, rMaxFront_, slopeT_);
0356         DDSolid solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rins, routs, 0.0, 2._pi);
0357         glog = DDLogicalPart(solid.ddname(), matter, solid);
0358 #ifdef EDM_ML_DEBUG
0359         edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " Tubs made of " << matName
0360                                       << " of dimensions " << rinB << ":" << rins << ", " << routF << ":" << routs
0361                                       << ", " << hthick << ", 0.0, 360.0 and positioned in: " << glog.name()
0362                                       << " number " << copy;
0363 #endif
0364         positionMix(glog, name, copy, thick_[ii], matter, rins, rMixLayer_[i], routs, zz, cpv);
0365       }
0366       DDTranslation r1(0, 0, zz);
0367       DDRotation rot;
0368       cpv.position(glog, module, copy, r1, rot);
0369       ++copyNumber_[ii];
0370 #ifdef EDM_ML_DEBUG
0371       edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << glog.name() << " number " << copy << " positioned in "
0372                                     << module.name() << " at " << r1 << " with no rotation";
0373 #endif
0374       zz += hthick;
0375     }  // End of loop over layers in a block
0376     zi = zo;
0377     laymin = laymax;
0378     if (std::abs(thickTot - layerThick_[i]) >= tol2_) {
0379       if (thickTot > layerThick_[i]) {
0380         edm::LogError("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " is smaller than " << thickTot
0381                                    << ": thickness of all its components **** ERROR ****";
0382       } else {
0383         edm::LogWarning("HGCalGeom") << "Thickness of the partition " << layerThick_[i] << " does not match with "
0384                                      << thickTot << " of the components";
0385       }
0386     }
0387   }  // End of loop over blocks
0388 }
0389 
0390 void DDHGCalHEAlgo::positionMix(const DDLogicalPart& glog,
0391                                 const std::string& nameM,
0392                                 int copyM,
0393                                 double thick,
0394                                 const DDMaterial& matter,
0395                                 double rin,
0396                                 double rmid,
0397                                 double rout,
0398                                 double zz,
0399                                 DDCompactView& cpv) {
0400   DDLogicalPart glog1;
0401   DDTranslation tran;
0402   DDRotation rot;
0403   for (unsigned int ly = 0; ly < layerTypeTop_.size(); ++ly) {
0404     int ii = layerTypeTop_[ly];
0405     copyNumberTop_[ii] = copyM;
0406   }
0407   for (unsigned int ly = 0; ly < layerTypeBot_.size(); ++ly) {
0408     int ii = layerTypeBot_[ly];
0409     copyNumberBot_[ii] = copyM;
0410   }
0411   double hthick = 0.5 * thick;
0412   // Make the top part first
0413   std::string name = nameM + "Top";
0414   DDSolid solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rmid, rout, 0.0, 2._pi);
0415   glog1 = DDLogicalPart(solid.ddname(), matter, solid);
0416 #ifdef EDM_ML_DEBUG
0417   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " Tubs made of " << matter.name()
0418                                 << " of dimensions " << rmid << ", " << rout << ", " << hthick << ", 0.0, 360.0";
0419 #endif
0420   cpv.position(glog1, glog, 1, tran, rot);
0421 #ifdef EDM_ML_DEBUG
0422   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << glog1.name() << " number 1 positioned in " << glog.name()
0423                                 << " at " << tran << " with no rotation";
0424 #endif
0425   double thickTot(0), zpos(-hthick);
0426   for (unsigned int ly = 0; ly < layerTypeTop_.size(); ++ly) {
0427     int ii = layerTypeTop_[ly];
0428     int copy = copyNumberTop_[ii];
0429     double hthickl = 0.5 * layerThickTop_[ii];
0430     thickTot += layerThickTop_[ii];
0431     name = namesTop_[ii] + std::to_string(copy);
0432 #ifdef EDM_ML_DEBUG
0433     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Layer " << ly << ":" << ii << " R " << rmid << ":" << rout
0434                                   << " Thick " << layerThickTop_[ii];
0435 #endif
0436     DDName matName(DDSplit(materialsTop_[ii]).first, DDSplit(materialsTop_[ii]).second);
0437     DDMaterial matter1(matName);
0438     solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthickl, rmid, rout, 0.0, 2._pi);
0439     DDLogicalPart glog2 = DDLogicalPart(solid.ddname(), matter1, solid);
0440 #ifdef EDM_ML_DEBUG
0441     double eta1 = -log(tan(0.5 * atan(rmid / zz)));
0442     double eta2 = -log(tan(0.5 * atan(rout / zz)));
0443     edm::LogVerbatim("HGCalGeom") << name << " z|rin|rout " << zz << ":" << rmid << ":" << rout << " eta " << eta1
0444                                   << ":" << eta2;
0445     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " Tubs made of " << matName
0446                                   << " of dimensions " << rmid << ", " << rout << ", " << hthickl << ", 0.0, 360.0";
0447 #endif
0448     zpos += hthickl;
0449     DDTranslation r1(0, 0, zpos);
0450     cpv.position(glog2, glog1, copy, r1, rot);
0451 #ifdef EDM_ML_DEBUG
0452     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Position " << glog2.name() << " number " << copy << " in "
0453                                   << glog1.name() << " at " << r1 << " with no rotation";
0454 #endif
0455     ++copyNumberTop_[ii];
0456     zpos += hthickl;
0457   }
0458   if (std::abs(thickTot - thick) >= tol2_) {
0459     if (thickTot > thick) {
0460       edm::LogError("HGCalGeom") << "Thickness of the partition " << thick << " is smaller than " << thickTot
0461                                  << ": thickness of all its components in the top part **** ERROR ****";
0462     } else {
0463       edm::LogWarning("HGCalGeom") << "Thickness of the partition " << thick << " does not match with " << thickTot
0464                                    << " of the components in top part";
0465     }
0466   }
0467 
0468   // Make the bottom part next
0469   name = nameM + "Bottom";
0470   solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthick, rin, rmid, 0.0, 2._pi);
0471   glog1 = DDLogicalPart(solid.ddname(), matter, solid);
0472 #ifdef EDM_ML_DEBUG
0473   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " Tubs made of " << matter.name()
0474                                 << " of dimensions " << rin << ", " << rmid << ", " << hthick << ", 0.0, 360.0";
0475 #endif
0476   cpv.position(glog1, glog, 1, tran, rot);
0477 #ifdef EDM_ML_DEBUG
0478   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << glog1.name() << " number 1 positioned in " << glog.name()
0479                                 << " at " << tran << " with no rotation";
0480 #endif
0481   thickTot = 0;
0482   zpos = -hthick;
0483   for (unsigned int ly = 0; ly < layerTypeBot_.size(); ++ly) {
0484     int ii = layerTypeBot_[ly];
0485     int copy = copyNumberBot_[ii];
0486     double hthickl = 0.5 * layerThickBot_[ii];
0487     thickTot += layerThickBot_[ii];
0488     name = namesBot_[ii] + std::to_string(copy);
0489 #ifdef EDM_ML_DEBUG
0490     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Layer " << ly << ":" << ii << " R " << rin << ":" << rmid
0491                                   << " Thick " << layerThickBot_[ii];
0492 #endif
0493     DDName matName(DDSplit(materialsBot_[ii]).first, DDSplit(materialsBot_[ii]).second);
0494     DDMaterial matter1(matName);
0495     solid = DDSolidFactory::tubs(DDName(name, nameSpace_), hthickl, rin, rmid, 0.0, 2._pi);
0496     DDLogicalPart glog2 = DDLogicalPart(solid.ddname(), matter1, solid);
0497 #ifdef EDM_ML_DEBUG
0498     double eta1 = -log(tan(0.5 * atan(rin / zz)));
0499     double eta2 = -log(tan(0.5 * atan(rmid / zz)));
0500     edm::LogVerbatim("HGCalGeom") << name << " z|rin|rout " << zz << ":" << rin << ":" << rmid << " eta " << eta1 << ":"
0501                                   << eta2;
0502     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << solid.name() << " Tubs made of " << matName
0503                                   << " of dimensions " << rin << ", " << rmid << ", " << hthickl << ", 0.0, 360.0";
0504 #endif
0505     zpos += hthickl;
0506     DDTranslation r1(0, 0, zpos);
0507     cpv.position(glog2, glog1, copy, r1, rot);
0508 #ifdef EDM_ML_DEBUG
0509     edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Position " << glog2.name() << " number " << copy << " in "
0510                                   << glog1.name() << " at " << r1 << " with no rotation";
0511 #endif
0512     if (layerSenseBot_[ly] != 0) {
0513 #ifdef EDM_ML_DEBUG
0514       edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: z " << (zz + zpos) << " Center " << copy << ":"
0515                                     << (copy - firstLayer_) << ":" << layerCenter_[copy - firstLayer_];
0516 #endif
0517       positionSensitive(glog2, rin, rmid, zz + zpos, layerSenseBot_[ly], layerCenter_[copy - firstLayer_], cpv);
0518     }
0519     zpos += hthickl;
0520     ++copyNumberBot_[ii];
0521   }
0522   if (std::abs(thickTot - thick) >= tol2_) {
0523     if (thickTot > thick) {
0524       edm::LogError("HGCalGeom") << "Thickness of the partition " << thick << " is smaller than " << thickTot
0525                                  << ": thickness of all its components in the top part **** ERROR ****";
0526     } else {
0527       edm::LogWarning("HGCalGeom") << "Thickness of the partition " << thick << " does not match with " << thickTot
0528                                    << " of the components in top part";
0529     }
0530   }
0531 }
0532 
0533 void DDHGCalHEAlgo::positionSensitive(const DDLogicalPart& glog,
0534                                       double rin,
0535                                       double rout,
0536                                       double zpos,
0537                                       int layertype,
0538                                       int layercenter,
0539                                       DDCompactView& cpv) {
0540   static const double sqrt3 = std::sqrt(3.0);
0541   double r = 0.5 * (waferSize_ + waferSepar_);
0542   double R = 2.0 * r / sqrt3;
0543   double dy = 0.75 * R;
0544   int N = (int)(0.5 * rout / r) + 2;
0545   std::pair<double, double> xyoff = geomTools_.shiftXY(layercenter, (waferSize_ + waferSepar_));
0546 #ifdef EDM_ML_DEBUG
0547   int ium(0), ivm(0), iumAll(0), ivmAll(0), kount(0), ntot(0), nin(0);
0548   std::vector<int> ntype(6, 0);
0549   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << glog.ddname() << " rin:rout " << rin << ":" << rout << " zpos "
0550                                 << zpos << " N " << N << " for maximum u, v Offset; Shift " << xyoff.first << ":"
0551                                 << xyoff.second << " WaferSize " << (waferSize_ + waferSepar_);
0552 #endif
0553   for (int u = -N; u <= N; ++u) {
0554     for (int v = -N; v <= N; ++v) {
0555       int nr = 2 * v;
0556       int nc = -2 * u + v;
0557       double xpos = xyoff.first + nc * r;
0558       double ypos = xyoff.second + nr * dy;
0559       std::pair<int, int> corner = HGCalGeomTools::waferCorner(xpos, ypos, r, R, rin, rout, false);
0560 #ifdef EDM_ML_DEBUG
0561       int iu = std::abs(u);
0562       int iv = std::abs(v);
0563       ++ntot;
0564 #endif
0565       if (corner.first > 0) {
0566         int type = waferType_->getType(xpos, ypos, zpos);
0567         int copy = HGCalTypes::packTypeUV(type, u, v);
0568         if (layertype > 1)
0569           type += 3;
0570 #ifdef EDM_ML_DEBUG
0571         if (iu > ium)
0572           ium = iu;
0573         if (iv > ivm)
0574           ivm = iv;
0575         kount++;
0576         if (copies_.count(copy) == 0)
0577           copies_.insert(copy);
0578 #endif
0579         if (corner.first == (int)(HGCalParameters::k_CornerSize)) {
0580 #ifdef EDM_ML_DEBUG
0581           if (iu > iumAll)
0582             iumAll = iu;
0583           if (iv > ivmAll)
0584             ivmAll = iv;
0585           ++nin;
0586 #endif
0587           DDTranslation tran(xpos, ypos, 0.0);
0588           DDRotation rotation;
0589           DDName name = DDName(DDSplit(wafers_[type]).first, DDSplit(wafers_[type]).second);
0590           cpv.position(name, glog.ddname(), copy, tran, rotation);
0591 #ifdef EDM_ML_DEBUG
0592           ++ntype[type];
0593           edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: " << name << " number " << copy << " positioned in "
0594                                         << glog.ddname() << " at " << tran << " with no rotation";
0595 #endif
0596         }
0597       }
0598     }
0599   }
0600 #ifdef EDM_ML_DEBUG
0601   edm::LogVerbatim("HGCalGeom") << "DDHGCalHEAlgo: Maximum # of u " << ium << ":" << iumAll << " # of v " << ivm << ":"
0602                                 << ivmAll << " and " << nin << ":" << kount << ":" << ntot << " wafers (" << ntype[0]
0603                                 << ":" << ntype[1] << ":" << ntype[2] << ":" << ntype[3] << ":" << ntype[4] << ":"
0604                                 << ntype[5] << ") for " << glog.ddname() << " R " << rin << ":" << rout;
0605 #endif
0606 }
0607 
0608 DEFINE_EDM_PLUGIN(DDAlgorithmFactory, DDHGCalHEAlgo, "hgcal:DDHGCalHEAlgo");