RomanPot_geo.cpp
Go to the documentation of this file.
1 #include "DD4hep/DetFactoryHelper.h"
2 #include "DD4hep/Printout.h"
3 #include "DD4hep/Shapes.h"
4 #include "TMath.h"
5 #include "DDRec/Surface.h"
6 #include "DDRec/DetectorData.h"
7 #include "XML/Layering.h"
8 
9 using namespace std;
10 using namespace dd4hep;
11 using namespace dd4hep::rec;
12 using namespace ROOT::Math;
13 
42 static Ref_t build_detector(Detector& dtor, xml_h e, SensitiveDetector sens) {
43  xml_det_t x_det = e;
44  //Layering layering (e);
45  int det_id = x_det.id();
46  string det_name = x_det.nameStr();
47  xml_dim_t dim = x_det.dimensions();
48  double pixel_x = dim.x();
49  double pixel_y = dim.y();
50  double delta = dim.delta();
51  xml_dim_t frame_dim = x_det.child(_U(frame));
52  double frame_x = frame_dim.x();
53  double frame_y = frame_dim.y();
54  double frame_z = frame_dim.z();
55  string frame_vis = dd4hep::getAttrOrDefault(frame_dim, _Unicode(vis), "RPLayerVis");
56  xml_dim_t pos = x_det.position();
57  double zoffset = pos.z_offset();
58  double x_pos = dd4hep::getAttrOrDefault(pos, _Unicode(x),0.0);
59  double y_pos = dd4hep::getAttrOrDefault(pos, _Unicode(y),0.0);
60  double z_pos = dd4hep::getAttrOrDefault(pos, _Unicode(z),0.0);
61  bool rotated = pos.attr<bool>(_Unicode(rotation));
62  double pos_out = pos.attr<double>(_Unicode(vmax));
63  double curr_pos = pos.attr<double>(_Unicode(v));
64  Material vacuum = dtor.material("Vacuum");
65  Material aluminum = dtor.material("Aluminum");
66  Material frame_mat = dtor.material("Aluminum");
67 
68  DetElement sdet(det_name, det_id);
69  Assembly assembly(det_name + "_assembly");
70 
71  sens.setType("tracker");
72  string module_name = "RomanPot";
73 
74  string vis0 = dd4hep::getAttrOrDefault(x_det, _Unicode(vis), "BlueVis");
75 
76  sdet.setAttributes(dtor, assembly, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());
77 
78  double rp_detector_tube_length = frame_x + pos_out; // TODO: Shortest possible tube given the RP dimensions. Will probably change in the future.
79  double rp_detector_tube_thickness = 5.0 * dd4hep::mm;
80  double rp_detector_tube_radius = 0.55*sqrt(frame_z * frame_z + frame_y * frame_y); // TODO: Tightest possible fit. Will probably change in the future.
81  Tube rp_detector_vacuum_tube(0.0, rp_detector_tube_radius + rp_detector_tube_thickness, rp_detector_tube_length);
82  Tube rp_detector_vacuum_tube2(0.0, rp_detector_tube_radius, rp_detector_tube_length);
83 
84  // These are currently hardcoded. Will probably change in the future.
85  double rp_chamber_thickness = 10.0 * dd4hep::mm;
86  double rp_chamber_radius = 5.0 * dd4hep::cm;
87  double rp_chamber_length = rp_detector_tube_radius*1.5;
88 
89  // All the rotation and translation shenanigans related to RP unit orientation happen here.
90  // `rot1` and `rot2` are used for in construction of the tubes, `shift_pos1` and `shift_pos2`
91  // properly transform the RP detector positions.
92  auto rot1 = rotated ? Rotation3D(RotationX(0.5*M_PI)) : Rotation3D(RotationY(0.5*M_PI));
93  auto rot2 = rotated ? Rotation3D(RotationX(-0.5*M_PI)) : Rotation3D(RotationY(-0.5*M_PI));
94  auto shift_pos1 = rotated ? Transform3D(RotationZ(0.5 * M_PI) * Translation3D(curr_pos + 0.5 * frame_x, 0.0, 0.0))
95  : Transform3D(Translation3D(curr_pos + 0.5 * frame_x, 0, 0));
96  auto shift_pos2 = rotated ? Transform3D(RotationZ(-0.5 * M_PI) * Translation3D(curr_pos + 0.5 * frame_x, 0.0, 0.0))
97  : Transform3D(RotationZ(M_PI) * Translation3D(curr_pos + 0.5 * frame_x, 0, 0));
98  auto det_offset = Position(0.5*(pixel_x + delta - frame_x), 0, 0);
99 
100  // Construct the RP enclosure shell.
101  Tube rp_beam_pipe_tube(rp_chamber_radius, rp_chamber_radius + rp_chamber_thickness, rp_chamber_length/2.0);
102  Tube rp_beam_vacuum_tube(0.0, rp_chamber_radius + rp_chamber_thickness, rp_chamber_length);
103  Tube rp_beam_vacuum_tube2(0.0, rp_chamber_radius, rp_chamber_length);
104  UnionSolid rp_chamber_tee_outer1(rp_beam_vacuum_tube, rp_detector_vacuum_tube, rot1);
105  UnionSolid rp_chamber_tee_outer(rp_chamber_tee_outer1, rp_detector_vacuum_tube, rot2);
106  UnionSolid rp_chamber_tee_inner1(rp_beam_vacuum_tube2, rp_detector_vacuum_tube2, rot1);
107  UnionSolid rp_chamber_tee_inner(rp_chamber_tee_inner1, rp_detector_vacuum_tube2, rot2);
108  //SubtractionSolid shell(rp_chamber_tee_outer,rp_chamber_tee_inner);
109  Volume rp_chamber_vol("rp_chamber_walls_vol", rp_chamber_tee_outer, aluminum);
110  Volume rp_vacuum_vol("rp_chamber_vacuum_vol", rp_chamber_tee_inner, vacuum);
111 
112  auto chamber_pv = assembly.placeVolume(rp_chamber_vol);
113  auto vacuum_pv = rp_chamber_vol.placeVolume(rp_vacuum_vol);
114  vacuum_pv.addPhysVolID("element", 1);
115  rp_chamber_vol.setVisAttributes(vis0.c_str());
116  rp_vacuum_vol.setVisAttributes(vis0.c_str());
117  //rp_vacuum_vol.setVisAttributes(dtor.invisible());
118 
119  DetElement vacuum_DE(sdet, "vacuum_DE", 1);
120  // Make the enclosure box for the detector frame
121  Box frame_box_full(0.5 * frame_x, 0.5 * frame_y, 0.5 * frame_z);
122  Box det_cutoff(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), 0.5 * frame_z);
123  // Make the cutout where whe detector will sit
124  // The frame with a hole for the detector
125  SubtractionSolid frame_box(frame_box_full, det_cutoff, det_offset);
126  // Here we assum there is only one unique layer which is repeated.
127  // Layer
128  xml_comp_t x_layer = x_det.child(_U(layer));
129  const int repeat = x_layer.repeat();
130  double layer_thickness = 0.0;
131  for (xml_coll_t si(x_layer, _U(slice)); si; ++si) {
132  xml_comp_t x_slice = si;
133  layer_thickness += x_slice.thickness();
134  }
135  std::cout << " layer thickness = " << layer_thickness << " \n";
136  double total_thickness = repeat*layer_thickness;
137  // --- d1 ---
138  //
139  Volume frame_box_full1_vol("frame_box_full1", frame_box_full, vacuum);
140  PlacedVolume frame_box_full1_pv = rp_vacuum_vol.placeVolume(frame_box_full1_vol, shift_pos1);
141  DetElement frame1_DE(vacuum_DE, "frame_DE1", 1);
142  frame_box_full1_pv.addPhysVolID("frame", 1);
143 
144  frame1_DE.setAttributes(dtor, frame_box_full1_vol, x_det.regionStr(), x_det.limitsStr(), frame_vis);
145 
146  Volume frame1_vol("xsensor_frame1", frame_box, frame_mat);
147  PlacedVolume frame1_pv = frame_box_full1_vol.placeVolume(frame1_vol);
148 
149  // Loop over layers
150  double layer_pos_offset = -0.5 * total_thickness;
151  int layer_num = 1;
152  // How many times does the layer repeat. Defined in the compact description.
153  for (int j = 0; j < repeat; j++) {
154  string layer_name = _toString(layer_num, "layer1_%d");
155  // double layer_thickness =
156  // layering.layer(layer_num - 1)->thickness(); // Get the thickness of the current layer.
157  Position layer_pos = Position(0, 0, layer_pos_offset + 0.5 * layer_thickness);
158  Box layer_box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), layer_thickness / 2.0);
159  Volume layer_vol(layer_name, layer_box, vacuum);
160  DetElement layer(frame1_DE, layer_name, layer_num);
161  layer.setAttributes(dtor, layer_vol, x_det.regionStr(), x_det.limitsStr(), frame_vis);
162 
163  // Loop over slices within the layer
164  double slice_pos_offset = -0.5 * layer_thickness;
165  int slice_num = 1;
166  for (xml_coll_t si(x_layer, _U(slice)); si; ++si) {
167  xml_comp_t x_slice = si;
168  string slice_name = _toString(slice_num, "slice1_%d");
169  double slice_thickness = x_slice.thickness();
170  Position slice_pos = Position(0, 0, slice_pos_offset + 0.5 * slice_thickness);
171  Box slice_box = x_slice.isSensitive()
172  ? Box(0.5 * pixel_x, 0.5 * pixel_y, slice_thickness / 2.0)
173  : Box(0.5 * (pixel_x), 0.5 * (pixel_y), slice_thickness / 2.0);
174  Volume slice_vol(slice_name, slice_box, dtor.material(x_slice.materialStr()));
175  string slice_vis = dd4hep::getAttrOrDefault(x_slice, _Unicode(vis), "BlueVis");
176  slice_vol.setVisAttributes(slice_vis.c_str());
177  DetElement slice(layer, slice_name, slice_num);
178  if (x_slice.isSensitive()) {
179  Vector3D u(1., 0., 0.);
180  Vector3D v(0., 1., 0.);
181  Vector3D n(0., 0., 1.);
182  Vector3D o(0., 0., 0.);
183  double inner_thickness = slice_pos_offset + slice_thickness / 2.0;
184  double outer_thickness = layer_thickness - inner_thickness;
185  SurfaceType type(SurfaceType::Sensitive);
186  VolPlane surf(slice_vol, type, inner_thickness, outer_thickness, u, v, n, o);
187  slice_vol.setSensitiveDetector(sens);
188  sens.setType("tracker");
189  }
190  slice.setAttributes(dtor, slice_vol, x_det.regionStr(), x_det.limitsStr(), slice_vis);
191  // Place the slice.
192  PlacedVolume slice_pv = layer_vol.placeVolume(slice_vol, slice_pos);
193  slice_pv.addPhysVolID("slice", slice_num );
194  slice.setPlacement(slice_pv);
195  slice_pos_offset += slice_thickness; // Move the position offset for the next slice.
196  ++slice_num;
197  }
198  // Place the layer.
199  PlacedVolume layer_pv = frame_box_full1_vol.placeVolume(layer_vol, layer_pos + det_offset);
200  layer_pv.addPhysVolID("layer", layer_num );
201  layer.setPlacement(layer_pv);
202  layer_pos_offset += layer_thickness;
203  ++layer_num;
204  }
205 
206  // --- d2 ---
207  //Volume frame2_vol("xsensor_frame2", frame_box, frame_mat);
208  //PlacedVolume frame2_pv = rp_vacuum_vol.placeVolume(frame2_vol, shift_pos2);
209  //DetElement frame2_DE(vacuum_DE, "frame_DE2", 2);
210  //frame2_pv.addPhysVolID("frame", 2);
211  //frame2_DE.setAttributes(dtor, frame2_vol, x_det.regionStr(), x_det.limitsStr(), "RedVis");
212 
213  Volume frame_box_full2_vol("frame_box_full2", frame_box_full, vacuum);
214  PlacedVolume frame_box_full2_pv = rp_vacuum_vol.placeVolume(frame_box_full2_vol, shift_pos2);
215  DetElement frame2_DE(vacuum_DE, "frame_DE2", 2);
216  frame_box_full2_pv.addPhysVolID("frame", 2);
217  frame2_DE.setAttributes(dtor, frame_box_full2_vol, x_det.regionStr(), x_det.limitsStr(), frame_vis);
218 
219  Volume frame2_vol("xsensor_frame2", frame_box, frame_mat);
220  PlacedVolume frame2_pv = frame_box_full2_vol.placeVolume(frame2_vol);
221 
222  // Loop over layers
223  layer_pos_offset = -0.5 * total_thickness;
224  layer_num = 1;
225  for (int j = 0; j < repeat; j++) {
226  string layer_name = _toString(layer_num, "layer2_%d");
227  // double layer_thickness = layering.layer(layer_num - 1)->thickness(); // Get the thickness
228  // of the current layer.
229  Position layer_pos = Position(0, 0, layer_pos_offset + 0.5 * layer_thickness);
230  Box layer_box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), layer_thickness / 2.0);
231  Volume layer_vol(layer_name, layer_box, vacuum);
232  DetElement layer(frame2_DE, layer_name, layer_num);
233  layer.setAttributes(dtor, layer_vol, x_det.regionStr(), x_det.limitsStr(),
234  frame_vis);
235 
236  // Loop over slices within the layer
237  double slice_pos_offset = -0.5 * layer_thickness;
238  int slice_num = 1;
239  for (xml_coll_t si(x_layer, _U(slice)); si; ++si) {
240  xml_comp_t x_slice = si;
241  string slice_name = _toString(slice_num, "slice2_%d");
242  double slice_thickness = x_slice.thickness();
243  Position slice_pos = Position(0, 0, slice_pos_offset + 0.5 * slice_thickness);
244  Box slice_box = x_slice.isSensitive()
245  ? Box(0.5 * pixel_x, 0.5 * pixel_y, slice_thickness / 2.0)
246  : Box(0.5 * pixel_x, 0.5 * pixel_y, slice_thickness / 2.0);
247  Volume slice_vol(slice_name, slice_box, dtor.material(x_slice.materialStr()));
248  string slice_vis = dd4hep::getAttrOrDefault(x_slice, _Unicode(vis), "BlueVis");
249  slice_vol.setVisAttributes(slice_vis.c_str());
250  DetElement slice(layer, slice_name, slice_num);
251  if (x_slice.isSensitive()) {
252  Vector3D u(1., 0., 0.);
253  Vector3D v(0., 1., 0.);
254  Vector3D n(0., 0., 1.);
255  Vector3D o(0., 0., 0.);
256  double inner_thickness = slice_pos_offset + slice_thickness / 2.0;
257  double outer_thickness = layer_thickness - inner_thickness;
258  SurfaceType type(SurfaceType::Sensitive);
259  VolPlane surf(slice_vol, type, inner_thickness, outer_thickness, u, v, n, o);
260  slice_vol.setSensitiveDetector(sens);
261  sens.setType("tracker");
262  }
263  slice.setAttributes(dtor, slice_vol, x_det.regionStr(), x_det.limitsStr(), slice_vis);
264  // Place the slice.
265  PlacedVolume slice_pv = layer_vol.placeVolume(slice_vol, slice_pos);
266  slice_pv.addPhysVolID("slice", slice_num );
267  slice.setPlacement(slice_pv);
268  slice_pos_offset += slice_thickness; // Move the position offset for the next slice.
269  ++slice_num;
270  }
271  // Place the layer.
272  PlacedVolume layer_pv = frame_box_full2_vol.placeVolume(layer_vol, layer_pos + det_offset);
273  layer_pv.addPhysVolID("layer", layer_num );
274  layer.setPlacement(layer_pv);
275  layer_pos_offset += layer_thickness;
276  ++layer_num;
277  }
278 
279  auto assembly_pv = dtor.pickMotherVolume(sdet).placeVolume(assembly, Position(x_pos, y_pos, z_pos + zoffset));
280  assembly_pv.addPhysVolID("system", det_id);
281  sdet.setPlacement(assembly_pv);
282 
283  assembly->GetShape()->ComputeBBox();
284  return sdet;
285 }
287 DECLARE_DETELEMENT(RomanPot, build_detector)
static Ref_t build_detector(Detector &dtor, xml_h e, SensitiveDetector sens)
Detector
Definition: DDG4.py:69
Namespace for the AIDA detector description toolkit.