Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VEmAdjointModel Class Referenceabstract

#include <G4VEmAdjointModel.hh>

Inheritance diagram for G4VEmAdjointModel:
Collaboration diagram for G4VEmAdjointModel:

Public Member Functions

 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)
 

Protected Member Functions

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

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 72 of file G4VEmAdjointModel.hh.

Constructor & Destructor Documentation

G4VEmAdjointModel::G4VEmAdjointModel ( const G4String nam)

Definition at line 41 of file G4VEmAdjointModel.cc.

41  :
42 name(nam)
43 // lowLimit(0.1*keV), highLimit(100.0*TeV), fluc(0), name(nam), pParticleChange(0)
44 {
50  currentCouple=0;
52 }
size_t RegisterEmAdjointModel(G4VEmAdjointModel *)
G4VEmModel * theDirectEMModel
G4MaterialCutsCouple * currentCouple
G4double additional_weight_correction_factor_for_post_step_outside_model
static G4AdjointCSManager * GetAdjointCSManager()

Here is the call graph for this function:

G4VEmAdjointModel::~G4VEmAdjointModel ( )
virtual

Definition at line 55 of file G4VEmAdjointModel.cc.

56 {;}

Member Function Documentation

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

Reimplemented in G4AdjointhIonisationModel, G4AdjointBremsstrahlungModel, G4AdjointComptonModel, and G4AdjointPhotoElectricModel.

Definition at line 59 of file G4VEmAdjointModel.cc.

62 {
63  DefineCurrentMaterial(aCouple);
64  preStepEnergy=primEnergy;
65 
66  std::vector<G4double>* CS_Vs_Element = &CS_Vs_ElementForProdToProjCase;
67  if (IsScatProjToProjCase) CS_Vs_Element = &CS_Vs_ElementForScatProjToProjCase;
69  this,
70  primEnergy,
72  IsScatProjToProjCase,
73  *CS_Vs_Element);
74  if (IsScatProjToProjCase) lastAdjointCSForScatProjToProjCase = lastCS;
76 
77 
78 
79  return lastCS;
80 
81 }
G4double lastAdjointCSForScatProjToProjCase
G4Material * currentMaterial
std::vector< G4double > CS_Vs_ElementForProdToProjCase
std::vector< G4double > CS_Vs_ElementForScatProjToProjCase
void DefineCurrentMaterial(const G4MaterialCutsCouple *couple)
G4double currentTcutForDirectSecond
G4double lastAdjointCSForProdToProjCase
static G4AdjointCSManager * GetAdjointCSManager()
G4double ComputeAdjointCS(G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase, std::vector< G4double > &AdjointCS_for_each_element)

Here is the call graph for this function:

Here is the caller graph for this function:

std::vector< std::vector< G4double > * > G4VEmAdjointModel::ComputeAdjointCrossSectionVectorPerAtomForScatProj ( G4double  kinEnergyProd,
G4double  Z,
G4double  A = 0.,
G4int  nbin_pro_decade = 10 
)

Definition at line 293 of file G4VEmAdjointModel.cc.

301  kinEnergyScatProjForIntegration = kinEnergyScatProj;
302 
303  //compute the vector of integrated cross sections
304  //-------------------
305 
306  G4double minEProj= GetSecondAdjEnergyMinForScatProjToProjCase(kinEnergyScatProj);
307  G4double maxEProj= GetSecondAdjEnergyMaxForScatProjToProjCase(kinEnergyScatProj);
308  G4double dEmax=maxEProj-kinEnergyScatProj;
309  G4double dEmin=GetLowEnergyLimit();
310  G4double dE1=dEmin;
311  G4double dE2=dEmin;
312 
313 
314  std::vector< double>* log_ESec_vector = new std::vector< double>();
315  std::vector< double>* log_Prob_vector = new std::vector< double>();
316  log_ESec_vector->push_back(std::log(dEmin));
317  log_Prob_vector->push_back(-50.);
318  G4int nbins=std::max( int(std::log10(dEmax/dEmin))*nbin_pro_decade,5);
319  G4double fE=std::pow(dEmax/dEmin,1./nbins);
320 
321 
322 
323 
324 
325  G4double int_cross_section=0.;
326 
327  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
328  while (dE1 <dEmax*0.9999999999999){
329  dE2=dE1*fE;
330  int_cross_section +=integral.Simpson(this,
331  &G4VEmAdjointModel::DiffCrossSectionFunction2,minEProj+dE1,std::min(minEProj+dE2,maxEProj), 5);
332  //G4cout<<"int_cross_section "<<minEProj+dE1<<'\t'<<int_cross_section<<G4endl;
333  log_ESec_vector->push_back(std::log(std::min(dE2,maxEProj-minEProj)));
334  log_Prob_vector->push_back(std::log(int_cross_section));
335  dE1=dE2;
336 
337  }
338 
339 
340  std::vector< std::vector<G4double> *> res_mat;
341  res_mat.clear();
342  if (int_cross_section >0.) {
343  res_mat.push_back(log_ESec_vector);
344  res_mat.push_back(log_Prob_vector);
345  }
346 
347  return res_mat;
348 }
G4double kinEnergyScatProjForIntegration
int G4int
Definition: G4Types.hh:78
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
double A(double temperature)
typedef int(XMLCALL *XML_NotStandaloneHandler)(void *userData)
G4double GetLowEnergyLimit()
G4double DiffCrossSectionFunction2(G4double kinEnergyProj)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double Simpson(T &typeT, F f, G4double a, G4double b, G4int n)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

std::vector< std::vector< G4double > * > G4VEmAdjointModel::ComputeAdjointCrossSectionVectorPerAtomForSecond ( G4double  kinEnergyProd,
G4double  Z,
G4double  A = 0.,
G4int  nbin_pro_decade = 10 
)

Definition at line 239 of file G4VEmAdjointModel.cc.

244 {
248  kinEnergyProdForIntegration = kinEnergyProd;
249 
250  //compute the vector of integrated cross sections
251  //-------------------
252 
253  G4double minEProj= GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
254  G4double maxEProj= GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
255  G4double E1=minEProj;
256  std::vector< double>* log_ESec_vector = new std::vector< double>();
257  std::vector< double>* log_Prob_vector = new std::vector< double>();
258  log_ESec_vector->clear();
259  log_Prob_vector->clear();
260  log_ESec_vector->push_back(std::log(E1));
261  log_Prob_vector->push_back(-50.);
262 
263  G4double E2=std::pow(10.,double( int(std::log10(minEProj)*nbin_pro_decade)+1)/nbin_pro_decade);
264  G4double fE=std::pow(10.,1./nbin_pro_decade);
265  G4double int_cross_section=0.;
266 
267  if (std::pow(fE,5.)>(maxEProj/minEProj)) fE = std::pow(maxEProj/minEProj,0.2);
268 
269  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
270  while (E1 <maxEProj*0.9999999){
271  //G4cout<<E1<<'\t'<<E2<<G4endl;
272 
273  int_cross_section +=integral.Simpson(this,
274  &G4VEmAdjointModel::DiffCrossSectionFunction1,E1,std::min(E2,maxEProj*0.99999999), 5);
275  log_ESec_vector->push_back(std::log(std::min(E2,maxEProj)));
276  log_Prob_vector->push_back(std::log(int_cross_section));
277  E1=E2;
278  E2*=fE;
279 
280  }
281  std::vector< std::vector<G4double>* > res_mat;
282  res_mat.clear();
283  if (int_cross_section >0.) {
284  res_mat.push_back(log_ESec_vector);
285  res_mat.push_back(log_Prob_vector);
286  }
287 
288  return res_mat;
289 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4double DiffCrossSectionFunction1(G4double kinEnergyProj)
double A(double temperature)
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
typedef int(XMLCALL *XML_NotStandaloneHandler)(void *userData)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double Simpson(T &typeT, F f, G4double a, G4double b, G4int n)
double G4double
Definition: G4Types.hh:76
G4double kinEnergyProdForIntegration

Here is the call graph for this function:

std::vector< std::vector< G4double > * > G4VEmAdjointModel::ComputeAdjointCrossSectionVectorPerVolumeForScatProj ( G4Material aMaterial,
G4double  kinEnergyProd,
G4int  nbin_pro_decade = 10 
)

Definition at line 403 of file G4VEmAdjointModel.cc.

408  SelectedMaterial= aMaterial;
409  kinEnergyScatProjForIntegration = kinEnergyScatProj;
410 
411  //compute the vector of integrated cross sections
412  //-------------------
413 
414  G4double minEProj= GetSecondAdjEnergyMinForScatProjToProjCase(kinEnergyScatProj);
415  G4double maxEProj= GetSecondAdjEnergyMaxForScatProjToProjCase(kinEnergyScatProj);
416 
417 
418  G4double dEmax=maxEProj-kinEnergyScatProj;
419  G4double dEmin=GetLowEnergyLimit();
420  G4double dE1=dEmin;
421  G4double dE2=dEmin;
422 
423 
424  std::vector< double>* log_ESec_vector = new std::vector< double>();
425  std::vector< double>* log_Prob_vector = new std::vector< double>();
426  log_ESec_vector->push_back(std::log(dEmin));
427  log_Prob_vector->push_back(-50.);
428  G4int nbins=std::max( int(std::log10(dEmax/dEmin))*nbin_pro_decade,5);
429  G4double fE=std::pow(dEmax/dEmin,1./nbins);
430 
431  G4double int_cross_section=0.;
432 
433  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
434  while (dE1 <dEmax*0.9999999999999){
435  dE2=dE1*fE;
436  int_cross_section +=integral.Simpson(this,
437  &G4VEmAdjointModel::DiffCrossSectionFunction2,minEProj+dE1,std::min(minEProj+dE2,maxEProj), 5);
438  log_ESec_vector->push_back(std::log(std::min(dE2,maxEProj-minEProj)));
439  log_Prob_vector->push_back(std::log(int_cross_section));
440  dE1=dE2;
441 
442  }
443 
444 
445 
446 
447 
448  std::vector< std::vector<G4double> *> res_mat;
449  res_mat.clear();
450  if (int_cross_section >0.) {
451  res_mat.push_back(log_ESec_vector);
452  res_mat.push_back(log_Prob_vector);
453  }
454 
455  return res_mat;
456 }
G4double kinEnergyScatProjForIntegration
G4Material * SelectedMaterial
int G4int
Definition: G4Types.hh:78
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
G4double GetLowEnergyLimit()
G4double DiffCrossSectionFunction2(G4double kinEnergyProj)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double Simpson(T &typeT, F f, G4double a, G4double b, G4int n)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

std::vector< std::vector< G4double > * > G4VEmAdjointModel::ComputeAdjointCrossSectionVectorPerVolumeForSecond ( G4Material aMaterial,
G4double  kinEnergyProd,
G4int  nbin_pro_decade = 10 
)

Definition at line 351 of file G4VEmAdjointModel.cc.

356  SelectedMaterial= aMaterial;
357  kinEnergyProdForIntegration = kinEnergyProd;
358  //compute the vector of integrated cross sections
359  //-------------------
360 
361  G4double minEProj= GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
362  G4double maxEProj= GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
363  G4double E1=minEProj;
364  std::vector< double>* log_ESec_vector = new std::vector< double>();
365  std::vector< double>* log_Prob_vector = new std::vector< double>();
366  log_ESec_vector->clear();
367  log_Prob_vector->clear();
368  log_ESec_vector->push_back(std::log(E1));
369  log_Prob_vector->push_back(-50.);
370 
371  G4double E2=std::pow(10.,double( int(std::log10(minEProj)*nbin_pro_decade)+1)/nbin_pro_decade);
372  G4double fE=std::pow(10.,1./nbin_pro_decade);
373  G4double int_cross_section=0.;
374 
375  if (std::pow(fE,5.)>(maxEProj/minEProj)) fE = std::pow(maxEProj/minEProj,0.2);
376 
377  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
378  while (E1 <maxEProj*0.9999999){
379 
380  int_cross_section +=integral.Simpson(this,
381  &G4VEmAdjointModel::DiffCrossSectionFunction1,E1,std::min(E2,maxEProj*0.99999999), 5);
382  log_ESec_vector->push_back(std::log(std::min(E2,maxEProj)));
383  log_Prob_vector->push_back(std::log(int_cross_section));
384  E1=E2;
385  E2*=fE;
386 
387  }
388  std::vector< std::vector<G4double>* > res_mat;
389  res_mat.clear();
390 
391  if (int_cross_section >0.) {
392  res_mat.push_back(log_ESec_vector);
393  res_mat.push_back(log_Prob_vector);
394  }
395 
396 
397 
398  return res_mat;
399 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4Material * SelectedMaterial
G4double DiffCrossSectionFunction1(G4double kinEnergyProj)
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double Simpson(T &typeT, F f, G4double a, G4double b, G4int n)
double G4double
Definition: G4Types.hh:76
G4double kinEnergyProdForIntegration

Here is the call graph for this function:

void G4VEmAdjointModel::CorrectPostStepWeight ( G4ParticleChange fParticleChange,
G4double  old_weight,
G4double  adjointPrimKinEnergy,
G4double  projectileKinEnergy,
G4bool  IsScatProjToProjCase 
)
protectedvirtual

Reimplemented in G4AdjointIonIonisationModel, and G4AdjointPhotoElectricModel.

Definition at line 634 of file G4VEmAdjointModel.cc.

639 {
640  G4double new_weight=old_weight;
641  G4double w_corr =1./CS_biasing_factor;
643 
644 
646  if ( !IsScatProjToProjCase) lastCS=lastAdjointCSForProdToProjCase;
647  if ((adjointPrimKinEnergy-preStepEnergy)/preStepEnergy>0.001){ //Is that in all cases needed???
648  G4double post_stepCS=AdjointCrossSection(currentCouple, adjointPrimKinEnergy
649  ,IsScatProjToProjCase );
650  if (post_stepCS>0 && lastCS>0) w_corr*=post_stepCS/lastCS;
651  }
652 
653  new_weight*=w_corr;
654 
655  //G4cout<<"Post step "<<new_weight<<'\t'<<w_corr<<'\t'<<old_weight<<G4endl;
656  new_weight*=projectileKinEnergy/adjointPrimKinEnergy;//This is needed due to the biasing of diff CS
657  //by the factor adjointPrimKinEnergy/projectileKinEnergy
658 
659 
660 
661  fParticleChange->SetParentWeightByProcess(false);
662  fParticleChange->SetSecondaryWeightByProcess(false);
663  fParticleChange->ProposeParentWeight(new_weight);
664 }
G4double lastAdjointCSForScatProjToProjCase
G4double GetPostStepWeightCorrection()
void ProposeParentWeight(G4double finalWeight)
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
G4MaterialCutsCouple * currentCouple
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
G4double lastAdjointCSForProdToProjCase
double G4double
Definition: G4Types.hh:76
static G4AdjointCSManager * GetAdjointCSManager()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEmAdjointModel::DefineCurrentMaterial ( const G4MaterialCutsCouple couple)

Definition at line 693 of file G4VEmAdjointModel.cc.

694 { if(couple != currentCouple) {
695  currentCouple = const_cast<G4MaterialCutsCouple*> (couple);
696  currentMaterial = const_cast<G4Material*> (couple->GetMaterial());
697  currentCoupleIndex = couple->GetIndex();
699  size_t idx=56;
700  currentTcutForDirectSecond =0.00000000001;
705  if (idx <56){
706  const std::vector<G4double>* aVec = G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(idx);
708  }
709  }
710 
711 
712  }
713 }
static G4AdjointGamma * AdjointGamma()
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
size_t GetIndex() const
Definition: G4Material.hh:262
static G4AdjointElectron * AdjointElectron()
G4Material * currentMaterial
G4MaterialCutsCouple * currentCouple
static G4ProductionCutsTable * GetProductionCutsTable()
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
G4double currentTcutForDirectSecond
static G4AdjointPositron * AdjointPositron()
const G4Material * GetMaterial() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEmAdjointModel::DefineDirectEMModel ( G4VEmModel aModel)
inline

Definition at line 193 of file G4VEmAdjointModel.hh.

193 {theDirectEMModel = aModel;}
G4VEmModel * theDirectEMModel

Here is the caller graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionFunction1 ( G4double  kinEnergyProj)
protected

Definition at line 202 of file G4VEmAdjointModel.cc.

202  {
203 
204 
205  G4double bias_factor = CS_biasing_factor*kinEnergyProdForIntegration/kinEnergyProj;
206 
207 
208  if (UseMatrixPerElement ) {
210  }
211  else {
213  }
214 }
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
virtual G4double DiffCrossSectionPerAtomPrimToSecond(G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
G4Material * SelectedMaterial
double G4double
Definition: G4Types.hh:76
G4double kinEnergyProdForIntegration

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionFunction2 ( G4double  kinEnergyProj)
protected

Definition at line 218 of file G4VEmAdjointModel.cc.

218  {
219 
220  G4double bias_factor = CS_biasing_factor*kinEnergyScatProjForIntegration/kinEnergyProj;
221  if (UseMatrixPerElement ) {
223  }
224  else {
226 
227  }
228 
229 }
G4double kinEnergyScatProjForIntegration
virtual G4double DiffCrossSectionPerVolumePrimToScatPrim(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyScatProj)
G4Material * SelectedMaterial
virtual G4double DiffCrossSectionPerAtomPrimToScatPrim(G4double kinEnergyProj, G4double kinEnergyScatProj, G4double Z, G4double A=0.)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionPerAtomPrimToScatPrim ( G4double  kinEnergyProj,
G4double  kinEnergyScatProj,
G4double  Z,
G4double  A = 0. 
)
virtual

Reimplemented in G4AdjointComptonModel.

Definition at line 145 of file G4VEmAdjointModel.cc.

150 { G4double kinEnergyProd = kinEnergyProj - kinEnergyScatProj;
151  G4double dSigmadEprod;
152  if (kinEnergyProd <=0) dSigmadEprod=0;
153  else dSigmadEprod=DiffCrossSectionPerAtomPrimToSecond(kinEnergyProj,kinEnergyProd,Z,A);
154  return dSigmadEprod;
155 
156 }
virtual G4double DiffCrossSectionPerAtomPrimToSecond(G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
double A(double temperature)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionPerAtomPrimToSecond ( G4double  kinEnergyProj,
G4double  kinEnergyProd,
G4double  Z,
G4double  A = 0. 
)
virtual

Reimplemented in G4AdjointhIonisationModel, G4AdjointIonIonisationModel, G4AdjointComptonModel, and G4AdjointeIonisationModel.

Definition at line 113 of file G4VEmAdjointModel.cc.

118 {
119  G4double dSigmadEprod=0;
120  G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
121  G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
122 
123 
124  if (kinEnergyProj>Emin_proj && kinEnergyProj<=Emax_proj){ //the produced particle should have a kinetic energy smaller than the projectile
125 
126  /*G4double Tmax=kinEnergyProj;
127  if (second_part_of_same_type) Tmax = kinEnergyProj/2.;*/
128 
129  G4double E1=kinEnergyProd;
130  G4double E2=kinEnergyProd*1.000001;
131  G4double dE=(E2-E1);
134 
135  dSigmadEprod=(sigma1-sigma2)/dE;
136  }
137  return dSigmadEprod;
138 
139 
140 
141 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4ParticleDefinition * theDirectPrimaryPartDef
static const G4double dE
G4VEmModel * theDirectEMModel
double A(double temperature)
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:321
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionPerVolumeFunctionForIntegrationOverEkinProj ( G4double  EkinProd)
protected

Definition at line 233 of file G4VEmAdjointModel.cc.

234 {
236 }
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
G4Material * SelectedMaterial
G4double kinEnergyProjForIntegration

Here is the call graph for this function:

G4double G4VEmAdjointModel::DiffCrossSectionPerVolumePrimToScatPrim ( const G4Material aMaterial,
G4double  kinEnergyProj,
G4double  kinEnergyScatProj 
)
virtual

Definition at line 189 of file G4VEmAdjointModel.cc.

193 { G4double kinEnergyProd = kinEnergyProj - kinEnergyScatProj;
194  G4double dSigmadEprod;
195  if (kinEnergyProd <=0) dSigmadEprod=0;
196  else dSigmadEprod=DiffCrossSectionPerVolumePrimToSecond(aMaterial,kinEnergyProj,kinEnergyProd);
197  return dSigmadEprod;
198 
199 }
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

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

Reimplemented in G4AdjointBremsstrahlungModel.

Definition at line 161 of file G4VEmAdjointModel.cc.

165 {
166  G4double dSigmadEprod=0;
167  G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
168  G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
169 
170 
171  if (kinEnergyProj>Emin_proj && kinEnergyProj<=Emax_proj){
172  /*G4double Tmax=kinEnergyProj;
173  if (second_part_of_same_type) Tmax = kinEnergyProj/2.;*/
174  G4double E1=kinEnergyProd;
175  G4double E2=kinEnergyProd*1.0001;
176  G4double dE=(E2-E1);
177  G4double sigma1=theDirectEMModel->CrossSectionPerVolume(aMaterial,theDirectPrimaryPartDef,kinEnergyProj,E1,1.e20);
178  G4double sigma2=theDirectEMModel->CrossSectionPerVolume(aMaterial,theDirectPrimaryPartDef,kinEnergyProj,E2,1.e20);
179  dSigmadEprod=(sigma1-sigma2)/dE;
180  }
181  return dSigmadEprod;
182 
183 
184 
185 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:257
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4ParticleDefinition * theDirectPrimaryPartDef
static const G4double dE
G4VEmModel * theDirectEMModel
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

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

Reimplemented in G4AdjointBremsstrahlungModel, G4AdjointComptonModel, and G4AdjointPhotoElectricModel.

Definition at line 84 of file G4VEmAdjointModel.cc.

87 {
88  return AdjointCrossSection(aCouple, primEnergy,
89  IsScatProjToProjCase);
90 
91  /*
92  //To continue
93  DefineCurrentMaterial(aCouple);
94  preStepEnergy=primEnergy;
95  if (IsScatProjToProjCase){
96  G4double ekin=primEnergy*mass_ratio_projectile;
97  lastCS = G4AdjointCSManager::GetAdjointCSManager()->GetAdjointSigma(ekin, model_index,true, aCouple);
98  lastAdjointCSForScatProjToProjCase = lastCS;
99  //G4cout<<ekin<<std::endl;
100  }
101  else {
102  G4double ekin=primEnergy*mass_ratio_product;
103  lastCS = G4AdjointCSManager::GetAdjointCSManager()->GetAdjointSigma(ekin, model_index,false, aCouple);
104  lastAdjointCSForProdToProjCase = lastCS;
105  //G4cout<<ekin<<std::endl;
106  }
107  return lastCS;
108  */
109 }
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)

Here is the call graph for this function:

Here is the caller graph for this function:

G4ParticleDefinition* G4VEmAdjointModel::GetAdjointEquivalentOfDirectPrimaryParticleDefinition ( )
inline

Definition at line 181 of file G4VEmAdjointModel.hh.

G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
G4ParticleDefinition* G4VEmAdjointModel::GetAdjointEquivalentOfDirectSecondaryParticleDefinition ( )
inline

Definition at line 183 of file G4VEmAdjointModel.hh.

G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
G4bool G4VEmAdjointModel::GetApplyCutInRange ( )
inline

Definition at line 214 of file G4VEmAdjointModel.hh.

214 { return ApplyCutInRange;}

Here is the caller graph for this function:

G4double G4VEmAdjointModel::GetHighEnergyLimit ( )
inline

Definition at line 185 of file G4VEmAdjointModel.hh.

185 {return HighEnergyLimit;}
G4double G4VEmAdjointModel::GetLowEnergyLimit ( )
inline

Definition at line 187 of file G4VEmAdjointModel.hh.

187 {return LowEnergyLimit;}

Here is the caller graph for this function:

G4String G4VEmAdjointModel::GetName ( void  )
inline

Definition at line 216 of file G4VEmAdjointModel.hh.

216 { return name;}

Here is the caller graph for this function:

G4double G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented in G4AdjointIonIonisationModel, and G4AdjointhIonisationModel.

Definition at line 681 of file G4VEmAdjointModel.cc.

682 { return HighEnergyLimit;
683 }

Here is the caller graph for this function:

G4double G4VEmAdjointModel::GetSecondAdjEnergyMaxForScatProjToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented in G4AdjointIonIonisationModel, G4AdjointhIonisationModel, and G4AdjointComptonModel.

Definition at line 667 of file G4VEmAdjointModel.cc.

668 { G4double maxEProj= HighEnergyLimit;
669  if (second_part_of_same_type) maxEProj=std::min(kinEnergyScatProj*2.,HighEnergyLimit);
670  return maxEProj;
671 }
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented in G4AdjointIonIonisationModel, G4AdjointhIonisationModel, and G4AdjointComptonModel.

Definition at line 686 of file G4VEmAdjointModel.cc.

687 { G4double minEProj=PrimAdjEnergy;
688  if (second_part_of_same_type) minEProj=PrimAdjEnergy*2.;
689  return minEProj;
690 }
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4VEmAdjointModel::GetSecondAdjEnergyMinForScatProjToProjCase ( G4double  PrimAdjEnergy,
G4double  Tcut = 0 
)
virtual

Reimplemented in G4AdjointIonIonisationModel, and G4AdjointhIonisationModel.

Definition at line 674 of file G4VEmAdjointModel.cc.

675 { G4double Emin=PrimAdjEnergy;
676  if (ApplyCutInRange) Emin=PrimAdjEnergy+Tcut;
677  return Emin;
678 }
static const G4double Emin
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4bool G4VEmAdjointModel::GetSecondPartOfSameType ( )
inline

Definition at line 203 of file G4VEmAdjointModel.hh.

203 {return second_part_of_same_type;}
G4bool G4VEmAdjointModel::GetUseMatrix ( )
inline

Definition at line 211 of file G4VEmAdjointModel.hh.

211 {return UseMatrix;}

Here is the caller graph for this function:

G4bool G4VEmAdjointModel::GetUseMatrixPerElement ( )
inline

Definition at line 212 of file G4VEmAdjointModel.hh.

212 { return UseMatrixPerElement;}

Here is the caller graph for this function:

G4bool G4VEmAdjointModel::GetUseOnlyOneMatrixForAllElements ( )
inline

Definition at line 213 of file G4VEmAdjointModel.hh.

Here is the caller graph for this function:

G4double G4VEmAdjointModel::SampleAdjSecEnergyFromCSMatrix ( size_t  MatrixIndex,
G4double  prim_energy,
G4bool  IsScatProjToProjCase 
)
protected

Definition at line 459 of file G4VEmAdjointModel.cc.

460 {
461 
462 
463  G4AdjointCSMatrix* theMatrix= (*pOnCSMatrixForProdToProjBackwardScattering)[MatrixIndex];
464  if (IsScatProjToProjCase) theMatrix= (*pOnCSMatrixForScatProjToProjBackwardScattering)[MatrixIndex];
465  std::vector< double>* theLogPrimEnergyVector = theMatrix->GetLogPrimEnergyVector();
466 
467  if (theLogPrimEnergyVector->size() ==0){
468  G4cout<<"No data are contained in the given AdjointCSMatrix!"<<G4endl;
469  G4cout<<"The sampling procedure will be stopped."<<G4endl;
470  return 0.;
471 
472  }
473 
475  G4double aLogPrimEnergy = std::log(aPrimEnergy);
476  size_t ind =theInterpolator->FindPositionForLogVector(aLogPrimEnergy,*theLogPrimEnergyVector);
477 
478 
479  G4double aLogPrimEnergy1,aLogPrimEnergy2;
480  G4double aLogCS1,aLogCS2;
481  G4double log01,log02;
482  std::vector< double>* aLogSecondEnergyVector1 =0;
483  std::vector< double>* aLogSecondEnergyVector2 =0;
484  std::vector< double>* aLogProbVector1=0;
485  std::vector< double>* aLogProbVector2=0;
486  std::vector< size_t>* aLogProbVectorIndex1=0;
487  std::vector< size_t>* aLogProbVectorIndex2=0;
488 
489  theMatrix->GetData(ind, aLogPrimEnergy1,aLogCS1,log01, aLogSecondEnergyVector1,aLogProbVector1,aLogProbVectorIndex1);
490  theMatrix->GetData(ind+1, aLogPrimEnergy2,aLogCS2,log02, aLogSecondEnergyVector2,aLogProbVector2,aLogProbVectorIndex2);
491 
492  G4double rand_var = G4UniformRand();
493  G4double log_rand_var= std::log(rand_var);
494  G4double log_Tcut =std::log(currentTcutForDirectSecond);
495  G4double Esec=0;
496  G4double log_dE1,log_dE2;
497  G4double log_rand_var1,log_rand_var2;
498  G4double log_E1,log_E2;
499  log_rand_var1=log_rand_var;
500  log_rand_var2=log_rand_var;
501 
502  G4double Emin=0.;
503  G4double Emax=0.;
504  if (theMatrix->IsScatProjToProjCase()){ //case where Tcut plays a role
507  G4double dE=0;
508  if (Emin < Emax ){
509  if (ApplyCutInRange) {
510  if (second_part_of_same_type && currentTcutForDirectSecond>aPrimEnergy) return aPrimEnergy;
511 
512  log_rand_var1=log_rand_var+theInterpolator->InterpolateForLogVector(log_Tcut,*aLogSecondEnergyVector1,*aLogProbVector1);
513  log_rand_var2=log_rand_var+theInterpolator->InterpolateForLogVector(log_Tcut,*aLogSecondEnergyVector2,*aLogProbVector2);
514 
515  }
516  log_dE1 = theInterpolator->Interpolate(log_rand_var1,*aLogProbVector1,*aLogSecondEnergyVector1,"Lin");
517  log_dE2 = theInterpolator->Interpolate(log_rand_var2,*aLogProbVector2,*aLogSecondEnergyVector2,"Lin");
518  dE=std::exp(theInterpolator->LinearInterpolation(aLogPrimEnergy,aLogPrimEnergy1,aLogPrimEnergy2,log_dE1,log_dE2));
519  }
520 
521  Esec = aPrimEnergy +dE;
522  Esec=std::max(Esec,Emin);
523  Esec=std::min(Esec,Emax);
524 
525  }
526  else { //Tcut condition is already full-filled
527 
528  log_E1 = theInterpolator->Interpolate(log_rand_var,*aLogProbVector1,*aLogSecondEnergyVector1,"Lin");
529  log_E2 = theInterpolator->Interpolate(log_rand_var,*aLogProbVector2,*aLogSecondEnergyVector2,"Lin");
530 
531  Esec = std::exp(theInterpolator->LinearInterpolation(aLogPrimEnergy,aLogPrimEnergy1,aLogPrimEnergy2,log_E1,log_E2));
532  Emin=GetSecondAdjEnergyMinForProdToProjCase(aPrimEnergy);
533  Emax=GetSecondAdjEnergyMaxForProdToProjCase(aPrimEnergy);
534  Esec=std::max(Esec,Emin);
535  Esec=std::min(Esec,Emax);
536 
537  }
538 
539  return Esec;
540 
541 
542 
543 
544 
545 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4bool GetData(unsigned int i, G4double &aPrimEnergy, G4double &aCS, G4double &log0, std::vector< double > *&aLogSecondEnergyVector, std::vector< double > *&aLogProbVector, std::vector< size_t > *&aLogProbVectorIndex)
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
G4double InterpolateForLogVector(G4double &x, std::vector< G4double > &x_vec, std::vector< G4double > &y_vec)
G4double Interpolate(G4double &x, std::vector< G4double > &x_vec, std::vector< G4double > &y_vec, G4String InterPolMethod="Log")
static const G4double dE
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
std::vector< double > * GetLogPrimEnergyVector()
size_t FindPositionForLogVector(G4double &x, std::vector< G4double > &x_vec)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
G4double LinearInterpolation(G4double &x, G4double &x1, G4double &x2, G4double &y1, G4double &y2)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static const G4double Emin
static G4AdjointInterpolator * GetInstance()
static const G4double Emax
#define G4endl
Definition: G4ios.hh:61
G4double currentTcutForDirectSecond
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VEmAdjointModel::SampleAdjSecEnergyFromCSMatrix ( G4double  prim_energy,
G4bool  IsScatProjToProjCase 
)
protected

Definition at line 549 of file G4VEmAdjointModel.cc.

550 { SelectCSMatrix(IsScatProjToProjCase);
551  return SampleAdjSecEnergyFromCSMatrix(indexOfUsedCrossSectionMatrix, aPrimEnergy, IsScatProjToProjCase);
552 }
G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
size_t indexOfUsedCrossSectionMatrix
void SelectCSMatrix(G4bool IsScatProjToProjCase)

Here is the call graph for this function:

G4double G4VEmAdjointModel::SampleAdjSecEnergyFromDiffCrossSectionPerAtom ( G4double  prim_energy,
G4bool  IsScatProjToProjCase 
)
protectedvirtual

Definition at line 581 of file G4VEmAdjointModel.cc.

582 {
583  // here we try to use the rejection method
584  //-----------------------------------------
585 
586  const G4int iimax = 1000;
587  G4double E=0;
588  G4double x,xmin,greject,q;
589  if ( IsScatProjToProjCase){
592  xmin=Emin/Emax;
593  G4double grejmax = DiffCrossSectionPerAtomPrimToScatPrim(Emin,prim_energy,1)*prim_energy;
594 
595  G4int ii =0;
596  do {
597  q = G4UniformRand();
598  x = 1./(q*(1./xmin -1.) +1.);
599  E=x*Emax;
600  greject = DiffCrossSectionPerAtomPrimToScatPrim( E,prim_energy ,1)*prim_energy;
601  ++ii;
602  if(ii >= iimax) { break; }
603  }
604  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
605  while( greject < G4UniformRand()*grejmax );
606 
607  }
608  else {
611  xmin=Emin/Emax;
612  G4double grejmax = DiffCrossSectionPerAtomPrimToSecond(Emin,prim_energy,1);
613  G4int ii =0;
614  do {
615  q = G4UniformRand();
616  x = std::pow(xmin, q);
617  E=x*Emax;
618  greject = DiffCrossSectionPerAtomPrimToSecond( E,prim_energy ,1);
619  ++ii;
620  if(ii >= iimax) { break; }
621  }
622  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
623  while( greject < G4UniformRand()*grejmax );
624 
625 
626 
627  }
628 
629  return E;
630 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
tuple x
Definition: test.py:50
virtual G4double DiffCrossSectionPerAtomPrimToSecond(G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
int G4int
Definition: G4Types.hh:78
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
#define G4UniformRand()
Definition: Randomize.hh:97
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
virtual G4double DiffCrossSectionPerAtomPrimToScatPrim(G4double kinEnergyProj, G4double kinEnergyScatProj, G4double Z, G4double A=0.)
static const G4double Emin
static const G4double Emax
G4double currentTcutForDirectSecond
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

virtual void G4VEmAdjointModel::SampleSecondaries ( const G4Track aTrack,
G4bool  IsScatProjToProjCase,
G4ParticleChange fParticleChange 
)
pure virtual
void G4VEmAdjointModel::SelectCSMatrix ( G4bool  IsScatProjToProjCase)
protected

Definition at line 555 of file G4VEmAdjointModel.cc.

556 {
559  else if (!UseOnlyOneMatrixForAllElements) { //Select Material
560  std::vector<G4double>* CS_Vs_Element = &CS_Vs_ElementForScatProjToProjCase;
562  if ( !IsScatProjToProjCase) {
563  CS_Vs_Element = &CS_Vs_ElementForProdToProjCase;
565  }
566  G4double rand_var= G4UniformRand();
567  G4double SumCS=0.;
568  size_t ind=0;
569  for (size_t i=0;i<CS_Vs_Element->size();i++){
570  SumCS+=(*CS_Vs_Element)[i];
571  if (rand_var<=SumCS/lastCS){
572  ind=i;
573  break;
574  }
575  }
577  }
578 }
G4double lastAdjointCSForScatProjToProjCase
G4Material * currentMaterial
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:202
size_t indexOfUsedCrossSectionMatrix
#define G4UniformRand()
Definition: Randomize.hh:97
std::vector< G4double > CS_Vs_ElementForProdToProjCase
std::vector< G4double > CS_Vs_ElementForScatProjToProjCase
size_t GetIndex() const
Definition: G4Element.hh:182
G4bool UseOnlyOneMatrixForAllElements
G4double lastAdjointCSForProdToProjCase
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEmAdjointModel::SetAdditionalWeightCorrectionFactorForPostStepOutsideModel ( G4double  factor)
inline

Definition at line 220 of file G4VEmAdjointModel.hh.

Here is the caller graph for this function:

void G4VEmAdjointModel::SetAdjointEquivalentOfDirectPrimaryParticleDefinition ( G4ParticleDefinition aPart)

Definition at line 729 of file G4VEmAdjointModel.cc.

730 {
734  if (theAdjEquivOfDirectPrimPartDef->GetParticleName() =="adj_gamma")
736 }
G4ParticleDefinition * theDirectPrimaryPartDef
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
const G4String & GetParticleName() const
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
static G4Electron * Electron()
Definition: G4Electron.cc:94

Here is the call graph for this function:

void G4VEmAdjointModel::SetAdjointEquivalentOfDirectSecondaryParticleDefinition ( G4ParticleDefinition aPart)
inline

Definition at line 197 of file G4VEmAdjointModel.hh.

197  {
199  }
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
void G4VEmAdjointModel::SetApplyCutInRange ( G4bool  aBool)
inline

Definition at line 210 of file G4VEmAdjointModel.hh.

210 { ApplyCutInRange = aBool;}

Here is the caller graph for this function:

void G4VEmAdjointModel::SetCorrectWeightForPostStepInModel ( G4bool  aBool)
inline

Definition at line 219 of file G4VEmAdjointModel.hh.

Here is the caller graph for this function:

virtual void G4VEmAdjointModel::SetCSBiasingFactor ( G4double  aVal)
inlinevirtual

Definition at line 217 of file G4VEmAdjointModel.hh.

217 {CS_biasing_factor = aVal;}

Here is the caller graph for this function:

void G4VEmAdjointModel::SetCSMatrices ( std::vector< G4AdjointCSMatrix * > *  Vec1CSMatrix,
std::vector< G4AdjointCSMatrix * > *  Vec2CSMatrix 
)
inline

Definition at line 174 of file G4VEmAdjointModel.hh.

174  {
177 
178 
179  };
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForProdToProjBackwardScattering
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForScatProjToProjBackwardScattering

Here is the caller graph for this function:

void G4VEmAdjointModel::SetHighEnergyLimit ( G4double  aVal)

Definition at line 716 of file G4VEmAdjointModel.cc.

717 { HighEnergyLimit=aVal;
719 }
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4VEmModel * theDirectEMModel

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEmAdjointModel::SetLowEnergyLimit ( G4double  aVal)

Definition at line 722 of file G4VEmAdjointModel.cc.

723 {
724  LowEnergyLimit=aVal;
726 }
G4VEmModel * theDirectEMModel
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:739

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VEmAdjointModel::SetSecondPartOfSameType ( G4bool  aBool)
inline

Definition at line 201 of file G4VEmAdjointModel.hh.

201 {second_part_of_same_type =aBool;}

Here is the caller graph for this function:

void G4VEmAdjointModel::SetUseMatrix ( G4bool  aBool)
inline

Definition at line 205 of file G4VEmAdjointModel.hh.

205 { UseMatrix = aBool;}

Here is the caller graph for this function:

void G4VEmAdjointModel::SetUseMatrixPerElement ( G4bool  aBool)
inline

Definition at line 207 of file G4VEmAdjointModel.hh.

207 { UseMatrixPerElement = aBool;}

Here is the caller graph for this function:

void G4VEmAdjointModel::SetUseOnlyOneMatrixForAllElements ( G4bool  aBool)
inline

Definition at line 208 of file G4VEmAdjointModel.hh.

G4bool UseOnlyOneMatrixForAllElements

Here is the caller graph for this function:

Member Data Documentation

G4double G4VEmAdjointModel::additional_weight_correction_factor_for_post_step_outside_model
protected

Definition at line 345 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::ApplyCutInRange
protected

Definition at line 312 of file G4VEmAdjointModel.hh.

G4int G4VEmAdjointModel::ASelectedNucleus
protected

Definition at line 272 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::correct_weight_for_post_step_in_model
protected

Definition at line 344 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::CS_biasing_factor
protected

Definition at line 327 of file G4VEmAdjointModel.hh.

std::vector<G4double> G4VEmAdjointModel::CS_Vs_ElementForProdToProjCase
protected

Definition at line 286 of file G4VEmAdjointModel.hh.

std::vector<G4double> G4VEmAdjointModel::CS_Vs_ElementForScatProjToProjCase
protected

Definition at line 285 of file G4VEmAdjointModel.hh.

G4MaterialCutsCouple* G4VEmAdjointModel::currentCouple
protected

Definition at line 307 of file G4VEmAdjointModel.hh.

size_t G4VEmAdjointModel::currentCoupleIndex
protected

Definition at line 309 of file G4VEmAdjointModel.hh.

G4Material* G4VEmAdjointModel::currentMaterial
protected

Definition at line 306 of file G4VEmAdjointModel.hh.

size_t G4VEmAdjointModel::currentMaterialIndex
protected

Definition at line 308 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::currentTcutForDirectPrim
protected

Definition at line 310 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::currentTcutForDirectSecond
protected

Definition at line 311 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::HighEnergyLimit
protected

Definition at line 322 of file G4VEmAdjointModel.hh.

size_t G4VEmAdjointModel::indexOfUsedCrossSectionMatrix
protected

Definition at line 337 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::kinEnergyProdForIntegration
protected

Definition at line 275 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::kinEnergyProjForIntegration
protected

Definition at line 277 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::kinEnergyScatProjForIntegration
protected

Definition at line 276 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::lastAdjointCSForProdToProjCase
protected

Definition at line 290 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::lastAdjointCSForScatProjToProjCase
protected

Definition at line 289 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::lastCS
protected

Definition at line 288 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::LowEnergyLimit
protected

Definition at line 323 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::mass_ratio_product
protected

Definition at line 316 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::mass_ratio_projectile
protected

Definition at line 317 of file G4VEmAdjointModel.hh.

size_t G4VEmAdjointModel::model_index
protected

Definition at line 339 of file G4VEmAdjointModel.hh.

const G4String G4VEmAdjointModel::name
protected

Definition at line 267 of file G4VEmAdjointModel.hh.

std::vector< G4AdjointCSMatrix* >* G4VEmAdjointModel::pOnCSMatrixForProdToProjBackwardScattering
protected

Definition at line 283 of file G4VEmAdjointModel.hh.

std::vector< G4AdjointCSMatrix* >* G4VEmAdjointModel::pOnCSMatrixForScatProjToProjBackwardScattering
protected

Definition at line 284 of file G4VEmAdjointModel.hh.

G4VParticleChange* G4VEmAdjointModel::pParticleChange
protected

Definition at line 259 of file G4VEmAdjointModel.hh.

G4double G4VEmAdjointModel::preStepEnergy
protected

Definition at line 302 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::second_part_of_same_type
protected

Definition at line 298 of file G4VEmAdjointModel.hh.

G4Material* G4VEmAdjointModel::SelectedMaterial
protected

Definition at line 274 of file G4VEmAdjointModel.hh.

G4ParticleDefinition* G4VEmAdjointModel::theAdjEquivOfDirectPrimPartDef
protected

Definition at line 295 of file G4VEmAdjointModel.hh.

G4ParticleDefinition* G4VEmAdjointModel::theAdjEquivOfDirectSecondPartDef
protected

Definition at line 296 of file G4VEmAdjointModel.hh.

G4VEmModel* G4VEmAdjointModel::theDirectEMModel
protected

Definition at line 258 of file G4VEmAdjointModel.hh.

G4ParticleDefinition* G4VEmAdjointModel::theDirectPrimaryPartDef
protected

Definition at line 297 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::UseMatrix
protected

Definition at line 331 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::UseMatrixPerElement
protected

Definition at line 332 of file G4VEmAdjointModel.hh.

G4bool G4VEmAdjointModel::UseOnlyOneMatrixForAllElements
protected

Definition at line 333 of file G4VEmAdjointModel.hh.

G4int G4VEmAdjointModel::ZSelectedNucleus
protected

Definition at line 273 of file G4VEmAdjointModel.hh.


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