From 2e5d19c44afe8f3020c538d3af936bbb4d5d1097 Mon Sep 17 00:00:00 2001 From: Wouter Deconinck Date: Tue, 7 May 2024 12:20:59 -0500 Subject: [PATCH] fix: clang-format --- src/BackwardsTaggers_geo.cpp | 197 +++++++------- src/BeamPipeChain_geo.cpp | 147 ++++++----- src/IP6BeamPipe.cpp | 289 ++++++++++----------- src/LumiWindow_geo.cpp | 47 ++-- src/magnetVacuumFF.cpp | 481 ++++++++++++++++++----------------- 5 files changed, 595 insertions(+), 566 deletions(-) diff --git a/src/BackwardsTaggers_geo.cpp b/src/BackwardsTaggers_geo.cpp index 692f44e74..d68882309 100644 --- a/src/BackwardsTaggers_geo.cpp +++ b/src/BackwardsTaggers_geo.cpp @@ -7,7 +7,6 @@ #include "DDRec/DetectorData.h" #include "DDRec/Surface.h" - ////////////////////////////////////////////////// // Far backwards vacuum drift volume // Low Q2 tagger trackers placed either in or out of vacuum @@ -18,36 +17,35 @@ using namespace dd4hep; using namespace dd4hep::rec; // Helper function to make the tagger tracker detectors -static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement, SensitiveDetector& sens); - +static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement, + SensitiveDetector& sens); -static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) -{ +static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) { - xml_det_t x_det = e; - string detName = x_det.nameStr(); - int detID = x_det.id(); + xml_det_t x_det = e; + string detName = x_det.nameStr(); + int detID = x_det.id(); DetElement det(detName, detID); string vis_name = dd4hep::getAttrOrDefault(x_det, _Unicode(vis), "BackwardsBox"); // Dimensions of main beamline pipe - xml::Component dim = x_det.child(_Unicode(dimensions)); - double WidthL = dim.attr(_Unicode(xL)); - double WidthR = dim.attr(_Unicode(xR)); + xml::Component dim = x_det.child(_Unicode(dimensions)); + double WidthL = dim.attr(_Unicode(xL)); + double WidthR = dim.attr(_Unicode(xR)); double Width = (WidthL + WidthR) / 2; double Height = dim.y(); double Thickness = dim.z(); // Materials - Material Vacuum = desc.material("Vacuum"); - Material Copper = desc.material("Copper"); + Material Vacuum = desc.material("Vacuum"); + Material Copper = desc.material("Copper"); // Central focal point of the geometry xml::Component pos = x_det.child(_Unicode(focus)); - double off = pos.z(); + double off = pos.z(); // Beamline rotation xml_dim_t rot = x_det.rotation(); @@ -56,13 +54,13 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) double wall = dd4hep::getAttrOrDefault(x_det, _Unicode(wall), 1 * mm); // Make bounding box to make IntersectionSolid with other components - xml::Component BB = x_det.child(_Unicode(bounding)); - double BB_MinX = BB.xmin(); - double BB_MinY = BB.ymin(); - double BB_MinZ = BB.zmin(); - double BB_MaxX = BB.xmax(); - double BB_MaxY = BB.ymax(); - double BB_MaxZ = BB.zmax(); + xml::Component BB = x_det.child(_Unicode(bounding)); + double BB_MinX = BB.xmin(); + double BB_MinY = BB.ymin(); + double BB_MinZ = BB.zmin(); + double BB_MaxX = BB.xmax(); + double BB_MaxY = BB.ymax(); + double BB_MaxZ = BB.zmax(); double BB_X = abs(BB_MaxX - BB_MinX); double BB_Y = abs(BB_MaxY - BB_MinY); @@ -71,11 +69,11 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) Box Far_Backwards_Box(BB_X, BB_Y, BB_Z); // Entry box geometry description joining magnet, taggers and lumi - xml::Component EB = x_det.child(_Unicode(exitdim)); - double ED_X = EB.x(); - double ED_Y = EB.y(); - double ED_Z = off - EB.attr(_Unicode(lumiZ)); - double Lumi_R = EB.attr(_Unicode(lumiR)); + xml::Component EB = x_det.child(_Unicode(exitdim)); + double ED_X = EB.x(); + double ED_Y = EB.y(); + double ED_Z = off - EB.attr(_Unicode(lumiZ)); + double Lumi_R = EB.attr(_Unicode(lumiR)); // Maximum theta to exit the dipole from double exitTheta = EB.attr(_Unicode(maxTheta)); @@ -100,26 +98,28 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) Assembly DetAssembly("Tagger_vacuum_assembly"); Assembly DetAssemblyAir("Tagger_air_assembly"); - int nVacuum = 0; - int nAir = 0; + int nVacuum = 0; + int nAir = 0; //----------------------------------------------------------------- // Add Tagger box containers and vacuum box extension for modules //----------------------------------------------------------------- for (xml_coll_t mod(x_det, _Unicode(module)); mod; ++mod) { - int moduleID = dd4hep::getAttrOrDefault(mod, _Unicode(id), 0); + int moduleID = dd4hep::getAttrOrDefault(mod, _Unicode(id), 0); string moduleName = dd4hep::getAttrOrDefault(mod, _Unicode(name), "Tagger0"); // Offset from the electron beam - double tagoff = dd4hep::getAttrOrDefault(mod, _Unicode(offset_min), 50.0 * mm); + double tagoff = dd4hep::getAttrOrDefault(mod, _Unicode(offset_min), 50.0 * mm); // Overlap left beyond theta setting double overlap = dd4hep::getAttrOrDefault(mod, _Unicode(overlap), 0.0 * mm); // Theta coverage expected - double thetamin = dd4hep::getAttrOrDefault(mod, _Unicode(theta_min), 0.030 * rad) - rot.theta(); - double thetamax = dd4hep::getAttrOrDefault(mod, _Unicode(theta_max), 0.030 * rad) - rot.theta(); + double thetamin = + dd4hep::getAttrOrDefault(mod, _Unicode(theta_min), 0.030 * rad) - rot.theta(); + double thetamax = + dd4hep::getAttrOrDefault(mod, _Unicode(theta_max), 0.030 * rad) - rot.theta(); // Align box to max or minimum theta expected at the tagger from focal point bool max_align = dd4hep::getAttrOrDefault(mod, _Unicode(max_align), false); @@ -128,10 +128,10 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) bool extend_vacuum = dd4hep::getAttrOrDefault(mod, _Unicode(extend_vacuum), true); // Size f the actual tagger box, replicated in BackwardsTagger - xml_dim_t moddim = mod.child(_Unicode(dimensions)); - double w = moddim.x(); - double h = moddim.y(); - double tagboxL = moddim.z(); + xml_dim_t moddim = mod.child(_Unicode(dimensions)); + double w = moddim.x(); + double h = moddim.y(); + double tagboxL = moddim.z(); // Width and height of vacuum volume auto vac_w = w; @@ -142,8 +142,7 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) auto box_h = h + wall; // Angle in relation to the main beam - auto theta = thetamin; - + auto theta = thetamin; auto offsetx = -(box_w - wall) * (cos(theta)); auto offsetz = (box_w - wall) * (sin(theta)); @@ -151,18 +150,18 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) auto vacoffsetz = vac_w * (sin(theta)); auto l = (tagoff) / (sin(theta)) + tagboxL; - auto tagoffsetx = vacoffsetx - (l) * sin(theta); - auto tagoffsetz = vacoffsetz - (l) * cos(theta); + auto tagoffsetx = vacoffsetx - (l)*sin(theta); + auto tagoffsetz = vacoffsetz - (l)*cos(theta); if (max_align) { theta = thetamax; - offsetx = (overlap+box_w - wall) * (cos(theta)); - offsetz = -(overlap+box_w - wall) * (sin(theta)); - vacoffsetx = (overlap+vac_w) * (cos(theta)); - vacoffsetz = -(overlap+vac_w) * (sin(theta)); + offsetx = (overlap + box_w - wall) * (cos(theta)); + offsetz = -(overlap + box_w - wall) * (sin(theta)); + vacoffsetx = (overlap + vac_w) * (cos(theta)); + vacoffsetz = -(overlap + vac_w) * (sin(theta)); l = (2 * offsetx + tagoff) / sin(theta); - tagoffsetx = vacoffsetx - (l) * sin(theta); - tagoffsetz = vacoffsetz - (l) * cos(theta); + tagoffsetx = vacoffsetx - (l)*sin(theta); + tagoffsetz = vacoffsetz - (l)*cos(theta); } Box TagWallBox(box_w, box_h, l + wall); @@ -173,8 +172,10 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) Volume mother = DetAssemblyAir; if (extend_vacuum) { - Wall_Box = UnionSolid(Wall_Box, TagWallBox, Transform3D(rotate, Position(offsetx, 0, offsetz))); - Vacuum_Box = UnionSolid(Vacuum_Box, TagVacBox, Transform3D(rotate, Position(vacoffsetx, 0, vacoffsetz))); + Wall_Box = + UnionSolid(Wall_Box, TagWallBox, Transform3D(rotate, Position(offsetx, 0, offsetz))); + Vacuum_Box = UnionSolid(Vacuum_Box, TagVacBox, + Transform3D(rotate, Position(vacoffsetx, 0, vacoffsetz))); mother = DetAssembly; nVacuum++; } else { @@ -185,14 +186,16 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) PlacedVolume pv_mod2 = mother.placeVolume( TaggerAssembly, - Transform3D(rotate, Position(tagoffsetx, 0, - tagoffsetz))); // Very strange y is not centered and offset needs correcting for... - DetElement moddet(det,moduleName, moduleID); + Transform3D( + rotate, + Position( + tagoffsetx, 0, + tagoffsetz))); // Very strange y is not centered and offset needs correcting for... + DetElement moddet(det, moduleName, moduleID); pv_mod2.addPhysVolID("module", moduleID); moddet.setPlacement(pv_mod2); Make_Tagger(desc, mod, TaggerAssembly, moddet, sens); - } //----------------------------------------------------------------- @@ -230,23 +233,27 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) //----------------------------------------------------------------- double exitDist = BB_MinZ - off; double cutX = (ED_X - exitDist * tan(-rot.theta())) * cos(rot.theta()); - double cutZ = (ED_X - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta()); + double cutZ = + (ED_X - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta()); double cutXwall = (ED_X - wall - exitDist * tan(-rot.theta())) * cos(rot.theta()); - double cutZwall = (ED_X - wall - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta()); + double cutZwall = + (ED_X - wall - exitDist * tan(-rot.theta())) * sin(rot.theta()) + exitDist * cos(rot.theta()); - Wall_Box = IntersectionSolid(Wall_Box, Cut_Box, Transform3D(RotationY(exitTheta), Position(xbox - cutX, 0, cutZ))); + Wall_Box = IntersectionSolid(Wall_Box, Cut_Box, + Transform3D(RotationY(exitTheta), Position(xbox - cutX, 0, cutZ))); Vacuum_Box = - IntersectionSolid(Vacuum_Box, Cut_Box, Transform3D(RotationY(exitTheta), Position(xbox - cutXwall, 0, cutZwall))); + IntersectionSolid(Vacuum_Box, Cut_Box, + Transform3D(RotationY(exitTheta), Position(xbox - cutXwall, 0, cutZwall))); //----------------------------------------------------------------- // Cut solids so they are only in the far backwards box //----------------------------------------------------------------- RotationY rotate2(-rot.theta()); - Position position(0, 0, (exitDist - BB_Z) / cos(rot.theta())); + Position position(0, 0, (exitDist - BB_Z) / cos(rot.theta())); IntersectionSolid Wall_Box_Sub(Wall_Box, Far_Backwards_Box, Transform3D(rotate2, position)); IntersectionSolid Vacuum_Box_Sub(Vacuum_Box, Far_Backwards_Box, Transform3D(rotate2, position)); - SubtractionSolid Wall_Box_Out(Wall_Box_Sub, Vacuum_Box_Sub); + SubtractionSolid Wall_Box_Out(Wall_Box_Sub, Vacuum_Box_Sub); Volume vacVol("TaggerStation_Vacuum", Vacuum_Box_Sub, Vacuum); vacVol.setVisAttributes(desc.visAttributes("BackwardsVac")); @@ -264,28 +271,27 @@ static Ref_t create_detector(Detector& desc, xml_h e, SensitiveDetector sens) backAssembly.placeVolume(DetAssemblyAir); // placement in mother volume - Transform3D tr(RotationY(rot.theta()), Position(pos.x(), pos.y(), pos.z())); + Transform3D tr(RotationY(rot.theta()), Position(pos.x(), pos.y(), pos.z())); PlacedVolume detPV = desc.pickMotherVolume(det).placeVolume(backAssembly, tr); detPV.addPhysVolID("system", detID); det.setPlacement(detPV); - return det; } -static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement, SensitiveDetector& sens) -{ +static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetElement modElement, + SensitiveDetector& sens) { sens.setType("tracker"); Material Air = desc.material("Air"); Material Silicon = desc.material("Silicon"); - xml_dim_t moddim = mod.child(_Unicode(dimensions)); - double tag_w = moddim.x(); - double tag_h = moddim.y(); - double tagboxL = moddim.z(); + xml_dim_t moddim = mod.child(_Unicode(dimensions)); + double tag_w = moddim.x(); + double tag_h = moddim.y(); + double tagboxL = moddim.z(); Volume Tagger_Air; @@ -295,12 +301,14 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme // Add window layer and air-vacuum boxes for (xml_coll_t lay(mod, _Unicode(foilLayer)); lay; ++lay) { - string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "foil"); - string layerVis = dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerShieldingVis"); - double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); - double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 45*deg); - double layerThickness = dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 100 * um); - string layerMaterial = dd4hep::getAttrOrDefault(lay, _Unicode(material), "Copper"); + string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "foil"); + string layerVis = + dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerShieldingVis"); + double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); + double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 45 * deg); + double layerThickness = + dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 100 * um); + string layerMaterial = dd4hep::getAttrOrDefault(lay, _Unicode(material), "Copper"); Material FoilMaterial = desc.material(layerMaterial); @@ -313,11 +321,12 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme Tagger_Air = Volume("AirVolume", Box_Air, Air); Tagger_Air.setVisAttributes(desc.visAttributes("BackwardsAir")); - Box Foil_Box(tag_w/cos(layerRot)-0.5*layerThickness*tan(layerRot), tag_h, layerThickness / 2); + Box Foil_Box(tag_w / cos(layerRot) - 0.5 * layerThickness * tan(layerRot), tag_h, + layerThickness / 2); Volume layVol("FoilVolume", Foil_Box, FoilMaterial); layVol.setVisAttributes(desc.visAttributes(layerVis)); - env.placeVolume(layVol, Transform3D(rotate,Position(0, 0, tagboxL+tag_w*tan(layerRot)))); + env.placeVolume(layVol, Transform3D(rotate, Position(0, 0, tagboxL + tag_w * tan(layerRot)))); // Currently only one "foil" layer implemented break; @@ -326,12 +335,14 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme // Add window layer and air-vacuum boxes for (xml_coll_t lay(mod, _Unicode(windowLayer)); lay; ++lay) { - string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "window"); - string layerVis = dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerShieldingVis"); - double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); - double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 0); - double layerThickness = dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 1 * mm); - string layerMaterial = dd4hep::getAttrOrDefault(lay, _Unicode(material), "Copper"); + string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "window"); + string layerVis = + dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerShieldingVis"); + double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); + double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 0); + double layerThickness = + dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 1 * mm); + string layerMaterial = dd4hep::getAttrOrDefault(lay, _Unicode(material), "Copper"); Material WindowMaterial = desc.material(layerMaterial); @@ -344,7 +355,7 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme Tagger_Air = Volume("AirVolume", Box_Air, Air); Tagger_Air.setVisAttributes(desc.visAttributes("BackwardsAir")); - Box Window_Box(tag_w, tag_h, layerThickness / 2); + Box Window_Box(tag_w, tag_h, layerThickness / 2); Volume layVol("WindowVolume", Window_Box, WindowMaterial); layVol.setVisAttributes(desc.visAttributes(layerVis)); @@ -360,12 +371,14 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme int N_layers = 0; for (xml_coll_t lay(mod, _Unicode(trackLayer)); lay; ++lay, ++N_layers) { - int layerID = dd4hep::getAttrOrDefault(lay, _Unicode(id), 0); - string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "timepix"); - string layerVis = dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerLayerVis"); - double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 0); - double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); - double layerThickness = dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 200 * um); + int layerID = dd4hep::getAttrOrDefault(lay, _Unicode(id), 0); + string layerType = dd4hep::getAttrOrDefault(lay, _Unicode(type), "timepix"); + string layerVis = + dd4hep::getAttrOrDefault(lay, _Unicode(vis), "FFTrackerLayerVis"); + double layerRot = dd4hep::getAttrOrDefault(lay, _Unicode(angle), 0); + double layerZ = dd4hep::getAttrOrDefault(lay, _Unicode(z), 0 * mm); + double layerThickness = + dd4hep::getAttrOrDefault(lay, _Unicode(sensor_thickness), 200 * um); Volume mother = env; double MotherThickness = tagboxL; @@ -378,19 +391,17 @@ static void Make_Tagger(Detector& desc, xml_coll_t& mod, Assembly& env, DetEleme MotherThickness = airThickness / 2; } - Box Layer_Box(tag_w, tag_h, layerThickness / 2); + Box Layer_Box(tag_w, tag_h, layerThickness / 2); Volume layVol("Tagger_tracker_layer", Layer_Box, Silicon); layVol.setSensitiveDetector(sens); layVol.setVisAttributes(desc.visAttributes(layerVis)); - - PlacedVolume pv_layer = mother.placeVolume(layVol, Transform3D(rotate, Position(0, 0, MotherThickness - layerZ + layerThickness / 2))); + PlacedVolume pv_layer = mother.placeVolume( + layVol, Transform3D(rotate, Position(0, 0, MotherThickness - layerZ + layerThickness / 2))); pv_layer.addPhysVolID("layer", layerID); - DetElement laydet(modElement,"layerName"+std::to_string(layerID), layerID); + DetElement laydet(modElement, "layerName" + std::to_string(layerID), layerID); laydet.setPlacement(pv_layer); - - } } diff --git a/src/BeamPipeChain_geo.cpp b/src/BeamPipeChain_geo.cpp index abec2b18b..7213cfec6 100644 --- a/src/BeamPipeChain_geo.cpp +++ b/src/BeamPipeChain_geo.cpp @@ -18,20 +18,19 @@ using namespace std; using namespace dd4hep; -static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector /* sens */) -{ +static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector /* sens */) { using namespace ROOT::Math; - xml_det_t x_det = e; - string det_name = x_det.nameStr(); + xml_det_t x_det = e; + string det_name = x_det.nameStr(); DetElement sdet(det_name, x_det.id()); - Assembly assembly(det_name + "_assembly"); - Material m_SS = description.material("StainlessSteel"); - Material m_Cu = description.material("Copper"); - Material m_Vacuum = description.material("Vacuum"); - string vis_name = dd4hep::getAttrOrDefault(x_det, _Unicode(vis), "BeamPipeVis"); - double thickness = getAttrOrDefault(x_det, _Unicode(wall_thickness), 0); - double coating = getAttrOrDefault(x_det, _Unicode(coating_thickness), 0); + Assembly assembly(det_name + "_assembly"); + Material m_SS = description.material("StainlessSteel"); + Material m_Cu = description.material("Copper"); + Material m_Vacuum = description.material("Vacuum"); + string vis_name = dd4hep::getAttrOrDefault(x_det, _Unicode(vis), "BeamPipeVis"); + double thickness = getAttrOrDefault(x_det, _Unicode(wall_thickness), 0); + double coating = getAttrOrDefault(x_det, _Unicode(coating_thickness), 0); // beam pipe vector names; @@ -43,89 +42,101 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector / vector rOuters2; // Grab info for beamline magnets - for( xml_coll_t pipe_coll(x_det, _Unicode(pipe)); pipe_coll; pipe_coll++ ) { // pipes + for (xml_coll_t pipe_coll(x_det, _Unicode(pipe)); pipe_coll; pipe_coll++) { // pipes - xml_comp_t pipe( pipe_coll ); + xml_comp_t pipe(pipe_coll); - names.push_back( getAttrOrDefault(pipe, _Unicode(name), "") ); + names.push_back(getAttrOrDefault(pipe, _Unicode(name), "")); // Vectors momentarily filled with zeros for pipes in between magnets - xCenters.push_back( getAttrOrDefault(pipe, _Unicode(xcenter), 0) ); - zCenters.push_back( getAttrOrDefault(pipe, _Unicode(zcenter), 0) ); - lengths.push_back( getAttrOrDefault(pipe, _Unicode(length), 0) ); - thetas.push_back( getAttrOrDefault(pipe, _Unicode(theta), 0) ); - rOuters1.push_back( getAttrOrDefault(pipe, _Unicode(rout1), 0) ); - rOuters2.push_back( getAttrOrDefault(pipe, _Unicode(rout2), 0) ); + xCenters.push_back(getAttrOrDefault(pipe, _Unicode(xcenter), 0)); + zCenters.push_back(getAttrOrDefault(pipe, _Unicode(zcenter), 0)); + lengths.push_back(getAttrOrDefault(pipe, _Unicode(length), 0)); + thetas.push_back(getAttrOrDefault(pipe, _Unicode(theta), 0)); + rOuters1.push_back(getAttrOrDefault(pipe, _Unicode(rout1), 0)); + rOuters2.push_back(getAttrOrDefault(pipe, _Unicode(rout2), 0)); } // Calculate parameters for connecting pipes in between magnets - for( uint pipeN = 0; pipeN < names.size(); pipeN++ ) { - - if( lengths[pipeN] > 0 ) { continue; } // pipe parameters already set to nonzero values - if( pipeN == 0 ) { continue; } // can't create pipe for an empty starting slot - if( (pipeN+1) == names.size() ) { continue; } // can't create pipe for an empty end slot - - double x = ( xCenters[pipeN-1] - lengths[pipeN-1]/2.*sin(thetas[pipeN-1]) + xCenters[pipeN+1] + lengths[pipeN+1]/2.*sin(thetas[pipeN+1]) ) / 2.; - double z = ( zCenters[pipeN-1] - lengths[pipeN-1]/2.*cos(thetas[pipeN-1]) + zCenters[pipeN+1] + lengths[pipeN+1]/2.*cos(thetas[pipeN+1]) ) / 2.; - double deltaX = (xCenters[pipeN-1] - lengths[pipeN-1]/2.*sin(thetas[pipeN-1])) - (xCenters[pipeN+1] + lengths[pipeN+1]/2.*sin(thetas[pipeN+1])); - double deltaZ = (zCenters[pipeN-1] - lengths[pipeN-1]/2.*cos(thetas[pipeN-1])) - (zCenters[pipeN+1] + lengths[pipeN+1]/2.*cos(thetas[pipeN+1])); - double l = sqrt( pow(deltaX, 2) + pow(deltaZ, 2) ); - double theta = atan( deltaX / deltaZ ); + for (uint pipeN = 0; pipeN < names.size(); pipeN++) { + + if (lengths[pipeN] > 0) { + continue; + } // pipe parameters already set to nonzero values + if (pipeN == 0) { + continue; + } // can't create pipe for an empty starting slot + if ((pipeN + 1) == names.size()) { + continue; + } // can't create pipe for an empty end slot + + double x = (xCenters[pipeN - 1] - lengths[pipeN - 1] / 2. * sin(thetas[pipeN - 1]) + + xCenters[pipeN + 1] + lengths[pipeN + 1] / 2. * sin(thetas[pipeN + 1])) / + 2.; + double z = (zCenters[pipeN - 1] - lengths[pipeN - 1] / 2. * cos(thetas[pipeN - 1]) + + zCenters[pipeN + 1] + lengths[pipeN + 1] / 2. * cos(thetas[pipeN + 1])) / + 2.; + double deltaX = (xCenters[pipeN - 1] - lengths[pipeN - 1] / 2. * sin(thetas[pipeN - 1])) - + (xCenters[pipeN + 1] + lengths[pipeN + 1] / 2. * sin(thetas[pipeN + 1])); + double deltaZ = (zCenters[pipeN - 1] - lengths[pipeN - 1] / 2. * cos(thetas[pipeN - 1])) - + (zCenters[pipeN + 1] + lengths[pipeN + 1] / 2. * cos(thetas[pipeN + 1])); + double l = sqrt(pow(deltaX, 2) + pow(deltaZ, 2)); + double theta = atan(deltaX / deltaZ); // Small air gap between connecting and magnet beam pipes to avoid G4 overlap errors - if( (theta != thetas[pipeN-1]) || (theta != thetas[pipeN+1]) ) { + if ((theta != thetas[pipeN - 1]) || (theta != thetas[pipeN + 1])) { l -= 0.5; } xCenters[pipeN] = x; zCenters[pipeN] = z; - lengths[pipeN] = l; - thetas[pipeN] = theta; - rOuters1[pipeN] = rOuters2[pipeN-1]; - rOuters2[pipeN] = rOuters1[pipeN+1]; + lengths[pipeN] = l; + thetas[pipeN] = theta; + rOuters1[pipeN] = rOuters2[pipeN - 1]; + rOuters2[pipeN] = rOuters1[pipeN + 1]; } // Add all pipes to the assembly - for( uint pipeN = 0; pipeN < xCenters.size(); pipeN++ ) { + for (uint pipeN = 0; pipeN < xCenters.size(); pipeN++) { // beam pipe matter - ConeSegment s_tube( - lengths[pipeN] / 2.0, // length - rOuters2[pipeN] - thickness, // r2 in - rOuters2[pipeN], // r2 out - rOuters1[pipeN] - thickness, // r1 in - rOuters1[pipeN] // r1 out - ); + ConeSegment s_tube(lengths[pipeN] / 2.0, // length + rOuters2[pipeN] - thickness, // r2 in + rOuters2[pipeN], // r2 out + rOuters1[pipeN] - thickness, // r1 in + rOuters1[pipeN] // r1 out + ); // beam pipe coating - ConeSegment s_coat( - lengths[pipeN] / 2.0, // length - rOuters2[pipeN] - thickness - coating, // r2 in - rOuters2[pipeN] - thickness, // r2 out - rOuters1[pipeN] - thickness - coating, // r1 in - rOuters1[pipeN] - thickness // r1 out - ); + ConeSegment s_coat(lengths[pipeN] / 2.0, // length + rOuters2[pipeN] - thickness - coating, // r2 in + rOuters2[pipeN] - thickness, // r2 out + rOuters1[pipeN] - thickness - coating, // r1 in + rOuters1[pipeN] - thickness // r1 out + ); // beam pipe vacuum - ConeSegment s_vacm( - lengths[pipeN] / 2.0, // length - 0, // r2 in - rOuters2[pipeN] - thickness - coating, // r2 out - 0, // r1 in - rOuters1[pipeN] - thickness - coating // r1 out - ); - - Volume v_tube("v_tube_" + names[pipeN], s_tube, m_SS); + ConeSegment s_vacm(lengths[pipeN] / 2.0, // length + 0, // r2 in + rOuters2[pipeN] - thickness - coating, // r2 out + 0, // r1 in + rOuters1[pipeN] - thickness - coating // r1 out + ); + + Volume v_tube("v_tube_" + names[pipeN], s_tube, m_SS); Volume v_coat("v_coating_" + names[pipeN], s_coat, m_Cu); - Volume v_vacm("v_vacuum_" + names[pipeN], s_vacm, m_Vacuum); + Volume v_vacm("v_vacuum_" + names[pipeN], s_vacm, m_Vacuum); - v_tube.setVisAttributes(description.visAttributes( vis_name ) ); - v_coat.setVisAttributes(description.visAttributes( vis_name + "Coating") ); + v_tube.setVisAttributes(description.visAttributes(vis_name)); + v_coat.setVisAttributes(description.visAttributes(vis_name + "Coating")); - assembly.placeVolume(v_tube, Transform3D( RotationY(thetas[pipeN]), Position(xCenters[pipeN], 0, zCenters[pipeN]))); - assembly.placeVolume(v_coat, Transform3D( RotationY(thetas[pipeN]), Position(xCenters[pipeN], 0, zCenters[pipeN]))); - assembly.placeVolume(v_vacm, Transform3D( RotationY(thetas[pipeN]), Position(xCenters[pipeN], 0, zCenters[pipeN]))); + assembly.placeVolume(v_tube, Transform3D(RotationY(thetas[pipeN]), + Position(xCenters[pipeN], 0, zCenters[pipeN]))); + assembly.placeVolume(v_coat, Transform3D(RotationY(thetas[pipeN]), + Position(xCenters[pipeN], 0, zCenters[pipeN]))); + assembly.placeVolume(v_vacm, Transform3D(RotationY(thetas[pipeN]), + Position(xCenters[pipeN], 0, zCenters[pipeN]))); } // Final placement auto pv_assembly = - description.pickMotherVolume(sdet).placeVolume( assembly, Position(0.0, 0.0, 0.0)); + description.pickMotherVolume(sdet).placeVolume(assembly, Position(0.0, 0.0, 0.0)); sdet.setPlacement(pv_assembly); diff --git a/src/IP6BeamPipe.cpp b/src/IP6BeamPipe.cpp index 55f118f49..b2e163fd7 100644 --- a/src/IP6BeamPipe.cpp +++ b/src/IP6BeamPipe.cpp @@ -33,22 +33,21 @@ using namespace dd4hep; * \endcode * */ -static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) -{ +static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) { using namespace ROOT::Math; - xml_det_t x_det = e; - string det_name = x_det.nameStr(); - xml_comp_t x_dettype = x_det.child(dd4hep::xml::Strng_t("type_flags")); + xml_det_t x_det = e; + string det_name = x_det.nameStr(); + xml_comp_t x_dettype = x_det.child(dd4hep::xml::Strng_t("type_flags")); unsigned int typeFlag = x_dettype.type(); DetElement sdet(det_name, x_det.id()); - Assembly assembly(det_name + "_assembly"); - Material m_SS = det.material("StainlessSteel"); - Material m_Cu = det.material("Copper"); - Material m_Be = det.material("Beryllium"); - Material m_Au = det.material("Gold"); - Material m_Vacuum = det.material("Vacuum"); - string vis_name = x_det.visStr(); + Assembly assembly(det_name + "_assembly"); + Material m_SS = det.material("StainlessSteel"); + Material m_Cu = det.material("Copper"); + Material m_Be = det.material("Beryllium"); + Material m_Au = det.material("Gold"); + Material m_Vacuum = det.material("Vacuum"); + string vis_name = x_det.visStr(); xml::Component IP_pipe_c = x_det.child(_Unicode(IP_pipe)); @@ -56,7 +55,8 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * double IP_beampipe_OD = IP_pipe_c.attr(_Unicode(OD)); double IP_beampipe_wall_thickness = IP_pipe_c.attr(_Unicode(wall_thickness)); double IP_beampipe_gold_thickness = IP_pipe_c.attr(_Unicode(gold_thickness)); - double IP_beampipe_ID = IP_beampipe_OD - 2.0 * IP_beampipe_gold_thickness - 2.0 * IP_beampipe_wall_thickness; + double IP_beampipe_ID = + IP_beampipe_OD - 2.0 * IP_beampipe_gold_thickness - 2.0 * IP_beampipe_wall_thickness; double IP_acts_beampipe_OD = IP_beampipe_ID - 5.0 * mm; double IP_acts_beampipe_ID = IP_acts_beampipe_OD - 1.0 * mm; @@ -64,17 +64,17 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * double downstream_straight_length = IP_pipe_c.attr(_Unicode(downstream_straight_length)); // central beampipe volume - Tube central_tube(0.5 * IP_acts_beampipe_ID, 0.5 * IP_acts_beampipe_OD, - 0.5 * (upstream_straight_length + downstream_straight_length)); - Volume central_volume("acts_central_beampipe_vol", central_tube, m_Vacuum); + Tube central_tube(0.5 * IP_acts_beampipe_ID, 0.5 * IP_acts_beampipe_OD, + 0.5 * (upstream_straight_length + downstream_straight_length)); + Volume central_volume("acts_central_beampipe_vol", central_tube, m_Vacuum); const double central_offset = -.5 * (upstream_straight_length - downstream_straight_length); - DetElement central_det(sdet, "acts_beampipe_central", 1); + DetElement central_det(sdet, "acts_beampipe_central", 1); // Set dd4hep variant parameters for conversion to ACTS tracking geometry central_det.setTypeFlag(typeFlag); - auto ¶ms = DD4hepDetectorHelper::ensureExtension(central_det); - int nBinPhi = 144; // fix later. Should take this from a xml tag - int nBinZ = 10; // fix later. Should take this from a xml tag + auto& params = DD4hepDetectorHelper::ensureExtension(central_det); + int nBinPhi = 144; // fix later. Should take this from a xml tag + int nBinZ = 10; // fix later. Should take this from a xml tag params.set("layer_material", true); params.set("layer_material_representing", true); params.set("layer_material_representing_binPhi", nBinPhi); @@ -94,28 +94,36 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * Tube downstream_IP_vacuum_fill(0.0, IP_acts_beampipe_ID / 2.0, downstream_straight_length / 2.0); Tube downstream_IP_acts_beampipe(IP_acts_beampipe_ID / 2.0, IP_acts_beampipe_OD / 2.0, downstream_straight_length / 2.0); - Tube downstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0, downstream_straight_length / 2.0); + Tube downstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0, + downstream_straight_length / 2.0); Tube downstream_IP_gold(IP_beampipe_ID / 2.0, IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, downstream_straight_length / 2.0); Tube downstream_IP_tube(IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, IP_beampipe_OD / 2.0, downstream_straight_length / 2.0); Tube upstream_IP_vacuum_fill(0.0, IP_acts_beampipe_ID / 2.0, upstream_straight_length / 2.0); - Tube upstream_IP_acts_beampipe(IP_acts_beampipe_ID / 2.0, IP_acts_beampipe_OD / 2.0, upstream_straight_length / 2.0); - Tube upstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0, upstream_straight_length / 2.0); + Tube upstream_IP_acts_beampipe(IP_acts_beampipe_ID / 2.0, IP_acts_beampipe_OD / 2.0, + upstream_straight_length / 2.0); + Tube upstream_IP_vacuum_padding(IP_acts_beampipe_OD / 2.0, IP_beampipe_ID / 2.0, + upstream_straight_length / 2.0); Tube upstream_IP_gold(IP_beampipe_ID / 2.0, IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, upstream_straight_length / 2.0); Tube upstream_IP_tube(IP_beampipe_ID / 2.0 + IP_beampipe_gold_thickness, IP_beampipe_OD / 2.0, upstream_straight_length / 2.0); - Volume v_downstream_IP_vacuum_fill("v_downstream_IP_vacuum_fill", downstream_IP_vacuum_fill, m_Vacuum); - Volume v_downstream_IP_acts_beampipe("v_downstream_IP_acts_beampipe", downstream_IP_acts_beampipe, m_Vacuum); - Volume v_downstream_IP_vacuum_padding("v_downstream_IP_vacuum_padding", downstream_IP_vacuum_padding, m_Vacuum); + Volume v_downstream_IP_vacuum_fill("v_downstream_IP_vacuum_fill", downstream_IP_vacuum_fill, + m_Vacuum); + Volume v_downstream_IP_acts_beampipe("v_downstream_IP_acts_beampipe", downstream_IP_acts_beampipe, + m_Vacuum); + Volume v_downstream_IP_vacuum_padding("v_downstream_IP_vacuum_padding", + downstream_IP_vacuum_padding, m_Vacuum); Volume v_downstream_IP_gold("v_downstream_IP_gold", downstream_IP_gold, m_Au); Volume v_downstream_IP_tube("v_downstream_IP_tube", downstream_IP_tube, m_Be); Volume v_upstream_IP_vacuum_fill("v_upstream_IP_vacuum_fill", upstream_IP_vacuum_fill, m_Vacuum); - Volume v_upstream_IP_acts_beampipe("v_upstream_IP_acts_beampipe", upstream_IP_acts_beampipe, m_Vacuum); - Volume v_upstream_IP_vacuum_padding("v_upstream_IP_vacuum_padding", upstream_IP_vacuum_padding, m_Vacuum); + Volume v_upstream_IP_acts_beampipe("v_upstream_IP_acts_beampipe", upstream_IP_acts_beampipe, + m_Vacuum); + Volume v_upstream_IP_vacuum_padding("v_upstream_IP_vacuum_padding", upstream_IP_vacuum_padding, + m_Vacuum); Volume v_upstream_IP_gold("v_upstream_IP_gold", upstream_IP_gold, m_Au); Volume v_upstream_IP_tube("v_upstream_IP_tube", upstream_IP_tube, m_Be); @@ -127,14 +135,17 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * assembly.placeVolume(v_upstream_IP_vacuum_fill, Position(0, 0, -upstream_straight_length / 2.0)); central_volume.placeVolume(v_upstream_IP_acts_beampipe, Position(0, 0, -upstream_straight_length / 2.0 - central_offset)); - assembly.placeVolume(v_upstream_IP_vacuum_padding, Position(0, 0, -upstream_straight_length / 2.0)); + assembly.placeVolume(v_upstream_IP_vacuum_padding, + Position(0, 0, -upstream_straight_length / 2.0)); assembly.placeVolume(v_upstream_IP_gold, Position(0, 0, -upstream_straight_length / 2.0)); assembly.placeVolume(v_upstream_IP_tube, Position(0, 0, -upstream_straight_length / 2.0)); - assembly.placeVolume(v_downstream_IP_vacuum_fill, Position(0, 0, downstream_straight_length / 2.0)); + assembly.placeVolume(v_downstream_IP_vacuum_fill, + Position(0, 0, downstream_straight_length / 2.0)); central_volume.placeVolume(v_downstream_IP_acts_beampipe, Position(0, 0, downstream_straight_length / 2.0 - central_offset)); - assembly.placeVolume(v_downstream_IP_vacuum_padding, Position(0, 0, downstream_straight_length / 2.0)); + assembly.placeVolume(v_downstream_IP_vacuum_padding, + Position(0, 0, downstream_straight_length / 2.0)); assembly.placeVolume(v_downstream_IP_gold, Position(0, 0, downstream_straight_length / 2.0)); assembly.placeVolume(v_downstream_IP_tube, Position(0, 0, downstream_straight_length / 2.0)); @@ -158,7 +169,7 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * std::vector rmax_vac, rmin_vac, z_vac; xml::Component vacuum = x_pipe.child(_Unicode(vacuum)); for (xml_coll_t x_zplane_i(vacuum, _Unicode(zplane)); x_zplane_i; ++x_zplane_i) { - xml_comp_t x_zplane = x_zplane_i; + xml_comp_t x_zplane = x_zplane_i; rmin_vac.push_back(0); rmax_vac.push_back(x_zplane.attr(_Unicode(OD)) / 2.0); z_vac.push_back(x_zplane.attr(_Unicode(z))); @@ -167,7 +178,7 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * std::vector rmax_coa, rmin_coa, z_coa; xml::Component coating = x_pipe.child(_Unicode(coating)); for (xml_coll_t x_zplane_i(coating, _Unicode(zplane)); x_zplane_i; ++x_zplane_i) { - xml_comp_t x_zplane = x_zplane_i; + xml_comp_t x_zplane = x_zplane_i; rmin_coa.push_back(0); rmax_coa.push_back(x_zplane.attr(_Unicode(OD)) / 2.0); z_coa.push_back(x_zplane.attr(_Unicode(z))); @@ -176,36 +187,31 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * std::vector rmax_mat, rmin_mat, z_mat; xml::Component matter = x_pipe.child(_Unicode(matter)); for (xml_coll_t x_zplane_i(matter, _Unicode(zplane)); x_zplane_i; ++x_zplane_i) { - xml_comp_t x_zplane = x_zplane_i; + xml_comp_t x_zplane = x_zplane_i; rmin_mat.push_back(0); rmax_mat.push_back(x_zplane.attr(_Unicode(OD)) / 2.0); z_mat.push_back(x_zplane.attr(_Unicode(z))); } - return std::tuple( - {0, 2.0 * M_PI, rmin_mat, rmax_mat, z_mat}, - {0, 2.0 * M_PI, rmin_coa, rmax_coa, z_coa}, - {0, 2.0 * M_PI, rmin_vac, rmax_vac, z_vac} - ); + return std::tuple({0, 2.0 * M_PI, rmin_mat, rmax_mat, z_mat}, + {0, 2.0 * M_PI, rmin_coa, rmax_coa, z_coa}, + {0, 2.0 * M_PI, rmin_vac, rmax_vac, z_vac}); }; //--------------------------------------------------------------------------------- // Create volumes - auto create_volumes = [&](const std::string& name, - xml::Component& x_pipe1, - xml::Component& x_pipe2, - xml_coll_t& x_additional_subtraction_i) { - + auto create_volumes = [&](const std::string& name, xml::Component& x_pipe1, + xml::Component& x_pipe2, xml_coll_t& x_additional_subtraction_i) { auto pipe1_polycones = zplane_to_polycones(x_pipe1); // lepton auto pipe2_polycones = zplane_to_polycones(x_pipe2); // hadron - auto crossing_angle = getAttrOrDefault(x_pipe2, _Unicode(crossing_angle), 0.0); - auto axis_intersection = getAttrOrDefault(x_pipe2, _Unicode(axis_intersection), 0.0); + auto crossing_angle = getAttrOrDefault(x_pipe2, _Unicode(crossing_angle), 0.0); + auto axis_intersection = getAttrOrDefault(x_pipe2, _Unicode(axis_intersection), 0.0); // transformation matrix: shift -> rotate -> shift auto tf = Transform3D(Position(0, 0, axis_intersection)); - if(name == "upstream") - tf *= Transform3D(RotationY( crossing_angle)); - else if(name == "downstream") - tf *= Transform3D(RotationY( std::abs(crossing_angle))); + if (name == "upstream") + tf *= Transform3D(RotationY(crossing_angle)); + else if (name == "downstream") + tf *= Transform3D(RotationY(std::abs(crossing_angle))); tf *= Transform3D(Position(0, 0, -axis_intersection)); // Global solids @@ -215,9 +221,9 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * if (name == "upstream" || name == "joint") { // union of all matter, coating, and vacuum - UnionSolid matter_union( std::get<0>(pipe1_polycones), std::get<0>(pipe2_polycones), tf); - UnionSolid coating_union( std::get<1>(pipe1_polycones), std::get<1>(pipe2_polycones), tf); - UnionSolid vacuum_union( std::get<2>(pipe1_polycones), std::get<2>(pipe2_polycones), tf); + UnionSolid matter_union(std::get<0>(pipe1_polycones), std::get<0>(pipe2_polycones), tf); + UnionSolid coating_union(std::get<1>(pipe1_polycones), std::get<1>(pipe2_polycones), tf); + UnionSolid vacuum_union(std::get<2>(pipe1_polycones), std::get<2>(pipe2_polycones), tf); // subtract vacuum from matter matter = SubtractionSolid(matter_union, vacuum_union); @@ -227,8 +233,7 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * vacuum = SubtractionSolid(vacuum_union, matter_union); // subtract coating from vacuum vacuum = SubtractionSolid(vacuum, coating_union); - } - else if (name == "downstream") { + } else if (name == "downstream") { auto thickness_pipe2 = getAttrOrDefault(x_pipe2, _Unicode(thickness), 0.0); auto thickness_coating = getAttrOrDefault(x_pipe2, _Unicode(coating), 0.0); auto horizontal_offset = getAttrOrDefault(x_pipe2, _Unicode(horizontal_offset), 0.0); @@ -243,125 +248,118 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * const double box_cut_sizex = 10000.0 * mm; const double box_cut_sizey = 10000.0 * mm; const double box_cut_sizez = 10000.0 * mm; - Box box_cut(box_cut_sizex,box_cut_sizey,box_cut_sizez); // virtual box to subtract from solids + Box box_cut(box_cut_sizex, box_cut_sizey, + box_cut_sizez); // virtual box to subtract from solids // horizontal offset of the outgoing h-beam pipe w.r.t. the incoming e-beam pipe (along X-axis) - const double p_mat_1[] = {horizontal_offset + thickness_coating + thickness_pipe2, 0, 0}; // matter + const double p_mat_1[] = {horizontal_offset + thickness_coating + thickness_pipe2, 0, + 0}; // matter const double p_coa_1[] = {horizontal_offset + thickness_coating, 0, 0}; // coating - const double p_vac_1[] = {horizontal_offset, 0, 0}; // vacuum - - pipe2_polycones_mat_cut = SubtractionSolid(std::get<0>(pipe2_polycones), box_cut, - tf * Transform3D(Position(p_mat_1[0] + box_cut_sizex, - p_mat_1[1], - p_mat_1[2]))); - pipe2_polycones_coa_cut = SubtractionSolid(std::get<1>(pipe2_polycones), box_cut, - tf * Transform3D(Position(p_coa_1[0] + box_cut_sizex, - p_coa_1[1], - p_coa_1[2]))); - pipe2_polycones_vac_cut = SubtractionSolid(std::get<2>(pipe2_polycones), box_cut, - tf * Transform3D(Position(p_vac_1[0] + box_cut_sizex, - p_vac_1[1], - p_vac_1[2]))); + const double p_vac_1[] = {horizontal_offset, 0, 0}; // vacuum + + pipe2_polycones_mat_cut = SubtractionSolid( + std::get<0>(pipe2_polycones), box_cut, + tf * Transform3D(Position(p_mat_1[0] + box_cut_sizex, p_mat_1[1], p_mat_1[2]))); + pipe2_polycones_coa_cut = SubtractionSolid( + std::get<1>(pipe2_polycones), box_cut, + tf * Transform3D(Position(p_coa_1[0] + box_cut_sizex, p_coa_1[1], p_coa_1[2]))); + pipe2_polycones_vac_cut = SubtractionSolid( + std::get<2>(pipe2_polycones), box_cut, + tf * Transform3D(Position(p_vac_1[0] + box_cut_sizex, p_vac_1[1], p_vac_1[2]))); // cut on the opposite side from the IP const double p_mat_2[] = {0, 0, cone_z_end + thickness_coating + thickness_pipe2}; // matter - const double p_coa_2[] = {0, 0, cone_z_end + thickness_coating}; // matter - const double p_vac_2[] = {0, 0, cone_z_end}; // vacuum - - pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut, box_cut, - tf * Transform3D(Position(p_mat_2[0], - p_mat_2[1], - p_mat_2[2] + box_cut_sizez))); - pipe2_polycones_coa_cut = SubtractionSolid(pipe2_polycones_coa_cut, box_cut, - tf * Transform3D(Position(p_coa_2[0], - p_coa_2[1], - p_coa_2[2] + box_cut_sizez))); - pipe2_polycones_vac_cut = SubtractionSolid(pipe2_polycones_vac_cut, box_cut, - tf * Transform3D(Position(p_vac_2[0], - p_vac_2[1], - p_vac_2[2] + box_cut_sizez))); + const double p_coa_2[] = {0, 0, cone_z_end + thickness_coating}; // matter + const double p_vac_2[] = {0, 0, cone_z_end}; // vacuum + + pipe2_polycones_mat_cut = SubtractionSolid( + pipe2_polycones_mat_cut, box_cut, + tf * Transform3D(Position(p_mat_2[0], p_mat_2[1], p_mat_2[2] + box_cut_sizez))); + pipe2_polycones_coa_cut = SubtractionSolid( + pipe2_polycones_coa_cut, box_cut, + tf * Transform3D(Position(p_coa_2[0], p_coa_2[1], p_coa_2[2] + box_cut_sizez))); + pipe2_polycones_vac_cut = SubtractionSolid( + pipe2_polycones_vac_cut, box_cut, + tf * Transform3D(Position(p_vac_2[0], p_vac_2[1], p_vac_2[2] + box_cut_sizez))); // cut on the IP side const double p_mat_3[] = {0, 0, cone_z_start - thickness_coating - thickness_pipe2}; // matter - const double p_coa_3[] = {0, 0, cone_z_start - thickness_coating}; // matter - const double p_vac_3[] = {0, 0, cone_z_start}; // vacuum - - pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut, box_cut, - tf * Transform3D(Position(p_mat_3[0], - p_mat_3[1], - p_mat_3[2] - box_cut_sizez))); - pipe2_polycones_coa_cut = SubtractionSolid(pipe2_polycones_coa_cut, box_cut, - tf * Transform3D(Position(p_coa_3[0], - p_coa_3[1], - p_coa_3[2] - box_cut_sizez))); - pipe2_polycones_vac_cut = SubtractionSolid(pipe2_polycones_vac_cut, box_cut, - tf * Transform3D(Position(p_vac_3[0], - p_vac_3[1], - p_vac_3[2] - box_cut_sizez))); + const double p_coa_3[] = {0, 0, cone_z_start - thickness_coating}; // matter + const double p_vac_3[] = {0, 0, cone_z_start}; // vacuum + + pipe2_polycones_mat_cut = SubtractionSolid( + pipe2_polycones_mat_cut, box_cut, + tf * Transform3D(Position(p_mat_3[0], p_mat_3[1], p_mat_3[2] - box_cut_sizez))); + pipe2_polycones_coa_cut = SubtractionSolid( + pipe2_polycones_coa_cut, box_cut, + tf * Transform3D(Position(p_coa_3[0], p_coa_3[1], p_coa_3[2] - box_cut_sizez))); + pipe2_polycones_vac_cut = SubtractionSolid( + pipe2_polycones_vac_cut, box_cut, + tf * Transform3D(Position(p_vac_3[0], p_vac_3[1], p_vac_3[2] - box_cut_sizez))); // add an extension to the h-beam pipe Tube tb_mat(0 * mm, extension_r + thickness_coating + extension_thickness, extension_z); Tube tb_coa(0 * mm, extension_r + thickness_coating, extension_z); Tube tb_vac(0 * mm, extension_r, extension_z); pipe2_polycones_mat_cut = UnionSolid(pipe2_polycones_mat_cut, tb_mat, - tf * Transform3D(RotationY(crossing_angle)) * - Transform3D(Position(0, 0, cone_z_end))); + tf * Transform3D(RotationY(crossing_angle)) * + Transform3D(Position(0, 0, cone_z_end))); pipe2_polycones_coa_cut = UnionSolid(pipe2_polycones_coa_cut, tb_coa, - tf * Transform3D(RotationY(crossing_angle)) * - Transform3D(Position(0, 0, cone_z_end))); + tf * Transform3D(RotationY(crossing_angle)) * + Transform3D(Position(0, 0, cone_z_end))); pipe2_polycones_vac_cut = UnionSolid(pipe2_polycones_vac_cut, tb_vac, - tf * Transform3D(RotationY(crossing_angle)) * - Transform3D(Position(0, 0, cone_z_end))); + tf * Transform3D(RotationY(crossing_angle)) * + Transform3D(Position(0, 0, cone_z_end))); // subtract h-vacuum from h-matter - pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut,pipe2_polycones_vac_cut); + pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut, pipe2_polycones_vac_cut); // subtract h-coating from h-matter - pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut,pipe2_polycones_coa_cut); + pipe2_polycones_mat_cut = SubtractionSolid(pipe2_polycones_mat_cut, pipe2_polycones_coa_cut); // subtract h-vacuum from h-coating - pipe2_polycones_coa_cut = SubtractionSolid(pipe2_polycones_coa_cut,pipe2_polycones_vac_cut); + pipe2_polycones_coa_cut = SubtractionSolid(pipe2_polycones_coa_cut, pipe2_polycones_vac_cut); // unite h-matter and e-matter - pipe2_polycones_mat_cut = UnionSolid(pipe2_polycones_mat_cut,std::get<0>(pipe1_polycones),tf); + pipe2_polycones_mat_cut = + UnionSolid(pipe2_polycones_mat_cut, std::get<0>(pipe1_polycones), tf); // unite h-coating and e-coating - pipe2_polycones_coa_cut = UnionSolid(pipe2_polycones_coa_cut,std::get<1>(pipe1_polycones),tf); + pipe2_polycones_coa_cut = + UnionSolid(pipe2_polycones_coa_cut, std::get<1>(pipe1_polycones), tf); // subtract e-vacuum from matter - matter = SubtractionSolid(pipe2_polycones_mat_cut,std::get<2>(pipe1_polycones),tf); + matter = SubtractionSolid(pipe2_polycones_mat_cut, std::get<2>(pipe1_polycones), tf); // subtract e-coating from matter - matter = SubtractionSolid(pipe2_polycones_mat_cut,std::get<1>(pipe1_polycones),tf); + matter = SubtractionSolid(pipe2_polycones_mat_cut, std::get<1>(pipe1_polycones), tf); // subtract e-vacuum from coating - coating = SubtractionSolid(pipe2_polycones_coa_cut,std::get<2>(pipe1_polycones),tf); + coating = SubtractionSolid(pipe2_polycones_coa_cut, std::get<2>(pipe1_polycones), tf); // unite h- and e-vacuum - vacuum = UnionSolid(pipe2_polycones_vac_cut,std::get<2>(pipe1_polycones),tf); + vacuum = UnionSolid(pipe2_polycones_vac_cut, std::get<2>(pipe1_polycones), tf); // subtract matter from vacuum - vacuum = SubtractionSolid(vacuum,matter); + vacuum = SubtractionSolid(vacuum, matter); // subtract coating from vacuum - vacuum = SubtractionSolid(vacuum,coating); - } - else { + vacuum = SubtractionSolid(vacuum, coating); + } else { throw std::runtime_error("Wrong side name (should be: upstream or downstream): " + name); } // subtract additional vacuum from matter and add it to the vacuum for (; x_additional_subtraction_i; ++x_additional_subtraction_i) { - xml_comp_t x_additional_subtraction = x_additional_subtraction_i; + xml_comp_t x_additional_subtraction = x_additional_subtraction_i; auto additional_polycones = zplane_to_polycones(x_additional_subtraction); - auto additional_crossing_angle = getAttrOrDefault(x_additional_subtraction, _Unicode(crossing_angle), 0.0); + auto additional_crossing_angle = + getAttrOrDefault(x_additional_subtraction, _Unicode(crossing_angle), 0.0); auto additional_tf = Transform3D(RotationY(additional_crossing_angle)); // final matter, coating, and vacuum - matter = SubtractionSolid(matter, std::get<2>(additional_polycones),tf * additional_tf); - coating = SubtractionSolid(coating,std::get<2>(additional_polycones),tf * additional_tf); - vacuum = UnionSolid(vacuum,std::get<2>(additional_polycones),tf * additional_tf); + matter = SubtractionSolid(matter, std::get<2>(additional_polycones), tf * additional_tf); + coating = SubtractionSolid(coating, std::get<2>(additional_polycones), tf * additional_tf); + vacuum = UnionSolid(vacuum, std::get<2>(additional_polycones), tf * additional_tf); } - return std::tuple( - {"v_" + name + "_matter", matter, m_SS}, - {"v_" + name + "_coating", coating, m_Cu}, - {"v_" + name + "_vacuum", vacuum, m_Vacuum} - ); + return std::tuple({"v_" + name + "_matter", matter, m_SS}, + {"v_" + name + "_coating", coating, m_Cu}, + {"v_" + name + "_vacuum", vacuum, m_Vacuum}); }; //--------------------------------------------------------------------------------- @@ -373,14 +371,15 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * xml::Component upstream_c = x_det.child(_Unicode(upstream)); xml::Component incoming_hadron_c = upstream_c.child(_Unicode(incoming_hadron)); xml::Component outgoing_lepton_c = upstream_c.child(_Unicode(outgoing_lepton)); - xml::Component upstream_lepton_hadron_c = upstream_c.child(_Unicode(upstream_lepton_hadron_joint)); - xml_coll_t additional_subtractions_upstream(upstream_c, _Unicode(additional_subtraction)); + xml::Component upstream_lepton_hadron_c = + upstream_c.child(_Unicode(upstream_lepton_hadron_joint)); + xml_coll_t additional_subtractions_upstream(upstream_c, _Unicode(additional_subtraction)); - auto volumes_upstream = create_volumes( - "upstream", outgoing_lepton_c, incoming_hadron_c,additional_subtractions_upstream); + auto volumes_upstream = create_volumes("upstream", outgoing_lepton_c, incoming_hadron_c, + additional_subtractions_upstream); - auto joint_upstream = create_volumes( - "joint", upstream_lepton_hadron_c, upstream_lepton_hadron_c,additional_subtractions_upstream); + auto joint_upstream = create_volumes("joint", upstream_lepton_hadron_c, upstream_lepton_hadron_c, + additional_subtractions_upstream); // reflect auto tf_upstream = Transform3D(RotationZYX(0, 0, 0)); @@ -408,19 +407,21 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * xml::Component downstream_c = x_det.child(_Unicode(downstream)); xml::Component incoming_lepton_c = downstream_c.child(_Unicode(incoming_lepton)); xml::Component outgoing_hadron_c = downstream_c.child(_Unicode(outgoing_hadron)); - xml_coll_t additional_subtractions_downstream(downstream_c, _Unicode(additional_subtraction)); + xml_coll_t additional_subtractions_downstream(downstream_c, _Unicode(additional_subtraction)); - auto volumes_downstream = create_volumes( - "downstream", incoming_lepton_c, outgoing_hadron_c,additional_subtractions_downstream); + auto volumes_downstream = create_volumes("downstream", incoming_lepton_c, outgoing_hadron_c, + additional_subtractions_downstream); // transform auto tf_downstream = - Transform3D(Position(0, 0, +getAttrOrDefault(outgoing_hadron_c, _Unicode(axis_intersection), 0.0))) * - Transform3D(RotationY(getAttrOrDefault(outgoing_hadron_c, _Unicode(crossing_angle), 0.0))) * - Transform3D(Position(0, 0, -getAttrOrDefault(outgoing_hadron_c, _Unicode(axis_intersection), 0.0))); + Transform3D( + Position(0, 0, +getAttrOrDefault(outgoing_hadron_c, _Unicode(axis_intersection), 0.0))) * + Transform3D(RotationY(getAttrOrDefault(outgoing_hadron_c, _Unicode(crossing_angle), 0.0))) * + Transform3D( + Position(0, 0, -getAttrOrDefault(outgoing_hadron_c, _Unicode(axis_intersection), 0.0))); // reflect - if(getAttrOrDefault(downstream_c, _Unicode(reflect), true)) + if (getAttrOrDefault(downstream_c, _Unicode(reflect), true)) tf_downstream = Transform3D(RotationZYX(0, M_PI, 0)); // add matter @@ -428,7 +429,7 @@ static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens * // add coating assembly.placeVolume(std::get<1>(volumes_downstream), tf_downstream); // add vacuum - if(getAttrOrDefault(downstream_c, _Unicode(place_vacuum), true)) + if (getAttrOrDefault(downstream_c, _Unicode(place_vacuum), true)) assembly.placeVolume(std::get<2>(volumes_downstream), tf_downstream); // ----------------------------- diff --git a/src/LumiWindow_geo.cpp b/src/LumiWindow_geo.cpp index 5f9446729..347a5bd06 100644 --- a/src/LumiWindow_geo.cpp +++ b/src/LumiWindow_geo.cpp @@ -12,38 +12,37 @@ using namespace std; using namespace dd4hep; -static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) -{ +static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) { sens.setType("calorimeter"); - xml_det_t x_det = e; - xml_comp_t x_dim = x_det.dimensions(); - xml_comp_t x_pos = x_det.position(); - xml_comp_t x_rot = x_det.rotation(); + xml_det_t x_det = e; + xml_comp_t x_dim = x_det.dimensions(); + xml_comp_t x_pos = x_det.position(); + xml_comp_t x_rot = x_det.rotation(); // - string det_name = x_det.nameStr(); - string mat_name = dd4hep::getAttrOrDefault( x_det, _U(material), "Aluminum" ); + string det_name = x_det.nameStr(); + string mat_name = dd4hep::getAttrOrDefault(x_det, _U(material), "Aluminum"); // - double sizeX = x_dim.x(); - double sizeY = x_dim.y(); - double sizeZ = x_dim.z(); - double posX = x_pos.x(); - double posY = x_pos.y(); - double posZ = x_pos.z(); - double rotX = x_rot.x(); - double rotY = x_rot.y(); - double rotZ = x_rot.z(); - - Box box( sizeX, sizeY, sizeZ ); - Volume vol( det_name + "_vol_ExitWindow", box, description.material( mat_name ) ); - vol.setVisAttributes( description.visAttributes(x_det.visStr()) ); + double sizeX = x_dim.x(); + double sizeY = x_dim.y(); + double sizeZ = x_dim.z(); + double posX = x_pos.x(); + double posY = x_pos.y(); + double posZ = x_pos.z(); + double rotX = x_rot.x(); + double rotY = x_rot.y(); + double rotZ = x_rot.z(); + + Box box(sizeX, sizeY, sizeZ); + Volume vol(det_name + "_vol_ExitWindow", box, description.material(mat_name)); + vol.setVisAttributes(description.visAttributes(x_det.visStr())); vol.setSensitiveDetector(sens); - Transform3D pos( RotationZYX(rotX, rotY, rotZ), Position(posX, posY, posZ) ); + Transform3D pos(RotationZYX(rotX, rotY, rotZ), Position(posX, posY, posZ)); DetElement det(det_name, x_det.id()); - Volume motherVol = description.pickMotherVolume( det ); - PlacedVolume phv = motherVol.placeVolume( vol, pos ); + Volume motherVol = description.pickMotherVolume(det); + PlacedVolume phv = motherVol.placeVolume(vol, pos); phv.addPhysVolID("system", x_det.id()); det.setPlacement(phv); diff --git a/src/magnetVacuumFF.cpp b/src/magnetVacuumFF.cpp index 76ec7d614..3a7f31a60 100644 --- a/src/magnetVacuumFF.cpp +++ b/src/magnetVacuumFF.cpp @@ -1,7 +1,6 @@ // SPDX-License-Identifier: LGPL-3.0-or-later // Copyright (C) 2023 Alex Jentsch - #include "DD4hep/DetFactoryHelper.h" #include "DD4hep/Printout.h" #include @@ -27,252 +26,260 @@ using namespace dd4hep; static double getRotatedZ(double z, double x, double angle); static double getRotatedX(double z, double x, double angle); -static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) -{ - - xml_det_t x_det = e; - string det_name = x_det.nameStr(); - DetElement sdet(det_name, x_det.id()); - Assembly assembly(det_name + "_assembly"); - Material m_Vac = det.material("Vacuum"); - string vis_name = x_det.visStr(); - double endOfCentralBeamPipe_z = dd4hep::getAttrOrDefault(x_det, _Unicode(endOfCentralBeamPipe_z), 0); - - PlacedVolume pv_assembly; - - //---------------------------------------------- - // Starting point is only the magnet centers, - // lengths, rotations, and radii --> - // everything else calculated internally to - // make it easier to update later. - //---------------------------------------------- - - bool makeIP_B0pfVacuum = true; //This is for the special gap location between IP and b0pf - - //information for actual FF magnets, with magnet centers as reference - vector radii_magnet; - vector lengths_magnet; - vector rotation_magnet; - vector x_elem_magnet; - vector y_elem_magnet; - vector z_elem_magnet; - - //calculated entrance/exit points of FF magnet - vector x_beg; - vector z_beg; - vector x_end; - vector z_end; - - //calculated center of gap regions between magnets, rotation, and length - vector angle_elem_gap; - vector z_gap; - vector x_gap; - vector length_gap; - - //storage elements for CutTube geometry element used for gaps - vector inRadius; - vector outRadius; - vector nxLow; - vector nyLow; - vector nzLow; - vector nxHigh; - vector nyHigh; - vector nzHigh; - vector phi_initial; - vector phi_final; - - for(xml_coll_t c(x_det,_U(element)); c; ++c){ - - xml_dim_t pos = c.child(_U(placement)); - double pos_x = pos.x(); - double pos_y = pos.y(); - double pos_z = pos.z(); - double pos_theta = pos.attr(_U(theta)); - xml_dim_t dims = c.child(_U(dimensions)); //dimensions(); - double dim_z = dims.z(); - xml_dim_t apperture = c.child(_Unicode(apperture)); - double app_r = apperture.r(); - - radii_magnet.push_back(app_r); // cm - lengths_magnet.push_back(dim_z); //cm - rotation_magnet.push_back(pos_theta); // radians - x_elem_magnet.push_back(pos_x*dd4hep::cm); - y_elem_magnet.push_back(pos_y*dd4hep::cm); - z_elem_magnet.push_back(pos_z*dd4hep::cm); - - } - - int numMagnets = radii_magnet.size(); //number of actual FF magnets between IP and FF detectors - int numGaps = numMagnets - 1; //number of gaps between magnets (excluding the IP to B0pf transition -- special case) - - //------------------------------------------- - // override numbers for the first element --> - // doesn't use the actual B0pf geometry!!! - // -->it's based on the B0 beam pipe - // this needs to be fixed later to read-in - // that beam pipe geometry - //------------------------------------------- - - radii_magnet[0] = 2.9; // cm - lengths_magnet[0] = 120.0; // cm - rotation_magnet[0] = -0.025; // radians - x_elem_magnet[0] = -16.5; // cm - y_elem_magnet[0] = 0.0; // cm - z_elem_magnet[0] = 640.0; // cm - - //------------------------------------------- - //calculate entrance/exit points of magnets - //------------------------------------------- - - for(int i = 0; i < numMagnets; i++){ - - // need to use the common coordinate system --> - // use x = z, and y = x to make things easier - - z_beg.push_back(getRotatedZ(-0.5*lengths_magnet[i], 0.0, rotation_magnet[i]) + z_elem_magnet[i]); - z_end.push_back(getRotatedZ( 0.5*lengths_magnet[i], 0.0, rotation_magnet[i]) + z_elem_magnet[i]); - x_beg.push_back(getRotatedX(-0.5*lengths_magnet[i], 0.0, rotation_magnet[i]) + x_elem_magnet[i]); - x_end.push_back(getRotatedX( 0.5*lengths_magnet[i], 0.0, rotation_magnet[i]) + x_elem_magnet[i]); - - } - - //------------------------------------------ - // this part is a bit ugly for now - - // it's to make the vacuum volume between the - // end of the IP beam pipe and the beginning of - // beginning of the B0pf magnet - // - // -->the volume will be calculated at the end - //------------------------------------------- - - double diameterReduce = 11.0*dd4hep::cm; //size reduction to avoid overlap with electron pipe - double vacuumDiameterEntrance = 25.792*dd4hep::cm - diameterReduce; //extracted from central_beampipe.xml, line 64 - double vacuumDiameterExit = 17.4*dd4hep::cm; //15mrad @ entrance to magnet to not overlap electron magnet - double crossingAngle = -0.025; //radians - double endOfCentralBeamPipe_x = endOfCentralBeamPipe_z*crossingAngle; - - //----------------------------------------------- - //calculate gap region center, length, and angle - //----------------------------------------------- - - for(int i = 1; i < numMagnets; i++){ - - angle_elem_gap.push_back((x_beg[i] - x_end[i-1])/(z_beg[i] - z_end[i-1])); - length_gap.push_back(sqrt(pow(z_beg[i] - z_end[i-1], 2) + pow(x_beg[i] - x_end[i-1], 2))); - z_gap.push_back(z_end[i-1] + 0.5*length_gap[i-1]*cos(angle_elem_gap[i-1])); - x_gap.push_back(x_end[i-1] + 0.5*length_gap[i-1]*sin(angle_elem_gap[i-1])); - - } - - //----------------------------------------------- - // fill CutTube storage elements - //----------------------------------------------- - - for(int gapIdx = 0; gapIdx < numGaps; gapIdx++){ - - inRadius.push_back(0.0); - outRadius.push_back(radii_magnet[gapIdx+1]); - phi_initial.push_back(0.0); - phi_final.push_back(2*M_PI); - nxLow.push_back(-(length_gap[gapIdx]/2.0)*sin(rotation_magnet[gapIdx]-angle_elem_gap[gapIdx])); - nyLow.push_back(0.0); - nzLow.push_back(-(length_gap[gapIdx]/2.0)*cos(rotation_magnet[gapIdx]-angle_elem_gap[gapIdx])); - nxHigh.push_back((length_gap[gapIdx]/2.0)*sin(rotation_magnet[gapIdx+1]-angle_elem_gap[gapIdx])); - nyHigh.push_back(0.0); - nzHigh.push_back((length_gap[gapIdx]/2.0)*cos(rotation_magnet[gapIdx+1]-angle_elem_gap[gapIdx])); - - } - - //----------------------- - // inside magnets - //----------------------- - - for(int pieceIdx = 0; pieceIdx < numMagnets; pieceIdx++){ - - std::string piece_name = Form("MagnetVacuum%d", pieceIdx); - - Tube magnetPiece(piece_name, 0.0, radii_magnet[pieceIdx], lengths_magnet[pieceIdx]/2); - Volume vpiece(piece_name, magnetPiece, m_Vac); - sdet.setAttributes(det, vpiece, x_det.regionStr(), x_det.limitsStr(), vis_name); - - auto pv = assembly.placeVolume(vpiece, Transform3D(RotationY(rotation_magnet[pieceIdx]), - Position(x_elem_magnet[pieceIdx], y_elem_magnet[pieceIdx], z_elem_magnet[pieceIdx]))); - pv.addPhysVolID("sector", 1); - - DetElement de(sdet, Form("sector%d_de", pieceIdx), 1); - de.setPlacement(pv); - - } - - //-------------------------- - //between magnets - //-------------------------- +static Ref_t create_detector(Detector& det, xml_h e, SensitiveDetector /* sens */) { + + xml_det_t x_det = e; + string det_name = x_det.nameStr(); + DetElement sdet(det_name, x_det.id()); + Assembly assembly(det_name + "_assembly"); + Material m_Vac = det.material("Vacuum"); + string vis_name = x_det.visStr(); + double endOfCentralBeamPipe_z = + dd4hep::getAttrOrDefault(x_det, _Unicode(endOfCentralBeamPipe_z), 0); + + PlacedVolume pv_assembly; + + //---------------------------------------------- + // Starting point is only the magnet centers, + // lengths, rotations, and radii --> + // everything else calculated internally to + // make it easier to update later. + //---------------------------------------------- + + bool makeIP_B0pfVacuum = true; //This is for the special gap location between IP and b0pf + + //information for actual FF magnets, with magnet centers as reference + vector radii_magnet; + vector lengths_magnet; + vector rotation_magnet; + vector x_elem_magnet; + vector y_elem_magnet; + vector z_elem_magnet; + + //calculated entrance/exit points of FF magnet + vector x_beg; + vector z_beg; + vector x_end; + vector z_end; + + //calculated center of gap regions between magnets, rotation, and length + vector angle_elem_gap; + vector z_gap; + vector x_gap; + vector length_gap; + + //storage elements for CutTube geometry element used for gaps + vector inRadius; + vector outRadius; + vector nxLow; + vector nyLow; + vector nzLow; + vector nxHigh; + vector nyHigh; + vector nzHigh; + vector phi_initial; + vector phi_final; + + for (xml_coll_t c(x_det, _U(element)); c; ++c) { + + xml_dim_t pos = c.child(_U(placement)); + double pos_x = pos.x(); + double pos_y = pos.y(); + double pos_z = pos.z(); + double pos_theta = pos.attr(_U(theta)); + xml_dim_t dims = c.child(_U(dimensions)); //dimensions(); + double dim_z = dims.z(); + xml_dim_t apperture = c.child(_Unicode(apperture)); + double app_r = apperture.r(); + + radii_magnet.push_back(app_r); // cm + lengths_magnet.push_back(dim_z); //cm + rotation_magnet.push_back(pos_theta); // radians + x_elem_magnet.push_back(pos_x * dd4hep::cm); + y_elem_magnet.push_back(pos_y * dd4hep::cm); + z_elem_magnet.push_back(pos_z * dd4hep::cm); + } + + int numMagnets = radii_magnet.size(); //number of actual FF magnets between IP and FF detectors + int numGaps = + numMagnets - + 1; //number of gaps between magnets (excluding the IP to B0pf transition -- special case) + + //------------------------------------------- + // override numbers for the first element --> + // doesn't use the actual B0pf geometry!!! + // -->it's based on the B0 beam pipe + // this needs to be fixed later to read-in + // that beam pipe geometry + //------------------------------------------- + + radii_magnet[0] = 2.9; // cm + lengths_magnet[0] = 120.0; // cm + rotation_magnet[0] = -0.025; // radians + x_elem_magnet[0] = -16.5; // cm + y_elem_magnet[0] = 0.0; // cm + z_elem_magnet[0] = 640.0; // cm + + //------------------------------------------- + //calculate entrance/exit points of magnets + //------------------------------------------- + + for (int i = 0; i < numMagnets; i++) { + + // need to use the common coordinate system --> + // use x = z, and y = x to make things easier + + z_beg.push_back(getRotatedZ(-0.5 * lengths_magnet[i], 0.0, rotation_magnet[i]) + + z_elem_magnet[i]); + z_end.push_back(getRotatedZ(0.5 * lengths_magnet[i], 0.0, rotation_magnet[i]) + + z_elem_magnet[i]); + x_beg.push_back(getRotatedX(-0.5 * lengths_magnet[i], 0.0, rotation_magnet[i]) + + x_elem_magnet[i]); + x_end.push_back(getRotatedX(0.5 * lengths_magnet[i], 0.0, rotation_magnet[i]) + + x_elem_magnet[i]); + } + + //------------------------------------------ + // this part is a bit ugly for now - + // it's to make the vacuum volume between the + // end of the IP beam pipe and the beginning of + // beginning of the B0pf magnet + // + // -->the volume will be calculated at the end + //------------------------------------------- + + double diameterReduce = 11.0 * dd4hep::cm; //size reduction to avoid overlap with electron pipe + double vacuumDiameterEntrance = + 25.792 * dd4hep::cm - diameterReduce; //extracted from central_beampipe.xml, line 64 + double vacuumDiameterExit = + 17.4 * dd4hep::cm; //15mrad @ entrance to magnet to not overlap electron magnet + double crossingAngle = -0.025; //radians + double endOfCentralBeamPipe_x = endOfCentralBeamPipe_z * crossingAngle; + + //----------------------------------------------- + //calculate gap region center, length, and angle + //----------------------------------------------- + + for (int i = 1; i < numMagnets; i++) { + + angle_elem_gap.push_back((x_beg[i] - x_end[i - 1]) / (z_beg[i] - z_end[i - 1])); + length_gap.push_back(sqrt(pow(z_beg[i] - z_end[i - 1], 2) + pow(x_beg[i] - x_end[i - 1], 2))); + z_gap.push_back(z_end[i - 1] + 0.5 * length_gap[i - 1] * cos(angle_elem_gap[i - 1])); + x_gap.push_back(x_end[i - 1] + 0.5 * length_gap[i - 1] * sin(angle_elem_gap[i - 1])); + } + + //----------------------------------------------- + // fill CutTube storage elements + //----------------------------------------------- + + for (int gapIdx = 0; gapIdx < numGaps; gapIdx++) { + + inRadius.push_back(0.0); + outRadius.push_back(radii_magnet[gapIdx + 1]); + phi_initial.push_back(0.0); + phi_final.push_back(2 * M_PI); + nxLow.push_back(-(length_gap[gapIdx] / 2.0) * + sin(rotation_magnet[gapIdx] - angle_elem_gap[gapIdx])); + nyLow.push_back(0.0); + nzLow.push_back(-(length_gap[gapIdx] / 2.0) * + cos(rotation_magnet[gapIdx] - angle_elem_gap[gapIdx])); + nxHigh.push_back((length_gap[gapIdx] / 2.0) * + sin(rotation_magnet[gapIdx + 1] - angle_elem_gap[gapIdx])); + nyHigh.push_back(0.0); + nzHigh.push_back((length_gap[gapIdx] / 2.0) * + cos(rotation_magnet[gapIdx + 1] - angle_elem_gap[gapIdx])); + } + + //----------------------- + // inside magnets + //----------------------- + + for (int pieceIdx = 0; pieceIdx < numMagnets; pieceIdx++) { + + std::string piece_name = Form("MagnetVacuum%d", pieceIdx); + + Tube magnetPiece(piece_name, 0.0, radii_magnet[pieceIdx], lengths_magnet[pieceIdx] / 2); + Volume vpiece(piece_name, magnetPiece, m_Vac); + sdet.setAttributes(det, vpiece, x_det.regionStr(), x_det.limitsStr(), vis_name); + + auto pv = assembly.placeVolume( + vpiece, Transform3D(RotationY(rotation_magnet[pieceIdx]), + Position(x_elem_magnet[pieceIdx], y_elem_magnet[pieceIdx], + z_elem_magnet[pieceIdx]))); + pv.addPhysVolID("sector", 1); + + DetElement de(sdet, Form("sector%d_de", pieceIdx), 1); + de.setPlacement(pv); + } + + //-------------------------- + //between magnets + //-------------------------- + + for (int pieceIdx = numMagnets; pieceIdx < numGaps + numMagnets; pieceIdx++) { + + int correctIdx = pieceIdx - numMagnets; + + std::string piece_name = Form("GapVacuum%d", correctIdx); + + CutTube gapPiece(piece_name, inRadius[correctIdx], outRadius[correctIdx], + length_gap[correctIdx] / 2, phi_initial[correctIdx], phi_final[correctIdx], + nxLow[correctIdx], nyLow[correctIdx], nzLow[correctIdx], nxHigh[correctIdx], + nyHigh[correctIdx], nzHigh[correctIdx]); + + Volume vpiece(piece_name, gapPiece, m_Vac); + sdet.setAttributes(det, vpiece, x_det.regionStr(), x_det.limitsStr(), vis_name); + + auto pv = assembly.placeVolume( + vpiece, Transform3D(RotationY(angle_elem_gap[correctIdx]), + Position(x_gap[correctIdx], 0.0, z_gap[correctIdx]))); + pv.addPhysVolID("sector", 1); + + DetElement de(sdet, Form("sector%d_de", pieceIdx), 1); + de.setPlacement(pv); + } + + //-------------------------------------------------------------- + //make and place vacuum volume to connect IP beam pipe to B0pf + //-------------------------------------------------------------- + + if (makeIP_B0pfVacuum) { + + double specialGapLength = sqrt(pow(z_beg[0] - endOfCentralBeamPipe_z, 2) + + pow(x_beg[0] - endOfCentralBeamPipe_x, 2)) - + 0.1; + double specialGap_z = 0.5 * specialGapLength * cos(crossingAngle) + endOfCentralBeamPipe_z; + double specialGap_x = 0.5 * specialGapLength * sin(crossingAngle) + endOfCentralBeamPipe_x; - for(int pieceIdx = numMagnets; pieceIdx < numGaps + numMagnets; pieceIdx++){ + std::string piece_name = Form("GapVacuum%d", numGaps + numMagnets); - int correctIdx = pieceIdx-numMagnets; + Cone specialGap(piece_name, specialGapLength / 2, 0.0, vacuumDiameterEntrance / 2, 0.0, + vacuumDiameterExit / 2); - std::string piece_name = Form("GapVacuum%d", correctIdx); + Volume specialGap_v(piece_name, specialGap, m_Vac); + sdet.setAttributes(det, specialGap_v, x_det.regionStr(), x_det.limitsStr(), vis_name); - CutTube gapPiece(piece_name, inRadius[correctIdx], outRadius[correctIdx], length_gap[correctIdx]/2, phi_initial[correctIdx], phi_final[correctIdx], - nxLow[correctIdx], nyLow[correctIdx], nzLow[correctIdx], nxHigh[correctIdx], nyHigh[correctIdx], nzHigh[correctIdx]); + auto pv = + assembly.placeVolume(specialGap_v, Transform3D(RotationY(crossingAngle), + Position(specialGap_x, 0.0, specialGap_z))); + pv.addPhysVolID("sector", 1); - Volume vpiece(piece_name, gapPiece, m_Vac); - sdet.setAttributes(det, vpiece, x_det.regionStr(), x_det.limitsStr(), vis_name); + DetElement de(sdet, Form("sector%d_de", numGaps + numMagnets), 1); + de.setPlacement(pv); + } - auto pv = assembly.placeVolume(vpiece, Transform3D(RotationY(angle_elem_gap[correctIdx]), - Position(x_gap[correctIdx], 0.0, z_gap[correctIdx]))); - pv.addPhysVolID("sector", 1); + //---------------------------------------------------- - DetElement de(sdet, Form("sector%d_de", pieceIdx), 1); - de.setPlacement(pv); + pv_assembly = det.pickMotherVolume(sdet).placeVolume(assembly); + pv_assembly.addPhysVolID("system", x_det.id()).addPhysVolID("barrel", 1); + sdet.setPlacement(pv_assembly); + assembly->GetShape()->ComputeBBox(); - - } - - //-------------------------------------------------------------- - //make and place vacuum volume to connect IP beam pipe to B0pf - //-------------------------------------------------------------- - - if(makeIP_B0pfVacuum){ - - double specialGapLength = sqrt(pow(z_beg[0] - endOfCentralBeamPipe_z, 2) + pow(x_beg[0] - endOfCentralBeamPipe_x, 2)) - 0.1; - double specialGap_z = 0.5*specialGapLength*cos(crossingAngle) + endOfCentralBeamPipe_z; - double specialGap_x = 0.5*specialGapLength*sin(crossingAngle) + endOfCentralBeamPipe_x; - - std::string piece_name = Form("GapVacuum%d", numGaps + numMagnets); - - Cone specialGap(piece_name, specialGapLength/2, 0.0, vacuumDiameterEntrance/2, 0.0, vacuumDiameterExit/2 ); - - Volume specialGap_v(piece_name, specialGap, m_Vac); - sdet.setAttributes(det, specialGap_v, x_det.regionStr(), x_det.limitsStr(), vis_name); - - auto pv = assembly.placeVolume(specialGap_v, Transform3D(RotationY(crossingAngle), Position(specialGap_x, 0.0, specialGap_z))); - pv.addPhysVolID("sector", 1); - - DetElement de(sdet, Form("sector%d_de", numGaps + numMagnets), 1); - de.setPlacement(pv); - - } - - //---------------------------------------------------- - - pv_assembly = det.pickMotherVolume(sdet).placeVolume(assembly); - pv_assembly.addPhysVolID("system", x_det.id()).addPhysVolID("barrel", 1); - sdet.setPlacement(pv_assembly); - assembly->GetShape()->ComputeBBox(); - - return sdet; + return sdet; } -double getRotatedZ(double z, double x, double angle){ +double getRotatedZ(double z, double x, double angle) { return z * cos(angle) - x * sin(angle); } - return z*cos(angle) - x*sin(angle); -} - -double getRotatedX(double z, double x, double angle){ - - return z*sin(angle) + x*cos(angle); -} +double getRotatedX(double z, double x, double angle) { return z * sin(angle) + x * cos(angle); } DECLARE_DETELEMENT(magnetElementInnerVacuum, create_detector)