Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
examples/extended/medical/dna/neuron/src/PrimaryGeneratorAction.cc
이 파일의 문서화 페이지로 가기
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 example is provided by the Geant4-DNA collaboration
27 // Any report or published results obtained using the Geant4-DNA software
28 // shall cite the following Geant4-DNA collaboration publication:
29 // Med. Phys. 37 (2010) 4692-4708
30 // and papers
31 // M. Batmunkh et al. J Radiat Res Appl Sci 8 (2015) 498-507
32 // O. Belov et al. Physica Medica 32 (2016) 1510-1520
33 // The Geant4-DNA web site is available at http://geant4-dna.org
34 //
35 // -------------------------------------------------------------------
36 // November 2016
37 // -------------------------------------------------------------------
38 //
39 // $ID$
42 
43 #include "PrimaryGeneratorAction.hh"
44 #include "G4SystemOfUnits.hh"
45 //
46 #include "CommandLineParser.hh"
47 #include "G4ParticleTable.hh"
48 #include "G4ParticleGun.hh"
49 #include "G4ParticleDefinition.hh"
50 #include "G4PhysicalConstants.hh"
51 #include "Randomize.hh"
52 #include "G4Event.hh"
53 #include "G4Box.hh"
54 #include "G4Orb.hh"
55 #include "G4LogicalVolume.hh"
56 #include "G4LogicalVolumeStore.hh"
57 #include "G4VPhysicalVolume.hh"
58 #include "G4PhysicalVolumeStore.hh"
59 
60 using namespace G4DNAPARSER ;
61 
62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
63 
66 {
67  G4int n_particle = 1;
68  fpParticleGun = new G4ParticleGun(n_particle);
69  //
71  G4ParticleDefinition* particle = particleTable->FindParticle("proton");
73  // default gun parameters
77 }
78 
79 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
80 
82 {
83  delete fpParticleGun;
84 }
85 
86 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
87 
89 {
90 
91  // Initial kinetic energy of particles beam as Gaussion distribution!
92 
93  // G4double Ekin = 10.*MeV;
94  // G4double deltaEkin = 9.0955e-5*MeV;
95  // fpParticleGun->SetParticleEnergy(G4RandGauss::shoot(Ekin,deltaEkin));
96 
97  // In order to avoid dependence of PrimaryGeneratorAction
98  // on DetectorConstruction class we get world volume
99  // from G4LogicalVolumeStore
100 
101  // We included three options for particles direction:
102 
103  G4double mediumRadius = 0.;
104  G4double boundingXHalfLength = 0.;
105  G4double boundingYHalfLength = 0.;
106  G4double boundingZHalfLength = 0.;
107  G4LogicalVolume* mediumLV
109  G4LogicalVolume* boundingLV
110  = G4LogicalVolumeStore::GetInstance()->GetVolume("BoundingSlice");
111  G4Orb* mediumSphere = 0;
112  G4Box* boundingSlice = 0;
113  if ( mediumLV && boundingLV)
114  {
115  mediumSphere = dynamic_cast< G4Orb*>(mediumLV->GetSolid());
116  boundingSlice = dynamic_cast< G4Box*>(boundingLV->GetSolid());
117  }
118  if ( mediumSphere && boundingSlice)
119  {
120  mediumRadius = mediumSphere->GetRadius();
121  boundingXHalfLength = boundingSlice->GetXHalfLength();
122  boundingYHalfLength = boundingSlice->GetYHalfLength();
123  boundingZHalfLength = boundingSlice->GetZHalfLength();
124 
127 
128  if ( CommandLineParser::GetParser()->GetCommandIfActive("-sXY"))
129  {
130  //G4cerr << " Initial beam position uniformly spread on a square! "
131  // << G4endl;
132  // INITIAL BEAM DIVERGENCE
134  // // INITIAL BEAM POSITION
136  CLHEP::RandFlat::shoot(-boundingXHalfLength,boundingXHalfLength),
137  CLHEP::RandFlat::shoot(-boundingYHalfLength,boundingYHalfLength),
138  -mediumRadius));
139  // Surface area of a square
140  // fGunArea = 4*boundingXHalfLength*boundingYHalfLength ;
141  //G4cerr << " Particle Fluence Area on a Square (um2) = "
142  // <<fGunArea / (um*um)<< G4endl;
143  }
144 
147 
148  else if ( CommandLineParser::GetParser()->GetCommandIfActive("-dXY"))
149  {
150  //G4cerr << " Initial beam position uniformly spread on a disk! "
151  // << G4endl;
153  G4double x0,y0,z0;
154  x0 = 100.*mm;
155  y0 = 100.*mm;
156  z0 = -mediumRadius; // mediumRadius;
157  while (! (std::sqrt(x0*x0+y0*y0)<= mediumRadius) )
158  {
159  x0 = CLHEP::RandFlat::shoot(-mediumRadius,mediumRadius);
160  y0 = CLHEP::RandFlat::shoot(-mediumRadius,mediumRadius);
161  }
163  // Surface area of a disk
164  // fGunArea = pi*mediumRadius*mediumRadius ;
165  // G4cerr << " Particle Fluence Area on a Disk (um2) = "
166  // <<fGunArea / (um*um)<< G4endl;
167 
168  }
169 
171  // Select a starting position on a sphere including the
172  // target volume and neuron morphology
173  else
174  {
175  //G4cerr << " Initial beam position uniformly spread on a sphere! "<< G4endl;
176  G4double cosTheta = 2.*G4UniformRand()-1;
177  G4double sinTheta = std::sqrt(1.-cosTheta*cosTheta);
178  G4double phi = twopi*G4UniformRand();
179  G4ThreeVector positionStart(mediumRadius*sinTheta*std::cos(phi),
180  mediumRadius*sinTheta*std::sin(phi),
181  mediumRadius*cosTheta);
182  fpParticleGun->SetParticlePosition(positionStart);
183  // To compute the direction, select a point inside the target volume
184  G4ThreeVector positionDir(
185  boundingXHalfLength*(2.*G4UniformRand()-1),
186  boundingYHalfLength*(2.*G4UniformRand()-1),
187  boundingZHalfLength*(2.*G4UniformRand()-1));
189  (positionDir-positionStart).unit());
190  // Surface area of sphere
191  fGunArea = 4.*pi*mediumRadius*mediumRadius ;
192  // G4cerr << " Particle Fluence Area on sphere (um2) = "
193  // <<fGunArea / (um*um)<< G4endl;
194  }
195 
196  }
197  else
198  {
199  G4cerr << "Bounding slice volume not found!" << G4endl;
200  G4cerr << "Default particle kinematic used" << G4endl;
201  }
202 
204 }
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4double GetYHalfLength() const
CLHEP::Hep3Vector G4ThreeVector
static G4ParticleTable * GetParticleTable()
static constexpr double MeV
Definition: G4SIunits.hh:214
G4double GetRadius() const
static constexpr double mm
Definition: G4SIunits.hh:115
#define G4endl
Definition: G4ios.hh:61
void SetParticlePosition(G4ThreeVector aPosition)
static constexpr double um
Definition: G4SIunits.hh:113
double G4double
Definition: G4Types.hh:76
void SetParticleDefinition(G4ParticleDefinition *aParticleDefinition)
Definition: G4Box.hh:64
#define G4UniformRand()
Definition: Randomize.hh:53
static constexpr double twopi
Definition: G4SIunits.hh:76
Definition: G4Orb.hh:62
G4double GetXHalfLength() const
G4GLOB_DLL std::ostream G4cerr
int G4int
Definition: G4Types.hh:78
virtual void GeneratePrimaryVertex(G4Event *evt)
static double shoot()
Definition: RandFlat.cc:59
G4VSolid * GetSolid() const
void SetParticleEnergy(G4double aKineticEnergy)
static constexpr double pi
Definition: G4SIunits.hh:75
G4double GetZHalfLength() const
void SetParticleMomentumDirection(G4ParticleMomentum aMomentumDirection)
static G4LogicalVolumeStore * GetInstance()
G4LogicalVolume * GetVolume(const G4String &name, G4bool verbose=true) const