292{
295
296 G4double maxKinEnergy = std::min(maxEnergy,tmax);
297 if(minKinEnergy >= maxKinEnergy) { return; }
298
299
300
301
302 G4double totEnergy = kineticEnergy + mass;
303 G4double etot2 = totEnergy*totEnergy;
304 G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/etot2;
305
309 if( 0.0 < spin ) { fmax += 0.5*maxKinEnergy*maxKinEnergy/etot2; }
310
313
314
315 do {
317 deltaKinEnergy = minKinEnergy*maxKinEnergy
318 /(minKinEnergy*(1.0 - rndm[0]) + maxKinEnergy*rndm[0]);
319
320 f = 1.0 - beta2*deltaKinEnergy/tmax;
321 if( 0.0 < spin ) {
322 f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
323 f += f1;
324 }
325
326
327 } while( fmax*rndm[1] > f);
328
329
330
331
332 G4double x = formfact*deltaKinEnergy*(deltaKinEnergy + 2*electron_mass_c2);
333 if(x > 1.e-6) {
334
337 if( 0.0 < spin ) {
338 G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
339 grej *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
340 }
341 if(grej > 1.1) {
342 G4cout <<
"### G4AtimaEnergyLossModel WARNING: grej= " << grej
344 << " Ekin(MeV)= " << kineticEnergy
345 << " delEkin(MeV)= " << deltaKinEnergy
347 }
348 if(rndmEngineMod->
flat() > grej) {
return; }
349 }
350
352
354
357
358 deltaDirection =
360
361 } else {
362
364 std::sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
365 G4double cost = deltaKinEnergy * (totEnergy + electron_mass_c2) /
367 if(cost > 1.0) { cost = 1.0; }
368 G4double sint = std::sqrt((1.0 - cost)*(1.0 + cost));
369
371
372 deltaDirection.
set(sint*cos(phi),sint*sin(phi), cost) ;
374 }
375
376
377
378
379
380
381
382
383
384
385
386
388
389 vdp->push_back(delta);
390
391
392 kineticEnergy -= deltaKinEnergy;
394 finalP = finalP.
unit();
395
398}
G4GLOB_DLL std::ostream G4cout
void set(double x, double y, double z)
Hep3Vector & rotateUz(const Hep3Vector &)
virtual void flatArray(const int size, double *vect)=0
const G4ThreeVector & GetMomentumDirection() const
G4ThreeVector GetMomentum() const
G4double GetTotalMomentum() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void SetProposedMomentumDirection(const G4ThreeVector &dir)
const G4String & GetParticleName() const
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4int SelectRandomAtomNumber(const G4Material *) const