Geant4  10.00.p01
G4VEmAdjointModel.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 // $Id: G4VEmAdjointModel.hh 66892 2013-01-17 10:57:59Z gunter $
27 //
29 // Module: G4VEMAdjointModel
30 // Author: L. Desorgher
31 // Organisation: SpaceIT GmbH
32 // Contract: ESA contract 21435/08/NL/AT
33 // Customer: ESA/ESTEC
35 //
36 // CHANGE HISTORY
37 // --------------
38 // ChangeHistory:
39 // 10 September 2009 Move to a virtual class. L. Desorgher
40 // 1st April 2007 creation by L. Desorgher
41 //
42 //-------------------------------------------------------------
43 // Documentation:
44 // Base class for Adjoint EM model. It is based on the use of direct G4VEmModel.
45 //
46 
47 
48 #ifndef G4VEmAdjointModel_h
49 #define G4VEmAdjointModel_h 1
50 
51 #include "globals.hh"
52 #include "G4DynamicParticle.hh"
53 #include "G4ParticleDefinition.hh"
54 #include "G4MaterialCutsCouple.hh"
55 #include "G4Material.hh"
56 #include "G4Element.hh"
57 #include "G4ElementVector.hh"
58 #include "Randomize.hh"
59 #include "G4ParticleDefinition.hh"
60 #include "G4VEmModel.hh"
61 #include "G4Electron.hh"
62 #include "G4Gamma.hh"
63 #include "G4ProductionCutsTable.hh"
64 
65 class G4PhysicsTable;
66 class G4Region;
67 class G4VParticleChange;
68 class G4ParticleChange;
69 class G4Track;
70 class G4AdjointCSMatrix;
71 
73 {
74 
75 public: // public methods
76 
77  G4VEmAdjointModel(const G4String& nam);
78 
79  virtual ~G4VEmAdjointModel();
80 
81  //------------------------------------------------------------------------
82  // Virtual methods to be implemented for the sample secondaries concrete model
83  //------------------------------------------------------------------------
84 
85  //virtual void Initialise()=0;
86 
87  virtual void SampleSecondaries(const G4Track& aTrack,
88  G4bool IsScatProjToProjCase,
89  G4ParticleChange* fParticleChange)=0;
90 
91 
92  //------------------------------------------------------------------------
93  // Methods for adjoint processes; may be overwritten if needed;
94  //------------------------------------------------------------------------
95 
96 
97  virtual G4double AdjointCrossSection(const G4MaterialCutsCouple* aCouple,
98  G4double primEnergy,
99  G4bool IsScatProjToProjCase);
100 
102  G4double primEnergy,
103  G4bool IsScatProjToProjCase);
104 
106  G4double kinEnergyProj, // kinetic energy of the primary particle before the interaction
107  G4double kinEnergyProd, // kinetic energy of the secondary particle
108  G4double Z,
109  G4double A = 0.);
110 
112  G4double kinEnergyProj, // kinetic energy of the primary particle before the interaction
113  G4double kinEnergyScatProj, // kinetic energy of the primary particle after the interaction
114  G4double Z,
115  G4double A = 0.);
116 
117 
118 
120  const G4Material* aMaterial,
121  G4double kinEnergyProj, // kinetic energy of the primary particle before the interaction
122  G4double kinEnergyProd // kinetic energy of the secondary particle
123  );
124 
126  const G4Material* aMaterial,
127  G4double kinEnergyProj, // kinetic energy of the primary particle before the interaction
128  G4double kinEnergyScatProj // kinetic energy of the primary particle after the interaction
129  );
130 
131 
132  //Energy limits of adjoint secondary
133  //------------------
134 
139 
140 
141 
142  //Other Methods
143  //---------------
144 
145  void DefineCurrentMaterial(const G4MaterialCutsCouple* couple);
146 
147 
148  std::vector< std::vector< double>* > ComputeAdjointCrossSectionVectorPerAtomForSecond(
149  G4double kinEnergyProd,
150  G4double Z,
151  G4double A = 0.,
152  G4int nbin_pro_decade=10
153  );
154  std::vector< std::vector< double>* > ComputeAdjointCrossSectionVectorPerAtomForScatProj(
155  G4double kinEnergyProd,
156  G4double Z,
157  G4double A = 0.,
158  G4int nbin_pro_decade=10
159  );
160 
161  std::vector< std::vector< double>* > ComputeAdjointCrossSectionVectorPerVolumeForSecond(
162  G4Material* aMaterial,
163  G4double kinEnergyProd,
164  G4int nbin_pro_decade=10
165  );
166  std::vector< std::vector< double>* > ComputeAdjointCrossSectionVectorPerVolumeForScatProj(
167  G4Material* aMaterial,
168  G4double kinEnergyProd,
169  G4int nbin_pro_decade=10
170  );
171 
172 
173 
174  inline void SetCSMatrices(std::vector< G4AdjointCSMatrix* >* Vec1CSMatrix, std::vector< G4AdjointCSMatrix* >* Vec2CSMatrix){
177 
178 
179  };
180 
182 
184 
186 
188 
189  void SetHighEnergyLimit(G4double aVal);
190 
191  void SetLowEnergyLimit(G4double aVal);
192 
193  inline void DefineDirectEMModel(G4VEmModel* aModel){theDirectEMModel = aModel;}
194 
196 
199  }
200 
202 
204 
205  inline void SetUseMatrix(G4bool aBool) { UseMatrix = aBool;}
206 
207  inline void SetUseMatrixPerElement(G4bool aBool){ UseMatrixPerElement = aBool;}
209 
210  inline void SetApplyCutInRange(G4bool aBool){ ApplyCutInRange = aBool;}
211  inline G4bool GetUseMatrix() {return UseMatrix;}
215 
216  inline G4String GetName(){ return name;}
217  inline virtual void SetCSBiasingFactor(G4double aVal) {CS_biasing_factor = aVal;}
218 
219 protected:
220 
221  //Some of them can be overriden by daughter classes
222 
223 
227 
228 
229 
230  //General methods to sample secondary energy
231  //--------------------------------------
232  G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex,G4double prim_energy,G4bool IsScatProjToProjCase);
233  G4double SampleAdjSecEnergyFromCSMatrix(G4double prim_energy,G4bool IsScatProjToProjCase);
234  void SelectCSMatrix(G4bool IsScatProjToProjCase);
235 
236  virtual G4double SampleAdjSecEnergyFromDiffCrossSectionPerAtom(G4double prim_energy,G4bool IsScatProjToProjCase);
237 
238 
239 
240  //Post Step weight correction
241  //----------------------------
242  virtual void CorrectPostStepWeight(G4ParticleChange* fParticleChange,
243  G4double old_weight,
244  G4double adjointPrimKinEnergy,
245  G4double projectileKinEnergy,
246  G4bool IsScatProjToProjCase);
247 
248 
249 
250 
251 
252 
253 protected: //attributes
254 
257 
258 
259 
260 
261  //Name
262  //-----
263 
264  const G4String name;
265 
266  //Needed for CS integration at the initialisation phase
267  //-----------------------------------------------------
268 
275 
276 
277  //for the adjoint simulation we need for each element or material:
278  //an adjoint CS Matrix
279  //-----------------------------
280 
281  std::vector< G4AdjointCSMatrix* >* pOnCSMatrixForProdToProjBackwardScattering;
282  std::vector< G4AdjointCSMatrix* >* pOnCSMatrixForScatProjToProjBackwardScattering;
283  std::vector<G4double> CS_Vs_ElementForScatProjToProjCase;
284  std::vector<G4double> CS_Vs_ElementForProdToProjCase;
285 
289 
290 
291 
292 
293  //particle definition
294  //------------------
295 
300 
301 
302  //Prestep energy
303  //-------------
305 
306  //Current couple material
307  //----------------------
315 
316 
317 
318 
319  //For ions
320  //---------
323 
324 
325  //Energy limits
326  //-------------
327 
330 
331 
332 
333  //Cross Section biasing factor
334  //---------------------------
336 
337 
338  //Type of Model with Matrix or not
339  //--------------------------------
341  G4bool UseMatrixPerElement; //other possibility is per Material
343 
344 
345  //Index of Cross section matrices to be used
346  //------------
348 
349  size_t model_index;
350 
351 
352 
353 
354 
355 
356 
357 
358 };
359 
360 
361 #endif
362 
G4double kinEnergyScatProjForIntegration
virtual ~G4VEmAdjointModel()
virtual G4double DiffCrossSectionPerVolumePrimToSecond(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
G4double lastAdjointCSForScatProjToProjCase
G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
G4ParticleDefinition * theDirectPrimaryPartDef
void SetAdjointEquivalentOfDirectSecondaryParticleDefinition(G4ParticleDefinition *aPart)
virtual G4double GetAdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
G4double currentTcutForDirectPrim
virtual G4double SampleAdjSecEnergyFromDiffCrossSectionPerAtom(G4double prim_energy, G4bool IsScatProjToProjCase)
G4double DiffCrossSectionPerVolumeFunctionForIntegrationOverEkinProj(G4double EkinProd)
virtual G4double DiffCrossSectionPerVolumePrimToScatPrim(const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyScatProj)
G4VParticleChange * pParticleChange
G4bool GetUseOnlyOneMatrixForAllElements()
G4ParticleDefinition * GetAdjointEquivalentOfDirectPrimaryParticleDefinition()
virtual G4double DiffCrossSectionPerAtomPrimToSecond(G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
G4Material * currentMaterial
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
void SetLowEnergyLimit(G4double aVal)
G4Material * SelectedMaterial
int G4int
Definition: G4Types.hh:78
size_t indexOfUsedCrossSectionMatrix
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
void SetHighEnergyLimit(G4double aVal)
virtual void SampleSecondaries(const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)=0
void SetUseMatrixPerElement(G4bool aBool)
G4double DiffCrossSectionFunction1(G4double kinEnergyProj)
G4VEmModel * theDirectEMModel
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
std::vector< G4double > CS_Vs_ElementForProdToProjCase
G4VEmAdjointModel(const G4String &nam)
void SetSecondPartOfSameType(G4bool aBool)
std::vector< G4double > CS_Vs_ElementForScatProjToProjCase
bool G4bool
Definition: G4Types.hh:79
virtual G4double DiffCrossSectionPerAtomPrimToScatPrim(G4double kinEnergyProj, G4double kinEnergyScatProj, G4double Z, G4double A=0.)
G4double GetLowEnergyLimit()
void SetUseMatrix(G4bool aBool)
void SetCSMatrices(std::vector< G4AdjointCSMatrix * > *Vec1CSMatrix, std::vector< G4AdjointCSMatrix * > *Vec2CSMatrix)
void SetAdjointEquivalentOfDirectPrimaryParticleDefinition(G4ParticleDefinition *aPart)
static const G4double A[nN]
G4double DiffCrossSectionFunction2(G4double kinEnergyProj)
G4MaterialCutsCouple * currentCouple
void DefineCurrentMaterial(const G4MaterialCutsCouple *couple)
G4bool UseOnlyOneMatrixForAllElements
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerAtomForScatProj(G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
void SetUseOnlyOneMatrixForAllElements(G4bool aBool)
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForProdToProjBackwardScattering
void DefineDirectEMModel(G4VEmModel *aModel)
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerAtomForSecond(G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
G4double GetHighEnergyLimit()
G4double kinEnergyProjForIntegration
G4ParticleDefinition * GetAdjointEquivalentOfDirectSecondaryParticleDefinition()
void SelectCSMatrix(G4bool IsScatProjToProjCase)
virtual void SetCSBiasingFactor(G4double aVal)
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForScatProjToProjBackwardScattering
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerVolumeForSecond(G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
G4bool GetSecondPartOfSameType()
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
G4double currentTcutForDirectSecond
std::vector< std::vector< double > * > ComputeAdjointCrossSectionVectorPerVolumeForScatProj(G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
G4double lastAdjointCSForProdToProjCase
void SetApplyCutInRange(G4bool aBool)
double G4double
Definition: G4Types.hh:76
virtual void CorrectPostStepWeight(G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
G4double kinEnergyProdForIntegration