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

#include <G4PolarizedCompton.hh>

Inheritance diagram for G4PolarizedCompton:
Collaboration diagram for G4PolarizedCompton:

Public Member Functions

 G4PolarizedCompton (const G4String &processName="pol-compt", G4ProcessType type=fElectromagnetic)
 
virtual ~G4PolarizedCompton ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &) override
 
virtual void PrintInfo () override
 
void SetModel (const G4String &name)
 
- Public Member Functions inherited from G4VEmProcess
 G4VEmProcess (const G4String &name, G4ProcessType type=fElectromagnetic)
 
virtual ~G4VEmProcess ()
 
virtual void ProcessDescription (std::ostream &outFile) const
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &) override
 
virtual void StartTracking (G4Track *) override
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &) override
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
 
G4double CrossSectionPerVolume (G4double kineticEnergy, const G4MaterialCutsCouple *couple)
 
G4double ComputeCrossSectionPerAtom (G4double kineticEnergy, G4double Z, G4double A=0., G4double cut=0.0)
 
G4double MeanFreePath (const G4Track &track)
 
G4double GetLambda (G4double &kinEnergy, const G4MaterialCutsCouple *couple)
 
void SetLambdaBinning (G4int nbins)
 
void SetMinKinEnergy (G4double e)
 
void SetMinKinEnergyPrim (G4double e)
 
void SetMaxKinEnergy (G4double e)
 
G4PhysicsTableLambdaTable () const
 
G4PhysicsTableLambdaTablePrim () const
 
const G4ParticleDefinitionParticle () const
 
const G4ParticleDefinitionSecondaryParticle () const
 
G4VEmModelSelectModelForMaterial (G4double kinEnergy, size_t &idxRegion) const
 
void AddEmModel (G4int, G4VEmModel *, const G4Region *region=nullptr)
 
G4VEmModelEmModel (G4int index=1) const
 
void SetEmModel (G4VEmModel *, G4int index=1)
 
void UpdateEmModel (const G4String &, G4double, G4double)
 
G4int GetNumberOfModels () const
 
G4int GetNumberOfRegionModels (size_t couple_index) const
 
G4VEmModelGetRegionModel (G4int idx, size_t couple_index) const
 
G4VEmModelGetModelByIndex (G4int idx=0, G4bool ver=false) const
 
const G4VEmModelGetCurrentModel () const
 
const G4ElementGetCurrentElement () const
 
void SetCrossSectionBiasingFactor (G4double f, G4bool flag=true)
 
G4double CrossSectionBiasingFactor () const
 
void ActivateForcedInteraction (G4double length=0.0, const G4String &r="", G4bool flag=true)
 
void ActivateSecondaryBiasing (const G4String &region, G4double factor, G4double energyLimit)
 
void SetIntegral (G4bool val)
 
void SetBuildTableFlag (G4bool val)
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual void InitialiseProcess (const G4ParticleDefinition *) override
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &) override
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition) override
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
 
- Protected Member Functions inherited from G4VEmProcess
virtual G4double MinPrimaryEnergy (const G4ParticleDefinition *, const G4Material *)
 
G4VEmModelSelectModel (G4double &kinEnergy, size_t index)
 
G4PhysicsVectorLambdaPhysicsVector (const G4MaterialCutsCouple *)
 
G4int LambdaBinning () const
 
G4double MinKinEnergy () const
 
G4double MaxKinEnergy () const
 
G4double PolarAngleLimit () const
 
G4bool IsIntegral () const
 
G4double RecalculateLambda (G4double kinEnergy, const G4MaterialCutsCouple *couple)
 
G4ParticleChangeForGammaGetParticleChange ()
 
void SetParticle (const G4ParticleDefinition *p)
 
void SetSecondaryParticle (const G4ParticleDefinition *p)
 
size_t CurrentMaterialCutsCoupleIndex () const
 
G4double GetGammaEnergyCut ()
 
G4double GetElectronEnergyCut ()
 
void SetStartFromNullFlag (G4bool val)
 
void SetSplineFlag (G4bool val)
 
const G4ElementGetTargetElement () const
 
const G4IsotopeGetTargetIsotope () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VEmProcess
G4ParticleChangeForGamma fParticleChange
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 72 of file G4PolarizedCompton.hh.

Constructor & Destructor Documentation

G4PolarizedCompton::G4PolarizedCompton ( const G4String processName = "pol-compt",
G4ProcessType  type = fElectromagnetic 
)
explicit

Definition at line 73 of file G4PolarizedCompton.cc.

74  :
75  G4VEmProcess (processName, type),
76  buildAsymmetryTable(true),
77  useAsymmetryTable(true),
78  isInitialised(false),
79  mType(10),
80  targetPolarization(0.0,0.0,0.0)
81 {
83  SetBuildTableFlag(true);
87  SetSplineFlag(true);
88  emModel = nullptr;
89 }
G4VEmProcess(const G4String &name, G4ProcessType type=fElectromagnetic)
Definition: G4VEmProcess.cc:91
void SetBuildTableFlag(G4bool val)
void SetSplineFlag(G4bool val)
void SetStartFromNullFlag(G4bool val)
void SetMinKinEnergyPrim(G4double e)
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
void SetSecondaryParticle(const G4ParticleDefinition *p)
static G4Electron * Electron()
Definition: G4Electron.cc:94
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

G4PolarizedCompton::~G4PolarizedCompton ( )
virtual

Definition at line 93 of file G4PolarizedCompton.cc.

94 {
95  CleanTable();
96 }

Member Function Documentation

void G4PolarizedCompton::BuildPhysicsTable ( const G4ParticleDefinition part)
overrideprotectedvirtual

Reimplemented from G4VEmProcess.

Definition at line 278 of file G4PolarizedCompton.cc.

279 {
280  // *** build (unpolarized) cross section tables (Lambda)
282  if(buildAsymmetryTable && emModel) {
283  G4bool isMaster = true;
284  const G4PolarizedCompton* masterProcess =
285  static_cast<const G4PolarizedCompton*>(GetMasterProcess());
286  if(masterProcess && masterProcess != this) { isMaster = false; }
287  if(isMaster) { BuildAsymmetryTable(part); }
288  }
289 }
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
virtual void BuildPhysicsTable(const G4ParticleDefinition &) override
bool G4bool
Definition: G4Types.hh:79

Here is the call graph for this function:

G4double G4PolarizedCompton::GetMeanFreePath ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
)
overrideprotectedvirtual

Reimplemented from G4VEmProcess.

Definition at line 155 of file G4PolarizedCompton.cc.

158 {
159  // *** get unploarised mean free path from lambda table ***
160  G4double mfp = G4VEmProcess::GetMeanFreePath(aTrack, previousStepSize, condition);
161 
162  if (theAsymmetryTable && useAsymmetryTable && mfp < DBL_MAX) {
163  mfp *= ComputeSaturationFactor(aTrack);
164  }
165  if (verboseLevel>=2) {
166  G4cout << "G4PolarizedCompton::MeanFreePath: " << mfp / mm << " mm " << G4endl;
167  }
168  return mfp;
169 }
G4double condition(const G4ErrorSymMatrix &m)
static constexpr double mm
Definition: G4SIunits.hh:115
G4int verboseLevel
Definition: G4VProcess.hh:368
virtual G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

void G4PolarizedCompton::InitialiseProcess ( const G4ParticleDefinition )
overrideprotectedvirtual

Implements G4VEmProcess.

Definition at line 118 of file G4PolarizedCompton.cc.

119 {
120  if(!isInitialised) {
121  isInitialised = true;
122  if(0 == mType) {
123  if(!EmModel(1)) { SetEmModel(new G4KleinNishinaCompton(), 1); }
124  } else {
125  emModel = new G4PolarizedComptonModel();
126  SetEmModel(emModel, 1);
127  }
129  EmModel(1)->SetLowEnergyLimit(param->MinKinEnergy());
131  AddEmModel(1, EmModel(1));
132  }
133 }
G4double MaxKinEnergy() const
G4VEmModel * EmModel(G4int index=1) const
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
void SetEmModel(G4VEmModel *, G4int index=1)
G4double MinKinEnergy() const
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
static G4EmParameters * Instance()
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:739

Here is the call graph for this function:

G4bool G4PolarizedCompton::IsApplicable ( const G4ParticleDefinition p)
overridevirtual

Implements G4VEmProcess.

Definition at line 111 of file G4PolarizedCompton.cc.

112 {
113  return (&p == G4Gamma::Gamma());
114 }
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86

Here is the call graph for this function:

G4double G4PolarizedCompton::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
overrideprotectedvirtual

Reimplemented from G4VEmProcess.

Definition at line 173 of file G4PolarizedCompton.cc.

177 {
178  // save previous values
181 
182  // *** compute unpolarized step limit ***
183  // this changes theNumberOfInteractionLengthLeft and currentInteractionLength
185  previousStepSize,
186  condition);
187  G4double x0 = x;
188  G4double satFact = 1.0;
189 
190  // *** add corrections on polarisation ***
191  if (theAsymmetryTable && useAsymmetryTable && x < DBL_MAX) {
192  satFact = ComputeSaturationFactor(aTrack);
193  G4double curLength = currentInteractionLength*satFact;
194  G4double prvLength = iLength*satFact;
195  if(nLength > 0.0) {
197  std::max(nLength - previousStepSize/prvLength, 0.0);
198  }
199  x = theNumberOfInteractionLengthLeft * curLength;
200  }
201  if (verboseLevel>=2) {
202  G4cout << "G4PolarizedCompton::PostStepGPIL: "
203  << std::setprecision(8) << x/mm << " mm;" << G4endl
204  << " unpolarized value: "
205  << std::setprecision(8) << x0/mm << " mm." << G4endl;
206  }
207  return x;
208 }
G4double condition(const G4ErrorSymMatrix &m)
static constexpr double mm
Definition: G4SIunits.hh:115
G4int verboseLevel
Definition: G4VProcess.hh:368
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
tuple x
Definition: test.py:50
G4GLOB_DLL std::ostream G4cout
G4double currentInteractionLength
Definition: G4VProcess.hh:297
T max(const T t1, const T t2)
brief Return the largest of the two arguments
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) override

Here is the call graph for this function:

void G4PolarizedCompton::PrintInfo ( )
overridevirtual

Implements G4VEmProcess.

Definition at line 137 of file G4PolarizedCompton.cc.

138 {
139  G4cout << " Total cross sections has a good parametrisation"
140  << " from 10 KeV to (100/Z) GeV"
141  << "\n Sampling according " << EmModel(1)->GetName() << " model"
142  << G4endl;
143 }
G4VEmModel * EmModel(G4int index=1) const
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
const G4String & GetName() const
Definition: G4VEmModel.hh:802

Here is the call graph for this function:

void G4PolarizedCompton::SetModel ( const G4String name)

Definition at line 147 of file G4PolarizedCompton.cc.

148 {
149  if(ss == "Klein-Nishina") { mType = 0; }
150  if(ss == "Polarized-Compton") { mType = 10; }
151 }

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