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 /// \file medical/GammaTherapy/src/PrimaryGeneratorAction.cc 28 /// \brief Implementation of the PrimaryGeneratorAction class 29 // 30 31 //--------------------------------------------------------------------------- 32 // 33 // ClassName: PrimaryGeneratorAction 34 // 35 // Description: Generate primary beam 36 // 37 // Authors: V.Grichine, V.Ivanchenko 38 // 39 // Modified: 40 // 41 //---------------------------------------------------------------------------- 42 // 43 44 #include "PrimaryGeneratorAction.hh" 45 46 #include "DetectorConstruction.hh" 47 #include "PrimaryGeneratorMessenger.hh" 48 49 #include "G4ParticleDefinition.hh" 50 #include "G4ParticleGun.hh" 51 #include "G4ParticleTable.hh" 52 #include "G4PhysicalConstants.hh" 53 #include "G4SystemOfUnits.hh" 54 #include "Randomize.hh" 55 56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 58 59 PrimaryGeneratorAction::PrimaryGeneratorAction(DetectorConstruction* pDet) : fDetector(pDet) 60 { 61 InitializeMe(); 62 } 63 64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 65 66 void PrimaryGeneratorAction::InitializeMe() 67 { 68 fVerbose = fDetector->GetVerbose(); 69 fMessenger = new PrimaryGeneratorMessenger(this); 70 fParticleGun = new G4ParticleGun(); 71 fCounter = 0; 72 fX0 = 0.0; 73 fY0 = 0.0; 74 fZ0 = 0.0; 75 fSigmaX = 1.5 * mm; 76 fSigmaY = 1.5 * mm; 77 fSigmaZ = 0.0; 78 fSigmaE = 0.0; 79 fRMax2 = 2.5 * 2.5 * mm * mm; 80 fSigmaTheta = 0.0; 81 // fSigmaTheta = 0.17*degree; 82 fMinCosTheta = 2.0; 83 SetBeamEnergy(50.0 * MeV); 84 fPosition = G4ThreeVector(fX0, fY0, fZ0); 85 fDirection = G4ThreeVector(0.0, 0.0, 1.0); 86 fGauss = true; 87 } 88 89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 90 91 PrimaryGeneratorAction::~PrimaryGeneratorAction() 92 { 93 delete fParticleGun; 94 delete fMessenger; 95 } 96 97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 98 99 void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) 100 { 101 fCounter++; 102 103 // Simulation of beam position 104 G4double x = fX0; 105 G4double y = fY0; 106 G4double z = fDetector->GetGeneratorPosZ(); 107 do { 108 if (0.0 < fSigmaX) { 109 x = G4RandGauss::shoot(fX0, fSigmaX); 110 } 111 if (0.0 < fSigmaY) { 112 y = G4RandGauss::shoot(fY0, fSigmaY); 113 } 114 } while (x * x + y * y > fRMax2); 115 116 fPosition = G4ThreeVector(x, y, z); 117 fParticleGun->SetParticlePosition(fPosition); 118 119 // Simulation of beam direction 120 G4double ux = fDirection.x(); 121 G4double uy = fDirection.y(); 122 G4double uz = fDirection.z(); 123 124 // Beam particles are uniformly distributed over phi, cosTheta 125 if (1.0 > fMinCosTheta) { 126 uz = fMinCosTheta + (1.0 - fMinCosTheta) * G4UniformRand(); 127 ux = std::sqrt((1.0 - uz) * (1.0 + uz)); 128 } 129 else if (fSigmaTheta > 0.0) { 130 ux = G4RandGauss::shoot(0.0, fSigmaTheta); 131 uz = std::sqrt((1.0 - ux) * (1.0 + ux)); 132 } 133 134 G4double phi = twopi * G4UniformRand(); 135 uy = ux; 136 ux *= std::cos(phi); 137 uy *= std::sin(phi); 138 fDirection = G4ThreeVector(ux, uy, uz); 139 140 fParticleGun->SetParticleMomentumDirection(fDirection); 141 142 // Simulation of beam kinetic energy 143 G4double kinEnergy = fEnergy; 144 145 if (fGauss == "flatE") { 146 kinEnergy = fEnergy - fSigmaE + 2. * fSigmaE * G4UniformRand(); 147 } 148 else if (0.0 < fSigmaE) { 149 kinEnergy = fEnergy + G4RandGauss::shoot(0.0, fSigmaE); 150 } 151 fParticleGun->SetParticleEnergy(kinEnergy); 152 153 if (fVerbose > 0) { 154 G4ParticleDefinition* particle = fParticleGun->GetParticleDefinition(); 155 G4String particleName = particle->GetParticleName(); 156 G4cout << "Event# " << fCounter << " Beam particle is generated by PrimaryGeneratorAction " 157 << G4endl; 158 G4cout << "ParticleName= " << particleName << " PDGcode= " << particle->GetPDGEncoding() 159 << std::setprecision(5) << " KinEnergy(GeV)= " << fEnergy / GeV 160 << " x(mm)= " << x / mm << " y(mm)= " << y / mm << " z(mm)= " << z / mm 161 << " ux= " << ux << " uy= " << uy << " uz= " << uz << G4endl; 162 } 163 164 fParticleGun->GeneratePrimaryVertex(anEvent); 165 } 166 167 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 168 169 void PrimaryGeneratorAction::SetBeamEnergy(G4double val) 170 { 171 fEnergy = val; 172 if (fEnergy < fDetector->GetMaxEnergy()) fDetector->SetMaxEnergy(fEnergy); 173 } 174 175 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 176