| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
26 // <<
27 /// \file medical/GammaTherapy/src/PrimaryGene <<
28 /// \brief Implementation of the PrimaryGenera <<
29 // <<
30 23
31 //-------------------------------------------- 24 //---------------------------------------------------------------------------
32 // 25 //
33 // ClassName: PrimaryGeneratorAction 26 // ClassName: PrimaryGeneratorAction
34 // 27 //
35 // Description: Generate primary beam 28 // Description: Generate primary beam
36 // 29 //
37 // Authors: V.Grichine, V.Ivanchenko 30 // Authors: V.Grichine, V.Ivanchenko
38 // 31 //
39 // Modified: 32 // Modified:
40 // 33 //
41 //-------------------------------------------- 34 //----------------------------------------------------------------------------
42 // 35 //
43 36
44 #include "PrimaryGeneratorAction.hh" 37 #include "PrimaryGeneratorAction.hh"
45 <<
46 #include "DetectorConstruction.hh" 38 #include "DetectorConstruction.hh"
47 #include "PrimaryGeneratorMessenger.hh" 39 #include "PrimaryGeneratorMessenger.hh"
48 << 40 #include "Randomize.hh"
49 #include "G4ParticleDefinition.hh" <<
50 #include "G4ParticleGun.hh" 41 #include "G4ParticleGun.hh"
51 #include "G4ParticleTable.hh" 42 #include "G4ParticleTable.hh"
52 #include "G4PhysicalConstants.hh" << 43 #include "G4ParticleDefinition.hh"
53 #include "G4SystemOfUnits.hh" << 44 #include "Histo.hh"
54 #include "Randomize.hh" <<
55 45
56 //....oooOO0OOooo........oooOO0OOooo........oo 46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
57 //....oooOO0OOooo........oooOO0OOooo........oo 47 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
58 48
59 PrimaryGeneratorAction::PrimaryGeneratorAction << 49 PrimaryGeneratorAction::PrimaryGeneratorAction(DetectorConstruction* pDet):
>> 50 fDetector(pDet)
60 { 51 {
61 InitializeMe(); 52 InitializeMe();
62 } 53 }
63 54
64 //....oooOO0OOooo........oooOO0OOooo........oo 55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
65 56
66 void PrimaryGeneratorAction::InitializeMe() 57 void PrimaryGeneratorAction::InitializeMe()
67 { 58 {
68 fVerbose = fDetector->GetVerbose(); << 59 theMessenger = new PrimaryGeneratorMessenger(this);
69 fMessenger = new PrimaryGeneratorMessenger(t << 60 particleGun = new G4ParticleGun();
70 fParticleGun = new G4ParticleGun(); << 61 counter = 0;
71 fCounter = 0; << 62 verbose = 0;
72 fX0 = 0.0; << 63 x0 = 0.0;
73 fY0 = 0.0; << 64 y0 = 0.0;
74 fZ0 = 0.0; << 65 z0 = 0.0;
75 fSigmaX = 1.5 * mm; << 66 sigmaX = 1.5*mm;
76 fSigmaY = 1.5 * mm; << 67 sigmaY = 1.5*mm;
77 fSigmaZ = 0.0; << 68 sigmaZ = 0.0;
78 fSigmaE = 0.0; << 69 sigmaE = 0.0;
79 fRMax2 = 2.5 * 2.5 * mm * mm; << 70 rMax2 = 2.5*2.5*mm*mm;
80 fSigmaTheta = 0.0; << 71 sigmaTheta = 0.0;
81 // fSigmaTheta = 0.17*degree; << 72 // sigmaTheta = 0.17*degree;
82 fMinCosTheta = 2.0; << 73 minCosTheta = 2.0;
83 SetBeamEnergy(50.0 * MeV); << 74 SetBeamEnergy(50.0*MeV);
84 fPosition = G4ThreeVector(fX0, fY0, fZ0); << 75 position = G4ThreeVector(x0,y0,z0);
85 fDirection = G4ThreeVector(0.0, 0.0, 1.0); << 76 direction = G4ThreeVector(0.0,0.0,1.0);
86 fGauss = true; << 77 m_gauss = true;
>> 78 if(energy < (Histo::GetPointer())->GetMaxEnergy())
>> 79 (Histo::GetPointer())->SetMaxEnergy(energy);
87 } 80 }
88 81
89 //....oooOO0OOooo........oooOO0OOooo........oo 82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
90 83
91 PrimaryGeneratorAction::~PrimaryGeneratorActio 84 PrimaryGeneratorAction::~PrimaryGeneratorAction()
92 { 85 {
93 delete fParticleGun; << 86 delete particleGun;
94 delete fMessenger; << 87 delete theMessenger;
95 } 88 }
96 89
97 //....oooOO0OOooo........oooOO0OOooo........oo 90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
98 91
99 void PrimaryGeneratorAction::GeneratePrimaries 92 void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
100 { 93 {
101 fCounter++; << 94 counter++ ;
>> 95 verbose = (Histo::GetPointer())->GetVerbose();
102 96
103 // Simulation of beam position 97 // Simulation of beam position
104 G4double x = fX0; << 98 G4double x = x0;
105 G4double y = fY0; << 99 G4double y = y0;
106 G4double z = fDetector->GetGeneratorPosZ(); 100 G4double z = fDetector->GetGeneratorPosZ();
107 do { 101 do {
108 if (0.0 < fSigmaX) { << 102 if(0.0 < sigmaX) x = G4RandGauss::shoot(x0,sigmaX);
109 x = G4RandGauss::shoot(fX0, fSigmaX); << 103 if(0.0 < sigmaY) y = G4RandGauss::shoot(y0,sigmaY);
110 } << 104 } while (x*x + y*y > rMax2);
111 if (0.0 < fSigmaY) { <<
112 y = G4RandGauss::shoot(fY0, fSigmaY); <<
113 } <<
114 } while (x * x + y * y > fRMax2); <<
115 105
116 fPosition = G4ThreeVector(x, y, z); << 106 position = G4ThreeVector(x,y,z);
117 fParticleGun->SetParticlePosition(fPosition) << 107 particleGun->SetParticlePosition(position);
118 108
119 // Simulation of beam direction 109 // Simulation of beam direction
120 G4double ux = fDirection.x(); << 110 G4double ux = direction.x();
121 G4double uy = fDirection.y(); << 111 G4double uy = direction.y();
122 G4double uz = fDirection.z(); << 112 G4double uz = direction.z();
123 113
124 // Beam particles are uniformly distributed 114 // Beam particles are uniformly distributed over phi, cosTheta
125 if (1.0 > fMinCosTheta) { << 115 if(1.0 > minCosTheta) {
126 uz = fMinCosTheta + (1.0 - fMinCosTheta) * << 116 uz = minCosTheta + (1.0 - minCosTheta)*G4UniformRand() ;
127 ux = std::sqrt((1.0 - uz) * (1.0 + uz)); << 117 ux = std::sqrt(1.0 - uz*uz) ;
>> 118 } else if (sigmaTheta > 0.0) {
>> 119 ux = G4RandGauss::shoot(0.0,sigmaTheta);
>> 120 uz = std::sqrt(1.0 - ux*ux);
128 } 121 }
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 122
140 fParticleGun->SetParticleMomentumDirection(f << 123 G4double phi = twopi*G4UniformRand() ;
>> 124 uy = ux ;
>> 125 ux *= std::cos(phi) ;
>> 126 uy *= std::sin(phi) ;
>> 127 direction = G4ThreeVector(ux,uy,uz) ;
>> 128
>> 129 direction = direction.unit();
>> 130 particleGun->SetParticleMomentumDirection(direction);
>> 131 G4ParticleDefinition* particle = particleGun->GetParticleDefinition();
141 132
142 // Simulation of beam kinetic energy 133 // Simulation of beam kinetic energy
143 G4double kinEnergy = fEnergy; << 134 G4double kinEnergy = energy;
144 135
145 if (fGauss == "flatE") { << 136 if(m_gauss == "flatE") kinEnergy = energy - sigmaE + 2.*sigmaE*G4UniformRand();
146 kinEnergy = fEnergy - fSigmaE + 2. * fSigm << 137 else if(0.0 < sigmaE) kinEnergy = energy + G4RandGauss::shoot(0.0,sigmaE);
147 } << 138
148 else if (0.0 < fSigmaE) { << 139 particleGun->SetParticleEnergy(kinEnergy);
149 kinEnergy = fEnergy + G4RandGauss::shoot(0 << 140
150 } << 141 G4String particleName = particle->GetParticleName() ;
151 fParticleGun->SetParticleEnergy(kinEnergy); <<
152 142
153 if (fVerbose > 0) { << 143 if(verbose > 0) {
154 G4ParticleDefinition* particle = fParticle << 144 G4cout << "Event# " << counter
155 G4String particleName = particle->GetParti << 145 << " Beam particle is generated by PrimaryGeneratorAction "
156 G4cout << "Event# " << fCounter << " Bea <<
157 << G4endl; 146 << G4endl;
158 G4cout << "ParticleName= " << particleName << 147 G4cout << "ParticleName= " << particleName
159 << std::setprecision(5) << " KinE << 148 << " PDGcode= " << particle->GetPDGEncoding()
160 << " x(mm)= " << x / mm << " y(mm << 149 << std::setprecision(5)
161 << " ux= " << ux << " uy= " << uy << 150 << " KinEnergy(GeV)= "
162 } << 151 << energy/GeV
>> 152 << " x(mm)= "
>> 153 << x/mm
>> 154 << " y(mm)= "
>> 155 << y/mm
>> 156 << " z(mm)= "
>> 157 << z/mm
>> 158 << " ux= "
>> 159 << ux
>> 160 << " uy= "
>> 161 << uy
>> 162 << " uz= "
>> 163 << uz
>> 164 << G4endl;
>> 165 }
163 166
164 fParticleGun->GeneratePrimaryVertex(anEvent) << 167 particleGun->GeneratePrimaryVertex(anEvent);
165 } 168 }
166 169
167 //....oooOO0OOooo........oooOO0OOooo........oo 170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
168 171
>> 172 void PrimaryGeneratorAction::SetBeamSigmaE(G4double val)
>> 173 {
>> 174 sigmaE = val;
>> 175 }
>> 176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 177
169 void PrimaryGeneratorAction::SetBeamEnergy(G4d 178 void PrimaryGeneratorAction::SetBeamEnergy(G4double val)
170 { 179 {
171 fEnergy = val; << 180 energy = val;
172 if (fEnergy < fDetector->GetMaxEnergy()) fDe << 181 if(energy < (Histo::GetPointer())->GetMaxEnergy())
>> 182 (Histo::GetPointer())->SetMaxEnergy(energy);
173 } 183 }
174 184
175 //....oooOO0OOooo........oooOO0OOooo........oo 185 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 186
>> 187
>> 188
176 189