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

#include <G4FieldManager.hh>

Public Member Functions

 G4FieldManager (G4Field *detectorField=nullptr, G4ChordFinder *pChordFinder=nullptr, G4bool b=true)
 
 G4FieldManager (G4MagneticField *detectorMagneticField)
 
virtual ~G4FieldManager ()
 
 G4FieldManager (const G4FieldManager &)=delete
 
G4FieldManageroperator= (const G4FieldManager &)=delete
 
G4bool SetDetectorField (G4Field *detectorField, G4int failMode=0)
 
void ProposeDetectorField (G4Field *detectorField)
 
void ChangeDetectorField (G4Field *detectorField)
 
const G4FieldGetDetectorField () const
 
G4bool DoesFieldExist () const
 
void CreateChordFinder (G4MagneticField *detectorMagField)
 
void SetChordFinder (G4ChordFinder *aChordFinder)
 
G4ChordFinderGetChordFinder ()
 
const G4ChordFinderGetChordFinder () const
 
virtual void ConfigureForTrack (const G4Track *)
 
G4double GetDeltaIntersection () const
 
G4double GetDeltaOneStep () const
 
void SetAccuraciesWithDeltaOneStep (G4double valDeltaOneStep)
 
void SetDeltaOneStep (G4double valueD1step)
 
void SetDeltaIntersection (G4double valueDintersection)
 
G4double GetMinimumEpsilonStep () const
 
void SetMinimumEpsilonStep (G4double newEpsMin)
 
G4double GetMaximumEpsilonStep () const
 
void SetMaximumEpsilonStep (G4double newEpsMax)
 
G4bool DoesFieldChangeEnergy () const
 
void SetFieldChangesEnergy (G4bool value)
 
virtual G4FieldManagerClone () const
 

Detailed Description

Definition at line 84 of file G4FieldManager.hh.

Constructor & Destructor Documentation

◆ G4FieldManager() [1/3]

G4FieldManager::G4FieldManager ( G4Field detectorField = nullptr,
G4ChordFinder pChordFinder = nullptr,
G4bool  b = true 
)

Definition at line 41 of file G4FieldManager.cc.

44 : fDetectorField(detectorField),
45 fChordFinder(pChordFinder),
46 fDelta_One_Step_Value( fDefault_Delta_One_Step_Value ),
47 fDelta_Intersection_Val( fDefault_Delta_Intersection_Val ),
48 fEpsilonMin( fEpsilonMinDefault ),
49 fEpsilonMax( fEpsilonMaxDefault)
50{
51 if ( detectorField != nullptr )
52 {
53 fFieldChangesEnergy = detectorField->DoesFieldChangeEnergy();
54 }
55 else
56 {
57 fFieldChangesEnergy = fieldChangesEnergy;
58 }
59
60 // Add to store
61 //
63}
static void Register(G4FieldManager *pVolume)
virtual G4bool DoesFieldChangeEnergy() const =0

◆ G4FieldManager() [2/3]

G4FieldManager::G4FieldManager ( G4MagneticField detectorMagneticField)

Definition at line 65 of file G4FieldManager.cc.

66 : fDetectorField(detectorField), fAllocatedChordFinder(true),
67 fDelta_One_Step_Value( fDefault_Delta_One_Step_Value ),
68 fDelta_Intersection_Val( fDefault_Delta_Intersection_Val ),
69 fEpsilonMin( fEpsilonMinDefault ),
70 fEpsilonMax( fEpsilonMaxDefault )
71{
72 fChordFinder = new G4ChordFinder( detectorField );
73
74 // Add to store
75 //
77}

◆ ~G4FieldManager()

G4FieldManager::~G4FieldManager ( )
virtual

Definition at line 139 of file G4FieldManager.cc.

140{
141 if( fAllocatedChordFinder )
142 {
143 delete fChordFinder;
144 }
146}
static void DeRegister(G4FieldManager *pVolume)

◆ G4FieldManager() [3/3]

G4FieldManager::G4FieldManager ( const G4FieldManager )
delete

Member Function Documentation

◆ ChangeDetectorField()

void G4FieldManager::ChangeDetectorField ( G4Field detectorField)
inline

◆ Clone()

G4FieldManager * G4FieldManager::Clone ( ) const
virtual

Definition at line 79 of file G4FieldManager.cc.

80{
81 G4Field* aField = nullptr;
82 G4FieldManager* aFM = nullptr;
83 G4ChordFinder* aCF = nullptr;
84 try {
85 if ( fDetectorField != nullptr )
86 {
87 aField = fDetectorField->Clone();
88 }
89
90 // Create a new field manager, note that we do not set
91 // any chordfinder now
92 //
93 aFM = new G4FieldManager( aField , nullptr , fFieldChangesEnergy );
94
95 // Check if originally we have the fAllocatedChordFinder variable
96 // set, in case, call chord constructor
97 //
98 if ( fAllocatedChordFinder )
99 {
100 aFM->CreateChordFinder( dynamic_cast<G4MagneticField*>(aField) );
101 }
102 else
103 {
104 // Chord was specified by user, should we clone?
105 // TODO: For the moment copy pointer, to be understood
106 // if cloning of ChordFinder is needed
107 //
108 aCF = fChordFinder; /*->Clone*/
109 aFM->fChordFinder = aCF;
110 }
111
112 // Copy values of other variables
113
114 aFM->fEpsilonMax = fEpsilonMax;
115 aFM->fEpsilonMin = fEpsilonMin;
116 aFM->fDelta_Intersection_Val = fDelta_Intersection_Val;
117 aFM->fDelta_One_Step_Value = fDelta_One_Step_Value;
118 // TODO: Should we really add to the store the cloned FM?
119 // Who will use this?
120 }
121 catch ( ... )
122 {
123 // Failed creating clone: probably user did not implement Clone method
124 // in derived classes?
125 // Perform clean-up after ourselves...
126 delete aField;
127 delete aFM;
128 delete aCF;
129 throw;
130 }
131 return aFM;
132}
void CreateChordFinder(G4MagneticField *detectorMagField)
virtual G4Field * Clone() const
Definition: G4Field.cc:54

Referenced by G4VUserDetectorConstruction::CloneF().

◆ ConfigureForTrack()

void G4FieldManager::ConfigureForTrack ( const G4Track )
virtual

Definition at line 134 of file G4FieldManager.cc.

135{
136 // Default is to do nothing!
137}

Referenced by G4CoupledTransportation::AlongStepGetPhysicalInteractionLength(), and G4ITTransportation::AlongStepGetPhysicalInteractionLength().

◆ CreateChordFinder()

void G4FieldManager::CreateChordFinder ( G4MagneticField detectorMagField)

Definition at line 149 of file G4FieldManager.cc.

150{
151 if ( fAllocatedChordFinder )
152 {
153 delete fChordFinder;
154 }
155 fAllocatedChordFinder = false;
156
157 if( detectorMagField != nullptr )
158 {
159 fChordFinder = new G4ChordFinder( detectorMagField );
160 fAllocatedChordFinder = true;
161 }
162 else
163 {
164 fChordFinder = nullptr;
165 }
166}

Referenced by Clone().

◆ DoesFieldChangeEnergy()

◆ DoesFieldExist()

G4bool G4FieldManager::DoesFieldExist ( ) const
inline

◆ GetChordFinder() [1/2]

G4ChordFinder * G4FieldManager::GetChordFinder ( )
inline

◆ GetChordFinder() [2/2]

const G4ChordFinder * G4FieldManager::GetChordFinder ( ) const
inline

◆ GetDeltaIntersection()

G4double G4FieldManager::GetDeltaIntersection ( ) const
inline

◆ GetDeltaOneStep()

G4double G4FieldManager::GetDeltaOneStep ( ) const
inline

◆ GetDetectorField()

◆ GetMaximumEpsilonStep()

G4double G4FieldManager::GetMaximumEpsilonStep ( ) const
inline

◆ GetMinimumEpsilonStep()

G4double G4FieldManager::GetMinimumEpsilonStep ( ) const
inline

◆ operator=()

G4FieldManager & G4FieldManager::operator= ( const G4FieldManager )
delete

◆ ProposeDetectorField()

void G4FieldManager::ProposeDetectorField ( G4Field detectorField)
inline

◆ SetAccuraciesWithDeltaOneStep()

void G4FieldManager::SetAccuraciesWithDeltaOneStep ( G4double  valDeltaOneStep)
inline

◆ SetChordFinder()

void G4FieldManager::SetChordFinder ( G4ChordFinder aChordFinder)
inline

◆ SetDeltaIntersection()

void G4FieldManager::SetDeltaIntersection ( G4double  valueDintersection)
inline

◆ SetDeltaOneStep()

void G4FieldManager::SetDeltaOneStep ( G4double  valueD1step)
inline

◆ SetDetectorField()

G4bool G4FieldManager::SetDetectorField ( G4Field detectorField,
G4int  failMode = 0 
)

Definition at line 180 of file G4FieldManager.cc.

182{
183 G4VIntegrationDriver* driver = nullptr;
184 G4EquationOfMotion* equation = nullptr;
185 // G4bool compatibleField = false;
186 G4bool ableToSet = false;
187
188 fDetectorField = pDetectorField;
189 InitialiseFieldChangesEnergy();
190
191 // Must 'propagate' the field to the dependent classes
192 //
193 if( fChordFinder != nullptr )
194 {
195 failMode= std::max( failMode, 1) ;
196 // If a chord finder exists, warn in case of error!
197
198 driver = fChordFinder->GetIntegrationDriver();
199 if( driver != nullptr )
200 {
201 equation = driver->GetEquationOfMotion();
202
203 // Should check the compatibility between the
204 // field and the equation HERE
205
206 if( equation != nullptr )
207 {
208 equation->SetFieldObj(pDetectorField);
209 ableToSet = true;
210 }
211 }
212 }
213
214 if( !ableToSet && (failMode > 0) )
215 {
216 // If this fails, report the issue !
217
219 msg << "Unable to set the field in the dependent objects of G4FieldManager"
220 << G4endl;
221 msg << "All the dependent classes must be fully initialised,"
222 << "before it is possible to call this method." << G4endl;
223 msg << "The problem encountered was the following: " << G4endl;
224 if( fChordFinder == nullptr ) { msg << " No ChordFinder. " ; }
225 else if( driver == nullptr ) { msg << " No Integration Driver set. ";}
226 else if( equation == nullptr ) { msg << " No Equation found. " ; }
227 // else if( !compatibleField ) { msg << " Field not compatible. ";}
228 else { msg << " Can NOT find reason for failure. ";}
229 msg << G4endl;
230 G4ExceptionSeverity severity = (failMode != 1)
232 G4Exception("G4FieldManager::SetDetectorField", "Geometry001",
233 severity, msg);
234 }
235 return ableToSet;
236}
G4ExceptionSeverity
@ JustWarning
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
bool G4bool
Definition: G4Types.hh:86
#define G4endl
Definition: G4ios.hh:57
G4VIntegrationDriver * GetIntegrationDriver()
void SetFieldObj(G4Field *pField)
virtual G4EquationOfMotion * GetEquationOfMotion()=0

◆ SetFieldChangesEnergy()

void G4FieldManager::SetFieldChangesEnergy ( G4bool  value)
inline

◆ SetMaximumEpsilonStep()

void G4FieldManager::SetMaximumEpsilonStep ( G4double  newEpsMax)
inline

◆ SetMinimumEpsilonStep()

void G4FieldManager::SetMinimumEpsilonStep ( G4double  newEpsMin)
inline

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