325{
327
329 const G4double minKinEnergy = std::min(cut, tmax);
330 const G4double maxKinEnergy = std::min(maxEnergy, tmax);
331 if(minKinEnergy >= maxKinEnergy) { return; }
332
333
334
335
336 G4double totEnergy = kineticEnergy + mass;
337 G4double etot2 = totEnergy*totEnergy;
338 G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/etot2;
339
343 if( 0.0 < spin ) { fmax += 0.5*maxKinEnergy*maxKinEnergy/etot2; }
344
345 CLHEP::HepRandomEngine* rndmEngineMod = G4Random::getTheEngine();
347
348
349 do {
351 deltaKinEnergy = minKinEnergy*maxKinEnergy
352 /(minKinEnergy*(1.0 - rndm[0]) + maxKinEnergy*rndm[0]);
353
354 f = 1.0 - beta2*deltaKinEnergy/tmax;
355 if( 0.0 < spin ) {
356 f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
357 f += f1;
358 }
359
360
361 } while( fmax*rndm[1] > f);
362
363
364
365
366 G4double x = formfact*deltaKinEnergy;
367 if(x > 1.e-6) {
368
371 if( 0.0 < spin ) {
372 G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
373 grej *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
374 }
375 if(grej > 1.1) {
376 G4cout <<
"### G4LindhardSorensenIonModel WARNING: grej= " << grej
378 << " Ekin(MeV)= " << kineticEnergy
379 << " delEkin(MeV)= " << deltaKinEnergy
381 }
382 if(rndmEngineMod->
flat() > grej) {
return; }
383 }
384
386
388
391
392 deltaDirection =
394
395 } else {
396
398 std::sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
399 G4double cost = deltaKinEnergy * (totEnergy + electron_mass_c2) /
401 cost = std::min(cost, 1.0);
402 G4double sint = std::sqrt((1.0 - cost)*(1.0 + cost));
403
405
406 deltaDirection.
set(sint*std::cos(phi),sint*std::sin(phi), cost) ;
408 }
409
410
411
412
413
414
415
416
417
418
419
420
421 auto delta = new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
422
423 vdp->push_back(delta);
424
425
426 kineticEnergy -= deltaKinEnergy;
428 finalP = finalP.
unit();
429
430 fParticleChange->SetProposedKineticEnergy(kineticEnergy);
431 fParticleChange->SetProposedMomentumDirection(finalP);
432}
CLHEP::Hep3Vector G4ThreeVector
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
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
G4double GetTotalMomentum() const
const G4String & GetParticleName() const
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4int SelectRandomAtomNumber(const G4Material *) const