ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4CascadParticle.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4CascadParticle.cc
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 // 20100112 M. Kelsey -- Remove G4CascadeMomentum, use G4LorentzVector directly
28 // 20100114 M. Kelsey -- Replace vector<G4Double> position with G4ThreeVector
29 // 20101012 M. Kelsey -- Check for negative d2 in getPathToTheNextZone()
30 // 20110806 M. Kelsey -- Add fill() function to replicate ctor/op=() action
31 // 20110922 M. Kelsey -- Follow G4InuclParticle::print(ostream&) migration,
32 // Add stream argument to print(), add operator<<().
33 // 20111017 M. Kelsey -- Add check for zero momentum in path calculation.
34 // 20130221 M. Kelsey -- Set verbosity using global envvar, for both ctors.
35 // 20130304 M. Kelsey -- Add index data member, for use with G4CascadeHistory,
36 // and explicit copy operations
37 
38 #include "G4CascadParticle.hh"
39 #include "G4CascadeParameters.hh"
40 #include "G4ios.hh"
41 #include <cmath>
42 
43 
44 // Default constructor produces non-functional object
45 
47  : verboseLevel(G4CascadeParameters::verbose()),
48  current_zone(-1), current_path(-1.), movingIn(false),
49  reflectionCounter(0), reflected(false), generation(-1), historyId(-1) {
50  if (verboseLevel > 3) {
51  G4cout << " >>> G4CascadParticle::G4CascadParticle" << G4endl;
52  }
53 }
54 
57  const G4ThreeVector& pos, G4int izone, G4double cpath,
58  G4int gen)
59  : verboseLevel(G4CascadeParameters::verbose()),
60  theParticle(particle), position(pos),
61  current_zone(izone), current_path(cpath), movingIn(true),
62  reflectionCounter(0), reflected(false), generation(gen), historyId(-1) {
63  if (verboseLevel > 3) {
64  G4cout << " >>> G4CascadParticle::G4CascadParticle "
65  << particle.getDefinition()->GetParticleName()
66  << " @ " << pos << G4endl;
67  }
68 }
69 
70 // Copy contents, including history information
71 
73  if (&cpart != this) { // Avoid unnecessary work
74  verboseLevel = cpart.verboseLevel;
75  theParticle = cpart.theParticle;
76  position = cpart.position;
77  current_zone = cpart.current_zone;
78  current_path = cpart.current_path;
79  movingIn = cpart.movingIn;
81  reflected = cpart.reflected;
82  generation = cpart.generation;
83  historyId = cpart.historyId;
84  }
85 
86  return *this;
87 }
88 
89 // Analogue to operator=() to support filling vectors w/o temporaries
90 
92  const G4ThreeVector& pos, G4int izone,
93  G4double cpath, G4int gen) {
94  if (verboseLevel > 3) G4cout << " >>> G4CascadParticle::fill" << G4endl;
95 
97  position = pos;
98  current_zone = izone;
99  current_path = cpath;
100  movingIn = true;
101  reflectionCounter = 0;
102  reflected = false;
103  generation = gen;
104  historyId = -1; // New particle, no history entry yet
105 }
106 
107 
109  G4double rz_out) {
110  if (verboseLevel > 3) {
111  G4cout << " >>> G4CascadParticle::getPathToTheNextZone rz_in " << rz_in
112  << " rz_out " << rz_out << G4endl;
113  }
114 
115  const G4LorentzVector& mom = getMomentum();
116 
117  G4double path = -1.0;
118  G4double rp = mom.vect().dot(position);
119  G4double rr = position.mag2();
120  G4double pp = mom.vect().mag2();
121 
122  if (std::abs(pp) < 1e-9) { // Cut-off for "at rest" is 1 eV momentum
123  if (verboseLevel > 3) G4cout << " at rest; path length is zero" << G4endl;
124 
125  if (current_zone == 0) movingIn = false; // Allow 'stuck' to escape
126  return 0.;
127  }
128 
129  G4double ra = rr - rp * rp / pp;
130  pp = std::sqrt(pp);
131  G4double ds;
132  G4double d2;
133 
134  if (verboseLevel > 3) {
135  G4cout << " current_zone " << current_zone << " rr " << rr
136  << " rp " << rp << " pp " << pp << " ra " << ra << G4endl;
137  }
138 
139  if (current_zone == 0 || rp > 0.0) {
140  d2 = rz_out * rz_out - ra;
141  if (d2 > 0.0) {
142  ds = 1.0;
143  movingIn = false;
144  } else {
145  d2 = rz_in * rz_in - ra;
146  ds = -1.0;
147  movingIn = true;
148  }
149  } else {
150  d2 = rz_in * rz_in - ra;
151  if (d2 > 0.0) {
152  ds = -1.0;
153  movingIn = true;
154  } else {
155  d2 = rz_out * rz_out - ra;
156  ds = 1.0;
157  movingIn = false;
158  }
159  }
160 
161  if (verboseLevel > 3) G4cout << " ds " << ds << " d2 " << d2 << G4endl;
162 
163  if (d2 < 0.0 && d2 > -1e-6) d2 = 0.; // Account for round-off
164 
165  if (d2 > 0.0) path = ds * std::sqrt(d2) - rp / pp; // Avoid FPE failure
166 
167  return path;
168 }
169 
171  if (verboseLevel > 3) {
172  G4cout << " >>> G4CascadParticle::propagateAlongThePath" << G4endl;
173  }
174 
175  position += getMomentum().vect().unit()*path;
176 }
177 
178 
179 // Proper stream output (just calls print())
180 
181 std::ostream& operator<<(std::ostream& os, const G4CascadParticle& part) {
182  part.print(os);
183  return os;
184 }
185 
186 void G4CascadParticle::print(std::ostream& os) const {
187  os << " pos " << position << " zone " << current_zone
188  << " current_path " << current_path
189  << " reflectionCounter " << reflectionCounter << G4endl
190  << theParticle << G4endl;
191 }
192