Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:15:57

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2022 Wouter Deconinck, Whitney Armstrong, Sylvester Joosten
0003 
0004 //==========================================================================
0005 //
0006 //      <detector name ="DetName" type="Beampipe" >
0007 //      <layer id="#(int)" inner_r="#(double)" outer_z="#(double)" >
0008 //      <slice material="string" thickness="#(double)" >
0009 //      </layer>
0010 //      </detector>
0011 //==========================================================================
0012 #include "DD4hep/DetFactoryHelper.h"
0013 #include "DD4hep/Printout.h"
0014 #include "TMath.h"
0015 #include <XML/Helper.h>
0016 #include "XML/Utilities.h"
0017 #include "DD4hepDetectorHelper.h"
0018 
0019 using namespace std;
0020 using namespace dd4hep;
0021 
0022 /** \addtogroup beamline Beamline Instrumentation
0023  */
0024 
0025 /** \addtogroup IRChamber Interaction Region Vacuum Chamber.
0026  * \brief Type: **IRChamber**.
0027  * \ingroup beamline
0028  *
0029  *
0030  * \code
0031  *   <detector>
0032  *   </detector>
0033  * \endcode
0034  *
0035  */
0036 static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) {
0037 
0038   using namespace ROOT::Math;
0039   xml_det_t x_det       = e;
0040   string det_name       = x_det.nameStr();
0041   xml_comp_t x_dettype  = x_det.child(dd4hep::xml::Strng_t("type_flags"));
0042   unsigned int typeFlag = x_dettype.type();
0043   DetElement sdet(det_name, x_det.id());
0044   Assembly assembly(det_name + "_assembly");
0045   Material m_Al     = det.material("Aluminum");
0046   Material m_Be     = det.material("Beryllium");
0047   Material m_Au     = det.material("Gold");
0048   Material m_Vacuum = det.material("Vacuum");
0049   string vis_name   = x_det.visStr();
0050 
0051   xml::Component IP_pipe_c = x_det.child(_Unicode(IP_pipe));
0052 
0053   // IP
0054   double IP_beampipe_OD             = IP_pipe_c.attr<double>(_Unicode(OD));
0055   double IP_beampipe_wall_thickness = IP_pipe_c.attr<double>(_Unicode(wall_thickness));
0056   double IP_beampipe_gold_thickness = IP_pipe_c.attr<double>(_Unicode(gold_thickness));
0057   double IP_beampipe_ID =
0058       IP_beampipe_OD - 2.0 * IP_beampipe_gold_thickness - 2.0 * IP_beampipe_wall_thickness;
0059   double IP_acts_beampipe_OD = IP_beampipe_ID - 5.0 * mm;
0060   double IP_acts_beampipe_ID = IP_acts_beampipe_OD - 1.0 * mm;
0061 
0062   double upstream_straight_length   = IP_pipe_c.attr<double>(_Unicode(upstream_straight_length));
0063   double downstream_straight_length = IP_pipe_c.attr<double>(_Unicode(downstream_straight_length));
0064 
0065   // central beampipe volume
0066   Tube central_tube(0.5 * IP_acts_beampipe_ID, 0.5 * IP_acts_beampipe_OD,
0067                     0.5 * (upstream_straight_length + downstream_straight_length));
0068   Volume central_volume("acts_central_beampipe_vol", central_tube, m_Vacuum);
0069   const double central_offset = -.5 * (upstream_straight_length - downstream_straight_length);
0070   DetElement central_det(sdet, "acts_beampipe_central", 1);
0071 
0072   // Set dd4hep variant parameters for conversion to ACTS tracking geometry
0073   central_det.setTypeFlag(typeFlag);
0074   auto& params = DD4hepDetectorHelper::ensureExtension<dd4hep::rec::VariantParameters>(central_det);
0075   int nBinPhi  = 144; // fix later. Should take this from a xml tag
0076   int nBinZ    = 10;  // fix later. Should take this from a xml tag
0077   params.set<bool>("layer_material", true);
0078   params.set<bool>("layer_material_representing", true);
0079   params.set<int>("layer_material_representing_binPhi", nBinPhi);
0080   params.set<int>("layer_material_representing_binZ", nBinZ);
0081 
0082   // -----------------------------
0083   // IP beampipe
0084   //
0085   // setup for the central IP beampipe:
0086   //
0087   // /-------\ Be wall
0088   //  /-----\  Au coating
0089   //   /---\   Vacuum padding (5mm)
0090   //    /-\    Fake vacuum beampipe (1mm)
0091   //     -     Vacuum filled inner beampipe
0092   //
0093   Tube downstream_IP_vacuum_fill(0.0, IP_acts_beampipe_ID / 2.0, downstream_straight_length / 2.0);
0094   Tube downstream_IP_acts_beampipe(IP_acts_beampipe_ID / 2.0, IP_acts_beampipe_OD / 2.0,
0095                                    downstream_straight_length / 2.0);
0096   Tube downstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0,
0097                                     downstream_straight_length / 2.0);
0098   Tube downstream_IP_gold(IP_beampipe_ID / 2.0, IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness,
0099                           downstream_straight_length / 2.0);
0100   Tube downstream_IP_tube(IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, IP_beampipe_OD / 2.0,
0101                           downstream_straight_length / 2.0);
0102 
0103   Tube upstream_IP_vacuum_fill(0.0, IP_acts_beampipe_ID / 2.0, upstream_straight_length / 2.0);
0104   Tube upstream_IP_acts_beampipe(IP_acts_beampipe_ID / 2.0, IP_acts_beampipe_OD / 2.0,
0105                                  upstream_straight_length / 2.0);
0106   Tube upstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0,
0107                                   upstream_straight_length / 2.0);
0108   Tube upstream_IP_gold(IP_beampipe_ID / 2.0, IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness,
0109                         upstream_straight_length / 2.0);
0110   Tube upstream_IP_tube(IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, IP_beampipe_OD / 2.0,
0111                         upstream_straight_length / 2.0);
0112 
0113   Volume v_downstream_IP_vacuum_fill("v_downstream_IP_vacuum_fill", downstream_IP_vacuum_fill,
0114                                      m_Vacuum);
0115   Volume v_downstream_IP_acts_beampipe("v_downstream_IP_acts_beampipe", downstream_IP_acts_beampipe,
0116                                        m_Vacuum);
0117   Volume v_downstream_IP_vacuum_padding("v_downstream_IP_vacuum_padding",
0118                                         downstream_IP_vacuum_padding, m_Vacuum);
0119   Volume v_downstream_IP_gold("v_downstream_IP_gold", downstream_IP_gold, m_Au);
0120   Volume v_downstream_IP_tube("v_downstream_IP_tube", downstream_IP_tube, m_Be);
0121   Volume v_upstream_IP_vacuum_fill("v_upstream_IP_vacuum_fill", upstream_IP_vacuum_fill, m_Vacuum);
0122   Volume v_upstream_IP_acts_beampipe("v_upstream_IP_acts_beampipe", upstream_IP_acts_beampipe,
0123                                      m_Vacuum);
0124   Volume v_upstream_IP_vacuum_padding("v_upstream_IP_vacuum_padding", upstream_IP_vacuum_padding,
0125                                       m_Vacuum);
0126   Volume v_upstream_IP_gold("v_upstream_IP_gold", upstream_IP_gold, m_Au);
0127   Volume v_upstream_IP_tube("v_upstream_IP_tube", upstream_IP_tube, m_Be);
0128 
0129   sdet.setAttributes(det, v_upstream_IP_gold, x_det.regionStr(), x_det.limitsStr(), vis_name);
0130   sdet.setAttributes(det, v_upstream_IP_tube, x_det.regionStr(), x_det.limitsStr(), vis_name);
0131   sdet.setAttributes(det, v_downstream_IP_gold, x_det.regionStr(), x_det.limitsStr(), vis_name);
0132   sdet.setAttributes(det, v_downstream_IP_tube, x_det.regionStr(), x_det.limitsStr(), vis_name);
0133 
0134   assembly.placeVolume(v_upstream_IP_vacuum_fill, Position(0, 0, -upstream_straight_length / 2.0));
0135   central_volume.placeVolume(v_upstream_IP_acts_beampipe,
0136                              Position(0, 0, -upstream_straight_length / 2.0 - central_offset));
0137   assembly.placeVolume(v_upstream_IP_vacuum_padding,
0138                        Position(0, 0, -upstream_straight_length / 2.0));
0139   assembly.placeVolume(v_upstream_IP_gold, Position(0, 0, -upstream_straight_length / 2.0));
0140   assembly.placeVolume(v_upstream_IP_tube, Position(0, 0, -upstream_straight_length / 2.0));
0141 
0142   assembly.placeVolume(v_downstream_IP_vacuum_fill,
0143                        Position(0, 0, downstream_straight_length / 2.0));
0144   central_volume.placeVolume(v_downstream_IP_acts_beampipe,
0145                              Position(0, 0, downstream_straight_length / 2.0 - central_offset));
0146   assembly.placeVolume(v_downstream_IP_vacuum_padding,
0147                        Position(0, 0, downstream_straight_length / 2.0));
0148   assembly.placeVolume(v_downstream_IP_gold, Position(0, 0, downstream_straight_length / 2.0));
0149   assembly.placeVolume(v_downstream_IP_tube, Position(0, 0, downstream_straight_length / 2.0));
0150 
0151   auto central_pv = assembly.placeVolume(central_volume, Position(0, 0, +central_offset));
0152   central_det.setPlacement(central_pv);
0153 
0154   // Helper function to create polycone pairs (shell and vacuum)
0155   auto zplane_to_polycones = [](xml::Component& x_pipe) {
0156     std::vector<double> zero, rmax, rmin, z;
0157     for (xml_coll_t x_zplane_i(x_pipe, _Unicode(zplane)); x_zplane_i; ++x_zplane_i) {
0158       xml_comp_t x_zplane = x_zplane_i;
0159       auto thickness      = getAttrOrDefault(x_zplane, _Unicode(thickness), x_pipe.thickness());
0160       thickness += getAttrOrDefault(x_zplane, _Unicode(extra_thickness), 0.0);
0161       zero.push_back(0);
0162       rmax.push_back(x_zplane.attr<double>(_Unicode(OD)) / 2.0);
0163       rmin.push_back(x_zplane.attr<double>(_Unicode(OD)) / 2.0 - thickness);
0164       z.push_back(x_zplane.attr<double>(_Unicode(z)));
0165     }
0166     return std::make_pair<Polycone, Polycone>({0, 2.0 * M_PI, rmin, rmax, z},
0167                                               {0, 2.0 * M_PI, zero, rmin, z});
0168   };
0169 
0170   auto create_volumes = [&](const std::string& name, xml::Component& x_pipe1,
0171                             xml::Component& x_pipe2, xml_coll_t& x_additional_subtraction_i,
0172                             bool subtract_vacuum_from_matter = true,
0173                             bool subtract_matter_from_vacuum = false) {
0174     auto pipe1_polycones = zplane_to_polycones(x_pipe1);
0175     auto pipe2_polycones = zplane_to_polycones(x_pipe2);
0176 
0177     auto crossing_angle    = getAttrOrDefault(x_pipe2, _Unicode(crossing_angle), 0.0);
0178     auto axis_intersection = getAttrOrDefault(x_pipe2, _Unicode(axis_intersection), 0.0);
0179 
0180     auto tf = Transform3D(Position(0, 0, axis_intersection)) *
0181               Transform3D(RotationY(crossing_angle)) *
0182               Transform3D(Position(0, 0, -axis_intersection));
0183 
0184     // union of all matter and vacuum
0185     UnionSolid matter_union(pipe1_polycones.first, pipe2_polycones.first, tf);
0186     UnionSolid vacuum_union(pipe1_polycones.second, pipe2_polycones.second, tf);
0187 
0188     // subtract vacuum from matter
0189     BooleanSolid matter;
0190     if (subtract_vacuum_from_matter) {
0191       matter = SubtractionSolid(matter_union, vacuum_union);
0192     } else {
0193       matter = matter_union;
0194     }
0195     // subtract matter from vacuum
0196     BooleanSolid vacuum;
0197     if (subtract_matter_from_vacuum) {
0198       vacuum = SubtractionSolid(vacuum_union, matter_union);
0199     } else {
0200       vacuum = vacuum_union;
0201     }
0202 
0203     // subtract additional vacuum from matter
0204     for (; x_additional_subtraction_i; ++x_additional_subtraction_i) {
0205       xml_comp_t x_additional_subtraction = x_additional_subtraction_i;
0206       auto additional_polycones           = zplane_to_polycones(x_additional_subtraction);
0207       auto additional_crossing_angle =
0208           getAttrOrDefault(x_additional_subtraction, _Unicode(crossing_angle), 0.0);
0209       auto additional_axis_intersection =
0210           getAttrOrDefault(x_additional_subtraction, _Unicode(axis_intersection), 0.0);
0211       auto additional_tf = Transform3D(Position(0, 0, additional_axis_intersection)) *
0212                            Transform3D(RotationY(additional_crossing_angle)) *
0213                            Transform3D(Position(0, 0, -additional_axis_intersection));
0214       matter = SubtractionSolid(matter, additional_polycones.second, additional_tf);
0215     }
0216 
0217     return std::make_pair<Volume, Volume>({"v_" + name + "_matter", matter, m_Al},
0218                                           {"v_" + name + "_vacuum", vacuum, m_Vacuum});
0219   };
0220 
0221   // -----------------------------
0222   // Upstream:
0223   // - incoming hadron tube: straight section, tapered section, straight section
0224   // - outgoing electron tube: tapered section, straight section
0225 
0226   xml::Component upstream_c        = x_det.child(_Unicode(upstream));
0227   xml::Component incoming_hadron_c = upstream_c.child(_Unicode(incoming_hadron));
0228   xml::Component outgoing_lepton_c = upstream_c.child(_Unicode(outgoing_lepton));
0229   xml_coll_t additional_subtractions_upstream(upstream_c, _Unicode(additional_subtraction));
0230   bool subtract_vacuum_upstream =
0231       getAttrOrDefault<bool>(upstream_c, _Unicode(subtract_vacuum), true);
0232   bool subtract_matter_upstream =
0233       getAttrOrDefault<bool>(upstream_c, _Unicode(subtract_matter), true);
0234   auto volumes_upstream = create_volumes("upstream", outgoing_lepton_c, incoming_hadron_c,
0235                                          additional_subtractions_upstream, subtract_vacuum_upstream,
0236                                          subtract_matter_upstream);
0237 
0238   auto tf_upstream = Transform3D(RotationZYX(0, 0, 0));
0239   if (getAttrOrDefault<bool>(upstream_c, _Unicode(reflect), true)) {
0240     tf_upstream = Transform3D(RotationZYX(0, M_PI, 0));
0241   }
0242   assembly.placeVolume(volumes_upstream.first, tf_upstream);
0243   if (getAttrOrDefault<bool>(upstream_c, _Unicode(place_vacuum), true)) {
0244     assembly.placeVolume(volumes_upstream.second, tf_upstream);
0245   }
0246 
0247   // -----------------------------
0248   // downstream:
0249   // - incoming electron tube: tube with tube cut out
0250   // - outgoing hadron tube: cone centered at scattering angle
0251   // (incoming electron tube internal touching to outgoing hadron tube)
0252 
0253   xml::Component downstream_c      = x_det.child(_Unicode(downstream));
0254   xml::Component incoming_lepton_c = downstream_c.child(_Unicode(incoming_lepton));
0255   xml::Component outgoing_hadron_c = downstream_c.child(_Unicode(outgoing_hadron));
0256   xml_coll_t additional_subtractions_downstream(downstream_c, _Unicode(additional_subtraction));
0257   bool subtract_vacuum_downstream =
0258       getAttrOrDefault<bool>(downstream_c, _Unicode(subtract_vacuum), true);
0259   bool subtract_matter_downstream =
0260       getAttrOrDefault<bool>(downstream_c, _Unicode(subtract_matter), true);
0261   auto volumes_downstream = create_volumes("downstream", incoming_lepton_c, outgoing_hadron_c,
0262                                            additional_subtractions_downstream,
0263                                            subtract_vacuum_downstream, subtract_matter_downstream);
0264 
0265   auto tf_downstream = Transform3D(RotationZYX(0, 0, 0));
0266   if (getAttrOrDefault<bool>(downstream_c, _Unicode(reflect), true)) {
0267     tf_downstream = Transform3D(RotationZYX(0, M_PI, 0));
0268   }
0269   assembly.placeVolume(volumes_downstream.first, tf_downstream);
0270   if (getAttrOrDefault<bool>(downstream_c, _Unicode(place_vacuum), true)) {
0271     assembly.placeVolume(volumes_downstream.second, tf_downstream);
0272   }
0273 
0274   // -----------------------------
0275   // final placement
0276   auto pv_assembly = det.pickMotherVolume(sdet).placeVolume(assembly);
0277   pv_assembly.addPhysVolID("system", sdet.id());
0278   sdet.setPlacement(pv_assembly);
0279   assembly->GetShape()->ComputeBBox();
0280   return sdet;
0281 }
0282 
0283 DECLARE_DETELEMENT(IP6BeamPipe, create_detector)