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

Last change on this file since 1340 was 1340, checked in by garnier, 14 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.