| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // Code developed by
27 // Silvia Pozzi (1), silvia.pozzi@iss.it
28 // Barbara Caccia (1), barbara.caccia@iss.it
29 // Carlo Mancini Terracciano (2), carlo.mancin
30 // (1) Istituto Superiore di Sanita' and INFN
31 // (2) Univ. La Sapienza and INFN Roma, Italy
32
33 #include "DetectorConstruction.hh"
34 #include "DetectorMessenger.hh"
35
36 #include <G4LogicalVolume.hh>
37 #include <G4PVPlacement.hh>
38 #include <G4PVReplica.hh>
39 #include <G4ProductionCuts.hh>
40 #include <G4NistManager.hh>
41 #include <G4SystemOfUnits.hh>
42 #include <G4VisAttributes.hh>
43 #include <G4Box.hh>
44 #include <G4VSolid.hh>
45 #include <G4Cons.hh>
46 #include <G4Tubs.hh>
47 #include <G4BooleanSolid.hh>
48 #include <G4SubtractionSolid.hh>
49 #include <G4SDManager.hh>
50 #include <G4GlobalMagFieldMessenger.hh>
51 #include <G4UnionSolid.hh>
52 #include <G4SubtractionSolid.hh>
53 #include <G4RunManager.hh>
54 #include <G4GeometryManager.hh>
55 #include <G4PhysicalVolumeStore.hh>
56 #include <G4LogicalVolumeStore.hh>
57 #include <G4SolidStore.hh>
58 #include <G4Transform3D.hh>
59 #include <CLHEP/Vector/Rotation.h>
60
61
62 G4VPhysicalVolume* DetectorConstruction::Const
63 {
64 detectorMessenger = new DetectorMessenger(th
65
66 //// Saturn 43
67 G4Material* air = G4NistManager::Instance()-
68
69 G4double worldSizeX = 3.5 * m;
70 G4double worldSizeY = 3.5 * m;
71 G4double worldSizeZ = 3.5 * m;
72
73 G4VSolid* worldBox = new G4Box("world", worl
74
75 G4LogicalVolume* worldLog = new G4LogicalVol
76
77 G4VisAttributes* visAttr = new G4VisAttribut
78 visAttr->SetVisibility(false);
79 worldLog->SetVisAttributes(visAttr);
80
81 world_phys = new G4PVPlacement(nullptr, {},
82
83 ConstructMaterials();
84 ConstructAccelerator();
85 ConstructPhantom();
86
87 return world_phys;
88 }
89
90 void DetectorConstruction::ConstructPhantom()
91 {
92 G4Material* water = G4NistManager::Instance(
93
94 G4ThreeVector phantomHalfDim = {phantomSideD
95
96 G4VSolid* Phantom = new G4Box("Phantom", pha
97 G4LogicalVolume* Phantom_log = new G4Logical
98
99 G4VSolid* spessore = new G4Box("Spessore", p
100 G4LogicalVolume* spessore_log = new G4Logica
101 new G4PVPlacement(nullptr, {0.,0., 1.25*mm},
102
103 // axes origin on the face of the phantom
104 G4ThreeVector phantomPos = {0., 0., phantomH
105
106 phantom_phys = new G4PVPlacement(nullptr, ph
107 "phantom_phys", Phantom_log, world_ph
108
109 PhantomSegmentation(Phantom_log);
110
111 // Region for cuts
112 G4Region *regVol= new G4Region("fullWaterPha
113 G4ProductionCuts* cuts = new G4ProductionCut
114 cuts->SetProductionCut(0.1*mm);
115 regVol->SetProductionCuts(cuts);
116
117 Phantom_log -> SetRegion(regVol);
118 regVol->AddRootLogicalVolume(Phantom_log);
119
120 //--------- Visualization attributes -------
121 /// Phantom
122 G4VisAttributes* simpleH2OVisAtt= new G4VisA
123 simpleH2OVisAtt->SetVisibility(true);
124 simpleH2OVisAtt->SetForceSolid(false);
125 Phantom_log->SetVisAttributes(simpleH2OVisAt
126
127 }
128
129 G4Material* DetectorConstruction::GetMaterial(
130 {
131 G4Material* material = nullptr;
132
133 if (materialName == "kapton")
134 {
135 material = mat_Kapton;
136 }
137 else if (materialName == "steel")
138 {
139 material = mat_Ssteel;
140 }
141 else if (materialName == "xc10")
142 {
143 material = mat_XC10;
144 }
145 else if (materialName == "wnicu")
146 {
147 material = mat_WNICU;
148 }
149 else
150 {
151 G4cerr << "* DetectorConstruction::GetMa
152 materialName << " not found" << G4en
153 }
154
155 return material;
156 }
157
158 void DetectorConstruction::ConstructMaterials(
159 {
160 G4NistManager* nist = G4NistManager::Insta
161
162 G4Material* el_H = nist -> FindOrBuildMater
163 G4Material* el_C = nist -> FindOrBuildMater
164 G4Material* el_N = nist -> FindOrBuildMater
165 G4Material* el_O = nist -> FindOrBuildMater
166 G4Material* el_Si = nist -> FindOrBuildMater
167 G4Material* el_Cr = nist -> FindOrBuildMater
168 G4Material* el_Mn = nist -> FindOrBuildMater
169 G4Material* el_Fe = nist -> FindOrBuildMater
170 G4Material* el_Ni = nist -> FindOrBuildMater
171 G4Material* el_Cu = nist -> FindOrBuildMater
172 G4Material* el_W = nist -> FindOrBuildMater
173
174 mat_Ssteel = new G4Material("StainlessSteel"
175 mat_Ssteel -> AddMaterial(el_Fe, 0.6898);
176 mat_Ssteel -> AddMaterial(el_Cr, 0.18);
177 mat_Ssteel -> AddMaterial(el_Ni, 0.10);
178 mat_Ssteel -> AddMaterial(el_Mn, 0.02);
179 mat_Ssteel -> AddMaterial(el_Si, 0.01);
180 mat_Ssteel -> AddMaterial(el_C, 0.0002);
181
182 mat_XC10 = new G4Material("CARBON_STEEL", 7.
183 mat_XC10 -> AddMaterial(el_Fe, 0.993);
184 mat_XC10 -> AddMaterial(el_Mn, 0.006);
185 mat_XC10 -> AddMaterial(el_C, 0.001);
186
187 mat_WNICU = new G4Material("Denal(WNICU)", 1
188 mat_WNICU -> AddMaterial(el_W, 0.905);
189 mat_WNICU -> AddMaterial(el_Ni, 0.07);
190 mat_WNICU -> AddMaterial(el_Cu, 0.025);
191
192 mat_Kapton = new G4Material("Kapton", 1.42*g
193 mat_Kapton -> AddMaterial(el_C, 69.1133 *per
194 mat_Kapton -> AddMaterial(el_O, 20.9235 *per
195 mat_Kapton -> AddMaterial(el_N, 7.3270 *perC
196 mat_Kapton -> AddMaterial(el_H, 2.6362 *perC
197
198 }
199
200 void DetectorConstruction::ConstructAccelerato
201 {
202 G4bool checkOverlap = true;
203
204 G4Material* Vacuum = G4NistManager::Instance
205
206 accHalfSize = {600.*mm, 600.*mm, 600.*mm};
207 G4ThreeVector initialCentre = {0.*mm, 0.*mm,
208
209 G4Box* accWorldBox = new G4Box("accWorld", a
210 G4LogicalVolume* accWorldLV = new G4LogicalV
211 accWorld_phys = new G4PVPlacement(nullptr, i
212
213 G4VisAttributes* simpleAlSVisAtt = new G4Vis
214 simpleAlSVisAtt -> SetVisibility(false);
215 accWorldLV -> SetVisAttributes(simpleAlSVisA
216
217 ConstructTarget();
218 ConstructPrimaryCollimator();
219 ConstructVacuumWindow();
220 ConstructFlatteningFilter();
221 ConstructIonizationChamber();
222 ConstructJawsX();
223 ConstructJawsY();
224 }
225
226 void DetectorConstruction::ConstructTarget()
227 {
228 G4bool checkOverlap = true;
229
230 G4ThreeVector targetCentre, boxHalfSide, tub
231 G4double tubeRadius, tubeHalfLengthZ;
232 targetCentre.set(0, 0, -3.5*mm);
233 boxHalfSide.set(5.*mm, 5.*mm, 7.5*mm);
234
235 tubeRadius = 3.*mm;
236 tubeHalfLengthZ = 5.5*mm; // 11 mm tot
237
238 G4Material* diskMaterial = G4NistManager::In
239 G4Material* targetMaterial = G4NistManager::
240
241 // Physical volumes
242 G4Box* box = new G4Box("targetBox", boxHalfS
243 G4Tubs* tube = new G4Tubs("targetTube", 0.*m
244
245 G4SubtractionSolid* BoxMinusTube = new G4Sub
246 G4LogicalVolume* BoxMinusTubeLV = new G4Logi
247 new G4PVPlacement(0, targetCentre, "TargetPV
248
249 G4Tubs* disk = new G4Tubs("targetDisk", 0.*m
250 G4LogicalVolume* diskLV = new G4LogicalVolum
251 new G4PVPlacement(0, {0.,0.,-15.*mm}, "diskP
252
253 // Visualization
254 G4VisAttributes* simpleAlSVisAtt;
255 simpleAlSVisAtt = new G4VisAttributes(G4Colo
256 simpleAlSVisAtt -> SetVisibility(true);
257 simpleAlSVisAtt -> SetForceSolid(true);
258 BoxMinusTubeLV -> SetVisAttributes(simpleAl
259 diskLV -> SetVisAttributes(simpleAlSVisAtt);
260
261 // Region for cuts
262 G4Region *regVol;
263 regVol = new G4Region("TargetR");
264 G4ProductionCuts* cuts = new G4ProductionCut
265 cuts -> SetProductionCut(1.*mm);
266 regVol -> SetProductionCuts(cuts);
267 BoxMinusTubeLV -> SetRegion(regVol);
268 regVol -> AddRootLogicalVolume(BoxMinusTubeL
269 diskLV -> SetRegion(regVol);
270 regVol -> AddRootLogicalVolume(diskLV);
271 }
272
273 void DetectorConstruction::ConstructPrimaryCol
274 {
275 G4bool checkOverlap = true;
276
277 G4double bigTube1HalfLengthZ, bigTube1Radius
278 G4double bigTube2HalfLengthZ, bigTube2Radius
279 G4double bigTube3HalfLengthZ, bigTube3Radius
280
281 G4double tube1Radius, tube1HalfLengthZ, tube
282
283 G4double cone1RadiusMin, cone1RadiusMax, con
284 G4double cone2RadiusMin, cone2RadiusMax, con
285 G4double cone3RadiusMin, cone3RadiusMax, con
286
287 // upper principal tube = bigTube1
288 bigTube1Radius = 141./2.*mm;
289 bigTube1HalfLengthZ = 79./2. * mm;
290 bigTube1Z = 5. * mm + bigTube1HalfLe
291 tube1Radius = 34./2. * mm;
292 tube1HalfLengthZ = 15.86/2. * mm;
293 cone1RadiusMin = 34./2.*mm; // upper cone
294 cone1RadiusMax = 54./2.*mm;
295 cone1HalfLengthZ = (bigTube1HalfLengthZ - t
296
297 // middle principal tube = bigTube2
298 bigTube2Radius = 141./2. * mm;
299 bigTube2HalfLengthZ = 57.5/2. * mm;
300 bigTube2Z = bigTube1Z + bigTube1Half
301 tube2Radius = 53./2. * mm;
302 tube2HalfLengthZ = 1.35/2. * mm;
303 cone2RadiusMin = 53./2.*mm; // middle con
304 cone2RadiusMax = 81./2.*mm;
305 cone2HalfLengthZ = bigTube2HalfLengthZ - tu
306
307 // lower principal tube = bigTube3
308 bigTube3Radius = 241./2.*mm;
309 bigTube3HalfLengthZ = 40.50/2.*mm;
310 bigTube3Z = bigTube2Z + bigTube2HalfLe
311 cone3RadiusMin = 88.06/2.*mm; // lower co
312 cone3RadiusMax = 109.76/2.*mm;
313 cone3HalfLengthZ = 40.50/2.*mm;
314
315 G4Material* upperTubeMaterial = GetMaterial(
316 G4Material* middleTubeMaterial = GetMaterial
317 G4Material* lowerTubeMaterial = G4NistManage
318
319 // Physical volumes
320 G4Tubs* bigTube1 = new G4Tubs("PriCollBigTub
321 G4Tubs* tube1 = new G4Tubs("PriCollTube1", 0
322 G4Cons* cone1 = new G4Cons("PriCollCone1", 0
323
324 G4Tubs* bigTube2 = new G4Tubs("PriCollBigTub
325 G4Tubs* tube2 = new G4Tubs("PriCollTube2", 0
326 G4Cons* cone2 = new G4Cons("PriCollCone2", 0
327
328 G4Tubs* bigTube3 = new G4Tubs("PriCollBigTub
329 G4Cons* cone3 = new G4Cons("PriCollCone3", 0
330
331 G4UnionSolid* tube1AndCone1 = new G4UnionSol
332 tube1, cone1, nullptr, {0., 0., (tube1
333 G4SubtractionSolid* bigTube1NotTube1AndCone1
334 "PriColltube1NotTube1AndCone1", bigTube1
335 {0., 0., (-bigTube1HalfLengthZ + tube1Ha
336 G4LogicalVolume* bigTube1NotTube1AndCone1LV
337 bigTube1NotTube1AndCone1, upperTubeMater
338
339 G4UnionSolid* tube2AndCone2 = new G4UnionSol
340 tube2, cone2, nullptr, {0., 0., (tube2
341 G4SubtractionSolid* bigTube2NotTube2AndCone2
342 "PriCollBigTube2NotTube2Cone2", bigTube2
343 {0., 0., (-bigTube2HalfLengthZ + tube2Ha
344 G4LogicalVolume* bigTube2NotTube2AndCone2LV
345 bigTube2NotTube2AndCone2, middleTubeMate
346
347 G4SubtractionSolid* bigTube3NotCone3 = new G
348 G4LogicalVolume* bigTube3NotCone3LV = new G
349 lowerTubeMaterial, "PriCollBigTube3NotCo
350
351 new G4PVPlacement(nullptr, {0, 0, bigTube1Z}
352 bigTube1NotTube1AndCone1LV, accWorld_phy
353 new G4PVPlacement(nullptr, {0, 0, bigTube2Z}
354 bigTube2NotTube2AndCone2LV, accWorld_phy
355 new G4PVPlacement(nullptr, {0, 0, bigTube3Z}
356 bigTube3NotCone3LV, accWorld_phys, false
357
358 // Visualization
359 G4VisAttributes* simpleAlSVisAtt1 = new G4Vi
360 G4VisAttributes* simpleAlSVisAtt2 = new G4Vi
361 G4VisAttributes* simpleAlSVisAtt3 = new G4Vi
362 simpleAlSVisAtt1 -> SetVisibility(true);
363 simpleAlSVisAtt1 -> SetForceSolid(false);
364 simpleAlSVisAtt1 -> SetForceWireframe(true);
365 bigTube1NotTube1AndCone1LV -> SetVisAttribut
366 simpleAlSVisAtt2 -> SetVisibility(true);
367 simpleAlSVisAtt2 -> SetForceSolid(false);
368 simpleAlSVisAtt2 -> SetForceWireframe(true);
369 bigTube2NotTube2AndCone2LV -> SetVisAttribut
370 simpleAlSVisAtt3 -> SetVisibility(true);
371 simpleAlSVisAtt3 -> SetForceSolid(false);
372 simpleAlSVisAtt3 -> SetForceWireframe(true);
373 bigTube3NotCone3LV -> SetVisAttributes(simpl
374
375 // Region for cuts
376 G4Region* regVol = new G4Region("primaryColl
377 G4ProductionCuts* cuts = new G4ProductionCut
378 cuts->SetProductionCut(1.*cm);
379 regVol->SetProductionCuts(cuts);
380 bigTube1NotTube1AndCone1LV -> SetRegion(regV
381 regVol -> AddRootLogicalVolume(bigTube1NotTu
382 bigTube2NotTube2AndCone2LV -> SetRegion(regV
383 regVol -> AddRootLogicalVolume(bigTube2NotTu
384 bigTube3NotCone3LV -> SetRegion(regVol);
385 regVol -> AddRootLogicalVolume(bigTube3NotCo
386 }
387
388 void DetectorConstruction::ConstructFlattening
389 {
390 G4double tube1HalfLengthZ; //tube1Radius,
391 G4double cone1RadiusMin, cone1RadiusMax, con
392 G4double cone2RadiusMin, cone2RadiusMax, con
393 G4double cone3RadiusMin, cone3RadiusMax, con
394
395 tubeFFRadius = 108./2.*mm; // top principa
396 tube1HalfLengthZ = 7.5/2.*mm;
397
398 cone1RadiusMin = 54./2.*mm; // top cone
399 cone1RadiusMax = 76./2.*mm;
400 cone1HalfLengthZ = 13.7/2.*mm;
401
402 cone2RadiusMin = 8./2.*mm; // mid cone
403 cone2RadiusMax = 54./2.*mm;
404 cone2HalfLengthZ = (44.3-13.7)/2.*mm;
405
406 cone3RadiusMin = 0.000001*mm; // bottom co
407 cone3RadiusMax = 8./2.*mm;
408 cone3HalfLengthZ = (46.8-44.3)/2.*mm;
409 tubeFFFirstFaceZ = 149.50*mm;
410
411 ffZ = tubeFFFirstFaceZ + tube1HalfLengthZ; /
412
413 G4Material* ffMaterial = GetMaterial("steel"
414
415 // Physical volumes
416 G4Tubs* tube1 = new G4Tubs("ffTube1", 0.*mm,
417 G4Cons* cone1 = new G4Cons("ffCone1", 0., co
418 G4Cons* cone2 = new G4Cons("ffCone2", 0., co
419 G4Cons* cone3 = new G4Cons("ffCone3", 0., co
420
421 G4double pos = (cone1HalfLengthZ - tube1Half
422 G4UnionSolid *tube1AndCone1 = new G4UnionSol
423 cone1, 0, {0., 0., pos});
424
425 pos += cone1HalfLengthZ + cone2HalfLengthZ;
426 G4UnionSolid* tubeCone1AndCone2 = new G4Unio
427 tube1AndCone1, cone2, 0, {0., 0., pos});
428
429 pos += cone2HalfLengthZ + cone3HalfLengthZ;
430 G4UnionSolid* tubeCone12AndCone3 = new G4Uni
431 tubeCone1AndCone2, cone3, 0, {0., 0., po
432
433 G4LogicalVolume* ffLV = new G4LogicalVolume(
434
435 new G4PVPlacement(0, {0,0,ffZ}, "ffPV", ffLV
436
437 // Visualization
438 G4VisAttributes* simpleAlSVisAtt= new G4VisA
439 simpleAlSVisAtt->SetVisibility(true);
440 simpleAlSVisAtt->SetForceSolid(true);
441 ffLV->SetVisAttributes(simpleAlSVisAtt);
442
443 // Region for cuts
444 G4Region* regVol;
445 regVol= new G4Region("ffR");
446 G4ProductionCuts* cuts = new G4ProductionCut
447 cuts->SetProductionCut(0.1*cm);
448 regVol->SetProductionCuts(cuts);
449 ffLV->SetRegion(regVol);
450 regVol->AddRootLogicalVolume(ffLV);
451 }
452
453 void DetectorConstruction::ConstructIonization
454 {
455 G4bool checkOverlap = true;
456
457 G4double tubeRadius, tubeA1Z, tubeA2Z, tubeA
458 G4double kaptonThickness, AlThickness1, AlTh
459
460 kaptonThickness = 0.025/2. *mm;
461 AlThickness1 = 1.6e-5/2. *cm;
462 AlThickness8 = 8e-6/2. *cm;
463
464 tubeRadius = 110./2.*mm; // upper principal
465 tubeA1Z = (202.5 + AlThickness8)*mm;
466 tubeA2Z = (204.5 + AlThickness1)*mm;
467 tubeA3Z = (206.5 + AlThickness8)*mm;
468 tubeA4Z = (207.5 + AlThickness8)*mm;
469 tubeA5Z = (209.5 + AlThickness1)*mm;
470 tubeA6Z = (211.5 + AlThickness8)*mm;
471
472 G4double pos1 = (AlThickness1 + kaptonThickn
473 G4double pos8 = (AlThickness8 + kaptonThickn
474
475 G4Material* el_Al = G4NistManager::Instance(
476 G4Material* Kapton = GetMaterial("kapton");
477
478 // Physical volumes
479 G4Tubs* tubeAl8 = new G4Tubs("ICTube1A", 0.*
480 G4Tubs* tubeK = new G4Tubs("ICTube1B", 0.*mm
481 G4Tubs* tubeAl1 = new G4Tubs("ICTube2A", 0.*
482
483 G4LogicalVolume* IC_Al8_LV = new G4LogicalVo
484 G4LogicalVolume* IC_Al1_LV = new G4LogicalVo
485 G4LogicalVolume* IC_K_LV = new G4LogicalVo
486
487 new G4PVPlacement(0, {0,0,tubeA1Z}, "IC_AL8_
488 new G4PVPlacement(0, {0,0,tubeA1Z + pos8}, "
489
490 new G4PVPlacement(0, {0,0,tubeA2Z}, "IC_AL1_
491 new G4PVPlacement(0, {0,0,tubeA2Z + pos1}, "
492
493 new G4PVPlacement(0, {0,0,tubeA3Z}, "IC_AL8_
494 new G4PVPlacement(0, {0,0,tubeA3Z + pos8}, "
495
496 new G4PVPlacement(0, {0,0,tubeA4Z}, "IC_AL8_
497 new G4PVPlacement(0, {0,0,tubeA4Z + pos8}, "
498
499 new G4PVPlacement(0, {0,0,tubeA5Z}, "IC_AL1_
500 new G4PVPlacement(0, {0,0,tubeA5Z + pos1}, "
501
502 new G4PVPlacement(0, {0,0,tubeA6Z}, "IC_AL8_
503 new G4PVPlacement(0, {0,0,tubeA6Z + pos8}, "
504
505 // Visualization
506 G4VisAttributes* simpleAlSVisAttAl1 = new G4
507 simpleAlSVisAttAl1 -> SetVisibility(true);
508 simpleAlSVisAttAl1 -> SetForceSolid(true);
509 IC_Al1_LV -> SetVisAttributes(simpleAlSVisAt
510
511 G4VisAttributes* simpleAlSVisAttAl8 = new G4
512 simpleAlSVisAttAl8 -> SetVisibility(true);
513 simpleAlSVisAttAl8 -> SetForceSolid(true);
514 IC_Al8_LV -> SetVisAttributes(simpleAlSVisAt
515
516 G4VisAttributes* simpleAlSVisAttK = new G4Vi
517 simpleAlSVisAttK -> SetVisibility(true);
518 simpleAlSVisAttK -> SetForceSolid(true);
519 IC_K_LV -> SetVisAttributes(simpleAlSVisAttK
520
521 // Region for cuts
522 G4Region * regVol = new G4Region("IonChamber
523 G4ProductionCuts* cuts = new G4ProductionCut
524 cuts -> SetProductionCut(0.1*mm);
525 regVol -> SetProductionCuts(cuts);
526
527 IC_Al1_LV -> SetRegion(regVol);
528 regVol -> AddRootLogicalVolume(IC_Al1_LV);
529 IC_Al8_LV -> SetRegion(regVol);
530 regVol -> AddRootLogicalVolume(IC_Al8_LV);
531 IC_K_LV -> SetRegion(regVol);
532 regVol -> AddRootLogicalVolume(IC_K_LV);
533 }
534
535 void DetectorConstruction::ConstructVacuumWind
536 {
537 G4double tubeRadius, tubeHalfLengthZ, tubeZ;
538
539 tubeRadius = 116.53/2.*mm; // upper princ
540 tubeHalfLengthZ = 2./2.*mm;
541 tubeZ = 215.75*mm + tubeHalfLengthZ;
542
543 G4Material* el_Al= G4NistManager::Instance()
544
545 // Physical volumes
546 G4Tubs* tube = new G4Tubs("VacWindowTube", 0
547 G4LogicalVolume* tubeLV = new G4LogicalVolum
548
549 new G4PVPlacement(0, {0,0,tubeZ}, "VacWindow
550
551 // Visualization
552 G4VisAttributes* simpleAlSVisAtt = new G4Vis
553 simpleAlSVisAtt -> SetVisibility(true);
554 simpleAlSVisAtt -> SetForceSolid(true);
555 tubeLV -> SetVisAttributes(simpleAlSVisAtt);
556
557 // Region for cuts
558 G4Region* regVol= new G4Region("VacWindowR")
559 G4ProductionCuts* cuts = new G4ProductionCut
560 cuts->SetProductionCut(0.1*cm);
561 regVol->SetProductionCuts(cuts);
562 tubeLV -> SetRegion(regVol);
563 regVol -> AddRootLogicalVolume(tubeLV);
564 }
565 void DetectorConstruction::ConstructJawsX()
566 {
567 G4bool checkOverlap = false;
568
569 G4double boxSide = 101.*mm;
570 G4double box1HalfLengthZ = 3.*mm;
571 G4double box2HalfLengthZ = 35.*mm;
572 G4double box3HalfLengthZ = 35.*mm;
573 G4double box4HalfLengthZ = 27.*mm;
574 G4double box5HalfLengthZ = 10.*mm;
575 G4double box6HalfLengthZ = 5.*mm;
576
577 G4Material* el_Pb = G4NistManager::Instance(
578 G4Material* XC10 = GetMaterial("xc10");
579 G4Material* WNICU = GetMaterial("wnicu");
580
581 G4ThreeVector boxJawXSide = {boxSide, boxSid
582
583 G4Box* boxJaw1X = new G4Box("Jaw1Xbox", boxJ
584 G4LogicalVolume* boxJaw1XLV = new G4LogicalV
585 G4NistManager::Instance()->FindOrBuildMa
586 G4Box* boxJaw2X = new G4Box("Jaw2Xbox", boxJ
587 G4LogicalVolume* boxJaw2XLV = new G4LogicalV
588 G4NistManager::Instance()->FindOrBuild
589
590 jaw1XInitialPos.setX(boxJawXSide.getX()/2.);
591 jaw1XInitialPos.setY(0.);
592 jaw1XInitialPos.setZ(275.5*mm + boxJawXSide.
593
594 G4RotationMatrix* cRotationJaw1X = new G4Rot
595 cRotationJaw1X -> rotateY(std::fabs(std::ata
596
597 G4ThreeVector position = *cRotationJaw1X * j
598 *cRotationJaw1X = cRotationJaw1X -> inverse(
599
600 boxJaw1X_phys =
601 new G4PVPlacement (cRotationJaw1X, positio
602 "Jaw1XPV", boxJaw1XLV, accWorld_phys,
603
604 jaw2XInitialPos.setX(-jaw1XInitialPos.getX()
605 jaw2XInitialPos.setY(jaw1XInitialPos.getY())
606 jaw2XInitialPos.setZ(jaw1XInitialPos.getZ())
607
608 G4RotationMatrix* cRotationJaw2X = new G4Rot
609 cRotationJaw2X -> rotateY(-std::fabs(std::at
610
611 position = *cRotationJaw2X * jaw2XInitialPos
612 *cRotationJaw2X = cRotationJaw2X -> inverse(
613
614 boxJaw2X_phys =
615 new G4PVPlacement (cRotationJaw2X, pos
616 "Jaw2XPV", boxJaw2XLV, accWorld_ph
617
618 // Physical volumes
619 G4Box* box1 = new G4Box("Jaws1XBox1", boxSid
620 G4Box* box2 = new G4Box("Jaws1XBox2", boxSid
621 G4Box* box3 = new G4Box("Jaws1XBox3", boxSid
622 G4Box* box4 = new G4Box("Jaws1XBox4", boxSid
623 G4Box* box5 = new G4Box("Jaws1XBox5", boxSid
624 G4Box* box6 = new G4Box("Jaws1XBox6", boxSid
625 G4LogicalVolume* box1LV1 = new G4LogicalVolu
626 G4LogicalVolume* box2LV1 = new G4LogicalVolu
627 G4LogicalVolume* box3LV1 = new G4LogicalVolu
628 G4LogicalVolume* box4LV1 = new G4LogicalVolu
629 G4LogicalVolume* box5LV1 = new G4LogicalVolu
630 G4LogicalVolume* box6LV1 = new G4LogicalVolu
631
632 G4LogicalVolume* box1LV2 = new G4LogicalVolu
633 G4LogicalVolume* box2LV2 = new G4LogicalVolu
634 G4LogicalVolume* box3LV2 = new G4LogicalVolu
635 G4LogicalVolume* box4LV2 = new G4LogicalVolu
636 G4LogicalVolume* box5LV2 = new G4LogicalVolu
637 G4LogicalVolume* box6LV2 = new G4LogicalVolu
638
639 G4ThreeVector centre = {0., 0., 0.};
640 G4double zCentreCurrentBox = -boxJawXSide.ge
641
642 centre.setZ(zCentreCurrentBox);
643 new G4PVPlacement(nullptr, centre, "Jaws1XPV
644 new G4PVPlacement(nullptr, centre, "Jaws2XPV
645
646 zCentreCurrentBox += (box1HalfLengthZ+box2Ha
647 centre.setZ(zCentreCurrentBox);
648 new G4PVPlacement(nullptr, centre, "Jaws1XPV
649 new G4PVPlacement(nullptr, centre, "Jaws2XPV
650
651 zCentreCurrentBox += (box2HalfLengthZ+box3Ha
652 centre.setZ(zCentreCurrentBox);
653 new G4PVPlacement(nullptr, centre, "Jaws1XPV
654 new G4PVPlacement(nullptr, centre, "Jaws2XPV
655
656 zCentreCurrentBox += (box3HalfLengthZ+box4Ha
657 centre.setZ(zCentreCurrentBox);
658 new G4PVPlacement(nullptr, centre, "Jaws1XPV
659 new G4PVPlacement(nullptr, centre, "Jaws2XPV
660
661 zCentreCurrentBox += (box4HalfLengthZ+box5Ha
662 centre.setZ(zCentreCurrentBox);
663 new G4PVPlacement(nullptr, centre, "Jaws1XPV
664 new G4PVPlacement(nullptr, centre, "Jaws2XPV
665
666 zCentreCurrentBox += (box5HalfLengthZ+box6Ha
667 centre.setZ(zCentreCurrentBox);
668 new G4PVPlacement(nullptr, centre, "Jaws1XPV
669 new G4PVPlacement(nullptr, centre, "Jaws2XPV
670
671 // Visualization
672 G4VisAttributes* visAtt = new G4VisAttribute
673 visAtt -> SetVisibility(false);
674 visAtt -> SetForceSolid(false);
675
676 G4VisAttributes* simpleAlSVisAttPb = new G4V
677 simpleAlSVisAttPb -> SetVisibility(true);
678 simpleAlSVisAttPb -> SetForceSolid(true);
679
680 G4VisAttributes* simpleAlSVisAttXC10 = new G
681 simpleAlSVisAttXC10 -> SetVisibility(true);
682 simpleAlSVisAttXC10 ->SetForceSolid(true);
683
684 G4VisAttributes* simpleAlSVisAttWNICU = new
685 simpleAlSVisAttWNICU ->SetVisibility(true);
686 simpleAlSVisAttWNICU ->SetForceSolid(true);
687
688 box1LV1 -> SetVisAttributes(simpleAlSVisAttX
689 box2LV1 -> SetVisAttributes(simpleAlSVisAttP
690 box3LV1 -> SetVisAttributes(simpleAlSVisAttW
691 box4LV1 -> SetVisAttributes(simpleAlSVisAttP
692 box5LV1 -> SetVisAttributes(simpleAlSVisAttP
693 box6LV1 -> SetVisAttributes(simpleAlSVisAttW
694
695 box1LV2 -> SetVisAttributes(simpleAlSVisAttX
696 box2LV2 -> SetVisAttributes(simpleAlSVisAttP
697 box3LV2 -> SetVisAttributes(simpleAlSVisAttW
698 box4LV2 -> SetVisAttributes(simpleAlSVisAttP
699 box5LV2 -> SetVisAttributes(simpleAlSVisAttP
700 box6LV2 -> SetVisAttributes(simpleAlSVisAttW
701
702 boxJaw1XLV -> SetVisAttributes(visAtt);
703 boxJaw2XLV -> SetVisAttributes(visAtt);
704
705 // Region for cuts
706 G4Region* regVol = new G4Region("JawsXR");
707 G4ProductionCuts* cuts = new G4ProductionCut
708 cuts -> SetProductionCut(2.*cm);
709 regVol -> SetProductionCuts(cuts);
710 box1LV1 -> SetRegion(regVol);
711 regVol -> AddRootLogicalVolume(box1LV1);
712 box2LV1 -> SetRegion(regVol);
713 regVol -> AddRootLogicalVolume(box2LV1);
714 box3LV1 -> SetRegion(regVol);
715 regVol -> AddRootLogicalVolume(box3LV1);
716 box4LV1 -> SetRegion(regVol);
717 regVol -> AddRootLogicalVolume(box4LV1);
718 box5LV1 -> SetRegion(regVol);
719 regVol -> AddRootLogicalVolume(box5LV1);
720 box6LV1 -> SetRegion(regVol);
721 regVol -> AddRootLogicalVolume(box6LV1);
722
723 box1LV2 -> SetRegion(regVol);
724 regVol -> AddRootLogicalVolume(box1LV2);
725 box2LV2 -> SetRegion(regVol);
726 regVol -> AddRootLogicalVolume(box2LV2);
727 box3LV2 -> SetRegion(regVol);
728 regVol -> AddRootLogicalVolume(box3LV2);
729 box4LV2 -> SetRegion(regVol);
730 regVol -> AddRootLogicalVolume(box4LV2);
731 box5LV2 -> SetRegion(regVol);
732 regVol -> AddRootLogicalVolume(box5LV2);
733 box6LV2 -> SetRegion(regVol);
734 regVol -> AddRootLogicalVolume(box6LV2);
735 }
736
737 void DetectorConstruction::ConstructJawsY()
738 {
739 G4bool checkOverlap = false;
740
741 G4double boxSide = 101.*mm;
742 G4double box1HalfLengthZ = 15.*mm;
743 G4double box2HalfLengthZ = 31.*mm;
744 G4double box3HalfLengthZ = 35.*mm;
745 G4double box4HalfLengthZ = 21.*mm;
746
747 G4Material* el_Pb = G4NistManager::Instance(
748 G4Material* XC10 = GetMaterial("xc10");
749 G4Material* WNICU = GetMaterial("wnicu");
750
751 G4ThreeVector boxJawYSide = {boxSide, boxSid
752
753 G4Box* boxJaw1Y = new G4Box("Jaw1Ybox", boxJ
754 G4LogicalVolume* boxJaw1YLV = new G4LogicalV
755 G4NistManager::Instance()->FindOrBuildMa
756 G4Box* boxJaw2Y = new G4Box("Jaw2Ybox", boxJ
757 G4LogicalVolume* boxJaw2YLV = new G4Logica
758 G4NistManager::Instance()->FindOrBuild
759
760 jaw1YInitialPos.setX(0.);
761 jaw1YInitialPos.setY(boxJawYSide.getX()/2.);
762 jaw1YInitialPos.setZ(248.5*mm + boxJawYSide.
763
764 G4RotationMatrix* cRotationJaw1Y = new G4Rot
765 cRotationJaw1Y -> rotateX(-std::fabs(std::at
766
767 G4ThreeVector position = *cRotationJaw1Y * j
768 *cRotationJaw1Y = cRotationJaw1Y -> inverse(
769
770 boxJaw1Y_phys =
771 new G4PVPlacement (cRotationJaw1Y, posit
772 "Jaw1YPV", boxJaw1YLV, accWorld_phys
773
774 jaw2YInitialPos.setX(jaw1YInitialPos.getX())
775 jaw2YInitialPos.setY(-jaw1YInitialPos.getY()
776 jaw2YInitialPos.setZ(jaw1YInitialPos.getZ())
777
778 G4RotationMatrix* cRotationJaw2Y = new G4Rot
779 cRotationJaw2Y -> rotateX(std::fabs(std::ata
780
781 position = *cRotationJaw2Y * jaw2YInitialPos
782 *cRotationJaw2Y = cRotationJaw2Y -> inverse(
783
784 boxJaw2Y_phys =
785 new G4PVPlacement (cRotationJaw2Y, posit
786 "Jaw2YPV", boxJaw2YLV, accWorld_phys
787
788 // Physical volumes
789 G4Box *box1 = new G4Box("Jaws1YBox1", boxSid
790 G4Box *box2 = new G4Box("Jaws1YBox2", boxSid
791 G4Box *box3 = new G4Box("Jaws1YBox3", boxSid
792 G4Box *box4 = new G4Box("Jaws1YBox4", boxSid
793 G4LogicalVolume* box1LV1 = new G4LogicalVolu
794 G4LogicalVolume* box2LV1 = new G4LogicalVolu
795 G4LogicalVolume* box3LV1 = new G4LogicalVolu
796 G4LogicalVolume* box4LV1 = new G4LogicalVolu
797
798 G4LogicalVolume* box1LV2 = new G4LogicalVolu
799 G4LogicalVolume* box2LV2 = new G4LogicalVolu
800 G4LogicalVolume* box3LV2 = new G4LogicalVolu
801 G4LogicalVolume* box4LV2 = new G4LogicalVolu
802
803 G4ThreeVector centre = {0., 0., 0.};
804 G4double zCentreCurrentBox = -boxJawYSide.ge
805
806 centre.setZ(zCentreCurrentBox);
807 new G4PVPlacement(nullptr, centre, "Jaws1YPV
808 new G4PVPlacement(nullptr, centre, "Jaws2YPV
809
810 zCentreCurrentBox += (box1HalfLengthZ+box2Ha
811 centre.setZ(zCentreCurrentBox);
812 new G4PVPlacement(nullptr, centre, "Jaws1YPV
813 new G4PVPlacement(nullptr, centre, "Jaws2YPV
814
815 zCentreCurrentBox += (box2HalfLengthZ+box3Ha
816 centre.setZ(zCentreCurrentBox);
817 new G4PVPlacement(nullptr, centre, "Jaws1YPV
818 new G4PVPlacement(nullptr, centre, "Jaws2YPV
819
820 zCentreCurrentBox += (box3HalfLengthZ+box4Ha
821 centre.setZ(zCentreCurrentBox);
822 new G4PVPlacement(nullptr, centre, "Jaws1YPV
823 new G4PVPlacement(nullptr, centre, "Jaws2YPV
824
825 // Visualization
826 G4VisAttributes* visAtt = new G4VisAttribute
827 visAtt -> SetVisibility(false);
828 visAtt -> SetForceSolid(false);
829
830 G4VisAttributes* simpleAlSVisAttPb = new G4V
831 simpleAlSVisAttPb -> SetVisibility(true);
832 simpleAlSVisAttPb -> SetForceSolid(true);
833
834 G4VisAttributes* simpleAlSVisAttXC10 = new G
835 simpleAlSVisAttXC10 -> SetVisibility(true);
836 simpleAlSVisAttXC10 -> SetForceSolid(true);
837
838 G4VisAttributes* simpleAlSVisAttWNICU = new
839 simpleAlSVisAttWNICU -> SetVisibility(true);
840 simpleAlSVisAttWNICU -> SetForceSolid(true);
841
842 box1LV1 -> SetVisAttributes(simpleAlSVisAttW
843 box2LV1 -> SetVisAttributes(simpleAlSVisAttP
844 box3LV1 -> SetVisAttributes(simpleAlSVisAttW
845 box4LV1 -> SetVisAttributes(simpleAlSVisAttP
846
847 box1LV2 -> SetVisAttributes(simpleAlSVisAttW
848 box2LV2 -> SetVisAttributes(simpleAlSVisAttP
849 box3LV2 -> SetVisAttributes(simpleAlSVisAttW
850 box4LV2 -> SetVisAttributes(simpleAlSVisAttP
851
852 boxJaw1YLV -> SetVisAttributes(visAtt);
853 boxJaw2YLV -> SetVisAttributes(visAtt);
854
855 // Region for cuts
856 G4Region* regVol = new G4Region("JawsYR");
857 G4ProductionCuts* cuts = new G4ProductionCut
858 cuts -> SetProductionCut(2.*cm);
859 regVol -> SetProductionCuts(cuts);
860
861 box1LV1 -> SetRegion(regVol);
862 regVol -> AddRootLogicalVolume(box1LV1);
863 box2LV1 -> SetRegion(regVol);
864 regVol -> AddRootLogicalVolume(box2LV1);
865 box3LV1 -> SetRegion(regVol);
866 regVol -> AddRootLogicalVolume(box3LV1);
867 box4LV1 -> SetRegion(regVol);
868 regVol -> AddRootLogicalVolume(box4LV1);
869
870 box1LV2 -> SetRegion(regVol);
871 regVol -> AddRootLogicalVolume(box1LV2);
872 box2LV2 -> SetRegion(regVol);
873 regVol -> AddRootLogicalVolume(box2LV2);
874 box3LV2 -> SetRegion(regVol);
875 regVol -> AddRootLogicalVolume(box3LV2);
876 box4LV2 -> SetRegion(regVol);
877 regVol -> AddRootLogicalVolume(box4LV2);
878 }
879
880 void DetectorConstruction::PhantomSegmentation
881 {
882
883 G4ThreeVector phantomHalfDim = {phantomSideD
884
885 nSideCells = CheckPhantomSegmentation(phanto
886 nDepthCells = CheckPhantomSegmentation(phant
887
888 G4Material* water = G4NistManager::Instance(
889
890 G4String yRepName("RepY");
891 G4VSolid* solYRep = new G4Box(yRepName, phan
892 phantomHalfDim.getY()/nSideCells, phanto
893 G4LogicalVolume* logYRep = new G4LogicalVolu
894 new G4PVReplica(yRepName, logYRep, Phantom_l
895 2.*phantomHalfDim.getY()/nSideCells);
896
897 G4String xRepName("RepX");
898 G4VSolid* solXRep = new G4Box(xRepName, phan
899 phantomHalfDim.getY()/nSideCells, phanto
900 G4LogicalVolume* logXRep = new G4LogicalVolu
901 new G4PVReplica(xRepName, logXRep, logYRep,
902 2.*phantomHalfDim.getX()/nSideCells);
903
904 G4String zVoxName("phantomSens");
905 G4VSolid* solVoxel = new G4Box(zVoxName, pha
906 phantomHalfDim.getY()/nSideCells, phanto
907
908 LVPhantomSens = new G4LogicalVolume(solVoxel
909 new G4PVReplica(zVoxName, LVPhantomSens, log
910 2.*phantomHalfDim.getZ()/nDepthCells);
911
912 }
913
914 G4int DetectorConstruction::CheckPhantomSegmen
915 {
916 // I want an odd number of voxels
917 // check whether nCells is integer
918 if ( nCells != (int)nCells )
919 {
920 nCells = std::ceil(nCells);
921 }
922 // check whether nCells is even
923 if ( std::fmod(nCells,2) == 0)
924 {
925 nCells += 1;
926 }
927
928 return nCells;
929 }
930
931 G4double DetectorConstruction::GetAccOriginPos
932 {
933 return sourceToSkinDistance;
934 }
935
936 G4double DetectorConstruction::GetFFilterRadiu
937 {
938 return tubeFFRadius;
939 }
940
941 G4double DetectorConstruction::GetFFilterZ()
942 {
943 return tubeFFFirstFaceZ;
944 }
945
946 G4double DetectorConstruction::GetVoxelDepthDi
947 {
948 return voxelDepthDim;
949 }
950
951 G4int DetectorConstruction::GetNumberSideCells
952 {
953 return nSideCells;
954 }
955
956 G4int DetectorConstruction::GetNumberDepthCell
957 {
958 return nDepthCells;
959 }
960
961 G4int DetectorConstruction::GetPhantomDepth()
962 {
963 return phantomSideDim;
964 }
965
966 void DetectorConstruction::SetJaws(G4double va
967 {
968 jawAperture = value;
969 }
970
971 void DetectorConstruction::SetTargetPosition(G
972 {
973 sourceToSkinDistance = value;
974 }
975
976 void DetectorConstruction::SetPhantomSide(G4do
977 {
978 phantomSideDim = value;
979 }
980
981 void DetectorConstruction::SetVoxelSide(G4doub
982 {
983 voxelSideDim = value;
984 }
985
986 void DetectorConstruction::SetVoxelDepth(G4dou
987 {
988 voxelDepthDim = value;
989 }
990
991 void DetectorConstruction::UpdateGeometry(G4St
992 {
993
994 if ( string == "fieldSide" )
995 {
996 //
997 G4double halfFieldSide = value/2;
998 G4ThreeVector position;
999 G4RotationMatrix* rMatrix1X = new G4Rotati
1000
1001 rMatrix1X -> rotateY(std::fabs(std::atan(
1002 position = *rMatrix1X * jaw1XInitialPos;
1003 boxJaw1X_phys -> SetTranslation(position)
1004 *rMatrix1X = rMatrix1X -> inverse();
1005 boxJaw1X_phys -> SetRotation(rMatrix1X);
1006
1007 G4RotationMatrix* rMatrix2X = new G4Rotat
1008 rMatrix2X -> rotateY(-std::fabs(std::atan
1009
1010 position = *rMatrix2X * jaw2XInitialPos;
1011 boxJaw2X_phys -> SetTranslation(position)
1012 *rMatrix2X = rMatrix2X -> inverse();
1013 boxJaw2X_phys -> SetRotation(rMatrix2X);
1014
1015 G4RotationMatrix* rMatrix1Y = new G4Rotat
1016 rMatrix1Y -> rotateX(-std::fabs(std::atan
1017 position = *rMatrix1Y * jaw1YInitialPos;
1018 boxJaw1Y_phys -> SetTranslation(position)
1019 *rMatrix1Y = rMatrix1Y -> inverse();
1020 boxJaw1Y_phys -> SetRotation(rMatrix1Y);
1021
1022 G4RotationMatrix* rMatrix2Y = new G4Rotat
1023 rMatrix2Y -> rotateX(std::fabs(std::atan(
1024
1025 position = *rMatrix2Y * jaw2YInitialPos;
1026 boxJaw2Y_phys -> SetTranslation(position)
1027 *rMatrix2Y = rMatrix2Y -> inverse();
1028 boxJaw2Y_phys -> SetRotation(rMatrix2Y);
1029 }
1030 else if ( string == "sourceToSkinDistance")
1031 {
1032 accWorld_phys -> SetTranslation({0., 0.,
1033 }
1034 else if (string == "phantomSide")
1035 {
1036 G4Box* box = (G4Box *) phantom_phys -> Ge
1037
1038 // change phantom side
1039 box -> SetXHalfLength(phantomSideDim/2.);
1040 box -> SetYHalfLength(phantomSideDim/2.);
1041 box -> SetZHalfLength(phantomSideDim/2.);
1042 phantom_phys -> SetTranslation({0., 0., p
1043
1044 // clear daughters and rebuild them
1045 phantom_phys -> GetLogicalVolume() -> Cle
1046 PhantomSegmentation(phantom_phys -> GetLo
1047 }
1048 else if (string == "voxelSide")
1049 {
1050 phantom_phys -> GetLogicalVolume() -> Cle
1051 PhantomSegmentation(phantom_phys -> GetLo
1052 }
1053 else if (string == "voxelDepth")
1054 {
1055 phantom_phys -> GetLogicalVolume() -> Cle
1056 PhantomSegmentation(phantom_phys -> GetLo
1057 }
1058 else
1059 G4cerr << "*** DetectorMessenger::UpdateG
1060
1061 G4RunManager::GetRunManager()->GeometryHasB
1062 }
1063
1064
1065