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