Back to home page

EIC code displayed by LXR

 
 

    


Warning, file /DD4hep/examples/DDCMS/src/plugins/DDEcalEndcapAlgo.cpp was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 #include "DD4hep/DetFactoryHelper.h"
0002 #include "DDCMS/DDCMSPlugins.h"
0003 #include <Math/AxisAngle.h>
0004 #include <Math/Rotation3D.h>
0005 #include <Math/Vector3D.h>
0006 
0007 #include <CLHEP/Geometry/Transform3D.h>
0008 #include <CLHEP/Units/PhysicalConstants.h>
0009 #include <CLHEP/Units/SystemOfUnits.h>
0010 #include <CLHEP/Geometry/Point3D.h>
0011 #include <CLHEP/Geometry/Plane3D.h>
0012 #include <CLHEP/Geometry/Vector3D.h>
0013 #include <CLHEP/Geometry/Transform3D.h>
0014 #include <CLHEP/Vector/EulerAngles.h>
0015 
0016 #include <array>
0017 #include <string>
0018 #include <vector>
0019 
0020 using namespace std;
0021 using namespace dd4hep;
0022 
0023 using DD3Vector = ROOT::Math::DisplacementVector3D<ROOT::Math::Cartesian3D<double>>;
0024 using DDTranslation = ROOT::Math::DisplacementVector3D<ROOT::Math::Cartesian3D<double> >;
0025 using DDRotation = ROOT::Math::Rotation3D;
0026 using DDRotationMatrix = ROOT::Math::Rotation3D;
0027 using DDAxisAngle = ROOT::Math::AxisAngle;
0028 
0029 namespace {
0030 
0031   constexpr long double piRadians(M_PI);
0032   constexpr long double degPerRad = 180. / piRadians;  // Degrees per radian
0033   constexpr double operator"" _mm(long double length) { return length * 0.1; }
0034   constexpr long double operator"" _deg(long double value) { return value / degPerRad; }
0035 
0036   // Define Endcap Supercrystal class
0037   class DDEcalEndcapTrap {
0038   public:
0039     DDEcalEndcapTrap(const int hand, const double front, const double rear, const double length);
0040     DDEcalEndcapTrap() = delete;
0041 
0042     void rotate(const DDRotationMatrix& rot);
0043     //void rotate(const DDTranslation& frontCentre, const DDTranslation& rearCentre);
0044     void translate(const DDTranslation& trans);
0045 
0046     void rotateX(const double angle);
0047     void rotateY(const double angle);
0048     void translate();
0049     void moveto(const DDTranslation& frontCentre, const DDTranslation& rearCentre);
0050     double elevationAngle(const DDTranslation& trans);
0051     double polarAngle(const DDTranslation& trans);
0052     double elevationAngle();
0053     double polarAngle();
0054     DDTranslation cornerPos(const int icorner);
0055     void cornerPos(const int icorner, const DDTranslation& cc);
0056     DDTranslation centrePos();
0057     DDTranslation fcentrePos();
0058     DDTranslation rcentrePos();
0059     void calculateCorners();
0060     void calculateCentres();
0061     DDRotationMatrix rotation() { return m_rotation; }
0062     //void print();
0063 
0064   private:
0065     DDRotationMatrix m_rotation;
0066     DDTranslation m_translation;
0067 
0068     double m_centre[4];
0069     double m_fcentre[4];
0070     double m_rcentre[4];
0071     double m_corners[25];
0072     double m_front;
0073     double m_rear;
0074     double m_length;
0075 
0076     int m_hand;
0077     int m_update;
0078   };
0079   // Implementation of DDEcalEndcapTrap class
0080 
0081   DDEcalEndcapTrap::DDEcalEndcapTrap(const int hand, const double front, const double rear, const double length) {
0082     //
0083     //  Initialise corners of supercrystal.
0084 
0085     // Start out with bottom surface on (x,z) plane, front face in (x,y) plane.
0086 
0087     double xsign;
0088 
0089     if (hand == 2) {
0090       xsign = -1.;
0091     } else {
0092       xsign = 1.;
0093     }
0094 
0095     m_hand = hand;
0096     m_front = front;
0097     m_rear = rear;
0098     m_length = length;
0099 
0100     int icorner;
0101     icorner = 1;
0102     m_corners[3 * icorner - 3] = xsign * front;
0103     m_corners[3 * icorner - 2] = front;
0104     m_corners[3 * icorner - 1] = 0.;
0105     icorner = 2;
0106     m_corners[3 * icorner - 3] = xsign * front;
0107     m_corners[3 * icorner - 2] = 0.;
0108     m_corners[3 * icorner - 1] = 0.;
0109     icorner = 3;
0110     m_corners[3 * icorner - 3] = 0.;
0111     m_corners[3 * icorner - 2] = 0.;
0112     m_corners[3 * icorner - 1] = 0.;
0113     icorner = 4;
0114     m_corners[3 * icorner - 3] = 0.;
0115     m_corners[3 * icorner - 2] = front;
0116     m_corners[3 * icorner - 1] = 0.;
0117 
0118     icorner = 5;
0119     m_corners[3 * icorner - 3] = xsign * rear;
0120     m_corners[3 * icorner - 2] = rear;
0121     m_corners[3 * icorner - 1] = length;
0122     icorner = 6;
0123     m_corners[3 * icorner - 3] = xsign * rear;
0124     m_corners[3 * icorner - 2] = 0.;
0125     m_corners[3 * icorner - 1] = length;
0126     icorner = 7;
0127     m_corners[3 * icorner - 3] = 0.;
0128     m_corners[3 * icorner - 2] = 0.;
0129     m_corners[3 * icorner - 1] = length;
0130     icorner = 8;
0131     m_corners[3 * icorner - 3] = 0.;
0132     m_corners[3 * icorner - 2] = rear;
0133     m_corners[3 * icorner - 1] = length;
0134 
0135     calculateCentres();
0136 
0137     // Move centre of SC to (0,0,0)
0138 
0139     translate();
0140 
0141     // Rotate into standard position (face centres on z axis)
0142 
0143     //  this->rotate();
0144 
0145     calculateCentres();
0146   }
0147 #if 0
0148   void DDEcalEndcapTrap::rotate(const DDTranslation& /* frontCentre */, const DDTranslation& /* rearCentre */) {
0149     //
0150     //  Rotate supercrystal to bring front and rear face centres to specified points
0151     //
0152   }
0153 #endif
0154   void DDEcalEndcapTrap::rotate(const DDRotationMatrix& rot) {
0155     //
0156     //  Rotate supercrystal by specified rotation about (0,0,0)
0157     //
0158 
0159     int icorner;
0160     DDTranslation cc;
0161     for (icorner = 1; icorner <= 8; icorner++) {
0162       cc = cornerPos(icorner);
0163       cc = rot * cc;
0164       cornerPos(icorner, cc);
0165     }
0166     m_rotation = rot * m_rotation;
0167     calculateCentres();
0168   }
0169 
0170   void DDEcalEndcapTrap::translate() {
0171     translate(-1. * centrePos());
0172   }
0173 
0174   void DDEcalEndcapTrap::translate(const DDTranslation& trans) {
0175     //
0176     //  Translate supercrystal by specified amount
0177     //
0178 
0179     DDTranslation tcorner;
0180     for (int icorner = 1; icorner <= 8; icorner++) {
0181       tcorner = cornerPos(icorner) + trans;
0182       cornerPos(icorner, tcorner);
0183     }
0184     calculateCentres();
0185     m_translation = trans + m_translation;
0186   }
0187 
0188   void DDEcalEndcapTrap::moveto(const DDTranslation& frontCentre, const DDTranslation& rearCentre) {
0189     //
0190     //  Rotate (about X then about Y) and translate supercrystal to bring axis joining front and rear face centres parallel to line connecting specified points
0191     //
0192 
0193     //  Get azimuthal and polar angles of current axis and target axis
0194     double currentTheta = elevationAngle();
0195     double currentPhi = polarAngle();
0196     double targetTheta = elevationAngle(frontCentre - rearCentre);
0197     double targetPhi = polarAngle(frontCentre - rearCentre);
0198 
0199     //  Rotate to correct angle (X then Y)
0200     rotateX(targetTheta - currentTheta);
0201     rotateY(targetPhi - currentPhi);
0202 
0203     //  Translate SC to final position
0204     DDTranslation targetCentre = 0.5 * (frontCentre + rearCentre);
0205     translate(targetCentre - centrePos());
0206   }
0207 
0208   void DDEcalEndcapTrap::rotateX(const double angle) {
0209     //
0210     //  Rotate SC through given angle about X axis
0211     //
0212 
0213     const CLHEP::HepRotation tmp(CLHEP::Hep3Vector(1., 0., 0.), angle);
0214 
0215     rotate(DDRotationMatrix(tmp.xx(), tmp.xy(), tmp.xz(), tmp.yx(), tmp.yy(), tmp.yz(), tmp.zx(), tmp.zy(), tmp.zz()));
0216   }
0217 
0218   void DDEcalEndcapTrap::rotateY(const double angle) {
0219     //
0220     //  Rotate SC through given angle about Y axis
0221     //
0222     const CLHEP::HepRotation tmp(CLHEP::Hep3Vector(0., 1., 0.), angle);
0223 
0224     rotate(DDRotationMatrix(tmp.xx(), tmp.xy(), tmp.xz(), tmp.yx(), tmp.yy(), tmp.yz(), tmp.zx(), tmp.zy(), tmp.zz()));
0225   }
0226 
0227   void DDEcalEndcapTrap::calculateCentres() {
0228     //
0229     //  Calculate crystal centre and front & rear face centres
0230     //
0231 
0232     int ixyz, icorner;
0233 
0234     for (ixyz = 0; ixyz < 3; ixyz++) {
0235       m_centre[ixyz] = 0;
0236       m_fcentre[ixyz] = 0;
0237       m_rcentre[ixyz] = 0;
0238     }
0239 
0240     for (icorner = 1; icorner <= 4; icorner++) {
0241       for (ixyz = 0; ixyz < 3; ixyz++) {
0242     m_centre[ixyz] = m_centre[ixyz] + 0.125 * m_corners[3 * icorner - 3 + ixyz];
0243     m_fcentre[ixyz] = m_fcentre[ixyz] + 0.25 * m_corners[3 * icorner - 3 + ixyz];
0244       }
0245     }
0246     for (icorner = 5; icorner <= 8; icorner++) {
0247       for (ixyz = 0; ixyz < 3; ixyz++) {
0248     m_centre[ixyz] = m_centre[ixyz] + 0.125 * m_corners[3 * icorner - 3 + ixyz];
0249     m_rcentre[ixyz] = m_rcentre[ixyz] + 0.25 * m_corners[3 * icorner - 3 + ixyz];
0250       }
0251     }
0252   }
0253 
0254   DDTranslation DDEcalEndcapTrap::cornerPos(const int icorner) {
0255     //
0256     //  Return specified corner as a DDTranslation
0257     //
0258     return DDTranslation(m_corners[3 * icorner - 3], m_corners[3 * icorner - 2], m_corners[3 * icorner - 1]);
0259   }
0260 
0261   void DDEcalEndcapTrap::cornerPos(const int icorner, const DDTranslation& cornerxyz) {
0262     //
0263     //  Save position of specified corner.
0264     //
0265     for (int ixyz = 0; ixyz < 3; ixyz++) {
0266       m_corners[3 * icorner - 3 + ixyz] = (0 == ixyz ? cornerxyz.x() : (1 == ixyz ? cornerxyz.y() : cornerxyz.z()));
0267       ;
0268     }
0269   }
0270 
0271   DDTranslation DDEcalEndcapTrap::centrePos() {
0272     //
0273     //  Return SC centre as a DDTranslation
0274     //
0275     return DDTranslation(m_centre[0], m_centre[1], m_centre[2]);
0276   }
0277 
0278   DDTranslation DDEcalEndcapTrap::fcentrePos() {
0279     //
0280     //  Return SC front face centre as a DDTranslation
0281     //
0282     return DDTranslation(m_fcentre[0], m_fcentre[1], m_fcentre[2]);
0283   }
0284 
0285   DDTranslation DDEcalEndcapTrap::rcentrePos() {
0286     //
0287     //  Return SC rear face centre as a DDTranslation
0288     //
0289     return DDTranslation(m_rcentre[0], m_rcentre[1], m_rcentre[2]);
0290   }
0291 
0292   double DDEcalEndcapTrap::elevationAngle(const DDTranslation& trans) {
0293     //
0294     //  Return elevation angle (out of x-z plane) of a given translation (seen as a vector from the origin).
0295     //
0296     double sintheta = trans.y() / trans.r();
0297     return asin(sintheta);
0298   }
0299 
0300   double DDEcalEndcapTrap::elevationAngle() {
0301     //
0302     //  Return elevation angle (out of x-z plane) of SC in current position.
0303     //
0304     DDTranslation current = fcentrePos() - rcentrePos();
0305     return elevationAngle(current);
0306   }
0307 
0308   double DDEcalEndcapTrap::polarAngle(const DDTranslation& trans) {
0309     //
0310     //  Return polar angle (from x to z) of a given translation (seen as a vector from the origin).
0311     //
0312     double tanphi = trans.x() / trans.z();
0313     return atan(tanphi);
0314   }
0315 
0316   double DDEcalEndcapTrap::polarAngle() {
0317     //
0318     //  Return elevation angle (out of x-z plane) of SC in current position.
0319     //
0320     DDTranslation current = fcentrePos() - rcentrePos();
0321     return polarAngle(current);
0322   }
0323 #if 0
0324   void DDEcalEndcapTrap::print() {
0325     //
0326     //  Print SC coordinates for debugging
0327     //
0328     for (int ic = 1; ic <= 8; ic++) {
0329       /* DDTranslation cc = */  cornerPos(ic);
0330     }
0331   }
0332 #endif
0333   namespace {
0334     struct Endcap {
0335       string mat;
0336       double zOff;
0337 
0338       string quaName;
0339       string quaMat;
0340 
0341       string crysMat;
0342       string wallMat;
0343 
0344       double crysLength;
0345       double crysRear;
0346       double crysFront;
0347       double sCELength;
0348       double sCERear;
0349       double sCEFront;
0350       double sCALength;
0351       double sCARear;
0352       double sCAFront;
0353       double sCAWall;
0354       double sCHLength;
0355       double sCHSide;
0356 
0357       double nSCTypes;
0358       vector<double> vecEESCProf;
0359       double nColumns;
0360       vector<double> vecEEShape;
0361       double nSCCutaway;
0362       vector<double> vecEESCCutaway;
0363       double nSCquad;
0364       vector<double> vecEESCCtrs;
0365       double nCRSC;
0366       vector<double> vecEECRCtrs;
0367 
0368       array<double, 3> cutParms;
0369       string cutBoxName;
0370 
0371       string envName;
0372       string alvName;
0373       string intName;
0374       string cryName;
0375 
0376       DDTranslation cryFCtr[5][5];
0377       DDTranslation cryRCtr[5][5];
0378       DDTranslation scrFCtr[10][10];
0379       DDTranslation scrRCtr[10][10];
0380 
0381       double pFHalf;
0382       double pFFifth;
0383       double pF45;
0384 
0385       vector<double> vecEESCLims;
0386 
0387       double iLength;
0388       double iXYOff;
0389       double cryZOff;
0390       double zFront;
0391     };
0392 
0393     const Rotation3D& myrot(dd4hep::cms::Namespace& ns, const string& nam, const Rotation3D& r) {
0394       ns.addRotation(nam, r);
0395       return ns.rotation(ns.prepend(nam));
0396     }
0397 
0398     string_view mynamespace(string_view input) {
0399       string_view v = input;
0400       auto trim_pos = v.find(':');
0401       if (trim_pos != v.npos)
0402     v.remove_suffix(v.size() - (trim_pos + 1));
0403       return v;
0404     }
0405   }  // namespace
0406 
0407   static long algorithm(dd4hep::Detector& /* description */, dd4hep::cms::ParsingContext& ctxt, xml_h e,
0408             SensitiveDetector& /* sens */) {
0409     dd4hep::cms::Namespace ns(ctxt, e, true);
0410     dd4hep::cms::AlgoArguments args(ctxt, e);
0411 
0412     // TRICK!
0413     string myns{mynamespace(args.parentName()).data(), mynamespace(args.parentName()).size()};
0414 
0415     Endcap ee;
0416     ee.mat = args.str("EEMat");
0417     ee.zOff = args.dble("EEzOff");
0418 
0419     ee.quaName = args.str("EEQuaName");
0420     ee.quaMat = args.str("EEQuaMat");
0421     ee.crysMat = args.str("EECrysMat");
0422     ee.wallMat = args.str("EEWallMat");
0423     ee.crysLength = args.dble("EECrysLength");
0424     ee.crysRear = args.dble("EECrysRear");
0425     ee.crysFront = args.dble("EECrysFront");
0426     ee.sCELength = args.dble("EESCELength");
0427     ee.sCERear = args.dble("EESCERear");
0428     ee.sCEFront = args.dble("EESCEFront");
0429     ee.sCALength = args.dble("EESCALength");
0430     ee.sCARear = args.dble("EESCARear");
0431     ee.sCAFront = args.dble("EESCAFront");
0432     ee.sCAWall = args.dble("EESCAWall");
0433     ee.sCHLength = args.dble("EESCHLength");
0434     ee.sCHSide = args.dble("EESCHSide");
0435     ee.nSCTypes = args.dble("EEnSCTypes");
0436     ee.nColumns = args.dble("EEnColumns");
0437     ee.nSCCutaway = args.dble("EEnSCCutaway");
0438     ee.nSCquad = args.dble("EEnSCquad");
0439     ee.nCRSC = args.dble("EEnCRSC");
0440     ee.vecEESCProf = args.vecDble("EESCProf");
0441     ee.vecEEShape = args.vecDble("EEShape");
0442     ee.vecEESCCutaway = args.vecDble("EESCCutaway");
0443     ee.vecEESCCtrs = args.vecDble("EESCCtrs");
0444     ee.vecEECRCtrs = args.vecDble("EECRCtrs");
0445 
0446     ee.cutBoxName = args.str("EECutBoxName");
0447 
0448     ee.envName = args.str("EEEnvName");
0449     ee.alvName = args.str("EEAlvName");
0450     ee.intName = args.str("EEIntName");
0451     ee.cryName = args.str("EECryName");
0452 
0453     ee.pFHalf = args.dble("EEPFHalf");
0454     ee.pFFifth = args.dble("EEPFFifth");
0455     ee.pF45 = args.dble("EEPF45");
0456 
0457     ee.vecEESCLims = args.vecDble("EESCLims");
0458     ee.iLength = args.dble("EEiLength");
0459     ee.iXYOff = args.dble("EEiXYOff");
0460     ee.cryZOff = args.dble("EECryZOff");
0461     ee.zFront = args.dble("EEzFront");
0462 
0463     //  Position supercrystals in EE Quadrant
0464 
0465     //********************************* cutbox for trimming edge SCs
0466     const double cutWid(ee.sCERear / sqrt(2.));
0467     ee.cutParms[0] = cutWid;
0468     ee.cutParms[1] = cutWid;
0469     ee.cutParms[2] = ee.sCELength / sqrt(2.);
0470     Solid eeCutBox = Box(ee.cutBoxName, ee.cutParms[0], ee.cutParms[1], ee.cutParms[2]);
0471     //**************************************************************
0472 
0473     const double zFix(ee.zFront - 3172.0_mm);  // fix for changing z offset
0474 
0475     //** fill supercrystal front and rear center positions from xml input
0476     for (unsigned int iC(0); iC != (unsigned int)ee.nSCquad; ++iC) {
0477       const unsigned int iOff(8 * iC);
0478       const unsigned int ix((unsigned int)ee.vecEESCCtrs[iOff + 0]);
0479       const unsigned int iy((unsigned int)ee.vecEESCCtrs[iOff + 1]);
0480 
0481       assert(ix > 0 && ix < 11 && iy > 0 && iy < 11);
0482 
0483       ee.scrFCtr[ix - 1][iy - 1] =
0484         DDTranslation(ee.vecEESCCtrs[iOff + 2], ee.vecEESCCtrs[iOff + 4], ee.vecEESCCtrs[iOff + 6] + zFix);
0485 
0486       ee.scrRCtr[ix - 1][iy - 1] =
0487         DDTranslation(ee.vecEESCCtrs[iOff + 3], ee.vecEESCCtrs[iOff + 5], ee.vecEESCCtrs[iOff + 7] + zFix);
0488     }
0489 
0490     //** fill crystal front and rear center positions from xml input
0491     for (unsigned int iC(0); iC != 25; ++iC) {
0492       const unsigned int iOff(8 * iC);
0493       const unsigned int ix((unsigned int)ee.vecEECRCtrs[iOff + 0]);
0494       const unsigned int iy((unsigned int)ee.vecEECRCtrs[iOff + 1]);
0495 
0496       assert(ix > 0 && ix < 6 && iy > 0 && iy < 6);
0497 
0498       ee.cryFCtr[ix - 1][iy - 1] =
0499         DDTranslation(ee.vecEECRCtrs[iOff + 2], ee.vecEECRCtrs[iOff + 4], ee.vecEECRCtrs[iOff + 6]);
0500 
0501       ee.cryRCtr[ix - 1][iy - 1] =
0502         DDTranslation(ee.vecEECRCtrs[iOff + 3], ee.vecEECRCtrs[iOff + 5], ee.vecEECRCtrs[iOff + 7]);
0503     }
0504 
0505     Solid eeCRSolid = Trap(ee.cryName,
0506                0.5 * ee.crysLength,
0507                atan((ee.crysRear - ee.crysFront) / (sqrt(2.) * ee.crysLength)),
0508                45._deg,
0509                0.5 * ee.crysFront,
0510                0.5 * ee.crysFront,
0511                0.5 * ee.crysFront,
0512                0._deg,
0513                0.5 * ee.crysRear,
0514                0.5 * ee.crysRear,
0515                0.5 * ee.crysRear,
0516                0._deg);
0517     Volume eeCRLog = Volume(ee.cryName, eeCRSolid, ns.material(ee.crysMat));
0518 
0519     for (unsigned int isc(0); isc < ee.nSCTypes; ++isc) {
0520       unsigned int iSCType = isc + 1;
0521       const string anum(std::to_string(iSCType));
0522       const double eFront(0.5 * ee.sCEFront);
0523       const double eRear(0.5 * ee.sCERear);
0524       const double eAng(atan((ee.sCERear - ee.sCEFront) / (sqrt(2.) * ee.sCELength)));
0525       const double ffived(45._deg);
0526       const double zerod(0._deg);
0527       string eeSCEnvName(1 == iSCType ? ee.envName + std::to_string(iSCType)
0528              : (ee.envName + std::to_string(iSCType) + "Tmp"));
0529       Solid eeSCEnv = ns.addSolidNS(
0530                     eeSCEnvName,
0531                     Trap(eeSCEnvName, 0.5 * ee.sCELength, eAng, ffived, eFront, eFront, eFront, zerod, eRear, eRear, eRear, zerod));
0532 
0533       const double aFront(0.5 * ee.sCAFront);
0534       const double aRear(0.5 * ee.sCARear);
0535       const double aAng(atan((ee.sCARear - ee.sCAFront) / (sqrt(2.) * ee.sCALength)));
0536       string eeSCAlvName(
0537              (1 == iSCType ? ee.alvName + std::to_string(iSCType) : (ee.alvName + std::to_string(iSCType) + "Tmp")));
0538       Solid eeSCAlv = ns.addSolidNS(
0539                     eeSCAlvName,
0540                     Trap(eeSCAlvName, 0.5 * ee.sCALength, aAng, ffived, aFront, aFront, aFront, zerod, aRear, aRear, aRear, zerod));
0541 
0542       const double dwall(ee.sCAWall);
0543       const double iFront(aFront - dwall);
0544       const double iRear(iFront);
0545       const double iLen(ee.iLength);
0546       string eeSCIntName(1 == iSCType ? ee.intName + std::to_string(iSCType)
0547              : (ee.intName + std::to_string(iSCType) + "Tmp"));
0548       Solid eeSCInt = ns.addSolidNS(eeSCIntName,
0549                     Trap(eeSCIntName,
0550                      iLen / 2.,
0551                      atan((ee.sCARear - ee.sCAFront) / (sqrt(2.) * ee.sCALength)),
0552                      ffived,
0553                      iFront,
0554                      iFront,
0555                      iFront,
0556                      zerod,
0557                      iRear,
0558                      iRear,
0559                      iRear,
0560                      zerod));
0561 
0562       const double dz(-0.5 * (ee.sCELength - ee.sCALength));
0563       const double dxy(0.5 * dz * (ee.sCERear - ee.sCEFront) / ee.sCELength);
0564       const double zIOff(-(ee.sCALength - iLen) / 2.);
0565       const double xyIOff(ee.iXYOff);
0566 
0567       Volume eeSCELog;
0568       Volume eeSCALog;
0569       Volume eeSCILog;
0570 
0571       if (1 == iSCType) {  // standard SC in this block
0572     eeSCELog = ns.addVolumeNS(Volume(myns + ee.envName + std::to_string(iSCType), eeSCEnv, ns.material(ee.mat)));
0573     eeSCALog = Volume(ee.alvName + std::to_string(iSCType), eeSCAlv, ns.material(ee.wallMat));
0574     eeSCILog = Volume(ee.intName + std::to_string(iSCType), eeSCInt, ns.material(ee.mat));
0575       } else {  // partial SCs this block: create subtraction volumes as appropriate
0576     const double half(ee.cutParms[0] - ee.pFHalf * ee.crysRear);
0577     const double fifth(ee.cutParms[0] + ee.pFFifth * ee.crysRear);
0578     const double fac(ee.pF45);
0579 
0580     const double zmm(0.0_mm);
0581 
0582     DDTranslation cutTra(
0583                  2 == iSCType ? DDTranslation(zmm, half, zmm)
0584                  : (3 == iSCType ? DDTranslation(half, zmm, zmm)
0585                 : (4 == iSCType ? DDTranslation(zmm, -fifth, zmm)
0586                    : (5 == iSCType ? DDTranslation(-half * fac, -half * fac, zmm)
0587                       : DDTranslation(-fifth, zmm, zmm)))));
0588 
0589     const CLHEP::HepRotationZ cutm(ffived);
0590 
0591     Rotation3D cutRot(5 != iSCType ? Rotation3D()
0592               : myrot(ns,
0593                   "EECry5Rot",
0594                   Rotation3D(cutm.xx(),
0595                          cutm.xy(),
0596                          cutm.xz(),
0597                          cutm.yx(),
0598                          cutm.yy(),
0599                          cutm.yz(),
0600                          cutm.zx(),
0601                          cutm.zy(),
0602                          cutm.zz())));
0603 
0604     Solid eeCutEnv = SubtractionSolid(ee.envName + std::to_string(iSCType),
0605                       ns.solid(ee.envName + std::to_string(iSCType) + "Tmp"),
0606                       eeCutBox,
0607                       Transform3D(cutRot, cutTra));
0608 
0609     const DDTranslation extra(dxy, dxy, dz);
0610 
0611     Solid eeCutAlv = SubtractionSolid(ee.alvName + std::to_string(iSCType),
0612                       ns.solid(ee.alvName + std::to_string(iSCType) + "Tmp"),
0613                       eeCutBox,
0614                       Transform3D(cutRot, cutTra - extra));
0615 
0616     const double mySign(iSCType < 4 ? +1. : -1.);
0617 
0618     const DDTranslation extraI(xyIOff + mySign * 2.0_mm, xyIOff + mySign * 2.0_mm, zIOff);
0619 
0620     Solid eeCutInt = SubtractionSolid(ee.intName + std::to_string(iSCType),
0621                       ns.solid(ee.intName + std::to_string(iSCType) + "Tmp"),
0622                       eeCutBox,
0623                       Transform3D(cutRot, cutTra - extraI));
0624 
0625     eeSCELog = ns.addVolumeNS(Volume(myns + ee.envName + std::to_string(iSCType), eeCutEnv, ns.material(ee.mat)));
0626     eeSCALog = Volume(ee.alvName + std::to_string(iSCType), eeCutAlv, ns.material(ee.wallMat));
0627     eeSCILog = Volume(ee.intName + std::to_string(iSCType), eeCutInt, ns.material(ee.mat));
0628       }
0629       eeSCELog.placeVolume(eeSCALog, iSCType * 100 + 1, Position(dxy, dxy, dz));
0630       eeSCALog.placeVolume(eeSCILog, iSCType * 100 + 1, Position(xyIOff, xyIOff, zIOff));
0631 
0632       DDTranslation croffset(0., 0., 0.);
0633 
0634       // Position crystals within parent supercrystal interior volume
0635       static const unsigned int ncol(5);
0636 
0637       if (iSCType > 0 && iSCType <= ee.nSCTypes) {
0638     const unsigned int icoffset((iSCType - 1) * ncol - 1);
0639 
0640     // Loop over columns of SC
0641     for (unsigned int icol(1); icol <= ncol; ++icol) {
0642       // Get column limits for this SC type from xml input
0643       const int ncrcol((int)ee.vecEESCProf[icoffset + icol]);
0644 
0645       const int imin(0 < ncrcol ? 1 : (0 > ncrcol ? ncol + ncrcol + 1 : 0));
0646       const int imax(0 < ncrcol ? ncrcol : (0 > ncrcol ? ncol : 0));
0647 
0648       if (imax > 0) {
0649         // Loop over crystals in this row
0650         for (int irow(imin); irow <= imax; ++irow) {
0651           // Create crystal as a DDEcalEndcapTrap object and calculate rotation and
0652           // translation required to position it in the SC.
0653           DDEcalEndcapTrap crystal(1, ee.crysFront, ee.crysRear, ee.crysLength);
0654 
0655           crystal.moveto(ee.cryFCtr[icol - 1][irow - 1], ee.cryRCtr[icol - 1][irow - 1]);
0656 
0657           string rname("EECrRoC" + std::to_string(icol) + "R" + std::to_string(irow));
0658 
0659           eeSCALog.placeVolume(
0660                    eeCRLog,
0661                    100 * iSCType + 10 * (icol - 1) + (irow - 1),
0662                    Transform3D(
0663                            myrot(ns, rname, crystal.rotation()),
0664                            Position(crystal.centrePos().x(), crystal.centrePos().y(), crystal.centrePos().z() - ee.cryZOff)));
0665         }
0666       }
0667     }
0668       }
0669     }
0670 
0671     //** Loop over endcap columns
0672     for (int icol = 1; icol <= int(ee.nColumns); icol++) {
0673       //**  Loop over SCs in column, using limits from xml input
0674       for (int irow = int(ee.vecEEShape[2 * icol - 2]); irow <= int(ee.vecEEShape[2 * icol - 1]); ++irow) {
0675     if (ee.vecEESCLims[0] <= icol && ee.vecEESCLims[1] >= icol && ee.vecEESCLims[2] <= irow &&
0676         ee.vecEESCLims[3] >= irow) {
0677       // Find SC type (complete or partial) for this location
0678       unsigned int isctype = 1;
0679 
0680       for (unsigned int ii = 0; ii < (unsigned int)(ee.nSCCutaway); ++ii) {
0681         if ((ee.vecEESCCutaway[3 * ii] == icol) && (ee.vecEESCCutaway[3 * ii + 1] == irow)) {
0682           isctype = int(ee.vecEESCCutaway[3 * ii + 2]);
0683         }
0684       }
0685 
0686       // Create SC as a DDEcalEndcapTrap object and calculate rotation and
0687       // translation required to position it in the endcap.
0688       DDEcalEndcapTrap scrys(1, ee.sCEFront, ee.sCERear, ee.sCELength);
0689       scrys.moveto(ee.scrFCtr[icol - 1][irow - 1], ee.scrRCtr[icol - 1][irow - 1]);
0690       scrys.translate(DDTranslation(0., 0., -ee.zOff));
0691 
0692       string rname(ee.envName + std::to_string(isctype) + std::to_string(icol) + "R" + std::to_string(irow));
0693       // Position SC in endcap
0694       Volume quaLog = ns.volume(ee.quaName);
0695       Volume childEnvLog = ns.volume(myns + ee.envName + std::to_string(isctype));
0696       quaLog.placeVolume(childEnvLog,
0697                  100 * isctype + 10 * (icol - 1) + (irow - 1),
0698                  Transform3D(scrys.rotation(), scrys.centrePos()));
0699     }
0700       }
0701     }
0702 
0703     return 1;
0704   }
0705 }
0706 
0707 DECLARE_DDCMS_DETELEMENT(DDCMS_ecal_DDEcalEndcapAlgo, algorithm)