Belle II Software development
GeoFarBeamLineCreator.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <ir/geometry/GeoFarBeamLineCreator.h>
10
11#include <geometry/Materials.h>
12#include <geometry/CreatorFactory.h>
13#include <geometry/utilities.h>
14#include <framework/gearbox/Unit.h>
15#include <ir/simulation/SensitiveDetector.h>
16#include <simulation/background/BkgSensitiveDetector.h>
17
18#include <cmath>
19#include <boost/format.hpp>
20#include <boost/foreach.hpp>
21#include <boost/algorithm/string.hpp>
22
23#include <G4LogicalVolume.hh>
24#include <G4PVPlacement.hh>
25
26//Shapes
27#include <G4Box.hh>
28#include <G4Tubs.hh>
29#include <G4Torus.hh>
30#include <G4Polycone.hh>
31#include <G4Trd.hh>
32#include <G4IntersectionSolid.hh>
33#include <G4SubtractionSolid.hh>
34#include <G4UserLimits.hh>
35
36using namespace std;
37using namespace boost;
38
39namespace Belle2 {
46 using namespace geometry;
47
48 namespace ir {
49
50 //-----------------------------------------------------------------
51 // Register the Creator
52 //-----------------------------------------------------------------
53
54 geometry::CreatorFactory<GeoFarBeamLineCreator> GeoFarBeamLineFactory("FarBeamLineCreator");
55
56 //-----------------------------------------------------------------
57 // Implementation
58 //-----------------------------------------------------------------
59
61 {
63 }
64
66 {
67 delete m_sensitive;
68 }
69
70 void GeoFarBeamLineCreator::createGeometry(G4LogicalVolume& topVolume, GeometryTypes)
71 {
72
73 const int N = 2;
74
75 double stepMax = 5.0 * Unit::mm;
76 int flag_limitStep = int(m_config.getParameter("LimitStepLength"));
77
78
79 //double unitFactor = 10.0;
80 const double unitFactor = Unit::cm / Unit::mm;
81
82 map<string, FarBeamLineElement> elements;
83
84
85 //--------------
86 //- limits
87
88 //--------------
89 //- TubeR
90
92
93 //get parameters from .xml file
94 std::string prep = "TubeR.";
95
96 int TubeR_N = int(m_config.getParameter(prep + "N"));
97
98 std::vector<double> TubeR_Z(TubeR_N);
99 std::vector<double> TubeR_R(TubeR_N);
100 std::vector<double> TubeR_r(TubeR_N);
101
102 for (int i = 0; i < TubeR_N; ++i) {
103 ostringstream ossZID;
104 ossZID << "Z" << i;
105
106 ostringstream ossRID;
107 ossRID << "R" << i;
108
109 ostringstream ossrID;
110 ossrID << "r" << i;
111
112 TubeR_Z[i] = m_config.getParameter(prep + ossZID.str()) * unitFactor;
113 TubeR_R[i] = m_config.getParameter(prep + ossRID.str()) * unitFactor;
114 TubeR_r[i] = m_config.getParameter(prep + ossrID.str(), 0.0) * unitFactor;
115 }
116
117 tubeR.transform = G4Translate3D(0.0, 0.0, 0.0);
118
119 //define geometry
120 tubeR.geo = new G4Polycone("geo_TubeR_name", 0.0, 2 * M_PI, TubeR_N, &(TubeR_Z[0]), &(TubeR_r[0]), &(TubeR_R[0]));
121
122 tubeR.logi = NULL;
123
124 elements["TubeR"] = tubeR;
125
126 //--------------
127 //- TubeL
128
129 FarBeamLineElement tubeL;
130
131 //get parameters from .xml file
132 prep = "TubeL.";
133
134 int TubeL_N = int(m_config.getParameter(prep + "N"));
135
136 std::vector<double> TubeL_Z(TubeL_N);
137 std::vector<double> TubeL_R(TubeL_N);
138 std::vector<double> TubeL_r(TubeL_N);
139
140 for (int i = 0; i < TubeL_N; ++i) {
141 ostringstream ossZID;
142 ossZID << "Z" << i;
143
144 ostringstream ossRID;
145 ossRID << "R" << i;
146
147 ostringstream ossrID;
148 ossrID << "r" << i;
149
150 TubeL_Z[i] = m_config.getParameter(prep + ossZID.str()) * unitFactor;
151 TubeL_R[i] = m_config.getParameter(prep + ossRID.str()) * unitFactor;
152 TubeL_r[i] = m_config.getParameter(prep + ossrID.str(), 0.0) * unitFactor;
153 }
154
155 tubeL.transform = G4Translate3D(0.0, 0.0, 0.0);
156
157 //define geometry
158 tubeL.geo = new G4Polycone("geo_TubeL_name", 0.0, 2 * M_PI, TubeL_N, &(TubeL_Z[0]), &(TubeL_r[0]), &(TubeL_R[0]));
159
160 tubeL.logi = NULL;
161
162 elements["TubeL"] = tubeL;
163
164 std::vector<double> zero_r(N, 0.);
165
166 std::vector<std::string> straightSections;
167 boost::split(straightSections, m_config.getParameterStr("Straight"), boost::is_any_of(" "));
168 for (const auto& name : straightSections) {
169 //--------------
170 //- Create straight element
171
172 prep = name + ".";
173 string type = m_config.getParameterStr(prep + "type");
174
175
176 FarBeamLineElement polycone;
177
178 std::vector<double> Polycone_Z(N);
179 std::vector<double> Polycone_R(N);
180 std::vector<double> Polycone_r(N);
181 Polycone_Z[0] = 0;
182 Polycone_Z[1] = m_config.getParameter(prep + "L") * unitFactor;
183 Polycone_R[0] = m_config.getParameter(prep + "R") * unitFactor;
184 Polycone_R[1] = m_config.getParameter(prep + "R") * unitFactor;
185 Polycone_r[0] = m_config.getParameter(prep + "r") * unitFactor;
186 Polycone_r[1] = m_config.getParameter(prep + "r") * unitFactor;
187
188 double Polycone_X0 = m_config.getParameter(prep + "X0") * unitFactor;
189 double Polycone_Y0 = m_config.getParameter(prep + "Y0", 0) * unitFactor;
190 double Polycone_Z0 = m_config.getParameter(prep + "Z0") * unitFactor;
191 double Polycone_PHI = m_config.getParameter(prep + "PHI");
192 double Polycone_PHIYZ = m_config.getParameter(prep + "PHIYZ", 0);
193
194 polycone.transform = G4Translate3D(Polycone_X0, Polycone_Y0, Polycone_Z0);
195 polycone.transform = polycone.transform * G4RotateY3D(Polycone_PHI / Unit::rad);
196 if (Polycone_PHIYZ != 0)
197 polycone.transform = polycone.transform * G4RotateX3D(Polycone_PHIYZ / Unit::rad);
198
199 //define geometry
200 string subtract = m_config.getParameterStr(prep + "Subtract", "");
201 string intersect = m_config.getParameterStr(prep + "Intersect", "");
202
203 string geo_polyconexx_name = "geo_" + name + "xx_name";
204 string geo_polyconex_name = "geo_" + name + "x_name";
205 string geo_polycone_name = "geo_" + name + "_name";
206
207 G4VSolid* geo_polyconexx(NULL), *geo_polycone(NULL);
208
209 if (subtract != "" || intersect != "")
210 if (type == "pipe") // for pipes inner space will be created as vacuum
211 geo_polyconexx = new G4Polycone(geo_polyconexx_name, 0.0, 2 * M_PI, N, &(Polycone_Z[0]), &(zero_r[0]), &(Polycone_R[0]));
212 else
213 geo_polyconexx = new G4Polycone(geo_polyconexx_name, 0.0, 2 * M_PI, N, &(Polycone_Z[0]), &(Polycone_r[0]), &(Polycone_R[0]));
214 else if (type == "pipe") // for pipes inner space will be created as vacuum
215 geo_polycone = new G4Polycone(geo_polycone_name, 0.0, 2 * M_PI, N, &(Polycone_Z[0]), &(zero_r[0]), &(Polycone_R[0]));
216 else
217 geo_polycone = new G4Polycone(geo_polycone_name, 0.0, 2 * M_PI, N, &(Polycone_Z[0]), &(Polycone_r[0]), &(Polycone_R[0]));
218
219
220 if (subtract != "" && intersect != "") {
221 G4VSolid* geo_polyconex = new G4SubtractionSolid(geo_polyconex_name, geo_polyconexx, elements[subtract].geo,
222 polycone.transform.inverse()*elements[subtract].transform);
223 geo_polycone = new G4IntersectionSolid(geo_polycone_name, geo_polyconex, elements[intersect].geo,
224 polycone.transform.inverse()*elements[intersect].transform);
225 } else if (subtract != "")
226 geo_polycone = new G4SubtractionSolid(geo_polycone_name, geo_polyconexx, elements[subtract].geo,
227 polycone.transform.inverse()*elements[subtract].transform);
228 else if (intersect != "")
229 geo_polycone = new G4IntersectionSolid(geo_polycone_name, geo_polyconexx, elements[intersect].geo,
230 polycone.transform.inverse()*elements[intersect].transform);
231
232 polycone.geo = geo_polycone;
233
234 // define logical volume
235 string strMat_polycone = m_config.getParameterStr(prep + "Material");
236 G4Material* mat_polycone = Materials::get(strMat_polycone);
237 string logi_polycone_name = "logi_" + name + "_name";
238 G4LogicalVolume* logi_polycone = new G4LogicalVolume(polycone.geo, mat_polycone, logi_polycone_name);
239 setColor(*logi_polycone, "#CC0000");
240 setVisibility(*logi_polycone, false);
241
242 polycone.logi = logi_polycone;
243
244 //put volume
245 string phys_polycone_name = "phys_" + name + "_name";
246 new G4PVPlacement(polycone.transform, logi_polycone, phys_polycone_name, &topVolume, false, 0);
247
248 elements[name] = polycone;
249
250 double sum = 0.0;
251 for (int i = 0; i < N; ++i)
252 sum += Polycone_r[i]; // check that there is a space inside a pipe
253 if (type == "pipe" && sum != 0) { // add vacuum inside a pipe
254
255 FarBeamLineElement vacuum;
256
257 string nameVac = name + "Vac";
258
259 //define geometry
260 string geo_vacuumxx_name = "geo_" + nameVac + "xx_name";
261 string geo_vacuum_name = "geo_" + nameVac + "_name";
262
263 G4VSolid* geo_vacuumxx, *geo_vacuum;
264
265 geo_vacuumxx = new G4Polycone(geo_vacuumxx_name, 0.0, 2 * M_PI, N, &(Polycone_Z[0]), &(zero_r[0]), &(Polycone_r[0]));
266 geo_vacuum = new G4IntersectionSolid(geo_vacuumxx_name, geo_vacuumxx, geo_polycone);
267
268 vacuum.geo = geo_vacuum;
269 vacuum.transform = polycone.transform;
270
271 // define logical volume
272 G4Material* mat_vacuum = Materials::get("Vacuum");
273 string logi_vacuum_name = "logi_" + nameVac + "_name";
274 G4LogicalVolume* logi_vacuum = new G4LogicalVolume(vacuum.geo, mat_vacuum, logi_vacuum_name);
275 if (flag_limitStep) logi_vacuum->SetUserLimits(new G4UserLimits(stepMax));
276 setColor(*logi_vacuum, "#000000");
277 setVisibility(*logi_vacuum, false);
278
279 vacuum.logi = logi_vacuum;
280
281 //put volume
282 string phys_vacuum_name = "phys_" + nameVac + "_name";
283 //new G4PVPlacement(vacuum.transform, logi_vacuum, phys_vacuum_name, &topVolume, false, 0);
284 new G4PVPlacement(0, G4ThreeVector(0, 0, 0), logi_vacuum, phys_vacuum_name, logi_polycone, false, 0);
285
286 elements[nameVac] = vacuum;
287 }
288 }
289
290
291 std::vector<std::string> bendingSections;
292 boost::split(bendingSections, m_config.getParameterStr("Bending"), boost::is_any_of(" "));
293 for (const auto& name : bendingSections) {
294 //--------------
295 //- Create torus element
296
297 prep = name + ".";
298 string type = m_config.getParameterStr(prep + "type");
299
300 FarBeamLineElement torus;
301
302 double torus_r = m_config.getParameter(prep + "r") * unitFactor;
303 double torus_R = m_config.getParameter(prep + "R") * unitFactor;
304 double torus_RT = m_config.getParameter(prep + "RT") * unitFactor;
305 double torus_X0 = m_config.getParameter(prep + "X0") * unitFactor;
306 double torus_Z0 = m_config.getParameter(prep + "Z0") * unitFactor;
307 double torus_SPHI = m_config.getParameter(prep + "SPHI");
308 double torus_DPHI = m_config.getParameter(prep + "DPHI");
309
310 torus.transform = G4Translate3D(torus_X0, 0.0, torus_Z0);
311 torus.transform = torus.transform * G4RotateX3D(M_PI / 2 / Unit::rad);
312
313 //define geometry
314 string subtract = m_config.getParameterStr(prep + "Subtract", "");
315 string intersect = m_config.getParameterStr(prep + "Intersect", "");
316
317 string geo_torusxx_name = "geo_" + name + "xx_name";
318 string geo_torusx_name = "geo_" + name + "x_name";
319 string geo_torus_name = "geo_" + name + "_name";
320
321 G4VSolid* geo_torus(NULL);
322
323 if (subtract != "" || intersect != "") {
324 G4VSolid* geo_torusxx(NULL);
325 if (type == "pipe") // for pipes inner space will be created as vacuum
326 geo_torusxx = new G4Torus(geo_torusxx_name, 0, torus_R, torus_RT, torus_SPHI, torus_DPHI);
327 else
328 geo_torusxx = new G4Torus(geo_torusxx_name, torus_r, torus_R, torus_RT, torus_SPHI, torus_DPHI);
329 if (subtract != "" && intersect != "") {
330 G4VSolid* geo_torusx = new G4SubtractionSolid(geo_torusx_name, geo_torusxx, elements[subtract].geo,
331 torus.transform.inverse()*elements[subtract].transform);
332 geo_torus = new G4IntersectionSolid(geo_torus_name, geo_torusx, elements[intersect].geo,
333 torus.transform.inverse()*elements[intersect].transform);
334 } else if (subtract != "")
335 geo_torus = new G4SubtractionSolid(geo_torus_name, geo_torusxx, elements[subtract].geo,
336 torus.transform.inverse()*elements[subtract].transform);
337 else
338 geo_torus = new G4IntersectionSolid(geo_torus_name, geo_torusxx, elements[intersect].geo,
339 torus.transform.inverse()*elements[intersect].transform);
340 } else if (type == "pipe") // for pipes inner space will be created as vacuum
341 geo_torus = new G4Torus(geo_torus_name, 0, torus_R, torus_RT, torus_SPHI, torus_DPHI);
342 else
343 geo_torus = new G4Torus(geo_torus_name, torus_r, torus_R, torus_RT, torus_SPHI, torus_DPHI);
344
345 torus.geo = geo_torus;
346
347 // define logical volume
348 string strMat_torus = m_config.getParameterStr(prep + "Material");
349 G4Material* mat_torus = Materials::get(strMat_torus);
350 string logi_torus_name = "logi_" + name + "_name";
351 G4LogicalVolume* logi_torus = new G4LogicalVolume(torus.geo, mat_torus, logi_torus_name);
352 setColor(*logi_torus, "#CC0000");
353 setVisibility(*logi_torus, false);
354
355 torus.logi = logi_torus;
356
357 //put volume
358 string phys_torus_name = "phys_" + name + "_name";
359 new G4PVPlacement(torus.transform, logi_torus, phys_torus_name, &topVolume, false, 0);
360
361 elements[name] = torus;
362
363 if (type == "pipe" && torus_r != 0) { // add vacuum inside a pipe
364
365 FarBeamLineElement vacuum;
366
367 string nameVac = name + "Vac";
368
369 //define geometry
370 string geo_vacuumxx_name = "geo_" + nameVac + "xx_name";
371 string geo_vacuum_name = "geo_" + nameVac + "_name";
372
373 G4VSolid* geo_vacuumxx, *geo_vacuum;
374
375 geo_vacuumxx = new G4Torus(geo_vacuumxx_name, 0.0, torus_r, torus_RT, torus_SPHI, torus_DPHI);
376 geo_vacuum = new G4IntersectionSolid(geo_vacuum_name, geo_vacuumxx, geo_torus);
377
378 vacuum.geo = geo_vacuum;
379 vacuum.transform = torus.transform;
380
381 // define logical volume
382 G4Material* mat_vacuum = Materials::get("Vacuum");
383 string logi_vacuum_name = "logi_" + nameVac + "_name";
384 G4LogicalVolume* logi_vacuum = new G4LogicalVolume(vacuum.geo, mat_vacuum, logi_vacuum_name);
385 if (flag_limitStep) logi_vacuum->SetUserLimits(new G4UserLimits(stepMax));
386 setColor(*logi_vacuum, "#000000");
387 setVisibility(*logi_vacuum, false);
388
389 vacuum.logi = logi_vacuum;
390
391 //put volume
392 string phys_vacuum_name = "phys_" + nameVac + "_name";
393 //new G4PVPlacement(vacuum.transform, logi_vacuum, phys_vacuum_name, &topVolume, false, 0);
394 new G4PVPlacement(0, G4ThreeVector(0, 0, 0), logi_vacuum, phys_vacuum_name, logi_torus, false, 0);
395
396 elements[nameVac] = vacuum;
397 }
398 }
399
400
401 //--------------------------------------------------------------------------------------------
402 //- Gate shields, end-of-tunnel concrete shields, polyethylene shields, collimator shields
403
404 std::vector<std::string> shields;
405 boost::split(shields, m_config.getParameterStr("Shield"), boost::is_any_of(" "));
406 for (const auto& name : shields) {
407 prep = name + ".";
408
409 //- Shield made as box with optional subtracted box-shaped inner space (hole)
410
411 double shield_W = m_config.getParameter(prep + "W") * unitFactor;
412 double shield_H = m_config.getParameter(prep + "H") * unitFactor;
413 double shield_L = m_config.getParameter(prep + "L") * unitFactor;
414 double shield_X0 = m_config.getParameter(prep + "X0") * unitFactor;
415 double shield_Y0 = m_config.getParameter(prep + "Y0") * unitFactor;
416 double shield_Z0 = m_config.getParameter(prep + "Z0") * unitFactor;
417
418 double shield_hole_W = m_config.getParameter(prep + "holeW", 0) * unitFactor;
419 double shield_hole_H = m_config.getParameter(prep + "holeH", 0) * unitFactor;
420 double shield_hole_L = m_config.getParameter(prep + "holeL", 0) * unitFactor;
421 double shield_hole_dX = m_config.getParameter(prep + "holeDX", 0) * unitFactor;
422 double shield_hole_dY = m_config.getParameter(prep + "holeDY", 0) * unitFactor;
423 double shield_hole_dZ = m_config.getParameter(prep + "holeDZ", 0) * unitFactor;
424
425 double shield_PHI = m_config.getParameter(prep + "PHI");
426
427 // storable element
428 FarBeamLineElement shield;
429
430 shield.transform = G4Translate3D(shield_X0, shield_Y0, shield_Z0);
431 shield.transform = shield.transform * G4RotateY3D(shield_PHI / Unit::rad);
432
433 G4Transform3D transform_shield_hole = G4Translate3D(shield_hole_dX, shield_hole_dY, shield_hole_dZ);
434
435 //define geometry
436 string geo_shieldx_name = "geo_" + name + "x_name";
437 string geo_shield_hole_name = "geo_" + name + "_hole_name";
438 string geo_shield_name = "geo_" + name + "_name";
439
440 if (shield_hole_W == 0 || shield_hole_H == 0 || shield_hole_L == 0) {
441 G4Box* geo_shield = new G4Box(geo_shield_name, shield_W / 2.0, shield_H / 2.0, shield_L / 2.0);
442
443 shield.geo = geo_shield;
444 } else {
445 G4Box* geo_shieldx = new G4Box(geo_shieldx_name, shield_W / 2.0, shield_H / 2.0, shield_L / 2.0);
446 G4Box* geo_shield_hole = new G4Box(geo_shield_hole_name, shield_hole_W / 2.0, shield_hole_H / 2.0, shield_hole_L / 2.0);
447 G4SubtractionSolid* geo_shield = new G4SubtractionSolid(geo_shield_name, geo_shieldx, geo_shield_hole,
448 transform_shield_hole);
449
450 shield.geo = geo_shield;
451 }
452
453 string strMat_shield = m_config.getParameterStr(prep + "Material");
454 G4Material* mat_shield = Materials::get(strMat_shield);
455
456 string logi_shield_name = "logi_" + name + "_name";
457 G4LogicalVolume* logi_shield = new G4LogicalVolume(shield.geo, mat_shield, logi_shield_name);
458
459 shield.logi = logi_shield;
460
461 //put volume
462 setColor(*logi_shield, "#0000CC");
463 //setVisibility(*logi_shield, false);
464 string phys_shield_name = "phys_" + name + "_name";
465 new G4PVPlacement(shield.transform, shield.logi, phys_shield_name, &topVolume, false, 0);
466
467 elements[name] = shield;
468 }
469
470
471 //--------------
472 //- Tube (virtual tube for radiation level study)
473
474 //define geometry
475 G4Tubs* geo_Tube = new G4Tubs("geo_Tube_name", 3995 * CLHEP::mm, 4000 * CLHEP::mm, 29 * CLHEP::m, 0. * CLHEP::deg, 360.*CLHEP::deg);
476 G4Material* mat_Tube = Materials::get("G4_Si");
477 G4LogicalVolume* logi_Tube = new G4LogicalVolume(geo_Tube, mat_Tube, "logi_Tube_name");
478
479 //put volume
480 setColor(*logi_Tube, "#CC0000");
481 setVisibility(*logi_Tube, false);
482 bool radiation_study = false;
483 // cppcheck-suppress knownConditionTrueFalse
484 if (radiation_study && elements.count("GateShieldL")) {
485 new G4PVPlacement(elements["GateShieldL"].transform, logi_Tube, "phys_Tube_name", &topVolume, false, 0);
486 }
487
488
489 //---------------------------
490 // for dose simulation
491 //---------------------------
492
493 // cppcheck-suppress knownConditionTrueFalse
494 if (radiation_study) {
495 //neutron shield (poly)
496 if (elements.count("PolyShieldL"))
497 elements["PolyShieldL"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1001));
498 if (elements.count("PolyShieldR"))
499 elements["PolyShieldR"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1002));
500
501 //additional neutron shield (concrete)
502 if (elements.count("ConcreteShieldL"))
503 elements["ConcreteShieldL"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1003));
504 if (elements.count("ConcreteShieldR"))
505 elements["ConcreteShieldR"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1004));
506
507 //gate shield (concrete)
508 if (elements.count("GateShieldL"))
509 elements["GateShieldL"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1005));
510 if (elements.count("GateShieldR"))
511 elements["GateShieldR"].logi->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1006));
512
513 //virtual material outside gate-shield
514 logi_Tube->SetSensitiveDetector(new BkgSensitiveDetector("IR", 1007));
515 }
516
517
518 //------------------
519 //- Collimators
520
521 std::vector<std::string> collimators;
522 boost::split(collimators, m_config.getParameterStr("Collimator"), boost::is_any_of(" "));
523 for (const auto& name : collimators) {
524 //- Collimators consist of two independent jaws (trapezoids), identical in shape, positioned opposite to each other
525 //- Each jaw consists of copper body and high Z head
526
527 prep = name + ".";
528
529 string type = m_config.getParameterStr(prep + "type");
530 string motherVolume = m_config.getParameterStr(prep + "MotherVolume");
531 string motherVolumeVacuum = motherVolume + "Vac";
532
533 // If zz < 0 (positioned at negative z) vertical collimator is flipped when rotated into Mother Volume system
534 G4Scale3D scale;
535 G4Rotate3D rotation;
536 G4Translate3D translation;
537 elements[motherVolumeVacuum].transform.getDecomposition(scale, rotation, translation);
538 double zz = rotation.zz();
539
540 // d1, d2 are collimator jaws displacements from beam center, d1<0, d2>0
541 // Z is collimator position inside its Mother Volume
542 double collimator_d1 = m_config.getParameter(prep + "d1") * unitFactor;
543 double collimator_d2 = m_config.getParameter(prep + "d2") * unitFactor;
544 double collimator_fullH = m_config.getParameter(prep + "fullH") * unitFactor;
545 double collimator_headH = m_config.getParameter(prep + "headH") * unitFactor;
546 double collimator_minW = m_config.getParameter(prep + "minW") * unitFactor;
547 double collimator_maxW = m_config.getParameter(prep + "maxW") * unitFactor;
548 double collimator_th = m_config.getParameter(prep + "th") * unitFactor;
549 double collimator_Z = m_config.getParameter(prep + "Z") * unitFactor;
550
551 B2WARNING("Collimator " << name << " displacement d1 is set to " << collimator_d1 << "mm (must be negative)");
552 B2WARNING("Collimator " << name << " displacement d2 is set to " << collimator_d2 << "mm (must be positive)");
553
554
555 // Collimator heads
556
557 // dx1,2 dy1,2 dz are trapezoid dimensions
558 double head_dx1;
559 double head_dx2;
560 double head_dy1;
561 double head_dy2;
562 double head_dz = collimator_headH / 2.0;
563 if (type == "vertical") {
564 head_dx1 = collimator_th / 2.0;
565 head_dx2 = collimator_th / 2.0;
566 head_dy1 = ((collimator_maxW - collimator_minW) * collimator_headH / collimator_fullH + collimator_minW) / 2.0;
567 head_dy2 = collimator_minW / 2.0;
568 } else {
569 head_dx1 = ((collimator_maxW - collimator_minW) * collimator_headH / collimator_fullH + collimator_minW) / 2.0;
570 head_dx2 = collimator_minW / 2.0;
571 head_dy1 = collimator_th / 2.0;
572 head_dy2 = collimator_th / 2.0;
573 }
574
575 // storable elements
576 FarBeamLineElement collimator_head1;
577 FarBeamLineElement collimator_head2;
578
579 // move collimator to position on beam line
580 G4Transform3D transform_head1 = G4Translate3D(0.0, 0.0, collimator_Z);
581 G4Transform3D transform_head2 = G4Translate3D(0.0, 0.0, collimator_Z);
582
583 // rotate and move collimator jaws to their relative positions
584 if (type == "vertical") {
585 transform_head1 = transform_head1 * G4Translate3D(0.0, -head_dz + collimator_d1, 0.0);
586 transform_head1 = transform_head1 * G4RotateX3D(-M_PI / 2 / Unit::rad);
587
588 transform_head2 = transform_head2 * G4Translate3D(0.0, head_dz + collimator_d2, 0.0);
589 transform_head2 = transform_head2 * G4RotateX3D(M_PI / 2 / Unit::rad);
590 } else {
591 if (zz > 0) {
592 transform_head1 = transform_head1 * G4Translate3D(-head_dz + collimator_d1, 0.0, 0.0);
593 transform_head1 = transform_head1 * G4RotateY3D(M_PI / 2 / Unit::rad);
594
595 transform_head2 = transform_head2 * G4Translate3D(head_dz + collimator_d2, 0.0, 0.0);
596 transform_head2 = transform_head2 * G4RotateY3D(-M_PI / 2 / Unit::rad);
597 } else {
598 transform_head1 = transform_head1 * G4Translate3D(head_dz - collimator_d1, 0.0, 0.0);
599 transform_head1 = transform_head1 * G4RotateY3D(-M_PI / 2 / Unit::rad);
600
601 transform_head2 = transform_head2 * G4Translate3D(-head_dz - collimator_d2, 0.0, 0.0);
602 transform_head2 = transform_head2 * G4RotateY3D(M_PI / 2 / Unit::rad);
603 }
604 }
605
606 collimator_head1.transform = transform_head1;
607 collimator_head2.transform = transform_head2;
608
609 // define geometry
610 string geo_headx_name = "geo_" + name + "_headx_name";
611
612 string geo_head1_name = "geo_" + name + "_head1_name";
613 string geo_head2_name = "geo_" + name + "_head2_name";
614
615 G4VSolid* geo_headx = new G4Trd(geo_headx_name, head_dx1, head_dx2, head_dy1, head_dy2, head_dz);
616
617 G4VSolid* geo_head1 = new G4IntersectionSolid(geo_head1_name, geo_headx, elements[motherVolumeVacuum].geo,
618 collimator_head1.transform.inverse());
619 G4VSolid* geo_head2 = new G4IntersectionSolid(geo_head2_name, geo_headx, elements[motherVolumeVacuum].geo,
620 collimator_head2.transform.inverse());
621
622 collimator_head1.geo = geo_head1;
623 collimator_head2.geo = geo_head2;
624
625 // define logical volume
626 string strMat_head = m_config.getParameterStr(prep + "HeadMaterial");
627 G4Material* mat_head = Materials::get(strMat_head);
628 string logi_head1_name = "logi_" + name + "_head1_name";
629 string logi_head2_name = "logi_" + name + "_head2_name";
630 G4LogicalVolume* logi_head1 = new G4LogicalVolume(geo_head1, mat_head, logi_head1_name);
631 G4LogicalVolume* logi_head2 = new G4LogicalVolume(geo_head2, mat_head, logi_head2_name);
632 setColor(*logi_head1, "#CC0000");
633 setColor(*logi_head2, "#CC0000");
634 setVisibility(*logi_head1, false);
635 setVisibility(*logi_head2, false);
636
637 // check if collimator is inside beam pipe
638 double volume_head1 = logi_head1->GetSolid()->GetCubicVolume();
639 double volume_head2 = logi_head2->GetSolid()->GetCubicVolume();
640
641 collimator_head1.logi = logi_head1;
642 collimator_head2.logi = logi_head2;
643
644 // put volume
645 string phys_head1_name = "phys_" + name + "_head1" + "_name";
646 string phys_head2_name = "phys_" + name + "_head2" + "_name";
647 if (volume_head1 != 0)
648 new G4PVPlacement(collimator_head1.transform, logi_head1, phys_head1_name, elements[motherVolumeVacuum].logi, false, 0);
649 if (volume_head2 != 0)
650 new G4PVPlacement(collimator_head2.transform, logi_head2, phys_head2_name, elements[motherVolumeVacuum].logi, false, 0);
651
652 // to use it later in "intersect" and "subtract"
653 collimator_head1.transform = collimator_head1.transform * elements[motherVolumeVacuum].transform;
654 collimator_head2.transform = collimator_head2.transform * elements[motherVolumeVacuum].transform;
655
656 string name_head1 = name + "_head1";
657 string name_head2 = name + "_head2";
658 elements[name_head1] = collimator_head1;
659 elements[name_head2] = collimator_head2;
660
661
662 // Collimator bodies
663
664 // dx1,2 dy1,2 dz are trapezoid dimensions
665 double body_dx1;
666 double body_dx2;
667 double body_dy1;
668 double body_dy2;
669 double body_dz = (collimator_fullH - collimator_headH) / 2.0;
670 if (type == "vertical") {
671 body_dx1 = collimator_th / 2.0;
672 body_dx2 = collimator_th / 2.0;
673 body_dy1 = collimator_maxW / 2.0;
674 body_dy2 = ((collimator_maxW - collimator_minW) * collimator_headH / collimator_fullH + collimator_minW) / 2.0;
675 } else {
676 body_dx1 = collimator_maxW / 2.0;
677 body_dx2 = ((collimator_maxW - collimator_minW) * collimator_headH / collimator_fullH + collimator_minW) / 2.0;
678 body_dy1 = collimator_th / 2.0;
679 body_dy2 = collimator_th / 2.0;
680 }
681
682 // storable elements
683 FarBeamLineElement collimator_body1;
684 FarBeamLineElement collimator_body2;
685
686 // reuse head transfomation with additional shift
687 if (type == "vertical") {
688 collimator_body1.transform = G4Translate3D(0.0, -head_dz - body_dz, 0.0) * transform_head1;
689 collimator_body2.transform = G4Translate3D(0.0, head_dz + body_dz, 0.0) * transform_head2;
690 } else {
691 if (zz > 0) {
692 collimator_body1.transform = G4Translate3D(-head_dz - body_dz, 0.0, 0.0) * transform_head1;
693 collimator_body2.transform = G4Translate3D(head_dz + body_dz, 0.0, 0.0) * transform_head2;
694 } else {
695 collimator_body1.transform = G4Translate3D(head_dz + body_dz, 0.0, 0.0) * transform_head1;
696 collimator_body2.transform = G4Translate3D(-head_dz - body_dz, 0.0, 0.0) * transform_head2;
697 }
698 }
699
700 // define geometry
701 string geo_bodyx_name = "geo_" + name + "_bodyx_name";
702
703 string geo_body1_name = "geo_" + name + "_body1_name";
704 string geo_body2_name = "geo_" + name + "_body2_name";
705
706 G4VSolid* geo_bodyx = new G4Trd(geo_bodyx_name, body_dx1, body_dx2, body_dy1, body_dy2, body_dz);
707
708 G4VSolid* geo_body1 = new G4IntersectionSolid(geo_body1_name, geo_bodyx, elements[motherVolumeVacuum].geo,
709 collimator_body1.transform.inverse());
710 G4VSolid* geo_body2 = new G4IntersectionSolid(geo_body2_name, geo_bodyx, elements[motherVolumeVacuum].geo,
711 collimator_body2.transform.inverse());
712
713 collimator_body1.geo = geo_body1;
714 collimator_body2.geo = geo_body2;
715
716 // define logical volume
717 string strMat_body = m_config.getParameterStr(prep + "Material");
718 G4Material* mat_body = Materials::get(strMat_body);
719 string logi_body1_name = "logi_" + name + "_body1_name";
720 string logi_body2_name = "logi_" + name + "_body2_name";
721 G4LogicalVolume* logi_body1 = new G4LogicalVolume(geo_body1, mat_body, logi_body1_name);
722 G4LogicalVolume* logi_body2 = new G4LogicalVolume(geo_body2, mat_body, logi_body2_name);
723 setColor(*logi_body1, "#CC0000");
724 setColor(*logi_body2, "#CC0000");
725 setVisibility(*logi_body1, false);
726 setVisibility(*logi_body2, false);
727
728 // check if collimator is inside beam pipe
729 double volume_body1 = logi_body1->GetSolid()->GetCubicVolume();
730 double volume_body2 = logi_body2->GetSolid()->GetCubicVolume();
731
732 collimator_body1.logi = logi_body1;
733 collimator_body2.logi = logi_body2;
734
735 // put volume
736 string phys_body1_name = "phys_" + name + "_body1" + "_name";
737 string phys_body2_name = "phys_" + name + "_body2" + "_name";
738 if (volume_body1 != 0)
739 new G4PVPlacement(collimator_body1.transform, logi_body1, phys_body1_name, elements[motherVolumeVacuum].logi, false, 0);
740 if (volume_body2 != 0)
741 new G4PVPlacement(collimator_body2.transform, logi_body2, phys_body2_name, elements[motherVolumeVacuum].logi, false, 0);
742
743 // to use it later in "intersect" and "subtract"
744 collimator_body1.transform = collimator_body1.transform * elements[motherVolumeVacuum].transform;
745 collimator_body2.transform = collimator_body2.transform * elements[motherVolumeVacuum].transform;
746
747 string name_body1 = name + "_body1";
748 string name_body2 = name + "_body2";
749 elements[name_body1] = collimator_body1;
750 elements[name_body2] = collimator_body2;
751 }
752 }
753 }
755}
The Class for BeamBackground Sensitive Detector.
const std::string & getParameterStr(const std::string &name) const
Get string parameter.
Definition: IRGeoBase.h:64
double getParameter(const std::string &name) const
Get parameter value.
Definition: IRGeoBase.h:41
static const double mm
[millimeters]
Definition: Unit.h:70
static const double rad
Standard of [angle].
Definition: Unit.h:50
static const double cm
Standard units with the value = 1.
Definition: Unit.h:47
static G4Material * get(const std::string &name)
Find given material.
Definition: Materials.h:63
FarBeamLineGeo m_config
geometry parameters object
GeoFarBeamLineCreator()
Constructor of the GeoFarBeamLineCreator class.
virtual ~GeoFarBeamLineCreator()
The destructor of the GeoFarBeamLineCreator class.
void createGeometry(G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Create detector geometry.
SensitiveDetector * m_sensitive
Sensitive detector.
The IR Sensitive Detector class.
int intersect(const TRGCDCLpar &lp1, const TRGCDCLpar &lp2, CLHEP::HepVector &v1, CLHEP::HepVector &v2)
intersection
Definition: Lpar.cc:249
void setVisibility(G4LogicalVolume &volume, bool visible)
Helper function to quickly set the visibility of a given volume.
Definition: utilities.cc:108
void setColor(G4LogicalVolume &volume, const std::string &color)
Set the color of a logical volume.
Definition: utilities.cc:100
GeometryTypes
Flag indiciating the type of geometry to be used.
Abstract base class for different kinds of events.
STL namespace.
The struct for FarBeamLineElement.
G4LogicalVolume * logi
Logical volume.
G4Transform3D transform
Transformation.