]>
Commit | Line | Data |
---|---|---|
2817d3e2 | 1 | // $Id$ |
0375ab40 | 2 | // Category: global |
2817d3e2 | 3 | // |
4 | // See the class description in the header file. | |
5 | ||
0375ab40 | 6 | #include "TG4G3CutVector.h" |
2817d3e2 | 7 | #include "TG4G3Defaults.h" |
8 | ||
9 | #include <G4Track.hh> | |
10 | #include <G4ParticleDefinition.hh> | |
11 | #include <G4VProcess.hh> | |
12 | ||
72095f7c | 13 | //_____________________________________________________________________________ |
0375ab40 | 14 | TG4G3CutVector::TG4G3CutVector() |
2817d3e2 | 15 | { |
16 | // initialize fCutVector | |
58c0119e | 17 | fCutVector = new TG4doubleVector(); |
2817d3e2 | 18 | for (G4int i=0; i<kNoG3Cuts; i++) fCutVector->insert(0.); |
19 | } | |
20 | ||
72095f7c | 21 | //_____________________________________________________________________________ |
0375ab40 | 22 | TG4G3CutVector::TG4G3CutVector(const TG4G3CutVector& right) |
2817d3e2 | 23 | { |
58c0119e | 24 | // allocation |
25 | fCutVector = new TG4doubleVector(); | |
26 | ||
27 | // copy stuff | |
28 | *this = right; | |
2817d3e2 | 29 | } |
30 | ||
72095f7c | 31 | //_____________________________________________________________________________ |
0375ab40 | 32 | TG4G3CutVector::~TG4G3CutVector() { |
2817d3e2 | 33 | // |
34 | delete fCutVector; | |
35 | } | |
36 | ||
37 | // operators | |
38 | ||
72095f7c | 39 | //_____________________________________________________________________________ |
0375ab40 | 40 | TG4G3CutVector& TG4G3CutVector::operator=(const TG4G3CutVector& right) |
2817d3e2 | 41 | { |
42 | // check assignement to self | |
43 | if (this == &right) return *this; | |
44 | ||
45 | // initialize fCutVector | |
46 | fCutVector->clear(); | |
47 | for (G4int i=0; i<kNoG3Cuts; i++) { | |
48 | fCutVector->insert((*right.fCutVector)[i]); | |
49 | } | |
50 | ||
51 | return *this; | |
52 | } | |
53 | ||
72095f7c | 54 | //_____________________________________________________________________________ |
0375ab40 | 55 | G4double TG4G3CutVector::operator[](G4int index) const |
2817d3e2 | 56 | { |
57 | // | |
58 | if (index < kNoG3Cuts) | |
59 | return (*fCutVector)[index]; | |
60 | else { | |
61 | TG4Globals::Exception( | |
0375ab40 | 62 | "TG4G3CutVector::operator[]: index out of the vector scope"); |
2817d3e2 | 63 | return 0.; |
64 | } | |
65 | } | |
66 | ||
67 | // public methods | |
68 | ||
72095f7c | 69 | //_____________________________________________________________________________ |
0375ab40 | 70 | void TG4G3CutVector::SetG3Cut(TG4G3Cut cut, G4double cutValue) |
2817d3e2 | 71 | { |
72 | // Sets the cutValue for the specified cut. | |
73 | // --- | |
74 | ||
0375ab40 | 75 | if (cut<kNoG3Cuts) { |
76 | (*fCutVector)[cut] = cutValue; | |
2817d3e2 | 77 | } |
78 | else { | |
79 | TG4Globals::Exception( | |
0375ab40 | 80 | "TG4G3CutVector::SetG3Cut: Inconsistent cut."); |
2817d3e2 | 81 | } |
82 | } | |
83 | ||
72095f7c | 84 | //_____________________________________________________________________________ |
0375ab40 | 85 | void TG4G3CutVector::SetG3Defaults() |
2817d3e2 | 86 | { |
87 | // Sets G3 default values for all cuts. | |
88 | // --- | |
89 | ||
90 | for (G4int i=0; i<kNoG3Cuts; i++) { | |
91 | (*fCutVector)[i] = TG4G3Defaults::CutValue(i); | |
92 | } | |
93 | } | |
94 | ||
72095f7c | 95 | //_____________________________________________________________________________ |
0375ab40 | 96 | G4double TG4G3CutVector::GetMinEkine(const G4Track& track) const |
2817d3e2 | 97 | { |
98 | // Returns the cut value for the particle associated with | |
99 | // specified track. | |
100 | // --- | |
101 | ||
102 | G4ParticleDefinition* particle = track.GetDefinition(); | |
103 | G4String particleName = particle->GetParticleName(); | |
104 | ||
105 | if (particleName == "gamma") { | |
106 | return GetMinEkineForGamma(track); | |
107 | } | |
108 | else if (particleName == "e-") { | |
109 | return GetMinEkineForElectron(track); | |
110 | } | |
111 | else if ((particle->GetParticleType() == "baryon") || | |
112 | (particle->GetParticleType() == "meson") || | |
113 | (particle->GetParticleType() == "nucleus")) { | |
114 | if (particle->GetPDGCharge() == 0) | |
115 | return GetMinEkineForHadron(track); | |
116 | else | |
117 | return GetMinEkineForNeutralHadron(track); | |
118 | } | |
119 | else if ((particleName == "mu-") || (particleName == "mu+")) { | |
120 | return GetMinEkineForMuon(track); | |
121 | } | |
122 | else { | |
0375ab40 | 123 | G4String text = "TG4G3CutVector::GetMinEkine: \n"; |
2817d3e2 | 124 | text = text + " The kinetic energy cut for " + particleName; |
125 | text = text + " is not defined."; | |
126 | TG4Globals::Warning(text); | |
127 | return 0.; | |
128 | } | |
129 | } | |
130 | ||
72095f7c | 131 | //_____________________________________________________________________________ |
0375ab40 | 132 | G4double TG4G3CutVector::GetMinEkineForGamma(const G4Track& track) const |
2817d3e2 | 133 | { |
134 | // Returns the cut value for gamma. | |
135 | // (Cut is not applied for "opticalphoton" | |
136 | // as it is treated in G4 as a particle different | |
137 | // from "gamma" in G4.) | |
138 | // --- | |
139 | ||
0375ab40 | 140 | G4cout << "TG4G3CutVector::GetMinEkineForGamma start" << G4endl; |
2817d3e2 | 141 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); |
142 | G4String processName = ""; | |
143 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
144 | ||
145 | if ((processName == "eBrem") || (processName == "IeBrem")) { | |
146 | return (*fCutVector)[kBCUTE]; | |
147 | } | |
148 | else if ((processName == "MuBrems") || (processName == "IMuBrems") || | |
149 | (processName == "//hBrems")|| (processName == "//IhBrems")) { | |
150 | // hadron Brehmstrahlung is not defined in G4 | |
151 | return (*fCutVector)[kBCUTM]; | |
152 | } | |
153 | else { | |
154 | return (*fCutVector)[kCUTGAM]; | |
155 | } | |
156 | } | |
157 | ||
72095f7c | 158 | //_____________________________________________________________________________ |
0375ab40 | 159 | G4double TG4G3CutVector::GetMinEkineForElectron(const G4Track& track) const |
2817d3e2 | 160 | { |
161 | // Returns the cut value for e-. | |
162 | // Should these cuts be applied to e+ too ?? | |
163 | // --- | |
164 | ||
165 | const G4VProcess* kpCreatorProcess = track.GetCreatorProcess(); | |
166 | G4String processName = ""; | |
167 | if (kpCreatorProcess) processName = kpCreatorProcess->GetProcessName(); | |
168 | ||
169 | if ((processName == "eIoni") || (processName == "IeIoni")) { | |
170 | // !! Geant4 treats delta rays + continuous energy loss | |
171 | // within one process | |
172 | return (*fCutVector)[kDCUTE]; | |
173 | } | |
174 | else if ((processName == "MuIoni") || (processName == "IMuIoni")) { | |
175 | // !! Geant4 treats delta rays + continuous energy loss | |
176 | // within one process | |
177 | return (*fCutVector)[kDCUTM]; | |
178 | } | |
179 | else { | |
180 | return (*fCutVector)[kCUTELE]; | |
181 | } | |
182 | } | |
183 | ||
72095f7c | 184 | //_____________________________________________________________________________ |
0375ab40 | 185 | G4double TG4G3CutVector::GetMinEkineForHadron(const G4Track& track) const |
2817d3e2 | 186 | { |
187 | // Returns the cut value for charged hadron. | |
188 | // --- | |
189 | ||
190 | return (*fCutVector)[kCUTHAD]; | |
191 | } | |
192 | ||
72095f7c | 193 | //_____________________________________________________________________________ |
0375ab40 | 194 | G4double TG4G3CutVector::GetMinEkineForNeutralHadron(const G4Track& track) const |
2817d3e2 | 195 | { |
196 | // Returns the cut value for neutral hadron. | |
197 | // --- | |
198 | ||
199 | return (*fCutVector)[kCUTNEU]; | |
200 | } | |
201 | ||
72095f7c | 202 | //_____________________________________________________________________________ |
0375ab40 | 203 | G4double TG4G3CutVector::GetMinEkineForMuon(const G4Track& track) const |
2817d3e2 | 204 | { |
205 | // Returns the cut value for neutral muon. | |
206 | // --- | |
207 | ||
208 | return (*fCutVector)[kCUTMUO]; | |
209 | } | |
210 | ||
72095f7c | 211 | //_____________________________________________________________________________ |
0375ab40 | 212 | G4double TG4G3CutVector::GetMinEkineForOther(const G4Track& track) const |
2817d3e2 | 213 | { |
214 | // Returns 0. | |
215 | // --- | |
216 | ||
217 | return 0.; | |
218 | } | |
219 |