Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:16:00

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2022 Whitney Armstrong
0003 
0004 /** \addtogroup Trackers Trackers
0005  * \brief Type: **BarrelTrackerWithFrame**.
0006  * \author W. Armstrong
0007  *
0008  * \ingroup trackers
0009  *
0010  * @{
0011  */
0012 #include "DD4hep/DetFactoryHelper.h"
0013 #include "DD4hep/Printout.h"
0014 #include "DD4hep/Shapes.h"
0015 #include "DDRec/DetectorData.h"
0016 #include "DDRec/Surface.h"
0017 #include "XML/Layering.h"
0018 #include "XML/Utilities.h"
0019 #include "DD4hepDetectorHelper.h"
0020 #include <array>
0021 #include <map>
0022 
0023 using namespace std;
0024 using namespace dd4hep;
0025 using namespace dd4hep::rec;
0026 using namespace dd4hep::detail;
0027 
0028 /** Endcap Trapezoidal Tracker.
0029  *
0030  * @author Whitney Armstrong
0031  *
0032  */
0033 static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) {
0034   typedef vector<PlacedVolume> Placements;
0035   xml_det_t x_det = e;
0036   Material vacuum = description.vacuum();
0037   int det_id      = x_det.id();
0038   string det_name = x_det.nameStr();
0039   bool reflect    = x_det.reflect(false);
0040   DetElement sdet(det_name, det_id);
0041   Assembly assembly(det_name);
0042 
0043   Material air     = description.material("Air");
0044   Volume motherVol = description.pickMotherVolume(sdet);
0045   int m_id = 0, c_id = 0, n_sensor = 0;
0046   map<string, Volume> modules;
0047   map<string, Placements> sensitives;
0048   map<string, std::vector<VolPlane>> volplane_surfaces;
0049   map<string, std::array<double, 2>> module_thicknesses;
0050   PlacedVolume pv;
0051 
0052   // Set detector type flag
0053   dd4hep::xml::setDetectorTypeFlag(x_det, sdet);
0054   auto& params = DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(sdet);
0055 
0056   // Add the volume boundary material if configured
0057   for (xml_coll_t bmat(x_det, _Unicode(boundary_material)); bmat; ++bmat) {
0058     xml_comp_t x_boundary_material = bmat;
0059     DD4hepDetectorHelper::xmlToProtoSurfaceMaterial(x_boundary_material, params,
0060                                                     "boundary_material");
0061   }
0062 
0063   assembly.setVisAttributes(description.invisible());
0064   sens.setType("tracker");
0065 
0066   for (xml_coll_t su(x_det, _U(support)); su; ++su) {
0067     xml_comp_t x_support     = su;
0068     double support_thickness = getAttrOrDefault(x_support, _U(thickness), 2.0 * mm);
0069     double support_length    = getAttrOrDefault(x_support, _U(length), 2.0 * mm);
0070     double support_rmin      = getAttrOrDefault(x_support, _U(rmin), 2.0 * mm);
0071     double support_zstart    = getAttrOrDefault(x_support, _U(zstart), 2.0 * mm);
0072     std::string support_name =
0073         getAttrOrDefault<std::string>(x_support, _Unicode(name), "support_tube");
0074     std::string support_vis = getAttrOrDefault<std::string>(x_support, _Unicode(vis), "AnlRed");
0075     xml_dim_t pos(x_support.child(_U(position), false));
0076     xml_dim_t rot(x_support.child(_U(rotation), false));
0077     Solid support_solid;
0078     if (x_support.hasChild(_U(shape))) {
0079       xml_comp_t shape(x_support.child(_U(shape)));
0080       string shape_type = shape.typeStr();
0081       support_solid     = xml::createShape(description, shape_type, shape);
0082     } else {
0083       support_solid = Tube(support_rmin, support_rmin + support_thickness, support_length / 2);
0084     }
0085     Transform3D tr =
0086         Transform3D(Rotation3D(),
0087                     Position(0, 0, (reflect ? -1.0 : 1.0) * (support_zstart + support_length / 2)));
0088     if (pos.ptr() && rot.ptr()) {
0089       Rotation3D rot3D(RotationZYX(rot.z(0), rot.y(0), rot.x(0)));
0090       Position pos3D(pos.x(0), pos.y(0), pos.z(0));
0091       tr = Transform3D(rot3D, pos3D);
0092     } else if (pos.ptr()) {
0093       tr = Transform3D(Rotation3D(), Position(pos.x(0), pos.y(0), pos.z(0)));
0094     } else if (rot.ptr()) {
0095       Rotation3D rot3D(RotationZYX(rot.z(0), rot.y(0), rot.x(0)));
0096       tr = Transform3D(rot3D, Position());
0097     }
0098     Material support_mat = description.material(x_support.materialStr());
0099     Volume support_vol(support_name, support_solid, support_mat);
0100     support_vol.setVisAttributes(description.visAttributes(support_vis));
0101     pv = assembly.placeVolume(support_vol, tr);
0102     // pv = assembly.placeVolume(support_vol, Position(0, 0, support_zstart + support_length / 2));
0103   }
0104 
0105   for (xml_coll_t mi(x_det, _U(module)); mi; ++mi, ++m_id) {
0106     xml_comp_t x_mod = mi;
0107     string m_nam     = x_mod.nameStr();
0108     xml_comp_t trd   = x_mod.trd();
0109 
0110     double posY;
0111     double x1              = trd.x1();
0112     double x2              = trd.x2();
0113     double z               = trd.z();
0114     double total_thickness = 0.;
0115     xml_coll_t ci(x_mod, _U(module_component));
0116     for (ci.reset(), total_thickness = 0.0; ci; ++ci)
0117       total_thickness += xml_comp_t(ci).thickness();
0118 
0119     double thickness_so_far = 0.0;
0120     double y1               = total_thickness / 2;
0121     double y2               = total_thickness / 2;
0122     Trapezoid m_solid(x1, x2, y1, y2, z);
0123     Volume m_volume(m_nam, m_solid, vacuum);
0124     m_volume.setVisAttributes(description.visAttributes(x_mod.visStr()));
0125 
0126     Solid frame_s;
0127     if (x_mod.hasChild(_U(frame))) {
0128       // build frame from trd (assumed to be smaller)
0129       xml_comp_t m_frame     = x_mod.child(_U(frame));
0130       xml_comp_t f_pos       = m_frame.child(_U(position));
0131       xml_comp_t frame_trd   = m_frame.trd();
0132       double frame_thickness = getAttrOrDefault(m_frame, _U(thickness), total_thickness);
0133       double frame_x1        = frame_trd.x1();
0134       double frame_x2        = frame_trd.x2();
0135       double frame_z         = frame_trd.z();
0136       // make the frame match the total thickness if thickness attribute is not given
0137       Trapezoid f_solid1(x1, x2, frame_thickness / 2.0, frame_thickness / 2.0, z);
0138       Trapezoid f_solid(frame_x1, frame_x2, frame_thickness / 2.0, frame_thickness / 2.0, frame_z);
0139       SubtractionSolid frame_shape(f_solid1, f_solid);
0140       frame_s = frame_shape;
0141 
0142       Material f_mat = description.material(m_frame.materialStr());
0143       Volume f_vol(m_nam + "_frame", frame_shape, f_mat);
0144       f_vol.setVisAttributes(description.visAttributes(m_frame.visStr()));
0145 
0146       // figure out how to best place
0147       pv = m_volume.placeVolume(f_vol, Position(f_pos.x(), f_pos.y(), f_pos.z()));
0148     }
0149 
0150     for (ci.reset(), n_sensor = 1, c_id = 0, posY = -y1; ci; ++ci, ++c_id) {
0151       xml_comp_t c     = ci;
0152       double c_thick   = c.thickness();
0153       auto comp_x1     = getAttrOrDefault(c, _Unicode(x1), x1);
0154       auto comp_x2     = getAttrOrDefault(c, _Unicode(x2), x2);
0155       auto comp_height = getAttrOrDefault(c, _Unicode(height), z);
0156 
0157       Material c_mat = description.material(c.materialStr());
0158       string c_name  = _toString(c_id, "component%d");
0159 
0160       Trapezoid comp_s1(comp_x1, comp_x2, c_thick / 2e0, c_thick / 2e0, comp_height);
0161       Solid comp_shape = comp_s1;
0162       if (frame_s.isValid()) {
0163         comp_shape = SubtractionSolid(comp_s1, frame_s);
0164       }
0165       Volume c_vol(c_name, comp_shape, c_mat);
0166 
0167       c_vol.setVisAttributes(description.visAttributes(c.visStr()));
0168       pv = m_volume.placeVolume(c_vol, Position(0, posY + c_thick / 2, 0));
0169       if (c.isSensitive()) {
0170         module_thicknesses[m_nam] = {thickness_so_far + c_thick / 2.0,
0171                                      total_thickness - thickness_so_far - c_thick / 2.0};
0172         // std::cout << " adding sensitive volume" << c_name << "\n";
0173         sdet.check(n_sensor > 2,
0174                    "SiTrackerEndcap2::fromCompact: " + c_name + " Max of 2 modules allowed!");
0175         pv.addPhysVolID("sensor", n_sensor);
0176         c_vol.setSensitiveDetector(sens);
0177         sensitives[m_nam].push_back(pv);
0178         ++n_sensor;
0179         // -------- create a measurement plane for the tracking surface attched to the sensitive volume -----
0180         Vector3D u(0., 0., -1.);
0181         Vector3D v(-1., 0., 0.);
0182         Vector3D n(0., 1., 0.);
0183         // Vector3D o( 0. , 0. , 0. ) ;
0184 
0185         // compute the inner and outer thicknesses that need to be assigned to the tracking surface
0186         // depending on wether the support is above or below the sensor
0187         double inner_thickness = module_thicknesses[m_nam][0];
0188         double outer_thickness = module_thicknesses[m_nam][1];
0189 
0190         SurfaceType type(SurfaceType::Sensitive);
0191 
0192         // if( isStripDetector )
0193         //  type.setProperty( SurfaceType::Measurement1D , true ) ;
0194 
0195         VolPlane surf(c_vol, type, inner_thickness, outer_thickness, u, v, n); //,o ) ;
0196         volplane_surfaces[m_nam].push_back(surf);
0197 
0198         //--------------------------------------------
0199       }
0200       posY += c_thick;
0201       thickness_so_far += c_thick;
0202     }
0203     modules[m_nam] = m_volume;
0204   }
0205 
0206   for (xml_coll_t li(x_det, _U(layer)); li; ++li) {
0207     xml_comp_t x_layer(li);
0208     int l_id    = x_layer.id();
0209     int mod_num = 1;
0210 
0211     xml_comp_t l_env  = x_layer.child(_U(envelope));
0212     string layer_name = det_name + std::string("_layer") + std::to_string(l_id);
0213 
0214     std::string layer_vis = l_env.attr<std::string>(_Unicode(vis));
0215     double layer_rmin     = l_env.attr<double>(_Unicode(rmin));
0216     double layer_rmax     = l_env.attr<double>(_Unicode(rmax));
0217     double layer_length   = l_env.attr<double>(_Unicode(length));
0218     double layer_zstart   = l_env.attr<double>(_Unicode(zstart));
0219     double layer_center_z = layer_zstart + layer_length / 2.0;
0220     // printout(INFO,"ROOTGDMLParse","+++ Read geometry from GDML file file:%s",input.c_str());
0221     // std::cout << "SiTracker Endcap layer " << l_id << " zstart = " << layer_zstart/dd4hep::mm << "mm ( " <<
0222     // layer_length/dd4hep::mm << " mm thick )\n";
0223 
0224     // Assembly    layer_assembly(layer_name);
0225     // assembly.placeVolume(layer_assembly);
0226     Tube layer_tub(layer_rmin, layer_rmax, layer_length / 2);
0227     Volume layer_vol(layer_name, layer_tub, air); // Create the layer envelope volume.
0228     layer_vol.setVisAttributes(description.visAttributes(layer_vis));
0229 
0230     PlacedVolume layer_pv;
0231     if (reflect) {
0232       layer_pv = assembly.placeVolume(
0233           layer_vol, Transform3D(RotationZYX(0.0, -M_PI, 0.0), Position(0, 0, -layer_center_z)));
0234       layer_pv.addPhysVolID("layer", l_id);
0235       layer_name += "_N";
0236     } else {
0237       layer_pv = assembly.placeVolume(layer_vol, Position(0, 0, layer_center_z));
0238       layer_pv.addPhysVolID("layer", l_id);
0239       layer_name += "_P";
0240     }
0241     DetElement layer_element(sdet, layer_name, l_id);
0242     layer_element.setPlacement(layer_pv);
0243 
0244     auto& layerParams =
0245         DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(layer_element);
0246 
0247     for (xml_coll_t lmat(x_layer, _Unicode(layer_material)); lmat; ++lmat) {
0248       xml_comp_t x_layer_material = lmat;
0249       DD4hepDetectorHelper::xmlToProtoSurfaceMaterial(x_layer_material, layerParams,
0250                                                       "layer_material");
0251     }
0252 
0253     for (xml_coll_t ri(x_layer, _U(ring)); ri; ++ri) {
0254       xml_comp_t x_ring    = ri;
0255       double r             = x_ring.r();
0256       double phi0          = x_ring.phi0(0);
0257       double zstart        = x_ring.zstart();
0258       double dz            = x_ring.dz(0);
0259       int nmodules         = x_ring.nmodules();
0260       string m_nam         = x_ring.moduleStr();
0261       Volume m_vol         = modules[m_nam];
0262       double iphi          = 2 * M_PI / nmodules;
0263       double phi           = phi0;
0264       Placements& sensVols = sensitives[m_nam];
0265 
0266       for (int k = 0; k < nmodules; ++k) {
0267         string m_base = _toString(l_id, "layer%d") + _toString(mod_num, "_module%d");
0268         double x      = -r * std::cos(phi);
0269         double y      = -r * std::sin(phi);
0270 
0271         if (!reflect) {
0272           DetElement module(layer_element, m_base + "_pos", det_id);
0273           pv = layer_vol.placeVolume(m_vol, Transform3D(RotationZYX(0, -M_PI / 2 - phi, -M_PI / 2),
0274                                                         Position(x, y, zstart + dz)));
0275           pv.addPhysVolID("module", mod_num);
0276           module.setPlacement(pv);
0277           for (size_t ic = 0; ic < sensVols.size(); ++ic) {
0278             PlacedVolume sens_pv = sensVols[ic];
0279             DetElement comp_elt(module, sens_pv.volume().name(), mod_num);
0280             auto& comp_elt_params =
0281                 DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(comp_elt);
0282             comp_elt_params.set<string>("axis_definitions", "XZY");
0283             comp_elt.setPlacement(sens_pv);
0284             volSurfaceList(comp_elt)->push_back(volplane_surfaces[m_nam][ic]);
0285           }
0286         } else {
0287           pv = layer_vol.placeVolume(m_vol, Transform3D(RotationZYX(0, -M_PI / 2 - phi, -M_PI / 2),
0288                                                         Position(x, y, -zstart - dz)));
0289           pv.addPhysVolID("module", mod_num);
0290           DetElement r_module(layer_element, m_base + "_neg", det_id);
0291           r_module.setPlacement(pv);
0292           for (size_t ic = 0; ic < sensVols.size(); ++ic) {
0293             PlacedVolume sens_pv = sensVols[ic];
0294             DetElement comp_elt(r_module, sens_pv.volume().name(), mod_num);
0295             auto& comp_elt_params =
0296                 DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(comp_elt);
0297             comp_elt_params.set<string>("axis_definitions", "XZY");
0298             comp_elt.setPlacement(sens_pv);
0299             volSurfaceList(comp_elt)->push_back(volplane_surfaces[m_nam][ic]);
0300           }
0301         }
0302         dz = -dz;
0303         phi += iphi;
0304         ++mod_num;
0305       }
0306     }
0307   }
0308   pv = motherVol.placeVolume(assembly, Position(0, 0, (reflect ? -1.0e-9 : 1.0e-9)));
0309   pv.addPhysVolID("system", det_id);
0310   sdet.setPlacement(pv);
0311   return sdet;
0312 }
0313 
0314 DECLARE_DETELEMENT(epic_TrapEndcapTracker, create_detector)