epic/src/hadronDownstreamBeamPipe.cpp

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2022 Alex Jentsch, Whitney Armstrong, Wouter Deconinck
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
0017 using namespace std;
0018 using namespace dd4hep;
0019
0020 /** \addtogroup beamline Beamline Instrumentation
0021  */
0022
0023 /** \addtogroup IRChamber Interaction Region Vacuum Chamber.
0024  * \brief Type: **IRChamber**.
0025  * \ingroup beamline
0026  *
0027  *
0028  * \code
0029  *   <detector>
0030  *   </detector>
0031  * \endcode
0032  *
0033  */
0034 static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) {
0035
0036   using namespace ROOT::Math;
0037   xml_det_t x_det = e;
0038   string det_name = x_det.nameStr();
0039   // Material   air       = det.air();
0040   DetElement sdet(det_name, x_det.id());
0041   Assembly assembly(det_name + "_assembly");
0042   // Material   m_Cu    = det.material("Copper");
0043   // Material   m_Al    = det.material("Aluminum");
0044   Material m_Be   = det.material("Beryllium");
0045   Material m_SS   = det.material("StainlessSteel");
0046   Material m_vac  = det.material("Vacuum");
0047   string vis_name = x_det.visStr();
0048
0049   PlacedVolume pv_assembly;
0050
0051   //xml::Component pos = x_det.position();
0052   // xml::Component rot   = x_det.rotation();
0053
0054   //int numPipePieces = 4; //number of individual pipe sections
0055
0056   double b0_hadron_tube_inner_r = 2.9;   // cm
0057   double b0_hadron_tube_outer_r = 3.1;   // cm
0058   double b0_hadron_tube_length  = 120.0; // cm
0059
0060   double drift_hadron_section_1_inner_r = 20.0;
0061   double drift_hadron_section_1_outer_r = 20.2;
0062
0063   double drift_hadron_section_3_inner_r_ent = 20.0;
0064   double drift_hadron_section_3_outer_r_ent = 20.2;
0065   double drift_hadron_section_3_inner_r_ex  = 5.0;
0066   double drift_hadron_section_3_outer_r_ex  = 5.2;
0067   double drift_hadron_section_3_length      = 150.0;
0068
0069   double drift_hadron_section_3_x = 0.0;
0070   double drift_hadron_section_3_z = 0.0;
0071
0072   double drift_hadron_section_4_inner_r = 5.0;
0073   double drift_hadron_section_4_outer_r = 5.2;
0074   double drift_hadron_section_4_length  = 850.0;
0075
0076   double drift_hadron_section_4_x = 0.0;
0077   double drift_hadron_section_4_z = 0.0;
0078
0079   //Calculatate full drift region from line formula for proton orbit
0080
0081   //(z, x)
0082   //double orbit_start[2] = {22.0623828, -0.6543372}; //meters!
0083   //double orbit_end[2] = {38.5445361, -1.4039456}; //meters!
0084   //22.07774534
0085   double orbit_start[2] = {22.07774534, -0.650777226}; //meters
0086   double orbit_end[2]   = {38.54362489, -1.436245325}; //meters
0087
0088   //calculate straight line formula x = slope*z + intercept
0089
0090   double slope     = (orbit_end[1] - orbit_start[1]) / (orbit_end[0] - orbit_start[0]);
0091   double intercept = orbit_start[1] - (slope * orbit_start[0]);
0092
0093   // This is the beam tube in the B0 magnet for the hadron beam
0094   // doesn't use the slope information calculated before - it stands alone
0095
0096   Tube b0_hadron_tube(b0_hadron_tube_inner_r, b0_hadron_tube_outer_r, b0_hadron_tube_length / 2.0);
0097   Volume v_b0_hadron_tube("v_b0_hadron_tube", b0_hadron_tube, m_Be);
0098   sdet.setAttributes(det, v_b0_hadron_tube, x_det.regionStr(), x_det.limitsStr(), vis_name);
0099
0100   //----------------------------------
0101   //    build drift beam pipe here
0102   //----------------------------------
0103
0104   double z_start_pipe[3] = {orbit_start[0], 30.000, 31.500};
0105   double z_end_pipe[3]   = {30.000, 31.500, 40.000};
0106
0107   for (int iSection = 0; iSection < 3; iSection++) {
0108
0109     double z_endpoint   = z_end_pipe[iSection]; //meters
0110     double x_endpoint   = (slope * z_endpoint) + intercept;
0111     double x_startpoint = (slope * z_start_pipe[iSection]) + intercept;
0112
0113     double length =
0114         sqrt(pow(z_endpoint - z_start_pipe[iSection], 2) + pow(x_endpoint - x_startpoint, 2));
0115     double z_center = (0.5 * length + z_start_pipe[iSection]) * cos(slope);
0116     double x_center = (slope * z_center) + intercept;
0117
0118     double entrance_r_inner  = 0.0; //drift_hadron_section_1_inner_r;
0119     double exit_radius_inner = 0.0;
0120     double entrance_r_outer  = 0.0; //drift_hadron_section_1_inner_r;
0121     double exit_radius_outer = 0.0;
0122
0123     if (iSection < 1) {
0124       entrance_r_inner  = drift_hadron_section_1_inner_r;
0125       exit_radius_inner = drift_hadron_section_1_inner_r;
0126       entrance_r_outer  = drift_hadron_section_1_outer_r;
0127       exit_radius_outer = drift_hadron_section_1_outer_r;
0128       length            = length - 0.04;
0129     }
0130     if (iSection == 1) {
0131       entrance_r_inner              = drift_hadron_section_3_inner_r_ent;
0132       exit_radius_inner             = drift_hadron_section_3_inner_r_ex;
0133       entrance_r_outer              = drift_hadron_section_3_outer_r_ent;
0134       exit_radius_outer             = drift_hadron_section_3_outer_r_ex;
0135       drift_hadron_section_3_x      = x_center;
0136       drift_hadron_section_3_z      = z_center;
0137       drift_hadron_section_3_length = length - 0.02;
0138       //old numbers commented out for reference - A. Jentsch
0139       //length = length - 0.02;
0140       //x_center = -99.25250431/100.0;
0141       //z_center = 2924.185347/100.0;
0142       //length = drift_hadron_section_3_length/100.0;
0143     }
0144     if (iSection == 2) {
0145       entrance_r_inner              = drift_hadron_section_4_inner_r;
0146       exit_radius_inner             = drift_hadron_section_4_inner_r;
0147       entrance_r_outer              = drift_hadron_section_4_outer_r;
0148       exit_radius_outer             = drift_hadron_section_4_outer_r;
0149       drift_hadron_section_4_x      = x_center;
0150       drift_hadron_section_4_z      = z_center;
0151       drift_hadron_section_4_length = length - 0.02;
0152       length                        = length - 0.02;
0153       //old numbers commented out for reference - A. Jentsch
0154       // x_center =  -123.076799/100.0;
0155       // z_center = 3423.617428/100.0;
0156       //length = drift_hadron_section_4_length/100.0;
0157     }
0158
0159     Cone drift_pipe((length * 100.0) / 2.0, entrance_r_inner, entrance_r_outer, exit_radius_inner,
0160                     exit_radius_outer);
0161
0162     Volume v_pipe(Form("v_drift_tube_pipe_%d", iSection), drift_pipe, m_SS);
0163     sdet.setAttributes(det, v_pipe, x_det.regionStr(), x_det.limitsStr(), vis_name);
0164
0165     auto pv_pipe = assembly.placeVolume(
0166         v_pipe, Transform3D(RotationY(slope),
0167                             Position(100.0 * x_center, 0.0, 100.0 * z_center))); // 2353.06094)));
0168     pv_pipe.addPhysVolID("sector", 1);
0169     DetElement pipe_de(sdet, Form("sector_pipe_%d_de", iSection), 1);
0170     pipe_de.setPlacement(pv_pipe);
0171   }
0172
0173   //------------------------------
0174   //  make vacuum volumes here
0175   //------------------------------
0176
0177   //last two entries are dummy numbers for now
0178   double z_start_points[7]    = {orbit_start[0] + 0.03, 22.590, 24.590, 26.055, 28.055, 20.0, 20.0};
0179   double z_endpoints_array[7] = {22.499, 24.499, 25.980, 27.980, z_start_pipe[1], 25.0, 25.0};
0180
0181   for (int iVac = 0; iVac < 7; iVac++) {
0182
0183     double z_endpoint   = z_endpoints_array[iVac]; //meters
0184     double x_endpoint   = (slope * z_endpoint) + intercept;
0185     double x_startpoint = (slope * z_start_points[iVac]) + intercept;
0186
0187     double length =
0188         sqrt(pow(z_endpoint - z_start_points[iVac], 2) + pow(x_endpoint - x_startpoint, 2));
0189     double z_center = (0.5 * length + z_start_points[iVac]) * cos(slope);
0190     double x_center = (slope * z_center) + intercept;
0191
0192     double entrance_r_inner  = 0.0; //drift_hadron_section_1_inner_r;
0193     double exit_radius_inner = 0.0;
0194
0195     if (iVac < 5) {
0196       entrance_r_inner  = drift_hadron_section_1_inner_r;
0197       exit_radius_inner = drift_hadron_section_1_inner_r;
0198     }
0199     if (iVac == 5) {
0200       entrance_r_inner  = drift_hadron_section_1_inner_r;
0201       exit_radius_inner = drift_hadron_section_3_inner_r_ex;
0202       x_center          = drift_hadron_section_3_x;             //-99.25250431/100.0;
0203       z_center          = drift_hadron_section_3_z;             //2924.185347/100.0;
0204       length            = drift_hadron_section_3_length - 0.02; ///100.0;
0205     }
0206     if (iVac == 6) {
0207       entrance_r_inner  = drift_hadron_section_4_inner_r;
0208       exit_radius_inner = drift_hadron_section_4_inner_r;
0209       //x_center =  -123.076799/100.0;
0210       //z_center = 3423.617428/100.0;
0211       //length = drift_hadron_section_4_length/100.0;
0212       x_center = drift_hadron_section_4_x;
0213       z_center = drift_hadron_section_4_z;
0214       length   = drift_hadron_section_4_length;
0215     }
0216
0217     Cone drift_vacuum((length * 100.0) / 2.0, 0.0, entrance_r_inner - 0.5, 0.0,
0218                       exit_radius_inner - 0.5);
0219
0220     Volume v_vacuum(Form("v_drift_tube_vacuum_%d", iVac), drift_vacuum, m_vac);
0221     sdet.setAttributes(det, v_vacuum, x_det.regionStr(), x_det.limitsStr(), "AnlBlue");
0222
0223     auto pv_vacuum = assembly.placeVolume(
0224         v_vacuum, Transform3D(RotationY(slope),
0225                               Position(100.0 * x_center, 0.0, 100.0 * z_center))); // 2353.06094)));
0226     pv_vacuum.addPhysVolID("sector", 1);
0227     DetElement vacuum_de(sdet, Form("sector_vac_%d_de", iVac), 1);
0228     vacuum_de.setPlacement(pv_vacuum);
0229   }
0230
0231   // Transform3D posAndRot(RotationZYX(rot.z(), rot.y(), rot.x()), Position(pos.x(), pos.y(), pos.z()));
0232   // Transform3D posAndRot(RotationZYX(rot.z(), rot.y(), rot.x()), Position(x_position, y_position, z_position));
0233
0234   pv_assembly = det.pickMotherVolume(sdet).placeVolume(assembly); //, posAndRot);
0235   pv_assembly.addPhysVolID("system", x_det.id()).addPhysVolID("barrel", 1);
0236   sdet.setPlacement(pv_assembly);
0237   assembly->GetShape()->ComputeBBox();
0238   return sdet;
0239 }
0240
0241 DECLARE_DETELEMENT(hadronDownstreamBeamPipe, create_detector)