epic/src/CurvedBarrelTracker_geo.cpp

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2022 - 2025 Whitney Armstrong, Jonathan Witte, Shujie Li
0003
0004 /** Curved Barrel tracker
0005  * - Derived from "BarrelTrackerWithFrame_geo.cpp".
0006  *
0007  * - Designed to process "vertex_barrel.xml":
0008  *       - When the upper and lower modules are provided, will use the
0009  *         Build-in EIC-LAS RSU structure with four sections and inactive areas
0010  *
0011  *
0012  * \code
0013  * \endcode
0014  *
0015  *
0016  * @author Whitney Armstrong, Jonathan Witte, Shujie Li
0017  */
0018
0019 #include "DD4hep/DetFactoryHelper.h"
0020 #include "DD4hep/Printout.h"
0021 #include "DD4hep/Shapes.h"
0022 #include "DDRec/DetectorData.h"
0023 #include "DDRec/Surface.h"
0024 #include "XML/Layering.h"
0025 #include "XML/Utilities.h"
0026 #include <array>
0027 #include "DD4hepDetectorHelper.h"
0028
0029 using namespace std;
0030 using namespace dd4hep;
0031 using namespace dd4hep::rec;
0032 using namespace dd4hep::detail;
0033
0034 static Ref_t create_CurvedBarrelTracker(Detector& description, xml_h e, SensitiveDetector sens) {
0035
0036   xml_det_t x_det = e;
0037   Material air    = description.air();
0038   int det_id      = x_det.id();
0039   string det_name = x_det.nameStr();
0040   DetElement sdet(det_name, det_id);
0041
0042   typedef vector<PlacedVolume> Placements;
0043   map<string, Volume> volumes;
0044   map<string, Placements> sensitives;
0045   map<string, std::vector<VolPlane>> volplane_surfaces;
0046   map<string, std::vector<double>> module_length;
0047
0048   PlacedVolume pv, pv_frame;
0049
0050   // Set detector type flag
0051   dd4hep::xml::setDetectorTypeFlag(x_det, sdet);
0052   auto& params = DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(sdet);
0053
0054   // Add the volume boundary material if configured
0055   for (xml_coll_t bmat(x_det, _Unicode(boundary_material)); bmat; ++bmat) {
0056     xml_comp_t x_boundary_material = bmat;
0057     DD4hepDetectorHelper::xmlToProtoSurfaceMaterial(x_boundary_material, params,
0058                                                     "boundary_material");
0059   }
0060
0061   Assembly assembly(det_name);
0062
0063   sens.setType("tracker");
0064
0065   // loop over the modules
0066   for (xml_coll_t mi(x_det, _U(module)); mi; ++mi) {
0067     xml_comp_t x_mod = mi;
0068     string m_nam     = x_mod.nameStr();
0069     double m_rmin    = x_mod.rmin();
0070     double m_length  = x_mod.length();
0071     double m_width   = x_mod.width();
0072     module_length[m_nam].push_back(m_length);
0073
0074     if (volumes.find(m_nam) != volumes.end()) {
0075       printout(ERROR, "CurvedBarrelTracker",
0076                string((string("Module with named ") + m_nam + string(" already exists."))).c_str());
0077       throw runtime_error("Logics error in building modules.");
0078     }
0079     int ncomponents        = 0;
0080     int sensor_number      = 1;
0081     double total_thickness = 0;
0082
0083     // Compute module total thickness from components
0084     xml_coll_t ci(x_mod, _U(module_component));
0085     for (ci.reset(), total_thickness = 0.0; ci; ++ci) {
0086       total_thickness += xml_comp_t(ci).thickness();
0087     }
0088     // the module assembly volume
0089     Assembly m_vol(m_nam);
0090     volumes[m_nam] = m_vol;
0091     m_vol.setVisAttributes(description.visAttributes(x_mod.visStr()));
0092
0093     double thickness_so_far = 0.0;
0094     for (xml_coll_t mci(x_mod, _U(module_component)); mci; ++mci, ++ncomponents) {
0095       Volume c_vol;
0096       xml_comp_t x_comp  = mci;
0097       const string c_nam = x_comp.nameStr(); // _toString(ncomponents, "component%d");
0098       const string c_mat = x_comp.materialStr();
0099       double c_thickness = x_comp.thickness();
0100       double c_width     = getAttrOrDefault(x_comp, _U(width), m_width);
0101       double c_length    = getAttrOrDefault(x_comp, _U(length), m_length);
0102       double c_rmin      = m_rmin + thickness_so_far;
0103       double c_dphi      = c_width / c_rmin;
0104
0105       if (c_nam == "RSU") { // for RSU, create ONE upper or lower 3-tile sections.
0106         // **** hard-coded RSU design with 12 tiles, plus backbones, readout pads, biasing
0107         // Having issue including multiple sensitive surfaces in one Tube module (worked with box).
0108         // Therefore use type "upper" and "lower" to create two mirrored tile sections. (left right is identical)
0109         //
0110         // "|" = backbone. Length alone z.
0111         //
0112         // | ------readout-------- | -------readout--------
0113         // | tilex3                | tilex3
0114         // | ------biasing-------- | -------biasing--------
0115         // | ------biasing-------- | -------biasing--------
0116         // | tilex3                | tilex3
0117         // | ------readout-------- | -------readout--------
0118         const string frame_vis        = "SVTReadoutVis";
0119         const string c_type           = x_comp.typeStr();
0120         const double BiasingWidth     = 0.06 * mm; // need to x2 for two sets
0121         const double ReadoutPadsWidth = m_width - BiasingWidth - c_width;
0122         const double BackboneLength   = m_length - c_length; //0.06*mm;
0123
0124         double px = 0, py = 0, pz = 0;
0125         double c_z0 = -m_length / 2;
0126         double c_z1 = c_z0 + BackboneLength;
0127         double c_z2 = m_length / 2;
0128         pz          = (c_z1 + c_z2) / 2; // section central z
0129
0130         double c_phi0 = 0;
0131         double c_phi1, c_phi2;
0132         double c_phi3 = m_width / m_rmin;
0133         string c_nam1, c_nam2;
0134         if (c_type == "upper") {
0135           c_phi1 = BiasingWidth / c_rmin;
0136           c_phi2 = c_phi1 + c_dphi;
0137           c_nam1 = "biasing";
0138           c_nam2 = "readout";
0139         } else if (c_type == "lower") {
0140           c_phi1 = ReadoutPadsWidth / c_rmin;
0141           c_phi2 = c_phi1 + c_dphi;
0142           c_nam1 = "readout";
0143           c_nam2 = "biasing";
0144         } else {
0145           printout(ERROR, "CurvedBarrelTracker",
0146                    string((string("Module ") + m_nam + string(": invalid RSU component type [") +
0147                            c_type + string("], should be upper or lower")))
0148                        .c_str());
0149           throw runtime_error("Logics error in building modules.");
0150         }
0151
0152         // *** inactive areas
0153         // biasing and readout (horizontal)
0154         Tube f_tube1(c_rmin, c_rmin + c_thickness, c_length / 2, c_phi0, c_phi1);
0155         Volume f_vol1(c_nam1, f_tube1, description.material(c_mat));
0156         pv_frame = m_vol.placeVolume(f_vol1, Position(px, py, pz));
0157         f_vol1.setVisAttributes(description, frame_vis);
0158
0159         Tube f_tube2(c_rmin, c_rmin + c_thickness, c_length / 2, c_phi2, c_phi3);
0160         Volume f_vol2(c_nam2, f_tube2, description.material(c_mat));
0161         pv_frame = m_vol.placeVolume(f_vol2, Position(px, py, pz));
0162         f_vol2.setVisAttributes(description, frame_vis);
0163
0164         // backbone (vertical)
0165         Tube f_tube3(c_rmin, c_rmin + c_thickness, BackboneLength / 2, c_phi0, c_phi3);
0166         Volume f_vol3("backbone", f_tube3, description.material(c_mat));
0167         pv_frame = m_vol.placeVolume(f_vol3, Position(px, py, (c_z0 + c_z1) / 2));
0168         f_vol3.setVisAttributes(description, frame_vis);
0169
0170         // *** sensitive tile
0171         Tube c_tube(c_rmin, c_rmin + c_thickness, c_length / 2, c_phi1, c_phi2);
0172         c_vol = Volume(c_nam + "_" + c_type, c_tube, description.material(c_mat));
0173         pv    = m_vol.placeVolume(c_vol, Position(px, py, pz));
0174       } else { // for regular component, no difference b/w upper/lower
0175         Tube c_tube(c_rmin, c_rmin + c_thickness, c_length / 2, 0, c_dphi);
0176         c_vol = Volume(c_nam, c_tube, description.material(c_mat));
0177         pv    = m_vol.placeVolume(c_vol, Position(0, 0, 0));
0178       }
0179       c_vol.setRegion(description, x_comp.regionStr());
0180       c_vol.setLimitSet(description, x_comp.limitsStr());
0181       c_vol.setVisAttributes(description, x_comp.visStr());
0182       if (x_comp.isSensitive()) {
0183         pv.addPhysVolID("sensor", sensor_number++);
0184         c_vol.setSensitiveDetector(sens);
0185         sensitives[m_nam].push_back(pv);
0186         // -------- create a measurement plane for the tracking surface attched to the sensitive volume -----
0187         Vector3D u(-1., 0., 0.);
0188         Vector3D v(0., -1., 0.);
0189         Vector3D n(0., 0., 1.);
0190         // compute the inner and outer thicknesses that need to be assigned to the tracking surface
0191         // depending on wether the support is above or below the sensor
0192         double inner_thickness = thickness_so_far + c_thickness / 2.0;
0193         double outer_thickness = total_thickness - inner_thickness;
0194
0195         SurfaceType type(SurfaceType::Sensitive);
0196         VolPlane surf(c_vol, type, inner_thickness, outer_thickness, u, v, n); //,o ) ;
0197         volplane_surfaces[m_nam].push_back(surf);
0198       }
0199       thickness_so_far += c_thickness;
0200     }
0201   }
0202
0203   // now build the layers by alternating upper and lower modules.
0204   // each layer has a gap b/w the upper and lower half (see phi0 and phi1)
0205   for (xml_coll_t li(x_det, _U(layer)); li; ++li) {
0206     xml_comp_t x_layer  = li;
0207     xml_comp_t x_barrel = x_layer.child(_U(barrel_envelope));
0208     xml_comp_t x_layout = x_layer.child(_U(rphi_layout));
0209     xml_comp_t z_layout = x_layer.child(_U(z_layout)); // Get the <z_layout> element.
0210     int lay_id          = x_layer.id();
0211     string m_nam        = x_layer.moduleStr();
0212     string lay_nam      = det_name + _toString(x_layer.id(), "_layer%d");
0213     Tube lay_tub(x_barrel.inner_r(), x_barrel.outer_r(), x_barrel.z_length() / 2.0);
0214
0215     Volume lay_vol(lay_nam, lay_tub, air); // Create the layer envelope volume.
0216     Position lay_pos(0, 0, getAttrOrDefault(x_barrel, _U(z0), 0.));
0217     lay_vol.setVisAttributes(description.visAttributes(x_layer.visStr()));
0218
0219     int l_nphi      = x_layout.nphi(); // Number of modules in phi.
0220     double phi0     = x_layout.phi0(); // Starting phi of first module.
0221     int nphi[2]     = {int((l_nphi + 1) / 2),
0222                        int(l_nphi / 2)}; // number of modules in uppper and lower modules.
0223     double phi_incr = (M_PI * 2 - phi0 * 2) / l_nphi; // Phi increment for one module.
0224     double z0       = z_layout.z0();                  // Z position of first module in phi.
0225     double nz       = z_layout.nz();                  // Number of modules to place in z.
0226
0227     Volume module_env[2];
0228     Placements sensVols[2];
0229     string m_nams[2] = {m_nam + "_upper", m_nam + "_lower"};
0230     // if both upper and lower modules are provided (for RSU tiles)
0231     if ((volumes.find(m_nams[0]) != volumes.end()) && (volumes.find(m_nams[1]) != volumes.end())) {
0232       if (nphi[0] != nphi[1]) {
0233         printout(ERROR, "CurvedBarrelTracker",
0234                  string((string("Layer ") + lay_nam +
0235                          string(": nphi must be even number to allow upper and lower modules"))
0236                             .c_str()));
0237         throw runtime_error("Logics error in building modules.");
0238       }
0239       module_env[0] = volumes[m_nams[0]];
0240       module_env[1] = volumes[m_nams[1]];
0241       sensVols[0]   = sensitives[m_nams[0]];
0242       sensVols[1]   = sensitives[m_nams[1]];
0243     }
0244     // for other regular modules
0245     else if (volumes.find(m_nam) != volumes.end()) {
0246       module_env[0] = volumes[m_nam];
0247       module_env[1] = volumes[m_nam];
0248       sensVols[0]   = sensitives[m_nam];
0249       sensVols[1]   = sensitives[m_nam];
0250       m_nams[0]     = m_nam;
0251       m_nams[1]     = m_nam;
0252     }
0253     DetElement lay_elt(sdet, lay_nam, lay_id);
0254
0255     // the local coordinate systems of modules in dd4hep and acts differ
0256     // see http://acts.web.cern.ch/ACTS/latest/doc/group__DD4hepPlugins.html
0257     auto& layerParams =
0258         DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(lay_elt);
0259
0260     for (xml_coll_t lmat(x_layer, _Unicode(layer_material)); lmat; ++lmat) {
0261       xml_comp_t x_layer_material = lmat;
0262       DD4hepDetectorHelper::xmlToProtoSurfaceMaterial(x_layer_material, layerParams,
0263                                                       "layer_material");
0264     }
0265
0266     // Z increment for module placement along Z axis.
0267     // Adjust for z0 at center of module rather than
0268     // the end of cylindrical envelope.
0269     // double z_incr = nz > 1 ? (2.0 * abs(z0)) / (nz - 1) : 0.0;
0270     // Starting z for module placement along Z axis.
0271     int module = 1;
0272     // Loop over the number of modules in phi.
0273     for (int kk = 0; kk < 2; kk++) {
0274       int iphi      = nphi[kk];
0275       double z_incr = module_length[m_nams[kk]][0];
0276       for (int ii = 0; ii < iphi; ii++) {
0277         double dphi = phi0;
0278         if (ii > (iphi / 2) - 1)
0279           dphi = 2 * phi0;
0280         // Loop over the number of modules in z.
0281         double module_z = z0;
0282         for (int j = 0; j < nz; j++, module_z += z_incr) {
0283           string module_name = _toString(module, "module%d");
0284           DetElement mod_elt(lay_elt, module_name, module);
0285           Transform3D tr(
0286               RotationZYX(dphi + phi_incr * kk + phi_incr * ii * 2, 0, 0),
0287               Position(0, 0, module_z)); // altering upper and lower module to fill every other row
0288           pv = lay_vol.placeVolume(module_env[kk], tr);
0289           pv.addPhysVolID("module", module);
0290           mod_elt.setPlacement(pv);
0291           for (size_t ic = 0; ic < sensVols[kk].size(); ++ic) {
0292             PlacedVolume sens_pv = sensVols[kk][ic];
0293             DetElement comp_de(mod_elt, std::string("de_") + sens_pv.volume().name(), module);
0294             comp_de.setPlacement(sens_pv);
0295             auto& comp_de_params =
0296                 DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(comp_de);
0297             comp_de_params.set<string>("axis_definitions", "XYZ");
0298             // note from Wouter: Region allows setting specific physics in regions of space. Limits similar but with production limits. Currently unused in simulations. But we'd sort of tried to make sure we don't prevent ourselves from using it entirely.
0299             // comp_de.setAttributes(description, sens_pv.volume(), x_layer.regionStr(), x_layer.limitsStr(),
0300             //                       xml_det_t(xmleles[m_nam]).visStr());
0301             //
0302             volSurfaceList(comp_de)->push_back(volplane_surfaces[m_nams[kk]][ic]);
0303           }
0304
0305           /// Increase counters etc.
0306           module++;
0307         }
0308       }
0309     }
0310     // Create the PhysicalVolume for the layer.
0311     pv = assembly.placeVolume(lay_vol, lay_pos); // Place layer in mother
0312     pv.addPhysVolID("layer", lay_id);            // Set the layer ID.
0313     lay_elt.setAttributes(description, lay_vol, x_layer.regionStr(), x_layer.limitsStr(),
0314                           x_layer.visStr());
0315     lay_elt.setPlacement(pv);
0316   }
0317   sdet.setAttributes(description, assembly, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());
0318   assembly.setVisAttributes(description.invisible());
0319   pv = description.pickMotherVolume(sdet).placeVolume(assembly);
0320   pv.addPhysVolID("system", det_id); // Set the subdetector system ID.
0321   sdet.setPlacement(pv);
0322   return sdet;
0323 }
0324
0325 //@}
0326 // clang-format off
0327 DECLARE_DETELEMENT(epic_CylinderSVTBarrel,create_CurvedBarrelTracker)