| Geant4 Cross Reference |
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 //
28 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
30
31 #include "DetectorConstruction.hh"
32
33 #include "G4AutoDelete.hh"
34 #include "G4Box.hh"
35 #include "G4Colour.hh"
36 #include "G4Ellipsoid.hh"
37 #include "G4GeometryManager.hh"
38 #include "G4GlobalMagFieldMessenger.hh"
39 #include "G4LogicalVolume.hh"
40 #include "G4Material.hh"
41 #include "G4NistManager.hh"
42 #include "G4PVPlacement.hh"
43 #include "G4PhysicalConstants.hh"
44 #include "G4PhysicalVolumeStore.hh"
45 #include "G4RunManager.hh"
46 #include "G4Sphere.hh"
47 #include "G4SystemOfUnits.hh"
48 #include "G4UniformMagField.hh"
49 #include "G4UnitsTable.hh"
50 #include "G4UserLimits.hh"
51 #include "G4VisAttributes.hh"
52
53 #include "DetectorMessenger.hh"
54
55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56
57 DetectorConstruction::DetectorConstruction()
58 : G4VUserDetectorConstruction(),
59 fAbsorberMaterial(nullptr),
60 fWorldMaterial(nullptr),
61 fSolidWorld(nullptr),
62 fLogicWorld(nullptr),
63 fPhysiWorld(nullptr),
64 fSolidAbsorber(nullptr),
65 fLogicAbsorber(nullptr),
66 fPhysiAbsorber(nullptr),
67 material1(nullptr),
68 material2(nullptr),
69 material3(nullptr),
70 material4(nullptr),
71 material5(nullptr),
72 material6(nullptr),
73 material_GDP(nullptr),
74 ellipse1(nullptr),
75 logicEllipse1(nullptr),
76 physiEllipse1(nullptr),
77 ellipse2(nullptr),
78 logicEllipse2(nullptr),
79 physiEllipse2(nullptr),
80 ellipse3(nullptr),
81 logicEllipse3(nullptr),
82 physiEllipse3(nullptr),
83 ellipse4(nullptr),
84 logicEllipse4(nullptr),
85 physiEllipse4(nullptr),
86 ellipse5(nullptr),
87 logicEllipse5(nullptr),
88 physiEllipse5(nullptr),
89 ellipse6(nullptr),
90 logicEllipse6(nullptr),
91 physiEllipse6(nullptr),
92 solid_GDP(nullptr),
93 logic_GDP(nullptr),
94 physi_GDP(nullptr),
95 fDetectorMessenger(nullptr)
96 {
97 phantom_type = 1;
98
99 DefineMaterials();
100
101 if (phantom_type == 1) {
102 // default parameter values of the box
103 fAbsorberThickness = 40. * um;
104 fAbsorberSizeYZ = 40. * um;
105
106 // SetAbsorberMaterial("G4_Cu");
107 SetAbsorberMaterial("Body_10times");
108 // SetAbsorberMaterial("Body_real");
109 }
110 fXposAbs = 0. * cm;
111 ComputeGeomParameters();
112
113 // materials
114 SetWorldMaterial("G4_Galactic");
115
116 // create commands for interactive definition of the box
117 fDetectorMessenger = new DetectorMessenger(this);
118 }
119
120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
121
122 DetectorConstruction::~DetectorConstruction()
123 {
124 delete fDetectorMessenger;
125 }
126
127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
128
129 void DetectorConstruction::DefineMaterials()
130 {
131 //
132 // define Elements
133 //
134 auto H = new G4Element("Hydrogen", "H", 1, 1.01 * g / mole);
135 auto C = new G4Element("Carbon", "C", 6, 12.01 * g / mole);
136 auto N = new G4Element("Nitrogen", "N", 7, 14.01 * g / mole);
137 auto O = new G4Element("Oxygen", "O", 8, 16.00 * g / mole);
138 auto P = new G4Element("Phosphorus", "P", 15, 30.97 * g / mole);
139 auto S = new G4Element("Sulphur", "S", 16, 32.07 * g / mole);
140 auto Cl = new G4Element("Chlorine", "Cl", 17, 35.45 * g / mole);
141 auto K = new G4Element("Potassium", "K", 19, 39.10 * g / mole);
142 auto Ca = new G4Element("Calcium", "Ca", 20, 40.08 * g / mole);
143 auto Ti = new G4Element("Titanium", "Ti", 22, 47.87 * g / mole);
144 auto Ge = new G4Element("Germanium", "Ge", 32., 72.59 * g / mole);
145
146 //
147 // define biological materials
148 //
149 auto BioMaterial = new G4Material("C10H17O3N2", 1.218 * g / cm3, 4);
150 BioMaterial->AddElement(C, 10);
151 BioMaterial->AddElement(H, 17);
152 BioMaterial->AddElement(O, 3);
153 BioMaterial->AddElement(N, 2);
154
155 // FIRST ELLIPSOID : external surface of the worm ("skin")
156 // Original (real) composition in dry mass (the worm is dehydrated by
157 // lyophilization) Each element content is defined by its mass fraction (must
158 // be normalised to 1)
159 auto Skin_real = new G4Material("Skin_real", 0.49729 * g / cm3, 5);
160 Skin_real->AddElement(P, 0.47 * perCent);
161 Skin_real->AddElement(S, 0.21 * perCent);
162 Skin_real->AddElement(Cl, 0.05 * perCent);
163 Skin_real->AddElement(K, 0.48 * perCent);
164 Skin_real->AddMaterial(BioMaterial, 98.79 * perCent);
165
166 // Skin1: mass fractions of mineral elements 1O times more than original
167 // (real) mass fractions
168 auto Skin_10times = new G4Material("Skin_10times", 0.49729 * g / cm3, 5);
169 Skin_10times->AddElement(P, 4.71 * perCent);
170 Skin_10times->AddElement(S, 2.10 * perCent);
171 Skin_10times->AddElement(Cl, 0.46 * perCent);
172 Skin_10times->AddElement(K, 4.77 * perCent);
173 Skin_10times->AddMaterial(BioMaterial, 87.96 * perCent);
174
175 // SECOND ELLIPSOID : "body" of the worm, limited by the internal surface of
176 // the "skin" Original (real) composition in dry mass (the worm is dehydrated
177 // by lyophilization) Each element content is defined by its mass fraction
178 // (must be normalised to 1)
179 auto Body_real = new G4Material("Body_real", 0.40853 * g / cm3, 2);
180 Body_real->AddElement(P, 0.72 * perCent);
181 Body_real->AddMaterial(BioMaterial, 99.28 * perCent);
182
183 // Body_10times: mass fractions of mineral elements 1O times more than
184 // original (real) mass fractions
185 auto Body_10times = new G4Material("Body_10times", 0.40853 * g / cm3, 2);
186 Body_10times->AddElement(P, 7.23 * perCent);
187 Body_10times->AddMaterial(BioMaterial, 92.77 * perCent);
188
189 // THIRD ELLIPSOID : Nucleus of cell 1 (contained inside ellipsoid 2)
190 // Original (real) composition in dry mass (the worm is dehydrated by
191 // lyophilization) Each element content is defined by its mass fraction (must
192 // be normalised to 1)
193 auto Cell1_real = new G4Material("Cell1_real", 0.66262 * g / cm3, 3);
194 Cell1_real->AddElement(P, 0.57 * perCent);
195 Cell1_real->AddElement(Ca, 1.55 * perCent);
196 Cell1_real->AddMaterial(BioMaterial, 97.88 * perCent);
197
198 // Cell1_10times: mass fractions of mineral elements 1O times more than
199 // original (real) mass fractions
200 auto Cell1_10times = new G4Material("Cell1_10times", 0.66262 * g / cm3, 3);
201 Cell1_10times->AddElement(P, 5.66 * perCent);
202 Cell1_10times->AddElement(Ca, 15.49 * perCent);
203 Cell1_10times->AddMaterial(BioMaterial, 78.85 * perCent);
204
205 // FOURTH ELLIPSOID : Nucleus of cell 2 (contained inside ellipsoid 2)
206 // Original (real) composition in dry mass (the worm is dehydrated by
207 // lyophilization) Each element content is defined by its mass fraction (must
208 // be normalised to 1)
209 auto Cell2_real = new G4Material("Cell2_real", 0.60227 * g / cm3, 3);
210 Cell2_real->AddElement(P, 0.61 * perCent);
211 Cell2_real->AddElement(Ca, 1.44 * perCent);
212 Cell2_real->AddMaterial(BioMaterial, 97.95 * perCent);
213
214 // Cell2_10times: mass fractions of mineral elements 1O times more than
215 // original (real) mass fractions
216 auto Cell2_10times = new G4Material("Cell2_10times", 0.60227 * g / cm3, 3);
217 Cell2_10times->AddElement(P, 6.10 * perCent);
218 Cell2_10times->AddElement(Ca, 14.45 * perCent);
219 Cell2_10times->AddMaterial(BioMaterial, 79.45 * perCent);
220
221 // FIFTH ELLIPSOID: intestine of the worm
222 // Original (real) composition in dry mass (the worm is dehydrated by
223 // lyophilization) Each element content is defined by its mass fraction (must
224 // be normalised to 1)
225 auto Intestine_real = new G4Material("Intestine_real", 0.54058 * g / cm3, 4);
226 Intestine_real->AddElement(S, 0.12 * perCent);
227 Intestine_real->AddElement(Cl, 0.02 * perCent);
228 Intestine_real->AddElement(K, 0.20 * perCent);
229 Intestine_real->AddMaterial(BioMaterial, 99.66 * perCent);
230
231 // Intestine_10times: mass fractions of mineral elements 1O times more than
232 // original (real) mass fractions
233 auto Intestine_10times = new G4Material("Intestine_10times", 0.54058 * g / cm3, 4);
234 Intestine_10times->AddElement(S, 1.17 * perCent);
235 Intestine_10times->AddElement(Cl, 0.23 * perCent);
236 Intestine_10times->AddElement(K, 1.99 * perCent);
237 Intestine_10times->AddMaterial(BioMaterial, 96.61 * perCent);
238
239 // SIXTH ELLIPSOID: Nucleus of cell Ti (contained inside ellipsoid 5)
240 // Original (real) composition in dry mass (the worm is dehydrated by
241 // lyophilization) Each element content is defined by its mass fraction (must
242 // be normalised to 1)
243 auto CellTi_real = new G4Material("CellTi_real", 0.75138 * g / cm3, 5);
244 CellTi_real->AddElement(S, 0.11 * perCent);
245 CellTi_real->AddElement(Cl, 0.02 * perCent);
246 CellTi_real->AddElement(K, 0.18 * perCent);
247 CellTi_real->AddElement(Ti, 0.05 * perCent);
248 CellTi_real->AddMaterial(BioMaterial, 99.64 * perCent);
249
250 // CellTi_200times: mass fractions of mineral elements 1O times more than
251 // original (real) mass fractions except for Ti which is multiplied by 200
252 // times
253 auto CellTi_200times = new G4Material("CellTi_200times", 0.75138 * g / cm3, 5);
254 CellTi_200times->AddElement(S, 1.05 * perCent);
255 CellTi_200times->AddElement(Cl, 0.22 * perCent);
256 CellTi_200times->AddElement(K, 1.79 * perCent);
257 CellTi_200times->AddElement(Ti, 10.89 * perCent); // Ti density 200 times
258 CellTi_200times->AddMaterial(BioMaterial, 86.05 * perCent);
259
260 // CellTi_1000times: mass fractions of mineral elements 10 times more than
261 // original (real) mass fractions except for Ti which is multiplied by 1000
262 // times
263 auto CellTi_1000times = new G4Material("CellTi_1000times", 0.75138 * g / cm3, 5);
264 CellTi_1000times->AddElement(S, 1.05 * perCent);
265 CellTi_1000times->AddElement(Cl, 0.22 * perCent);
266 CellTi_1000times->AddElement(K, 1.79 * perCent);
267 CellTi_1000times->AddElement(Ti, 54.47 * perCent); // Ti density 1000 times
268 CellTi_1000times->AddMaterial(BioMaterial, 42.47 * perCent);
269
270 //
271 // define simple materials Cl-doped polystyrene capsules (PS)
272 //
273 // Polystyrene
274 auto Cl_Polystyrene = new G4Material("Cl-doped Polystyrene", 1.18425 * g / cm3, 3);
275 Cl_Polystyrene->AddElement(C, 97);
276 Cl_Polystyrene->AddElement(H, 97);
277 Cl_Polystyrene->AddElement(Cl, 6);
278
279 // Ge-doped glow discharge polymer (GDP)
280 auto GDP = new G4Material("GDP", 1.08 * g / cm3, 3);
281 GDP->AddElement(C, 76430);
282 GDP->AddElement(H, 107002);
283 GDP->AddElement(Ge, 744);
284
285 //
286 // define a material from elements.
287 //
288 auto H2O = new G4Material("Water", 1.000 * g / cm3, 2);
289 H2O->AddElement(H, 2);
290 H2O->AddElement(O, 1);
291 H2O->GetIonisation()->SetMeanExcitationEnergy(78 * eV);
292
293 auto Air = new G4Material("Air", 1.290 * mg / cm3, 2);
294 Air->AddElement(N, 0.7);
295 Air->AddElement(O, 0.3);
296
297 //
298 // example of vacuum
299 //
300 G4double density = universe_mean_density; // from PhysicalConstants.h
301 G4double pressure = 3.e-18 * pascal;
302 G4double temperature = 2.73 * kelvin;
303 new G4Material("Galactic", 1, 1.01 * g / mole, density, kStateGas, temperature, pressure);
304 }
305
306 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
307
308 void DetectorConstruction::ComputeGeomParameters()
309 {
310 if (phantom_type == 1) {
311 // Compute derived parameters of the box
312 fXstartAbs = fXposAbs - 0.5 * fAbsorberThickness;
313 fXendAbs = fXposAbs + 0.5 * fAbsorberThickness;
314
315 G4double xmax = std::max(std::abs(fXstartAbs), std::abs(fXendAbs));
316 fWorldSizeX = 4 * xmax;
317 fWorldSizeYZ = 2 * fAbsorberSizeYZ;
318
319 if (nullptr != fPhysiWorld) {
320 ChangeGeometry();
321 }
322 }
323 else if (phantom_type == 2) {
324 fWorldSizeX = 1 * mm;
325 fWorldSizeYZ = 1 * mm;
326 }
327 else if (phantom_type == 3) {
328 fWorldSizeX = 1.5 * mm;
329 fWorldSizeYZ = 1.5 * mm;
330 }
331 }
332
333 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
334
335 G4VPhysicalVolume* DetectorConstruction::Construct()
336 {
337 if (nullptr != fPhysiWorld) {
338 return fPhysiWorld;
339 }
340 // World
341 //
342 fSolidWorld = new G4Box("World", // its name
343 fWorldSizeX / 2, fWorldSizeYZ / 2, fWorldSizeYZ / 2); // its size
344
345 fLogicWorld = new G4LogicalVolume(fSolidWorld, // its solid
346 fWorldMaterial, // its material
347 "World"); // its name
348
349 fPhysiWorld = new G4PVPlacement(nullptr, // no rotation
350 G4ThreeVector(0., 0., 0.), // at (0,0,0)
351 fLogicWorld, // its logical volume
352 "World", // its name
353 nullptr, // its mother volume
354 false, // no boolean operation
355 0); // copy number
356 if (phantom_type == 1) {
357 Construct_Phantom1();
358 }
359 else if (phantom_type == 2) {
360 Construct_Phantom2();
361 }
362 else if (phantom_type == 3) {
363 Construct_Phantom3();
364 }
365
366 // fLogicWorld->SetVisAttributes (G4VisAttributes::GetInvisible());
367
368 // always return the physical World
369 //
370 return fPhysiWorld;
371 }
372
373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
374
375 void DetectorConstruction::PrintGeomParameters()
376 {
377 if (phantom_type == 1) {
378 G4cout << "\n" << fWorldMaterial << G4endl;
379 G4cout << "\n" << fAbsorberMaterial << G4endl;
380
381 G4cout << "\n The WORLD is made of " << G4BestUnit(fWorldSizeX, "Length") << " of "
382 << fWorldMaterial->GetName();
383 G4cout << ". The transverse size (YZ) of the world is " << G4BestUnit(fWorldSizeYZ, "Length")
384 << G4endl;
385 G4cout << " The ABSORBER is made of " << G4BestUnit(fAbsorberThickness, "Length") << " of "
386 << fAbsorberMaterial->GetName();
387 G4cout << ". The transverse size (YZ) is " << G4BestUnit(fAbsorberSizeYZ, "Length") << G4endl;
388 G4cout << " X position of the middle of the absorber " << G4BestUnit(fXposAbs, "Length");
389 G4cout << G4endl;
390 }
391 else if (phantom_type == 2) {
392 G4cout << "The upper part of C.elegans is made of 6 ellipsoids" << G4endl;
393 }
394 else if (phantom_type == 3) {
395 G4cout << "A microsphere inertial confinement fusion target is contructed, made of: " << G4endl;
396 G4cout << "\n" << material_GDP << G4endl;
397 }
398 }
399
400 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
401
402 void DetectorConstruction::SetAbsorberMaterial(const G4String& materialChoice)
403 {
404 // search the material by its name
405 G4Material* pttoMaterial = G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
406
407 if (pttoMaterial && fAbsorberMaterial != pttoMaterial) {
408 fAbsorberMaterial = pttoMaterial;
409 if (fLogicAbsorber) {
410 fLogicAbsorber->SetMaterial(fAbsorberMaterial);
411 }
412 G4RunManager::GetRunManager()->PhysicsHasBeenModified();
413 }
414 }
415
416 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
417
418 void DetectorConstruction::SetWorldMaterial(const G4String& materialChoice)
419 {
420 // search the material by its name
421 G4Material* pttoMaterial = G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
422
423 if (pttoMaterial && fWorldMaterial != pttoMaterial) {
424 fWorldMaterial = pttoMaterial;
425 if (fLogicWorld) {
426 fLogicWorld->SetMaterial(fWorldMaterial);
427 }
428 G4RunManager::GetRunManager()->PhysicsHasBeenModified();
429 }
430 }
431
432 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
433
434 void DetectorConstruction::SetAbsorberThickness(G4double val)
435 {
436 fAbsorberThickness = val;
437 ComputeGeomParameters();
438 }
439
440 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
441
442 void DetectorConstruction::SetAbsorberSizeYZ(G4double val)
443 {
444 fAbsorberSizeYZ = val;
445 ComputeGeomParameters();
446 }
447
448 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
449
450 void DetectorConstruction::SetWorldSizeX(G4double val)
451 {
452 fWorldSizeX = val;
453 ComputeGeomParameters();
454 }
455
456 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
457
458 void DetectorConstruction::SetWorldSizeYZ(G4double val)
459 {
460 fWorldSizeYZ = val;
461 ComputeGeomParameters();
462 }
463
464 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
465
466 void DetectorConstruction::SetAbsorberXpos(G4double val)
467 {
468 fXposAbs = val;
469 }
470
471 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
472
473 void DetectorConstruction::ConstructSDandField()
474 {
475 if (fFieldMessenger.Get() == nullptr) {
476 // Create global magnetic field messenger.
477 // Uniform magnetic field is then created automatically if
478 // the field value is not zero.
479 G4ThreeVector fieldValue = G4ThreeVector();
480 auto msg = new G4GlobalMagFieldMessenger(fieldValue);
481 // msg->SetVerboseLevel(1);
482 G4AutoDelete::Register(msg);
483 fFieldMessenger.Put(msg);
484 }
485 }
486
487 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
488
489 void DetectorConstruction::ChangeGeometry()
490 {
491 if (phantom_type == 1) {
492 fSolidWorld->SetXHalfLength(fWorldSizeX * 0.5);
493 fSolidWorld->SetYHalfLength(fWorldSizeYZ * 0.5);
494 fSolidWorld->SetZHalfLength(fWorldSizeYZ * 0.5);
495
496 fSolidAbsorber->SetXHalfLength(fAbsorberThickness * 0.5);
497 fSolidAbsorber->SetYHalfLength(fAbsorberSizeYZ * 0.5);
498 fSolidAbsorber->SetZHalfLength(fAbsorberSizeYZ * 0.5);
499 }
500 }
501
502 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
503
504 void DetectorConstruction::SetPhantomType(G4int value)
505 {
506 phantom_type = value;
507 ComputeGeomParameters();
508 }
509
510 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
511
512 void DetectorConstruction::Construct_Phantom1()
513 {
514 fSolidAbsorber =
515 new G4Box("Absorber", fAbsorberThickness / 2, fAbsorberSizeYZ / 2, fAbsorberSizeYZ / 2);
516
517 fLogicAbsorber = new G4LogicalVolume(fSolidAbsorber, // its solid
518 fAbsorberMaterial, // its material
519 "Absorber"); // its name
520
521 G4RotationMatrix rm;
522 G4double angle = 0 * deg;
523 rm.rotateZ(angle);
524
525 fPhysiAbsorber = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(fXposAbs * um, 0., 0.)),
526 fLogicAbsorber, // its logical volume
527 "Absorber", // its name
528 fLogicWorld, // its mother
529 false, // no boolean operation
530 0); // copy number
531
532 auto logVisAtt = new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745));
533 logVisAtt->SetForceSolid(true);
534 fLogicAbsorber->SetVisAttributes(logVisAtt);
535
536 PrintGeomParameters();
537 }
538
539 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
540
541 void DetectorConstruction::Construct_Phantom2()
542 {
543 G4NistManager* nist = G4NistManager::Instance();
544 // PHANTOM 3: C. ELEGANS WORM **********************************
545 // made of 6 ellipsoids
546
547 // FIRST ELLIPSOID : external surface of the worm ("skin")
548
549 material1 = nist->FindOrBuildMaterial("Skin_10times");
550 G4ThreeVector center1 =
551 G4ThreeVector(fXposAbs, 0,
552 0); // position of the center of the first ellipsoid is (fXposAbs,0,0)
553 // if we change fXposAbs, all the phantom will be translated along X axis
554 // as the position of the center of all other ellipsoids is defined according
555 // to fXposAbs
556 G4double halfx1 = 20.61 * um; // Semiaxis in X
557 G4double halfy1 = 21.42 * um; // Semiaxis in Y
558 G4double halfz1 = 187.82 * um; // Semiaxis in Z
559 G4double pzTopCut1 = 187.82 * um; // here the top of the ellipsoid is not cut away
560 G4double pzBottomCut1 = 0 * um; // here the bottom of the ellipsoid is cut
561 // away below z = 0 (original configuration)
562 // G4double pzBottomCut1 = 18.07*um; // test with position of the slice of
563 // interest for vis.mac visualisation of the slice z = 18.07 µm
564
565 ellipse1 = new G4Ellipsoid("Ellipse1", halfx1, halfy1, halfz1, pzBottomCut1,
566 pzTopCut1); // cut ellipsoid (original configuration)
567
568 logicEllipse1 = new G4LogicalVolume(ellipse1, // its solid
569 material1, // its material
570 "Ellipse1"); // its name
571
572 physiEllipse1 = new G4PVPlacement(nullptr, // no rotation
573 center1, // its position
574 logicEllipse1, // its logical volume
575 "Ellipse1", // its name
576 fLogicWorld, // its mother
577 false, // no boolean operation
578 0); // copy number
579
580 // SECOND ELLIPSOID : "body" of the worm, limited by the internal surface of
581 // the "skin"
582 material2 = nist->FindOrBuildMaterial("Body_10times");
583 G4ThreeVector center2 = G4ThreeVector(-0.39 * um, 0, 0); // position according to ellipsoid 1
584 G4double halfx2 = 18.61 * um; // Semiaxis in X
585 G4double halfy2 = 19.01 * um; // Semiaxis in Y
586 G4double halfz2 = 186.64 * um; // Semiaxis in Z
587
588 G4double pzTopCut2 = 186.64 * um; // here the top of the ellipsoid is not cut away
589 G4double pzBottomCut2 = 0 * um; // here the bottom of the ellipsoid is cut
590 // away below z = 0 (original configuration)
591 // G4double pzBottomCut2 = 18.07*um; // test with position of the slice of
592 // interest for vis.mac visualisation of the slice z = 18.07 µm
593
594 // Ellipse2 = new G4Ellipsoid("Ellipse2",halfx2,halfy2,halfz2); // test with
595 // entire ellipsoid
596 ellipse2 = new G4Ellipsoid("Ellipse2", halfx2, halfy2, halfz2, pzBottomCut2,
597 pzTopCut2); // cut ellipsoid
598
599 logicEllipse2 = new G4LogicalVolume(ellipse2, // its solid
600 material2, // its material
601 "Ellipse2"); // its name
602
603 physiEllipse2 = new G4PVPlacement(nullptr, // no rotation
604 center2, // its position
605 logicEllipse2, // its logical volume
606 "Ellipse2", // its name
607 logicEllipse1, // its mother
608 false, // no boolean operation
609 0); // copy number
610
611 // THIRD ELLIPSOID : Nucleus of cell 1 (contained inside ellipsoid 2)
612 material3 = nist->FindOrBuildMaterial("Cell1_10times");
613 G4ThreeVector center3 =
614 G4ThreeVector(2.36 * um, -7.09 * um, 18.07 * um); // position according to ellipsoid 2
615 G4double halfx3 = 1.95 * um; // Semiaxis in X
616 G4double halfy3 = 3.23 * um; // Semiaxis in Y
617 G4double halfz3 = 4.32 * um; // Semiaxis in Z
618
619 // if we want to visualize slice at z = 18.07 µm using vis.mac,
620 // we should cut all the ellipsoids (cut away all parts below 18.07 µm)
621 // Visualization test: the lower part of ellipse 3 is cut off
622 // G4double pzTopCut3 =4.32 *um; // Visualization test
623 // G4double pzBottomCut3 = 0*um; // Visualization test: cut at zero before
624 // moving ellipsoid 3
625 // // to position defined just above => it will be cut at 18.07 µm after
626 // moving
627
628 ellipse3 = new G4Ellipsoid("Ellipse3", halfx3, halfy3, halfz3); // entire ellipsoid
629 // ellipse3 = new
630 // G4Ellipsoid("Ellipse3",halfx3,halfy3,halfz3,pzBottomCut3,pzTopCut3); //
631 // Visualization test
632
633 logicEllipse3 = new G4LogicalVolume(ellipse3, // its solid
634 material3, // its material
635 "Ellipse3"); // its name
636
637 physiEllipse3 = new G4PVPlacement(nullptr, // no rotation
638 center3, // its position
639 logicEllipse3, // its logical volume
640 "Ellipse3", // its name
641 logicEllipse2, // its mother
642 false, // no boolean operation
643 0); // copy number
644
645 // FOURTH ELLIPSOID : Nucleus of cell 2 (contained inside ellipsoid 2)
646 material4 = nist->FindOrBuildMaterial("Cell2_10times");
647 G4ThreeVector center4 =
648 G4ThreeVector(8.66 * um, -3.15 * um, 18.07 * um); // position according to ellipsoid 2
649 G4double halfx4 = 2.08 * um; // Semiaxis in X
650 G4double halfy4 = 2.46 * um; // Semiaxis in Y
651 G4double halfz4 = 4.32 * um; // Semiaxis in Z
652
653 // means the lower part of ellipse 4 was cut off
654 // G4double pzTopCut4 = 4.32*um;
655 // G4double pzBottomCut4 = 0*um;// Visualization test: cut at zero before
656 // moving ellipsoid 4
657 // // to position defined just above => it will be cut at 18.07 µm after
658 // moving
659
660 ellipse4 = new G4Ellipsoid("Ellipse4", halfx4, halfy4, halfz4); // entire ellipsoid
661 // Ellipse4 = new
662 // G4Ellipsoid("Ellipse4",halfx4,halfy4,halfz4,pzBottomCut4,pzTopCut4); //
663 // Visualization test
664
665 logicEllipse4 = new G4LogicalVolume(ellipse4, // its solid
666 material4, // its material
667 "Ellipse4"); // its name
668
669 physiEllipse4 = new G4PVPlacement(nullptr, // no rotation
670 center4, // its position
671 logicEllipse4, // its logical volume
672 "Ellipse4", // its name
673 logicEllipse2, // its mother
674 false, // no boolean operation
675 0); // copy number
676
677 // FIFTH ELLIPSOID: intestine of the worm
678 // Original (real) composition in dry mass (the worm is dehydrated by
679 // lyophilization) Each element content is defined by its mass fraction (must
680 // be normalised to 1)
681 material5 = nist->FindOrBuildMaterial("Intestine_10times");
682 G4ThreeVector center5 =
683 G4ThreeVector(1.64 * um, 0.61 * um, 0 * um); // position according to ellipsoid 2
684 G4RotationMatrix rm5;
685 G4double angle5 = -58.986 * deg; // ellipse 5 is rotated around z axis
686 // (rotation defined according to ellipse 2)
687 // G4double angle5=0*deg;
688 rm5.rotateZ(angle5);
689
690 G4double halfx5 = 3.67 * um; // Semiaxis in X
691 G4double halfy5 = 16.48 * um; // Semiaxis in Y
692 G4double halfz5 = 28.68 * um; // Semiaxis in Z
693
694 G4double pzTopCut5 = 28.68 * um; // here the top of the ellipsoid is not cut away
695 G4double pzBottomCut5 = 0 * um; // here the bottom of the ellipsoid is cut
696 // away below z = 0 // original configuration
697 // G4double pzBottomCut5 =18.07*um; // position of the slice of interest //
698 // Visualization test
699
700 // ellipse5 = new G4Ellipsoid("Ellipse5",halfx5,halfy5,halfz5); // entire
701 // ellipsoid // Visualization test
702 ellipse5 = new G4Ellipsoid("Ellipse5", halfx5, halfy5, halfz5, pzBottomCut5,
703 pzTopCut5); // original configuration
704
705 logicEllipse5 = new G4LogicalVolume(ellipse5, // its solid
706 material5, // its material
707 "Ellipse5"); // its name
708
709 physiEllipse5 = new G4PVPlacement(G4Transform3D(rm5, center5), // rotation
710 logicEllipse5, // its logical volume
711 "Ellipse5", // its name
712 logicEllipse2, // its mother
713 false, // no boolean operation
714 0); // copy number
715
716 // SIXTH ELLIPSOID: Nucleus of cell Ti (contained inside ellipsoid 5)
717 // Original (real) composition in dry mass (the worm is dehydrated by
718 // lyophilization) Each element content is defined by its mass fraction (must
719 // be normalised to 1)
720
721 material6 = nist->FindOrBuildMaterial("CellTi_1000times");
722 G4ThreeVector center6 =
723 G4ThreeVector(0.005 * um, 5.83 * um, 18.07 * um); // position according to ellipsoid 5
724 G4double halfx6 = 1.62 * um; // Semiaxis in X
725 G4double halfy6 = 1.95 * um; // Semiaxis in Y
726 G4double halfz6 = 1.62 * um; // Semiaxis in Z
727
728 // means the lower part of ellipse 6 was cut off
729 // G4double pzTopCut6 = 1.62*um;
730 // G4double pzBottomCut6 = 0*um; // Visualization test: cut at zero before
731 // moving ellipsoid 6
732 // // to position defined just above => it will be cut at 18.07 µm after
733 // moving. As ellipsoid 6 is contained in ellipsoid 5, it has been rotated at
734 // the same time as ellipsoid 5. We want to counterbalance this rotation and
735 // put ellipsoid 6 in the right direction (not rotated at all). So we must
736 // rotate ellipsoid 6 of the opposite angle as ellipsoid 5
737 G4RotationMatrix rm6;
738 G4double angle6 = 58.986 * deg; // ellipsoid 6 is rotated around z axis
739 // (rotation defined according to ellipse 5)
740 // G4double angle6= 0*deg;
741 rm6.rotateZ(angle6);
742 ellipse6 = new G4Ellipsoid("Ellipse6", halfx6, halfy6, halfz6); // entire ellipsoid
743 // ellipse6 = new
744 // G4Ellipsoid("Ellipse6",halfx6,halfy6,halfz6,pzBottomCut6,pzTopCut6); //
745 // Visualization test
746 logicEllipse6 = new G4LogicalVolume(ellipse6, // its solid
747 material6, // its material
748 "Ellipse6"); // its name
749
750 physiEllipse6 =
751 new G4PVPlacement(G4Transform3D(rm6, center6), // So now ellipsoid 6 is vertical again
752 // (appears not rotated at all)
753 logicEllipse6, // its logical volume
754 "Ellipse6", // its name
755 logicEllipse5, // its mother
756 false, // no boolean operation
757 0); // copy number
758
759 /* logicEllipse1->SetUserLimits(new G4UserLimits(0.5*um,DBL_MAX,DBL_MAX,0));
760
761 logicEllipse2->SetUserLimits(new G4UserLimits(0.5*um,DBL_MAX,DBL_MAX,0)); */
762
763 // As the material composing the C. elegans worm is very transparent to
764 // photons, the default steps in Geant4 (distance between pre-step and
765 // post-step point) are too long. To perform a more precise simulation,
766 // maximal values are defined for step length in each volume, according to its
767 // size and configuration.
768
769 logicEllipse1->SetUserLimits(new G4UserLimits(0.2 * um));
770 logicEllipse2->SetUserLimits(new G4UserLimits(0.3 * um));
771 logicEllipse3->SetUserLimits(new G4UserLimits(0.4 * um));
772 logicEllipse4->SetUserLimits(new G4UserLimits(0.4 * um));
773 logicEllipse5->SetUserLimits(new G4UserLimits(0.7 * um));
774 logicEllipse6->SetUserLimits(new G4UserLimits(0.3 * um));
775 // fLogicWorld->SetUserLimits(new G4UserLimits(1*um)); // world does not have
776 // to be simulated as precisely
777
778 // fLogicWorld ->SetVisAttributes(new
779 // G4VisAttributes(G4Colour(1.0,1.0,1.0)));
780 logicEllipse1->SetVisAttributes(new G4VisAttributes(G4Colour(1.0, 0.8, 1.0, 0.5)));
781 logicEllipse2->SetVisAttributes(new G4VisAttributes(G4Colour(1.0, 0.5, 1.0, 0.5)));
782 // logicEllipse3 ->SetVisAttributes(new
783 // G4VisAttributes(G4Colour(0.0,1.0,1.0))); logicEllipse4
784 // ->SetVisAttributes(new G4VisAttributes(G4Colour(0.5,0.3,0.5)));
785 logicEllipse5->SetVisAttributes(new G4VisAttributes(G4Colour(0.5, 0.5, 1.0, 0.5)));
786 // logicEllipse6 ->SetVisAttributes(new
787 // G4VisAttributes(G4Colour(0.5,0.8,0.2)));
788
789 auto logVisAtt3 = new G4VisAttributes(G4Colour(0.0, 1.0, 1.0));
790 logVisAtt3->SetForceSolid(true);
791 logicEllipse3->SetVisAttributes(logVisAtt3);
792
793 auto logVisAtt4 = new G4VisAttributes(G4Colour(0.5, 0.3, 0.5));
794 logVisAtt4->SetForceSolid(true);
795 logicEllipse4->SetVisAttributes(logVisAtt4);
796
797 auto logVisAtt6 = new G4VisAttributes(G4Colour(0.5, 0.8, 0.2));
798 logVisAtt6->SetForceSolid(true);
799 logicEllipse6->SetVisAttributes(logVisAtt6);
800
801 PrintGeomParameters();
802 }
803
804 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
805
806 void DetectorConstruction::Construct_Phantom3()
807 {
808 G4NistManager* nist = G4NistManager::Instance();
809 material_GDP = nist->FindOrBuildMaterial("GDP");
810
811 G4double pRmin = 171 * um; // thickness = 25 um
812 G4double pRmax = 196 * um;
813
814 // G4double pRmin = 146.8*um; // thickness = 10 um
815 // G4double pRmax = 156.8*um;
816
817 // G4double pRmin = 293.6*um; // thickness = 20 um
818 // G4double pRmax = 313.6*um;
819
820 solid_GDP = new G4Sphere("GDP", // name
821 pRmin, pRmax, 0., 2 * pi, 0., pi);
822
823 logic_GDP = new G4LogicalVolume(solid_GDP, // its solid
824 material_GDP, // its material
825 "GDP"); // its name
826
827 physi_GDP = new G4PVPlacement(nullptr, // no rotation
828 G4ThreeVector(fXposAbs, 0., 0.), // its position
829 logic_GDP, // its logical volume
830 "GDP", // its name
831 fLogicWorld, // its mother
832 false, // no boolean operation
833 0); // copy number
834
835 logic_GDP->SetUserLimits(new G4UserLimits(1 * um));
836
837 auto logVisAtt = new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745, 0.5));
838 // new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745));
839 logVisAtt->SetForceSolid(true);
840 // logVisAtt->SetVisibility(true);
841
842 logic_GDP->SetVisAttributes(logVisAtt);
843
844 PrintGeomParameters();
845 }
846
847 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
848