Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-27 09:18:47

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // This is the *BASIC* version of IORT, a Geant4-based application
0027 //
0028 // Main Authors: G.Russo(a,b), C.Casarino*(c), G.C. Candiano(c), G.A.P. Cirrone(d), F.Romano(d)
0029 // Contributor Authors: S.Guatelli(e)
0030 // Past Authors: G.Arnetta(c), S.E.Mazzaglia(d)
0031 //    
0032 //   (a) Fondazione Istituto San Raffaele G.Giglio, Cefalù, Italy
0033 //   (b) IBFM-CNR , Segrate (Milano), Italy
0034 //   (c) LATO (Laboratorio di Tecnologie Oncologiche), Cefalù, Italy
0035 //   (d) Laboratori Nazionali del Sud of the INFN, Catania, Italy
0036 //   (e) University of Wollongong, Australia
0037 //
0038 //   *Corresponding author, email to carlo.casarino@polooncologicocefalu.it
0039 //////////////////////////////////////////////////////////////////////////////////////////////
0040 
0041 #include "G4SystemOfUnits.hh"
0042 #include "IORTPrimaryGeneratorAction.hh"
0043 #include "IORTPrimaryGeneratorMessenger.hh"
0044 
0045 #include "globals.hh"
0046 #include "G4Event.hh"
0047 #include "G4ParticleGun.hh"
0048 #include "G4ParticleTable.hh"
0049 #include "G4ParticleDefinition.hh"
0050 #include "Randomize.hh"
0051   
0052 IORTPrimaryGeneratorAction::IORTPrimaryGeneratorAction()
0053 {
0054   // Define the messenger
0055   gunMessenger = new IORTPrimaryGeneratorMessenger(this);
0056 
0057   particleGun  = new G4ParticleGun();
0058 
0059   SetDefaultPrimaryParticle();  
0060 }  
0061 
0062 IORTPrimaryGeneratorAction::~IORTPrimaryGeneratorAction()
0063 {
0064   delete particleGun;
0065 
0066   delete gunMessenger;
0067 }
0068   
0069 void IORTPrimaryGeneratorAction::SetDefaultPrimaryParticle()
0070 {    
0071   // ****************************
0072   // Default primary particle
0073   // ****************************
0074   
0075   // Define primary particles: electrons // protons
0076   G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
0077   G4ParticleDefinition* particle = particleTable -> FindParticle("e-");   // ("proton")
0078   particleGun -> SetParticleDefinition(particle); 
0079 
0080   // Define the energy of primary particles:
0081   // gaussian distribution with mean energy = 10.0 *MeV   
0082   // and sigma = 400.0 *keV
0083   G4double defaultMeanKineticEnergy = 10.0 *CLHEP::MeV;   
0084   meanKineticEnergy = defaultMeanKineticEnergy;
0085 
0086   G4double defaultsigmaEnergy = 100.0 *CLHEP::keV;   
0087   sigmaEnergy = defaultsigmaEnergy;
0088  
0089   // Define the parameters of the initial position: 
0090   // the y, z coordinates have a gaussian distribution 
0091   
0092   G4double defaultX0 = -862.817 *CLHEP::mm;                 
0093   X0 = defaultX0;
0094 
0095   G4double defaultY0 = 0.0 *CLHEP::mm;  
0096   Y0 = defaultY0;
0097 
0098   G4double defaultZ0 = 0.0 *CLHEP::mm;  
0099   Z0 = defaultZ0;
0100 
0101   G4double defaultsigmaY = 1. *CLHEP::mm;  
0102   sigmaY = defaultsigmaY;
0103 
0104   G4double defaultsigmaZ = 1. *CLHEP::mm;  
0105   sigmaZ = defaultsigmaZ;
0106 
0107   // Define the parameters of the momentum of primary particles: 
0108   // The momentum along the y and z axis has a gaussian distribution
0109  
0110 /* 
0111   G4double defaultsigmaMomentumY = 0.0;  
0112   sigmaMomentumY = defaultsigmaMomentumY;
0113 
0114   G4double defaultsigmaMomentumZ = 0.0;  
0115   sigmaMomentumZ = defaultsigmaMomentumZ;
0116 */
0117 
0118   G4double defaultTheta = 6.0 *CLHEP::deg;  
0119   Theta = defaultTheta;
0120 
0121 }
0122 
0123 void IORTPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
0124 {
0125   // ****************************************
0126   // Set the beam angular apread 
0127   // and spot size
0128   // beam spot size
0129   // ****************************************
0130 
0131   // Set the position of the primary particles
0132   G4double x = X0;
0133   G4double y = Y0;
0134   G4double z = Z0;
0135 
0136   if ( sigmaY > 0.0 )
0137     {
0138       y += G4RandGauss::shoot( Y0, sigmaY );
0139     }
0140  
0141   if ( sigmaZ > 0.0 )
0142     {
0143       z += G4RandGauss::shoot( Z0, sigmaZ );
0144     }
0145 
0146   particleGun -> SetParticlePosition(G4ThreeVector( x , y , z ) );
0147  
0148   // ********************************************
0149   // Set the beam energy and energy spread
0150   // ********************************************
0151 
0152   G4double kineticEnergy = G4RandGauss::shoot( meanKineticEnergy, sigmaEnergy );
0153   particleGun -> SetParticleEnergy ( kineticEnergy );
0154 
0155   // Set the direction of the primary particles  
0156   
0157 
0158   /*                          
0159   G4double momentumX = 1.0;
0160   G4double momentumY = 0.0;
0161   G4double momentumZ = 0.0;
0162 
0163   if ( sigmaMomentumY  > 0.0 )
0164     {
0165       momentumY += G4RandGauss::shoot( 0., sigmaMomentumY );
0166     }
0167   if ( sigmaMomentumZ  > 0.0 )
0168     {
0169       momentumZ += G4RandGauss::shoot( 0., sigmaMomentumZ );
0170     }
0171  
0172   particleGun -> SetParticleMomentumDirection( G4ThreeVector(momentumX,momentumY,momentumZ) );
0173 
0174   */
0175   
0176                
0177   G4double Mx;    
0178   G4double My;
0179   G4double Mz;
0180   G4double condizione;
0181   
0182 while (true)  {
0183 
0184   //Mx =  CLHEP::RandFlat::shoot(0.9,1);
0185   //My =  CLHEP::RandFlat::shoot(-0.1,0.1);
0186   //Mz =  CLHEP::RandFlat::shoot(-0.1,0.1);
0187 
0188   Mx =  CLHEP::RandFlat::shoot(0.7,1);
0189   My =  CLHEP::RandFlat::shoot(-0.3,0.3); // ranges good for 0<Theta<20
0190   Mz =  CLHEP::RandFlat::shoot(-0.3,0.3);
0191 
0192   condizione = std::sqrt(Mx*Mx + My*My + Mz*Mz);
0193 
0194  
0195   if (condizione < 1)  {
0196     Mx = Mx/condizione;
0197     My = My/condizione;
0198     Mz = Mz/condizione;
0199 
0200 
0201     if (Mx > std::cos(Theta)) { 
0202       break;
0203         }
0204     }
0205 }
0206   
0207  
0208   particleGun -> SetParticleMomentumDirection( G4ThreeVector(Mx,My,Mz) );
0209   
0210 
0211   // Generate a primary particle
0212   particleGun -> GeneratePrimaryVertex( anEvent ); 
0213 } 
0214 
0215 void IORTPrimaryGeneratorAction::SetmeanKineticEnergy (G4double val )  
0216 {
0217   meanKineticEnergy = val;
0218   G4cout << "The mean Kinetic energy of the incident beam has been changed to (MeV):" 
0219          << meanKineticEnergy/MeV << G4endl; 
0220 } 
0221 
0222 void IORTPrimaryGeneratorAction::SetsigmaEnergy (G4double val )  
0223 { 
0224  sigmaEnergy = val;
0225  G4cout << "The sigma of the kinetic energy of the incident beam has been changed to (MeV):"
0226         << sigmaEnergy/MeV << G4endl; 
0227 }
0228 
0229 void IORTPrimaryGeneratorAction::SetXposition (G4double val )  
0230 { X0 = val;}
0231 
0232 void IORTPrimaryGeneratorAction::SetYposition (G4double val )  
0233 { Y0 = val;}
0234 
0235 void IORTPrimaryGeneratorAction::SetZposition (G4double val )  
0236 { Z0 = val;}
0237 
0238 void IORTPrimaryGeneratorAction::SetsigmaY (G4double val )  
0239 { sigmaY = val;}
0240 
0241 void IORTPrimaryGeneratorAction::SetsigmaZ (G4double val )  
0242 { sigmaZ = val;}
0243 
0244 /*
0245 void IORTPrimaryGeneratorAction::SetsigmaMomentumY (G4double val )  
0246 { sigmaMomentumY = val;}
0247 
0248 void IORTPrimaryGeneratorAction::SetsigmaMomentumZ (G4double val )  
0249 { sigmaMomentumZ = val;}
0250 */
0251 
0252 void IORTPrimaryGeneratorAction::SetTheta (G4double val ) 
0253 { Theta = val;}
0254 
0255 
0256 G4double IORTPrimaryGeneratorAction::GetmeanKineticEnergy(void)
0257 { return meanKineticEnergy;}
0258