source: trunk/source/track/src/G4VParticleChange.cc@ 1345

Last change on this file since 1345 was 1340, checked in by garnier, 15 years ago

update ti head
File size: 12.6 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26//
27// $Id: G4VParticleChange.cc,v 1.22 2010/07/21 09:30:15 gcosmo Exp $
28// GEANT4 tag $Name: track-V09-03-09 $
29//
30//
31// --------------------------------------------------------------
32// GEANT 4 class implementation file
33//
34//
35// ------------------------------------------------------------
36// Implemented for the new scheme 23 Mar. 1998 H.Kurahige
37// --------------------------------------------------------------
38
39#include "G4VParticleChange.hh"
40#include "G4Track.hh"
41#include "G4Step.hh"
42#include "G4TrackFastVector.hh"
43#include "G4ExceptionSeverity.hh"
44
45const G4double G4VParticleChange::accuracyForWarning = 1.0e-9;
46const G4double G4VParticleChange::accuracyForException = 0.001;
47
48G4VParticleChange::G4VParticleChange():
49 theNumberOfSecondaries(0),
50 theSizeOftheListOfSecondaries(G4TrackFastVectorSize),
51 theStatusChange(fAlive),
52 theSteppingControlFlag(NormalCondition),
53 theLocalEnergyDeposit(0.0),
54 theNonIonizingEnergyDeposit(0.0),
55 theTrueStepLength(0.0),
56 theFirstStepInVolume(false),
57 theLastStepInVolume(false),
58 verboseLevel(1),
59 theParentWeight(1.0),
60 fSetSecondaryWeightByProcess(false),
61 fSetParentWeightByProcess(true)
62{
63 debugFlag = false;
64#ifdef G4VERBOSE
65 // activate CHeckIt if in VERBOSE mode
66 debugFlag = true;
67#endif
68 theListOfSecondaries = new G4TrackFastVector();
69}
70
71G4VParticleChange::~G4VParticleChange() {
72 // check if tracks still exist in theListOfSecondaries
73 if (theNumberOfSecondaries>0) {
74#ifdef G4VERBOSE
75 if (verboseLevel>0) {
76 G4cerr << "G4VParticleChange::~G4VParticleChange() Warning ";
77 G4cerr << "theListOfSecondaries is not empty ";
78 }
79#endif
80 for (G4int index= 0; index<theNumberOfSecondaries; index++){
81 if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
82 }
83 }
84 delete theListOfSecondaries;
85}
86
87// copy and assignment operators are implemented as "shallow copy"
88G4VParticleChange::G4VParticleChange(const G4VParticleChange &right):
89 theNumberOfSecondaries(0),
90 theSizeOftheListOfSecondaries(G4TrackFastVectorSize),
91 theStatusChange(fAlive),
92 theSteppingControlFlag(NormalCondition),
93 theLocalEnergyDeposit(0.0),
94 theNonIonizingEnergyDeposit(0.0),
95 theFirstStepInVolume(false),
96 theLastStepInVolume(false),
97 verboseLevel(1),
98 theParentWeight(1.0),
99 fSetSecondaryWeightByProcess(false)
100{
101 debugFlag = false;
102#ifdef G4VERBOSE
103 // activate CHeckIt if in VERBOSE mode
104 debugFlag = true;
105#endif
106
107 theListOfSecondaries = right.theListOfSecondaries;
108 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries;
109 theNumberOfSecondaries = right.theNumberOfSecondaries;
110 theStatusChange = right.theStatusChange;
111 theTrueStepLength = right.theTrueStepLength;
112 theLocalEnergyDeposit = right.theLocalEnergyDeposit;
113 theNonIonizingEnergyDeposit = right.theNonIonizingEnergyDeposit;
114 theSteppingControlFlag = right.theSteppingControlFlag;
115 fSetParentWeightByProcess = right.fSetParentWeightByProcess;
116
117 theFirstStepInVolume = right.theFirstStepInVolume;
118 theLastStepInVolume = right.theLastStepInVolume;
119}
120
121
122G4VParticleChange & G4VParticleChange::operator=(const G4VParticleChange &right)
123{
124 debugFlag = false;
125#ifdef G4VERBOSE
126 // activate CHeckIt if in VERBOSE mode
127 debugFlag = true;
128#endif
129 if (this != &right)
130 {
131 theListOfSecondaries = right.theListOfSecondaries;
132 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries;
133 theNumberOfSecondaries = right.theNumberOfSecondaries;
134 theStatusChange = right.theStatusChange;
135 theTrueStepLength = right.theTrueStepLength;
136 theLocalEnergyDeposit = right.theLocalEnergyDeposit;
137 theNonIonizingEnergyDeposit = right.theNonIonizingEnergyDeposit;
138 theSteppingControlFlag = right.theSteppingControlFlag;
139 fSetParentWeightByProcess = right.fSetParentWeightByProcess;
140
141 theFirstStepInVolume = right.theFirstStepInVolume;
142 theLastStepInVolume = right.theLastStepInVolume;
143 }
144 return *this;
145}
146
147G4bool G4VParticleChange::operator==(const G4VParticleChange &right) const
148{
149 return (this == (G4VParticleChange *) &right);
150}
151
152
153G4bool G4VParticleChange::operator!=(const G4VParticleChange &right) const
154{
155 return (this != (G4VParticleChange *) &right);
156}
157
158//----------------------------------------------------------------
159// methods for printing messages
160//
161
162void G4VParticleChange::DumpInfo() const
163{
164
165// Show header
166 G4int olprc = G4cout.precision(3);
167 G4cout << " -----------------------------------------------"
168 << G4endl;
169 G4cout << " G4ParticleChange Information " << std::setw(20) << G4endl;
170 G4cout << " -----------------------------------------------"
171 << G4endl;
172
173 G4cout << " # of 2ndaries : "
174 << std::setw(20) << theNumberOfSecondaries
175 << G4endl;
176
177 if (theNumberOfSecondaries >0) {
178 G4cout << " Pointer to 2ndaries : "
179 << std::setw(20) << GetSecondary(0)
180 << G4endl;
181 G4cout << " (Showed only 1st one)"
182 << G4endl;
183 }
184 G4cout << " -----------------------------------------------"
185 << G4endl;
186
187 G4cout << " Energy Deposit (MeV): "
188 << std::setw(20) << theLocalEnergyDeposit/MeV
189 << G4endl;
190
191 G4cout << " Non-ionizing Energy Deposit (MeV): "
192 << std::setw(20) << theNonIonizingEnergyDeposit/MeV
193 << G4endl;
194
195 G4cout << " Track Status : "
196 << std::setw(20);
197 if( theStatusChange == fAlive ){
198 G4cout << " Alive";
199 } else if( theStatusChange == fStopButAlive ){
200 G4cout << " StopButAlive";
201 } else if( theStatusChange == fStopAndKill ){
202 G4cout << " StopAndKill";
203 } else if( theStatusChange == fKillTrackAndSecondaries ){
204 G4cout << " KillTrackAndSecondaries";
205 } else if( theStatusChange == fSuspend ){
206 G4cout << " Suspend";
207 } else if( theStatusChange == fPostponeToNextEvent ){
208 G4cout << " PostponeToNextEvent";
209 }
210 G4cout << G4endl;
211 G4cout << " True Path Length (mm) : "
212 << std::setw(20) << theTrueStepLength/mm
213 << G4endl;
214 G4cout << " Stepping Control : "
215 << std::setw(20) << theSteppingControlFlag
216 << G4endl;
217 if (theFirstStepInVolume) {
218 G4cout << " First Step In the voulme : " << G4endl;
219 }
220 if (theLastStepInVolume) {
221 G4cout << " Last Step In the voulme : " << G4endl;
222 }
223 G4cout.precision(olprc);
224}
225
226G4bool G4VParticleChange::CheckIt(const G4Track& )
227{
228
229 G4bool exitWithError = false;
230 G4double accuracy;
231
232 // Energy deposit should not be negative
233 G4bool itsOKforEnergy = true;
234 accuracy = -1.0*theLocalEnergyDeposit/MeV;
235 if (accuracy > accuracyForWarning) {
236#ifdef G4VERBOSE
237 G4cout << " G4VParticleChange::CheckIt : ";
238 G4cout << "the energy deposit is negative !!" << G4endl;
239 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl;
240#endif
241 itsOKforEnergy = false;
242 if (accuracy > accuracyForException) exitWithError = true;
243 }
244
245 // true path length should not be negative
246 G4bool itsOKforStepLength = true;
247 accuracy = -1.0*theTrueStepLength/mm;
248 if (accuracy > accuracyForWarning) {
249#ifdef G4VERBOSE
250 G4cout << " G4VParticleChange::CheckIt : ";
251 G4cout << "the true step length is negative !!" << G4endl;
252 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl;
253#endif
254 itsOKforStepLength = false;
255 if (accuracy > accuracyForException) exitWithError = true;
256 }
257
258 G4bool itsOK = itsOKforStepLength && itsOKforEnergy ;
259 // dump out information of this particle change
260#ifdef G4VERBOSE
261 if (! itsOK ){
262 G4cout << " G4VParticleChange::CheckIt " <<G4endl;
263 DumpInfo();
264 }
265#endif
266
267 // Exit with error
268 if (exitWithError) {
269 G4Exception("G4VParticleChange::CheckIt",
270 "100",
271 EventMustBeAborted,
272 "step length and/or energy deposit was illegal");
273 }
274
275 // correction
276 if ( !itsOKforStepLength ) {
277 theTrueStepLength = (1.e-12)*mm;
278 }
279 if ( !itsOKforEnergy ) {
280 theLocalEnergyDeposit = 0.0;
281 }
282 return itsOK;
283}
284
285G4bool G4VParticleChange::CheckSecondary(G4Track& aTrack)
286{
287 G4bool exitWithError = false;
288 G4double accuracy;
289
290 // MomentumDirection should be unit vector
291 G4bool itsOKforMomentum = true;
292 accuracy = std::fabs((aTrack.GetMomentumDirection()).mag2()-1.0);
293 if (accuracy > accuracyForWarning) {
294#ifdef G4VERBOSE
295 G4cout << " G4VParticleChange::CheckSecondary : ";
296 G4cout << "the Momentum direction is not unit vector !!" << G4endl;
297 G4cout << " Difference: " << accuracy << G4endl;
298#endif
299 itsOKforMomentum = false;
300 if (accuracy > accuracyForException) exitWithError = true;
301 }
302
303 // Kinetic Energy should not be negative
304 G4bool itsOKforEnergy;
305 accuracy = -1.0*(aTrack.GetKineticEnergy())/MeV;
306 if (accuracy < accuracyForWarning) {
307 itsOKforEnergy = true;
308 } else {
309#ifdef G4VERBOSE
310 G4cout << " G4VParticleChange::CheckSecondary : ";
311 G4cout << "the kinetic energy is negative !!" << G4endl;
312 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl;
313#endif
314 itsOKforEnergy = false;
315 if (accuracy < accuracyForException) { exitWithError = false;}
316 else { exitWithError = true; }
317 }
318 G4bool itsOKforProperTime = true;
319 //accuracy = (aTrack.GetProperTime())/ns;
320 // if (accuracy > accuracyForWarning) {
321#ifdef G4VERBOSE
322 // G4cout << " G4VParticleChange::CheckSecondary : ";
323 // G4cout << "the proper time goes back !!" << G4endl;
324 // G4cout << " Difference: " << accuracy << "[ns] " <<G4endl;
325#endif
326 // itsOKforProperTime = false;
327 // if (accuracy > accuracyForException) exitWithError = true;
328 //}
329
330 // Exit with error
331 if (exitWithError) {
332 G4Exception("G4VParticleChange::CheckSecondary",
333 "10",
334 EventMustBeAborted,
335 "momentum, energy and/or proper time was illegal");
336 }
337
338 G4bool itsOK = itsOKforMomentum && itsOKforEnergy && itsOKforProperTime;
339
340 //correction
341 if (!itsOKforMomentum) {
342 G4double vmag = (aTrack.GetMomentumDirection()).mag();
343 aTrack.SetMomentumDirection((1./vmag)*aTrack.GetMomentumDirection());
344 }
345 //if (!itsOKforProperTime) {
346 // aTrack.SetProperTime(0.0);
347 //}
348 if (!itsOKforEnergy) {
349 aTrack.SetKineticEnergy(0.0);
350 }
351
352 return itsOK;
353}
354
355
356G4double G4VParticleChange::GetAccuracyForWarning() const
357{
358 return accuracyForWarning;
359}
360
361G4double G4VParticleChange::GetAccuracyForException() const
362{
363 return accuracyForException;
364}
365
366void G4VParticleChange::AddSecondary(G4Track *aTrack)
367{
368 if (debugFlag) CheckSecondary(*aTrack);
369
370 if (!fSetSecondaryWeightByProcess){
371 // pass the weight of parent track
372 aTrack->SetWeight(theParentWeight);
373 }
374
375 // add a secondary after size check
376 if (theSizeOftheListOfSecondaries > theNumberOfSecondaries) {
377 theListOfSecondaries->SetElement(theNumberOfSecondaries, aTrack);
378 theNumberOfSecondaries++;
379 } else {
380#ifdef G4VERBOSE
381 if (verboseLevel>0) {
382 G4cerr << "G4VParticleChange::AddSecondary() Warning ";
383 G4cerr << "theListOfSecondaries is full !! " << G4endl;
384 G4cerr << " The object will not be added in theListOfSecondaries" << G4endl;
385 }
386#endif
387 }
388}
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
Note: See TracBrowser for help on using the repository browser.