]>
Commit | Line | Data |
---|---|---|
2817d3e2 | 1 | // $Id$ |
0375ab40 | 2 | // Category: global |
2817d3e2 | 3 | // |
e5967ab3 | 4 | // Author: I. Hrivnacova |
5 | // | |
6 | // Class TG4G3CutVector | |
7 | // -------------------- | |
2817d3e2 | 8 | // See the class description in the header file. |
9 | ||
0375ab40 | 10 | #include "TG4G3CutVector.h" |
2817d3e2 | 11 | #include "TG4G3Defaults.h" |
12 | ||
13 | #include <G4Track.hh> | |
14 | #include <G4ParticleDefinition.hh> | |
15 | #include <G4VProcess.hh> | |
16 | ||
2e385ef6 | 17 | #include <g4std/strstream> |
18 | ||
e5967ab3 | 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 | ||
72095f7c | 24 | //_____________________________________________________________________________ |
0375ab40 | 25 | TG4G3CutVector::TG4G3CutVector() |
e5967ab3 | 26 | : fDeltaRaysOn(true) |
2817d3e2 | 27 | { |
e5967ab3 | 28 | // fill name vector |
29 | if (fgCutNameVector.size() == 0) FillCutNameVector(); | |
30 | ||
2817d3e2 | 31 | // initialize fCutVector |
e5967ab3 | 32 | for (G4int i=0; i<kTOFMAX; i++) fCutVector.push_back(0.); |
33 | fCutVector.push_back(DBL_MAX); | |
2817d3e2 | 34 | } |
35 | ||
72095f7c | 36 | //_____________________________________________________________________________ |
0375ab40 | 37 | TG4G3CutVector::TG4G3CutVector(const TG4G3CutVector& right) |
e5967ab3 | 38 | : fCutVector(right.fCutVector.size()) |
2817d3e2 | 39 | { |
58c0119e | 40 | // copy stuff |
41 | *this = right; | |
2817d3e2 | 42 | } |
43 | ||
72095f7c | 44 | //_____________________________________________________________________________ |
0375ab40 | 45 | TG4G3CutVector::~TG4G3CutVector() { |
2817d3e2 | 46 | // |
2817d3e2 | 47 | } |
48 | ||
49 | // operators | |
50 | ||
72095f7c | 51 | //_____________________________________________________________________________ |
0375ab40 | 52 | TG4G3CutVector& TG4G3CutVector::operator=(const TG4G3CutVector& right) |
2817d3e2 | 53 | { |
54 | // check assignement to self | |
55 | if (this == &right) return *this; | |
56 | ||
e5967ab3 | 57 | // copy fCutVector |
58 | for (G4int i=0; i<kNoG3Cuts; i++) fCutVector[i] = right.fCutVector[i]; | |
59 | ||
60 | fDeltaRaysOn = right.fDeltaRaysOn; | |
2817d3e2 | 61 | |
62 | return *this; | |
63 | } | |
64 | ||
72095f7c | 65 | //_____________________________________________________________________________ |
0375ab40 | 66 | G4double TG4G3CutVector::operator[](G4int index) const |
2817d3e2 | 67 | { |
68 | // | |
69 | if (index < kNoG3Cuts) | |
e5967ab3 | 70 | return fCutVector[index]; |
2817d3e2 | 71 | else { |
72 | TG4Globals::Exception( | |
0375ab40 | 73 | "TG4G3CutVector::operator[]: index out of the vector scope"); |
2817d3e2 | 74 | return 0.; |
75 | } | |
76 | } | |
77 | ||
e5967ab3 | 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 | ||
2817d3e2 | 100 | // public methods |
101 | ||
72095f7c | 102 | //_____________________________________________________________________________ |
e5967ab3 | 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) | |
2817d3e2 | 136 | { |
137 | // Sets the cutValue for the specified cut. | |
138 | // --- | |
139 | ||
e5967ab3 | 140 | if (cut>=kNoG3Cuts) { |
2817d3e2 | 141 | TG4Globals::Exception( |
0375ab40 | 142 | "TG4G3CutVector::SetG3Cut: Inconsistent cut."); |
2817d3e2 | 143 | } |
e5967ab3 | 144 | |
145 | fCutVector[cut] = cutValue; | |
2817d3e2 | 146 | } |
147 | ||
72095f7c | 148 | //_____________________________________________________________________________ |
0375ab40 | 149 | void TG4G3CutVector::SetG3Defaults() |
2817d3e2 | 150 | { |
151 | // Sets G3 default values for all cuts. | |
152 | // --- | |
153 | ||
e5967ab3 | 154 | for (G4int i=0; i<kNoG3Cuts; i++) |
155 | fCutVector[i] = TG4G3Defaults::Instance()->CutValue(i); | |
2817d3e2 | 156 | } |
157 | ||
72095f7c | 158 | //_____________________________________________________________________________ |
e5967ab3 | 159 | G4bool TG4G3CutVector::Update(const TG4G3CutVector& vector) |
2817d3e2 | 160 | { |
e5967ab3 | 161 | // Update only those values that have not yet been set ( =0.) |
162 | // with given vector. | |
163 | // Returns true if some value was modified. | |
2817d3e2 | 164 | // --- |
165 | ||
e5967ab3 | 166 | G4bool result = false; |
2817d3e2 | 167 | |
e5967ab3 | 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; | |
2817d3e2 | 175 | } |
176 | ||
e5967ab3 | 177 | //_____________________________________________________________________________ |
2e385ef6 | 178 | G4String TG4G3CutVector::Format() const |
e5967ab3 | 179 | { |
2e385ef6 | 180 | // Formats the output into a string. |
e5967ab3 | 181 | // --- |
182 | ||
2e385ef6 | 183 | strstream tmpStream; |
184 | tmpStream << " G3 cut vector:" << G4endl; | |
e5967ab3 | 185 | if (fDeltaRaysOn) |
2e385ef6 | 186 | tmpStream << " Delta rays On" << G4endl; |
e5967ab3 | 187 | else |
2e385ef6 | 188 | tmpStream << " Delta rays Off" << G4endl; |
e5967ab3 | 189 | |
190 | // energy cuts | |
191 | for (G4int i=0; i<kTOFMAX; i++) { | |
192 | ||
2e385ef6 | 193 | tmpStream << " " << fgCutNameVector[i] << " cut value (MeV): "; |
e5967ab3 | 194 | |
195 | if (i == kDCUTE && !fDeltaRaysOn) | |
2e385ef6 | 196 | tmpStream << fgkDCUTEOff/MeV << G4endl; |
e5967ab3 | 197 | else if (i == kDCUTM && !fDeltaRaysOn) |
2e385ef6 | 198 | tmpStream << fgkDCUTMOff/MeV << G4endl; |
e5967ab3 | 199 | else |
2e385ef6 | 200 | tmpStream << fCutVector[i]/MeV << G4endl; |
e5967ab3 | 201 | } |
202 | ||
203 | // time cut | |
2e385ef6 | 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(); | |
e5967ab3 | 217 | } |
218 | ||
72095f7c | 219 | //_____________________________________________________________________________ |
0375ab40 | 220 | G4double TG4G3CutVector::GetMinEkineForGamma(const G4Track& track) const |
2817d3e2 | 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 | ||
2817d3e2 | 228 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); |
229 | G4String processName = ""; | |
230 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
231 | ||
232 | if ((processName == "eBrem") || (processName == "IeBrem")) { | |
e5967ab3 | 233 | return fCutVector[kBCUTE]; |
2817d3e2 | 234 | } |
e5967ab3 | 235 | else if ((processName == "MuBrems") || (processName == "IMuBrems")) { |
236 | //(processName == "//hBrems")|| (processName == "//IhBrems") | |
2817d3e2 | 237 | // hadron Brehmstrahlung is not defined in G4 |
e5967ab3 | 238 | return fCutVector[kBCUTM]; |
2817d3e2 | 239 | } |
240 | else { | |
e5967ab3 | 241 | return fCutVector[kCUTGAM]; |
2817d3e2 | 242 | } |
243 | } | |
244 | ||
72095f7c | 245 | //_____________________________________________________________________________ |
0375ab40 | 246 | G4double TG4G3CutVector::GetMinEkineForElectron(const G4Track& track) const |
2817d3e2 | 247 | { |
e5967ab3 | 248 | // Returns the cut value for e-. |
2817d3e2 | 249 | // --- |
250 | ||
251 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); | |
252 | G4String processName = ""; | |
253 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
254 | ||
255 | if ((processName == "eIoni") || (processName == "IeIoni")) { | |
e5967ab3 | 256 | // delta rays by e-, e+ |
257 | ||
258 | if (fDeltaRaysOn) | |
259 | return fCutVector[kDCUTE]; | |
260 | else | |
261 | return fgkDCUTEOff; | |
2817d3e2 | 262 | } |
e5967ab3 | 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; | |
2817d3e2 | 271 | } |
e5967ab3 | 272 | else if (processName == "MuPairProd") { |
273 | //direct pair production by muons | |
274 | ||
275 | return fCutVector[kPPCUTM]; | |
276 | } | |
2817d3e2 | 277 | else { |
e5967ab3 | 278 | return fCutVector[kCUTELE]; |
2817d3e2 | 279 | } |
280 | } | |
281 | ||
72095f7c | 282 | //_____________________________________________________________________________ |
e5967ab3 | 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 | |
2817d3e2 | 302 | { |
303 | // Returns the cut value for charged hadron. | |
304 | // --- | |
305 | ||
e5967ab3 | 306 | return fCutVector[kCUTHAD]; |
2817d3e2 | 307 | } |
308 | ||
72095f7c | 309 | //_____________________________________________________________________________ |
0375ab40 | 310 | G4double TG4G3CutVector::GetMinEkineForNeutralHadron(const G4Track& track) const |
2817d3e2 | 311 | { |
312 | // Returns the cut value for neutral hadron. | |
313 | // --- | |
314 | ||
e5967ab3 | 315 | return fCutVector[kCUTNEU]; |
2817d3e2 | 316 | } |
317 | ||
72095f7c | 318 | //_____________________________________________________________________________ |
0375ab40 | 319 | G4double TG4G3CutVector::GetMinEkineForMuon(const G4Track& track) const |
2817d3e2 | 320 | { |
321 | // Returns the cut value for neutral muon. | |
322 | // --- | |
323 | ||
e5967ab3 | 324 | return fCutVector[kCUTMUO]; |
2817d3e2 | 325 | } |
326 | ||
72095f7c | 327 | //_____________________________________________________________________________ |
0375ab40 | 328 | G4double TG4G3CutVector::GetMinEkineForOther(const G4Track& track) const |
2817d3e2 | 329 | { |
330 | // Returns 0. | |
331 | // --- | |
332 | ||
333 | return 0.; | |
334 | } | |
e5967ab3 | 335 | //_____________________________________________________________________________ |
336 | G4bool TG4G3CutVector::IsCut() const | |
337 | { | |
338 | // Returns true if any of cuts is set. | |
339 | // --- | |
2817d3e2 | 340 | |
e5967ab3 | 341 | for (G4int i=0; i<kNoG3Cuts; i++) |
342 | if (fCutVector[i] > fgkTolerance ) return true; | |
343 | ||
344 | return false; | |
345 | } |