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

Last change on this file since 826 was 826, checked in by garnier, 16 years ago

import all except CVS

File size: 12.5 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.20 2007/03/25 22:54:52 kurasige Exp $
28// GEANT4 tag $Name:  $
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   theFirstStepInVolume(false),
56   theLastStepInVolume(false),
57   verboseLevel(1),
58   theParentWeight(1.0),
59   fSetSecondaryWeightByProcess(false), 
60   fSetParentWeightByProcess(true) 
61{
62   debugFlag = false;
63#ifdef G4VERBOSE
64  // activate CHeckIt if in VERBOSE mode
65  debugFlag = true;
66#endif
67   theListOfSecondaries = new G4TrackFastVector();
68}
69
70G4VParticleChange::~G4VParticleChange() {
71  // check if tracks still exist in theListOfSecondaries
72  if (theNumberOfSecondaries>0) {
73#ifdef G4VERBOSE
74    if (verboseLevel>0) {
75      G4cerr << "G4VParticleChange::~G4VParticleChange() Warning  ";
76      G4cerr << "theListOfSecondaries is not empty ";
77    }
78#endif
79    for (G4int index= 0; index<theNumberOfSecondaries; index++){
80      if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
81    }
82  }
83  delete theListOfSecondaries; 
84}
85
86// copy and assignment operators are implemented as "shallow copy"
87G4VParticleChange::G4VParticleChange(const G4VParticleChange &right):
88   theNumberOfSecondaries(0),
89   theSizeOftheListOfSecondaries(G4TrackFastVectorSize),
90   theStatusChange(fAlive),
91   theSteppingControlFlag(NormalCondition),     
92   theLocalEnergyDeposit(0.0),
93   theNonIonizingEnergyDeposit(0.0),
94   theFirstStepInVolume(false),
95   theLastStepInVolume(false),
96   verboseLevel(1),
97   theParentWeight(1.0),
98   fSetSecondaryWeightByProcess(false)
99{
100   debugFlag = false;
101#ifdef G4VERBOSE
102  // activate CHeckIt if in VERBOSE mode
103  debugFlag = true;
104#endif
105
106  theListOfSecondaries = right.theListOfSecondaries;
107  theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries;
108  theNumberOfSecondaries = right.theNumberOfSecondaries;
109  theStatusChange = right.theStatusChange;
110  theTrueStepLength = right.theTrueStepLength;
111  theLocalEnergyDeposit = right.theLocalEnergyDeposit;
112  theNonIonizingEnergyDeposit = right.theNonIonizingEnergyDeposit;
113  theSteppingControlFlag = right.theSteppingControlFlag;
114  fSetParentWeightByProcess = right.fSetParentWeightByProcess;   
115
116  theFirstStepInVolume = right.theFirstStepInVolume;
117  theLastStepInVolume =  right.theLastStepInVolume;
118}
119
120
121G4VParticleChange & G4VParticleChange::operator=(const G4VParticleChange &right)
122{
123   debugFlag = false;
124#ifdef G4VERBOSE
125  // activate CHeckIt if in VERBOSE mode
126  debugFlag = true;
127#endif
128   if (this != &right)
129   {
130      theListOfSecondaries = right.theListOfSecondaries;
131      theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries;
132      theNumberOfSecondaries = right.theNumberOfSecondaries;
133      theStatusChange = right.theStatusChange;
134      theTrueStepLength = right.theTrueStepLength;
135      theLocalEnergyDeposit = right.theLocalEnergyDeposit;
136      theNonIonizingEnergyDeposit = right.theNonIonizingEnergyDeposit;
137      theSteppingControlFlag = right.theSteppingControlFlag;
138      fSetParentWeightByProcess = right.fSetParentWeightByProcess;   
139
140      theFirstStepInVolume = right.theFirstStepInVolume;
141      theLastStepInVolume =  right.theLastStepInVolume;
142   }
143   return *this;
144}
145
146G4bool G4VParticleChange::operator==(const G4VParticleChange &right) const
147{
148   return (this == (G4VParticleChange *) &right);
149}
150
151
152G4bool G4VParticleChange::operator!=(const G4VParticleChange &right) const
153{
154   return (this != (G4VParticleChange *) &right);
155}
156
157//----------------------------------------------------------------
158// methods for printing messages 
159//
160 
161void G4VParticleChange::DumpInfo() const
162{
163
164// Show header
165  G4cout.precision(3);
166  G4cout << "      -----------------------------------------------" 
167       << G4endl;
168  G4cout << "        G4ParticleChange Information  " << std::setw(20) << G4endl;
169  G4cout << "      -----------------------------------------------" 
170       << G4endl;
171
172  G4cout << "        # of 2ndaries       : " 
173       << std::setw(20) << theNumberOfSecondaries
174       << G4endl;
175
176  if (theNumberOfSecondaries >0) {
177    G4cout << "        Pointer to 2ndaries : " 
178         << std::setw(20) << GetSecondary(0)
179         << G4endl;
180    G4cout << "        (Showed only 1st one)"
181         << G4endl;
182  }
183  G4cout << "      -----------------------------------------------" 
184       << G4endl;
185
186  G4cout << "        Energy Deposit (MeV): " 
187       << std::setw(20) << theLocalEnergyDeposit/MeV
188       << G4endl;
189
190  G4cout << "        Non-ionizing Energy Deposit (MeV): " 
191       << std::setw(20) << theNonIonizingEnergyDeposit/MeV
192       << G4endl;
193
194  G4cout << "        Track Status        : " 
195       << std::setw(20);
196       if( theStatusChange == fAlive ){
197         G4cout << " Alive";
198       } else if( theStatusChange == fStopButAlive ){
199           G4cout << " StopButAlive";
200       } else if( theStatusChange == fStopAndKill ){
201           G4cout << " StopAndKill";
202       } else if( theStatusChange  == fKillTrackAndSecondaries ){
203           G4cout << " KillTrackAndSecondaries";
204       } else if( theStatusChange  == fSuspend ){
205           G4cout << " Suspend";
206       } else if( theStatusChange == fPostponeToNextEvent ){
207           G4cout << " PostponeToNextEvent";
208       }
209       G4cout << G4endl;
210  G4cout << "        True Path Length (mm) : " 
211       << std::setw(20) << theTrueStepLength/mm
212       << G4endl;
213  G4cout << "        Stepping Control      : " 
214       << std::setw(20) << theSteppingControlFlag
215       << G4endl;   
216  if (theFirstStepInVolume) {
217    G4cout << "    First Step In the voulme  : "  << G4endl;
218  }
219  if (theLastStepInVolume) {
220    G4cout << "    Last Step In the voulme  : "  << G4endl;
221  }
222}
223
224G4bool G4VParticleChange::CheckIt(const G4Track& )
225{
226
227  G4bool    exitWithError = false;
228  G4double  accuracy;
229
230  // Energy deposit should not be negative
231  G4bool itsOKforEnergy = true;
232  accuracy = -1.0*theLocalEnergyDeposit/MeV;
233  if (accuracy > accuracyForWarning) {
234#ifdef G4VERBOSE
235    G4cout << "  G4VParticleChange::CheckIt    : ";
236    G4cout << "the energy deposit  is negative  !!" << G4endl;
237    G4cout << "  Difference:  " << accuracy  << "[MeV] " <<G4endl;
238#endif
239    itsOKforEnergy = false;
240    if (accuracy > accuracyForException) exitWithError = true;
241  }
242 
243  // true path length should not be negative
244  G4bool itsOKforStepLength = true;
245  accuracy = -1.0*theTrueStepLength/mm;
246  if (accuracy > accuracyForWarning) {
247#ifdef G4VERBOSE
248    G4cout << "  G4VParticleChange::CheckIt    : ";
249    G4cout << "the true step length is negative  !!" << G4endl;
250    G4cout << "  Difference:  " << accuracy  << "[MeV] " <<G4endl;
251#endif
252    itsOKforStepLength = false;
253    if (accuracy > accuracyForException) exitWithError = true;
254  }
255
256  G4bool itsOK = itsOKforStepLength && itsOKforEnergy ;
257  // dump out information of this particle change
258#ifdef G4VERBOSE
259  if (! itsOK ){
260    G4cout << " G4VParticleChange::CheckIt " <<G4endl;
261    DumpInfo();
262  }
263#endif
264
265  // Exit with error
266  if (exitWithError) {
267    G4Exception("G4VParticleChange::CheckIt",
268                "100",
269                EventMustBeAborted,
270                "step length and/or energy deposit was  illegal");
271  }
272
273  // correction
274  if ( !itsOKforStepLength ) {
275    theTrueStepLength =  (1.e-12)*mm;
276  } 
277  if ( !itsOKforEnergy ) {
278    theLocalEnergyDeposit = 0.0;
279  }
280  return itsOK;
281}
282
283G4bool G4VParticleChange::CheckSecondary(G4Track& aTrack)
284{
285  G4bool    exitWithError = false;
286  G4double  accuracy;
287
288  // MomentumDirection should be unit vector
289  G4bool itsOKforMomentum = true; 
290  accuracy = std::abs((aTrack.GetMomentumDirection()).mag2()-1.0);
291  if (accuracy > accuracyForWarning) {
292#ifdef G4VERBOSE
293    G4cout << " G4VParticleChange::CheckSecondary  :   ";
294    G4cout << "the Momentum direction is not unit vector !!" << G4endl;
295    G4cout << "  Difference:  " << accuracy << G4endl;
296#endif
297    itsOKforMomentum = false;
298    if (accuracy > accuracyForException) exitWithError = true;
299  }
300 
301  // Kinetic Energy should not be negative
302  G4bool itsOKforEnergy;
303  accuracy = -1.0*(aTrack.GetKineticEnergy())/MeV;
304  if (accuracy < accuracyForWarning) {
305    itsOKforEnergy = true;
306  } else {
307#ifdef G4VERBOSE
308    G4cout << " G4VParticleChange::CheckSecondary  :   ";
309    G4cout << "the kinetic energy is negative  !!" << G4endl;
310    G4cout << "  Difference:  " << accuracy  << "[MeV] " <<G4endl;
311#endif
312    itsOKforEnergy = false;
313    if (accuracy < accuracyForException) { exitWithError = false;}
314    else { exitWithError = true; }
315  }
316  G4bool itsOKforProperTime = true;
317  //accuracy = (aTrack.GetProperTime())/ns;
318  //  if (accuracy > accuracyForWarning) {
319#ifdef G4VERBOSE
320  //  G4cout << "  G4VParticleChange::CheckSecondary    : ";
321  //  G4cout << "the proper time goes back  !!" << G4endl;
322  //  G4cout << "  Difference:  " << accuracy  << "[ns] " <<G4endl;
323#endif
324  //  itsOKforProperTime = false;
325  //  if (accuracy > accuracyForException) exitWithError = true;
326  //}
327 
328  // Exit with error
329  if (exitWithError) {
330    G4Exception("G4VParticleChange::CheckSecondary",
331                "10",
332                EventMustBeAborted, 
333                "momentum, energy and/or proper time was illegal");
334  }
335
336  G4bool itsOK = itsOKforMomentum && itsOKforEnergy && itsOKforProperTime;
337
338  //correction
339  if (!itsOKforMomentum) {
340    G4double vmag = (aTrack.GetMomentumDirection()).mag();
341    aTrack.SetMomentumDirection((1./vmag)*aTrack.GetMomentumDirection());
342  }
343  //if (!itsOKforProperTime) {
344  //  aTrack.SetProperTime(0.0);
345  //}
346  if (!itsOKforEnergy) {
347    aTrack.SetKineticEnergy(0.0);
348  }
349
350  return itsOK;
351}
352
353
354G4double G4VParticleChange::GetAccuracyForWarning() const
355{
356  return accuracyForWarning;
357}
358
359G4double G4VParticleChange::GetAccuracyForException() const
360{
361  return accuracyForException;
362}
363
364void G4VParticleChange::AddSecondary(G4Track *aTrack)
365{
366  if (debugFlag) CheckSecondary(*aTrack);
367
368  if (!fSetSecondaryWeightByProcess){
369    // pass the weight of parent track
370    aTrack->SetWeight(theParentWeight);
371  }
372
373  // add a secondary after size check
374  if (theSizeOftheListOfSecondaries > theNumberOfSecondaries) {
375    theListOfSecondaries->SetElement(theNumberOfSecondaries, aTrack);
376    theNumberOfSecondaries++;
377  } else {
378#ifdef G4VERBOSE
379    if (verboseLevel>0) {
380      G4cerr << "G4VParticleChange::AddSecondary() Warning  ";
381      G4cerr << "theListOfSecondaries is full !! " << G4endl;
382      G4cerr << " The object will not be added in theListOfSecondaries" << G4endl;
383    }
384#endif
385  }
386}
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
Note: See TracBrowser for help on using the repository browser.