Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4AdjointBremsstrahlungModel Class Reference

#include <G4AdjointBremsstrahlungModel.hh>

+ Inheritance diagram for G4AdjointBremsstrahlungModel:

Public Member Functions

 G4AdjointBremsstrahlungModel (G4VEmModel *aModel)
 
 G4AdjointBremsstrahlungModel ()
 
 ~G4AdjointBremsstrahlungModel ()
 
virtual void SampleSecondaries (const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
 
void RapidSampleSecondaries (const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
 
virtual G4double DiffCrossSectionPerVolumePrimToSecond (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
 
G4double DiffCrossSectionPerVolumePrimToSecondApproximated1 (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
 
G4double DiffCrossSectionPerVolumePrimToSecondApproximated2 (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
 
virtual G4double AdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double GetAdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
- Public Member Functions inherited from G4VEmAdjointModel
 G4VEmAdjointModel (const G4String &nam)
 
virtual ~G4VEmAdjointModel ()
 
virtual void SampleSecondaries (const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)=0
 
virtual G4double AdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double GetAdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double DiffCrossSectionPerAtomPrimToSecond (G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
 
virtual G4double DiffCrossSectionPerAtomPrimToScatPrim (G4double kinEnergyProj, G4double kinEnergyScatProj, G4double Z, G4double A=0.)
 
virtual G4double DiffCrossSectionPerVolumePrimToSecond (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
 
virtual G4double DiffCrossSectionPerVolumePrimToScatPrim (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyScatProj)
 
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase (G4double PrimAdjEnergy)
 
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase (G4double PrimAdjEnergy, G4double Tcut=0)
 
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase (G4double PrimAdjEnergy)
 
virtual G4double GetSecondAdjEnergyMinForProdToProjCase (G4double PrimAdjEnergy)
 
void DefineCurrentMaterial (const G4MaterialCutsCouple *couple)
 
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerAtomForSecond (G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
 
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerAtomForScatProj (G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
 
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerVolumeForSecond (G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
 
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerVolumeForScatProj (G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
 
void SetCSMatrices (std::vector< G4AdjointCSMatrix * > *Vec1CSMatrix, std::vector< G4AdjointCSMatrix * > *Vec2CSMatrix)
 
G4ParticleDefinitionGetAdjointEquivalentOfDirectPrimaryParticleDefinition ()
 
G4ParticleDefinitionGetAdjointEquivalentOfDirectSecondaryParticleDefinition ()
 
G4double GetHighEnergyLimit ()
 
G4double GetLowEnergyLimit ()
 
void SetHighEnergyLimit (G4double aVal)
 
void SetLowEnergyLimit (G4double aVal)
 
void DefineDirectEMModel (G4VEmModel *aModel)
 
void SetAdjointEquivalentOfDirectPrimaryParticleDefinition (G4ParticleDefinition *aPart)
 
void SetAdjointEquivalentOfDirectSecondaryParticleDefinition (G4ParticleDefinition *aPart)
 
void SetSecondPartOfSameType (G4bool aBool)
 
G4bool GetSecondPartOfSameType ()
 
void SetUseMatrix (G4bool aBool)
 
void SetUseMatrixPerElement (G4bool aBool)
 
void SetUseOnlyOneMatrixForAllElements (G4bool aBool)
 
void SetApplyCutInRange (G4bool aBool)
 
G4bool GetUseMatrix ()
 
G4bool GetUseMatrixPerElement ()
 
G4bool GetUseOnlyOneMatrixForAllElements ()
 
G4bool GetApplyCutInRange ()
 
G4String GetName ()
 
virtual void SetCSBiasingFactor (G4double aVal)
 
void SetCorrectWeightForPostStepInModel (G4bool aBool)
 
void SetAdditionalWeightCorrectionFactorForPostStepOutsideModel (G4double factor)
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmAdjointModel
G4double DiffCrossSectionFunction1 (G4double kinEnergyProj)
 
G4double DiffCrossSectionFunction2 (G4double kinEnergyProj)
 
G4double DiffCrossSectionPerVolumeFunctionForIntegrationOverEkinProj (G4double EkinProd)
 
G4double SampleAdjSecEnergyFromCSMatrix (size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
 
G4double SampleAdjSecEnergyFromCSMatrix (G4double prim_energy, G4bool IsScatProjToProjCase)
 
void SelectCSMatrix (G4bool IsScatProjToProjCase)
 
virtual G4double SampleAdjSecEnergyFromDiffCrossSectionPerAtom (G4double prim_energy, G4bool IsScatProjToProjCase)
 
virtual void CorrectPostStepWeight (G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
 
- Protected Attributes inherited from G4VEmAdjointModel
G4VEmModeltheDirectEMModel
 
G4VParticleChangepParticleChange
 
const G4String name
 
G4int ASelectedNucleus
 
G4int ZSelectedNucleus
 
G4MaterialSelectedMaterial
 
G4double kinEnergyProdForIntegration
 
G4double kinEnergyScatProjForIntegration
 
G4double kinEnergyProjForIntegration
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForProdToProjBackwardScattering
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForScatProjToProjBackwardScattering
 
std::vector< G4doubleCS_Vs_ElementForScatProjToProjCase
 
std::vector< G4doubleCS_Vs_ElementForProdToProjCase
 
G4double lastCS
 
G4double lastAdjointCSForScatProjToProjCase
 
G4double lastAdjointCSForProdToProjCase
 
G4ParticleDefinitiontheAdjEquivOfDirectPrimPartDef
 
G4ParticleDefinitiontheAdjEquivOfDirectSecondPartDef
 
G4ParticleDefinitiontheDirectPrimaryPartDef
 
G4bool second_part_of_same_type
 
G4double preStepEnergy
 
G4MaterialcurrentMaterial
 
G4MaterialCutsCouplecurrentCouple
 
size_t currentMaterialIndex
 
size_t currentCoupleIndex
 
G4double currentTcutForDirectPrim
 
G4double currentTcutForDirectSecond
 
G4bool ApplyCutInRange
 
G4double mass_ratio_product
 
G4double mass_ratio_projectile
 
G4double HighEnergyLimit
 
G4double LowEnergyLimit
 
G4double CS_biasing_factor
 
G4bool UseMatrix
 
G4bool UseMatrixPerElement
 
G4bool UseOnlyOneMatrixForAllElements
 
size_t indexOfUsedCrossSectionMatrix
 
size_t model_index
 
G4bool correct_weight_for_post_step_in_model
 
G4double additional_weight_correction_factor_for_post_step_outside_model
 

Detailed Description

Definition at line 61 of file G4AdjointBremsstrahlungModel.hh.

Constructor & Destructor Documentation

◆ G4AdjointBremsstrahlungModel() [1/2]

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel ( G4VEmModel aModel)

Definition at line 45 of file G4AdjointBremsstrahlungModel.cc.

45 :
46 G4VEmAdjointModel("AdjointeBremModel")
47{
48 SetUseMatrix(false);
50
51 theDirectStdBremModel = aModel;
52 theDirectEMModel=theDirectStdBremModel;
53 theEmModelManagerForFwdModels = new G4EmModelManager();
54 isDirectModelInitialised = false;
56 G4Region* r=0;
57 theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);
58
60 highKinEnergy= 1.*GeV;
61 lowKinEnergy = 1.0*keV;
62
63 lastCZ =0.;
64
65
70
71
73
74
75}
static G4AdjointElectron * AdjointElectron()
static G4AdjointGamma * AdjointGamma()
static G4Electron * Electron()
Definition: G4Electron.cc:93
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
void SetUseMatrixPerElement(G4bool aBool)
G4VEmModel * theDirectEMModel
void SetUseMatrix(G4bool aBool)
G4ParticleDefinition * theDirectPrimaryPartDef
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
void SetApplyCutInRange(G4bool aBool)
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef

◆ G4AdjointBremsstrahlungModel() [2/2]

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel ( )

Definition at line 78 of file G4AdjointBremsstrahlungModel.cc.

78 :
79 G4VEmAdjointModel("AdjointeBremModel")
80{
81 SetUseMatrix(false);
83
84 theDirectStdBremModel = new G4SeltzerBergerModel();
85 theDirectEMModel=theDirectStdBremModel;
86 theEmModelManagerForFwdModels = new G4EmModelManager();
87 isDirectModelInitialised = false;
89 G4Region* r=0;
90 theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);
91 // theDirectPenelopeBremModel =0;
93 highKinEnergy= 1.*GeV;
94 lowKinEnergy = 1.0*keV;
95 lastCZ =0.;
100}

◆ ~G4AdjointBremsstrahlungModel()

G4AdjointBremsstrahlungModel::~G4AdjointBremsstrahlungModel ( )

Definition at line 103 of file G4AdjointBremsstrahlungModel.cc.

104{if (theDirectStdBremModel) delete theDirectStdBremModel;
105 if (theEmModelManagerForFwdModels) delete theEmModelManagerForFwdModels;
106}

Member Function Documentation

◆ AdjointCrossSection()

G4double G4AdjointBremsstrahlungModel::AdjointCrossSection ( const G4MaterialCutsCouple aCouple,
G4double  primEnergy,
G4bool  IsScatProjToProjCase 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 399 of file G4AdjointBremsstrahlungModel.cc.

402{ if (!isDirectModelInitialised) {
403 theEmModelManagerForFwdModels->Initialise(G4Electron::Electron(),G4Gamma::Gamma(),1.,0);
404 isDirectModelInitialised =true;
405 }
406 if (UseMatrix) return G4VEmAdjointModel::AdjointCrossSection(aCouple,primEnergy,IsScatProjToProjCase);
407 DefineCurrentMaterial(aCouple);
408 G4double Cross=0.;
409 lastCZ=theDirectEMModel->CrossSectionPerVolume(aCouple->GetMaterial(),theDirectPrimaryPartDef,100.*MeV,100.*MeV/std::exp(1.));//this give the constant above
410
411 if (!IsScatProjToProjCase ){
412 G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(primEnergy);
413 G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(primEnergy);
414 if (Emax_proj>Emin_proj && primEnergy > currentTcutForDirectSecond) Cross= CS_biasing_factor*lastCZ*std::log(Emax_proj/Emin_proj);
415 }
416 else {
419 if (Emax_proj>Emin_proj) Cross= lastCZ*std::log((Emax_proj-primEnergy)*Emin_proj/Emax_proj/(Emin_proj-primEnergy));
420
421 }
422 return Cross;
423}
double G4double
Definition: G4Types.hh:83
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
const G4Material * GetMaterial() const
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
void DefineCurrentMaterial(const G4MaterialCutsCouple *couple)
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
G4double currentTcutForDirectSecond
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:254

Referenced by GetAdjointCrossSection().

◆ DiffCrossSectionPerVolumePrimToSecond()

G4double G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecond ( const G4Material aMaterial,
G4double  kinEnergyProj,
G4double  kinEnergyProd 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 323 of file G4AdjointBremsstrahlungModel.cc.

327{if (!isDirectModelInitialised) {
328 theEmModelManagerForFwdModels->Initialise(G4Electron::Electron(),G4Gamma::Gamma(),1.,0);
329 isDirectModelInitialised =true;
330 }
331/*
332 return DiffCrossSectionPerVolumePrimToSecondApproximated2(aMaterial,
333 kinEnergyProj,
334 kinEnergyProd);
335*/
337 kinEnergyProj,
338 kinEnergyProd);
339}
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)

Referenced by RapidSampleSecondaries().

◆ DiffCrossSectionPerVolumePrimToSecondApproximated1()

G4double G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecondApproximated1 ( const G4Material aMaterial,
G4double  kinEnergyProj,
G4double  kinEnergyProd 
)

Definition at line 343 of file G4AdjointBremsstrahlungModel.cc.

348{
349 G4double dCrossEprod=0.;
350 G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
351 G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
352
353
354 //In this approximation we consider that the secondary gammas are sampled with 1/Egamma energy distribution
355 //This is what is applied in the discrete standard model before the rejection test that make a correction
356 //The application of the same rejection function is not possible here.
357 //The differentiation of the CS over Ecut does not produce neither a good differential CS. That is due to the
358 // fact that in the discrete model the differential CS and the integrated CS are both fitted but separatly and
359 // therefore do not allow a correct numerical differentiation of the integrated CS to get the differential one.
360 // In the future we plan to use the brem secondary spectra from the G4Penelope implementation
361
362 if (kinEnergyProj>Emin_proj && kinEnergyProj<=Emax_proj){
363 G4double sigma=theDirectEMModel->CrossSectionPerVolume(aMaterial,theDirectPrimaryPartDef,kinEnergyProj,1.*keV);
364 dCrossEprod=sigma/kinEnergyProd/std::log(kinEnergyProj/keV);
365 }
366 return dCrossEprod;
367
368}

◆ DiffCrossSectionPerVolumePrimToSecondApproximated2()

G4double G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecondApproximated2 ( const G4Material aMaterial,
G4double  kinEnergyProj,
G4double  kinEnergyProd 
)

Definition at line 372 of file G4AdjointBremsstrahlungModel.cc.

377{
378 //In this approximation we derive the direct cross section over Tcut=gamma energy, en after apply the Migdla correction factor
379 //used in the direct model
380
381 G4double dCrossEprod=0.;
382
383 const G4ElementVector* theElementVector = material->GetElementVector();
384 const double* theAtomNumDensityVector = material->GetAtomicNumDensityVector();
385 G4double dum=0.;
386 G4double E1=kinEnergyProd,E2=kinEnergyProd*1.001;
387 G4double dE=E2-E1;
388 for (size_t i=0; i<material->GetNumberOfElements(); i++) {
389 G4double C1=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum ,E1);
390 G4double C2=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum,E2);
391 dCrossEprod += theAtomNumDensityVector[i] * (C1-C2)/dE;
392 }
393
394 return dCrossEprod;
395
396}
std::vector< G4Element * > G4ElementVector
const double C2
#define C1
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:359

◆ GetAdjointCrossSection()

G4double G4AdjointBremsstrahlungModel::GetAdjointCrossSection ( const G4MaterialCutsCouple aCouple,
G4double  primEnergy,
G4bool  IsScatProjToProjCase 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 425 of file G4AdjointBremsstrahlungModel.cc.

428{
429 return AdjointCrossSection(aCouple, primEnergy,IsScatProjToProjCase);
430 lastCZ=theDirectEMModel->CrossSectionPerVolume(aCouple->GetMaterial(),theDirectPrimaryPartDef,100.*MeV,100.*MeV/std::exp(1.));//this give the constant above
431 return G4VEmAdjointModel::GetAdjointCrossSection(aCouple, primEnergy,IsScatProjToProjCase);
432
433}
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
virtual G4double GetAdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)

◆ RapidSampleSecondaries()

void G4AdjointBremsstrahlungModel::RapidSampleSecondaries ( const G4Track aTrack,
G4bool  IsScatProjToProjCase,
G4ParticleChange fParticleChange 
)

Definition at line 188 of file G4AdjointBremsstrahlungModel.cc.

191{
192
193 const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
195
196
197 G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
198 G4double adjointPrimTotalEnergy = theAdjointPrimary->GetTotalEnergy();
199
200 if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
201 return;
202 }
203
204 G4double projectileKinEnergy =0.;
205 G4double gammaEnergy=0.;
206 G4double diffCSUsed=0.;
207 if (!IsScatProjToProjCase){
208 gammaEnergy=adjointPrimKinEnergy;
209 G4double Emax = GetSecondAdjEnergyMaxForProdToProjCase(adjointPrimKinEnergy);
210 G4double Emin= GetSecondAdjEnergyMinForProdToProjCase(adjointPrimKinEnergy);;
211 if (Emin>=Emax) return;
212 projectileKinEnergy=Emin*std::pow(Emax/Emin,G4UniformRand());
213 diffCSUsed=CS_biasing_factor*lastCZ/projectileKinEnergy;
214
215 }
216 else { G4double Emax = GetSecondAdjEnergyMaxForScatProjToProjCase(adjointPrimKinEnergy);
218 if (Emin>=Emax) return;
219 G4double f1=(Emin-adjointPrimKinEnergy)/Emin;
220 G4double f2=(Emax-adjointPrimKinEnergy)/Emax/f1;
221 projectileKinEnergy=adjointPrimKinEnergy/(1.-f1*std::pow(f2,G4UniformRand()));
222 gammaEnergy=projectileKinEnergy-adjointPrimKinEnergy;
223 diffCSUsed=lastCZ*adjointPrimKinEnergy/projectileKinEnergy/gammaEnergy;
224
225 }
226
227
228
229
230 //Weight correction
231 //-----------------------
232 //First w_corr is set to the ratio between adjoint total CS and fwd total CS
233 //if this has to be done in the model
234 //For the case of forced interaction this will be done in the PostStepDoIt of the
235 //forced interaction
236 //It is important to set the weight before the vreation of the secondary
237 //
241 }
242 //G4cout<<"Correction factor start in brem model "<<w_corr<<std::endl;
243
244
245 //Then another correction is needed due to the fact that a biaised differential CS has been used rather than the one consistent with the direct model
246 //Here we consider the true diffCS as the one obtained by the numericla differentiation over Tcut of the direct CS, corrected by the Migdal term.
247 //Basically any other differential CS diffCS could be used here (example Penelope).
248
249 G4double diffCS = DiffCrossSectionPerVolumePrimToSecond(currentMaterial, projectileKinEnergy, gammaEnergy);
250 /*G4cout<<"diffCS "<<diffCS <<std::endl;
251 G4cout<<"diffCS_Used "<<diffCSUsed <<std::endl;*/
252 w_corr*=diffCS/diffCSUsed;
253
254
255 G4double new_weight = aTrack.GetWeight()*w_corr;
256 /*G4cout<<"New weight brem "<<new_weight<<std::endl;
257 G4cout<<"Weight correction brem "<<w_corr<<std::endl;*/
258 fParticleChange->SetParentWeightByProcess(false);
259 fParticleChange->SetSecondaryWeightByProcess(false);
260 fParticleChange->ProposeParentWeight(new_weight);
261
262 //Kinematic
263 //---------
265 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
266 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
267 G4double projectileP = std::sqrt(projectileP2);
268
269
270 //Use the angular model of the forward model to generate the gamma direction
271 //---------------------------------------------------------------------------
272//Dum dynamic particle to use the model
273 G4DynamicParticle * aDynPart = new G4DynamicParticle(G4Electron::Electron(),G4ThreeVector(0.,0.,1.)*projectileP);
274
275 //Get the element from the direct model
277 projectileKinEnergy,currentTcutForDirectSecond);
278 G4int Z=elm->GetZasInt();
279 G4double energy = aDynPart->GetTotalEnergy()-gammaEnergy;
280 G4ThreeVector projectileMomentum =
282 G4double phi = projectileMomentum.getPhi();
283
284/*
285 //Angle of the gamma direction with the projectile taken from G4eBremsstrahlungModel
286 //------------------------------------------------
287 G4double u;
288 const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ;
289
290 if (9./(9.+d) > G4UniformRand()) u = - std::log(G4UniformRand()*G4UniformRand())/a1;
291 else u = - std::log(G4UniformRand()*G4UniformRand())/a2;
292
293 G4double theta = u*electron_mass_c2/projectileTotalEnergy;
294
295 G4double sint = std::sin(theta);
296 G4double cost = std::cos(theta);
297
298 G4double phi = twopi * G4UniformRand() ;
299 G4ThreeVector projectileMomentum;
300 projectileMomentum=G4ThreeVector(std::cos(phi)*sint,std::sin(phi)*sint,cost)*projectileP; //gamma frame
301*/
302 if (IsScatProjToProjCase) {//the adjoint primary is the scattered e-
303 G4ThreeVector gammaMomentum = (projectileTotalEnergy-adjointPrimTotalEnergy)*G4ThreeVector(0.,0.,1.);
304 G4ThreeVector dirProd=projectileMomentum-gammaMomentum;
305 G4double cost1 = std::cos(dirProd.angle(projectileMomentum));
306 G4double sint1 = std::sqrt(1.-cost1*cost1);
307 projectileMomentum=G4ThreeVector(std::cos(phi)*sint1,std::sin(phi)*sint1,cost1)*projectileP;
308 }
309
310 projectileMomentum.rotateUz(theAdjointPrimary->GetMomentumDirection());
311
312 if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
313 fParticleChange->ProposeTrackStatus(fStopAndKill);
314 fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
315 }
316 else {
317 fParticleChange->ProposeEnergy(projectileKinEnergy);
318 fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
319 }
320}
CLHEP::Hep3Vector G4ThreeVector
@ fStopAndKill
int G4int
Definition: G4Types.hh:85
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector unit() const
double angle(const Hep3Vector &) const
double getPhi() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
G4double GetPostStepWeightCorrection()
static G4AdjointCSManager * GetAdjointCSManager()
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4double GetTotalEnergy() const
G4int GetZasInt() const
Definition: G4Element.hh:131
void AddSecondary(G4Track *aSecondary)
void ProposeEnergy(G4double finalEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4double GetWeight() const
const G4DynamicParticle * GetDynamicParticle() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4double additional_weight_correction_factor_for_post_step_outside_model
G4Material * currentMaterial
G4bool correct_weight_for_post_step_in_model
G4MaterialCutsCouple * currentCouple
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:611
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:570
void ProposeTrackStatus(G4TrackStatus status)
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
void ProposeParentWeight(G4double finalWeight)
G4double energy(const ThreeVector &p, const G4double m)

Referenced by SampleSecondaries().

◆ SampleSecondaries()

void G4AdjointBremsstrahlungModel::SampleSecondaries ( const G4Track aTrack,
G4bool  IsScatProjToProjCase,
G4ParticleChange fParticleChange 
)
virtual

Implements G4VEmAdjointModel.

Definition at line 110 of file G4AdjointBremsstrahlungModel.cc.

113{
114 if (!UseMatrix) return RapidSampleSecondaries(aTrack,IsScatProjToProjCase,fParticleChange);
115
116 const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
118
119
120 G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
121 G4double adjointPrimTotalEnergy = theAdjointPrimary->GetTotalEnergy();
122
123 if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
124 return;
125 }
126
127 G4double projectileKinEnergy = SampleAdjSecEnergyFromCSMatrix(adjointPrimKinEnergy,
128 IsScatProjToProjCase);
129 //Weight correction
130 //-----------------------
131 CorrectPostStepWeight(fParticleChange,
132 aTrack.GetWeight(),
133 adjointPrimKinEnergy,
134 projectileKinEnergy,
135 IsScatProjToProjCase);
136
137
138 //Kinematic
139 //---------
141 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
142 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
143 G4double projectileP = std::sqrt(projectileP2);
144
145
146 //Angle of the gamma direction with the projectile taken from G4eBremsstrahlungModel
147 //------------------------------------------------
148 G4double u;
149 const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ;
150
151 if (9./(9.+d) > G4UniformRand()) u = - std::log(G4UniformRand()*G4UniformRand())/a1;
152 else u = - std::log(G4UniformRand()*G4UniformRand())/a2;
153
154 G4double theta = u*electron_mass_c2/projectileTotalEnergy;
155
156 G4double sint = std::sin(theta);
157 G4double cost = std::cos(theta);
158
159 G4double phi = twopi * G4UniformRand() ;
160
161 G4ThreeVector projectileMomentum;
162 projectileMomentum=G4ThreeVector(std::cos(phi)*sint,std::sin(phi)*sint,cost)*projectileP; //gamma frame
163 if (IsScatProjToProjCase) {//the adjoint primary is the scattered e-
164 G4ThreeVector gammaMomentum = (projectileTotalEnergy-adjointPrimTotalEnergy)*G4ThreeVector(0.,0.,1.);
165 G4ThreeVector dirProd=projectileMomentum-gammaMomentum;
166 G4double cost1 = std::cos(dirProd.angle(projectileMomentum));
167 G4double sint1 = std::sqrt(1.-cost1*cost1);
168 projectileMomentum=G4ThreeVector(std::cos(phi)*sint1,std::sin(phi)*sint1,cost1)*projectileP;
169
170 }
171
172 projectileMomentum.rotateUz(theAdjointPrimary->GetMomentumDirection());
173
174
175
176 if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
177 fParticleChange->ProposeTrackStatus(fStopAndKill);
178 fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
179 }
180 else {
181 fParticleChange->ProposeEnergy(projectileKinEnergy);
182 fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
183
184 }
185}
void RapidSampleSecondaries(const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
virtual void CorrectPostStepWeight(G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)

The documentation for this class was generated from the following files: