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