| 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 // Hadrontherapy advanced example for Geant4
27 // See more at: https://twiki.cern.ch/twiki/bi
28
29 #include "G4UnitsTable.hh"
30 #include "G4SDManager.hh"
31 #include "G4RunManager.hh"
32 #include "G4GeometryManager.hh"
33 #include "G4SolidStore.hh"
34 #include "G4PhysicalVolumeStore.hh"
35 #include "G4LogicalVolumeStore.hh"
36 #include "G4Box.hh"
37 #include "G4LogicalVolume.hh"
38 #include "G4ThreeVector.hh"
39 #include "G4PVPlacement.hh"
40 #include "globals.hh"
41 #include "G4Transform3D.hh"
42 #include "G4RotationMatrix.hh"
43 #include "G4Colour.hh"
44 #include "G4UserLimits.hh"
45 #include "G4UnitsTable.hh"
46 #include "G4VisAttributes.hh"
47 #include "G4NistManager.hh"
48 #include "HadrontherapyDetectorConstruction.hh
49 #include "HadrontherapyDetectorROGeometry.hh"
50 #include "HadrontherapyDetectorMessenger.hh"
51 #include "HadrontherapyDetectorSD.hh"
52 #include "HadrontherapyMatrix.hh"
53 #include "HadrontherapyLet.hh"
54 #include "PassiveProtonBeamLine.hh"
55 #include "BESTPassiveProtonBeamLine.hh"
56 #include "HadrontherapyMatrix.hh"
57
58 #include "HadrontherapyRBE.hh"
59 #include "G4SystemOfUnits.hh"
60
61 #include <cmath>
62
63
64
65 HadrontherapyDetectorConstruction* Hadronthera
66 //////////////////////////////////////////////
67 HadrontherapyDetectorConstruction::Hadronthera
68 : motherPhys(physicalTreatmentRoom), // pointe
69 detectorSD(0), detectorROGeometry(0), matrix(0
70 phantom(0), detector(0),
71 phantomLogicalVolume(0), detectorLogicalVolume
72 phantomPhysicalVolume(0), detectorPhysicalVolu
73 aRegion(0)
74 {
75
76
77 /* NOTE! that the HadrontherapyDetectorCon
78 * does NOT inherit from G4VUserDetectorCo
79 * So the Construct() mandatory virtual me
80 * like the passiveProtonBeamLIne, ...
81 */
82
83 // Messenger to change parameters of the p
84 detectorMessenger = new HadrontherapyDetec
85
86 // Default detector voxels size
87 // 200 slabs along the beam direction (X)
88 sizeOfVoxelAlongX = 200 *um;
89 sizeOfVoxelAlongY = 4 *cm;
90 sizeOfVoxelAlongZ = 4 *cm;
91
92 // Define here the material of the water p
93 SetPhantomMaterial("G4_WATER");
94 // Construct geometry (messenger commands)
95 // SetDetectorSize(4.*cm, 4.*cm, 4.*cm);
96 SetDetectorSize(4. *cm, 4. *cm, 4. *cm);
97 SetPhantomSize(40. *cm, 40. *cm, 40. *cm);
98
99 SetPhantomPosition(G4ThreeVector(20. *cm,
100 SetDetectorToPhantomPosition(G4ThreeVector
101 SetDetectorPosition();
102 //GetDetectorToWorldPosition();
103
104 // Write virtual parameters to the real on
105 UpdateGeometry();
106
107
108
109 }
110
111 //////////////////////////////////////////////
112 HadrontherapyDetectorConstruction::~Hadronther
113 {
114 delete detectorROGeometry;
115 delete matrix;
116 delete detectorMessenger;
117 }
118
119 //////////////////////////////////////////////
120 HadrontherapyDetectorConstruction* Hadronthera
121 {
122 return instance;
123 }
124
125 //////////////////////////////////////////////
126 // ConstructPhantom() is the method that const
127 // (or water phantom) in the usual Medical phy
128 // A water phantom can be considered a good ap
129 void HadrontherapyDetectorConstruction::Constr
130 {
131 // Definition of the solid volume of the P
132 phantom = new G4Box("Phantom",
133 phantomSizeX/2,
134 phantomSizeY/2,
135 phantomSizeZ/2);
136
137 // Definition of the logical volume of the
138 phantomLogicalVolume = new G4LogicalVolume
139
140
141
142 // Definition of the physics volume of the
143 phantomPhysicalVolume = new G4PVPlacement(
144
145
146
147
148
149
150
151 // Visualisation attributes of the phantom
152 red = new G4VisAttributes(G4Colour(255/255
153 red -> SetVisibility(true);
154 red -> SetForceSolid(true);
155 red -> SetForceWireframe(true);
156 phantomLogicalVolume -> SetVisAttributes(r
157 }
158
159 //////////////////////////////////////////////
160 // ConstructDetector() is the method the recon
161 // inside the water phantom. It is a volume, l
162 //
163 // **************************
164 // * water phantom *
165 // * *
166 // * *
167 // *--------------- *
168 // Beam * - *
169 // -----> * detector - *
170 // * - *
171 // *--------------- *
172 // * *
173 // * *
174 // * *
175 // **************************
176 //
177 // The detector can be dived in slices or voxe
178 // and inside it different quantities (dose di
179 // can be stored.
180 void HadrontherapyDetectorConstruction::Constr
181
182 {
183 // Definition of the solid volume of the D
184 detector = new G4Box("Detector",
185
186 phantomSizeX/2,
187
188 phantomSizeY/2,
189
190 phantomSizeZ/2);
191
192 // Definition of the logic volume of the P
193 detectorLogicalVolume = new G4LogicalVolum
194
195
196
197 // Definition of the physical volume of th
198 detectorPhysicalVolume = new G4PVPlacement
199
200
201
202
203
204
205 // Visualisation attributes of the detecto
206 skyBlue = new G4VisAttributes( G4Colour(13
207 skyBlue -> SetVisibility(true);
208 skyBlue -> SetForceSolid(true);
209 //skyBlue -> SetForceWireframe(true);
210 detectorLogicalVolume -> SetVisAttributes(
211
212 // **************
213 // **************
214 // Cut per Region
215 // **************
216 //
217 // A smaller cut is fixed in the phantom t
218 // required accuracy
219 if (!aRegion)
220 {
221 aRegion = new G4Region("DetectorLog");
222 detectorLogicalVolume -> SetRegion(aRe
223 aRegion->AddRootLogicalVolume( detecto
224 }
225 }
226
227 //////////////////////////////////////////////
228 void HadrontherapyDetectorConstruction::Initia
229
230
231
232 {
233 RO->Initialize(detectorToWorldPosition,
234 detectorSizeX/2,
235 detectorSizeY/2,
236 detectorSizeZ/2,
237 numberOfVoxelsAlongX,
238 numberOfVoxelsAlongY,
239 numberOfVoxelsAlongZ);
240 }
241 void HadrontherapyDetectorConstruction::Virtua
242 {
243
244 //Virtual plane
245 VirtualLayerPosition = G4ThreeVector(0*cm,
246 NewSource= Varbool;
247 if(NewSource == true)
248 {
249 // std::cout<<"trr"<<std::endl;
250 G4Material* airNist = G4NistManager::
251
252 solidVirtualLayer = new G4Box("Virtual
253 1.*um,
254 20.*cm,
255 40.*cm);
256
257 logicVirtualLayer = new G4LogicalVolum
258
259
260
261
262 physVirtualLayer= new G4PVPlacement(0,
263 "V
264 lo
265 mo
266 fa
267 0)
268
269 logicVirtualLayer -> SetVisAttributes(
270 }
271
272
273
274
275 }
276
277
278 //////////////////////////////////////////////
279 void HadrontherapyDetectorConstruction::Param
280 {
281 // Check phantom/detector sizes & relative
282 if (!IsInside(detectorSizeX,
283 detectorSizeY,
284 detectorSizeZ,
285 phantomSizeX,
286 phantomSizeY,
287 phantomSizeZ,
288 detectorToPhantomPosition
289 ))
290 G4Exception("HadrontherapyDetectorCons
291
292 // Check Detector sizes respect to the vox
293
294 if ( detectorSizeX < sizeOfVoxelAlongX) {
295 G4Exception("HadrontherapyDetectorCons
296 }
297 if ( detectorSizeY < sizeOfVoxelAlongY) {
298 G4Exception(" HadrontherapyDetectorCon
299 }
300 if ( detectorSizeZ < sizeOfVoxelAlongZ) {
301 G4Exception(" HadrontherapyDetectorCon
302 }
303 }
304
305 //////////////////////////////////////////////
306 G4bool HadrontherapyDetectorConstruction::SetP
307 {
308
309 if (G4Material* pMat = G4NistManager::Inst
310 {
311 phantomMaterial = pMat;
312 detectorMaterial = pMat;
313 if (detectorLogicalVolume && phantomLo
314 {
315 detectorLogicalVolume -> SetMateri
316 phantomLogicalVolume -> SetMateri
317
318 G4RunManager::GetRunManager() -> P
319 G4RunManager::GetRunManager() -> G
320 G4cout << "The material of Phantom
321 }
322 }
323 else
324 {
325 G4cout << "WARNING: material \"" << ma
326 " table [located in $G4INSTALL/source/
327 G4cout << "Use command \"/parameter/ni
328 return false;
329 }
330
331 return true;
332 }
333 //////////////////////////////////////////////
334 void HadrontherapyDetectorConstruction::SetPha
335 {
336 if (sizeX > 0.) phantomSizeX = sizeX;
337 if (sizeY > 0.) phantomSizeY = sizeY;
338 if (sizeZ > 0.) phantomSizeZ = sizeZ;
339 }
340
341 //////////////////////////////////////////////
342 void HadrontherapyDetectorConstruction::SetDet
343 {
344 if (sizeX > 0.) {detectorSizeX = sizeX;}
345 if (sizeY > 0.) {detectorSizeY = sizeY;}
346 if (sizeZ > 0.) {detectorSizeZ = sizeZ;}
347 SetVoxelSize(sizeOfVoxelAlongX, sizeOfVoxe
348 }
349
350 //////////////////////////////////////////////
351 void HadrontherapyDetectorConstruction::SetVox
352 {
353 if (sizeX > 0.) {sizeOfVoxelAlongX = sizeX
354 if (sizeY > 0.) {sizeOfVoxelAlongY = sizeY
355 if (sizeZ > 0.) {sizeOfVoxelAlongZ = sizeZ
356 }
357
358 //////////////////////////////////////////////
359 void HadrontherapyDetectorConstruction::SetPha
360 {
361 phantomPosition = pos;
362 }
363
364 //////////////////////////////////////////////
365 void HadrontherapyDetectorConstruction::SetDet
366 {
367 detectorToPhantomPosition = displ;
368 }
369
370 void HadrontherapyDetectorConstruction::SetVir
371 {
372
373 VirtualLayerPosition = position;
374 physVirtualLayer->SetTranslation(VirtualLa
375
376 }
377 //////////////////////////////////////////////
378 void HadrontherapyDetectorConstruction::Update
379 {
380 /*
381 * Check parameters consistency
382 */
383 ParametersCheck();
384
385 G4GeometryManager::GetInstance() -> OpenGe
386 if (phantom)
387 {
388 phantom -> SetXHalfLength(phantomSizeX
389 phantom -> SetYHalfLength(phantomSizeY
390 phantom -> SetZHalfLength(phantomSizeZ
391
392 phantomPhysicalVolume -> SetTranslatio
393 }
394 else ConstructPhantom();
395
396
397 // Get the center of the detector
398 SetDetectorPosition();
399 if (detector)
400 {
401
402 detector -> SetXHalfLength(detectorSiz
403 detector -> SetYHalfLength(detectorSiz
404 detector -> SetZHalfLength(detectorSiz
405
406 detectorPhysicalVolume -> SetTranslati
407 }
408 else ConstructDetector();
409
410 //std::cout<<NewSource<<std::endl;
411 /*if(NewSource)
412 {
413 std::cout<<"via"<<std::endl;
414 }*/
415
416
417 // std::cout<<"i"<<std::endl;
418 // std::cout<<VirtualLayerPosition<<std::e
419 // physVirtualLayer->SetTranslation(Virtua
420
421
422
423
424
425 // Round to nearest integer number of voxe
426
427 numberOfVoxelsAlongX = G4lrint(detectorSiz
428 sizeOfVoxelAlongX = ( detectorSizeX / numb
429 numberOfVoxelsAlongY = G4lrint(detectorSiz
430 sizeOfVoxelAlongY = ( detectorSizeY / numb
431 numberOfVoxelsAlongZ = G4lrint(detectorSiz
432 sizeOfVoxelAlongZ = ( detectorSizeZ / numb
433 PassiveProtonBeamLine *ppbl= (PassiveProto
434
435 G4RunManager::GetRunManager()->GetUserDete
436
437 HadrontherapyDetectorROGeometry* RO = (Had
438
439 //Set parameters, either for the Construct
440 RO->Initialize(GetDetectorToWorldPosition(
441 detectorSizeX/2,
442 detectorSizeY/2,
443 detectorSizeZ/2,
444 numberOfVoxelsAlongX,
445 numberOfVoxelsAlongY,
446 numberOfVoxelsAlongZ);
447
448 //This method below has an effect only if
449 RO->UpdateROGeometry();
450
451
452
453 volumeOfVoxel = sizeOfVoxelAlongX * sizeOf
454 massOfVoxel = detectorMaterial -> GetDensi
455 // This will clear the existing matrix (t
456 matrix = HadrontherapyMatrix::GetInstance(
457
458
459
460
461
462 // Initialize RBE
463 HadrontherapyRBE::CreateInstance(numberOfV
464
465 // Comment out the line below if let calcu
466 // Initialize LET with energy of primaries
467 if ( (let = HadrontherapyLet::GetInstance(
468 {
469 HadrontherapyLet::GetInstance() -> Ini
470 }
471
472
473 // Initialize analysis
474 // Inform the kernel about the new geometr
475 G4RunManager::GetRunManager() -> GeometryH
476 G4RunManager::GetRunManager() -> PhysicsHa
477
478 PrintParameters();
479
480 // CheckOverlaps();
481 }
482
483 //////////////////////////////////////////////
484 //Check of the geometry
485 //////////////////////////////////////////////
486 void HadrontherapyDetectorConstruction::CheckO
487 {
488 G4PhysicalVolumeStore* thePVStore = G4Phys
489 G4cout << thePVStore->size() << " physical
490 G4bool overlapFlag = false;
491 G4int res=1000;
492 G4double tol=0.; //tolerance
493 for (size_t i=0;i<thePVStore->size();i++)
494 {
495 //overlapFlag = (*thePVStore)[i]->Chec
496 overlapFlag = (*thePVStore)[i]->CheckO
497 if (overlapFlag)
498 G4cout << "Check: there are overlappin
499 }
500
501 //////////////////////////////////////////////
502 void HadrontherapyDetectorConstruction::PrintP
503 {
504
505 G4cout << "The (X,Y,Z) dimensions of the p
506 G4BestUnit( phantom -> GetXHalfLength()*2.
507 G4BestUnit( phantom -> GetYHalfLength()*2.
508 G4BestUnit( phantom -> GetZHalfLength()*2.
509
510 G4cout << "The (X,Y,Z) dimensions of the d
511 G4BestUnit( detector -> GetXHalfLength()*2
512 G4BestUnit( detector -> GetYHalfLength()*2
513 G4BestUnit( detector -> GetZHalfLength()*2
514
515 G4cout << "Displacement between Phantom an
516 G4cout << "DX= "<< G4BestUnit(phantomPosit
517 "DY= "<< G4BestUnit(phantomPosition.getY()
518 "DZ= "<< G4BestUnit(phantomPosition.getZ()
519
520 G4cout << "The (X,Y,Z) sizes of the Voxels
521 G4BestUnit(sizeOfVoxelAlongX, "Length") <
522 G4BestUnit(sizeOfVoxelAlongY, "Length") <
523 G4BestUnit(sizeOfVoxelAlongZ, "Length") <<
524
525 G4cout << "The number of Voxels along (X,Y
526 numberOfVoxelsAlongX << ',' <<
527 numberOfVoxelsAlongY <<',' <<
528 numberOfVoxelsAlongZ << ')' << G4endl;
529 }
530
531
532