| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 /// \file electromagnetic/TestEm5/src/Tracking
27 /// \brief Implementation of the TrackingActio
28 //
29 //
30 //
31 //....oooOO0OOooo........oooOO0OOooo........oo
32 //....oooOO0OOooo........oooOO0OOooo........oo
33
34 #include "TrackingAction.hh"
35
36 #include "DetectorConstruction.hh"
37 #include "EventAction.hh"
38 #include "HistoManager.hh"
39 #include "Run.hh"
40
41 #include "G4PhysicalConstants.hh"
42 #include "G4RunManager.hh"
43 #include "G4SystemOfUnits.hh"
44 #include "G4Track.hh"
45
46 //....oooOO0OOooo........oooOO0OOooo........oo
47
48 TrackingAction::TrackingAction(DetectorConstru
49 : fDetector(det), fEventAction(event)
50 {}
51
52 //....oooOO0OOooo........oooOO0OOooo........oo
53
54 void TrackingAction::PreUserTrackingAction(con
55 {
56 // few initialisations
57 //
58 if (aTrack->GetTrackID() == 1) {
59 fXstartAbs = fDetector->GetxstartAbs();
60 fXendAbs = fDetector->GetxendAbs();
61 fPrimaryCharge = aTrack->GetDefinition()->
62 fDirX = aTrack->GetMomentumDirection().x()
63 }
64 }
65
66 //....oooOO0OOooo........oooOO0OOooo........oo
67
68 void TrackingAction::PostUserTrackingAction(co
69 {
70 G4AnalysisManager* analysisManager = G4Analy
71
72 Run* run = static_cast<Run*>(G4RunManager::G
73
74 G4ThreeVector position = aTrack->GetPosition
75 G4ThreeVector vertex = aTrack->GetVertexPosi
76 G4double charge = aTrack->GetDefinition()->G
77 G4double energy = aTrack->GetKineticEnergy()
78
79 G4bool transmit =
80 (((fDirX >= 0.0 && position.x() >= fXendAb
81 && energy > 0.0);
82 G4bool reflect =
83 (((fDirX >= 0.0 && position.x() <= fXstart
84 && energy > 0.0);
85 G4bool notabsor = (transmit || reflect);
86
87 // transmitted + reflected particles counter
88 //
89 G4int flag = 0;
90 if (charge == fPrimaryCharge) flag = 1;
91 if (aTrack->GetTrackID() == 1) flag = 2;
92 if (transmit) fEventAction->SetTransmitFlag(
93 if (reflect) fEventAction->SetReflectFlag(fl
94
95 //
96 // histograms
97 //
98 G4bool charged = (charge != 0.);
99 G4bool neutral = !charged;
100
101 // energy spectrum at exit
102 // zero energy charged particle are not take
103 //
104 G4int id = 0;
105 if (transmit && charged)
106 id = 10;
107 else if (transmit && neutral)
108 id = 20;
109 else if (reflect && charged)
110 id = 30;
111 else if (reflect && neutral)
112 id = 40;
113
114 if (id > 0) {
115 analysisManager->FillH1(id, energy);
116 }
117
118 // energy leakage
119 //
120 if (notabsor) {
121 G4int trackID = aTrack->GetTrackID();
122 G4int index = 0;
123 if (trackID > 1) index = 1; // primary=0,
124 G4double eleak = aTrack->GetKineticEnergy(
125 if ((aTrack->GetDefinition() == G4Positron
126 eleak += 2 * electron_mass_c2;
127 run->AddEnergyLeak(eleak, index);
128 }
129
130 // space angle distribution at exit : dN/dOm
131 //
132 G4ThreeVector direction = aTrack->GetMomentu
133 id = 0;
134 if (transmit && charged)
135 id = 12;
136 else if (transmit && neutral)
137 id = 22;
138 else if (reflect && charged)
139 id = 32;
140 else if (reflect && neutral)
141 id = 42;
142
143 if (id > 0 && std::abs(direction.x()) < 1.0)
144 G4double theta = std::acos(direction.x());
145 if (theta > 0.0) {
146 G4double dteta = analysisManager->GetH1W
147 G4double unit = analysisManager->GetH1Un
148 G4double weight = unit / (twopi * std::s
149 /*
150 G4cout << "theta, dteta, unit, weight: "
151 << theta << " "
152 << dteta << " "
153 << unit << " "
154 << weight << G4endl;
155 */
156 analysisManager->FillH1(id, theta, weigh
157 }
158 }
159
160 // energy fluence at exit : dE(MeV)/dOmega
161 //
162 id = 0;
163 if (transmit && charged)
164 id = 11;
165 else if (reflect && charged)
166 id = 31;
167 else if (transmit && neutral)
168 id = 21;
169 else if (reflect && neutral)
170 id = 41;
171
172 if (id > 0 && std::abs(direction.x()) < 1.0)
173 G4double theta = std::acos(direction.x());
174 if (theta > 0.0) {
175 G4double dteta = analysisManager->GetH1W
176 G4double unit = analysisManager->GetH1Un
177 G4double weight = unit / (twopi * std::s
178 weight *= (aTrack->GetKineticEnergy() /
179 analysisManager->FillH1(id, theta, weigh
180 }
181 }
182
183 // projected angles distribution at exit
184 //
185 id = 0;
186 if (transmit && charged)
187 id = 13;
188 else if (transmit && neutral)
189 id = 23;
190 else if (reflect && charged)
191 id = 33;
192 else if (reflect && neutral)
193 id = 43;
194
195 if (id > 0) {
196 if (direction.x() != 0.0) {
197 G4double tet = std::atan(direction.y() /
198 analysisManager->FillH1(id, tet);
199 if (transmit && (flag == 2)) run->AddMsc
200
201 tet = std::atan(direction.z() / std::fab
202 analysisManager->FillH1(id, tet);
203 if (transmit && (flag == 2)) run->AddMsc
204 }
205 }
206
207 // projected position and radius at exit
208 //
209 if (transmit && energy > 0.0) {
210 G4double y = position.y(), z = position.z(
211 G4double r = std::sqrt(y * y + z * z);
212 analysisManager->FillH1(14, y);
213 analysisManager->FillH1(14, z);
214 analysisManager->FillH1(15, r);
215 }
216
217 // x-vertex of charged secondaries
218 //
219 if ((aTrack->GetParentID() == 1) && charged)
220 G4double xVertex = (aTrack->GetVertexPosit
221 analysisManager->FillH1(6, xVertex);
222 if (notabsor) analysisManager->FillH1(7, x
223 }
224 }
225
226 //....oooOO0OOooo........oooOO0OOooo........oo
227