Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • EIC/detectors/athena
  • zwzhao/athena
  • FernandoTA/athena
  • palspeic/athena
4 results
Show changes
Hide whitespace changes
Inline Side-by-side
Showing
with 814 additions and 289 deletions
compact/mrich_alt.xml → compact/unused/mrich_alt.xml
View file @ 1008a184
......@@ -2,9 +2,9 @@
<comment> MRICH (alternative design) </comment>
<define>
<constant name="MRICH_rmin" value="10*cm"/>
<constant name="MRICH_rmax" value="BackwardPID_rmax"/>
<constant name="MRICH_length" value="BackwardPID_length"/>
<constant name="MRICH_zmin" value="BackwardPID_zmin"/>
<constant name="MRICH_rmax" value="BackwardPIDRegion_rmax"/>
<constant name="MRICH_length" value="BackwardPIDRegion_length"/>
<constant name="MRICH_zmin" value="BackwardPIDRegion_zmin"/>
<constant name="MRICHAerogel_thickness" value="30.0*mm"/>
<constant name="MRICHAerogel_width" value="126.5*mm"/>
......
compact/vertex_tracker.xml
View file @ 1008a184
......@@ -5,8 +5,8 @@
THis value probably can live in the file that includes this one.
</comment>
<constant name="ITS3Thickness" value="40*um"/>
<constant name="TrackerCarbon_thickness" value="0.12*mm"/>
<constant name="TrackerEndcapAluminum_thickness" value="0.15*mm"/>
<constant name="VertexCarbon_thickness" value="0.12*mm"/>
<constant name="VertexEndcapAluminumThickness" value="0.15*mm"/>
<documentation>
#### Vertex Tracker Barrel Parameters
......@@ -27,7 +27,14 @@ Simple carbon fiber support shell.
</documentation>
<constant name="VertexBarrel_length" value="VertexTrackerBarrel_length"/>
<constant name="VertexBarrel_length" value="300.0*mm"/>
<constant name="VertexTrackerEndcapP_rmin" value="VertexTrackingRegion_rmin"/>
<constant name="VertexTrackerEndcapN_rmin" value="VertexTrackingRegion_rmin"/>
<constant name="VertexTrackerEndcapP_rmax" value="VertexTrackingRegion_rmax"/>
<constant name="VertexTrackerEndcapN_rmax" value="VertexTrackingRegion_rmax"/>
<constant name="VertexTrackerEndcapP_zmax" value="VertexTrackingRegionP_zmax"/>
<constant name="VertexTrackerEndcapN_zmax" value="VertexTrackingRegionN_zmax"/>
<constant name="VertexBarrelLayer_length" value="VertexBarrel_length - 1*mm"/>
<constant name="VertexBarrelMod_length" value="VertexBarrel_length - 2*mm"/>
......@@ -98,8 +105,11 @@ Simple carbon fiber support shell.
</display>
<detectors>
<documentation level="5">
### Actual detectors
</documentation>
<detector
id="VertexBarrel_ID"
id="VertexBarrel_0_ID"
name="VertexBarrel"
type="athena_VertexBarrel"
readout="VertexBarrelHits"
......@@ -177,7 +187,7 @@ Simple carbon fiber support shell.
</detector>
<detector
id="VertexEndcapP_ID"
id="VertexEndcapP_0_ID"
name="VertexEndcapP"
type="athena_TrapEndcapTracker"
readout="VertexEndcapHits"
......@@ -201,8 +211,8 @@ Simple carbon fiber support shell.
<module name="Module1" vis="AnlProcess_Blue">
<trd x1="VertexEndcapMod1_x1/2.0" x2="VertexEndcapMod1_x2/2.0" z="VertexEndcapMod1_y/2"/>
<module_component thickness="ITS3Thickness" material="Silicon" sensitive="true"/>
<module_component thickness="TrackerEndcapAluminum_thickness" material="Aluminum"/>
<module_component thickness="TrackerCarbon_thickness" material="CarbonFiber"/>
<module_component thickness="VertexEndcapAluminumThickness" material="Aluminum"/>
<module_component thickness="VertexCarbon_thickness" material="CarbonFiber"/>
</module>
<layer id="1">
<envelope vis="TrackerVis"
......@@ -220,7 +230,7 @@ Simple carbon fiber support shell.
<detector
id="VertexEndcapN_ID"
id="VertexEndcapN_0_ID"
name="VertexEndcapN"
type="athena_TrapEndcapTracker"
readout="VertexEndcapHits"
......@@ -243,8 +253,8 @@ Simple carbon fiber support shell.
<module name="Module1" vis="AnlProcess_Blue">
<trd x1="VertexEndcapMod1_x1/2.0" x2="VertexEndcapMod1_x2/2.0" z="VertexEndcapMod1_y/2"/>
<module_component thickness="ITS3Thickness" material="Silicon" sensitive="true"/>
<module_component thickness="TrackerEndcapAluminum_thickness" material="Aluminum"/>
<module_component thickness="TrackerCarbon_thickness" material="CarbonFiber"/>
<module_component thickness="VertexEndcapAluminumThickness" material="Aluminum"/>
<module_component thickness="VertexCarbon_thickness" material="CarbonFiber"/>
</module>
<layer id="1">
<envelope vis="TrackerVis"
......@@ -260,58 +270,6 @@ Simple carbon fiber support shell.
</layer>
</detector>
<!--
<detector id="VertexEndcapP_ID"
name="VertexEndcapP"
type="athena_SimpleDiskTracker"
readout="VertexEndcapHits"
insideTrackingVolume="true"
reflect="false" vis="TrackerVis">
<position x="0" y="0" z="0.0*mm"/>
<layer id="1" vis="AnlOrange"
inner_z="VertexTrackerEndcapP_zmin + 0.5*VertexTrackerEndcap_delta"
inner_r="VertexTrackerEndcapP_rmin-3*mm"
outer_r="VertexTrackerEndcapP_rmax">
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
<slice material="Silicon" thickness="1.0*mm" vis="AnlOrange" sensitive="true"/>
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
</layer>
<layer id="2" vis="AnlOrange"
inner_z="VertexTrackerEndcapP_zmin + 1.5*VertexTrackerEndcap_delta"
inner_r="VertexTrackerEndcapP_rmin"
outer_r="VertexTrackerEndcapP_rmax">
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
<slice material="Silicon" thickness="1.0*mm" vis="AnlOrange" sensitive="true"/>
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
</layer>
</detector>
<detector id="VertexEndcapN_ID"
name="VertexEndcapN"
type="athena_SimpleDiskTracker"
readout="VertexEndcapHits"
insideTrackingVolume="true"
reflect="true" vis="TrackerVis">
<position x="0" y="0" z="-0.0*mm-1.0e-9*mm"/>
<layer id="1" vis="AnlOrange"
inner_z="VertexTrackerEndcapN_zmin + 0.5*VertexTrackerEndcap_delta"
inner_r="VertexTrackerEndcapN_rmin"
outer_r="VertexTrackerEndcapN_rmax">
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
<slice material="Silicon" thickness="1.0*mm" vis="AnlOrange" sensitive="true"/>
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
</layer>
<layer id="2" vis="AnlOrange"
inner_z="VertexTrackerEndcapN_zmin + 1.5*VertexTrackerEndcap_delta"
inner_r="VertexTrackerEndcapN_rmin"
outer_r="VertexTrackerEndcapN_rmax">
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
<slice material="Silicon" thickness="1.0*mm" vis="AnlOrange" sensitive="true"/>
<slice material="Air" thickness="1.0*mm" vis="AnlOrange" />
</layer>
</detector>
-->
</detectors>
<readouts>
......
compact/vertex_tracker_3layers.xml 0 → 100644
View file @ 1008a184
<?xml version="1.0" encoding="UTF-8"?>
<lccdd>
<define>
<comment>
Main parameters
</comment>
<constant name="SiVertexSensor_thickness" value="40*um"/>
<constant name="VertexBarrelMod_length" value="280.0*mm"/>
<constant name="VertexBarrelMod_rmin" value="3.3*cm"/>
<constant name="VertexBarrelMod_offset" value="1.05*cm"/>
<constant name="VertexBarrelMod_count" value="3"/>
<constant name="VertexBarrelShell_thickness" value="300.0*um"/>
<documentation>
#### Vertex Tracker Barrel Parameters
- The sensor modules are 2 half-cylinders.
- There are 2 sensitive layers
- Each sensor has a thickness is 40um
- There is an outer shell for structural support 300um thick.
- The ID of this shell is set (arbitrarily) to 10 cm.
##### Sensor layers
Currently there are 2 sensor layers. Each is composed of 2 half-cylinders modules with only 40um of silicon thickness.
##### Support shell
Simple carbon fiber support shell.
</documentation>
<constant name="VertexBarrelEnvelope_length" value="VertexTrackingRegion_length"/>
<constant name="VertexBarrelLayer_length" value="VertexBarrelMod_length + 1*um"/>
<constant name="VertexBarrelLayer_thickness" value="0.2*cm"/>
<constant name="VertexBarrelMod_thickness" value="0.1*cm"/>
<constant name="VertexBarrelMod1_rmin" value="VertexBarrelMod_rmin"/>
<constant name="VertexBarrelMod2_rmin" value="VertexBarrelMod_rmin + 1 * VertexBarrelMod_offset"/>
<constant name="VertexBarrelMod3_rmin" value="VertexBarrelMod_rmin + 2 * VertexBarrelMod_offset"/>
<constant name="VertexBarrelLayer1_rmin" value="VertexBarrelMod_rmin - VertexBarrelLayer_thickness/2.0"/>
<constant name="VertexBarrelLayer1_rmax" value="VertexBarrelLayer1_rmin + VertexBarrelLayer_thickness"/>
<constant name="VertexBarrelLayer2_rmin" value="VertexBarrelLayer1_rmin + 1 * VertexBarrelMod_offset"/>
<constant name="VertexBarrelLayer2_rmax" value="VertexBarrelLayer2_rmin + VertexBarrelLayer_thickness"/>
<constant name="VertexBarrelLayer3_rmin" value="VertexBarrelLayer1_rmin + 2 * VertexBarrelMod_offset"/>
<constant name="VertexBarrelLayer3_rmax" value="VertexBarrelLayer3_rmin + VertexBarrelLayer_thickness"/>
<comment>
"Support" is to "shell" like "layer" is to "module", and is need for the flat stave barrel implementation.
</comment>
<constant name="VertexBarrelShell_rmin" value="VertexBarrelMod_rmin + VertexBarrelMod_count * VertexBarrelMod_offset"/>
<constant name="VertexBarrelShell_rmax" value="VertexBarrelShell_rmin + VertexBarrelShell_thickness"/>
<constant name="VertexBarrelShell_length" value="VertexBarrelMod_length"/>
<constant name="VertexBarrelSupport_thickness" value="0.1*cm"/>
<constant name="VertexBarrelSupport_rmin" value="VertexBarrelShell_rmin-VertexBarrelSupport_thickness/2.0"/>
<constant name="VertexBarrelSupport_rmax" value="VertexBarrelSupport_rmin + VertexBarrelSupport_thickness"/>
<constant name="VertexBarrelSupport_length" value="VertexBarrelLayer_length"/>
<comment>
Extra parameters to approximate a cylinder as a set of skinny staves
due to ACTS limitations.
</comment>
<constant name="VertexBarrelStave_count" value="128"/>
<constant name="VertexBarrelStave1_width" value="2*VertexBarrelMod1_rmin * tan(180*degree/VertexBarrelStave_count)"/>
<constant name="VertexBarrelStave2_width" value="2*VertexBarrelMod2_rmin * tan(180*degree/VertexBarrelStave_count)"/>
<constant name="VertexBarrelStave3_width" value="2*VertexBarrelMod3_rmin * tan(180*degree/VertexBarrelStave_count)"/>
<constant name="VertexBarrelShellStave_width" value="2*VertexBarrelShell_rmin * tan(180*degree/VertexBarrelStave_count)"/>
</define>
<display>
</display>
<detectors>
<documentation level="5">
### Actual detectors
</documentation>
<detector
id="VertexBarrel_0_ID"
name="VertexBarrel"
type="athena_VertexBarrel"
readout="VertexBarrelHits"
insideTrackingVolume="true">
<dimensions
rmin="VertexBarrelLayer1_rmin"
rmax="VertexBarrelSupport_rmax"
length="VertexBarrelEnvelope_length" />
<comment>Vertex Barrel Modules</comment>
<module name="Module1" vis="VertexLayerVis">
<module_component name="ITS3"
material="Silicon"
sensitive="true"
width="VertexBarrelStave1_width"
length="VertexBarrelMod_length"
thickness="SiVertexSensor_thickness"
vis="VertexLayerVis" />
</module>
<module name="Module2" vis="VertexLayerVis">
<module_component name="ITS3"
material="Silicon"
sensitive="true"
width="VertexBarrelStave2_width"
length="VertexBarrelMod_length"
thickness="SiVertexSensor_thickness"
vis="VertexLayerVis" />
</module>
<module name="Module3" vis="VertexLayerVis">
<module_component name="ITS3"
material="Silicon"
sensitive="true"
width="VertexBarrelStave3_width"
length="VertexBarrelMod_length"
thickness="SiVertexSensor_thickness"
vis="VertexLayerVis" />
</module>
<module name="SupportShell" vis="VertexSupportVis">
<module_component name="CF Shell"
material="CarbonFiber"
sensitive="false"
width="VertexBarrelShellStave_width"
length="VertexBarrelShell_length"
thickness="VertexBarrelShell_thickness"
vis="VertexSupportVis" />
</module>
<comment> Layers composed of many arrayed modules </comment>
<layer module="Module1" id="1" vis="VertexLayerVis">
<barrel_envelope
inner_r="VertexBarrelLayer1_rmin"
outer_r="VertexBarrelLayer1_rmax"
z_length="VertexBarrelLayer_length" />
<layer_material surface="outer" binning="binPhi,binZ" bins0="VertexBarrelStave_count" bins1="100" />
<comment>
phi0 : Starting phi of first module.
phi_tilt : Phi tilt of a module.
rc : Radius of the module center.
nphi : Number of modules in phi.
rphi_dr : The delta radius of every other module.
z0 : Z position of first module in phi.
nz : Number of modules to place in z.
dr : Radial displacement parameter, of every other module.
</comment>
<rphi_layout phi_tilt="0.0*degree" nphi="VertexBarrelStave_count" phi0="0.0" rc="VertexBarrelMod1_rmin" dr="0.0 * mm"/>
<z_layout dr="0.0 * mm" z0="0.0 * mm" nz="1"/>
</layer>
<layer module="Module2" id="2" vis="VertexLayerVis">
<barrel_envelope
inner_r="VertexBarrelLayer2_rmin"
outer_r="VertexBarrelLayer2_rmax"
z_length="VertexBarrelLayer_length" />
<layer_material surface="outer" binning="binPhi,binZ" bins0="VertexBarrelStave_count" bins1="100" />
<rphi_layout phi_tilt="0.0*degree" nphi="VertexBarrelStave_count" phi0="0.0" rc="VertexBarrelMod2_rmin" dr="0.0 * mm"/>
<z_layout dr="0.0 * mm" z0="0.0 * mm" nz="1"/>
</layer>
<layer module="Module3" id="3" vis="VertexLayerVis">
<barrel_envelope
inner_r="VertexBarrelLayer3_rmin"
outer_r="VertexBarrelLayer3_rmax"
z_length="VertexBarrelLayer_length" />
<layer_material surface="outer" binning="binPhi,binZ" bins0="VertexBarrelStave_count" bins1="100" />
<rphi_layout phi_tilt="0.0*degree" nphi="VertexBarrelStave_count" phi0="0.0" rc="VertexBarrelMod3_rmin" dr="0.0 * mm"/>
<z_layout dr="0.0 * mm" z0="0.0 * mm" nz="1"/>
</layer>
<layer module="SupportShell" id="VertexBarrelMod_count + 1" vis="VertexSupportLayerVis">
<barrel_envelope
inner_r="VertexBarrelSupport_rmin"
outer_r="VertexBarrelSupport_rmax"
z_length="VertexBarrelSupport_length" />
<layer_material surface="outer" binning="binPhi,binZ" bins0="VertexBarrelStave_count" bins1="100" />
<rphi_layout phi_tilt="0.0*degree" nphi="VertexBarrelStave_count" phi0="0.0" rc="VertexBarrelShell_rmin" dr="0.0 * mm"/>
<z_layout dr="0.0 * mm" z0="0.0 * mm" nz="1"/>
</layer>
</detector>
</detectors>
<readouts>
<readout name="VertexBarrelHits">
<segmentation type="CartesianGridXY" grid_size_x="0.010*mm" grid_size_y="0.010*mm" />
<id>system:8,layer:4,module:12,sensor:2,x:32:-16,y:-16</id>
</readout>
</readouts>
</lccdd>
src/BarrelBarDetectorWithSideFrame_geo.cpp
View file @ 1008a184
......@@ -16,6 +16,7 @@
#include "DDRec/DetectorData.h"
#include "DDRec/Surface.h"
#include "XML/Layering.h"
#include <array>
using namespace std;
using namespace dd4hep;
......
src/BarrelCalorimeterInterlayers_geo.cpp
View file @ 1008a184
......@@ -21,7 +21,7 @@ using namespace dd4hep;
using namespace dd4hep::detail;
typedef ROOT::Math::XYPoint Point;
// fiber placement helpers, defined in BarrelCalorimeterHybrid_geo
// fiber placement helpers, defined below
vector<vector<Point>> fiberPositions(double radius, double x_spacing, double z_spacing,
double x, double z, double phi, double spacing_tol = 1e-2);
std::pair<int, int> getNdivisions(double x, double z, double dx, double dz);
......@@ -277,23 +277,22 @@ void buildSupport(Detector& desc, Volume &mod_vol, xml_comp_t x_support,
double beam_width = 2. * trd_x1_support / (n_beams + 1); // quick hack to make some gap between T beams
double beam_gap = getAttrOrDefault(x_support, _Unicode(beam_gap), 3.*cm);
// build H-shape beam
// build T-shape beam
double beam_space_x = beam_width + beam_gap;
double beam_space_z = support_thickness - beam_thickness;
double cross_thickness = support_thickness - 2.*beam_thickness;
double beam_pos_z = -beam_thickness / 2.;
double beam_center_z = support_thickness / 2. + beam_pos_z;
double cross_thickness = support_thickness - beam_thickness;
double beam_pos_z = beam_thickness / 2.;
double beam_center_z = support_thickness / 2. - beam_pos_z;
Box beam_vert_s(beam_thickness / 2. , trd_y, cross_thickness / 2.);
Box beam_vert_s(beam_thickness / 2., trd_y, cross_thickness / 2.);
Box beam_hori_s(beam_width / 2., trd_y, beam_thickness / 2.);
UnionSolid T_beam_s(beam_hori_s, beam_vert_s, Position(0., 0., beam_space_z / 2.));
UnionSolid H_beam_s(T_beam_s, beam_hori_s, Position(0., 0., support_thickness - beam_thickness));
Volume H_beam_vol("H_beam", H_beam_s, desc.material(x_support.materialStr()));
UnionSolid T_beam_s(beam_hori_s, beam_vert_s, Position(0., 0., support_thickness / 2.));
Volume H_beam_vol("H_beam", T_beam_s, desc.material(x_support.materialStr()));
H_beam_vol.setVisAttributes(desc, x_support.visStr());
// place H beams first
double beam_start_x = - (n_beams - 1) * (beam_width + beam_gap) / 2.;
for (int i = 0; i < n_beams; ++i) {
Position beam_pos(beam_start_x + i * (beam_width + beam_gap), 0., - support_thickness / 2. - beam_pos_z);
Position beam_pos(beam_start_x + i * (beam_width + beam_gap), 0., - support_thickness / 2. + beam_pos_z);
env_vol.placeVolume(H_beam_vol, beam_pos);
}
......@@ -303,10 +302,10 @@ void buildSupport(Detector& desc, Volume &mod_vol, xml_comp_t x_support,
Volume cross_vol("cross_center_beam", cross_s, desc.material(x_support.materialStr()));
cross_vol.setVisAttributes(desc, x_support.visStr());
for (int i = 0; i < n_beams - 1; ++i) {
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), 0., 0.));
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), 0., beam_pos_z));
for (int j = 1; j < n_cross_supports; j++) {
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), -j * grid_size, 0.));
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), j * grid_size, 0.));
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), -j * grid_size, beam_pos_z));
env_vol.placeVolume(cross_vol, Position(beam_start_x + beam_space_x * (i + 0.5), j * grid_size, beam_pos_z));
}
}
......@@ -315,24 +314,135 @@ void buildSupport(Detector& desc, Volume &mod_vol, xml_comp_t x_support,
double cross_edge_x = trd_x1_support + beam_start_x - beam_thickness / 2.;
double cross_trd_x1 = cross_edge_x + std::tan(hphi) * beam_thickness;
double cross_trd_x2 = cross_trd_x1 + 2.* std::tan(hphi) * cross_thickness;
double edge_pos_x = beam_start_x - beam_thickness / 2. - cross_trd_x1 / 2.;
double edge_pos_x = beam_start_x - cross_trd_x1 / 2. - beam_thickness / 2;
Trapezoid cross_s2_trd (cross_trd_x1 / 2., cross_trd_x2 / 2.,
beam_thickness / 2., beam_thickness / 2., cross_thickness / 2.);
Box cross_s2_box ((cross_trd_x2 - cross_trd_x1)/2., beam_thickness / 2., cross_thickness / 2.);
SubtractionSolid cross_s2(cross_s2_trd, cross_s2_box, Position((cross_trd_x1 + cross_trd_x2) / 4., 0., 0.));
Box cross_s2_box ((cross_trd_x2 - cross_trd_x1)/4., beam_thickness / 2., cross_thickness / 2.);
SubtractionSolid cross_s2(cross_s2_trd, cross_s2_box, Position((cross_trd_x2 + cross_trd_x1)/4., 0., 0.));
Volume cross_vol2("cross_edge_beam", cross_s2, desc.material(x_support.materialStr()));
cross_vol2.setVisAttributes(desc, x_support.visStr());
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, 0., 0.));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, 0., 0.) * RotationZ(M_PI)));
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, 0., beam_pos_z));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, 0., beam_pos_z) * RotationZ(M_PI)));
for (int j = 1; j < n_cross_supports; j++) {
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, -j * grid_size, 0.));
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, j * grid_size, 0.));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, -j * grid_size, 0.) * RotationZ(M_PI)));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, j * grid_size, 0.) * RotationZ(M_PI)));
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, -j * grid_size, beam_pos_z));
env_vol.placeVolume(cross_vol2, Position(edge_pos_x, j * grid_size, beam_pos_z));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, -j * grid_size, beam_pos_z) * RotationZ(M_PI)));
env_vol.placeVolume(cross_vol2, Transform3D(Translation3D(-edge_pos_x, j * grid_size, beam_pos_z) * RotationZ(M_PI)));
}
mod_vol.placeVolume(env_vol, Position(0.0, 0.0, l_pos_z + support_thickness/2.));
}
// Fill fiber lattice into trapezoid starting from position (0,0) in x-z coordinate system
vector<vector<Point>> fiberPositions(double radius, double x_spacing, double z_spacing, double x, double z, double phi,
double spacing_tol)
{
// z_spacing - distance between fiber layers in z
// x_spacing - distance between fiber centers in x
// x - half-length of the shorter (bottom) base of the trapezoid
// z - height of the trapezoid
// phi - angle between z and trapezoid arm
vector<vector<Point>> positions;
int z_layers = floor((z / 2 - radius - spacing_tol) / z_spacing); // number of layers that fit in z/2
double z_pos = 0.;
double x_pos = 0.;
for (int l = -z_layers; l < z_layers + 1; l++) {
vector<Point> xline;
z_pos = l * z_spacing;
double x_max = x + (z / 2. + z_pos) * tan(phi) - spacing_tol; // calculate max x at particular z_pos
(l % 2 == 0) ? x_pos = 0. : x_pos = x_spacing / 2; // account for spacing/2 shift
while (x_pos < (x_max - radius)) {
xline.push_back(Point(x_pos, z_pos));
if (x_pos != 0.)
xline.push_back(Point(-x_pos, z_pos)); // using symmetry around x=0
x_pos += x_spacing;
}
// Sort fiber IDs for a better organization
sort(xline.begin(), xline.end(), [](const Point& p1, const Point& p2) { return p1.x() < p2.x(); });
positions.emplace_back(std::move(xline));
}
return positions;
}
// Calculate number of divisions for the readout grid for the fiber layers
std::pair<int, int> getNdivisions(double x, double z, double dx, double dz)
{
// x and z defined as in vector<Point> fiberPositions
// dx, dz - size of the grid in x and z we want to get close to with the polygons
// See also descripltion when the function is called
double SiPMsize = 13.0 * mm;
double grid_min = SiPMsize + 3.0 * mm;
if (dz < grid_min) {
dz = grid_min;
}
if (dx < grid_min) {
dx = grid_min;
}
int nfit_cells_z = floor(z / dz);
int n_cells_z = nfit_cells_z;
if (nfit_cells_z == 0)
n_cells_z++;
int nfit_cells_x = floor((2 * x) / dx);
int n_cells_x = nfit_cells_x;
if (nfit_cells_x == 0)
n_cells_x++;
return std::make_pair(n_cells_x, n_cells_z);
}
// Calculate dimensions of the polygonal grid in the cartesian coordinate system x-z
vector<tuple<int, Point, Point, Point, Point>> gridPoints(int div_x, int div_z, double x, double z, double phi)
{
// x, z and phi defined as in vector<Point> fiberPositions
// div_x, div_z - number of divisions in x and z
double dz = z / div_z;
std::vector<std::tuple<int, Point, Point, Point, Point>> points;
for (int iz = 0; iz < div_z + 1; iz++) {
for (int ix = 0; ix < div_x + 1; ix++) {
double A_z = -z / 2 + iz * dz;
double B_z = -z / 2 + (iz + 1) * dz;
double len_x_for_z = 2 * (x + iz * dz * tan(phi));
double len_x_for_z_plus_1 = 2 * (x + (iz + 1) * dz * tan(phi));
double dx_for_z = len_x_for_z / div_x;
double dx_for_z_plus_1 = len_x_for_z_plus_1 / div_x;
double A_x = -len_x_for_z / 2. + ix * dx_for_z;
double B_x = -len_x_for_z_plus_1 / 2. + ix * dx_for_z_plus_1;
double C_z = B_z;
double D_z = A_z;
double C_x = B_x + dx_for_z_plus_1;
double D_x = A_x + dx_for_z;
int id = ix + div_x * iz;
auto A = Point(A_x, A_z);
auto B = Point(B_x, B_z);
auto C = Point(C_x, C_z);
auto D = Point(D_x, D_z);
// vertex points filled in the clock-wise direction
points.push_back(make_tuple(id, A, B, C, D));
}
}
return points;
}
DECLARE_DETELEMENT(athena_EcalBarrelInterlayers, create_detector)
src/BarrelTrackerWithFrame_geo.cpp
View file @ 1008a184
......@@ -13,10 +13,15 @@
#include "DDRec/DetectorData.h"
#include "XML/Layering.h"
#include "XML/Utilities.h"
#include <array>
#if defined(USE_ACTSDD4HEP)
#include "ActsDD4hep/ActsExtension.hpp"
#include "ActsDD4hep/ConvertMaterial.hpp"
#else
#include "Acts/Plugins/DD4hep/ActsExtension.hpp"
#include "Acts/Plugins/DD4hep/ConvertDD4hepMaterial.hpp"
#endif
using namespace std;
using namespace dd4hep;
......@@ -166,24 +171,26 @@ static Ref_t create_BarrelTrackerWithFrame(Detector& description, xml_h e, Sensi
double thickness_so_far = 0.0;
double thickness_sum = -total_thickness/2.0;
for (xml_coll_t ci(x_mod, _U(module_component)); ci; ++ci, ++ncomponents) {
xml_comp_t x_comp = ci;
xml_comp_t x_pos = x_comp.position(false);
xml_comp_t x_rot = x_comp.rotation(false);
string c_nam = _toString(ncomponents, "component%d");
Box c_box(x_comp.width() / 2, x_comp.length() / 2, x_comp.thickness() / 2);
Volume c_vol(c_nam, c_box, description.material(x_comp.materialStr()));
xml_comp_t x_comp = ci;
xml_comp_t x_pos = x_comp.position(false);
xml_comp_t x_rot = x_comp.rotation(false);
const string c_nam = _toString(ncomponents, "component%d");
Box c_box(x_comp.width() / 2, x_comp.length() / 2, x_comp.thickness() / 2);
Volume c_vol(c_nam, c_box, description.material(x_comp.materialStr()));
// Utility variable for the relative z-offset based off the previous components
const double zoff = thickness_sum+x_comp.thickness() / 2.0;
if (x_pos && x_rot) {
Position c_pos(x_pos.x(0), x_pos.y(0), x_pos.z(0));
Position c_pos(x_pos.x(0), x_pos.y(0), x_pos.z(0) + zoff);
RotationZYX c_rot(x_rot.z(0), x_rot.y(0), x_rot.x(0));
pv = m_vol.placeVolume(c_vol, Transform3D(c_rot, c_pos));
} else if (x_rot) {
Position c_pos(0, 0, thickness_sum + x_comp.thickness() / 2.0);
pv = m_vol.placeVolume(c_vol, Transform3D(RotationZYX(x_rot.z(0), x_rot.y(0), x_rot.x(0)),c_pos));
Position c_pos(0, 0, zoff);
pv = m_vol.placeVolume(c_vol, Transform3D(RotationZYX(x_rot.z(0), x_rot.y(0), x_rot.x(0)), c_pos));
} else if (x_pos) {
pv = m_vol.placeVolume(c_vol, Position(x_pos.x(0), x_pos.y(0), x_pos.z(0)));
pv = m_vol.placeVolume(c_vol, Position(x_pos.x(0), x_pos.y(0), x_pos.z(0) + zoff));
} else {
pv = m_vol.placeVolume(c_vol, Position(0,0,thickness_sum+x_comp.thickness()/2.0));
pv = m_vol.placeVolume(c_vol, Position(0, 0, zoff));
}
c_vol.setRegion(description, x_comp.regionStr());
c_vol.setLimitSet(description, x_comp.limitsStr());
......@@ -217,6 +224,11 @@ static Ref_t create_BarrelTrackerWithFrame(Detector& description, xml_h e, Sensi
}
thickness_sum += x_comp.thickness();
thickness_so_far += x_comp.thickness();
// apply relative offsets in z-position used to stack components side-by-side
if (x_pos) {
thickness_sum += x_pos.z(0);
thickness_so_far += x_pos.z(0);
}
}
}
......@@ -341,3 +353,4 @@ static Ref_t create_BarrelTrackerWithFrame(Detector& description, xml_h e, Sensi
DECLARE_DETELEMENT(BarrelTrackerWithFrame, create_BarrelTrackerWithFrame)
DECLARE_DETELEMENT(athena_TrackerBarrel, create_BarrelTrackerWithFrame)
DECLARE_DETELEMENT(athena_VertexBarrel, create_BarrelTrackerWithFrame)
DECLARE_DETELEMENT(athena_TOFBarrel, create_BarrelTrackerWithFrame)
src/CompositeTracker_geo.cpp 0 → 100644
View file @ 1008a184
//==========================================================================
// Based off DD4hep_SubDetectorAssembly
//
// This is a simple plugin to allow compositing different detectors
// into a single barrel or endcap for ACTS
// Note: positive/negative position strings differentiate between
// positive/negative endcaps
//--------------------------------------------------------------------------
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include "XML/Utilities.h"
#if defined(USE_ACTSDD4HEP)
#include "ActsDD4hep/ActsExtension.hpp"
#else
#include "Acts/Plugins/DD4hep/ActsExtension.hpp"
#endif
using namespace dd4hep;
using namespace dd4hep::detail;
static Ref_t create_element(Detector& description, xml_h e, Ref_t)
{
xml_det_t x_det(e);
const std::string det_name = x_det.nameStr();
DetElement sdet(det_name, x_det.id());
Volume vol;
Position pos;
const bool usePos = x_det.hasChild(_U(position));
sdet.setType("compound");
xml::setDetectorTypeFlag(e, sdet);
const std::string actsType = getAttrOrDefault(x_det, _Unicode(actsType), "endcap");
printout(DEBUG, det_name, "+++ Creating composite tracking detector (type: " + actsType + ")");
assert(actsType == "barrel" || actsType == "endcap");
// ACTS extension
{
Acts::ActsExtension* detWorldExt = new Acts::ActsExtension();
detWorldExt->addType(actsType, "detector");
sdet.addExtension<Acts::ActsExtension>(detWorldExt);
}
if (usePos) {
pos = Position(x_det.position().x(), x_det.position().y(), x_det.position().z());
}
vol = Assembly(det_name);
vol.setAttributes(description, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());
Volume mother = description.pickMotherVolume(sdet);
PlacedVolume pv;
if (usePos) {
pv = mother.placeVolume(vol, pos);
} else {
pv = mother.placeVolume(vol);
}
sdet.setPlacement(pv);
for (xml_coll_t c(x_det, _U(composite)); c; ++c) {
xml_dim_t component = c;
const std::string nam = component.nameStr();
description.declareParent(nam, sdet);
}
return sdet;
}
DECLARE_DETELEMENT(athena_CompositeTracker, create_element)
src/CylinderTrackerBarrel_geo.cpp
View file @ 1008a184
......@@ -2,11 +2,17 @@
// Specialized generic detector constructor
//==========================================================================
#include "Acts/Plugins/DD4hep/ActsExtension.hpp"
#include "Acts/Plugins/DD4hep/ConvertDD4hepMaterial.hpp"
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#if defined(USE_ACTSDD4HEP)
#include "ActsDD4hep/ActsExtension.hpp"
#include "ActsDD4hep/ConvertMaterial.hpp"
#else
#include "Acts/Plugins/DD4hep/ActsExtension.hpp"
#include "Acts/Plugins/DD4hep/ConvertDD4hepMaterial.hpp"
#endif
using namespace std;
using namespace dd4hep;
using namespace dd4hep::detail;
......
src/DIRC_geo.cpp
View file @ 1008a184
......@@ -162,7 +162,7 @@ static Ref_t createDetector(Detector& desc, xml_h e, SensitiveDetector sens)
Volume lGlue("lGlue", gGlue, desc.material(xml_glue.materialStr()));
lGlue.setVisAttributes(desc.visAttributes(xml_glue.visStr()));
sens.setType("photoncounter");
sens.setType("tracker");
lBar.setSensitiveDetector(sens);
int bars_repeat_z = 4; // TODO parametrize!
......@@ -320,7 +320,7 @@ static Ref_t createDetector(Detector& desc, xml_h e, SensitiveDetector sens)
// Box bar_geo(bar_thicknes / 2., bar_width / 2., bar_length / 2.);
// Volume bar_volume("cb_DIRC_bars_Logix", bar_geo, bar_material);
// bar_volume.setVisAttributes(desc.visAttributes(xml_det.visStr()));
// sens.setType("photoncounter");
// sens.setType("tracker");
// bar_volume.setSensitiveDetector(sens);
// for (int ia = 0; ia < bar_count; ia++) {
......
src/DRich_geo.cpp → src/DRICH_geo.cpp
View file @ 1008a184
......@@ -77,8 +77,8 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
double mirrorRmin = mirrorElem.attr<double>(_Unicode(rmin));
double mirrorRmax = mirrorElem.attr<double>(_Unicode(rmax));
double mirrorPhiw = mirrorElem.attr<double>(_Unicode(phiw));
double focusTuneLong = mirrorElem.attr<double>(_Unicode(focus_tune_long));
double focusTunePerp = mirrorElem.attr<double>(_Unicode(focus_tune_perp));
double focusTuneZ = mirrorElem.attr<double>(_Unicode(focus_tune_z));
double focusTuneX = mirrorElem.attr<double>(_Unicode(focus_tune_x));
// - sensor module
auto sensorElem = detElem.child(_Unicode(sensors)).child(_Unicode(module));
auto sensorMat = desc.material(sensorElem.attr<std::string>(_Unicode(material)));
......@@ -132,6 +132,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
// snout solids
double boreDelta = vesselRmin1 - vesselRmin0;
double snoutDelta = vesselRmax1 - vesselRmax0;
Cone vesselSnout(
snoutLength/2.0,
vesselRmin0,
......@@ -216,7 +217,60 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
// sensitive detector type
sens.setType("photoncounter");
sens.setType("tracker");
// BUILD RADIATOR ====================================================================
// attributes
double airGap = 0.01*mm; // air gap between aerogel and filter (FIXME? actually it's currently a gas gap)
// solid and volume: create aerogel and filter
Cone aerogelSolid(
aerogelThickness/2,
radiatorRmin,
radiatorRmax,
radiatorRmin + boreDelta * aerogelThickness / vesselLength,
radiatorRmax + snoutDelta * aerogelThickness / snoutLength
);
Cone filterSolid(
filterThickness/2,
radiatorRmin + boreDelta * (aerogelThickness + airGap) / vesselLength,
radiatorRmax + snoutDelta * (aerogelThickness + airGap) / snoutLength,
radiatorRmin + boreDelta * (aerogelThickness + airGap + filterThickness) / vesselLength,
radiatorRmax + snoutDelta * (aerogelThickness + airGap + filterThickness) / snoutLength
);
Volume aerogelVol( detName+"_aerogel", aerogelSolid, aerogelMat );
Volume filterVol( detName+"_filter", filterSolid, filterMat );
aerogelVol.setVisAttributes(aerogelVis);
filterVol.setVisAttributes(filterVis);
// aerogel placement and surface properties
// TODO [low-priority]: define skin properties for aerogel and filter
auto radiatorPos = Position(0., 0., radiatorFrontplane) + originFront;
auto aerogelPV = gasvolVol.placeVolume(aerogelVol,
Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) // re-center to originFront
* RotationY(radiatorPitch) // change polar angle to specified pitch // FIXME: probably broken (not yet in use anyway)
);
DetElement aerogelDE(det, "aerogel_de", 0);
aerogelDE.setPlacement(aerogelPV);
//SkinSurface aerogelSkin(desc, aerogelDE, "mirror_optical_surface", aerogelSurf, aerogelVol);
//aerogelSkin.isValid();
// filter placement and surface properties
if(!debug_optics) {
auto filterPV = gasvolVol.placeVolume(filterVol,
Translation3D(0., 0., airGap) // add an air gap
* Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) // re-center to originFront
* RotationY(radiatorPitch) // change polar angle
* Translation3D(0., 0., (aerogelThickness+filterThickness)/2.) // move to aerogel backplane
);
DetElement filterDE(det, "filter_de", 0);
filterDE.setPlacement(filterPV);
//SkinSurface filterSkin(desc, filterDE, "mirror_optical_surface", filterSurf, filterVol);
//filterSkin.isValid();
};
// SECTOR LOOP //////////////////////////////////
......@@ -230,43 +284,6 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
std::string secName = "sec" + std::to_string(isec);
// BUILD RADIATOR ====================================================================
// solid and volume: create aerogel and filter sectors
Tube aerogelSolid(radiatorRmin, radiatorRmax, aerogelThickness/2, -radiatorPhiw/2.0, radiatorPhiw/2.0);
Tube filterSolid( radiatorRmin, radiatorRmax, filterThickness/2, -radiatorPhiw/2.0, radiatorPhiw/2.0);
Volume aerogelVol( detName+"_aerogel_"+secName, aerogelSolid, aerogelMat );
Volume filterVol( detName+"_filter_"+secName, filterSolid, filterMat );
aerogelVol.setVisAttributes(aerogelVis);
filterVol.setVisAttributes(filterVis);
// aerogel placement and surface properties
// TODO [low-priority]: define skin properties for aerogel and filter
auto radiatorPos = Position(0., 0., radiatorFrontplane) + originFront;
auto aerogelPV = gasvolVol.placeVolume(aerogelVol,
RotationZ(sectorRotation) // rotate about beam axis to sector
* Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) // re-center to originFront
* RotationY(radiatorPitch) // change polar angle to specified pitch
);
DetElement aerogelDE(det, Form("aerogel_de%d", isec), isec);
aerogelDE.setPlacement(aerogelPV);
//SkinSurface aerogelSkin(desc, aerogelDE, Form("mirror_optical_surface%d", isec), aerogelSurf, aerogelVol);
//aerogelSkin.isValid();
// filter placement and surface properties
if(!debug_optics) {
auto filterPV = gasvolVol.placeVolume(filterVol,
RotationZ(sectorRotation) // rotate about beam axis to sector
* Translation3D(radiatorPos.x(), radiatorPos.y(), radiatorPos.z()) // re-center to originFront
* RotationY(radiatorPitch) // change polar angle
* Translation3D(0., 0., (aerogelThickness+filterThickness)/2.) // move to aerogel backplane
);
DetElement filterDE(det, Form("filter_de%d", isec), isec);
filterDE.setPlacement(filterPV);
//SkinSurface filterSkin(desc, filterDE, Form("mirror_optical_surface%d", isec), filterSurf, filterVol);
//filterSkin.isValid();
};
// BUILD SENSORS ====================================================================
......@@ -370,7 +387,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
// properties
sensorPV.addPhysVolID("sector", isec).addPhysVolID("module", imod);
DetElement sensorDE(det, Form("sensor_de%d_%d", isec, imod), 10000*isec+imod);
DetElement sensorDE(det, Form("sensor_de%d_%d", isec, imod), (imod<<3)|isec ); // id must match IRTAlgorithm usage
sensorDE.setPlacement(sensorPV);
if(!debug_optics) {
SkinSurface sensorSkin(desc, sensorDE, Form("sensor_optical_surface%d", isec), sensorSurf, sensorVol);
......@@ -395,16 +412,12 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
// BUILD MIRRORS ====================================================================
// attributes, re-defined w.r.t. IP, needed for mirror positioning
double zS = sensorSphCenterZ + vesselZmin; // sensor sphere attributes
double xS = sensorSphCenterX;
double rS = sensorSphRadius;
double B = vesselZmax - mirrorBackplane; // distance between IP and mirror back plane
// derive spherical mirror parameters `(zM,xM,rM)`, for given image point
// coordinates `(zI,xI)` and `dO`, defined as the z-distance between the
// object and the mirror surface; all coordinates are specified w.r.t. the
// object point coordinates
// object and the mirror surface
// - all coordinates are specified w.r.t. the object point coordinates
// - this is point-to-point focusing, but it can be used to effectively steer
// parallel-to-point focusing
double zM,xM,rM;
auto FocusMirror = [&zM,&xM,&rM](double zI, double xI, double dO) {
zM = dO*zI / (2*dO-zI);
......@@ -412,6 +425,12 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
rM = dO - zM;
};
// attributes, re-defined w.r.t. IP, needed for mirror positioning
double zS = sensorSphCenterZ + vesselZmin; // sensor sphere attributes
double xS = sensorSphCenterX;
double rS = sensorSphRadius;
double B = vesselZmax - mirrorBackplane; // distance between IP and mirror back plane
// focus 1: set mirror to focus IP on center of sensor sphere `(zS,xS)`
/*double zF = zS;
double xF = xS;
......@@ -423,6 +442,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
// - `focusTuneLong` is the distance to move, given as a fraction of `sensorSphRadius`
// - `focusTuneLong==0` means `(zF,xF)==(zS,xS)`
// - `focusTuneLong==1` means `(zF,xF)` will be on the sensor sphere, near the centroid
/*
double zC = sensorCentroidZ + vesselZmin;
double xC = sensorCentroidX;
double slopeF = (xC-xS) / (zC-zS);
......@@ -436,17 +456,13 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
zF += focusTunePerp * sensorSphRadius * std::cos(M_PI/2-thetaF);
xF += focusTunePerp * sensorSphRadius * std::sin(M_PI/2-thetaF);
FocusMirror(zF,xF,B);
/*
printf("(zC,xC) = ( %.2f, %.2f )\n",zC,xC);
printf("zS = %f\n",zS);
printf("xS = %f\n",xS);
printf("B = %f\n",B);
printf("zM = %f\n",zM);
printf("xM = %f\n",xM);
printf("rM = %f\n",rM);
*/
// focus 4: use (z,x) coordinates for tune parameters
double zF = zS + focusTuneZ;
double xF = xS + focusTuneX;
FocusMirror(zF,xF,B);
// re-define mirror attributes to be w.r.t vessel front plane
double mirrorCenterZ = zM - vesselZmin;
double mirrorCenterX = xM;
......@@ -471,6 +487,18 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
40*degree
);
/* CAUTION: if any of the relative placements or boolean operations below
* are changed, you MUST make sure this does not break access to the sphere
* primitive and positioning in Juggler `IRTAlgorithm`; cross check the
* mirror sphere attributes carefully!
*/
/*
// PRINT MIRROR ATTRIBUTES (before any sector z-rotation)
printf("zM = %f\n",zM); // sphere centerZ, w.r.t. IP
printf("xM = %f\n",xM); // sphere centerX, w.r.t. IP
printf("rM = %f\n",rM); // sphere radius
*/
// mirror placement transformation (note: transformations are in reverse order)
auto mirrorPos = Position(mirrorCenterX, 0., mirrorCenterZ) + originFront;
Transform3D mirrorPlacement(
......
src/FieldMapBrBz.cpp
View file @ 1008a184
#include <DD4hep/DetFactoryHelper.h>
#include <DD4hep/FieldTypes.h>
#include <DD4hep/Printout.h>
#include <XML/Utilities.h>
#include <cstdlib>
......@@ -12,6 +13,8 @@
#include <tuple>
namespace fs = std::filesystem;
#include "FileLoaderHelper.h"
using namespace dd4hep;
......@@ -19,7 +22,7 @@ using namespace dd4hep;
class FieldMapBrBz : public dd4hep::CartesianField::Object
{
public:
FieldMapBrBz(const std::string &field_type);
FieldMapBrBz(const std::string &field_type = "magnetic");
void Configure(double rmin, double rmax, double rstep, double zmin, double zmax, double zstep);
void LoadMap(const std::string &map_file, double scale);
void GetIndices(double r, double z, int &ir, int &iz, double &dr, double &dz);
......@@ -35,7 +38,7 @@ private:
};
// constructor
FieldMapBrBz::FieldMapBrBz(const std::string &field_type = "magnetic")
FieldMapBrBz::FieldMapBrBz(const std::string &field_type)
{
std::string ftype = field_type;
for (auto &c : ftype) { c = tolower(c); }
......@@ -185,20 +188,16 @@ static Ref_t create_field_map_brbz(Detector & /*lcdd*/, xml::Handle_t handle)
std::string field_map_file = x_par.attr<std::string>(_Unicode(field_map));
std::string field_map_url = x_par.attr<std::string>(_Unicode(url));
EnsureFileFromURLExists(field_map_url,field_map_file);
double field_map_scale = x_par.attr<double>(_Unicode(scale));
if( !fs::exists(fs::path(field_map_file)) ) {
auto ret = std::system(("mkdir -p fieldmaps && "
"wget " +
field_map_url + " -O " + field_map_file).c_str());
if (!fs::exists(fs::path(field_map_file))) {
std::cerr << "ERROR: file, " << field_map_file << ", does not exist\n";
std::quick_exit(1);
}
printout(ERROR, "FieldMapBrBz", "file " + field_map_file + " does not exist");
printout(ERROR, "FieldMapBrBz", "use a FileLoader plugin before the field element");
std::_Exit(EXIT_FAILURE);
}
auto map = new FieldMapBrBz(field_type);
map->Configure(r_dim.rmin(), r_dim.rmax(), r_dim.step(), z_dim.zmin(), z_dim.zmax(), z_dim.step());
......
src/FileLoader.cpp 0 → 100644
View file @ 1008a184
#include <DD4hep/DetFactoryHelper.h>
#include <DD4hep/Primitives.h>
#include <DD4hep/Factories.h>
#include <DD4hep/Printout.h>
#include <XML/Utilities.h>
#include <fmt/core.h>
#include <filesystem>
#include <iostream>
#include <cstdlib>
#include <string>
#include "FileLoaderHelper.h"
using namespace dd4hep;
void usage(int argc, char** argv) {
std::cout <<
"Usage: -plugin <name> -arg [-arg] \n"
" name: factory name FileLoader \n"
" cache:<string> cache location (may be read-only) \n"
" file:<string> file location \n"
" url:<string> url location \n"
" cmd:<string> download command with {0} for url, {1} for output \n"
"\tArguments given: " << arguments(argc,argv) << std::endl;
std::exit(EINVAL);
}
// Plugin to download files
long load_file(
Detector& /* desc */,
int argc,
char** argv
) {
// argument parsing
std::string cache, file, url;
std::string cmd("curl --retry 5 -f {0} -o {1}");
for (int i = 0; i < argc && argv[i]; ++i) {
if (0 == std::strncmp("cache:", argv[i], 6)) cache = (argv[i] + 6);
else if (0 == std::strncmp("file:", argv[i], 5)) file = (argv[i] + 5);
else if (0 == std::strncmp("url:", argv[i], 4)) url = (argv[i] + 4);
else if (0 == std::strncmp("cmd:", argv[i], 4)) cmd = (argv[i] + 4);
else usage(argc, argv);
}
printout(DEBUG, "FileLoader", "arg cache: " + cache);
printout(DEBUG, "FileLoader", "arg file: " + file);
printout(DEBUG, "FileLoader", "arg url: " + url);
printout(DEBUG, "FileLoader", "arg cmd: " + cmd);
// if file or url is empty, do nothing
if (file.empty()) {
printout(WARNING, "FileLoader", "no file specified");
return 0;
}
if (url.empty()) {
printout(WARNING, "FileLoader", "no url specified");
return 0;
}
EnsureFileFromURLExists(url, file, cache, cmd);
return 0;
}
DECLARE_APPLY(FileLoader, load_file)
src/FileLoaderHelper.h 0 → 100644
View file @ 1008a184
#pragma once
#include <DD4hep/DetFactoryHelper.h>
#include <DD4hep/Primitives.h>
#include <DD4hep/Factories.h>
#include <DD4hep/Printout.h>
#include <fmt/core.h>
#include <filesystem>
#include <iostream>
#include <cstdlib>
#include <string>
namespace fs = std::filesystem;
using namespace dd4hep;
// Function to download files
inline void
EnsureFileFromURLExists(
std::string url,
std::string file,
std::string cache = "",
std::string cmd = "curl --retry 5 -f {0} -o {1}"
) {
// parse cache for environment variables
auto pos = std::string::npos;
while ((pos = cache.find('$')) != std::string::npos) {
auto after = cache.find_first_not_of(
"ABCDEFGHIJKLMNOPQRSTUVWXYZ"
"abcdefghijklmnopqrstuvwxyz"
"0123456789"
"_",
pos + 1);
if (after == std::string::npos) after = cache.size(); // cache ends on env var
const std::string env_name(cache.substr(pos + 1, after - pos - 1));
auto env_ptr = std::getenv(env_name.c_str());
const std::string env_value(env_ptr != nullptr ? env_ptr : "");
cache.erase(pos, after - pos);
cache.insert(pos, env_value);
printout(INFO, "FileLoader", "$" + env_name + " -> " + env_value);
}
// create file path
fs::path file_path(file);
// create hash from url, hex of unsigned long long
std::string hash = fmt::format("{:016x}", dd4hep::detail::hash64(url)); // TODO: Use c++20 std::fmt
// create file parent path, if not exists
fs::path parent_path = file_path.parent_path();
if (!fs::exists(parent_path)) {
if (fs::create_directories(parent_path) == false) {
printout(ERROR, "FileLoader", "parent path " + parent_path.string() + " cannot be created");
printout(ERROR, "FileLoader", "check permissions and retry");
std::_Exit(EXIT_FAILURE);
}
}
// if file exists and is symlink to correct hash
fs::path hash_path(parent_path / hash);
if (fs::exists(file_path)
&& fs::equivalent(file_path, hash_path)) {
printout(INFO, "FileLoader", "Link " + file + " -> hash " + hash + " already exists");
return;
}
// if hash does not exist, we try to retrieve file from cache
if (!fs::exists(hash_path)) {
// recursive loop into cache directory
fs::path cache_path(cache);
printout(INFO, "FileLoader", "Cache " + cache_path.string());
if (fs::exists(cache_path)) {
for (auto const& dir_entry: fs::recursive_directory_iterator(cache_path)) {
if (!dir_entry.is_directory()) continue;
fs::path cache_dir_path = cache_path / dir_entry;
printout(INFO, "FileLoader", "Checking " + cache_dir_path.string());
fs::path cache_hash_path = cache_dir_path / hash;
if (fs::exists(cache_hash_path)) {
// symlink hash to cache/.../hash
printout(INFO, "FileLoader", "File " + file + " with hash " + hash + " found in " + cache_hash_path.string());
try {
fs::create_symlink(cache_hash_path, hash_path);
} catch (const fs::filesystem_error&) {
printout(ERROR, "FileLoader", "unable to link from " + hash_path.string() + " to " + cache_hash_path.string());
printout(ERROR, "FileLoader", "check permissions and retry");
std::_Exit(EXIT_FAILURE);
}
break;
}
}
}
}
// if hash does not exist, we try to retrieve file from url
if (!fs::exists(hash_path)) {
cmd = fmt::format(cmd, url, hash_path.c_str()); // TODO: Use c++20 std::fmt
printout(INFO, "FileLoader", "Downloading " + file + " as hash " + hash + " with " + cmd);
// run cmd
auto ret = std::system(cmd.c_str());
if (!fs::exists(hash_path)) {
printout(ERROR, "FileLoader", "unable to run cmd " + cmd);
printout(ERROR, "FileLoader", "check command and retry");
printout(ERROR, "FileLoader", "hint: allow insecure connections with -k");
std::_Exit(EXIT_FAILURE);
}
}
// check if file already exists
if (fs::exists(file_path)) {
// file already exists
if (fs::is_symlink(file_path)) {
// file is symlink
if (fs::equivalent(hash_path, fs::read_symlink(file_path))) {
// link points to correct path
return;
} else {
// link points to incorrect path
if (fs::remove(file_path) == false) {
printout(ERROR, "FileLoader", "unable to remove symlink " + file_path.string());
printout(ERROR, "FileLoader", "check permissions or remove manually");
std::_Exit(EXIT_FAILURE);
}
}
} else {
// file exists but not symlink
printout(ERROR, "FileLoader", "will not remove actual file " + file_path.string());
printout(ERROR, "FileLoader", "check content, remove manually, and retry");
std::_Exit(EXIT_FAILURE);
}
}
// file_path now does not exist
// symlink file_path to hash_path
try {
// use new path from hash so file link is local
fs::create_symlink(fs::path(hash), file_path);
} catch (const fs::filesystem_error&) {
printout(ERROR, "FileLoader", "unable to link from " + file_path.string() + " to " + hash_path.string());
printout(ERROR, "FileLoader", "check permissions and retry");
std::_Exit(EXIT_FAILURE);
}
}
src/GaseousRICH_geo.cpp
View file @ 1008a184
......@@ -75,7 +75,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
envVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
// sensitive detector type
sens.setType("photoncounter");
sens.setType("tracker");
// @TODO: place a radiator
// build_radiator(desc, envVol, detElement.child(_Unicode(radiator)), snout_front);
......
src/HybridCalorimeter_geo.cpp
View file @ 1008a184
......@@ -40,7 +40,9 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete
auto air_material = desc.material("Air");
double ROut = desc.constantAsDouble("EcalEndcapN_rmax");
double RIn = desc.constantAsDouble("EcalEndcapN_rmin");
double RIn_el = desc.constantAsDouble("EcalEndcapN_rmin1");
double ionCutout_dphi = desc.constantAsDouble("EcalEndcapNIonCutout_dphi");
double RIn = desc.constantAsDouble("EcalEndcapN_rmin2");
double SizeZ = desc.constantAsDouble("EcalEndcapN_thickness");
double thickness = desc.constantAsDouble("EcalEndcapN_thickness");
double trans_radius = desc.constantAsDouble("EcalEndcapNCrystal_rmax");
......@@ -57,23 +59,35 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete
double Crystal_z0 = desc.constantAsDouble("CrystalModule_z0");
// RIn and ROut will define outer tube embedding the calorimeter
// centers_rin/out define the maximum radius of module centers
// centers_rmin/out define the maximum radius of module centers
// so that modules are not overlapping with mother tube volume
double hypotenuse = sqrt(0.5 * glass_distance * glass_distance);
double centers_rin = RIn + hypotenuse + 1*mm;
double centers_rout = ROut - hypotenuse - 1*mm;
const double Crystal_offset = -0.5*(Crystal_thickness - thickness);
const double glassHypotenuse = std::hypot(glass_distance, glass_distance)/2;
const double crystalHypotenuse = glassHypotenuse/2;
// Offset these values a bit so we don't miss edge-blocks
const double glassCenters_rmax = ROut - glassHypotenuse + 1 * mm;
const double crystalCenters_rmin = RIn + crystalHypotenuse - 1 * mm;
// Also limits of the inner crystal blocks fill
const double cutout_tan = tan(ionCutout_dphi/2);
const double cutout_rmin = RIn_el + crystalHypotenuse - 1 * mm;
// Offset to align the modules at the zmin of the endcap,
const double Crystal_offset = -0.5 * (Crystal_thickness - thickness);
const double Glass_offset = -0.5*(Glass_thickness - thickness);
// envelope
// Assembly assembly(detName);
Tube outer_tube(RIn, ROut, SizeZ / 2.0, 0., 360.0 * deg);
Volume ecal_vol("negative_ecal", outer_tube, air_material);
// consists of an glass tube of the full thickness, and a crystal inner tube
// for the crystal that allows us to get closet to the beampipe without
// overlaps, and a partial electron tube that allows us to get closer to the
// electron beampipe in the region where there is no ion beampipe
Tube glass_tube(min(RIn + glassHypotenuse*2, trans_radius), ROut, SizeZ / 2.0, 0., 360.0 * deg);
Tube crystal_tube(RIn, min(RIn + glassHypotenuse*2, trans_radius), Crystal_thickness/ 2.0, 0., 360.0 * deg);
Tube electron_tube(RIn_el, RIn, Crystal_thickness / 2., ionCutout_dphi / 2., 360.0 * deg - ionCutout_dphi / 2);
UnionSolid outer_envelope(glass_tube, crystal_tube, Position(0,0,Crystal_offset));
UnionSolid envelope(outer_envelope, electron_tube, Position(0,0,Crystal_offset));
Volume ecal_vol("negative_ecal", envelope, air_material);
ecal_vol.setVisAttributes(desc.visAttributes("HybridEcalOuterVis"));
// TODO why 1cm and not something else?
double Glass_OuterR = ROut - 1 * cm ;
double Glass_InnerR = trans_radius;
......@@ -91,9 +105,11 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete
// GLASS
double diameter = 2 * Glass_OuterR;
// How many towers do we have per row/columnt?
// Add a gap + diameter as if we have N towers, we have N-1 gaps;
int towersInRow = std::ceil((diameter + Glass_Gap) / (Glass_Width + Glass_Gap));
// Can we fit an even or odd amount of glass blocks within our rmax?
// This determines the transition points between crystal and glass as we need the
// outer crystal arangement to be in multiples of 2 (aligned with glass)
const int towersInRow = std::ceil((diameter + Glass_Gap) / (Glass_Width + Glass_Gap));
const bool align_even = (towersInRow % 2);
// Is it odd or even number of towersInRow
double leftTowerPos, topTowerPos;
......@@ -114,109 +130,39 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete
int moduleIndex = 0;
// fmt::print("\nCE EMCAL GLASS SQUARE START\n");
// fmt::print("Glass_thickness = {} cm;\n", Glass_thickness / cm);
// fmt::print("Glass_Width = {} cm;\n", Glass_Width / cm);
// fmt::print("Glass_Gap = {} cm;\n", Glass_Gap / cm);
// fmt::print("Glass_InnerR = {} cm;\n", Glass_InnerR / cm);
// fmt::print("Glass_OuterR = {} cm;\n", Glass_OuterR / cm);
// fmt::print("Glass_PosZ = {} cm;\n", glass_shift_z / cm);
// fmt::print("Towers in Row/Col = {} cm;\n", glass_shift_z / cm);
// fmt::print("Top left tower pos = {:<10} {:<10} cm;\n", -leftTowerPos / cm, topTowerPos / cm);
// fmt::print("#Towers info:\n");
// fmt::print("#{:<5} {:<6} {:<3} {:<3} {:>10} {:>10} {}\n", "idx", "code", "col", "row", "x", "y", "name");
// We first do a "dry run", not really placing modules,
// but figuring out the ID scheme, number of modules, etc.
int glassModuleCount = 0;
int crystalModuleCount = 0;
int firstCrystRow = 1000000; // The first row, where crystals are started
int firstCrystCol = 1000000; // The fist column, where crystals are started
for(int rowIndex=0; rowIndex < towersInRow; rowIndex++) {
for(int colIndex=0; colIndex < towersInRow; colIndex++) {
double glass_x = leftTowerPos + colIndex * glass_distance;
double glass_y = topTowerPos + rowIndex * glass_distance;
double r = sqrt(glass_x * glass_x + glass_y * glass_y);
if (r < centers_rout && r > centers_rin) {
// we are placing something
if(r<trans_radius) {
// 4 crystal modules will be placed
crystalModuleCount+=4;
// Finding the first col and row where crystals start
// is the same algorithm as finding a minimum in array
if(colIndex<firstCrystCol) {
firstCrystCol = colIndex;
}
if(rowIndex<firstCrystRow) {
firstCrystRow = rowIndex;
}
}
else
{
// glass module will be places
glassModuleCount++;
}
}
}
}
// fmt::print("#Towers info:\n");
// fmt::print("#{:<5} {:<6} {:<3} {:<3} {:>10} {:>10} {}\n", "idx", "code", "col", "row", "x", "y", "name");
int glass_module_index = 0;
int cryst_module_index = 0;
for(int rowIndex=0; rowIndex < towersInRow; rowIndex++) {
for(int colIndex=0; colIndex < towersInRow; colIndex++) {
double glass_x = leftTowerPos + colIndex * (Glass_Width + Glass_Gap);
double glass_y = topTowerPos + rowIndex * (Glass_Width + Glass_Gap);
double r = sqrt(glass_x * glass_x + glass_y * glass_y);
if (r < centers_rout && r > centers_rin) {
int code = 1000 * rowIndex + colIndex;
std::string name = fmt::format("ce_EMCAL_glass_phys_{}", code);
if(r<trans_radius) {
// first crystal module
double crystal_x = glass_x - crystal_distance / 2;
double crystal_y = glass_y - crystal_distance / 2;
auto placement = ecal_vol.placeVolume(crystal_module, Position(crystal_x, crystal_y, Crystal_z0 + Crystal_offset));
placement.addPhysVolID("sector", 1);
placement.addPhysVolID("module", cryst_module_index++);
// second crystal module
crystal_x = glass_x + crystal_distance / 2;
crystal_y = glass_y - crystal_distance / 2;
placement = ecal_vol.placeVolume(crystal_module, Position(crystal_x, crystal_y, Crystal_z0 + Crystal_offset));
placement.addPhysVolID("sector", 1);
placement.addPhysVolID("module", cryst_module_index++);
// third crystal module
crystal_x = glass_x - crystal_distance / 2;
crystal_y = glass_y + crystal_distance / 2;
placement = ecal_vol.placeVolume(crystal_module, Position(crystal_x, crystal_y, Crystal_z0 + Crystal_offset));
placement.addPhysVolID("sector", 1);
placement.addPhysVolID("module", cryst_module_index++);
// forth crystal module
crystal_x = glass_x + crystal_distance / 2;
crystal_y = glass_y + crystal_distance / 2;
placement = ecal_vol.placeVolume(crystal_module, Position(crystal_x, crystal_y, Crystal_z0 + Crystal_offset));
placement.addPhysVolID("sector", 1);
placement.addPhysVolID("module", cryst_module_index++);
}
else
{
// glass module
auto placement = ecal_vol.placeVolume(glass_module, Position(glass_x, glass_y, Glass_z0 + Glass_offset));
placement.addPhysVolID("sector", 2);
placement.addPhysVolID("module", glass_module_index++);
const double glass_x = leftTowerPos + colIndex * (Glass_Width + Glass_Gap);
const double glass_y = topTowerPos + rowIndex * (Glass_Width + Glass_Gap);
const double r = std::hypot(glass_x, glass_y);
// crystal if within the transition radius (as defined by the equivalent glass
// block)
if (r < trans_radius) {
for (const auto dx : {-1, 1}) {
for (const auto& dy : {-1, 1}) {
const double crystal_x = glass_x + dx * crystal_distance / 2;
const double crystal_y = glass_y + dy * crystal_distance / 2;
const double crystal_r = std::hypot(crystal_x, crystal_y);
// check if crystal in the main crystal ring?
const bool mainRing = (crystal_r > crystalCenters_rmin);
const bool innerRing = !mainRing && crystal_r > cutout_rmin;
const bool ionCutout = crystal_x > 0 && fabs(crystal_y / crystal_x) < cutout_tan;
if (mainRing || (innerRing && !ionCutout)) {
auto placement =
ecal_vol.placeVolume(crystal_module, Position(crystal_x, crystal_y, Crystal_z0 + Crystal_offset));
placement.addPhysVolID("sector", 1);
placement.addPhysVolID("module", cryst_module_index++);
}
}
}
// fmt::print(" {:<5} {:<6} {:<3} {:<3} {:>10.4f} {:>10.4f} {}\n", towerIndex, code, colIndex, rowIndex, x / cm, y / cm, name);
//glass_module_index++;
// Glass block if within the rmax
} else if (r < glassCenters_rmax) {
// glass module
auto placement = ecal_vol.placeVolume(glass_module, Position(glass_x, glass_y, Glass_z0 + Glass_offset));
placement.addPhysVolID("sector", 2);
placement.addPhysVolID("module", glass_module_index++);
}
}
}
......
src/MRich_geo.cpp
View file @ 1008a184
......@@ -29,7 +29,7 @@ static Ref_t createDetector(Detector& description, xml::Handle_t e, SensitiveDet
string det_name = x_det.nameStr();
DetElement sdet(det_name, x_det.id());
Assembly assembly(det_name);
sens.setType("photoncounter");
sens.setType("tracker");
OpticalSurfaceManager surfMgr = description.surfaceManager();
bool projective = getAttrOrDefault(x_det, _Unicode(projective), false);
......@@ -437,7 +437,7 @@ static Ref_t createDetector(Detector& description, xml::Handle_t e, SensitiveDet
// modVol.placeVolume(mVol, Position(0., 0., -0.1*mm));
//
// modVol.setVisAttributes(description.visAttributes(mods.visStr()));
// sens.setType("photoncounter");
// sens.setType("tracker");
// modVol.setSensitiveDetector(sens);
//
// // place modules in the sectors (disk)
......