69 G4cout <<
"G4LEpp:ApplyYourself: incident particle: "
72 <<
", Px = " << Px/
GeV <<
" GeV/c"
73 <<
", Py = " << Py/
GeV <<
" GeV/c"
74 <<
", Pz = " << Pz/
GeV <<
" GeV/c" <<
G4endl;
76 <<
", kinetic energy = " << ek/
GeV <<
" GeV"
77 <<
", mass = " << E0/
GeV <<
" GeV"
78 <<
", charge = " << Q <<
G4endl;
90 E0 = std::sqrt(std::fabs(E02));
95 <<
", mass = " << E0/
GeV <<
" GeV"
96 <<
", charge = " << Q <<
G4endl;
100 if(cost > 1.0) { cost = 1.0; }
101 if(cost <-1.0) { cost =-1.0; }
102 G4double sint = std::sqrt((1.0 - cost)*(1.0 + cost));
112 G4double pseudoMass = std::sqrt(totalEnergy*totalEnergy - P*P);
119 G4double p = std::sqrt(px*px + py*py + pz*pz);
124 G4cout <<
" particle 1 momentum in CM " << px/
GeV
125 <<
" " << py/
GeV <<
" "
130 G4double pxnew = p*sint*std::cos(phi);
131 G4double pynew = p*sint*std::sin(phi);
135 if (px*px + py*py > 0) {
138 sint = (std::sqrt((1-cost)*(1+cost)) + std::sqrt(px*px+py*py)/p)/2;
140 if (std::fabs(px) > 0.000001*
GeV) ph = std::atan2(py,px);
143 px = (cost*cosp*pxnew - sinp*pynew + sint*cosp*pznew);
144 py = (cost*sinp*pxnew + cosp*pynew + sint*sinp*pznew);
145 pz = (-sint*pxnew + cost*pznew);
155 G4cout <<
" particle 1 momentum in CM " << px/
GeV <<
" " << py/
GeV <<
" "
166 G4double gammaCM = E1pM2/std::sqrt(E1pM2*E1pM2 - P*P);
169 G4cout <<
" betaCM " << betaCMx <<
" " << betaCMy <<
" "
170 << betaCMz <<
" " << betaCM <<
G4endl;
187 PA[4] = std::sqrt(M1*M1 + p*p);
189 G4double BETPA = BETA[1]*PA[1] + BETA[2]*PA[2] + BETA[3]*PA[3];
190 G4double BPGAM = (BETPA * BETA[4]/(BETA[4] + 1.) - PA[4]) * BETA[4];
192 PB[1] = PA[1] + BPGAM * BETA[1];
193 PB[2] = PA[2] + BPGAM * BETA[2];
194 PB[3] = PA[3] + BPGAM * BETA[3];
195 PB[4] = (PA[4] - BETPA) * BETA[4];
206 PA[4] = std::sqrt(M2*M2 + p*p);
208 BETPA = BETA[1]*PA[1] + BETA[2]*PA[2] + BETA[3]*PA[3];
209 BPGAM = (BETPA * BETA[4]/(BETA[4] + 1.) - PA[4]) * BETA[4];
211 PB[1] = PA[1] + BPGAM * BETA[1];
212 PB[2] = PA[2] + BPGAM * BETA[2];
213 PB[3] = PA[3] + BPGAM * BETA[3];
214 PB[4] = (PA[4] - BETPA) * BETA[4];
219 G4cout <<
" particle 1 momentum in LAB "
222 G4cout <<
" particle 2 momentum in LAB "
225 G4cout <<
" TOTAL momentum in LAB "
249 G4double ek = std::sqrt(plab*plab+nMass*nMass) - nMass;
259 G4int midBin = (je1 + je2)/2;
261 if (ek <
elab[midBin]) je2 = midBin;
264 while (je2 - je1 > 1);
275 dsig =
Sig[je2][0] -
Sig[je1][0];
277 b = Sig[je1][0] - rc*elab[je1];
284 G4int midBin = (ke1 + ke2)/2;
285 dsig = Sig[je2][midBin] - Sig[je1][midBin];
287 b = Sig[je1][midBin] - rc*elab[je1];
301 while (ke2 - ke1 > 1);
303 dsig = sigint2 - sigint1;
305 b = ke1 - rc*sigint1;
309 G4double t = 0.5*plab*plab*(1 - std::cos(theta));
static const G4float Sig[NENERGY][NANGLE]
static const double halfpi
void SetMomentum(const G4ThreeVector &momentum)
G4double GetKineticEnergy() const
static const G4float elab[NENERGY]
CLHEP::Hep3Vector G4ThreeVector
G4double GetTotalEnergy() const
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
G4ParticleDefinition * GetDefinition() const
G4DynamicParticle * ReturnTargetParticle() const
const G4String & GetParticleName() const
G4double GetTotalMomentum() const
void SetMinEnergy(G4double anEnergy)
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
const G4ParticleDefinition * GetDefinition() const
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
static const double twopi
static G4Proton * Proton()
const G4LorentzVector & Get4Momentum() const
void SetEnergyChange(G4double anEnergy)
G4double GetPDGMass() const
void SetMaxEnergy(const G4double anEnergy)
G4HadFinalState theParticleChange
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
G4double GetPDGCharge() const
void SetMomentumChange(const G4ThreeVector &aV)
G4ThreeVector GetMomentum() const
G4double GetTotalMomentum() const
G4double GetTotalEnergy() const