Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 23:30:05

0001 /////////////////////////////////////////////////////////////////////////////
0002 // File: DDHCalBarrelAlgo.cc
0003 //   adapted from CCal(G4)HcalBarrel.cc
0004 // Description: Geometry factory class for Hcal Barrel
0005 ///////////////////////////////////////////////////////////////////////////////
0006 
0007 #include <cmath>
0008 #include <algorithm>
0009 #include <map>
0010 #include <string>
0011 #include <vector>
0012 
0013 #include "DataFormats/Math/interface/angle_units.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "DetectorDescription/Core/interface/DDSplit.h"
0016 #include "DetectorDescription/DDCMS/interface/DDPlugins.h"
0017 #include "DetectorDescription/DDCMS/interface/DDutils.h"
0018 #include "DD4hep/DetFactoryHelper.h"
0019 
0020 //#define EDM_ML_DEBUG
0021 using namespace angle_units::operators;
0022 
0023 struct HcalBarrelAlgo {
0024   //General Volume
0025   //      <----- Zmax ------>
0026   //Router************************-------
0027   //      *                      *Rstep2|        Theta angle w.r.t. vertical
0028   //      *                      *---------------
0029   //      *                     *               |
0030   //      *                    *Theta[i]        Rmax[i]
0031   //      *                   *---------------  |
0032   //                        *Theta[0] Rmax[0]|  |
0033   //Rinner*****************----------------------
0034 
0035   std::string genMaterial;     //General material
0036   int nsectors;                //Number of potenital straight edges
0037   int nsectortot;              //Number of straight edges (actual)
0038   int nhalf;                   //Number of half modules
0039   double rinner, router;       //See picture
0040   int rzones;                  //  ....
0041   std::vector<double> theta;   //  .... (in degrees)
0042   std::vector<double> rmax;    //  ....
0043   std::vector<double> zoff;    //  ....
0044   std::vector<double> ttheta;  //tan(theta)
0045   std::string rotHalf;         //Rotation matrix of the second half
0046   std::string rotns;           //Name space for Rotation matrices
0047 
0048   //Upper layers inside general volume
0049   //     <---- Zout ---->
0050   //  |  ****************     |
0051   //  |  *              *     Wstep
0052   //  W  *              ***** |
0053   //  |  *                  *
0054   //  |  ********************
0055   //     <------ Zin ------->
0056   // Middle layers inside general volume
0057   //     <------ Zout ------>         Zout = Full sector Z at position
0058   //  |  ********************         Zin  = Full sector Z at position
0059   //  |  *                 *
0060   //  W  *                * Angle = Theta sector
0061   //  |  *               *  )
0062   //  |  ****************--------
0063   //     <------ Zin ------->
0064 
0065   // Lower layers
0066   //     <------ Zout ------>         Zin(i)=Zout(i-1)
0067   //  |  ********************         Zout(i)=Zin(i)+W(i)/tan(Theta(i))
0068   //  |  *                 *
0069   //  W  *                *  Theta
0070   //  |  *               *
0071   //  |  ****************--------
0072   //     <--- Zin ------>
0073 
0074   int nLayers;                          //Number of layers
0075   std::vector<int> layerId;             //Number identification
0076   std::vector<std::string> layerLabel;  //String identification
0077   std::vector<std::string> layerMat;    //Material
0078   std::vector<double> layerWidth;       //W in picture
0079   std::vector<double> layerD1;          //d1 in front picture
0080   std::vector<double> layerD2;          //d2 in front picture
0081   std::vector<double> layerAlpha;       //Angular width of the middle tiles
0082   std::vector<double> layerT1;          //t in front picture (side)
0083   std::vector<double> layerT2;          //t in front picture (front)
0084   std::vector<int> layerAbsorb;         //Absorber flag
0085   std::vector<double> layerGap;         //Gap at the edge
0086 
0087   int nAbsorber;                        //Number of absorber layers in middle
0088   std::vector<std::string> absorbName;  //Absorber name
0089   std::vector<std::string> absorbMat;   //Absorber material
0090   std::vector<double> absorbD;          //Distance from the bottom surface
0091   std::vector<double> absorbT;          //Thickness
0092   std::string middleMat;                //Material of the detector layer
0093   double middleD;                       //Distance from the bottom surface
0094   double middleW;                       //Half width
0095   int nMidAbs;                          //Number of absorbers in front part
0096   std::vector<std::string> midName;     //Absorber names in the front part
0097   std::vector<std::string> midMat;      //Absorber material
0098   std::vector<double> midW;             //Half width
0099   std::vector<double> midT;             //Thickness
0100 
0101   std::vector<std::string> sideMat;      //Material for special side layers
0102   std::vector<double> sideD;             //Depth from bottom surface
0103   std::vector<double> sideT;             //Thickness
0104   int nSideAbs;                          //Number of absorbers in special side
0105   std::vector<std::string> sideAbsName;  //Absorber name
0106   std::vector<std::string> sideAbsMat;   //Absorber material
0107   std::vector<double> sideAbsW;          //Half width
0108 
0109   // Detectors. Each volume inside the layer has the shape:
0110   //
0111   // ******************************* |
0112   // *\\\\\\\Plastic\\\\\\\\\\\\\\\* T2
0113   // ******************************* |
0114   // *////Scintillator/////////////* Tsc
0115   // ******************************* |
0116   // *\\\\\\\Plastic\\\\\\\\\\\\\\\* T1
0117   // ******************************* |   |
0118   // *         Air                 *     dP1
0119   // *******************************     |
0120   //
0121   std::string detMat;    //fill material
0122   std::string detRot;    //Rotation matrix for the 2nd
0123   std::string detMatPl;  //Plastic material
0124   std::string detMatSc;  //Scintillator material
0125   std::vector<int> detType;
0126   std::vector<double> detdP1;     //Air gap (side)
0127   std::vector<double> detdP2;     //Air gap (centre)
0128   std::vector<double> detT11;     //Back plastic thickness (side)
0129   std::vector<double> detT12;     //Back plastic thickness (centre)
0130   std::vector<double> detTsc;     //Scintillator
0131   std::vector<double> detT21;     //Front plastic thickness (side)
0132   std::vector<double> detT22;     //Front plastic thickness (centre)
0133   std::vector<double> detWidth1;  //Width of phi(1,4) megatiles
0134   std::vector<double> detWidth2;  //Width of phi(2,3) megatiles
0135   std::vector<int> detPosY;       //Positioning of phi(1,4) tiles - 0 centre
0136 
0137   std::string idName;       //Name of the "parent" volume.
0138   std::string idNameSpace;  //Namespace of this and ALL sub-parts
0139   int idOffset;             // Geant4 ID's...    = 3000;
0140 
0141   HcalBarrelAlgo() = delete;
0142 
0143   HcalBarrelAlgo(cms::DDParsingContext& ctxt, xml_h& e) {
0144     cms::DDNamespace ns(ctxt, e, true);
0145     cms::DDAlgoArguments args(ctxt, e);
0146 #ifdef EDM_ML_DEBUG
0147     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Creating an instance";
0148 #endif
0149 
0150     genMaterial = args.value<std::string>("MaterialName");
0151     nsectors = args.value<int>("NSector");
0152     nsectortot = args.value<int>("NSectorTot");
0153     nhalf = args.value<int>("NHalf");
0154     rinner = args.value<double>("RIn");
0155     router = args.value<double>("ROut");
0156     rzones = args.value<int>("RZones");
0157     rotHalf = args.value<std::string>("RotHalf");
0158     rotns = args.value<std::string>("RotNameSpace");
0159 
0160     theta = args.value<std::vector<double> >("Theta");
0161     rmax = args.value<std::vector<double> >("RMax");
0162     zoff = args.value<std::vector<double> >("ZOff");
0163     for (int i = 0; i < rzones; i++) {
0164       ttheta.emplace_back(tan(theta[i]));  //*deg already done in XML
0165     }
0166     if (rzones > 3)
0167       rmax[2] = (zoff[3] - zoff[2]) / ttheta[2];
0168 
0169 #ifdef EDM_ML_DEBUG
0170     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: General material " << genMaterial << "\tSectors " << nsectors
0171                                  << ", " << nsectortot << "\tHalves " << nhalf << "\tRotation matrix " << rotns << ":"
0172                                  << rotHalf << "\n\t\t" << cms::convert2mm(rinner) << "\t" << cms::convert2mm(router)
0173                                  << "\t" << rzones;
0174     for (int i = 0; i < rzones; i++)
0175       edm::LogVerbatim("HCalGeom") << "\tTheta[" << i << "] = " << theta[i] << "\trmax[" << i
0176                                    << "] = " << cms::convert2mm(rmax[i]) << "\tzoff[" << i
0177                                    << "] = " << cms::convert2mm(zoff[i]);
0178 #endif
0179     ///////////////////////////////////////////////////////////////
0180     //Layers
0181     nLayers = args.value<int>("NLayers");
0182 #ifdef EDM_ML_DEBUG
0183     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Layer\t" << nLayers;
0184 #endif
0185     layerId = args.value<std::vector<int> >("Id");
0186     layerLabel = args.value<std::vector<std::string> >("LayerLabel");
0187     layerMat = args.value<std::vector<std::string> >("LayerMat");
0188     layerWidth = args.value<std::vector<double> >("LayerWidth");
0189     layerD1 = args.value<std::vector<double> >("D1");
0190     layerD2 = args.value<std::vector<double> >("D2");
0191     layerAlpha = args.value<std::vector<double> >("Alpha2");
0192     layerT1 = args.value<std::vector<double> >("T1");
0193     layerT2 = args.value<std::vector<double> >("T2");
0194     layerAbsorb = args.value<std::vector<int> >("AbsL");
0195     layerGap = args.value<std::vector<double> >("Gap");
0196 #ifdef EDM_ML_DEBUG
0197     for (int i = 0; i < nLayers; i++)
0198       edm::LogVerbatim("HCalGeom") << layerLabel[i] << "\t" << layerId[i] << "\t" << layerMat[i] << "\t"
0199                                    << cms::convert2mm(layerWidth[i]) << "\t" << cms::convert2mm(layerD1[i]) << "\t"
0200                                    << cms::convert2mm(layerD2[i]) << "\t" << layerAlpha[i] << "\t"
0201                                    << cms::convert2mm(layerT1[i]) << "\t" << cms::convert2mm(layerT2[i]) << "\t"
0202                                    << layerAbsorb[i] << "\t" << cms::convert2mm(layerGap[i]);
0203 #endif
0204 
0205     ///////////////////////////////////////////////////////////////
0206     //Absorber Layers and middle part
0207     absorbName = args.value<std::vector<std::string> >("AbsorbName");
0208     absorbMat = args.value<std::vector<std::string> >("AbsorbMat");
0209     absorbD = args.value<std::vector<double> >("AbsorbD");
0210     absorbT = args.value<std::vector<double> >("AbsorbT");
0211     nAbsorber = absorbName.size();
0212 #ifdef EDM_ML_DEBUG
0213     for (int i = 0; i < nAbsorber; i++)
0214       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << absorbName[i] << " Material " << absorbMat[i] << " d "
0215                                    << cms::convert2mm(absorbD[i]) << " t " << cms::convert2mm(absorbT[i]);
0216 #endif
0217     middleMat = args.value<std::string>("MiddleMat");
0218     middleD = args.value<double>("MiddleD");
0219     middleW = args.value<double>("MiddleW");
0220 #ifdef EDM_ML_DEBUG
0221     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Middle material " << middleMat << " d "
0222                                  << cms::convert2mm(middleD) << " w " << cms::convert2mm(middleW);
0223 #endif
0224     midName = args.value<std::vector<std::string> >("MidAbsName");
0225     midMat = args.value<std::vector<std::string> >("MidAbsMat");
0226     midW = args.value<std::vector<double> >("MidAbsW");
0227     midT = args.value<std::vector<double> >("MidAbsT");
0228     nMidAbs = midName.size();
0229 #ifdef EDM_ML_DEBUG
0230     for (int i = 0; i < nMidAbs; i++)
0231       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << midName[i] << " Material " << midMat[i] << " W "
0232                                    << cms::convert2mm(midW[i]) << " T " << cms::convert2mm(midT[i]);
0233 #endif
0234 
0235     //Absorber layers in the side part
0236     sideMat = args.value<std::vector<std::string> >("SideMat");
0237     sideD = args.value<std::vector<double> >("SideD");
0238     sideT = args.value<std::vector<double> >("SideT");
0239 #ifdef EDM_ML_DEBUG
0240     for (unsigned int i = 0; i < sideMat.size(); i++)
0241       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Side material " << sideMat[i] << " d "
0242                                    << cms::convert2mm(sideD[i]) << " t " << cms::convert2mm(sideT[i]);
0243 #endif
0244     sideAbsName = args.value<std::vector<std::string> >("SideAbsName");
0245     sideAbsMat = args.value<std::vector<std::string> >("SideAbsMat");
0246     sideAbsW = args.value<std::vector<double> >("SideAbsW");
0247     nSideAbs = sideAbsName.size();
0248 #ifdef EDM_ML_DEBUG
0249     for (int i = 0; i < nSideAbs; i++)
0250       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << sideAbsName[i] << " Material " << sideAbsMat[i] << " W "
0251                                    << cms::convert2mm(sideAbsW[i]);
0252 #endif
0253 
0254     ///////////////////////////////////////////////////////////////
0255     // Detectors
0256 
0257     detMat = args.value<std::string>("DetMat");
0258     detRot = args.value<std::string>("DetRot");
0259     detMatPl = args.value<std::string>("DetMatPl");
0260     detMatSc = args.value<std::string>("DetMatSc");
0261 #ifdef EDM_ML_DEBUG
0262     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Detector (" << nLayers << ") Rotation matrix " << rotns << ":"
0263                                  << detRot << "\n\t\t" << detMat << "\t" << detMatPl << "\t" << detMatSc;
0264 #endif
0265     detType = args.value<std::vector<int> >("DetType");
0266     detdP1 = args.value<std::vector<double> >("DetdP1");
0267     detdP2 = args.value<std::vector<double> >("DetdP2");
0268     detT11 = args.value<std::vector<double> >("DetT11");
0269     detT12 = args.value<std::vector<double> >("DetT12");
0270     detTsc = args.value<std::vector<double> >("DetTsc");
0271     detT21 = args.value<std::vector<double> >("DetT21");
0272     detT22 = args.value<std::vector<double> >("DetT22");
0273     detWidth1 = args.value<std::vector<double> >("DetWidth1");
0274     detWidth2 = args.value<std::vector<double> >("DetWidth2");
0275     detPosY = args.value<std::vector<int> >("DetPosY");
0276 #ifdef EDM_ML_DEBUG
0277     for (int i = 0; i < nLayers; i++)
0278       edm::LogVerbatim("HCalGeom") << i + 1 << "\t" << detType[i] << "\t" << cms::convert2mm(detdP1[i]) << ", "
0279                                    << cms::convert2mm(detdP2[i]) << "\t" << cms::convert2mm(detT11[i]) << ", "
0280                                    << cms::convert2mm(detT12[i]) << "\t" << cms::convert2mm(detTsc[i]) << "\t"
0281                                    << cms::convert2mm(detT21[i]) << ", " << cms::convert2mm(detT22[i]) << "\t"
0282                                    << cms::convert2mm(detWidth1[i]) << "\t" << cms::convert2mm(detWidth2[i]) << "\t"
0283                                    << detPosY[i];
0284 #endif
0285 
0286     //  idName = parentName.name();
0287     idName = args.value<std::string>("MotherName");
0288     idNameSpace = ns.name();
0289     idOffset = args.value<int>("IdOffset");
0290 #ifdef EDM_ML_DEBUG
0291     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Parent " << args.parentName() << " idName " << idName
0292                                  << " NameSpace " << idNameSpace << " Offset " << idOffset;
0293     edm::LogVerbatim("HCalGeom") << "==>> Constructing DDHCalBarrelAlgo...\n"
0294                                  << "DDHCalBarrelAlgo: General volume...";
0295 #endif
0296 
0297     double alpha = (1._pi) / nsectors;
0298     double dphi = nsectortot * (2._pi) / nsectors;
0299     int nsec, ntot(15);
0300     if (nhalf == 1)
0301       nsec = 8;
0302     else
0303       nsec = 15;
0304     int nf = ntot - nsec;
0305 
0306     //Calculate zmin... see HCalBarrel.hh picture. For polyhedra
0307     //Rmin and Rmax are distances to vertex
0308     double zmax = zoff[3];
0309     double zstep5 = zoff[4];
0310     double zstep4 = (zoff[1] + rmax[1] * ttheta[1]);
0311     if ((zoff[2] + rmax[1] * ttheta[2]) > zstep4)
0312       zstep4 = (zoff[2] + rmax[1] * ttheta[2]);
0313     double zstep3 = (zoff[1] + rmax[0] * ttheta[1]);
0314     double zstep2 = (zoff[0] + rmax[0] * ttheta[0]);
0315     double zstep1 = (zoff[0] + rinner * ttheta[0]);
0316     double rout = router;
0317     double rout1 = rmax[3];
0318     double rin = rinner;
0319     double rmid1 = rmax[0];
0320     double rmid2 = rmax[1];
0321     double rmid3 = (zoff[4] - zoff[2]) / ttheta[2];
0322     double rmid4 = rmax[2];
0323 
0324     std::vector<double> pgonZ = {-zmax,
0325                                  -zstep5,
0326                                  -zstep5,
0327                                  -zstep4,
0328                                  -zstep3,
0329                                  -zstep2,
0330                                  -zstep1,
0331                                  0,
0332                                  zstep1,
0333                                  zstep2,
0334                                  zstep3,
0335                                  zstep4,
0336                                  zstep5,
0337                                  zstep5,
0338                                  zmax};
0339 
0340     std::vector<double> pgonRmin = {
0341         rmid4, rmid3, rmid3, rmid2, rmid1, rmid1, rin, rin, rin, rmid1, rmid1, rmid2, rmid3, rmid3, rmid4};
0342 
0343     std::vector<double> pgonRmax = {
0344         rout1, rout1, rout, rout, rout, rout, rout, rout, rout, rout, rout, rout, rout, rout1, rout1};
0345 
0346     std::vector<double> pgonZHalf = {0, zstep1, zstep2, zstep3, zstep4, zstep5, zstep5, zmax};
0347 
0348     std::vector<double> pgonRminHalf = {rin, rin, rmid1, rmid1, rmid2, rmid3, rmid3, rmid4};
0349 
0350     std::vector<double> pgonRmaxHalf = {rout, rout, rout, rout, rout, rout, rout1, rout1};
0351 
0352     std::string name("Null");
0353     dd4hep::Solid solid;
0354     if (nf == 0) {
0355       solid = dd4hep::Polyhedra(ns.prepend(idName), nsectortot, -alpha, dphi, pgonZ, pgonRmin, pgonRmax);
0356 #ifdef EDM_ML_DEBUG
0357       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << genMaterial
0358                                    << " with " << nsectortot << " sectors from " << convertRadToDeg(-alpha) << " to "
0359                                    << convertRadToDeg(-alpha + dphi) << " and with " << nsec << " sections ";
0360       for (unsigned int i = 0; i < pgonZ.size(); i++)
0361         edm::LogVerbatim("HCalGeom") << "\t"
0362                                      << "\tZ = " << cms::convert2mm(pgonZ[i])
0363                                      << "\tRmin = " << cms::convert2mm(pgonRmin[i])
0364                                      << "\tRmax = " << cms::convert2mm(pgonRmax[i]);
0365 #endif
0366     } else {
0367       solid = dd4hep::Polyhedra(ns.prepend(idName), nsectortot, -alpha, dphi, pgonZHalf, pgonRminHalf, pgonRmaxHalf);
0368 #ifdef EDM_ML_DEBUG
0369       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << genMaterial
0370                                    << " with " << nsectortot << " sectors from " << convertRadToDeg(-alpha) << " to "
0371                                    << convertRadToDeg(-alpha + dphi) << " and with " << nsec << " sections ";
0372       for (unsigned int i = 0; i < pgonZHalf.size(); i++)
0373         edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZHalf[i])
0374                                      << "\tRmin = " << cms::convert2mm(pgonRminHalf[i])
0375                                      << "\tRmax = " << cms::convert2mm(pgonRmaxHalf[i]);
0376 #endif
0377     }
0378 
0379     dd4hep::Material matter = ns.material(genMaterial);
0380     dd4hep::Volume genlogic(solid.name(), solid, matter);
0381     dd4hep::Volume parentName = ns.volume(args.parentName());
0382     parentName.placeVolume(genlogic, 1);
0383 #ifdef EDM_ML_DEBUG
0384     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << genlogic.name() << " number 1 positioned in "
0385                                  << parentName.name() << " at (0, 0, 0) with no rotation";
0386 #endif
0387     //Forward and backwards halfs
0388     name = idName + "Half";
0389     nf = (ntot + 1) / 2;
0390     solid = dd4hep::Polyhedra(ns.prepend(name), nsectortot, -alpha, dphi, pgonZHalf, pgonRminHalf, pgonRmaxHalf);
0391 #ifdef EDM_ML_DEBUG
0392     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << genMaterial
0393                                  << " with " << nsectortot << " sectors from " << convertRadToDeg(-alpha) << " to "
0394                                  << convertRadToDeg(-alpha + dphi) << " and with " << nf << " sections ";
0395     for (unsigned int i = 0; i < pgonZHalf.size(); i++)
0396       edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZHalf[i])
0397                                    << "\tRmin = " << cms::convert2mm(pgonRminHalf[i])
0398                                    << "\tRmax = " << cms::convert2mm(pgonRmaxHalf[i]);
0399 #endif
0400     dd4hep::Volume genlogich(solid.name(), solid, matter);
0401     genlogic.placeVolume(genlogich, 1);
0402 #ifdef EDM_ML_DEBUG
0403     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << genlogich.name() << " number 1 positioned in "
0404                                  << genlogic.name() << " at (0, 0, 0) with no rotation";
0405 #endif
0406     if (nhalf != 1) {
0407       dd4hep::Rotation3D rot = getRotation(rotHalf, rotns, ns);
0408       genlogic.placeVolume(genlogich, 2, rot);
0409 #ifdef EDM_ML_DEBUG
0410       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << genlogich.name() << " number 2 positioned in "
0411                                    << genlogic.name() << " at (0, 0, 0) with " << rot;
0412 #endif
0413     }  //end if (getNhalf...
0414 
0415     //Construct sector (from -alpha to +alpha)
0416     name = idName + "Module";
0417     solid = dd4hep::Polyhedra(ns.prepend(name), 1, -alpha, 2 * alpha, pgonZHalf, pgonRminHalf, pgonRmaxHalf);
0418     dd4hep::Volume seclogic(solid.name(), solid, matter);
0419 #ifdef EDM_ML_DEBUG
0420     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << matter.name()
0421                                  << " with 1 sector from " << convertRadToDeg(-alpha) << " to "
0422                                  << convertRadToDeg(alpha) << " and with " << nf << " sections";
0423     for (unsigned int i = 0; i < pgonZHalf.size(); i++)
0424       edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZHalf[i])
0425                                    << "\tRmin = " << cms::convert2mm(pgonRminHalf[i])
0426                                    << "\tRmax = " << cms::convert2mm(pgonRmaxHalf[i]);
0427 #endif
0428 
0429     for (int ii = 0; ii < nsectortot; ii++) {
0430       double phi = ii * 2 * alpha;
0431       dd4hep::Rotation3D rotation;
0432       if (phi != 0) {
0433         rotation = dd4hep::RotationZ(phi);
0434 #ifdef EDM_ML_DEBUG
0435         edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Creating a new "
0436                                      << "rotation around Z " << convertRadToDeg(phi);
0437 #endif
0438       }  //if phideg!=0
0439       genlogich.placeVolume(seclogic, ii + 1, rotation);
0440 #ifdef EDM_ML_DEBUG
0441       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << seclogic.name() << " number " << ii + 1
0442                                    << " positioned in " << genlogich.name() << " at (0, 0, 0) with " << rotation;
0443 #endif
0444     }
0445 
0446     //Construct the things inside the sector
0447 #ifdef EDM_ML_DEBUG
0448     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: Layers (" << nLayers << ") ...";
0449 #endif
0450     rin = rinner;
0451     for (int i = 0; i < nLayers; i++) {
0452       std::string name = idName + layerLabel[i];
0453       dd4hep::Material matter = ns.material(layerMat[i]);
0454       double width = layerWidth[i];
0455       double rout = rin + width;
0456 
0457       int in = 0, out = 0;
0458       for (int j = 0; j < rzones - 1; j++) {
0459         if (rin >= rmax[j])
0460           in = j + 1;
0461         if (rout > rmax[j])
0462           out = j + 1;
0463       }
0464       double zout = zoff[in] + rin * ttheta[in];
0465 
0466       //!!!!!!!!!!!!!!!!!Should be zero. And removed as soon as
0467       //vertical walls are allowed in SolidPolyhedra
0468       int nsec = 2;
0469       std::vector<double> pgonZ, pgonRmin, pgonRmax;
0470       // index 0
0471       pgonZ.emplace_back(0);
0472       pgonRmin.emplace_back(rin);
0473       pgonRmax.emplace_back(rout);
0474       // index 1
0475       pgonZ.emplace_back(zout);
0476       pgonRmin.emplace_back(rin);
0477       pgonRmax.emplace_back(rout);
0478       if (in == out) {
0479         if (in <= 3) {
0480           //index 2
0481           pgonZ.emplace_back(zoff[in] + rout * ttheta[in]);
0482           pgonRmin.emplace_back(pgonRmax[1]);
0483           pgonRmax.emplace_back(pgonRmax[1]);
0484           nsec++;
0485         }
0486       } else {
0487         if (in == 3) {
0488           //redo index 1, add index 2
0489           pgonZ[1] = (zoff[out] + rmax[out] * ttheta[out]);
0490           pgonZ.emplace_back(pgonZ[1]);
0491           pgonRmin.emplace_back(pgonRmin[1]);
0492           pgonRmax.emplace_back(rmax[in]);
0493           //index 3
0494           pgonZ.emplace_back(zoff[in] + rmax[in] * ttheta[in]);
0495           pgonRmin.emplace_back(pgonRmin[2]);
0496           pgonRmax.emplace_back(pgonRmax[2]);
0497           nsec += 2;
0498         } else {
0499           //index 2
0500           pgonZ.emplace_back(zoff[in] + rmax[in] * ttheta[in]);
0501           pgonRmin.emplace_back(rmax[in]);
0502           pgonRmax.emplace_back(pgonRmax[1]);
0503           nsec++;
0504           if (in == 0) {
0505             pgonZ.emplace_back(zoff[out] + rmax[in] * ttheta[out]);
0506             pgonRmin.emplace_back(pgonRmin[2]);
0507             pgonRmax.emplace_back(pgonRmax[2]);
0508             nsec++;
0509           }
0510           if (in <= 1) {
0511             pgonZ.emplace_back(zoff[out] + rout * ttheta[out]);
0512             pgonRmin.emplace_back(rout);
0513             pgonRmax.emplace_back(rout);
0514             nsec++;
0515           }
0516         }
0517       }
0518       //Solid & volume
0519       dd4hep::Solid solid;
0520       double alpha1 = alpha;
0521       if (layerGap[i] > 1.e-6) {
0522         double rmid = 0.5 * (rin + rout);
0523         double width = rmid * tan(alpha) - layerGap[i];
0524         alpha1 = atan(width / rmid);
0525 #ifdef EDM_ML_DEBUG
0526         edm::LogVerbatim("HCalGeom") << "\tAlpha_1 modified from " << convertRadToDeg(alpha) << " to "
0527                                      << convertRadToDeg(alpha1) << " Rmid " << cms::convert2mm(rmid)
0528                                      << " Reduced width " << cms::convert2mm(width);
0529 #endif
0530       }
0531       solid = dd4hep::Polyhedra(ns.prepend(name), 1, -alpha1, 2 * alpha1, pgonZ, pgonRmin, pgonRmax);
0532       dd4hep::Volume glog(solid.name(), solid, matter);
0533 #ifdef EDM_ML_DEBUG
0534       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " (Layer " << i << ") Polyhedra made of "
0535                                    << matter.name() << " with 1 sector from " << convertRadToDeg(-alpha1) << " to "
0536                                    << convertRadToDeg(alpha1) << " and with " << nsec << " sections";
0537       for (unsigned int k = 0; k < pgonZ.size(); k++)
0538         edm::LogVerbatim("HCalGeom") << "\t\t" << cms::convert2mm(pgonZ[k]) << "\t" << cms::convert2mm(pgonRmin[k])
0539                                      << "\t" << cms::convert2mm(pgonRmax[k]);
0540 #endif
0541 
0542       seclogic.placeVolume(glog, layerId[i]);
0543 #ifdef EDM_ML_DEBUG
0544       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " number " << layerId[i]
0545                                    << " positioned in " << seclogic.name() << " at (0,0,0) with no rotation";
0546 #endif
0547       constructInsideLayers(glog,
0548                             layerLabel[i],
0549                             layerId[i],
0550                             layerAbsorb[i],
0551                             rin,
0552                             layerD1[i],
0553                             alpha1,
0554                             layerD2[i],
0555                             layerAlpha[i],
0556                             layerT1[i],
0557                             layerT2[i],
0558                             ns);
0559       rin = rout;
0560     }
0561   }
0562 
0563   void constructInsideLayers(dd4hep::Volume& laylog,
0564                              const std::string& nm,
0565                              int id,
0566                              int nAbs,
0567                              double rin,
0568                              double d1,
0569                              double alpha1,
0570                              double d2,
0571                              double alpha2,
0572                              double t1,
0573                              double t2,
0574                              cms::DDNamespace& ns) {
0575 #ifdef EDM_ML_DEBUG
0576     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: \t\tInside layer " << id << "...";
0577 #endif
0578 
0579     ///////////////////////////////////////////////////////////////
0580     //Pointers to the Rotation Matrices and to the Materials
0581     dd4hep::Rotation3D rot = getRotation(detRot, rotns, ns);
0582 
0583     std::string nam0 = nm + "In";
0584     std::string name = idName + nam0;
0585     dd4hep::Material matter = ns.material(detMat);
0586 
0587     dd4hep::Solid solid;
0588     dd4hep::Volume glog, mother;
0589     double rsi, dx, dy, dz, x, y;
0590     int i, in;
0591     //Two lower volumes
0592     if (alpha1 > 0) {
0593       rsi = rin + d1;
0594       in = 0;
0595       for (i = 0; i < rzones - 1; i++) {
0596         if (rsi >= rmax[i])
0597           in = i + 1;
0598       }
0599       dx = 0.5 * t1;
0600       dy = 0.5 * rsi * (tan(alpha1) - tan(alpha2));
0601       dz = 0.5 * (zoff[in] + rsi * ttheta[in]);
0602       x = rsi + dx;
0603       y = 0.5 * rsi * (tan(alpha1) + tan(alpha2));
0604       dd4hep::Position r11(x, y, dz);
0605       dd4hep::Position r12(x, -y, dz);
0606 
0607       solid = dd4hep::Box(ns.prepend(name + "1"), dx, dy, dz);
0608 #ifdef EDM_ML_DEBUG
0609       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Box made of " << matter.name()
0610                                    << " of dimensions " << cms::convert2mm(dx) << ", " << cms::convert2mm(dy) << ", "
0611                                    << cms::convert2mm(dz);
0612 #endif
0613       glog = dd4hep::Volume(solid.name(), solid, matter);
0614 
0615       if (nAbs != 0) {
0616         mother = constructSideLayer(laylog, name, nAbs, rin, alpha1, ns);
0617       } else {
0618         mother = laylog;
0619       }
0620       mother.placeVolume(glog, idOffset + 1, r11);
0621       mother.placeVolume(glog, idOffset + 2, dd4hep::Transform3D(rot, r12));
0622 #ifdef EDM_ML_DEBUG
0623       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number " << (idOffset + 1)
0624                                    << " positioned in " << mother.name() << " at (" << cms::convert2mm(x) << ","
0625                                    << cms::convert2mm(y) << "," << cms::convert2mm(dz) << ") with no rotation\n"
0626                                    << "DDHCalBarrelAlgo: " << glog.name() << " Number " << (idOffset + 2)
0627                                    << " positioned in " << mother.name() << " at (" << cms::convert2mm(x) << ","
0628                                    << -cms::convert2mm(y) << "," << cms::convert2mm(dz) << ") with " << rot;
0629 #endif
0630       //Constructin the plastics and scintillators inside
0631       constructInsideDetectors(glog, nam0 + "1", id, dx, dy, dz, 1, ns);
0632     }
0633 
0634     //Upper volume
0635     rsi = rin + d2;
0636     in = 0;
0637     for (i = 0; i < rzones - 1; i++) {
0638       if (rsi >= rmax[i])
0639         in = i + 1;
0640     }
0641     dx = 0.5 * t2;
0642     dy = 0.5 * rsi * tan(alpha2);
0643     dz = 0.5 * (zoff[in] + rsi * ttheta[in]);
0644     x = rsi + dx;
0645     dd4hep::Position r21(x, dy, dz);
0646     dd4hep::Position r22(x, -dy, dz);
0647 
0648     solid = dd4hep::Box(ns.prepend(name + "2"), dx, dy, dz);
0649 #ifdef EDM_ML_DEBUG
0650     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Box made of " << matter.name()
0651                                  << " of dimensions " << cms::convert2mm(dx) << ", " << cms::convert2mm(dy) << ", "
0652                                  << cms::convert2mm(dz);
0653 #endif
0654     glog = dd4hep::Volume(solid.name(), solid, matter);
0655 
0656     if (nAbs < 0) {
0657       mother = constructMidLayer(laylog, name, rin, alpha1, ns);
0658     } else {
0659       mother = laylog;
0660     }
0661     mother.placeVolume(glog, idOffset + 3, r21);
0662     mother.placeVolume(glog, idOffset + 4, dd4hep::Transform3D(rot, r22));
0663 #ifdef EDM_ML_DEBUG
0664     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number " << (idOffset + 3)
0665                                  << " positioned in " << mother.name() << " at (" << cms::convert2mm(x) << ","
0666                                  << cms::convert2mm(dy) << "," << cms::convert2mm(dz)
0667                                  << ") with no rotation\nDDHCalBarrelAlgo: " << glog.name() << " Number "
0668                                  << (idOffset + 4) << " positioned in " << mother.name() << " at ("
0669                                  << cms::convert2mm(x) << "," << -cms::convert2mm(dy) << "," << cms::convert2mm(dz)
0670                                  << ") with " << rot;
0671 #endif
0672     //Constructin the plastics and scintillators inside
0673     constructInsideDetectors(glog, nam0 + "2", id, dx, dy, dz, 2, ns);
0674   }
0675 
0676   dd4hep::Volume constructSideLayer(
0677       dd4hep::Volume& laylog, const std::string& nm, int nAbs, double rin, double alpha, cms::DDNamespace& ns) {
0678     //Extra absorber layer
0679     int k = abs(nAbs) - 1;
0680     std::string namek = nm + "Side";
0681     double rsi = rin + sideD[k];
0682     int in = 0;
0683     for (int i = 0; i < rzones - 1; i++) {
0684       if (rsi >= rmax[i])
0685         in = i + 1;
0686     }
0687     std::vector<double> pgonZ, pgonRmin, pgonRmax;
0688     // index 0
0689     pgonZ.emplace_back(0.0);
0690     pgonRmin.emplace_back(rsi);
0691     pgonRmax.emplace_back(rsi + sideT[k]);
0692     // index 1
0693     pgonZ.emplace_back(zoff[in] + rsi * ttheta[in]);
0694     pgonRmin.emplace_back(rsi);
0695     pgonRmax.emplace_back(pgonRmax[0]);
0696     // index 2
0697     pgonZ.emplace_back(zoff[in] + pgonRmax[0] * ttheta[in]);
0698     pgonRmin.emplace_back(pgonRmax[1]);
0699     pgonRmax.emplace_back(pgonRmax[1]);
0700     dd4hep::Solid solid = dd4hep::Polyhedra(ns.prepend(namek), 1, -alpha, 2 * alpha, pgonZ, pgonRmin, pgonRmax);
0701     dd4hep::Material matter = ns.material(sideMat[k]);
0702     dd4hep::Volume glog(solid.name(), solid, matter);
0703 #ifdef EDM_ML_DEBUG
0704     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << sideMat[k]
0705                                  << " with 1 sector from " << convertRadToDeg(-alpha) << " to "
0706                                  << convertRadToDeg(alpha) << " and with " << pgonZ.size() << " sections";
0707     for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0708       edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0709                                    << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0710                                    << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0711 #endif
0712 
0713     laylog.placeVolume(glog, 1);
0714 #ifdef EDM_ML_DEBUG
0715     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number 1 positioned in " << laylog.name()
0716                                  << " at (0,0,0) with no rotation";
0717 #endif
0718     if (nAbs < 0) {
0719       dd4hep::Volume mother = glog;
0720       double rmid = pgonRmax[0];
0721       for (int i = 0; i < nSideAbs; i++) {
0722         double alpha1 = atan(sideAbsW[i] / rmid);
0723         if (alpha1 > 0) {
0724           std::string name = namek + sideAbsName[i];
0725           solid = dd4hep::Polyhedra(ns.prepend(name), 1, -alpha1, 2 * alpha1, pgonZ, pgonRmin, pgonRmax);
0726           dd4hep::Material matter = ns.material(sideAbsMat[i]);
0727           dd4hep::Volume log(solid.name(), solid, matter);
0728 #ifdef EDM_ML_DEBUG
0729           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << sideAbsMat[i]
0730                                        << " with 1 sector from " << convertRadToDeg(-alpha1) << " to "
0731                                        << convertRadToDeg(alpha1) << " and with " << pgonZ.size() << " sections";
0732           for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0733             edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0734                                          << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0735                                          << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0736 #endif
0737 
0738           mother.placeVolume(log, 1);
0739 #ifdef EDM_ML_DEBUG
0740           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << log.name() << " Number 1 positioned in "
0741                                        << mother.name() << " at (0,0,0) with no rotation";
0742 #endif
0743           mother = log;
0744         }
0745       }
0746     }
0747     return glog;
0748   }
0749 
0750   dd4hep::Volume constructMidLayer(
0751       dd4hep::Volume laylog, const std::string& nm, double rin, double alpha, cms::DDNamespace& ns) {
0752     dd4hep::Solid solid;
0753     dd4hep::Volume log, glog;
0754     std::string name = nm + "Mid";
0755     for (int k = 0; k < nAbsorber; k++) {
0756       std::string namek = name + absorbName[k];
0757       double rsi = rin + absorbD[k];
0758       int in = 0;
0759       for (int i = 0; i < rzones - 1; i++) {
0760         if (rsi >= rmax[i])
0761           in = i + 1;
0762       }
0763       std::vector<double> pgonZ, pgonRmin, pgonRmax;
0764       // index 0
0765       pgonZ.emplace_back(0.0);
0766       pgonRmin.emplace_back(rsi);
0767       pgonRmax.emplace_back(rsi + absorbT[k]);
0768       // index 1
0769       pgonZ.emplace_back(zoff[in] + rsi * ttheta[in]);
0770       pgonRmin.emplace_back(rsi);
0771       pgonRmax.emplace_back(pgonRmax[0]);
0772       // index 2
0773       pgonZ.emplace_back(zoff[in] + pgonRmax[0] * ttheta[in]);
0774       pgonRmin.emplace_back(pgonRmax[1]);
0775       pgonRmax.emplace_back(pgonRmax[1]);
0776       solid = dd4hep::Polyhedra(ns.prepend(namek), 1, -alpha, 2 * alpha, pgonZ, pgonRmin, pgonRmax);
0777       dd4hep::Material matter = ns.material(absorbMat[k]);
0778       log = dd4hep::Volume(solid.name(), solid, matter);
0779 #ifdef EDM_ML_DEBUG
0780       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << matter.name()
0781                                    << " with 1 sector from " << convertRadToDeg(-alpha) << " to "
0782                                    << convertRadToDeg(alpha) << " and with " << pgonZ.size() << " sections";
0783       for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0784         edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0785                                      << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0786                                      << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0787 #endif
0788 
0789       laylog.placeVolume(log, 1);
0790 #ifdef EDM_ML_DEBUG
0791       edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << log.name() << " Number 1 positioned in " << laylog.name()
0792                                    << " at (0,0,0) with no rotation";
0793 #endif
0794       if (k == 0) {
0795         double rmin = pgonRmin[0];
0796         double rmax = pgonRmax[0];
0797         dd4hep::Volume mother = log;
0798         for (int i = 0; i < 1; i++) {
0799           double alpha1 = atan(midW[i] / rmin);
0800           std::string namek = name + midName[i];
0801           solid = dd4hep::Polyhedra(ns.prepend(namek), 1, -alpha1, 2 * alpha1, pgonZ, pgonRmin, pgonRmax);
0802           dd4hep::Material matter1 = ns.material(midMat[i]);
0803           log = dd4hep::Volume(solid.name(), solid, matter1);
0804 #ifdef EDM_ML_DEBUG
0805           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of "
0806                                        << matter1.name() << " with 1 sector from " << convertRadToDeg(-alpha1) << " to "
0807                                        << convertRadToDeg(alpha1) << " and with " << pgonZ.size() << " sections";
0808           for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0809             edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0810                                          << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0811                                          << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0812 #endif
0813 
0814           mother.placeVolume(log, 1);
0815 #ifdef EDM_ML_DEBUG
0816           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << log.name() << " Number 1 positioned in "
0817                                        << mother.name() << " at (0,0,0) with no rotation";
0818 #endif
0819           mother = log;
0820         }
0821 
0822         // Now the layer with detectors
0823         double rmid = rmin + middleD;
0824         pgonRmin[0] = rmid;
0825         pgonRmax[0] = rmax;
0826         pgonRmin[1] = rmid;
0827         pgonRmax[1] = rmax;
0828         pgonZ[1] = zoff[in] + rmid * ttheta[in];
0829         pgonRmin[2] = rmax;
0830         pgonRmax[2] = rmax;
0831         pgonZ[2] = zoff[in] + rmax * ttheta[in];
0832         double alpha1 = atan(middleW / rmin);
0833         solid = dd4hep::Polyhedra(ns.prepend(name), 1, -alpha1, 2 * alpha1, pgonZ, pgonRmin, pgonRmax);
0834         dd4hep::Material matter1 = ns.material(middleMat);
0835         glog = dd4hep::Volume(solid.name(), solid, matter1);
0836 #ifdef EDM_ML_DEBUG
0837         edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of " << matter1.name()
0838                                      << " with 1 sector from " << convertRadToDeg(-alpha1) << " to "
0839                                      << convertRadToDeg(alpha1) << " and with " << pgonZ.size() << " sections";
0840         for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0841           edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0842                                        << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0843                                        << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0844 #endif
0845 
0846         mother.placeVolume(glog, 1);
0847 #ifdef EDM_ML_DEBUG
0848         edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number 1 positioned in "
0849                                      << mother.name() << " at (0,0,0) with no rotation";
0850 #endif
0851         // Now the remaining absorber layers
0852         for (int i = 1; i < nMidAbs; i++) {
0853           namek = name + midName[i];
0854           rmid = rmin + midT[i];
0855           pgonRmin[0] = rmin;
0856           pgonRmax[0] = rmid;
0857           pgonRmin[1] = rmin;
0858           pgonRmax[1] = rmid;
0859           pgonZ[1] = zoff[in] + rmin * ttheta[in];
0860           pgonRmin[2] = rmid;
0861           pgonRmax[2] = rmid;
0862           pgonZ[2] = zoff[in] + rmid * ttheta[in];
0863           alpha1 = atan(midW[i] / rmin);
0864           solid = dd4hep::Polyhedra(ns.prepend(namek), 1, -alpha1, 2 * alpha1, pgonZ, pgonRmin, pgonRmax);
0865           dd4hep::Material matter2 = ns.material(midMat[i]);
0866           log = dd4hep::Volume(solid.name(), solid, matter2);
0867 #ifdef EDM_ML_DEBUG
0868           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Polyhedra made of "
0869                                        << matter2.name() << " with 1 sector from " << convertRadToDeg(-alpha1) << " to "
0870                                        << convertRadToDeg(alpha1) << " and with " << pgonZ.size() << " sections";
0871           for (unsigned int ii = 0; ii < pgonZ.size(); ii++)
0872             edm::LogVerbatim("HCalGeom") << "\t\tZ = " << cms::convert2mm(pgonZ[ii])
0873                                          << "\tRmin = " << cms::convert2mm(pgonRmin[ii])
0874                                          << "\tRmax = " << cms::convert2mm(pgonRmax[ii]);
0875 #endif
0876 
0877           mother.placeVolume(log, i);
0878 #ifdef EDM_ML_DEBUG
0879           edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << log.name() << " Number " << i << " positioned in "
0880                                        << mother.name() << " at (0,0,0) with no rotation";
0881 #endif
0882           mother = log;
0883         }
0884       }
0885     }
0886     return glog;
0887   }
0888 
0889   void constructInsideDetectors(dd4hep::Volume& detector,
0890                                 const std::string& name,
0891                                 int id,
0892                                 double dx,
0893                                 double dy,
0894                                 double dz,
0895                                 int type,
0896                                 cms::DDNamespace& ns) {
0897 #ifdef EDM_ML_DEBUG
0898     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: \t\tInside detector " << id << "...";
0899 #endif
0900 
0901     dd4hep::Material plmatter = ns.material(detMatPl);
0902     dd4hep::Material scmatter = ns.material(detMatSc);
0903     std::string plname = DDSplit(detector.name()).first + "Plastic_";
0904     std::string scname = idName + "Scintillator" + name;
0905 
0906     id--;
0907     dd4hep::Solid solid;
0908     dd4hep::Volume glog;
0909     double wid, y = 0;
0910     double dx1, dx2, shiftX;
0911 
0912     if (type == 1) {
0913       wid = 0.5 * detWidth1[id];
0914       dx1 = 0.5 * detT11[id];
0915       dx2 = 0.5 * detT21[id];
0916       shiftX = detdP1[id];
0917       if (detPosY[id] > 0)
0918         y = -dy + wid;
0919     } else {
0920       wid = 0.5 * detWidth2[id];
0921       dx1 = 0.5 * detT12[id];
0922       dx2 = 0.5 * detT22[id];
0923       shiftX = detdP2[id];
0924     }
0925 
0926     solid = dd4hep::Box(ns.prepend(plname + "1"), dx1, wid, dz);
0927     glog = dd4hep::Volume(solid.name(), solid, plmatter);
0928 #ifdef EDM_ML_DEBUG
0929     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Box made of " << plmatter.name()
0930                                  << " of dimensions " << cms::convert2mm(dx1) << ", " << cms::convert2mm(wid) << ", "
0931                                  << cms::convert2mm(dz);
0932 #endif
0933 
0934     double x = shiftX + dx1 - dx;
0935     detector.placeVolume(glog, 1, dd4hep::Position(x, y, 0));
0936 #ifdef EDM_ML_DEBUG
0937     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number 1 positioned in " << detector.name()
0938                                  << " at (" << cms::convert2mm(x) << "," << cms::convert2mm(y)
0939                                  << ",0) with no rotation";
0940 #endif
0941     solid = dd4hep::Box(ns.prepend(scname), 0.5 * detTsc[id], wid, dz);
0942     glog = dd4hep::Volume(solid.name(), solid, scmatter);
0943 #ifdef EDM_ML_DEBUG
0944     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Box made of " << scmatter.name()
0945                                  << " of dimensions " << cms::convert2mm(0.5 * detTsc[id]) << ", "
0946                                  << cms::convert2mm(wid) << ", " << cms::convert2mm(dz);
0947 #endif
0948 
0949     x += dx1 + 0.5 * detTsc[id];
0950     int copyNo = id * 10 + detType[id];
0951     detector.placeVolume(glog, copyNo, dd4hep::Position(x, y, 0));
0952 #ifdef EDM_ML_DEBUG
0953     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number " << copyNo << " positioned in "
0954                                  << detector.name() << " at (" << cms::convert2mm(x) << "," << cms::convert2mm(y)
0955                                  << ",0) with no rotation";
0956 #endif
0957     solid = dd4hep::Box(ns.prepend(plname + "2"), dx2, wid, dz);
0958     glog = dd4hep::Volume(solid.name(), solid, plmatter);
0959 #ifdef EDM_ML_DEBUG
0960     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << solid.name() << " Box made of " << plmatter.name()
0961                                  << " of dimensions " << cms::convert2mm(dx2) << ", " << cms::convert2mm(wid) << ", "
0962                                  << cms::convert2mm(dz);
0963 #endif
0964 
0965     x += 0.5 * detTsc[id] + dx2;
0966     detector.placeVolume(glog, 1, dd4hep::Position(x, y, 0));
0967 #ifdef EDM_ML_DEBUG
0968     edm::LogVerbatim("HCalGeom") << "DDHCalBarrelAlgo: " << glog.name() << " Number 1 positioned in " << detector.name()
0969                                  << " at (" << cms::convert2mm(x) << "," << cms::convert2mm(y)
0970                                  << ",0) with no rotation";
0971 #endif
0972   }
0973 
0974   dd4hep::Rotation3D getRotation(std::string& rotation, std::string& rotns, cms::DDNamespace& ns) {
0975     std::string rot = (strchr(rotation.c_str(), NAMESPACE_SEP) == nullptr) ? (rotns + ":" + rotation) : rotation;
0976 #ifdef EDM_ML_DEBUG
0977     edm::LogVerbatim("HCalGeom") << "getRotation: " << rotation << ":" << rot << ":" << ns.rotation(rot);
0978 #endif
0979     return ns.rotation(rot);
0980   }
0981 };
0982 
0983 static long algorithm(dd4hep::Detector& /* description */, cms::DDParsingContext& ctxt, xml_h e) {
0984   HcalBarrelAlgo hcalbarrelalgo(ctxt, e);
0985 #ifdef EDM_ML_DEBUG
0986   edm::LogVerbatim("HCalGeom") << "<<== End of DDHCalBarrelAlgo construction";
0987 #endif
0988   return cms::s_executed;
0989 }
0990 
0991 // first argument is the type from the xml file
0992 DECLARE_DDCMS_DETELEMENT(DDCMS_hcal_DDHCalBarrelAlgo, algorithm)