Geant4  10.01.p03
G4VBiasingOperator.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: $
28 //
29 // --------------------------------------------------------------------
30 // GEANT 4 class header file
31 //
32 // Class Description:
33 //
34 // An abstract class to pilot the biasing in a logical volume. This
35 // class is for *making decisions* on biasing operations to be applied.
36 // These ones are represented by the G4VBiasingOperation class.
37 // The volume in which biasing is applied is specified by the
38 // AttachTo(const G4LogicalVolume *) method. This has to be specified
39 // at detector construction time in the method ConstructSDandField() of
40 // G4VUsedDetectorConstruction.
41 //
42 // At tracking time the biasing operator is messaged by each
43 // G4BiasingProcessInterface object attached to the current track. For
44 // example, if three physics processes are under biasing, and if an
45 // additional G4BiasingProcessInterface is present to handle non-physics
46 // based biasing (splitting, killing), the operator will be messaged by
47 // these four G4BiasingProcessInterface objects.
48 // The idendity of the calling G4BiasingProcessInterface is known
49 // to the G4VBiasingOperator by passing this process pointer to the
50 // operator.
51 //
52 // ** Mandatory methods: **
53 //
54 // Three types of biasing are to be decided by the G4VBiasingOperator:
55 //
56 // 1) non-physics-based biasing:
57 // -----------------------------
58 // Meant for pure killing/splitting/etc. biasing operations, not
59 // associated to a physics process:
60 //
61 // virtual G4VBiasingOperation* ProposeNonPhysicsBiasingOperation( const G4Track* track,
62 // const G4BiasingProcessInterface* callingProcess ) = 0;
63 //
64 // Arguments are the current track, and the G4BiasingProcessInterface
65 // pointer making the call to the operator. In this case, this process
66 // does not wrap a physics process and
67 // callingProcess->GetWrappedProcess() == 0.
68 //
69 // The G4VBiasingOperation pointer returned is the operation to be
70 // applied. Zero can be returned. This operation will limit the
71 // step and propose a final state.
72 //
73 // This method is the first operator method called, it is called at the
74 // by the PostStepGetPhysicalInterationLenght(...) method of the
75 // G4BiasingProcessInterface.
76 //
77 // 2) physics-based biasing:
78 // -------------------------
79 // Physics-based biasing operations are of two types:
80 // - biasing of the physics process occurence interaction law
81 // - biasing of the physics process final state production
82 //
83 // a) The biasing of the occurence interaction law is proposed by:
84 //
85 // virtual G4VBiasingOperation* ProposeOccurenceBiasingOperation( const G4Track* track,
86 // const G4BiasingProcessInterface* callingProcess ) = 0;
87 // The current G4Track pointer and the G4BiasingProcessInterface
88 // pointer of the process calling the operator are passed. The
89 // G4BiasingProcessInterface process wraps an actual physics process
90 // which pointer can be obtained with
91 // callingProcess->GetWrappedProcess() .
92 //
93 // The biasing operation returned will be asked for its biasing
94 // interaction by the calling process, which will be a const object
95 // for the process. All setup and sampling regarding this law should be done
96 // in the operator before returning the related operation to the process.
97 //
98 // This method is the second operator one called in a step, it is called by
99 // the PostStepGetPhysicalInterationLenght(...) method of the
100 // G4BiasingProcessInterface.
101 //
102 // b) The biasing of the physics process final state is proposed by:
103 //
104 // virtual G4VBiasingOperation* ProposeFinalStateBiasingOperation( const G4Track* track,
105 // const G4BiasingProcessInterface* callingProcess ) = 0;
106 //
107 // The operator can propose a biasing operation that will handle the
108 // physic procsse final state biasing. As in previous case a) the
109 // G4BiasingProcessInterface process wraps an actual physics process
110 // which pointer can be obtained with:
111 // callingProcess->GetWrappedProcess() .
112 //
113 // Cases a) and b) are handled independently, and one or two of these
114 // biasing types can be provided in the same step.
115 //
116 // This method is the last operator one called in a step, it is called
117 // by the PostStepDoIt(...) method of the G4BiasingProcessInterface.
118 //
119 //
120 // ** Optional methods: **
121 //
122 // At the end of the step, the operator is messaged by the G4BiasingProcessInterface
123 // for operation(s) which have been applied during the step. One of the two following
124 // methods is called:
125 //
126 // - In case of at most a single biasing operation was applied by the process, report in cases of:
127 // - a non-physics biasing operation applied, biasingCase == BAC_NonPhysics ;
128 // - physics-based biasing:
129 // - the operator requested no occurence, nor final state biasing, and did let the
130 // physics process go : biasingCase == BAC_None;
131 // - an occurence biasing, where operation denied the application of the PostStepDoIt(..)
132 // of the physics process (proposing a weight for this),
133 // biasingCase == BAC_DenyInteraction ;
134 // -
135 //
136 // virtual void OperationApplied( const G4BiasingProcessInterface* callingProcess,
137 // G4BiasingAppliedCase biasingCase,
138 // G4VBiasingOperation* operationApplied,
139 // const G4VParticleChange* particleChangeProduced );
140 // At most a single biasing operation was applied by the process:
141 // - a non-physics biasing operation was applied, biasingCase == BAC_NonPhysics ;
142 // - physics-based biasing:
143 // - the operator requested no biasing operations, and did let the physics
144 // process go : biasingCase == BAC_None;
145 // - an occurence biasing was proposed, which operation purpose was to deny the
146 // application of the PostStepDoIt(..) of the physics process (proposing a
147 // weight for this) : biasingCase == BAC_DenyInteraction ;
148 // - a single final state biasing was proposed, with no concomittant occurence:
149 // biasingCase == BAC_FinalState;
150 // The operation applied and final state passed to the tracking (particleChangeProduced) are
151 // passed as information to the operator.
152 //
153 // virtual void OperationApplied( const G4BiasingProcessInterface* callingProcess,
154 // G4BiasingAppliedCase biasingCase,
155 // G4VBiasingOperation* occurenceOperationApplied,
156 // G4double weightForOccurenceInteraction,
157 // G4VBiasingOperation* finalStateOperationApplied,
158 // const G4VParticleChange* particleChangeProduced );
159 // This method is called in case an occurence biasing operation has been applied during the step.
160 // Depending on if the occurence operation was applied alone and together with a final state
161 // operation, the biasingCase will take values:
162 // - occurence biasing alone : biasingCase == BAC_None ;
163 // in which case finalStateOperationApplied == 0;
164 // - occurence biasing + final state biasing : biasingCase == BAC_FinalState;
165 // The particleChangeProduced is the one *before* application of the weight for occurence : hence
166 // either the particle change of the physics process, or the physics process one, biased the final
167 // state biasing operation.
168 //
169 //
170 // ----------------G4VBiasingOperation ----------------
171 //
172 // Author: M.Verderi (LLR), November 2013
173 //
174 // --------------------------------------------------------------------
175 
176 #ifndef G4VBiasingOperator_hh
177 #define G4VBiasingOperator_hh 1
178 
179 #include "globals.hh"
180 
181 class G4VBiasingOperation;
182 class G4Track;
184 class G4LogicalVolume;
185 class G4VParticleChange;
187 #include <map>
188 #include <vector>
189 #include "G4BiasingAppliedCase.hh"
190 #include "G4Cache.hh"
191 
192 
194 
195  // -- State machine used to inform operators
196  // -- about run starting.
197  // -- Defined at the end of this file.
199 
200 public:
201  // ---------------
202  // -- Constructor:
203  // ---------------
205  virtual ~G4VBiasingOperator();
206 
207  // ----------------------------------------------
208  // -- abstract and user interface to sub-classes:
209  // ----------------------------------------------
210 protected:
211  // -- mandatory methods to let the operator tell about biasing operations to be applied:
212  // -------------------------------------------------------------------------------------
213  // -- These three methods have the same arguments passed : the current G4Track pointer, and the pointer of the
214  // -- G4BiasingProcessInterface instance calling this biasing operator. This same biasing operator will be called by each
215  // -- of the G4BiasingProcessInterface instances, meaning for example that:
216  // -- - if one G4BiasingProcessInterface with no wrapped physics process exits, ProposeNonPhysicsBiasingOperation(...)
217  // -- will be called one time at the beginning of the step,
218  // -- - if three G4BiasingProcessInterface instances exist, each of these one wrapping a physics process (eg
219  // -- conversion, Compton, photo-electric), ProposeOccurenceBiasingOperation(...) will be called three times,
220  // -- by each of these instances, at the beginning of the step and ProposeFinalStateBiasingOperation(...) will
221  // -- also be called by each of these instances, at the PostStepDoIt level.
222  // -- If a null pointer is returned, the analog -unbiased- behavior is adopted.
223  // -- non-physics-based biasing:
224  // -----------------------------
225  // -- [ First operator method called, at the PostStepGetPhysicalInterationLenght(...) level. ]
226  virtual G4VBiasingOperation* ProposeNonPhysicsBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
227  // -- physics-based biasing:
228  // -------------------------
229  // -- Method to propose an occurence biasing operation : ie a change of the interaction length distribution. The proposed
230  // -- biasing operation will then be asked for its interaction law.
231  // -- Note that *** all sanity checks regarding the operation and its interaction law will have to have been performed
232  // -- before returning the biasing operation pointer *** as no corrective/aborting actions will be possible beyond this point.
233  // -- The informations provided by the G4BiasingProcessInterface calling process (previous occurence operation, previous step length,
234  // -- etc.) might be useful for doing this. They will be useful also to decide with continuing with a same operation proposed
235  // -- in the previous step, updating the interaction law taking into account the new G4Track state and the previous step size.
236  // -- [ Second operator method called, at the PostStepGetPhysicalInterationLenght(...) level. ]
237  virtual G4VBiasingOperation* ProposeOccurenceBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
238  // -- [ Third operator method called, at the PostStepDoIt(...) level. ]
239  virtual G4VBiasingOperation* ProposeFinalStateBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
240 
241 protected:
242  // -- optionnal methods for further information passed to the operator:
243  // --------------------------------------------------------------------
244  // ---- report to operator about the operation applied, the biasingCase value provides what sort of biasing applied:
245  virtual void OperationApplied( const G4BiasingProcessInterface* callingProcess, G4BiasingAppliedCase biasingCase,
246  G4VBiasingOperation* operationApplied, const G4VParticleChange* particleChangeProduced );
247  // ---- same as above, report about the operation applied, for the case an occurence biasing was applied, together or not with a final state biasing.
248  // ---- The variable biasingCase tells if the final state is a biased one or not. **But in all cases**, this call happens only
249  // ---- for an occurence biaising : ie the occurence weight is applied on top of the particleChangeProduced, which is the particle
250  // ---- *before* the weight application for occurence biasing.
251  virtual void OperationApplied( const G4BiasingProcessInterface* callingProcess, G4BiasingAppliedCase biasingCase,
252  G4VBiasingOperation* occurenceOperationApplied, G4double weightForOccurenceInteraction,
253  G4VBiasingOperation* finalStateOperationApplied, const G4VParticleChange* particleChangeProduced );
254 protected:
255  // ---- method to inform operator that its biasing control is over (exit volume, or end of tracking):
256  // ---- [Called at the beginning of next step, or at the end of tracking.]
257  virtual void ExitBiasing( const G4Track* track, const G4BiasingProcessInterface* callingProcess );
258 
259 
260 protected:
261  // -- Utility:
262  // -----------
263  // ---- A utility method allowing to store secondaries produced by an operation.
264  // ---- The stored secondaries are the ones in the particle change
265  void RememberSecondaries( const G4BiasingProcessInterface* callingProcess,
266  const G4VBiasingOperation* operationApplied,
267  const G4VParticleChange* particleChangeProduced );
268  // ---- Informations about track is erased:
269  void ForgetTrack( const G4Track* track );
270 
271 public:
272  // ---- inform the operator of the start of the run:
273  virtual void StartRun() {}
274  // ---- inform the operator of the start (end) of the tracking of a new track:
275  virtual void StartTracking( const G4Track* /* track */ ) {}
276  virtual void EndTracking() {}
277 
278 
279 
280  // --------------------
281  // -- public interface:
282  // --------------------
283  // -- needed by user:
284 public:
285  const G4String GetName() const {return fName;}
286  void AttachTo( const G4LogicalVolume* ); // -- attach to single volume
287 
289  // -- all operators (might got to a manager):
290  static const std::vector < G4VBiasingOperator* >& GetBiasingOperators() {return fOperators.Get();}
291  // -- get operator associated to a logical volume:
292  static G4VBiasingOperator* GetBiasingOperator( const G4LogicalVolume* ); // -- might go to a manager ; or moved to volume
293 
294 
295 
296  // -- used by biasing process interface, or used by an other operator (not expected to be invoked differently than with these two cases):
297 public:
301  void ExitingBiasing( const G4Track* track, const G4BiasingProcessInterface* callingProcess );
302 
303 public:
304  void ReportOperationApplied( const G4BiasingProcessInterface* callingProcess, G4BiasingAppliedCase biasingCase,
305  G4VBiasingOperation* operationApplied, const G4VParticleChange* particleChangeProduced );
306  void ReportOperationApplied( const G4BiasingProcessInterface* callingProcess, G4BiasingAppliedCase biasingCase,
307  G4VBiasingOperation* occurenceOperationApplied, G4double weightForOccurenceInteraction,
308  G4VBiasingOperation* finalStateOperationApplied, const G4VParticleChange* particleChangeProduced );
309 
310 
311 public:
314 
315 
316 private:
318  // -- thread local:
319  // static std::map< const G4LogicalVolume*, G4VBiasingOperator* > fLogicalToSetupMap;
321  // -- thread local:
323  // static std::vector < G4VBiasingOperator* > fOperators;
324 
325  // -- thread local:
327 
328 
329  // -- For this operator:
330  std::vector< const G4LogicalVolume* > fRootVolumes;
331  std::map < const G4LogicalVolume*, G4int > fDepthInTree;
332 
333  // -- current operation:
337 
338  // -- previous operations:
346 
347 };
348 
349 // -- state machine to get biasing operators
350 // -- messaged at the beginning of runs:
351 #include "G4VStateDependent.hh"
353 public:
356 public:
357  G4bool Notify(G4ApplicationState requestedState);
358 private:
360 };
361 
362 #endif
const G4VBiasingOperation * fPreviousProposedOccurenceBiasingOperation
G4VBiasingOperation * GetProposedFinalStateBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
virtual void OperationApplied(const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
std::map< const G4LogicalVolume *, G4int > fDepthInTree
G4BiasingAppliedCase GetPreviousBiasingAppliedCase() const
G4String name
Definition: TRTMaterials.hh:40
value_type & Get() const
Definition: G4Cache.hh:282
std::vector< const G4LogicalVolume * > fRootVolumes
G4BiasingAppliedCase fPreviousBiasingAppliedCase
G4VBiasingOperation * GetProposedOccurenceBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
virtual G4VBiasingOperation * ProposeFinalStateBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
G4bool Notify(G4ApplicationState requestedState)
virtual void ExitBiasing(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
G4BiasingAppliedCase
const G4VBiasingOperation * fPreviousAppliedOccurenceBiasingOperation
static G4VectorCache< G4VBiasingOperator * > fOperators
static G4MapCache< const G4LogicalVolume *, G4VBiasingOperator * > fLogicalToSetupMap
bool G4bool
Definition: G4Types.hh:79
virtual G4VBiasingOperation * ProposeOccurenceBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
const G4VBiasingOperation * fPreviousAppliedNonPhysicsBiasingOperation
const G4VBiasingOperation * fPreviousProposedFinalStateBiasingOperation
void AttachTo(const G4LogicalVolume *)
const G4VBiasingOperation * fPreviousAppliedFinalStateBiasingOperation
virtual G4VBiasingOperation * ProposeNonPhysicsBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
static G4VBiasingOperator * GetBiasingOperator(const G4LogicalVolume *)
void ForgetTrack(const G4Track *track)
G4VBiasingOperator(G4String name)
void RememberSecondaries(const G4BiasingProcessInterface *callingProcess, const G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
virtual void StartTracking(const G4Track *)
G4VBiasingOperation * fFinalStateBiasingOperation
static const std::vector< G4VBiasingOperator * > & GetBiasingOperators()
static G4Cache< G4BiasingOperatorStateNotifier * > fStateNotifier
G4VBiasingOperation * fOccurenceBiasingOperation
double G4double
Definition: G4Types.hh:76
void ReportOperationApplied(const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
const G4VBiasingOperation * fPreviousProposedNonPhysicsBiasingOperation
G4VBiasingOperation * GetProposedNonPhysicsBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
void ExitingBiasing(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
const G4VBiasingOperation * GetBirthOperation(const G4Track *)
G4VBiasingOperation * fNonPhysicsBiasingOperation
const G4VBiasingOperation * GetPreviousNonPhysicsAppliedOperation()
virtual void EndTracking()
G4ApplicationState
const G4String GetName() const