]>
Commit | Line | Data |
---|---|---|
1 | // $Id$ | |
2 | // Category: global | |
3 | // | |
4 | // Author: I. Hrivnacova | |
5 | // | |
6 | // Class TG4G3CutVector | |
7 | // -------------------- | |
8 | // See the class description in the header file. | |
9 | ||
10 | #include "TG4G3CutVector.h" | |
11 | #include "TG4G3Defaults.h" | |
12 | ||
13 | #include <G4Track.hh> | |
14 | #include <G4ParticleDefinition.hh> | |
15 | #include <G4VProcess.hh> | |
16 | ||
17 | #include <g4std/strstream> | |
18 | ||
19 | const G4double TG4G3CutVector::fgkDCUTEOff = 10. * TeV; | |
20 | const G4double TG4G3CutVector::fgkDCUTMOff = 10. * TeV; | |
21 | const G4double TG4G3CutVector::fgkTolerance = 1. * keV; | |
22 | TG4StringVector TG4G3CutVector::fgCutNameVector; | |
23 | ||
24 | //_____________________________________________________________________________ | |
25 | TG4G3CutVector::TG4G3CutVector() | |
26 | : fDeltaRaysOn(true) | |
27 | { | |
28 | // fill name vector | |
29 | if (fgCutNameVector.size() == 0) FillCutNameVector(); | |
30 | ||
31 | // initialize fCutVector | |
32 | for (G4int i=0; i<kTOFMAX; i++) fCutVector.push_back(0.); | |
33 | fCutVector.push_back(DBL_MAX); | |
34 | } | |
35 | ||
36 | //_____________________________________________________________________________ | |
37 | TG4G3CutVector::TG4G3CutVector(const TG4G3CutVector& right) | |
38 | : fCutVector(right.fCutVector.size()) | |
39 | { | |
40 | // copy stuff | |
41 | *this = right; | |
42 | } | |
43 | ||
44 | //_____________________________________________________________________________ | |
45 | TG4G3CutVector::~TG4G3CutVector() { | |
46 | // | |
47 | } | |
48 | ||
49 | // operators | |
50 | ||
51 | //_____________________________________________________________________________ | |
52 | TG4G3CutVector& TG4G3CutVector::operator=(const TG4G3CutVector& right) | |
53 | { | |
54 | // check assignement to self | |
55 | if (this == &right) return *this; | |
56 | ||
57 | // copy fCutVector | |
58 | for (G4int i=0; i<kNoG3Cuts; i++) fCutVector[i] = right.fCutVector[i]; | |
59 | ||
60 | fDeltaRaysOn = right.fDeltaRaysOn; | |
61 | ||
62 | return *this; | |
63 | } | |
64 | ||
65 | //_____________________________________________________________________________ | |
66 | G4double TG4G3CutVector::operator[](G4int index) const | |
67 | { | |
68 | // | |
69 | if (index < kNoG3Cuts) | |
70 | return fCutVector[index]; | |
71 | else { | |
72 | TG4Globals::Exception( | |
73 | "TG4G3CutVector::operator[]: index out of the vector scope"); | |
74 | return 0.; | |
75 | } | |
76 | } | |
77 | ||
78 | // private methods | |
79 | ||
80 | //_____________________________________________________________________________ | |
81 | void TG4G3CutVector::FillCutNameVector() | |
82 | { | |
83 | // Defines fCutNameVector. | |
84 | // --- | |
85 | ||
86 | fgCutNameVector.push_back("CUTGAM"); | |
87 | fgCutNameVector.push_back("CUTELE"); | |
88 | fgCutNameVector.push_back("CUTNEU"); | |
89 | fgCutNameVector.push_back("CUTHAD"); | |
90 | fgCutNameVector.push_back("CUTMUO"); | |
91 | fgCutNameVector.push_back("BCUTE"); | |
92 | fgCutNameVector.push_back("BCUTM"); | |
93 | fgCutNameVector.push_back("DCUTE"); | |
94 | fgCutNameVector.push_back("DCUTM"); | |
95 | fgCutNameVector.push_back("PPCUTM"); | |
96 | fgCutNameVector.push_back("TOFMAX"); | |
97 | fgCutNameVector.push_back("NONE"); | |
98 | } | |
99 | ||
100 | // public methods | |
101 | ||
102 | //_____________________________________________________________________________ | |
103 | TG4G3Cut TG4G3CutVector::GetCut(const G4String& cutName) | |
104 | { | |
105 | // Retrieves corresponding TG4G3Cut constant from the cutName. | |
106 | // --- | |
107 | ||
108 | if (cutName == fgCutNameVector[kCUTGAM]) return kCUTGAM; | |
109 | else if (cutName == fgCutNameVector[kCUTELE]) return kCUTELE; | |
110 | else if (cutName == fgCutNameVector[kCUTNEU]) return kCUTNEU; | |
111 | else if (cutName == fgCutNameVector[kCUTHAD]) return kCUTHAD; | |
112 | else if (cutName == fgCutNameVector[kCUTMUO]) return kCUTMUO; | |
113 | else if (cutName == fgCutNameVector[kBCUTE]) return kBCUTE; | |
114 | else if (cutName == fgCutNameVector[kBCUTM]) return kBCUTM; | |
115 | else if (cutName == fgCutNameVector[kDCUTE]) return kDCUTE; | |
116 | else if (cutName == fgCutNameVector[kDCUTM]) return kDCUTM; | |
117 | else if (cutName == fgCutNameVector[kPPCUTM]) return kPPCUTM; | |
118 | else if (cutName == fgCutNameVector[kTOFMAX]) return kTOFMAX; | |
119 | else return kNoG3Cuts; | |
120 | } | |
121 | ||
122 | //_____________________________________________________________________________ | |
123 | const G4String& TG4G3CutVector::GetCutName(TG4G3Cut cut) | |
124 | { | |
125 | // Returns name of a specified cut. | |
126 | // --- | |
127 | ||
128 | // fill name vector | |
129 | if (fgCutNameVector.size() == 0) TG4G3CutVector::FillCutNameVector(); | |
130 | ||
131 | return fgCutNameVector[cut]; | |
132 | } | |
133 | ||
134 | //_____________________________________________________________________________ | |
135 | void TG4G3CutVector::SetCut(TG4G3Cut cut, G4double cutValue) | |
136 | { | |
137 | // Sets the cutValue for the specified cut. | |
138 | // --- | |
139 | ||
140 | if (cut>=kNoG3Cuts) { | |
141 | TG4Globals::Exception( | |
142 | "TG4G3CutVector::SetG3Cut: Inconsistent cut."); | |
143 | } | |
144 | ||
145 | fCutVector[cut] = cutValue; | |
146 | } | |
147 | ||
148 | //_____________________________________________________________________________ | |
149 | void TG4G3CutVector::SetG3Defaults() | |
150 | { | |
151 | // Sets G3 default values for all cuts. | |
152 | // --- | |
153 | ||
154 | for (G4int i=0; i<kNoG3Cuts; i++) | |
155 | fCutVector[i] = TG4G3Defaults::Instance()->CutValue(i); | |
156 | } | |
157 | ||
158 | //_____________________________________________________________________________ | |
159 | G4bool TG4G3CutVector::Update(const TG4G3CutVector& vector) | |
160 | { | |
161 | // Update only those values that have not yet been set ( =0.) | |
162 | // with given vector. | |
163 | // Returns true if some value was modified. | |
164 | // --- | |
165 | ||
166 | G4bool result = false; | |
167 | ||
168 | for (G4int i=0; i<kNoG3Cuts; i++) | |
169 | if (fCutVector[i] < fgkTolerance ) { | |
170 | fCutVector[i] = vector[i]; | |
171 | result = true; | |
172 | } | |
173 | ||
174 | return result; | |
175 | } | |
176 | ||
177 | //_____________________________________________________________________________ | |
178 | G4String TG4G3CutVector::Format() const | |
179 | { | |
180 | // Formats the output into a string. | |
181 | // --- | |
182 | ||
183 | strstream tmpStream; | |
184 | tmpStream << " G3 cut vector:" << G4endl; | |
185 | if (fDeltaRaysOn) | |
186 | tmpStream << " Delta rays On" << G4endl; | |
187 | else | |
188 | tmpStream << " Delta rays Off" << G4endl; | |
189 | ||
190 | // energy cuts | |
191 | for (G4int i=0; i<kTOFMAX; i++) { | |
192 | ||
193 | tmpStream << " " << fgCutNameVector[i] << " cut value (MeV): "; | |
194 | ||
195 | if (i == kDCUTE && !fDeltaRaysOn) | |
196 | tmpStream << fgkDCUTEOff/MeV << G4endl; | |
197 | else if (i == kDCUTM && !fDeltaRaysOn) | |
198 | tmpStream << fgkDCUTMOff/MeV << G4endl; | |
199 | else | |
200 | tmpStream << fCutVector[i]/MeV << G4endl; | |
201 | } | |
202 | ||
203 | // time cut | |
204 | tmpStream << " " << fgCutNameVector[kTOFMAX] << " cut value (s): " | |
205 | << fCutVector[kTOFMAX]/s << G4endl; | |
206 | ||
207 | return tmpStream.str(); | |
208 | } | |
209 | ||
210 | //_____________________________________________________________________________ | |
211 | void TG4G3CutVector::Print() const | |
212 | { | |
213 | // Prints the cuts. | |
214 | // --- | |
215 | ||
216 | G4cout << Format(); | |
217 | } | |
218 | ||
219 | //_____________________________________________________________________________ | |
220 | G4double TG4G3CutVector::GetMinEkineForGamma(const G4Track& track) const | |
221 | { | |
222 | // Returns the cut value for gamma. | |
223 | // (Cut is not applied for "opticalphoton" | |
224 | // as it is treated in G4 as a particle different | |
225 | // from "gamma" in G4.) | |
226 | // --- | |
227 | ||
228 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); | |
229 | G4String processName = ""; | |
230 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
231 | ||
232 | if ((processName == "eBrem") || (processName == "IeBrem")) { | |
233 | return fCutVector[kBCUTE]; | |
234 | } | |
235 | else if ((processName == "MuBrems") || (processName == "IMuBrems")) { | |
236 | //(processName == "//hBrems")|| (processName == "//IhBrems") | |
237 | // hadron Brehmstrahlung is not defined in G4 | |
238 | return fCutVector[kBCUTM]; | |
239 | } | |
240 | else { | |
241 | return fCutVector[kCUTGAM]; | |
242 | } | |
243 | } | |
244 | ||
245 | //_____________________________________________________________________________ | |
246 | G4double TG4G3CutVector::GetMinEkineForElectron(const G4Track& track) const | |
247 | { | |
248 | // Returns the cut value for e-. | |
249 | // --- | |
250 | ||
251 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); | |
252 | G4String processName = ""; | |
253 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
254 | ||
255 | if ((processName == "eIoni") || (processName == "IeIoni")) { | |
256 | // delta rays by e-, e+ | |
257 | ||
258 | if (fDeltaRaysOn) | |
259 | return fCutVector[kDCUTE]; | |
260 | else | |
261 | return fgkDCUTEOff; | |
262 | } | |
263 | else if ((processName == "MuIoni") || (processName == "IMuIoni") || | |
264 | (processName == "hIoni") || (processName == "IhIoni")) { | |
265 | // delta rays by other particles (mu, hadron) | |
266 | ||
267 | if (fDeltaRaysOn) | |
268 | return fCutVector[kDCUTM]; | |
269 | else | |
270 | return fgkDCUTMOff; | |
271 | } | |
272 | else if (processName == "MuPairProd") { | |
273 | //direct pair production by muons | |
274 | ||
275 | return fCutVector[kPPCUTM]; | |
276 | } | |
277 | else { | |
278 | return fCutVector[kCUTELE]; | |
279 | } | |
280 | } | |
281 | ||
282 | //_____________________________________________________________________________ | |
283 | G4double TG4G3CutVector::GetMinEkineForEplus(const G4Track& track) const | |
284 | { | |
285 | // Returns the cut value for e+. | |
286 | // --- | |
287 | ||
288 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); | |
289 | G4String processName = ""; | |
290 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
291 | ||
292 | if (processName == "MuPairProd") { | |
293 | //direct pair production by muons | |
294 | return fCutVector[kPPCUTM]; | |
295 | } | |
296 | else | |
297 | return fCutVector[kCUTELE]; | |
298 | } | |
299 | ||
300 | //_____________________________________________________________________________ | |
301 | G4double TG4G3CutVector::GetMinEkineForChargedHadron(const G4Track& track) const | |
302 | { | |
303 | // Returns the cut value for charged hadron. | |
304 | // --- | |
305 | ||
306 | return fCutVector[kCUTHAD]; | |
307 | } | |
308 | ||
309 | //_____________________________________________________________________________ | |
310 | G4double TG4G3CutVector::GetMinEkineForNeutralHadron(const G4Track& track) const | |
311 | { | |
312 | // Returns the cut value for neutral hadron. | |
313 | // --- | |
314 | ||
315 | return fCutVector[kCUTNEU]; | |
316 | } | |
317 | ||
318 | //_____________________________________________________________________________ | |
319 | G4double TG4G3CutVector::GetMinEkineForMuon(const G4Track& track) const | |
320 | { | |
321 | // Returns the cut value for neutral muon. | |
322 | // --- | |
323 | ||
324 | return fCutVector[kCUTMUO]; | |
325 | } | |
326 | ||
327 | //_____________________________________________________________________________ | |
328 | G4double TG4G3CutVector::GetMinEkineForOther(const G4Track& track) const | |
329 | { | |
330 | // Returns 0. | |
331 | // --- | |
332 | ||
333 | return 0.; | |
334 | } | |
335 | //_____________________________________________________________________________ | |
336 | G4bool TG4G3CutVector::IsCut() const | |
337 | { | |
338 | // Returns true if any of cuts is set. | |
339 | // --- | |
340 | ||
341 | for (G4int i=0; i<kNoG3Cuts; i++) | |
342 | if (fCutVector[i] > fgkTolerance ) return true; | |
343 | ||
344 | return false; | |
345 | } |