00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 #include "G4CascadParticle.hh"
00037 #include "G4ios.hh"
00038 #include <cmath>
00039
00040
00041
00042
00043 G4CascadParticle::G4CascadParticle()
00044 : verboseLevel(0), current_zone(-1), current_path(-1.), movingIn(false),
00045 reflectionCounter(0), reflected(false), generation(-1) {
00046 if (verboseLevel > 3) {
00047 G4cout << " >>> G4CascadParticle::G4CascadParticle" << G4endl;
00048 }
00049 }
00050
00051
00052
00053 void G4CascadParticle::fill(const G4InuclElementaryParticle& particle,
00054 const G4ThreeVector& pos, G4int izone,
00055 G4double cpath, G4int gen) {
00056 if (verboseLevel > 3) G4cout << " >>> G4CascadParticle::fill" << G4endl;
00057
00058 theParticle = particle;
00059 position = pos;
00060 current_zone = izone;
00061 current_path = cpath;
00062 movingIn = true;
00063 reflectionCounter = 0;
00064 reflected = false;
00065 generation = gen;
00066 }
00067
00068
00069 G4double G4CascadParticle::getPathToTheNextZone(G4double rz_in,
00070 G4double rz_out) {
00071 if (verboseLevel > 3) {
00072 G4cout << " >>> G4CascadParticle::getPathToTheNextZone rz_in " << rz_in
00073 << " rz_out " << rz_out << G4endl;
00074 }
00075
00076 const G4LorentzVector& mom = getMomentum();
00077
00078 G4double path = -1.0;
00079 G4double rp = mom.vect().dot(position);
00080 G4double rr = position.mag2();
00081 G4double pp = mom.vect().mag2();
00082
00083 if (std::abs(pp) < 1e-9) {
00084 if (verboseLevel > 3) G4cout << " at rest; path length is zero" << G4endl;
00085 return 0.;
00086 }
00087
00088 G4double ra = rr - rp * rp / pp;
00089 pp = std::sqrt(pp);
00090 G4double ds;
00091 G4double d2;
00092
00093 if (verboseLevel > 3) {
00094 G4cout << " current_zone " << current_zone << " rr " << rr
00095 << " rp " << rp << " pp " << pp << " ra " << ra << G4endl;
00096 }
00097
00098 if (current_zone == 0 || rp > 0.0) {
00099 d2 = rz_out * rz_out - ra;
00100 if (d2 > 0.0) {
00101 ds = 1.0;
00102 movingIn = false;
00103 } else {
00104 d2 = rz_in * rz_in - ra;
00105 ds = -1.0;
00106 movingIn = true;
00107 }
00108 } else {
00109 d2 = rz_in * rz_in - ra;
00110 if (d2 > 0.0) {
00111 ds = -1.0;
00112 movingIn = true;
00113 } else {
00114 d2 = rz_out * rz_out - ra;
00115 ds = 1.0;
00116 movingIn = false;
00117 }
00118 }
00119
00120 if (verboseLevel > 3) G4cout << " ds " << ds << " d2 " << d2 << G4endl;
00121
00122 if (d2 < 0.0 && d2 > -1e-6) d2 = 0.;
00123
00124 if (d2 > 0.0) path = ds * std::sqrt(d2) - rp / pp;
00125
00126 return path;
00127 }
00128
00129 void G4CascadParticle::propagateAlongThePath(G4double path) {
00130 if (verboseLevel > 3) {
00131 G4cout << " >>> G4CascadParticle::propagateAlongThePath" << G4endl;
00132 }
00133
00134 position += getMomentum().vect().unit()*path;
00135 }
00136
00137
00138
00139
00140 std::ostream& operator<<(std::ostream& os, const G4CascadParticle& part) {
00141 part.print(os);
00142 return os;
00143 }
00144
00145 void G4CascadParticle::print(std::ostream& os) const {
00146 os << theParticle << G4endl
00147 << " zone " << current_zone << " current_path " << current_path
00148 << " reflectionCounter " << reflectionCounter << G4endl
00149 << " x " << position.x() << " y " << position.y()
00150 << " z " << position.z() << G4endl;
00151 }
00152