Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
Public Member Functions | List of all members
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)
 
G4ParticleDefinition * GetAdjointEquivalentOfDirectPrimaryParticleDefinition ()
 
G4ParticleDefinition * GetAdjointEquivalentOfDirectSecondaryParticleDefinition ()
 
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)
 

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
G4VEmModel * theDirectEMModel
 
G4VParticleChange * pParticleChange
 
const G4String name
 
G4int ASelectedNucleus
 
G4int ZSelectedNucleus
 
G4Material * SelectedMaterial
 
G4double kinEnergyProdForIntegration
 
G4double kinEnergyScatProjForIntegration
 
G4double kinEnergyProjForIntegration
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForProdToProjBackwardScattering
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForScatProjToProjBackwardScattering
 
std::vector< G4double > CS_Vs_ElementForScatProjToProjCase
 
std::vector< G4double > CS_Vs_ElementForProdToProjCase
 
G4double lastCS
 
G4double lastAdjointCSForScatProjToProjCase
 
G4double lastAdjointCSForProdToProjCase
 
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
 
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
 
G4ParticleDefinition * theDirectPrimaryPartDef
 
G4bool second_part_of_same_type
 
G4double preStepEnergy
 
G4Material * currentMaterial
 
G4MaterialCutsCouple * currentCouple
 
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
 

Detailed Description

Definition at line 63 of file G4AdjointBremsstrahlungModel.hh.

Constructor & Destructor Documentation

◆ G4AdjointBremsstrahlungModel() [1/2]

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel ( G4VEmModel *  aModel)

Definition at line 46 of file G4AdjointBremsstrahlungModel.cc.

46 :
47 G4VEmAdjointModel("AdjointeBremModel")
48{
49 SetUseMatrix(false);
50 SetUseMatrixPerElement(false);
51
52 theDirectStdBremModel = aModel;
53 theDirectEMModel=theDirectStdBremModel;
54 theEmModelManagerForFwdModels = new G4EmModelManager();
55 isDirectModelInitialised = false;
56 G4VEmFluctuationModel* f=0;
57 G4Region* r=0;
58 theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);
59
60 SetApplyCutInRange(true);
61 highKinEnergy= 100.*TeV;
62 lowKinEnergy = 1.0*keV;
63
64 lastCZ =0.;
65
66
67 theAdjEquivOfDirectPrimPartDef =G4AdjointElectron::AdjointElectron();
68 theAdjEquivOfDirectSecondPartDef=G4AdjointGamma::AdjointGamma();
69 theDirectPrimaryPartDef=G4Electron::Electron();
70 second_part_of_same_type=false;
71
72 /*UsePenelopeModel=false;
73 if (UsePenelopeModel) {
74 G4PenelopeBremsstrahlungModel* thePenelopeModel = new G4PenelopeBremsstrahlungModel(G4Electron::Electron(),"PenelopeBrem");
75 theEmModelManagerForFwdModels = new G4EmModelManager();
76 isPenelopeModelInitialised = false;
77 G4VEmFluctuationModel* f=0;
78 G4Region* r=0;
79 theDirectEMModel=thePenelopeModel;
80 theEmModelManagerForFwdModels->AddEmModel(1, thePenelopeModel, f, r);
81 }
82 */
83
84
85
86}
G4AdjointElectron::AdjointElectron
static G4AdjointElectron * AdjointElectron()
Definition: G4AdjointElectron.cc:100
G4AdjointGamma::AdjointGamma
static G4AdjointGamma * AdjointGamma()
Definition: G4AdjointGamma.cc:84
G4Electron::Electron
static G4Electron * Electron()
Definition: G4Electron.cc:94
G4EmModelManager::AddEmModel
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
Definition: G4EmModelManager.cc:171
G4VEmAdjointModel::SetUseMatrixPerElement
void SetUseMatrixPerElement(G4bool aBool)
Definition: G4VEmAdjointModel.hh:207
G4VEmAdjointModel::theDirectEMModel
G4VEmModel * theDirectEMModel
Definition: G4VEmAdjointModel.hh:255
G4VEmAdjointModel::SetUseMatrix
void SetUseMatrix(G4bool aBool)
Definition: G4VEmAdjointModel.hh:205
G4VEmAdjointModel::theDirectPrimaryPartDef
G4ParticleDefinition * theDirectPrimaryPartDef
Definition: G4VEmAdjointModel.hh:298
G4VEmAdjointModel::theAdjEquivOfDirectSecondPartDef
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
Definition: G4VEmAdjointModel.hh:297
G4VEmAdjointModel::SetApplyCutInRange
void SetApplyCutInRange(G4bool aBool)
Definition: G4VEmAdjointModel.hh:210
G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
Definition: G4VEmAdjointModel.hh:296

◆ G4AdjointBremsstrahlungModel() [2/2]

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel ( )

Definition at line 89 of file G4AdjointBremsstrahlungModel.cc.

89 :
90 G4VEmAdjointModel("AdjointeBremModel")
91{
92 SetUseMatrix(false);
93 SetUseMatrixPerElement(false);
94
95 theDirectStdBremModel = new G4SeltzerBergerModel();
96 theDirectEMModel=theDirectStdBremModel;
97 theEmModelManagerForFwdModels = new G4EmModelManager();
98 isDirectModelInitialised = false;
99 G4VEmFluctuationModel* f=0;
100 G4Region* r=0;
101 theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);
102 // theDirectPenelopeBremModel =0;
103 SetApplyCutInRange(true);
104 highKinEnergy= 1.*GeV;
105 lowKinEnergy = 1.0*keV;
106 lastCZ =0.;
107 theAdjEquivOfDirectPrimPartDef =G4AdjointElectron::AdjointElectron();
108 theAdjEquivOfDirectSecondPartDef=G4AdjointGamma::AdjointGamma();
109 theDirectPrimaryPartDef=G4Electron::Electron();
110 second_part_of_same_type=false;
111
112}

◆ ~G4AdjointBremsstrahlungModel()

G4AdjointBremsstrahlungModel::~G4AdjointBremsstrahlungModel ( )

Definition at line 115 of file G4AdjointBremsstrahlungModel.cc.

116{if (theDirectStdBremModel) delete theDirectStdBremModel;
117 if (theEmModelManagerForFwdModels) delete theEmModelManagerForFwdModels;
118}

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= lastCZ*std::log(Emax_proj/Emin_proj);
415 }
416 else {
417 G4double Emax_proj = GetSecondAdjEnergyMaxForScatProjToProjCase(primEnergy);
418 G4double Emin_proj = GetSecondAdjEnergyMinForScatProjToProjCase(primEnergy,currentTcutForDirectSecond);
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}
G4double
double G4double
Definition: G4Types.hh:64
G4EmModelManager::Initialise
const G4DataVector * Initialise(const G4ParticleDefinition *, const G4ParticleDefinition *, G4double, G4int)
Definition: G4EmModelManager.cc:226
G4Gamma::Gamma
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4MaterialCutsCouple::GetMaterial
const G4Material * GetMaterial() const
Definition: G4MaterialCutsCouple.hh:150
G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
Definition: G4VEmAdjointModel.cc:671
G4VEmAdjointModel::DefineCurrentMaterial
void DefineCurrentMaterial(const G4MaterialCutsCouple *couple)
Definition: G4VEmAdjointModel.cc:683
G4VEmAdjointModel::GetSecondAdjEnergyMaxForScatProjToProjCase
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
Definition: G4VEmAdjointModel.cc:657
G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
Definition: G4VEmAdjointModel.cc:676
G4VEmAdjointModel::currentTcutForDirectSecond
G4double currentTcutForDirectSecond
Definition: G4VEmAdjointModel.hh:313
G4VEmAdjointModel::AdjointCrossSection
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
Definition: G4VEmAdjointModel.cc:58
G4VEmAdjointModel::GetSecondAdjEnergyMinForScatProjToProjCase
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
Definition: G4VEmAdjointModel.cc:664
G4VEmModel::CrossSectionPerVolume
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:186

Referenced by GetAdjointCrossSection().

◆ DiffCrossSectionPerVolumePrimToSecond()

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

Reimplemented from G4VEmAdjointModel.

Definition at line 310 of file G4AdjointBremsstrahlungModel.cc.

314{if (!isDirectModelInitialised) {
315 theEmModelManagerForFwdModels->Initialise(G4Electron::Electron(),G4Gamma::Gamma(),1.,0);
316 isDirectModelInitialised =true;
317 }
318
319 return DiffCrossSectionPerVolumePrimToSecondApproximated2(aMaterial,
320 kinEnergyProj,
321 kinEnergyProd);
322 /*return G4VEmAdjointModel::DiffCrossSectionPerVolumePrimToSecond(aMaterial,
323 kinEnergyProj,
324 kinEnergyProd);*/
325}
G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecondApproximated2
G4double DiffCrossSectionPerVolumePrimToSecondApproximated2(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
Definition: G4AdjointBremsstrahlungModel.cc:358

Referenced by RapidSampleSecondaries().

◆ DiffCrossSectionPerVolumePrimToSecondApproximated1()

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

Definition at line 329 of file G4AdjointBremsstrahlungModel.cc.

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

◆ DiffCrossSectionPerVolumePrimToSecondApproximated2()

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

Definition at line 358 of file G4AdjointBremsstrahlungModel.cc.

363{
364 //In this approximation we derive the direct cross section over Tcut=gamma energy, en after apply the Migdla correction factor
365 //used in the direct model
366
367 G4double dCrossEprod=0.;
368
369 const G4ElementVector* theElementVector = material->GetElementVector();
370 const double* theAtomNumDensityVector = material->GetAtomicNumDensityVector();
371 G4double dum=0.;
372 G4double E1=kinEnergyProd,E2=kinEnergyProd*1.001;
373 G4double dE=E2-E1;
374 for (size_t i=0; i<material->GetNumberOfElements(); i++) {
375 G4double C1=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum ,E1);
376 G4double C2=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum,E2);
377 dCrossEprod += theAtomNumDensityVector[i] * (C1-C2)/dE;
378
379 }
380
381 //Now the Migdal correction
382/*
383 G4double totalEnergy = kinEnergyProj+electron_mass_c2 ;
384 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy
385 *(material->GetElectronDensity());
386
387
388 G4double MigdalFactor = 1./(1.+kp2/(kinEnergyProd*kinEnergyProd)); // its seems that the factor used in the CS compuation i the direct
389 //model is different than the one used in the secondary sampling by a
390 //factor (1.+kp2) To be checked!
391
392 dCrossEprod*=MigdalFactor;
393 */
394 return dCrossEprod;
395
396}
G4ElementVector
std::vector< G4Element * > G4ElementVector
Definition: G4ElementVector.hh:44
C1
#define C1
G4VEmModel::ComputeCrossSectionPerAtom
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:240

Referenced by DiffCrossSectionPerVolumePrimToSecond().

◆ 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}
G4AdjointBremsstrahlungModel::AdjointCrossSection
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
Definition: G4AdjointBremsstrahlungModel.cc:399
G4VEmAdjointModel::GetAdjointCrossSection
virtual G4double GetAdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
Definition: G4VEmAdjointModel.cc:83

◆ RapidSampleSecondaries()

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

Definition at line 200 of file G4AdjointBremsstrahlungModel.cc.

203{
204
205 const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
206 DefineCurrentMaterial(aTrack.GetMaterialCutsCouple());
207
208
209 G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
210 G4double adjointPrimTotalEnergy = theAdjointPrimary->GetTotalEnergy();
211
212 if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
213 return;
214 }
215
216 G4double projectileKinEnergy =0.;
217 G4double gammaEnergy=0.;
218 G4double diffCSUsed=0.;
219 if (!IsScatProjToProjCase){
220 gammaEnergy=adjointPrimKinEnergy;
221 G4double Emax = GetSecondAdjEnergyMaxForProdToProjCase(adjointPrimKinEnergy);
222 G4double Emin= GetSecondAdjEnergyMinForProdToProjCase(adjointPrimKinEnergy);;
223 if (Emin>=Emax) return;
224 projectileKinEnergy=Emin*std::pow(Emax/Emin,G4UniformRand());
225 diffCSUsed=lastCZ/projectileKinEnergy;
226
227 }
228 else { G4double Emax = GetSecondAdjEnergyMaxForScatProjToProjCase(adjointPrimKinEnergy);
229 G4double Emin = GetSecondAdjEnergyMinForScatProjToProjCase(adjointPrimKinEnergy,currentTcutForDirectSecond);
230 if (Emin>=Emax) return;
231 G4double f1=(Emin-adjointPrimKinEnergy)/Emin;
232 G4double f2=(Emax-adjointPrimKinEnergy)/Emax/f1;
233 //G4cout<<"f1 and f2 "<<f1<<'\t'<<f2<<G4endl;
234 projectileKinEnergy=adjointPrimKinEnergy/(1.-f1*std::pow(f2,G4UniformRand()));
235 gammaEnergy=projectileKinEnergy-adjointPrimKinEnergy;
236 diffCSUsed=lastCZ*adjointPrimKinEnergy/projectileKinEnergy/gammaEnergy;
237
238 }
239
240
241
242
243 //Weight correction
244 //-----------------------
245 //First w_corr is set to the ratio between adjoint total CS and fwd total CS
246 G4double w_corr=G4AdjointCSManager::GetAdjointCSManager()->GetPostStepWeightCorrection();
247
248 //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
249 //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.
250 //Basically any other differential CS diffCS could be used here (example Penelope).
251
252 G4double diffCS = DiffCrossSectionPerVolumePrimToSecond(currentMaterial, projectileKinEnergy, gammaEnergy);
253 w_corr*=diffCS/diffCSUsed;
254
255 G4double new_weight = aTrack.GetWeight()*w_corr;
256 fParticleChange->SetParentWeightByProcess(false);
257 fParticleChange->SetSecondaryWeightByProcess(false);
258 fParticleChange->ProposeParentWeight(new_weight);
259
260 //Kinematic
261 //---------
262 G4double projectileM0 = theAdjEquivOfDirectPrimPartDef->GetPDGMass();
263 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
264 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
265 G4double projectileP = std::sqrt(projectileP2);
266
267
268 //Angle of the gamma direction with the projectile taken from G4eBremsstrahlungModel
269 //------------------------------------------------
270 G4double u;
271 const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ;
272
273 if (9./(9.+d) > G4UniformRand()) u = - std::log(G4UniformRand()*G4UniformRand())/a1;
274 else u = - std::log(G4UniformRand()*G4UniformRand())/a2;
275
276 G4double theta = u*electron_mass_c2/projectileTotalEnergy;
277
278 G4double sint = std::sin(theta);
279 G4double cost = std::cos(theta);
280
281 G4double phi = twopi * G4UniformRand() ;
282
283 G4ThreeVector projectileMomentum;
284 projectileMomentum=G4ThreeVector(std::cos(phi)*sint,std::sin(phi)*sint,cost)*projectileP; //gamma frame
285 if (IsScatProjToProjCase) {//the adjoint primary is the scattered e-
286 G4ThreeVector gammaMomentum = (projectileTotalEnergy-adjointPrimTotalEnergy)*G4ThreeVector(0.,0.,1.);
287 G4ThreeVector dirProd=projectileMomentum-gammaMomentum;
288 G4double cost1 = std::cos(dirProd.angle(projectileMomentum));
289 G4double sint1 = std::sqrt(1.-cost1*cost1);
290 projectileMomentum=G4ThreeVector(std::cos(phi)*sint1,std::sin(phi)*sint1,cost1)*projectileP;
291
292 }
293
294 projectileMomentum.rotateUz(theAdjointPrimary->GetMomentumDirection());
295
296
297
298 if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
299 fParticleChange->ProposeTrackStatus(fStopAndKill);
300 fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
301 }
302 else {
303 fParticleChange->ProposeEnergy(projectileKinEnergy);
304 fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
305
306 }
307}
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
fStopAndKill
@ fStopAndKill
Definition: G4TrackStatus.hh:56
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:53
CLHEP::Hep3Vector::unit
Hep3Vector unit() const
CLHEP::Hep3Vector::angle
double angle(const Hep3Vector &) const
CLHEP::Hep3Vector::rotateUz
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
G4AdjointBremsstrahlungModel::DiffCrossSectionPerVolumePrimToSecond
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
Definition: G4AdjointBremsstrahlungModel.cc:310
G4AdjointCSManager::GetPostStepWeightCorrection
G4double GetPostStepWeightCorrection()
Definition: G4AdjointCSManager.cc:525
G4AdjointCSManager::GetAdjointCSManager
static G4AdjointCSManager * GetAdjointCSManager()
Definition: G4AdjointCSManager.cc:61
G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4DynamicParticle::GetTotalEnergy
G4double GetTotalEnergy() const
G4ParticleChange::AddSecondary
void AddSecondary(G4Track *aSecondary)
Definition: G4ParticleChange.cc:218
G4ParticleChange::ProposeEnergy
void ProposeEnergy(G4double finalEnergy)
G4ParticleChange::ProposeMomentumDirection
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4Track::GetWeight
G4double GetWeight() const
G4Track::GetDynamicParticle
const G4DynamicParticle * GetDynamicParticle() const
G4Track::GetMaterialCutsCouple
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4VEmAdjointModel::currentMaterial
G4Material * currentMaterial
Definition: G4VEmAdjointModel.hh:308
G4VParticleChange::ProposeTrackStatus
void ProposeTrackStatus(G4TrackStatus status)
G4VParticleChange::SetSecondaryWeightByProcess
void SetSecondaryWeightByProcess(G4bool)
G4VParticleChange::SetParentWeightByProcess
void SetParentWeightByProcess(G4bool)
Definition: G4VParticleChange.cc:516
G4VParticleChange::ProposeParentWeight
void ProposeParentWeight(G4double finalWeight)

Referenced by SampleSecondaries().

◆ SampleSecondaries()

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

Implements G4VEmAdjointModel.

Definition at line 122 of file G4AdjointBremsstrahlungModel.cc.

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