00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048 #ifndef G4EmBiasingManager_h
00049 #define G4EmBiasingManager_h 1
00050
00051 #include "globals.hh"
00052 #include "G4ParticleDefinition.hh"
00053 #include "G4DynamicParticle.hh"
00054 #include "Randomize.hh"
00055 #include <vector>
00056
00057
00058
00059 class G4Region;
00060 class G4Track;
00061 class G4VEnergyLossProcess;
00062 class G4VEmModel;
00063 class G4MaterialCutsCouple;
00064 class G4ParticleChangeForLoss;
00065 class G4ParticleChangeForGamma;
00066
00067 class G4EmBiasingManager
00068 {
00069 public:
00070
00071 G4EmBiasingManager();
00072
00073 ~G4EmBiasingManager();
00074
00075 void Initialise(const G4ParticleDefinition& part,
00076 const G4String& procName, G4int verbose);
00077
00078
00079 void ActivateForcedInteraction(G4double length = 0.0,
00080 const G4String& r = "");
00081
00082
00083 void ActivateSecondaryBiasing(const G4String& region, G4double factor,
00084 G4double energyLimit);
00085
00086
00087 G4double GetStepLimit(G4int coupleIdx, G4double previousStep);
00088
00089
00090
00091
00092
00093
00094
00095 G4double ApplySecondaryBiasing(std::vector<G4DynamicParticle*>&,
00096 const G4Track& track,
00097 G4VEmModel* currentModel,
00098 G4ParticleChangeForGamma* pParticleChange,
00099 G4double& eloss,
00100 G4int coupleIdx,
00101 G4double tcut,
00102 G4double safety = 0.0);
00103
00104
00105 G4double ApplySecondaryBiasing(std::vector<G4DynamicParticle*>&,
00106 const G4Track& track,
00107 G4VEmModel* currentModel,
00108 G4ParticleChangeForLoss* pParticleChange,
00109 G4double& eloss,
00110 G4int coupleIdx,
00111 G4double tcut,
00112 G4double safety = 0.0);
00113
00114
00115 G4double ApplySecondaryBiasing(std::vector<G4Track*>&,
00116 G4int coupleIdx);
00117
00118 inline G4bool SecondaryBiasingRegion(G4int coupleIdx);
00119
00120 inline G4bool ForcedInteractionRegion(G4int coupleIdx);
00121
00122 inline void ResetForcedInteraction();
00123
00124 private:
00125
00126 void ApplyRangeCut(std::vector<G4DynamicParticle*>& vd,
00127 const G4Track& track,
00128 G4double& eloss,
00129 G4double safety);
00130
00131 G4double ApplySplitting(std::vector<G4DynamicParticle*>& vd,
00132 const G4Track& track,
00133 G4VEmModel* currentModel,
00134 G4int index,
00135 G4double tcut);
00136
00137 inline G4double ApplyRussianRoulette(std::vector<G4DynamicParticle*>& vd,
00138 G4int index);
00139
00140
00141 G4EmBiasingManager(G4EmBiasingManager &);
00142 G4EmBiasingManager & operator=(const G4EmBiasingManager &right);
00143
00144 G4int nForcedRegions;
00145 G4int nSecBiasedRegions;
00146 std::vector<const G4Region*> forcedRegions;
00147 std::vector<G4double> lengthForRegion;
00148 std::vector<const G4Region*> secBiasedRegions;
00149 std::vector<G4double> secBiasedWeight;
00150 std::vector<G4double> secBiasedEnegryLimit;
00151 std::vector<G4int> nBremSplitting;
00152
00153 std::vector<G4int> idxForcedCouple;
00154 std::vector<G4int> idxSecBiasedCouple;
00155
00156 std::vector<G4DynamicParticle*> tmpSecondaries;
00157
00158 G4VEnergyLossProcess* eIonisation;
00159
00160 const G4ParticleDefinition* theElectron;
00161
00162 G4double fSafetyMin;
00163 G4double currentStepLimit;
00164 G4bool startTracking;
00165 };
00166
00167 inline G4bool
00168 G4EmBiasingManager::SecondaryBiasingRegion(G4int coupleIdx)
00169 {
00170 G4bool res = false;
00171 if(nSecBiasedRegions > 0) {
00172 if(idxSecBiasedCouple[coupleIdx] >= 0) { res = true; }
00173 }
00174 return res;
00175 }
00176
00177 inline G4bool G4EmBiasingManager::ForcedInteractionRegion(G4int coupleIdx)
00178 {
00179 G4bool res = false;
00180 if(nForcedRegions > 0) {
00181 if(idxForcedCouple[coupleIdx] >= 0) { res = true; }
00182 }
00183 return res;
00184 }
00185
00186 inline void G4EmBiasingManager::ResetForcedInteraction()
00187 {
00188 startTracking = true;
00189 }
00190
00191
00192
00193 inline G4double
00194 G4EmBiasingManager::ApplyRussianRoulette(std::vector<G4DynamicParticle*>& vd,
00195 G4int index)
00196 {
00197 size_t n = vd.size();
00198 G4double weight = secBiasedWeight[index];
00199 for(size_t k=0; k<n; ++k) {
00200 if(G4UniformRand()*weight > 1.0) {
00201 const G4DynamicParticle* dp = vd[k];
00202 delete dp;
00203 vd[k] = 0;
00204 }
00205 }
00206 return weight;
00207 }
00208
00209
00210
00211 #endif