athena/src/BarrelCalorimeter_geo.cpp

0001 //==========================================================================
0002 //  AIDA Detector description implementation
0003 //--------------------------------------------------------------------------
0004 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
0005 // All rights reserved.
0006 //
0007 // For the licensing terms see $DD4hepINSTALL/LICENSE.
0008 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
0009 //
0010 // Author     : M.Frank
0011 //
0012 //==========================================================================
0013 //
0014 // Specialized generic detector constructor
0015 //
0016 //==========================================================================
0017 #include "DD4hep/DetFactoryHelper.h"
0018 #include "XML/Layering.h"
0019
0020 using namespace std;
0021 using namespace dd4hep;
0022 using namespace dd4hep::detail;
0023
0024 static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens)  {
0025   static double tolerance = 0e0;
0026   Layering      layering (e);
0027   xml_det_t     x_det     = e;
0028   Material      air       = description.air();
0029   int           det_id    = x_det.id();
0030   string        det_name  = x_det.nameStr();
0031   xml_comp_t    x_staves  = x_det.staves();
0032   xml_comp_t    x_dim     = x_det.dimensions();
0033   int           nsides    = x_dim.numsides();
0034   double        inner_r   = x_dim.rmin();
0035   double        dphi      = (2*M_PI/nsides);
0036   double        hphi      = dphi/2;
0037   double        support_thickness = 0.0;
0038   if(x_staves.hasChild("support")){
0039     support_thickness = getAttrOrDefault(x_staves.child(_U(support)), _U(thickness), 5.0 * cm);
0040   }
0041   double        mod_z     = layering.totalThickness() + support_thickness;
0042   double        outer_r   = inner_r + mod_z;
0043   double        totThick  = mod_z;
0044   double        offset    = x_det.attr<double>(_Unicode(offset));
0045   DetElement    sdet      (det_name,det_id);
0046   Volume        motherVol = description.pickMotherVolume(sdet);
0047   PolyhedraRegular hedra  (nsides,inner_r,inner_r+totThick+tolerance*2e0,x_dim.z());
0048   Volume        envelope  (det_name,hedra,air);
0049   PlacedVolume  env_phv   = motherVol.placeVolume(envelope,Transform3D(Translation3D(0,0,offset)*RotationZ(M_PI/nsides)));
0050
0051   env_phv.addPhysVolID("system",det_id);
0052   env_phv.addPhysVolID("barrel",0);
0053   sdet.setPlacement(env_phv);
0054
0055   DetElement    stave_det("stave0",det_id);
0056   double dx = 0.0; //mod_z / std::sin(dphi); // dx per layer
0057
0058   // Compute the top and bottom face measurements.
0059   double trd_x2 = (2 * std::tan(hphi) * outer_r - dx)/2 - tolerance;
0060   double trd_x1 = (2 * std::tan(hphi) * inner_r + dx)/2 - tolerance;
0061   double trd_y1 = x_dim.z()/2 - tolerance;
0062   double trd_y2 = trd_y1;
0063   double trd_z  = mod_z/2 - tolerance;
0064
0065   // Create the trapezoid for the stave.
0066   Trapezoid trd(trd_x1, // Outer side, i.e. the "short" X side.
0067                 trd_x2, // Inner side, i.e. the "long"  X side.
0068                 trd_y1, // Corresponds to subdetector (or module) Z.
0069                 trd_y2, //
0070                 trd_z); // Thickness, in Y for top stave, when rotated.
0071
0072   Volume mod_vol("stave",trd,air);
0073   double l_pos_z = -(layering.totalThickness() / 2) - support_thickness/2.0;
0074
0075   //double trd_x2_support = trd_x1;
0076   double trd_x1_support = (2 * std::tan(hphi) * outer_r - dx- support_thickness)/2 - tolerance;
0077
0078   Solid  support_frame_s;
0079   // optional stave support
0080   if(x_staves.hasChild("support")){
0081     xml_comp_t x_support         = x_staves.child(_U(support));
0082     // is the support on the inside surface?
0083     bool       is_inside_support = getAttrOrDefault<bool>(x_support, _Unicode(inside), true);
0084     // number of "beams" running the length of the stave.
0085     int    n_beams        = getAttrOrDefault<int>(x_support, _Unicode(n_beams), 3);
0086     double beam_thickness = support_thickness / 4.0; // maybe a parameter later...
0087     trd_x1_support        = (2 * std::tan(hphi) * (outer_r - support_thickness + beam_thickness)) / 2 - tolerance;
0088     double grid_size      = getAttrOrDefault(x_support, _Unicode(grid_size), 25.0 * cm);
0089     double beam_width     = 2.0 * trd_x1_support / (n_beams + 1); // quick hack to make some gap between T beams
0090
0091     double cross_beam_thickness    = support_thickness/4.0;
0092     //double trd_x1_support    = (2 * std::tan(hphi) * (inner_r + beam_thickness)) / 2 - tolerance;
0093     double trd_x2_support = trd_x2;
0094
0095     int n_cross_supports = std::floor((trd_y1-cross_beam_thickness)/grid_size);
0096
0097     Box        beam_vert_s(beam_thickness / 2.0 - tolerance, trd_y1, support_thickness / 2.0 - tolerance);
0098     Box        beam_hori_s(beam_width / 2.0 - tolerance, trd_y1, beam_thickness / 2.0 - tolerance);
0099     UnionSolid T_beam_s(beam_vert_s, beam_hori_s, Position(0, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
0100
0101     // cross supports
0102     Trapezoid  trd_support(trd_x1_support,trd_x2_support,
0103                            beam_thickness / 2.0 - tolerance, beam_thickness / 2.0 - tolerance,
0104                           support_thickness / 2.0 - tolerance - cross_beam_thickness/2.0);
0105     UnionSolid support_array_start_s(T_beam_s,trd_support,Position(0,0,cross_beam_thickness/2.0));
0106     for (int isup = 0; isup < n_cross_supports; isup++) {
0107       support_array_start_s = UnionSolid(support_array_start_s, trd_support, Position(0, -1.0 * isup * grid_size, cross_beam_thickness/2.0));
0108       support_array_start_s = UnionSolid(support_array_start_s, trd_support, Position(0, 1.0 * isup * grid_size, cross_beam_thickness/2.0));
0109     }
0110     support_array_start_s =
0111         UnionSolid(support_array_start_s, beam_hori_s,
0112                    Position(-1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
0113     support_array_start_s =
0114         UnionSolid(support_array_start_s, beam_hori_s,
0115                    Position(1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
0116     support_array_start_s =
0117         UnionSolid(support_array_start_s, beam_vert_s, Position(-1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, 0));
0118     support_array_start_s =
0119         UnionSolid(support_array_start_s, beam_vert_s, Position(1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, 0));
0120
0121     support_frame_s = support_array_start_s;
0122
0123     Material support_mat = description.material(x_support.materialStr());
0124     Volume   support_vol("support_frame_v", support_frame_s, support_mat);
0125     support_vol.setVisAttributes(description,x_support.visStr());
0126
0127     // figure out how to best place
0128     //auto pv = mod_vol.placeVolume(support_vol, Position(0.0, 0.0, l_pos_z + support_thickness / 2.0));
0129     auto pv = mod_vol.placeVolume(support_vol, Position(0.0, 0.0, -l_pos_z - support_thickness / 2.0));
0130   }
0131   //l_pos_z += support_thickness;
0132
0133   sens.setType("calorimeter");
0134   { // =====  buildBarrelStave(description, sens, module_volume) =====
0135     // Parameters for computing the layer X dimension:
0136     double stave_z  = trd_y1;
0137     double tan_hphi = std::tan(hphi);
0138     double l_dim_x  = trd_x1; // Starting X dimension for the layer.
0139
0140     // Loop over the sets of layer elements in the detector.
0141     int l_num = 1;
0142     for(xml_coll_t li(x_det,_U(layer)); li; ++li)  {
0143       xml_comp_t x_layer = li;
0144       int repeat = x_layer.repeat();
0145       // Loop over number of repeats for this layer.
0146       for (int j=0; j<repeat; j++)    {
0147         string l_name = _toString(l_num,"layer%d");
0148         double l_thickness = layering.layer(l_num-1)->thickness();  // Layer's thickness.
0149
0150         Position   l_pos(0,0,l_pos_z+l_thickness/2);      // Position of the layer.
0151         Box        l_box(l_dim_x-tolerance,stave_z-tolerance,l_thickness / 2-tolerance);
0152         Volume     l_vol(l_name,l_box,air);
0153         DetElement layer(stave_det, l_name, det_id);
0154
0155         // Loop over the sublayers or slices for this layer.
0156         int s_num = 1;
0157         double s_pos_z = -(l_thickness / 2);
0158         for(xml_coll_t si(x_layer,_U(slice)); si; ++si)  {
0159           xml_comp_t x_slice = si;
0160           string     s_name  = _toString(s_num,"slice%d");
0161           double     s_thick = x_slice.thickness();
0162           Box        s_box(l_dim_x-tolerance,stave_z-tolerance,s_thick / 2-tolerance);
0163           Volume     s_vol(s_name,s_box,description.material(x_slice.materialStr()));
0164           DetElement slice(layer,s_name,det_id);
0165
0166           if ( x_slice.isSensitive() ) {
0167             s_vol.setSensitiveDetector(sens);
0168           }
0169           slice.setAttributes(description,s_vol,x_slice.regionStr(),x_slice.limitsStr(),x_slice.visStr());
0170
0171           // Slice placement.
0172           PlacedVolume slice_phv = l_vol.placeVolume(s_vol,Position(0,0,s_pos_z+s_thick/2));
0173           slice_phv.addPhysVolID("slice", s_num);
0174           slice.setPlacement(slice_phv);
0175           // Increment Z position of slice.
0176           s_pos_z += s_thick;
0177
0178           // Increment slice number.
0179           ++s_num;
0180         }
0181
0182         // Set region, limitset, and vis of layer.
0183         layer.setAttributes(description,l_vol,x_layer.regionStr(),x_layer.limitsStr(),x_layer.visStr());
0184
0185         PlacedVolume layer_phv = mod_vol.placeVolume(l_vol,l_pos);
0186         layer_phv.addPhysVolID("layer", l_num);
0187         layer.setPlacement(layer_phv);
0188         // Increment to next layer Z position.
0189         double xcut = l_thickness * tan_hphi;
0190         l_dim_x += xcut;
0191         l_pos_z += l_thickness;
0192         ++l_num;
0193       }
0194     }
0195   }
0196
0197   // Set stave visualization.
0198   if ( x_staves )   {
0199     mod_vol.setVisAttributes(description.visAttributes(x_staves.visStr()));
0200
0201   }
0202   // Phi start for a stave.
0203   double phi = M_PI / nsides;
0204   double mod_x_off = dx / 2;             // Stave X offset, derived from the dx.
0205   double mod_y_off = inner_r + mod_z/2;  // Stave Y offset
0206
0207   // Create nsides staves.
0208   for (int i = 0; i < nsides; i++, phi -= dphi)      { // i is module number
0209     // Compute the stave position
0210     double m_pos_x = mod_x_off * std::cos(phi) - mod_y_off * std::sin(phi);
0211     double m_pos_y = mod_x_off * std::sin(phi) + mod_y_off * std::cos(phi);
0212     Transform3D tr(RotationZYX(0,phi,M_PI*0.5),Translation3D(-m_pos_x,-m_pos_y,0));
0213     PlacedVolume pv = envelope.placeVolume(mod_vol,tr);
0214     pv.addPhysVolID("system",det_id);
0215     pv.addPhysVolID("barrel",0);
0216     pv.addPhysVolID("module",i+1);
0217     DetElement sd = i==0 ? stave_det : stave_det.clone(_toString(i,"stave%d"));
0218     sd.setPlacement(pv);
0219     sdet.add(sd);
0220   }
0221
0222   // Set envelope volume attributes.
0223   envelope.setAttributes(description,x_det.regionStr(),x_det.limitsStr(),x_det.visStr());
0224   return sdet;
0225 }
0226
0227 DECLARE_DETELEMENT(athena_EcalBarrel,create_detector)
0228 DECLARE_DETELEMENT(athena_HcalBarrel,create_detector)