Geant4  10.00.p03
G4FieldManager.hh
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4FieldManager.hh 68055 2013-03-13 14:43:28Z gcosmo $
28 //
29 //
30 // class G4FieldManager
31 //
32 // Class description:
33 //
34 // A class to manage (Store) a pointer to the Field subclass that
35 // describes the field of a detector (magnetic, electric or other).
36 // Also stores a reference to the chord finder.
37 //
38 // The G4FieldManager class exists to allow the user program to specify
39 // the electric, magnetic and/or other field(s) of the detector.
40 //
41 // A field manager can be set to a logical volume (or to more than one),
42 // in order to vary its field from that of the world. In this manner
43 // a zero or constant field can override a global field, a more or
44 // less exact version can override the external approximation, lower
45 // or higher precision for tracking can be specified, a different
46 // stepper can be chosen for different volumes, ...
47 //
48 // It also stores a pointer to the ChordFinder object that can do the
49 // propagation in this field. All geometrical track "advancement"
50 // in the field is handled by this ChordFinder object.
51 //
52 // G4FieldManager allows the other classes/object (of the MagneticField
53 // & other class categories) to find out whether a detector field object
54 // exists and what that object is.
55 //
56 // The Chord Finder must be created either by calling CreateChordFinder
57 // for a Magnetic Field or by the user creating a a Chord Finder object
58 // "manually" and setting this pointer.
59 //
60 // A default FieldManager is created by the singleton class
61 // G4NavigatorForTracking and exists before main is called.
62 // However a new one can be created and given to G4NavigatorForTracking.
63 //
64 // Our current design envisions that one Field manager is
65 // valid for each region detector.
66 
67 // History:
68 // - 05.11.03 John Apostolakis, Added Min/MaximumEpsilonStep
69 // - 20.06.03 John Apostolakis, Abstract & ability to ConfigureForTrack
70 // - 10.03.97 John Apostolakis, design and implementation.
71 // -------------------------------------------------------------------
72 
73 #ifndef G4FIELDMANAGER_HH
74 #define G4FIELDMANAGER_HH 1
75 
76 #include "globals.hh"
77 
78 class G4Field;
79 class G4MagneticField;
80 class G4ChordFinder;
81 class G4Track; // Forward reference for parameter configuration
82 
84 {
85  public: // with description
86  G4FieldManager(G4Field *detectorField=0,
87  G4ChordFinder *pChordFinder=0,
88  G4bool b=true ); // fieldChangesEnergy is taken from field
89  // General constructor for any field.
90  // -> Must be set with field and chordfinder for use.
91  G4FieldManager(G4MagneticField *detectorMagneticField);
92  // Creates ChordFinder
93  // - assumes pure magnetic field (so Energy constant)
94  virtual ~G4FieldManager();
95 
96  G4bool SetDetectorField(G4Field *detectorField);
97  inline const G4Field* GetDetectorField() const;
98  inline G4bool DoesFieldExist() const;
99  // Set, get and check the field object
100 
101  void CreateChordFinder(G4MagneticField *detectorMagField);
102  inline void SetChordFinder(G4ChordFinder *aChordFinder);
103  inline G4ChordFinder* GetChordFinder();
104  inline const G4ChordFinder* GetChordFinder() const;
105  // Create, set or get the associated Chord Finder
106 
107  virtual void ConfigureForTrack( const G4Track * );
108  // Setup the choice of the configurable parameters
109  // relying on the current track's energy, particle identity, ..
110  // Note: In addition to the values of member variables,
111  // a user can use this to change the ChordFinder, the field, ...
112 
113  public: // with description
114 
115  inline G4double GetDeltaIntersection() const; // virtual ?
116  // Accuracy for boundary intersection.
117 
118  inline G4double GetDeltaOneStep() const; // virtual ?
119  // Accuracy for one tracking/physics step.
120 
121  inline void SetAccuraciesWithDeltaOneStep(G4double valDeltaOneStep);
122  // Sets both accuracies, maintaining a fixed ratio for accuracties
123  // of volume Intersection and Integration (in One Step)
124 
125  inline void SetDeltaOneStep(G4double valueD1step);
126  // Set accuracy for integration of one step. (only)
127  inline void SetDeltaIntersection(G4double valueDintersection);
128  // Set accuracy of intersection of a volume. (only)
129 
130  inline G4double GetMinimumEpsilonStep() const;
131  inline void SetMinimumEpsilonStep( G4double newEpsMin );
132  // Minimum for Relative accuracy of a Step
133 
134  inline G4double GetMaximumEpsilonStep() const;
135  inline void SetMaximumEpsilonStep( G4double newEpsMax );
136  // Maximum for Relative accuracy of a Step
137 
138  inline G4bool DoesFieldChangeEnergy() const;
139  inline void SetFieldChangesEnergy(G4bool value);
140  // For electric field this should be true
141  // For magnetic field this should be false
142 
143  virtual G4FieldManager* Clone() const;
144  //Needed for multi-threading, create a clone of this object
145 
146  private:
147 
150  // Private copy constructor and assignment operator.
151 
152  private:
153  // Dependent objects -- with state that depends on tracking
156 
157  G4bool fAllocatedChordFinder; // Did we used "new" to
158  // create fChordFinder ?
159  // INVARIANTS of tracking ---------------------------------------
160  //
161  // 1. CONSTANTS
162  const G4double fEpsilonMinDefault; // Can be 1.0e-5 to 1.0e-10 ...
163  const G4double fEpsilonMaxDefault; // Can be 1.0e-3 to 1.0e-8 ...
164 
165  // 2. CHARACTERISTIC of field
167 
168  // 3. PARAMETERS
169  //
170  // Values for the required accuracies
171  G4double fDelta_One_Step_Value; // for one tracking/physics step
172  G4double fDelta_Intersection_Val; // for boundary intersection
173 
176 
177  // Values for the small possible relative accuracy of a step
178  // (corresponding to the greatest possible integration accuracy)
181 
182 };
183 
184 // Our current design and implementation expect that a particular
185 // geometrical region has a Field manager.
186 // By default a Field Manager is created for the world volume, and
187 // will be utilised for all volumes unless it is overridden by a 'local'
188 // field manager.
189 
190 // Note also that a region with both electric E and magnetic B field will
191 // have these treated as one field.
192 // Similarly it could be extended to treat other fields as additional components
193 // of a single field type.
194 
195 
196 // Implementation of inline functions
197 
198 #include "G4FieldManager.icc"
199 
200 #endif /* G4FIELDMANAGER_HH */
virtual G4FieldManager * Clone() const
G4ChordFinder * GetChordFinder()
G4bool SetDetectorField(G4Field *detectorField)
G4double fDelta_One_Step_Value
G4double GetDeltaOneStep() const
G4double fEpsilonMin
G4double fDefault_Delta_Intersection_Val
void SetChordFinder(G4ChordFinder *aChordFinder)
G4double fDefault_Delta_One_Step_Value
void SetAccuraciesWithDeltaOneStep(G4double valDeltaOneStep)
virtual void ConfigureForTrack(const G4Track *)
void SetMinimumEpsilonStep(G4double newEpsMin)
const G4double fEpsilonMaxDefault
G4FieldManager(G4Field *detectorField=0, G4ChordFinder *pChordFinder=0, G4bool b=true)
bool G4bool
Definition: G4Types.hh:79
G4double GetMaximumEpsilonStep() const
G4bool DoesFieldExist() const
G4bool DoesFieldChangeEnergy() const
G4double fEpsilonMax
G4double GetDeltaIntersection() const
G4ChordFinder * fChordFinder
virtual ~G4FieldManager()
G4bool fFieldChangesEnergy
void SetDeltaIntersection(G4double valueDintersection)
void SetFieldChangesEnergy(G4bool value)
G4Field * fDetectorField
void SetMaximumEpsilonStep(G4double newEpsMax)
double G4double
Definition: G4Types.hh:76
void CreateChordFinder(G4MagneticField *detectorMagField)
const G4Field * GetDetectorField() const
G4double fDelta_Intersection_Val
G4bool fAllocatedChordFinder
void SetDeltaOneStep(G4double valueD1step)
const G4double fEpsilonMinDefault
G4FieldManager & operator=(const G4FieldManager &)
G4double GetMinimumEpsilonStep() const