85 G4bool IsScatProjToProjCase,
106 adjointPrimKinEnergy,
108 IsScatProjToProjCase);
117 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
118 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
125 if (IsScatProjToProjCase) {
128 G4double companionTotalEnergy =companionM0+ projectileKinEnergy-adjointPrimKinEnergy;
129 G4double companionP2 = companionTotalEnergy*companionTotalEnergy - companionM0*companionM0;
134 G4double P_parallel = (adjointPrimP*adjointPrimP + projectileP2 - companionP2)/(2.*adjointPrimP);
135 G4double P_perp = std::sqrt( projectileP2 - P_parallel*P_parallel);
139 projectileMomentum.rotateUz(dir_parallel);
143 if (!IsScatProjToProjCase ){
161 G4bool IsScatProjToProjCase,
179 if (!IsScatProjToProjCase){
181 eEnergy=adjointPrimKinEnergy;
184 if (Emin>=Emax)
return;
187 newCS=newCS*(b-
a)/eEnergy;
195 if (Emin>=Emax)
return;
196 G4double diff1=Emin-adjointPrimKinEnergy;
197 G4double diff2=Emax-adjointPrimKinEnergy;
199 G4double t1=adjointPrimKinEnergy*(1./diff1-1./diff2);
205 newCS=newCS*sum_t/adjointPrimKinEnergy/adjointPrimKinEnergy;
209 projectileKinEnergy =adjointPrimKinEnergy +1./(1./diff1-q);
214 projectileKinEnergy =1./(1./Emin-q);
217 projectileKinEnergy=Emin*std::pow(Emax/Emin,
G4UniformRand());
220 eEnergy=projectileKinEnergy-adjointPrimKinEnergy;
227 G4double diffCS_perAtom_Used=twopi_mc2_rcl2*mass*adjointPrimKinEnergy/projectileKinEnergy/projectileKinEnergy/eEnergy/eEnergy;
243 w_corr*=diffCS/diffCS_perAtom_Used;
260 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
261 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
268 if (IsScatProjToProjCase) {
271 G4double companionTotalEnergy =companionM0+ projectileKinEnergy-adjointPrimKinEnergy;
272 G4double companionP2 = companionTotalEnergy*companionTotalEnergy - companionM0*companionM0;
277 G4double P_parallel = (adjointPrimP*adjointPrimP + projectileP2 - companionP2)/(2.*adjointPrimP);
278 G4double P_perp = std::sqrt( projectileP2 - P_parallel*P_parallel);
282 projectileMomentum.rotateUz(dir_parallel);
286 if (!IsScatProjToProjCase ){
320 if (kinEnergyProj>Emin_proj && kinEnergyProj<=Emax_proj){
327 if (kinEnergyProj >2.*
MeV){
337 dSigmadEprod=(sigma1-sigma2)/dE;
338 if (dSigmadEprod>1.) {
341 G4cout<<
"dsigma "<<kinEnergyProj/
MeV<<
'\t'<<kinEnergyProd/
MeV<<
'\t'<<dSigmadEprod<<
G4endl;
354 G4double deltaKinEnergy = kinEnergyProd;
368 G4double etot2 = totEnergy*totEnergy;
369 G4double beta2 = kinEnergyProj*(kinEnergyProj + 2.0*
mass)/etot2;
372 f = 1.0 - beta2*deltaKinEnergy/Tmax;
374 f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
380 G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*
mass);
381 gg *= (1.0 +
magMoment2*(x2 - f1/f)/(1.0 + x2));
384 G4cout <<
"### G4BetheBlochModel in Adjoint Sim WARNING: g= " <<
g
407 pname !=
"deuteron" && pname !=
"triton") {
440 G4bool IsScatProjToProjCase)
448 if (!IsScatProjToProjCase ){
452 Cross*=(1./Emin_proj -1./Emax_proj)/primEnergy;
465 G4double diff1=Emin_proj-primEnergy;
466 G4double diff2=Emax_proj-primEnergy;
467 G4double t1=(1./diff1+1./Emin_proj-1./diff2-1./Emax_proj)/primEnergy;
469 G4double t2=2.*std::log(Emax_proj/Emin_proj)/primEnergy/primEnergy;
487 {
return PrimAdjEnergy+Tcut;
G4double one_minus_ratio_2
G4double GetKineticEnergy() const
CLHEP::Hep3Vector G4ThreeVector
G4double GetPostStepWeightCorrection()
const G4DynamicParticle * GetDynamicParticle() const
G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
G4ParticleDefinition * theDirectPrimaryPartDef
void ProposeParentWeight(G4double finalWeight)
G4double GetZ13(G4double Z)
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
void RapidSampleSecondaries(const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
static G4AdjointElectron * AdjointElectron()
G4double CS_biasing_factor
G4Material * currentMaterial
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
static G4NistManager * Instance()
const G4String & GetParticleName() const
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
virtual ~G4AdjointhIonisationModel()
G4double GetTotalMomentum() const
void DefineProjectileProperty()
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
G4VEmModel * theDirectEMModel
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
G4double GetElectronDensity() const
const G4ThreeVector & GetMomentumDirection() const
static G4Proton * Proton()
const G4String & GetParticleType() const
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
void DefineCurrentMaterial(const G4MaterialCutsCouple *couple)
G4bool UseOnlyOneMatrixForAllElements
G4VEmModel * theBraggDirectEMModel
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
G4double GetPDGMass() const
G4double mass_ratio_projectile
static G4AdjointProton * AdjointProton()
virtual void SampleSecondaries(const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
const G4double x[NPOINTSGL]
void ProposeEnergy(G4double finalEnergy)
G4double one_plus_ratio_2
void AddSecondary(G4Track *aSecondary)
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
G4double GetWeight() const
static const G4double Emin
G4double GetPDGSpin() const
static const G4double Emax
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
G4double currentTcutForDirectSecond
G4bool UseMatrixPerElement
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4double GetPDGMagneticMoment() const
void ProposeTrackStatus(G4TrackStatus status)
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
static const double eplus
G4double GetPDGCharge() const
static G4AdjointCSManager * GetAdjointCSManager()
G4int GetLeptonNumber() const
virtual void CorrectPostStepWeight(G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
virtual G4double DiffCrossSectionPerAtomPrimToSecond(G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
G4AdjointhIonisationModel(G4ParticleDefinition *projectileDefinition)
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4bool second_part_of_same_type