Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
HadrontherapyModulator.cc
Go to the documentation of this file.
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 // This is the *BASIC* version of Hadrontherapy, a Geant4-based application
27 // See more at: http://g4advancedexamples.lngs.infn.it/Examples/hadrontherapy
28 //
29 // Visit the Hadrontherapy web site (http://www.lns.infn.it/link/Hadrontherapy) to request
30 // the *COMPLETE* version of this program, together with its documentation;
31 // Hadrontherapy (both basic and full version) are supported by the Italian INFN
32 // Institute in the framework of the MC-INFN Group
33 //
34 
35 #include <fstream>
36 
37 #include "globals.hh"
38 #include "G4SystemOfUnits.hh"
39 #include "G4Material.hh"
40 #include "G4Tubs.hh"
41 #include "G4Box.hh"
42 #include "G4LogicalVolume.hh"
43 #include "G4VPhysicalVolume.hh"
44 #include "G4ThreeVector.hh"
45 #include "G4PVPlacement.hh"
46 #include "G4Transform3D.hh"
47 #include "G4RotationMatrix.hh"
48 #include "G4VisAttributes.hh"
49 #include "G4Colour.hh"
51 #include "G4Transform3D.hh"
52 #include "G4ios.hh"
53 #include "G4RunManager.hh"
54 #include "G4NistManager.hh"
55 #include <iostream>
56 
57 using namespace std;
58 
60  solidMod1(0), logicMod1(0), physiMod1(0),
61  solidMod2(0), logicMod2(0), physiMod2(0),
62  solidMod3(0), logicMod3(0), physiMod3(0),
63  solidMod4(0), logicMod4(0), physiMod4(0),
64  FileName("ModoulWeight.txt")
65 {
66  pi=4*atan(1.);
67  StepNumbers=22;
68  Weight=new G4double[StepNumbers];
69  StepThickness=new G4double[StepNumbers];
70  StartingAngle=new G4double[StepNumbers];
71  SpanningAngle=new G4double[StepNumbers];
72  PositionMod=new G4ThreeVector[StepNumbers];
73 
74 
75  solidMod=new G4Tubs *[StepNumbers];
76  logicMod=new G4LogicalVolume *[StepNumbers];
77  physiMod=new G4VPhysicalVolume *[(4*(StepNumbers-1)+1)];
78 
79  for (G4int i=0;i<StepNumbers;i++)
80  {
81  Weight[i]=0;
82  StepThickness[i]=0;
83  StartingAngle[i]=0;
84  SpanningAngle[i]=0;
85  PositionMod[i]=G4ThreeVector(0,0,0);
86  solidMod[i]=0;
87  logicMod[i]=0;
88 
89  }
90 
91  for (G4int i=0;i<4*(StepNumbers-1)+1;i++)
92  {
93  physiMod[i]=0;
94  }
95 
96 
97  ModulatorMessenger = new HadrontherapyModulatorMessenger(this);
99  rm = new G4RotationMatrix();
100  G4double phi = 270. *deg;
101  rm -> rotateY(phi);
102 }
105 {
106  delete rm;
107  delete [] Weight;
108  delete [] StepThickness;
109  delete [] StartingAngle;
110  delete [] SpanningAngle;
111  delete [] PositionMod;
112  delete [] solidMod;
113  delete [] logicMod;
114  delete [] physiMod;
115  delete ModulatorMessenger;
116 }
117 
120 {
121 /* Here we initialize the step properties of Modulator wheel, you can create your
122 specific modulator by changing the values in this class or writing them in an external
123 file and activate reading from file via a macrofile. */
124 
125  StepThickness[0]=0; Weight[0]=.14445;
126  StepThickness[1]=.8; Weight[1]=.05665;
127  StepThickness[2]=1.6; Weight[2]=.05049;
128  StepThickness[3]=2.4; Weight[3]=.04239;
129  StepThickness[4]=3.2; Weight[4]=.04313;
130  StepThickness[5]=4.0; Weight[5]=.03879;
131  StepThickness[6]=4.8; Weight[6]=.04182;
132  StepThickness[7]=5.6; Weight[7]=.03422;
133  StepThickness[8]=6.4; Weight[8]=.03469;
134  StepThickness[9]=7.2; Weight[9]=.03589;
135  StepThickness[10]=8.0; Weight[10]=.03633;
136  StepThickness[11]=8.8; Weight[11]=.03842;
137  StepThickness[12]=9.6; Weight[12]=.03688;
138  StepThickness[13]=10.4; Weight[13]=.03705;
139  StepThickness[14]=11.2; Weight[14]=.03773;
140  StepThickness[15]=12.0; Weight[15]=.03968;
141  StepThickness[16]=12.8; Weight[16]=.04058;
142  StepThickness[17]=13.6; Weight[17]=.03903;
143  StepThickness[18]=14.4; Weight[18]=.04370;
144  StepThickness[19]=15.2; Weight[19]=.03981;
145  StepThickness[20]=16.0; Weight[20]=.05226;
146  StepThickness[21]=16.8; Weight[21]=.03603;
148 
149 }
152 {
153  delete [] Weight;
154  delete [] StepThickness;
155  delete [] StartingAngle;
156  delete [] SpanningAngle;
157  delete [] PositionMod;
158  delete [] solidMod;
159  delete [] logicMod;
160  delete [] physiMod;
161  delete solidMod1;
162  delete logicMod1;
163  delete physiMod1;
164  delete solidMod2;
165  delete logicMod2;
166  delete physiMod2;
167  delete solidMod3;
168  delete logicMod3;
169  delete physiMod3;
170  delete solidMod4;
171  delete logicMod4;
172  delete physiMod4;
173 // The Modulator wheel properties is getting form an external file "ModoulWeight.txt"
174  File.open(Name, std::ios::in);
175  if(!File.is_open())
176  {
177  G4cout<<" WARNING: The File with name of "<<Name<<
178  " doesn't exist to get modulator step properties. please modify it and try again"<<G4endl;
179 
180  G4Exception("HadrontherapyModulator::ModulatorPropertiesFromFile( )", "Hadrontherapy0009"
181  , FatalException, "Error: No available external file for reading from");
182  }
183 
184  G4String string;
185  File >>string>> StepNumbers;
186  File >>string>>string>>string;
187 
188 
189  Weight=new G4double[StepNumbers];
190  StepThickness=new G4double[StepNumbers];
191  StartingAngle=new G4double[StepNumbers];
192  SpanningAngle=new G4double[StepNumbers];
193  PositionMod=new G4ThreeVector[StepNumbers];
194 
195 
196  solidMod=new G4Tubs *[StepNumbers];
197  logicMod=new G4LogicalVolume *[StepNumbers];
198  physiMod=new G4VPhysicalVolume *[(4*(StepNumbers-1)+1)];
199 
200  for(G4int i=0;i<StepNumbers;i++)
201  {
202  G4String stringX;
203  File>>stringX>> StepThickness[i]>>Weight[i];
204  }
205 
207  BuildSteps();
208 
209 
210 
211 }
214 {
215 
216  G4double TotalWeight=0;
217  // convert the absolute weight values to relative ones
218  for(G4int i=0;i<StepNumbers;i++)
219  {
220  TotalWeight+=Weight[i];
221  }
222 
223  for(G4int i=0;i<StepNumbers;i++)
224  {
225  Weight[i]=Weight[i]/TotalWeight;
226  }
227 
228  // To build the RMW step layers will be put one after another
229 
230  StartingAngle[0]=0 *deg;
231  SpanningAngle[0]=90 *deg;
232  G4double PositionModx;
233  G4double WholeStartingAngle=0 *deg;
234  G4double WholeThickness=0;
235  for(G4int i=1;i<StepNumbers;i++)
236  {
237  StartingAngle[i]=WholeStartingAngle+(Weight[i-1]*(2*pi))/8;
238  SpanningAngle[i]=90* deg -2*StartingAngle[i];
239  StepThickness[i]=StepThickness[i]-WholeThickness;
240  PositionModx=WholeThickness+StepThickness[i]/2.;
241  PositionMod[i]=G4ThreeVector(0,0,PositionModx);
242  WholeThickness+=StepThickness[i];
243  WholeStartingAngle=StartingAngle[i];
244  }
245 
246 
247 }
250 {
251  G4bool isotopes = false;
252  G4Material* airNist = G4NistManager::Instance()->FindOrBuildMaterial("G4_AIR", isotopes);
253 
254 
255  Mod0Mater = airNist;
256  ModMater = airNist; // You have to change modulator material via a macrofile (default is air)
257 
258  innerRadiusOfTheTube = 2.5 *cm;
259  outerRadiusOfTheTube = 9.5 *cm;
260 
261  // Mother of the modulator wheel
262  G4ThreeVector positionMotherMod = G4ThreeVector(-2160.50 *mm, 30 *mm, 50 *mm);
263 
264  G4Box* solidMotherMod = new G4Box("MotherMod", 12 *cm, 12 *cm, 12 *cm);
265 
266  logicMotherMod = new G4LogicalVolume(solidMotherMod, Mod0Mater,"MotherMod",0,0,0);
267 
268  physiMotherMod = new G4PVPlacement(rm,positionMotherMod, "MotherMod",
269  logicMotherMod,
270  motherVolume,
271  false,
272  0);
273  BuildSteps();
274 
275 
276 
277  }
280  {
281  //----------------------------------------------------------
282  // Mother volume of first quarter of the modulator
283  //----------------------------------------------------------
284 
285  G4double hightOfTheTube0 = 10.0 *cm;
286  G4double startAngleOfTheTube0 = 0 *deg;
287  G4double spanningAngleOfTheTube0 = 90 *deg;
288 
289  G4RotationMatrix rm1;
290  rm1.rotateZ(0 *deg);
291 
292  G4ThreeVector positionMod1 = G4ThreeVector(0*cm,0*cm,0*cm);
293 
294  solidMod1 = new G4Tubs("Mod1",
295  innerRadiusOfTheTube,
296  outerRadiusOfTheTube,
297  hightOfTheTube0/2.,
298  startAngleOfTheTube0,
299  spanningAngleOfTheTube0);
300 
301  logicMod1 = new G4LogicalVolume(solidMod1, Mod0Mater, "Mod1",0,0,0);
302 
303  physiMod1 = new G4PVPlacement(G4Transform3D(rm1, positionMod1),
304  logicMod1,
305  "Mod1",
306  logicMotherMod,
307  false,
308  0);
309 
310 
311  //----------------------------------------------------------
312  // modulator steps
313  //----------------------------------------------------------
314  for (G4int i=1;i<StepNumbers;i++)
315  {
316 
317  solidMod[i] = new G4Tubs("Modstep",
318  innerRadiusOfTheTube,
319  outerRadiusOfTheTube,
320  StepThickness[i]/2.,
321  StartingAngle[i],
322  SpanningAngle[i]);
323 
324  logicMod[i] = new G4LogicalVolume(solidMod[i],
325  ModMater, "Modstep",0,0,0);
326 
327  physiMod[i] = new G4PVPlacement(0,
328  PositionMod[i],
329  logicMod[i],
330  "Modstep",
331  logicMod1,
332  false,
333  0);
334 
335 
336  }
337 
339  //----------------------------------------------------------
340  // Mother volume of the second modulator quarter
341  //----------------------------------------------------------
342 
343 
344  G4RotationMatrix rm2;
345  rm2.rotateZ(90 *deg);
346 
347  G4ThreeVector positionMod2 = G4ThreeVector(0*cm,0*cm,0*cm);
348 
349  solidMod2 = new G4Tubs("Mod2",
350  innerRadiusOfTheTube,
351  outerRadiusOfTheTube,
352  hightOfTheTube0/2.,
353  startAngleOfTheTube0,
354  spanningAngleOfTheTube0);
355 
356  logicMod2 = new G4LogicalVolume(solidMod2,
357  Mod0Mater, "Mod2",0,0,0);
358 
359 
360  physiMod2 = new G4PVPlacement(G4Transform3D(rm2, positionMod2),
361  logicMod2,
362  "Mod2",
363  logicMotherMod,
364  false,
365  0);
366 
367 
368  for (G4int i=1;i<StepNumbers;i++)
369  {
370 
371  physiMod[StepNumbers+i-1] = new G4PVPlacement(0,
372  PositionMod[i],
373  logicMod[i],
374  "Modstep",
375  logicMod2,
376  false,
377  0);
378 
379  }
380 
381 
382 
383  //----------------------------------------------------------
384  // Mother volume of the third modulator quarter
385  //----------------------------------------------------------
386 
387 
389  rm3.rotateZ(180 *deg);
390 
391  G4ThreeVector positionMod3 = G4ThreeVector(0*cm,0*cm,0*cm);
392 
393  solidMod3 = new G4Tubs("Mod3",
394  innerRadiusOfTheTube,
395  outerRadiusOfTheTube,
396  hightOfTheTube0,
397  startAngleOfTheTube0/2.,
398  spanningAngleOfTheTube0);
399 
400  logicMod3 = new G4LogicalVolume(solidMod3,
401  Mod0Mater, "Mod3",0,0,0);
402 
403 
404  physiMod3 = new G4PVPlacement(G4Transform3D(rm3, positionMod3),
405  logicMod3, // its logical volume
406  "Mod3", // its name
407  logicMotherMod, // its mother volume
408  false, // no boolean operations
409  0); // no particular field
410 
411 
412 
413 
414  for (G4int i=1;i<StepNumbers;i++)
415  {
416 
417  physiMod[2*(StepNumbers-1)+i] = new G4PVPlacement(0,
418  PositionMod[i],
419  logicMod[i],
420  "Modstep",
421  logicMod3,
422  false,
423  0);
424 
425  }
426 
427  //----------------------------------------------------------
428  // Mother volume of the fourth modulator quarter
429  //----------------------------------------------------------
430 
431 
432  G4RotationMatrix rm4;
433  rm4.rotateZ(270 *deg);
434 
435  G4ThreeVector positionMod4 = G4ThreeVector(0*cm,0*cm,0*cm);
436 
437  solidMod4 = new G4Tubs("Mod4",
438  innerRadiusOfTheTube,
439  outerRadiusOfTheTube,
440  hightOfTheTube0,
441  startAngleOfTheTube0/2.,
442  spanningAngleOfTheTube0);
443 
444  logicMod4 = new G4LogicalVolume(solidMod4,
445  Mod0Mater, "Mod4",0,0,0);
446 
447 
448  physiMod4 = new G4PVPlacement(G4Transform3D(rm4, positionMod4),
449  logicMod4,
450  "Mod4",
451  logicMotherMod,
452  false,
453  0);
454 
455 
456 for (G4int i=1;i<StepNumbers;i++)
457  {
458  physiMod[3*(StepNumbers-1)+i] = new G4PVPlacement(0,
459  PositionMod[i],
460  logicMod[i],
461  "Modstep",
462  logicMod4,
463  false,
464  0);
465  }
466  // Inform the kernel about the new geometry
468  G4RunManager::GetRunManager() -> PhysicsHasBeenModified();
469  G4VisAttributes * red = new G4VisAttributes( G4Colour(1. ,0. ,0.));
470  red-> SetVisibility(true);
471  red-> SetForceSolid(true);
472  logicMotherMod -> SetVisAttributes(G4VisAttributes::Invisible);
473 
478 
479  for (G4int i=1;i<StepNumbers;i++)
480  {
481  logicMod[i] -> SetVisAttributes(red);
482  }
483 
484  }
485 
487 // Messenger values
490 {
491  G4double rotationAngle = angle;
492  rm -> rotateZ(rotationAngle);
493  physiMotherMod -> SetRotation(rm);
494  G4cout << "MODULATOR HAS BEEN ROTATED OF " << rotationAngle/deg
495  << " deg" << G4endl;
497 }
499 // Change modulator material
501 {
502  if (G4Material* NewMaterial = G4NistManager::Instance()->FindOrBuildMaterial(Material, false) )
503  {
504  if (NewMaterial)
505  {
506  for(G4int i=1;i<StepNumbers;i++)
507  {
508  logicMod[i] -> SetMaterial(NewMaterial);
509  // G4RunManager::GetRunManager() -> PhysicsHasBeenModified();
511 
512  // G4cout<<(logicMod[i]->GetMaterial()->GetName())<<G4endl;
513  }
514  G4cout << "The material of the Modulator wheel has been changed to " << Material << G4endl;
515  }
516  }
517  else
518  {
519  G4cout << "WARNING: material \"" << Material << "\" doesn't exist in NIST elements/materials"
520  " table [located in $G4INSTALL/source/materials/src/G4NistMaterialBuilder.cc]" << G4endl;
521  G4cout << "Use command \"/parameter/nist\" to see full materials list!" << G4endl;
522 
523 
524  }
525 }
526 
528 // Change modulator position in the beam line
530 {
531  G4ThreeVector NewModulatorPos=Pos;
532  physiMotherMod -> SetTranslation( NewModulatorPos);
534  G4cout << "The modulator wheel is translated to"<< NewModulatorPos/mm <<"mm " <<G4endl;
535 
536 }
538 //change modulator inner raduis
540 {
541 solidMod1 -> SetInnerRadius(newvalue);
542 solidMod2 -> SetInnerRadius(newvalue);
543 solidMod3 -> SetInnerRadius(newvalue);
544 solidMod4 -> SetInnerRadius(newvalue);
545  for(G4int i=1;i<StepNumbers;i++)
546  {
547  solidMod[i] -> SetInnerRadius(newvalue);}
549  G4cout << "InnerRadius of the Modulator Wheel has been changed to :"
550  << newvalue/mm<<" mm"<< G4endl;
551 }
553 //change modulator outer raduis
555 {
556 solidMod1 -> SetOuterRadius(newvalue);
557 solidMod2 -> SetOuterRadius(newvalue);
558 solidMod3 -> SetOuterRadius(newvalue);
559 solidMod4 -> SetOuterRadius(newvalue);
560  for(G4int i=1;i<StepNumbers;i++)
561  {
562  solidMod[i] -> SetOuterRadius(newvalue);}
564  G4cout << "OuterRadius of the Modulator Wheel has been changed to :"
565  << newvalue/mm<<" mm"<<G4endl;
566 }
569 
570 {
571 G4String Name=value;
572 if(value=="default" )
573 {
574 Name=FileName;
575 }
576 G4cout<<" Step properties of modulator will be get out from the external file "
577  <<Name<<G4endl;
579 }
void SetModulatorPosition(G4ThreeVector)
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
void GetDataFromFile(G4String value)
static constexpr double mm
Definition: G4SIunits.hh:115
Definition: test07.cc:36
CLHEP::Hep3Vector G4ThreeVector
CLHEP::HepRotation G4RotationMatrix
Definition: G4Box.hh:64
ush Pos
Definition: deflate.h:89
Definition: G4Tubs.hh:85
static G4double angle[DIM]
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
void BuildModulator(G4VPhysicalVolume *)
bool G4bool
Definition: G4Types.hh:79
static constexpr double cm
Definition: G4SIunits.hh:119
HepGeom::Transform3D G4Transform3D
def SetMaterial
Definition: EmPlot.py:25
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:79
static const G4VisAttributes Invisible
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
#define G4endl
Definition: G4ios.hh:61
void ModulatorPropertiesFromFile(G4String)
double G4double
Definition: G4Types.hh:76
static constexpr double deg
Definition: G4SIunits.hh:152
void SetVisAttributes(const G4VisAttributes *pVA)