80 G4cout <<
"G4LEnp:ApplyYourself: incident particle: "
83 <<
", Px = " << Px/
GeV <<
" GeV/c"
84 <<
", Py = " << Py/
GeV <<
" GeV/c"
85 <<
", Pz = " << Pz/
GeV <<
" GeV/c" <<
G4endl;
87 <<
", kinetic energy = " << ek/
GeV <<
" GeV"
88 <<
", mass = " << E0/
GeV <<
" GeV"
89 <<
", charge = " << Q <<
G4endl;
101 E0 = std::sqrt(std::abs(E02));
102 if (E02 < 0)E0 *= -1;
106 <<
", mass = " << E0/
GeV <<
" GeV"
107 <<
", charge = " << Q <<
G4endl;
116 G4int midBin = (je1 + je2)/2;
117 if (ek <
elab[midBin])
121 }
while (je2 - je1 > 1);
131 G4double b = sig[je1][0] - rc*elab[je1];
137 << ke1 <<
" " << ke2 <<
" "
138 << sigint1 <<
" " << sigint2 <<
G4endl;
141 G4int midBin = (ke1 + ke2)/2;
142 dsig = sig[je2][midBin] - sig[je1][midBin];
144 b = sig[je1][midBin] - rc*elab[je1];
146 if (sample < sigint) {
155 G4cout << ke1 <<
" " << ke2 <<
" "
156 << sigint1 <<
" " << sigint2 <<
G4endl;
158 }
while (ke2 - ke1 > 1);
160 dsig = sigint2 - sigint1;
162 b = ke1 - rc*sigint1;
167 G4cout <<
" energy bin " << je1 <<
" energy=" << elab[je1] <<
G4endl;
180 G4double pseudoMass = std::sqrt(totalEnergy*totalEnergy - P*P);
187 G4double p = std::sqrt(px*px + py*py + pz*pz);
192 G4cout <<
" particle 1 momentum in CM " << px/
GeV <<
" " << py/
GeV <<
" "
198 G4double pxnew = p*std::sin(theta)*std::cos(phi);
199 G4double pynew = p*std::sin(theta)*std::sin(phi);
203 if (px*px + py*py > 0) {
204 G4double cost, sint, ph, cosp, sinp;
206 sint = (std::sqrt(std::fabs((1-cost)*(1+cost))) + std::sqrt(px*px+py*py)/
p)/2;
208 if (std::abs(px) > 0.000001*
GeV) ph = std::atan2(py,px);
211 px = (cost*cosp*pxnew - sinp*pynew + sint*cosp*pznew);
212 py = (cost*sinp*pxnew + cosp*pynew + sint*sinp*pznew);
213 pz = (-sint*pxnew + cost*pznew);
223 G4cout <<
" particle 1 momentum in CM " << px/
GeV <<
" " << py/
GeV <<
" "
234 G4double gammaCM = E1pM2/std::sqrt(E1pM2*E1pM2 - P*P);
237 G4cout <<
" betaCM " << betaCMx <<
" " << betaCMy <<
" "
238 << betaCMz <<
" " << betaCM <<
G4endl;
255 PA[4] = std::sqrt(M1*M1 + p*p);
257 G4double BETPA = BETA[1]*PA[1] + BETA[2]*PA[2] + BETA[3]*PA[3];
258 G4double BPGAM = (BETPA * BETA[4]/(BETA[4] + 1.) - PA[4]) * BETA[4];
260 PB[1] = PA[1] + BPGAM * BETA[1];
261 PB[2] = PA[2] + BPGAM * BETA[2];
262 PB[3] = PA[3] + BPGAM * BETA[3];
263 PB[4] = (PA[4] - BETPA) * BETA[4];
274 PA[4] = std::sqrt(M2*M2 + p*p);
276 BETPA = BETA[1]*PA[1] + BETA[2]*PA[2] + BETA[3]*PA[3];
277 BPGAM = (BETPA * BETA[4]/(BETA[4] + 1.) - PA[4]) * BETA[4];
279 PB[1] = PA[1] + BPGAM * BETA[1];
280 PB[2] = PA[2] + BPGAM * BETA[2];
281 PB[3] = PA[3] + BPGAM * BETA[3];
282 PB[4] = (PA[4] - BETPA) * BETA[4];
287 G4cout <<
" particle 1 momentum in LAB "
290 G4cout <<
" particle 2 momentum in LAB "
293 G4cout <<
" TOTAL momentum in LAB "
316 G4double ek = std::sqrt(plab*plab+nMass*nMass) - nMass;
326 G4int midBin = (je1 + je2)/2;
327 if (ek <
elab[midBin])
331 }
while (je2 - je1 > 1);
342 G4double b = sig[je1][0] - rc*elab[je1];
348 G4int midBin = (ke1 + ke2)/2;
349 dsig = sig[je2][midBin] - sig[je1][midBin];
351 b = sig[je1][midBin] - rc*elab[je1];
364 }
while (ke2 - ke1 > 1);
366 dsig = sigint2 - sigint1;
368 b = ke1 - rc*sigint1;
372 G4double t = 0.5*plab*plab*(1-std::cos(theta));
static const G4float elab[NENERGY]
CLHEP::Hep3Vector G4ThreeVector
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
void SetMomentumChange(const G4ThreeVector &aV)
const G4ThreeVector & GetMomentumDirection() const
const G4String & GetParticleName() const
G4double GetPDGCharge() const
G4double GetTotalMomentum() const
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
static G4Proton * Proton()
G4double GetTotalEnergy() const
G4double GetPDGMass() const
G4DynamicParticle * ReturnTargetParticle() const
static constexpr double proton_mass_c2
void SetEnergyChange(G4double anEnergy)
G4ParticleDefinition * GetDefinition() const
static const G4float sig[NENERGY][NANGLE]
static constexpr double twopi
double A(double temperature)
static constexpr double halfpi
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
const G4LorentzVector & Get4Momentum() const
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
const G4ParticleDefinition * GetDefinition() const
G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
void SetMomentum(const G4ThreeVector &momentum)
static constexpr double degree
G4double GetKineticEnergy() const
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4GLOB_DLL std::ostream G4cout
G4double GetTotalMomentum() const
static constexpr double pi
G4double GetTotalEnergy() const
static constexpr double GeV
G4HadFinalState theParticleChange