Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/advanced/iort_therapy/src/IORTPrimaryGeneratorAction.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /examples/advanced/iort_therapy/src/IORTPrimaryGeneratorAction.cc (Version 11.3.0) and /examples/advanced/iort_therapy/src/IORTPrimaryGeneratorAction.cc (Version 9.5.p2)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 // This is the *BASIC* version of IORT, a Gean     26 // This is the *BASIC* version of IORT, a Geant4-based application
 27 //                                                 27 //
 28 // Main Authors: G.Russo(a,b), C.Casarino*(c),     28 // Main Authors: G.Russo(a,b), C.Casarino*(c), G.C. Candiano(c), G.A.P. Cirrone(d), F.Romano(d)
 29 // Contributor Authors: S.Guatelli(e)              29 // Contributor Authors: S.Guatelli(e)
 30 // Past Authors: G.Arnetta(c), S.E.Mazzaglia(d     30 // Past Authors: G.Arnetta(c), S.E.Mazzaglia(d)
 31 //                                                 31 //    
 32 //   (a) Fondazione Istituto San Raffaele G.Gi     32 //   (a) Fondazione Istituto San Raffaele G.Giglio, Cefalù, Italy
 33 //   (b) IBFM-CNR , Segrate (Milano), Italy        33 //   (b) IBFM-CNR , Segrate (Milano), Italy
 34 //   (c) LATO (Laboratorio di Tecnologie Oncol     34 //   (c) LATO (Laboratorio di Tecnologie Oncologiche), Cefalù, Italy
 35 //   (d) Laboratori Nazionali del Sud of the I     35 //   (d) Laboratori Nazionali del Sud of the INFN, Catania, Italy
 36 //   (e) University of Wollongong, Australia   <<  36 //   (e) University of Wallongong, Australia
 37 //                                                 37 //
 38 //   *Corresponding author, email to carlo.cas     38 //   *Corresponding author, email to carlo.casarino@polooncologicocefalu.it
 39 //////////////////////////////////////////////     39 //////////////////////////////////////////////////////////////////////////////////////////////
 40                                                    40 
 41 #include "G4SystemOfUnits.hh"                  << 
 42 #include "IORTPrimaryGeneratorAction.hh"           41 #include "IORTPrimaryGeneratorAction.hh"
 43 #include "IORTPrimaryGeneratorMessenger.hh"        42 #include "IORTPrimaryGeneratorMessenger.hh"
 44                                                << 
 45 #include "globals.hh"                          << 
 46 #include "G4Event.hh"                              43 #include "G4Event.hh"
 47 #include "G4ParticleGun.hh"                        44 #include "G4ParticleGun.hh"
 48 #include "G4ParticleTable.hh"                      45 #include "G4ParticleTable.hh"
 49 #include "G4ParticleDefinition.hh"                 46 #include "G4ParticleDefinition.hh"
 50 #include "Randomize.hh"                            47 #include "Randomize.hh"
 51                                                <<  48 #include "IORTAnalysisManager.hh"      
                                                   >>  49 
 52 IORTPrimaryGeneratorAction::IORTPrimaryGenerat     50 IORTPrimaryGeneratorAction::IORTPrimaryGeneratorAction()
 53 {                                                  51 {
 54   // Define the messenger                          52   // Define the messenger
 55   gunMessenger = new IORTPrimaryGeneratorMesse     53   gunMessenger = new IORTPrimaryGeneratorMessenger(this);
 56                                                    54 
 57   particleGun  = new G4ParticleGun();              55   particleGun  = new G4ParticleGun();
 58                                                    56 
 59   SetDefaultPrimaryParticle();                     57   SetDefaultPrimaryParticle();  
 60 }                                                  58 }  
 61                                                    59 
 62 IORTPrimaryGeneratorAction::~IORTPrimaryGenera     60 IORTPrimaryGeneratorAction::~IORTPrimaryGeneratorAction()
 63 {                                                  61 {
 64   delete particleGun;                              62   delete particleGun;
 65                                                    63 
 66   delete gunMessenger;                             64   delete gunMessenger;
 67 }                                                  65 }
 68                                                    66   
 69 void IORTPrimaryGeneratorAction::SetDefaultPri     67 void IORTPrimaryGeneratorAction::SetDefaultPrimaryParticle()
 70 {                                                  68 {    
 71   // ****************************                  69   // ****************************
 72   // Default primary particle                      70   // Default primary particle
 73   // ****************************                  71   // ****************************
 74                                                    72   
 75   // Define primary particles: electrons // pr     73   // Define primary particles: electrons // protons
 76   G4ParticleTable* particleTable = G4ParticleT     74   G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
 77   G4ParticleDefinition* particle = particleTab     75   G4ParticleDefinition* particle = particleTable -> FindParticle("e-");   // ("proton")
 78   particleGun -> SetParticleDefinition(particl     76   particleGun -> SetParticleDefinition(particle); 
 79                                                    77 
 80   // Define the energy of primary particles:       78   // Define the energy of primary particles:
 81   // gaussian distribution with mean energy =      79   // gaussian distribution with mean energy = 10.0 *MeV   
 82   // and sigma = 400.0 *keV                        80   // and sigma = 400.0 *keV
 83   G4double defaultMeanKineticEnergy = 10.0 *CL <<  81   G4double defaultMeanKineticEnergy = 10.0 *MeV;   
 84   meanKineticEnergy = defaultMeanKineticEnergy     82   meanKineticEnergy = defaultMeanKineticEnergy;
 85                                                    83 
 86   G4double defaultsigmaEnergy = 100.0 *CLHEP:: <<  84   G4double defaultsigmaEnergy = 100.0 *keV;   
 87   sigmaEnergy = defaultsigmaEnergy;                85   sigmaEnergy = defaultsigmaEnergy;
 88                                                <<  86   
                                                   >>  87 #ifdef G4ANALYSIS_USE_ROOT    
                                                   >>  88   // Write these values into the analysis if needed. Have to be written separately on change.
                                                   >>  89   IORTAnalysisManager::GetInstance()->setBeamMetaData(meanKineticEnergy, sigmaEnergy);
                                                   >>  90 #endif                       
                                                   >>  91 
 89   // Define the parameters of the initial posi     92   // Define the parameters of the initial position: 
 90   // the y, z coordinates have a gaussian dist <<  93   // the y, z coordinates have a gaussian distribution
                                                   >>  94   
 91                                                    95   
 92   G4double defaultX0 = -862.817 *CLHEP::mm;    <<  96   G4double defaultX0 = -862.817 *mm;                 
 93   X0 = defaultX0;                                  97   X0 = defaultX0;
 94                                                    98 
 95   G4double defaultY0 = 0.0 *CLHEP::mm;         <<  99   G4double defaultY0 = 0.0 *mm;  
 96   Y0 = defaultY0;                                 100   Y0 = defaultY0;
 97                                                   101 
 98   G4double defaultZ0 = 0.0 *CLHEP::mm;         << 102   G4double defaultZ0 = 0.0 *mm;  
 99   Z0 = defaultZ0;                                 103   Z0 = defaultZ0;
100                                                   104 
101   G4double defaultsigmaY = 1. *CLHEP::mm;      << 105   G4double defaultsigmaY = 1. *mm;  
102   sigmaY = defaultsigmaY;                         106   sigmaY = defaultsigmaY;
103                                                   107 
104   G4double defaultsigmaZ = 1. *CLHEP::mm;      << 108   G4double defaultsigmaZ = 1. *mm;  
105   sigmaZ = defaultsigmaZ;                         109   sigmaZ = defaultsigmaZ;
106                                                   110 
107   // Define the parameters of the momentum of     111   // Define the parameters of the momentum of primary particles: 
108   // The momentum along the y and z axis has a    112   // The momentum along the y and z axis has a gaussian distribution
109                                                   113  
110 /*                                                114 /* 
111   G4double defaultsigmaMomentumY = 0.0;           115   G4double defaultsigmaMomentumY = 0.0;  
112   sigmaMomentumY = defaultsigmaMomentumY;         116   sigmaMomentumY = defaultsigmaMomentumY;
113                                                   117 
114   G4double defaultsigmaMomentumZ = 0.0;           118   G4double defaultsigmaMomentumZ = 0.0;  
115   sigmaMomentumZ = defaultsigmaMomentumZ;         119   sigmaMomentumZ = defaultsigmaMomentumZ;
116 */                                                120 */
117                                                   121 
118   G4double defaultTheta = 6.0 *CLHEP::deg;     << 122   G4double defaultTheta = 6.0 *deg;  
119   Theta = defaultTheta;                           123   Theta = defaultTheta;
120                                                   124 
121 }                                                 125 }
122                                                   126 
123 void IORTPrimaryGeneratorAction::GeneratePrima    127 void IORTPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
124 {                                                 128 {
                                                   >> 129 #ifdef G4ANALYSIS_USE_ROOT
                                                   >> 130   // Increment the event counter
                                                   >> 131   IORTAnalysisManager::GetInstance()->startNewEvent();
                                                   >> 132 #endif
                                                   >> 133 
125   // ****************************************     134   // ****************************************
126   // Set the beam angular apread                  135   // Set the beam angular apread 
127   // and spot size                                136   // and spot size
128   // beam spot size                               137   // beam spot size
129   // ****************************************     138   // ****************************************
130                                                   139 
131   // Set the position of the primary particles    140   // Set the position of the primary particles
132   G4double x = X0;                                141   G4double x = X0;
133   G4double y = Y0;                                142   G4double y = Y0;
134   G4double z = Z0;                                143   G4double z = Z0;
135                                                   144 
136   if ( sigmaY > 0.0 )                             145   if ( sigmaY > 0.0 )
137     {                                             146     {
138       y += G4RandGauss::shoot( Y0, sigmaY );      147       y += G4RandGauss::shoot( Y0, sigmaY );
139     }                                             148     }
140                                                   149  
141   if ( sigmaZ > 0.0 )                             150   if ( sigmaZ > 0.0 )
142     {                                             151     {
143       z += G4RandGauss::shoot( Z0, sigmaZ );      152       z += G4RandGauss::shoot( Z0, sigmaZ );
144     }                                             153     }
145                                                   154 
146   particleGun -> SetParticlePosition(G4ThreeVe    155   particleGun -> SetParticlePosition(G4ThreeVector( x , y , z ) );
147                                                   156  
148   // *****************************************    157   // ********************************************
149   // Set the beam energy and energy spread        158   // Set the beam energy and energy spread
150   // *****************************************    159   // ********************************************
151                                                   160 
152   G4double kineticEnergy = G4RandGauss::shoot(    161   G4double kineticEnergy = G4RandGauss::shoot( meanKineticEnergy, sigmaEnergy );
153   particleGun -> SetParticleEnergy ( kineticEn    162   particleGun -> SetParticleEnergy ( kineticEnergy );
154                                                   163 
155   // Set the direction of the primary particle    164   // Set the direction of the primary particles  
156                                                   165   
157                                                   166 
158   /*                                              167   /*                          
159   G4double momentumX = 1.0;                       168   G4double momentumX = 1.0;
160   G4double momentumY = 0.0;                       169   G4double momentumY = 0.0;
161   G4double momentumZ = 0.0;                       170   G4double momentumZ = 0.0;
162                                                   171 
163   if ( sigmaMomentumY  > 0.0 )                    172   if ( sigmaMomentumY  > 0.0 )
164     {                                             173     {
165       momentumY += G4RandGauss::shoot( 0., sig    174       momentumY += G4RandGauss::shoot( 0., sigmaMomentumY );
166     }                                             175     }
167   if ( sigmaMomentumZ  > 0.0 )                    176   if ( sigmaMomentumZ  > 0.0 )
168     {                                             177     {
169       momentumZ += G4RandGauss::shoot( 0., sig    178       momentumZ += G4RandGauss::shoot( 0., sigmaMomentumZ );
170     }                                             179     }
171                                                   180  
172   particleGun -> SetParticleMomentumDirection(    181   particleGun -> SetParticleMomentumDirection( G4ThreeVector(momentumX,momentumY,momentumZ) );
173                                                   182 
174   */                                              183   */
175                                                   184   
176                                                   185                
177   G4double Mx;                                    186   G4double Mx;    
178   G4double My;                                    187   G4double My;
179   G4double Mz;                                    188   G4double Mz;
180   G4double condizione;                            189   G4double condizione;
181                                                   190   
182 while (true)  {                                   191 while (true)  {
183                                                   192 
184   //Mx =  CLHEP::RandFlat::shoot(0.9,1);       << 193   Mx =  CLHEP::RandFlat::shoot(0.9,1);
185   //My =  CLHEP::RandFlat::shoot(-0.1,0.1);    << 194   My =  CLHEP::RandFlat::shoot(-0.1,0.1);
186   //Mz =  CLHEP::RandFlat::shoot(-0.1,0.1);    << 195   Mz =  CLHEP::RandFlat::shoot(-0.1,0.1);
187                                                << 196   condizione = sqrt(Mx*Mx + My*My + Mz*Mz);
188   Mx =  CLHEP::RandFlat::shoot(0.7,1);         << 
189   My =  CLHEP::RandFlat::shoot(-0.3,0.3); // r << 
190   Mz =  CLHEP::RandFlat::shoot(-0.3,0.3);      << 
191                                                << 
192   condizione = std::sqrt(Mx*Mx + My*My + Mz*Mz << 
193                                                   197 
194                                                   198  
195   if (condizione < 1)  {                          199   if (condizione < 1)  {
196     Mx = Mx/condizione;                           200     Mx = Mx/condizione;
197     My = My/condizione;                           201     My = My/condizione;
198     Mz = Mz/condizione;                           202     Mz = Mz/condizione;
199                                                   203 
200                                                   204 
201     if (Mx > std::cos(Theta)) {                << 205     if (Mx > cos(Theta)) { 
202       break;                                      206       break;
203         }                                         207         }
204     }                                             208     }
205 }                                                 209 }
206                                                   210   
207                                                   211  
208   particleGun -> SetParticleMomentumDirection(    212   particleGun -> SetParticleMomentumDirection( G4ThreeVector(Mx,My,Mz) );
209                                                   213   
210                                                   214 
211   // Generate a primary particle                  215   // Generate a primary particle
212   particleGun -> GeneratePrimaryVertex( anEven    216   particleGun -> GeneratePrimaryVertex( anEvent ); 
213 }                                                 217 } 
214                                                   218 
215 void IORTPrimaryGeneratorAction::SetmeanKineti    219 void IORTPrimaryGeneratorAction::SetmeanKineticEnergy (G4double val )  
216 {                                                 220 {
217   meanKineticEnergy = val;                     << 221   meanKineticEnergy = val;
218   G4cout << "The mean Kinetic energy of the in << 222 #ifdef G4ANALYSIS_USE_ROOT
219          << meanKineticEnergy/MeV << G4endl;   << 223   // Update the beam-data in the analysis manager
                                                   >> 224   IORTAnalysisManager::GetInstance()->setBeamMetaData(meanKineticEnergy, sigmaEnergy);
                                                   >> 225 #endif
                                                   >> 226 
220 }                                                 227 } 
221                                                   228 
222 void IORTPrimaryGeneratorAction::SetsigmaEnerg    229 void IORTPrimaryGeneratorAction::SetsigmaEnergy (G4double val )  
223 {                                                 230 { 
224  sigmaEnergy = val;                            << 231   sigmaEnergy = val;
225  G4cout << "The sigma of the kinetic energy of << 232 #ifdef G4ANALYSIS_USE_ROOT
226         << sigmaEnergy/MeV << G4endl;          << 233   // Update the sigmaenergy in the metadata.
                                                   >> 234   IORTAnalysisManager::GetInstance()->setBeamMetaData(meanKineticEnergy, sigmaEnergy);
                                                   >> 235 #endif
227 }                                                 236 }
228                                                   237 
229 void IORTPrimaryGeneratorAction::SetXposition     238 void IORTPrimaryGeneratorAction::SetXposition (G4double val )  
230 { X0 = val;}                                      239 { X0 = val;}
231                                                   240 
232 void IORTPrimaryGeneratorAction::SetYposition     241 void IORTPrimaryGeneratorAction::SetYposition (G4double val )  
233 { Y0 = val;}                                      242 { Y0 = val;}
234                                                   243 
235 void IORTPrimaryGeneratorAction::SetZposition     244 void IORTPrimaryGeneratorAction::SetZposition (G4double val )  
236 { Z0 = val;}                                      245 { Z0 = val;}
237                                                   246 
238 void IORTPrimaryGeneratorAction::SetsigmaY (G4    247 void IORTPrimaryGeneratorAction::SetsigmaY (G4double val )  
239 { sigmaY = val;}                                  248 { sigmaY = val;}
240                                                   249 
241 void IORTPrimaryGeneratorAction::SetsigmaZ (G4    250 void IORTPrimaryGeneratorAction::SetsigmaZ (G4double val )  
242 { sigmaZ = val;}                                  251 { sigmaZ = val;}
243                                                   252 
244 /*                                                253 /*
245 void IORTPrimaryGeneratorAction::SetsigmaMomen    254 void IORTPrimaryGeneratorAction::SetsigmaMomentumY (G4double val )  
246 { sigmaMomentumY = val;}                          255 { sigmaMomentumY = val;}
247                                                   256 
248 void IORTPrimaryGeneratorAction::SetsigmaMomen    257 void IORTPrimaryGeneratorAction::SetsigmaMomentumZ (G4double val )  
249 { sigmaMomentumZ = val;}                          258 { sigmaMomentumZ = val;}
250 */                                                259 */
251                                                   260 
252 void IORTPrimaryGeneratorAction::SetTheta (G4d    261 void IORTPrimaryGeneratorAction::SetTheta (G4double val ) 
253 { Theta = val;}                                   262 { Theta = val;}
254                                                   263 
255                                                   264 
256 G4double IORTPrimaryGeneratorAction::GetmeanKi    265 G4double IORTPrimaryGeneratorAction::GetmeanKineticEnergy(void)
257 { return meanKineticEnergy;}                      266 { return meanKineticEnergy;}
258                                                   267 
259                                                   268