Back to home page

EIC code displayed by LXR

 
 

    


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

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2022 Simon Gardner
0003 
0004 #include "DD4hep/DetFactoryHelper.h"
0005 #include "DD4hep/OpticalSurfaces.h"
0006 #include "DD4hep/Printout.h"
0007 #include "DDRec/DetectorData.h"
0008 #include "DDRec/Surface.h"
0009 
0010 //////////////////////////////////////////////////
0011 // Far backwards vacuum drift volume
0012 // Low Q2 tagger trackers placed either in or out of vacuum
0013 //////////////////////////////////////////////////
0014 
0015 using namespace std;
0016 using namespace dd4hep;
0017 using namespace dd4hep::rec;
0018 
0019 // Helper function to make the tagger tracker detectors
0020 static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement,
0021                         SensitiveDetector& sens);
0022 
0023 static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) {
0024 
0025   xml_det_t x_det = e;
0026   string detName  = x_det.nameStr();
0027   int detID       = x_det.id();
0028 
0029   DetElement det(detName, detID);
0030 
0031   string vis_name = dd4hep::getAttrOrDefault<std::string>(x_det, _Unicode(vis), "BackwardsBox");
0032 
0033   // Dimensions of main beamline pipe
0034   xml::Component dim = x_det.child(_Unicode(dimensions));
0035   double WidthL      = dim.attr<double>(_Unicode(xL));
0036   double WidthR      = dim.attr<double>(_Unicode(xR));
0037 
0038   double Width     = (WidthL + WidthR) / 2;
0039   double Height    = dim.y();
0040   double Thickness = dim.z();
0041 
0042   // Materials
0043   Material Vacuum = desc.material("Vacuum");
0044   Material Steel  = desc.material("StainlessSteel");
0045 
0046   // Central focal point of the geometry
0047   xml::Component pos = x_det.child(_Unicode(focus));
0048   double off         = pos.z();
0049 
0050   // Beamline rotation
0051   xml_dim_t rot = x_det.rotation();
0052 
0053   // Beampipe thickness
0054   double wall = dd4hep::getAttrOrDefault<double>(x_det, _Unicode(wall), 1 * mm);
0055 
0056   // Make bounding box to make IntersectionSolid with other components
0057   xml::Component BB = x_det.child(_Unicode(bounding));
0058   double BB_MinX    = BB.xmin();
0059   double BB_MinY    = BB.ymin();
0060   double BB_MinZ    = BB.zmin();
0061   double BB_MaxX    = BB.xmax();
0062   double BB_MaxY    = BB.ymax();
0063   double BB_MaxZ    = BB.zmax();
0064 
0065   double BB_X = abs(BB_MaxX - BB_MinX);
0066   double BB_Y = abs(BB_MaxY - BB_MinY);
0067   double BB_Z = abs(BB_MaxZ - BB_MinZ);
0068 
0069   Box Far_Backwards_Box(BB_X, BB_Y, BB_Z);
0070 
0071   // Entry box geometry description joining magnet, taggers and lumi
0072   xml::Component EB = x_det.child(_Unicode(exitdim));
0073   double ED_X       = EB.x();
0074   double ED_Y       = EB.y();
0075   double ED_Z       = off - EB.attr<double>(_Unicode(lumiZ));
0076   double Lumi_R     = EB.attr<double>(_Unicode(lumiR));
0077 
0078   // Maximum theta to exit the dipole from
0079   double exitTheta = EB.attr<double>(_Unicode(maxTheta));
0080 
0081   // Generic box for making intersection solid with
0082   double xbox = 10 * m;
0083   double ybox = 10 * m;
0084   double zbox = 50 * m;
0085 
0086   Box Cut_Box(xbox, ybox, zbox);
0087 
0088   // Central pipe box
0089   //Tube Extended_Beam_Box(Width,Width+wall,Thickness); // More realistic tube pipe
0090   Box Extended_Beam_Box(Width + wall, Height + wall, Thickness); // Simpler box pipe
0091 
0092   // Central vacuum box
0093   //Tube Extended_Vacuum_Box(0,Width,Thickness); // More realistic tube pipe
0094   Box Extended_Vacuum_Box(Width, Height, Thickness); // Simpler box pipe
0095 
0096   Solid Wall_Box   = Extended_Beam_Box;
0097   Solid Vacuum_Box = Extended_Vacuum_Box;
0098 
0099   Assembly DetAssembly("Tagger_vacuum_assembly");
0100   Assembly DetAssemblyAir("Tagger_air_assembly");
0101   int nVacuum = 0;
0102   int nAir    = 0;
0103 
0104   //-----------------------------------------------------------------
0105   // Add Tagger box containers and vacuum box extension for modules
0106   //-----------------------------------------------------------------
0107   for (xml_coll_t mod(x_det, _Unicode(module)); mod; ++mod) {
0108 
0109     int moduleID      = dd4hep::getAttrOrDefault<int>(mod, _Unicode(id), 0);
0110     string moduleName = dd4hep::getAttrOrDefault<std::string>(mod, _Unicode(name), "Tagger0");
0111 
0112     // Offset from the electron beam
0113     double tagoff = dd4hep::getAttrOrDefault<double>(mod, _Unicode(offset_min), 50.0 * mm);
0114 
0115     // Overlap left beyond theta setting
0116     double overlap = dd4hep::getAttrOrDefault<double>(mod, _Unicode(overlap), 0.0 * mm);
0117 
0118     // Theta coverage expected
0119     double thetamin =
0120         dd4hep::getAttrOrDefault<double>(mod, _Unicode(theta_min), 0.030 * rad) - rot.theta();
0121     double thetamax =
0122         dd4hep::getAttrOrDefault<double>(mod, _Unicode(theta_max), 0.030 * rad) - rot.theta();
0123 
0124     // Align box to max or minimum theta expected at the tagger from focal point
0125     bool max_align = dd4hep::getAttrOrDefault<bool>(mod, _Unicode(max_align), false);
0126 
0127     // Extend the beam vacuum around the taggers
0128     bool extend_vacuum = dd4hep::getAttrOrDefault<bool>(mod, _Unicode(extend_vacuum), true);
0129 
0130     // Size f the actual tagger box, replicated in BackwardsTagger
0131     xml_dim_t moddim = mod.child(_Unicode(dimensions));
0132     double w         = moddim.x();
0133     double h         = moddim.y();
0134     double tagboxL   = moddim.z();
0135 
0136     // Width and height of vacuum volume
0137     auto vac_w = w;
0138     auto vac_h = h;
0139 
0140     // Width and height of box volume
0141     auto box_w = w + wall;
0142     auto box_h = h + wall;
0143 
0144     // Angle in relation to the main beam
0145     auto theta = thetamin;
0146 
0147     auto offsetx    = -(box_w - wall) * (cos(theta));
0148     auto offsetz    = (box_w - wall) * (sin(theta));
0149     auto vacoffsetx = -vac_w * (cos(theta));
0150     auto vacoffsetz = vac_w * (sin(theta));
0151     auto l          = (tagoff) / (sin(theta)) + tagboxL;
0152 
0153     auto tagoffsetx = vacoffsetx - (l)*sin(theta);
0154     auto tagoffsetz = vacoffsetz - (l)*cos(theta);
0155 
0156     if (max_align) {
0157       theta      = thetamax;
0158       offsetx    = (overlap + box_w - wall) * (cos(theta));
0159       offsetz    = -(overlap + box_w - wall) * (sin(theta));
0160       vacoffsetx = (overlap + vac_w) * (cos(theta));
0161       vacoffsetz = -(overlap + vac_w) * (sin(theta));
0162       l          = (2 * offsetx + tagoff) / sin(theta);
0163       tagoffsetx = vacoffsetx - (l)*sin(theta);
0164       tagoffsetz = vacoffsetz - (l)*cos(theta);
0165     }
0166 
0167     Box TagWallBox(box_w, box_h, l + wall);
0168     Box TagVacBox(vac_w, vac_h, l);
0169 
0170     RotationY rotate(theta);
0171 
0172     Volume mother = DetAssemblyAir;
0173 
0174     if (extend_vacuum) {
0175       Wall_Box =
0176           UnionSolid(Wall_Box, TagWallBox, Transform3D(rotate, Position(offsetx, 0, offsetz)));
0177       Vacuum_Box = UnionSolid(Vacuum_Box, TagVacBox,
0178                               Transform3D(rotate, Position(vacoffsetx, 0, vacoffsetz)));
0179       mother     = DetAssembly;
0180       nVacuum++;
0181     } else {
0182       nAir++;
0183     }
0184 
0185     Assembly TaggerAssembly("Tagger_module_assembly");
0186 
0187     PlacedVolume pv_mod2 = mother.placeVolume(
0188         TaggerAssembly,
0189         Transform3D(
0190             rotate,
0191             Position(
0192                 tagoffsetx, 0,
0193                 tagoffsetz))); // Very strange y is not centered and offset needs correcting for...
0194     DetElement moddet(det, moduleName, moduleID);
0195     pv_mod2.addPhysVolID("module", moduleID);
0196     moddet.setPlacement(pv_mod2);
0197 
0198     Make_Tagger(desc, mod, TaggerAssembly, moddet, sens);
0199   }
0200 
0201   //-----------------------------------------------------------------
0202   // Cut off any vacuum right of the main beamline
0203   //-----------------------------------------------------------------
0204 
0205   Wall_Box   = IntersectionSolid(Wall_Box, Cut_Box, Position(-xbox + Width + wall, 0, 0));
0206   Vacuum_Box = IntersectionSolid(Vacuum_Box, Cut_Box, Position(-xbox + Width, 0, 0));
0207 
0208   //-----------------------------------------------------------------
0209   // Luminosity connecting box
0210   //-----------------------------------------------------------------
0211   bool addLumi = dd4hep::getAttrOrDefault<bool>(x_det, _Unicode(lumi), true);
0212 
0213   if (addLumi) {
0214 
0215     Box Entry_Beam_Box(ED_X + wall, ED_Y + wall, ED_Z);
0216     Box Entry_Vacuum_Box(ED_X, ED_Y, ED_Z - wall);
0217     Tube Lumi_Exit(0, Lumi_R, ED_Z);
0218 
0219     // Future angled exit window and more realistic tube shaped pipe.
0220     // double angle = -pi/4;
0221     // CutTube Entry_Beam_Box  (ED_X, ED_X + wall, ED_Z,        0,2*pi, sin(angle),0,cos(angle), 0,0,1);
0222     // CutTube Entry_Vacuum_Box(0,    ED_X,        ED_Z - wall, 0,2*pi, sin(angle),0,cos(angle), 0,0,1);
0223     // CutTube Lumi_Exit       (0,    Lumi_R,      ED_Z,        0,2*pi, sin(angle),0,cos(angle), 0,0,1);
0224 
0225     // Add entry boxes to main beamline volume
0226     Wall_Box   = UnionSolid(Wall_Box, Entry_Beam_Box, Transform3D(RotationY(-rot.theta())));
0227     Vacuum_Box = UnionSolid(Vacuum_Box, Entry_Vacuum_Box, Transform3D(RotationY(-rot.theta())));
0228     Vacuum_Box = UnionSolid(Vacuum_Box, Lumi_Exit, Transform3D(RotationY(-rot.theta())));
0229   }
0230 
0231   //-----------------------------------------------------------------
0232   // Restrict tagger boxes into region defined by exitTheta from the dipole magnet
0233   //-----------------------------------------------------------------
0234   double exitDist = BB_MinZ - off;
0235   double cutX     = (ED_X - exitDist * tan(-rot.theta())) * cos(rot.theta());
0236   double cutZ =
0237       (ED_X - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta());
0238   double cutXwall = (ED_X - wall - exitDist * tan(-rot.theta())) * cos(rot.theta());
0239   double cutZwall =
0240       (ED_X - wall - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta());
0241 
0242   Wall_Box = IntersectionSolid(Wall_Box, Cut_Box,
0243                                Transform3D(RotationY(exitTheta), Position(xbox - cutX, 0, cutZ)));
0244   Vacuum_Box =
0245       IntersectionSolid(Vacuum_Box, Cut_Box,
0246                         Transform3D(RotationY(exitTheta), Position(xbox - cutXwall, 0, cutZwall)));
0247 
0248   //-----------------------------------------------------------------
0249   // Cut solids so they are only in the far backwards box
0250   //-----------------------------------------------------------------
0251   RotationY rotate2(-rot.theta());
0252   Position position(0, 0, (exitDist - BB_Z) / cos(rot.theta()));
0253 
0254   IntersectionSolid Wall_Box_Sub(Wall_Box, Far_Backwards_Box, Transform3D(rotate2, position));
0255   IntersectionSolid Vacuum_Box_Sub(Vacuum_Box, Far_Backwards_Box, Transform3D(rotate2, position));
0256   SubtractionSolid Wall_Box_Out(Wall_Box_Sub, Vacuum_Box_Sub);
0257 
0258   Volume vacVol("TaggerStation_Vacuum", Vacuum_Box_Sub, Vacuum);
0259   vacVol.setVisAttributes(desc.visAttributes("BackwardsVac"));
0260   if (nVacuum > 0)
0261     vacVol.placeVolume(DetAssembly);
0262 
0263   Volume wallVol("TaggerStation_Container", Wall_Box_Out, Steel);
0264   wallVol.setVisAttributes(desc.visAttributes(vis_name));
0265 
0266   Assembly backAssembly(detName + "_assembly");
0267   backAssembly.placeVolume(wallVol);
0268   backAssembly.placeVolume(vacVol);
0269 
0270   if (nAir > 0)
0271     backAssembly.placeVolume(DetAssemblyAir);
0272 
0273   // placement in mother volume
0274   Transform3D tr(RotationY(rot.theta()), Position(pos.x(), pos.y(), pos.z()));
0275   PlacedVolume detPV = desc.pickMotherVolume(det).placeVolume(backAssembly, tr);
0276   detPV.addPhysVolID("system", detID);
0277 
0278   det.setPlacement(detPV);
0279 
0280   return det;
0281 }
0282 
0283 static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement,
0284                         SensitiveDetector& sens) {
0285 
0286   sens.setType("tracker");
0287 
0288   Material Air     = desc.material("Air");
0289   Material Silicon = desc.material("Silicon");
0290 
0291   xml_dim_t moddim = mod.child(_Unicode(dimensions));
0292   double tag_w     = moddim.x();
0293   double tag_h     = moddim.y();
0294   double tagboxL   = moddim.z();
0295 
0296   Volume Tagger_Air;
0297 
0298   double airThickness    = 0;
0299   double vacuumThickness = tagboxL;
0300 
0301   // Add window layer and air-vacuum boxes
0302   for (xml_coll_t lay(mod, _Unicode(foilLayer)); lay; ++lay) {
0303 
0304     string layerType = dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(type), "foil");
0305     string layerVis =
0306         dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(vis), "FFTrackerShieldingVis");
0307     double layerZ   = dd4hep::getAttrOrDefault<double>(lay, _Unicode(z), 0 * mm);
0308     double layerRot = dd4hep::getAttrOrDefault<double>(lay, _Unicode(angle), 45 * deg);
0309     double layerThickness =
0310         dd4hep::getAttrOrDefault<double>(lay, _Unicode(sensor_thickness), 100 * um);
0311     string layerMaterial = dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(material), "Copper");
0312 
0313     Material FoilMaterial = desc.material(layerMaterial);
0314 
0315     RotationY rotate(layerRot);
0316 
0317     airThickness    = tagboxL - layerZ;
0318     vacuumThickness = tagboxL - airThickness;
0319 
0320     Box Box_Air(tag_w, tag_h, airThickness / 2);
0321     Tagger_Air = Volume("AirVolume", Box_Air, Air);
0322     Tagger_Air.setVisAttributes(desc.visAttributes("BackwardsAir"));
0323 
0324     Box Foil_Box(tag_w / cos(layerRot) - 0.5 * layerThickness * tan(layerRot), tag_h,
0325                  layerThickness / 2);
0326     Volume layVol("FoilVolume", Foil_Box, FoilMaterial);
0327     layVol.setVisAttributes(desc.visAttributes(layerVis));
0328 
0329     env.placeVolume(layVol, Transform3D(rotate, Position(0, 0, tagboxL + tag_w * tan(layerRot))));
0330 
0331     // Currently only one "foil" layer implemented
0332     break;
0333   }
0334 
0335   // Add window layer and air-vacuum boxes
0336   for (xml_coll_t lay(mod, _Unicode(windowLayer)); lay; ++lay) {
0337 
0338     string layerType = dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(type), "window");
0339     string layerVis =
0340         dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(vis), "FFTrackerShieldingVis");
0341     double layerZ   = dd4hep::getAttrOrDefault<double>(lay, _Unicode(z), 0 * mm);
0342     double layerRot = dd4hep::getAttrOrDefault<double>(lay, _Unicode(angle), 0);
0343     double layerThickness =
0344         dd4hep::getAttrOrDefault<double>(lay, _Unicode(sensor_thickness), 1 * mm);
0345     string layerMaterial = dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(material), "Copper");
0346 
0347     Material WindowMaterial = desc.material(layerMaterial);
0348 
0349     RotationY rotate(layerRot);
0350 
0351     airThickness    = tagboxL - layerZ;
0352     vacuumThickness = tagboxL - airThickness;
0353 
0354     Box Box_Air(tag_w, tag_h, airThickness / 2);
0355     Tagger_Air = Volume("AirVolume", Box_Air, Air);
0356     Tagger_Air.setVisAttributes(desc.visAttributes("BackwardsAir"));
0357 
0358     Box Window_Box(tag_w, tag_h, layerThickness / 2);
0359     Volume layVol("WindowVolume", Window_Box, WindowMaterial);
0360     layVol.setVisAttributes(desc.visAttributes(layerVis));
0361 
0362     Tagger_Air.placeVolume(layVol, Position(0, 0, airThickness / 2 - layerThickness / 2));
0363 
0364     env.placeVolume(Tagger_Air, Position(0, 0, tagboxL - airThickness / 2));
0365 
0366     // Currently only one "window" layer implemented
0367     break;
0368   }
0369 
0370   // Add Hodoscope layers
0371   int N_layers = 0;
0372   for (xml_coll_t lay(mod, _Unicode(trackLayer)); lay; ++lay, ++N_layers) {
0373 
0374     int layerID      = dd4hep::getAttrOrDefault<int>(lay, _Unicode(id), 0);
0375     string layerType = dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(type), "timepix");
0376     string layerVis =
0377         dd4hep::getAttrOrDefault<std::string>(lay, _Unicode(vis), "FFTrackerLayerVis");
0378     double layerRot = dd4hep::getAttrOrDefault<double>(lay, _Unicode(angle), 0);
0379     double layerZ   = dd4hep::getAttrOrDefault<double>(lay, _Unicode(z), 0 * mm);
0380     double layerThickness =
0381         dd4hep::getAttrOrDefault<double>(lay, _Unicode(sensor_thickness), 200 * um);
0382 
0383     Volume mother          = env;
0384     double MotherThickness = tagboxL;
0385 
0386     RotationY rotate(layerRot);
0387 
0388     if (layerZ > vacuumThickness) {
0389       mother = Tagger_Air;
0390       layerZ -= vacuumThickness;
0391       MotherThickness = airThickness / 2;
0392     }
0393 
0394     Box Layer_Box(tag_w, tag_h, layerThickness / 2);
0395     Volume layVol("Tagger_tracker_layer", Layer_Box, Silicon);
0396     layVol.setSensitiveDetector(sens);
0397     layVol.setVisAttributes(desc.visAttributes(layerVis));
0398 
0399     PlacedVolume pv_layer = mother.placeVolume(
0400         layVol, Transform3D(rotate, Position(0, 0, MotherThickness - layerZ + layerThickness / 2)));
0401     pv_layer.addPhysVolID("layer", layerID);
0402 
0403     DetElement laydet(modElement, "layerName" + std::to_string(layerID), layerID);
0404     laydet.setPlacement(pv_layer);
0405   }
0406 }
0407 
0408 DECLARE_DETELEMENT(BackwardsTagger, create_detector)