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 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	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 61 of file G4AdjointBremsstrahlungModel.hh.

Constructor & Destructor Documentation

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel

(

G4VEmModel *

aModel

)

Definition at line 46 of file G4AdjointBremsstrahlungModel.cc.

References G4EmModelManager::AddEmModel(), G4AdjointElectron::AdjointElectron(), G4AdjointGamma::AdjointGamma(), G4Electron::Electron(), python.hepunit::GeV, python.hepunit::keV, G4VEmAdjointModel::second_part_of_same_type, G4VEmAdjointModel::SetApplyCutInRange(), G4VEmAdjointModel::SetUseMatrix(), G4VEmAdjointModel::SetUseMatrixPerElement(), G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef, G4VEmAdjointModel::theAdjEquivOfDirectSecondPartDef, G4VEmAdjointModel::theDirectEMModel, and G4VEmAdjointModel::theDirectPrimaryPartDef.

                                                                            :
 G4VEmAdjointModel("AdjointeBremModel")
{ 
  SetUseMatrix(false);
  SetUseMatrixPerElement(false);
  
  theDirectStdBremModel = aModel;
  theDirectEMModel=theDirectStdBremModel;
  theEmModelManagerForFwdModels = new G4EmModelManager();
  isDirectModelInitialised = false;
  G4VEmFluctuationModel* f=0;
  G4Region* r=0;
  theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);

  SetApplyCutInRange(true);
  highKinEnergy= 1.*GeV;
  lowKinEnergy = 1.0*keV;

  lastCZ =0.;

  
  theAdjEquivOfDirectPrimPartDef =G4AdjointElectron::AdjointElectron();
  theAdjEquivOfDirectSecondPartDef=G4AdjointGamma::AdjointGamma();
  theDirectPrimaryPartDef=G4Electron::Electron();
  second_part_of_same_type=false;

  
  /*UsePenelopeModel=false;
  if (UsePenelopeModel) {
    G4PenelopeBremsstrahlungModel* thePenelopeModel = new G4PenelopeBremsstrahlungModel(G4Electron::Electron(),"PenelopeBrem");
    theEmModelManagerForFwdModels = new G4EmModelManager();
    isPenelopeModelInitialised = false;
    G4VEmFluctuationModel* f=0;
    G4Region* r=0;
    theDirectEMModel=thePenelopeModel;
    theEmModelManagerForFwdModels->AddEmModel(1, thePenelopeModel, f, r);
  }
  */    
  

  
}

G4AdjointBremsstrahlungModel::G4AdjointBremsstrahlungModel

(

)

Definition at line 90 of file G4AdjointBremsstrahlungModel.cc.

References G4EmModelManager::AddEmModel(), G4AdjointElectron::AdjointElectron(), G4AdjointGamma::AdjointGamma(), G4Electron::Electron(), python.hepunit::GeV, python.hepunit::keV, G4VEmAdjointModel::second_part_of_same_type, G4VEmAdjointModel::SetApplyCutInRange(), G4VEmAdjointModel::SetUseMatrix(), G4VEmAdjointModel::SetUseMatrixPerElement(), G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef, G4VEmAdjointModel::theAdjEquivOfDirectSecondPartDef, G4VEmAdjointModel::theDirectEMModel, and G4VEmAdjointModel::theDirectPrimaryPartDef.

                                                          :
 G4VEmAdjointModel("AdjointeBremModel")
{
  SetUseMatrix(false);
  SetUseMatrixPerElement(false);

  theDirectStdBremModel = new G4SeltzerBergerModel();
  theDirectEMModel=theDirectStdBremModel;
  theEmModelManagerForFwdModels = new G4EmModelManager();
  isDirectModelInitialised = false;
  G4VEmFluctuationModel* f=0;
  G4Region* r=0;
  theEmModelManagerForFwdModels->AddEmModel(1, theDirectStdBremModel, f, r);
 // theDirectPenelopeBremModel =0;
  SetApplyCutInRange(true);
  highKinEnergy= 1.*GeV;
  lowKinEnergy = 1.0*keV;
  lastCZ =0.;
  theAdjEquivOfDirectPrimPartDef =G4AdjointElectron::AdjointElectron();
  theAdjEquivOfDirectSecondPartDef=G4AdjointGamma::AdjointGamma();
  theDirectPrimaryPartDef=G4Electron::Electron();
  second_part_of_same_type=false;

}

G4AdjointBremsstrahlungModel::~G4AdjointBremsstrahlungModel

(

)

Definition at line 116 of file G4AdjointBremsstrahlungModel.cc.

{if (theDirectStdBremModel) delete theDirectStdBremModel;
 if (theEmModelManagerForFwdModels) delete theEmModelManagerForFwdModels;
}

Member Function Documentation

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

virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 400 of file G4AdjointBremsstrahlungModel.cc.

References G4VEmAdjointModel::AdjointCrossSection(), G4VEmModel::CrossSectionPerVolume(), G4VEmAdjointModel::CS_biasing_factor, G4VEmAdjointModel::currentTcutForDirectSecond, G4VEmAdjointModel::DefineCurrentMaterial(), G4Electron::Electron(), G4Gamma::Gamma(), G4MaterialCutsCouple::GetMaterial(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForScatProjToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForScatProjToProjCase(), G4EmModelManager::Initialise(), python.hepunit::MeV, G4VEmAdjointModel::theDirectEMModel, G4VEmAdjointModel::theDirectPrimaryPartDef, and G4VEmAdjointModel::UseMatrix.

Referenced by GetAdjointCrossSection().

{ if (!isDirectModelInitialised) {
    theEmModelManagerForFwdModels->Initialise(G4Electron::Electron(),G4Gamma::Gamma(),1.,0);
    isDirectModelInitialised =true;
  }
  if (UseMatrix) return G4VEmAdjointModel::AdjointCrossSection(aCouple,primEnergy,IsScatProjToProjCase);
  DefineCurrentMaterial(aCouple);
  G4double Cross=0.;
  lastCZ=theDirectEMModel->CrossSectionPerVolume(aCouple->GetMaterial(),theDirectPrimaryPartDef,100.*MeV,100.*MeV/std::exp(1.));//this give the constant above
  
  if (!IsScatProjToProjCase ){
    G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(primEnergy);
    G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(primEnergy);
    if (Emax_proj>Emin_proj && primEnergy > currentTcutForDirectSecond) Cross= 100.*CS_biasing_factor*lastCZ*std::log(Emax_proj/Emin_proj);
  }
  else {
    G4double Emax_proj = GetSecondAdjEnergyMaxForScatProjToProjCase(primEnergy);
    G4double Emin_proj = GetSecondAdjEnergyMinForScatProjToProjCase(primEnergy,currentTcutForDirectSecond);
    if (Emax_proj>Emin_proj) Cross= lastCZ*std::log((Emax_proj-primEnergy)*Emin_proj/Emax_proj/(Emin_proj-primEnergy));
    
  }
  return Cross; 
}                        

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

virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 311 of file G4AdjointBremsstrahlungModel.cc.

References DiffCrossSectionPerVolumePrimToSecondApproximated2(), G4Electron::Electron(), G4Gamma::Gamma(), and G4EmModelManager::Initialise().

Referenced by RapidSampleSecondaries().

{if (!isDirectModelInitialised) {
    theEmModelManagerForFwdModels->Initialise(G4Electron::Electron(),G4Gamma::Gamma(),1.,0);
    isDirectModelInitialised =true;
 }

 return  DiffCrossSectionPerVolumePrimToSecondApproximated2(aMaterial,
                                                         kinEnergyProj, 
                                                         kinEnergyProd);
 /*return G4VEmAdjointModel::DiffCrossSectionPerVolumePrimToSecond(aMaterial,
                                                         kinEnergyProj, 
                                                         kinEnergyProd);*/                                                               
}                     

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

Definition at line 330 of file G4AdjointBremsstrahlungModel.cc.

References G4VEmModel::CrossSectionPerVolume(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase(), python.hepunit::keV, G4VEmAdjointModel::theDirectEMModel, and G4VEmAdjointModel::theDirectPrimaryPartDef.

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

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

Definition at line 359 of file G4AdjointBremsstrahlungModel.cc.

References C1, G4VEmModel::ComputeCrossSectionPerAtom(), G4Material::GetAtomicNumDensityVector(), G4Material::GetElementVector(), G4Material::GetNumberOfElements(), G4VEmAdjointModel::theDirectEMModel, and G4VEmAdjointModel::theDirectPrimaryPartDef.

Referenced by DiffCrossSectionPerVolumePrimToSecond().

{
 //In this approximation we derive the direct cross section over Tcut=gamma energy, en after apply the Migdla correction factor 
  //used in the direct model
 
 G4double dCrossEprod=0.;
 
 const G4ElementVector* theElementVector = material->GetElementVector();
 const double* theAtomNumDensityVector = material->GetAtomicNumDensityVector();
 G4double dum=0.;
 G4double E1=kinEnergyProd,E2=kinEnergyProd*1.001;
 G4double dE=E2-E1;
 for (size_t i=0; i<material->GetNumberOfElements(); i++) { 
    G4double C1=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum ,E1);
    G4double C2=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,(*theElementVector)[i]->GetZ(),dum,E2);
    dCrossEprod += theAtomNumDensityVector[i] * (C1-C2)/dE;
   
 }
 
 //Now the Migdal correction
/*
 G4double totalEnergy = kinEnergyProj+electron_mass_c2 ;
 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy
                                             *(material->GetElectronDensity());
                         
 
 G4double MigdalFactor = 1./(1.+kp2/(kinEnergyProd*kinEnergyProd)); // its seems that the factor used in the CS compuation i the direct
                                    //model is different than the one used in the secondary sampling by a
                                    //factor (1.+kp2) To be checked!
 
 dCrossEprod*=MigdalFactor;
 */
 return dCrossEprod;
  
}

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

virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 426 of file G4AdjointBremsstrahlungModel.cc.

References AdjointCrossSection(), G4VEmModel::CrossSectionPerVolume(), G4VEmAdjointModel::GetAdjointCrossSection(), G4MaterialCutsCouple::GetMaterial(), python.hepunit::MeV, G4VEmAdjointModel::theDirectEMModel, and G4VEmAdjointModel::theDirectPrimaryPartDef.

{ 
  return AdjointCrossSection(aCouple, primEnergy,IsScatProjToProjCase);
  lastCZ=theDirectEMModel->CrossSectionPerVolume(aCouple->GetMaterial(),theDirectPrimaryPartDef,100.*MeV,100.*MeV/std::exp(1.));//this give the constant above
  return G4VEmAdjointModel::GetAdjointCrossSection(aCouple, primEnergy,IsScatProjToProjCase);
    
}

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

Definition at line 201 of file G4AdjointBremsstrahlungModel.cc.

References G4ParticleChange::AddSecondary(), CLHEP::Hep3Vector::angle(), G4VEmAdjointModel::CS_biasing_factor, G4VEmAdjointModel::currentMaterial, G4VEmAdjointModel::currentTcutForDirectSecond, G4VEmAdjointModel::DefineCurrentMaterial(), DiffCrossSectionPerVolumePrimToSecond(), python.hepunit::electron_mass_c2, fStopAndKill, G4UniformRand, G4AdjointCSManager::GetAdjointCSManager(), G4Track::GetDynamicParticle(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetMaterialCutsCouple(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGMass(), G4AdjointCSManager::GetPostStepWeightCorrection(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForScatProjToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForScatProjToProjCase(), G4DynamicParticle::GetTotalEnergy(), G4Track::GetWeight(), G4VEmAdjointModel::HighEnergyLimit, G4ParticleChange::ProposeEnergy(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeParentWeight(), G4VParticleChange::ProposeTrackStatus(), CLHEP::Hep3Vector::rotateUz(), G4VParticleChange::SetParentWeightByProcess(), G4VParticleChange::SetSecondaryWeightByProcess(), G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef, python.hepunit::twopi, and CLHEP::Hep3Vector::unit().

Referenced by SampleSecondaries().

{ 

 const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
 DefineCurrentMaterial(aTrack.GetMaterialCutsCouple());
 
 
 G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
 G4double adjointPrimTotalEnergy = theAdjointPrimary->GetTotalEnergy();
 
 if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
    return;
 }
  
 G4double projectileKinEnergy =0.;
 G4double gammaEnergy=0.;
 G4double diffCSUsed=0.; 
 if (!IsScatProjToProjCase){
    gammaEnergy=adjointPrimKinEnergy;
    G4double Emax = GetSecondAdjEnergyMaxForProdToProjCase(adjointPrimKinEnergy);
        G4double Emin=  GetSecondAdjEnergyMinForProdToProjCase(adjointPrimKinEnergy);;
    if (Emin>=Emax) return;
    projectileKinEnergy=Emin*std::pow(Emax/Emin,G4UniformRand());
    diffCSUsed=100.*CS_biasing_factor*lastCZ/projectileKinEnergy;
    
 }
 else { G4double Emax = GetSecondAdjEnergyMaxForScatProjToProjCase(adjointPrimKinEnergy);
    G4double Emin = GetSecondAdjEnergyMinForScatProjToProjCase(adjointPrimKinEnergy,currentTcutForDirectSecond);
    if (Emin>=Emax) return;
    G4double f1=(Emin-adjointPrimKinEnergy)/Emin;
    G4double f2=(Emax-adjointPrimKinEnergy)/Emax/f1;
    //G4cout<<"f1 and f2 "<<f1<<'\t'<<f2<<G4endl;
    projectileKinEnergy=adjointPrimKinEnergy/(1.-f1*std::pow(f2,G4UniformRand()));
    gammaEnergy=projectileKinEnergy-adjointPrimKinEnergy;
    diffCSUsed=lastCZ*adjointPrimKinEnergy/projectileKinEnergy/gammaEnergy;
    
 }
  
  
  
                            
 //Weight correction
 //-----------------------
 //First w_corr is set to the ratio between adjoint total CS and fwd total CS
 G4double w_corr=G4AdjointCSManager::GetAdjointCSManager()->GetPostStepWeightCorrection();

 //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
 //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.
 //Basically any other differential CS   diffCS could be used here (example Penelope). 
 
 G4double diffCS = DiffCrossSectionPerVolumePrimToSecond(currentMaterial, projectileKinEnergy, gammaEnergy);
 w_corr*=diffCS/diffCSUsed;
       
 G4double new_weight = aTrack.GetWeight()*w_corr;
 fParticleChange->SetParentWeightByProcess(false);
 fParticleChange->SetSecondaryWeightByProcess(false);
 fParticleChange->ProposeParentWeight(new_weight);
 
 //Kinematic
 //---------
 G4double projectileM0 = theAdjEquivOfDirectPrimPartDef->GetPDGMass();
 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;   
 G4double projectileP = std::sqrt(projectileP2);
 
 
 //Angle of the gamma direction with the projectile taken from G4eBremsstrahlungModel
 //------------------------------------------------
  G4double u;
  const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ;

  if (9./(9.+d) > G4UniformRand()) u = - std::log(G4UniformRand()*G4UniformRand())/a1;
     else                          u = - std::log(G4UniformRand()*G4UniformRand())/a2;

  G4double theta = u*electron_mass_c2/projectileTotalEnergy;

  G4double sint = std::sin(theta);
  G4double cost = std::cos(theta);

  G4double phi = twopi * G4UniformRand() ;
  
  G4ThreeVector projectileMomentum;
  projectileMomentum=G4ThreeVector(std::cos(phi)*sint,std::sin(phi)*sint,cost)*projectileP; //gamma frame
  if (IsScatProjToProjCase) {//the adjoint primary is the scattered e-
    G4ThreeVector gammaMomentum = (projectileTotalEnergy-adjointPrimTotalEnergy)*G4ThreeVector(0.,0.,1.);
    G4ThreeVector dirProd=projectileMomentum-gammaMomentum;
    G4double cost1 = std::cos(dirProd.angle(projectileMomentum));
    G4double sint1 =  std::sqrt(1.-cost1*cost1);
    projectileMomentum=G4ThreeVector(std::cos(phi)*sint1,std::sin(phi)*sint1,cost1)*projectileP;
  
  }
  
  projectileMomentum.rotateUz(theAdjointPrimary->GetMomentumDirection());
 
 
 
  if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
    fParticleChange->ProposeTrackStatus(fStopAndKill);
    fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
  }
  else {
    fParticleChange->ProposeEnergy(projectileKinEnergy);
    fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
    
  } 
} 

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

virtual

Implements G4VEmAdjointModel.

Definition at line 123 of file G4AdjointBremsstrahlungModel.cc.

References G4ParticleChange::AddSecondary(), CLHEP::Hep3Vector::angle(), G4VEmAdjointModel::CorrectPostStepWeight(), G4VEmAdjointModel::DefineCurrentMaterial(), python.hepunit::electron_mass_c2, fStopAndKill, G4UniformRand, G4Track::GetDynamicParticle(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetMaterialCutsCouple(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGMass(), G4DynamicParticle::GetTotalEnergy(), G4Track::GetWeight(), G4VEmAdjointModel::HighEnergyLimit, G4ParticleChange::ProposeEnergy(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), RapidSampleSecondaries(), CLHEP::Hep3Vector::rotateUz(), G4VEmAdjointModel::SampleAdjSecEnergyFromCSMatrix(), G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef, python.hepunit::twopi, CLHEP::Hep3Vector::unit(), and G4VEmAdjointModel::UseMatrix.

{
 if (!UseMatrix) return RapidSampleSecondaries(aTrack,IsScatProjToProjCase,fParticleChange); 

 const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
 DefineCurrentMaterial(aTrack.GetMaterialCutsCouple());
 
 
 G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
 G4double adjointPrimTotalEnergy = theAdjointPrimary->GetTotalEnergy();
 
 if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
    return;
 }
  
  G4double projectileKinEnergy = SampleAdjSecEnergyFromCSMatrix(adjointPrimKinEnergy,
                            IsScatProjToProjCase);
 //Weight correction
 //-----------------------                     
 CorrectPostStepWeight(fParticleChange, 
               aTrack.GetWeight(), 
               adjointPrimKinEnergy,
               projectileKinEnergy,
               IsScatProjToProjCase);   
 
 
 //Kinematic
 //---------
 G4double projectileM0 = theAdjEquivOfDirectPrimPartDef->GetPDGMass();
 G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
 G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;   
 G4double projectileP = std::sqrt(projectileP2);
 
 
 //Angle of the gamma direction with the projectile taken from G4eBremsstrahlungModel
 //------------------------------------------------
  G4double u;
  const G4double a1 = 0.625 , a2 = 3.*a1 , d = 27. ;

  if (9./(9.+d) > G4UniformRand()) u = - std::log(G4UniformRand()*G4UniformRand())/a1;
     else                          u = - std::log(G4UniformRand()*G4UniformRand())/a2;

  G4double theta = u*electron_mass_c2/projectileTotalEnergy;

  G4double sint = std::sin(theta);
  G4double cost = std::cos(theta);

  G4double phi = twopi * G4UniformRand() ;
  
  G4ThreeVector projectileMomentum;
  projectileMomentum=G4ThreeVector(std::cos(phi)*sint,std::sin(phi)*sint,cost)*projectileP; //gamma frame
  if (IsScatProjToProjCase) {//the adjoint primary is the scattered e-
    G4ThreeVector gammaMomentum = (projectileTotalEnergy-adjointPrimTotalEnergy)*G4ThreeVector(0.,0.,1.);
    G4ThreeVector dirProd=projectileMomentum-gammaMomentum;
    G4double cost1 = std::cos(dirProd.angle(projectileMomentum));
    G4double sint1 =  std::sqrt(1.-cost1*cost1);
    projectileMomentum=G4ThreeVector(std::cos(phi)*sint1,std::sin(phi)*sint1,cost1)*projectileP;
  
  }
  
  projectileMomentum.rotateUz(theAdjointPrimary->GetMomentumDirection());
 
 
 
  if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
    fParticleChange->ProposeTrackStatus(fStopAndKill);
    fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
  }
  else {
    fParticleChange->ProposeEnergy(projectileKinEnergy);
    fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
    
  } 
} 

---

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

• G4AdjointBremsstrahlungModel.hh
• G4AdjointBremsstrahlungModel.cc