Geant4_10
G4VBiasingOperator.cc
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 #include "G4VBiasingOperator.hh"
27 #include "G4VBiasingOperation.hh"
28 #include "G4VParticleChange.hh"
29 #include "G4BiasingTrackData.hh"
31 
32 
33 G4MapCache< const G4LogicalVolume*, G4VBiasingOperator* > G4VBiasingOperator::fLogicalToSetupMap;
34 G4VectorCache< G4VBiasingOperator* > G4VBiasingOperator::fOperators;
35 
36 
38  : fName(name)
39 {
40  fOperators.Push_back(this);
41 }
42 
44 {
45 }
46 
48 {
50  it = fLogicalToSetupMap.Find(logical);
51  if ( it == fLogicalToSetupMap.End() ) fLogicalToSetupMap[logical] = this;
52  else if ( (*it).second != this )
53  {
55  ed << "Biasing operator `" << GetName()
56  << "' can not be attached to Logical volume `"
57  << logical->GetName() << "' which is already used by an other operator !" << G4endl;
58  G4Exception("G4VBiasingOperator::AttachTo(...)",
59  "BIAS.MNG.01",
61  ed);
62  }
63 }
64 
65 
67 {
69  it = fLogicalToSetupMap.Find(logical);
70  if ( it == fLogicalToSetupMap.End() ) return 0;
71  else return (*it).second;
72 }
73 
75 {
76  fOccurenceBiasingOperation = ProposeOccurenceBiasingOperation(track, callingProcess);
77  // G4cout << GetName() << " Occurence : ";
78  // if ( fOccurenceBiasingOperation ) G4cout << fOccurenceBiasingOperation->GetName() << G4endl;
79  // else G4cout << " (none) " << G4endl;
80  return fOccurenceBiasingOperation;
81 }
82 
84 {
85  fFinalStateBiasingOperation = ProposeFinalStateBiasingOperation(track, callingProcess);
86  // G4cout << GetName() << " Final State : ";
87  // if ( fFinalStateBiasingOperation ) G4cout << fFinalStateBiasingOperation->GetName() << G4endl;
88  // else G4cout << " (none) " << G4endl;
89  return fFinalStateBiasingOperation;
90 }
91 
93 {
94  fNonPhysicsBiasingOperation = ProposeNonPhysicsBiasingOperation(track, callingProcess);
95  // G4cout << GetName() << " Non Physics : ";
96  // if ( fNonPhysicsBiasingOperation ) G4cout << fNonPhysicsBiasingOperation->GetName() << G4endl;
97  // else G4cout << " (none) " << G4endl;
98  return fNonPhysicsBiasingOperation;
99 }
100 
101 // void G4VBiasingOperator::RegisterSecondariesBirth( G4VBiasingOperation* operation, const G4VParticleChange* pChange )
102 // {
103 // for (G4int i2nd = 0; i2nd < pChange->GetNumberOfSecondaries (); i2nd++)
104 // RegisterSecondaryBirth( operation, pChange->GetSecondary (i2nd) );
105 // }
106 
107 // void G4VBiasingOperator::RegisterSecondaryBirth( G4VBiasingOperation* operation, const G4Track* track )
108 // {
109 // // G4cout << " --> Registering operation, track = " << operation->GetName() << " " << track << G4endl;
110 // fSecondaryBirthOperation[track] = operation;
111 // }
112 
114 {
116  if ( biasingData != 0 ) return biasingData->GetBirthOperation();
117  else return 0;
118 }
119 
120 
122  G4BiasingAppliedCase biasingCase,
123  G4VBiasingOperation* operationApplied,
124  const G4VParticleChange* particleChangeProduced )
125 {
126  fPreviousBiasingAppliedCase = biasingCase;
127  // if ( operationApplied != 0 )
128  // {
129  // if ( operationApplied->GetRememberSecondaries() )
130  // {
131  // for (G4int i2nd = 0; i2nd < particleChangeProduced->GetNumberOfSecondaries (); i2nd++)
132  // new G4BiasingTrackData( particleChangeProduced->GetSecondary (i2nd),
133  // operationApplied,
134  // this,
135  // callingProcess);
136  // }
137  // }
138  fPreviousAppliedOccurenceBiasingOperation = 0;
139  fPreviousAppliedFinalStateBiasingOperation = 0;
140  fPreviousAppliedNonPhysicsBiasingOperation = 0;
141  switch ( biasingCase )
142  {
143  case BAC_None:
144  break;
145  case BAC_NonPhysics:
146  fPreviousAppliedNonPhysicsBiasingOperation = operationApplied ;
147  break;
148  case BAC_DenyInteraction:
149  fPreviousAppliedOccurenceBiasingOperation = operationApplied;
150  break;
151  case BAC_FinalState:
152  fPreviousAppliedFinalStateBiasingOperation = operationApplied;
153  break;
154  case BAC_Occurence:
155  G4Exception("G4VBiasingOperator::ReportOperationApplied(...)",
156  "BIAS.MNG.02",
157  JustWarning,
158  "Internal logic error, please report !");
159  break;
160  default:
161  G4Exception("G4VBiasingOperator::ReportOperationApplied(...)",
162  "BIAS.MNG.03",
163  JustWarning,
164  "Internal logic error, please report !");
165  }
166  OperationApplied( callingProcess, biasingCase, operationApplied, particleChangeProduced );
167 }
168 
170  G4BiasingAppliedCase biasingCase,
171  G4VBiasingOperation* occurenceOperationApplied,
172  G4double weightForOccurenceInteraction,
173  G4VBiasingOperation* finalStateOperationApplied,
174  const G4VParticleChange* particleChangeProduced )
175 {
176  fPreviousBiasingAppliedCase = biasingCase;
177  // G4VBiasingOperation* operation(finalStateOperationApplied != 0 ? finalStateOperationApplied : occurenceOperationApplied);
178  // G4bool remember = occurenceOperationApplied->GetRememberSecondaries();
179  // if ( finalStateOperationApplied != 0 ) remember = remember || finalStateOperationApplied->GetRememberSecondaries();
180  // if ( remember )
181  // {
182  // for (G4int i2nd = 0; i2nd < particleChangeProduced->GetNumberOfSecondaries (); i2nd++)
183  // new G4BiasingTrackData( particleChangeProduced->GetSecondary (i2nd),
184  // operation,
185  // this,
186  // callingProcess);
187  // }
188  fPreviousAppliedOccurenceBiasingOperation = occurenceOperationApplied;
189  fPreviousAppliedFinalStateBiasingOperation = finalStateOperationApplied;
190  OperationApplied( callingProcess, biasingCase, occurenceOperationApplied, weightForOccurenceInteraction, finalStateOperationApplied, particleChangeProduced );
191 }
192 
193 
194 void G4VBiasingOperator::ExitingBiasing( const G4Track* track, const G4BiasingProcessInterface* callingProcess )
195 {
196  ExitBiasing( track, callingProcess );
197 
198  // -- reset all data members:
199  fOccurenceBiasingOperation = 0 ;
200  fFinalStateBiasingOperation = 0 ;
201  fNonPhysicsBiasingOperation = 0 ;
202  fPreviousProposedOccurenceBiasingOperation = 0 ;
203  fPreviousProposedFinalStateBiasingOperation = 0 ;
204  fPreviousProposedNonPhysicsBiasingOperation = 0 ;
205  fPreviousAppliedOccurenceBiasingOperation = 0 ;
206  fPreviousAppliedFinalStateBiasingOperation = 0 ;
207  fPreviousAppliedNonPhysicsBiasingOperation = 0 ;
208  fPreviousBiasingAppliedCase = BAC_None ;
209 }
210 
211 
213  const G4VBiasingOperation* operationApplied,
214  const G4VParticleChange* particleChangeProduced)
215 {
216  for (G4int i2nd = 0; i2nd < particleChangeProduced->GetNumberOfSecondaries (); i2nd++)
217  new G4BiasingTrackData( particleChangeProduced->GetSecondary (i2nd),
218  operationApplied,
219  this,
220  callingProcess);
221 }
222 
224 {
226  if ( biasingData != 0 ) delete biasingData;
227 }
228 
229 
230 // -- dummy empty implementations to allow letting arguments visible in the .hh
231 // -- but avoiding annoying warning messages about unused variables
232 // -- methods to inform operator that its biasing control is over:
234 {}
237 {}
241 {}
G4VBiasingOperation * GetProposedFinalStateBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
virtual void OperationApplied(const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
G4int GetNumberOfSecondaries() const
G4Track * GetSecondary(G4int anIndex) const
static G4BiasingTrackDataStore * GetInstance()
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4String GetName() const
iterator End()
Definition: G4Cache.hh:382
G4String fName
Definition: G4AttUtils.hh:55
const G4VBiasingOperation * GetBirthOperation() const
const XML_Char * name
Definition: expat.h:151
G4VBiasingOperation * GetProposedOccurenceBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
virtual G4VBiasingOperation * ProposeFinalStateBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
int G4int
Definition: G4Types.hh:78
G4BiasingTrackData * GetBiasingTrackData(const G4Track *track)
virtual void ExitBiasing(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
G4BiasingAppliedCase
virtual G4VBiasingOperation * ProposeOccurenceBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
void AttachTo(const G4LogicalVolume *)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual G4VBiasingOperation * ProposeNonPhysicsBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)=0
static G4VBiasingOperator * GetBiasingOperator(const G4LogicalVolume *)
void ForgetTrack(const G4Track *track)
void Push_back(const value_type &val)
Definition: G4Cache.hh:304
iterator Find(const key_type &k)
Definition: G4Cache.hh:388
G4VBiasingOperator(G4String name)
void RememberSecondaries(const G4BiasingProcessInterface *callingProcess, const G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void ReportOperationApplied(const G4BiasingProcessInterface *callingProcess, G4BiasingAppliedCase biasingCase, G4VBiasingOperation *operationApplied, const G4VParticleChange *particleChangeProduced)
G4VBiasingOperation * GetProposedNonPhysicsBiasingOperation(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
void ExitingBiasing(const G4Track *track, const G4BiasingProcessInterface *callingProcess)
const G4VBiasingOperation * GetBirthOperation(const G4Track *)
const G4String GetName() const