Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // >> 26 // $Id: G4eBremsstrahlungRelModel.hh 98737 2016-08-09 12:51:38Z gcosmo $ 26 // 27 // 27 // ------------------------------------------- 28 // ------------------------------------------------------------------- 28 // 29 // 29 // GEANT4 Class header file 30 // GEANT4 Class header file 30 // 31 // 31 // 32 // 32 // File name: G4eBremsstrahlungRelModel 33 // File name: G4eBremsstrahlungRelModel 33 // extention of standard G4eBre 34 // extention of standard G4eBremsstrahlungModel 34 // 35 // 35 // Author: Andreas Schaelicke << 36 // Author: Andreas Schaelicke 36 // 37 // 37 // Creation date: 28.03.2008 38 // Creation date: 28.03.2008 38 // 39 // 39 // Modifications: 40 // Modifications: 40 // 41 // 41 // 15.07.18 introduced data structures to sto << 42 // dent data for faster run-time com << 43 // 42 // 44 // Class Description: 43 // Class Description: 45 // 44 // 46 // Implementation of energy loss for gamma emi 45 // Implementation of energy loss for gamma emission by electrons and 47 // positrons including an improved version of 46 // positrons including an improved version of the LPM effect 48 47 49 // ------------------------------------------- 48 // ------------------------------------------------------------------- 50 // 49 // 51 50 52 #ifndef G4eBremsstrahlungRelModel_h 51 #ifndef G4eBremsstrahlungRelModel_h 53 #define G4eBremsstrahlungRelModel_h 1 52 #define G4eBremsstrahlungRelModel_h 1 54 53 55 #include "G4VEmModel.hh" 54 #include "G4VEmModel.hh" 56 #include <memory> << 55 #include "G4NistManager.hh" >> 56 #include "G4Exp.hh" >> 57 #include "G4Log.hh" 57 58 58 class G4ParticleChangeForLoss; 59 class G4ParticleChangeForLoss; >> 60 class G4PhysicsVector; 59 61 60 class G4eBremsstrahlungRelModel : public G4VEm << 62 class G4eBremsstrahlungRelModel : public G4VEmModel >> 63 { 61 64 62 public: 65 public: 63 66 64 explicit G4eBremsstrahlungRelModel(const G4P << 67 explicit G4eBremsstrahlungRelModel(const G4ParticleDefinition* p = 0, 65 const G4S << 68 const G4String& nam = "eBremLPM"); 66 69 67 ~G4eBremsstrahlungRelModel() override; << 70 virtual ~G4eBremsstrahlungRelModel(); 68 71 69 void Initialise(const G4ParticleDefinition*, << 72 virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&) override; >> 73 >> 74 virtual void InitialiseLocal(const G4ParticleDefinition*, >> 75 G4VEmModel* masterModel) override; >> 76 >> 77 virtual G4double ComputeDEDXPerVolume(const G4Material*, >> 78 const G4ParticleDefinition*, >> 79 G4double kineticEnergy, >> 80 G4double cutEnergy) override; >> 81 >> 82 virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition*, >> 83 G4double tkin, >> 84 G4double Z, G4double, >> 85 G4double cutEnergy, >> 86 G4double maxEnergy = DBL_MAX) override; >> 87 >> 88 virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*, >> 89 const G4MaterialCutsCouple*, >> 90 const G4DynamicParticle*, >> 91 G4double cutEnergy, >> 92 G4double maxEnergy) override; >> 93 >> 94 virtual void SetupForMaterial(const G4ParticleDefinition*, >> 95 const G4Material*,G4double) override; >> 96 >> 97 virtual G4double MinPrimaryEnergy(const G4Material*, >> 98 const G4ParticleDefinition*, >> 99 G4double cut) override; >> 100 >> 101 inline void SetLPMconstant(G4double val); >> 102 inline G4double LPMconstant() const; >> 103 >> 104 inline void SetLowestKinEnergy(G4double); >> 105 inline G4double LowestKinEnergy() const; 70 106 71 void InitialiseLocal(const G4ParticleDefinit << 72 G4VEmModel* masterModel << 73 << 74 G4double ComputeDEDXPerVolume(const G4Materi << 75 const G4Partic << 76 G4double ekin, << 77 G4double cutEn << 78 << 79 G4double ComputeCrossSectionPerAtom(const G4 << 80 G4double << 81 G4double << 82 G4double << 83 G4double << 84 G4double << 85 << 86 void SampleSecondaries(std::vector<G4Dynamic << 87 const G4MaterialCutsC << 88 const G4DynamicPartic << 89 G4double cutEnergy, << 90 G4double maxEnergy) o << 91 << 92 void SetupForMaterial(const G4ParticleDefini << 93 const G4Material*, G4d << 94 << 95 G4double MinPrimaryEnergy(const G4Material*, << 96 const G4ParticleDe << 97 G4double cutEnergy << 98 107 99 protected: 108 protected: 100 109 101 virtual G4double ComputeDXSectionPerAtom(G4d 110 virtual G4double ComputeDXSectionPerAtom(G4double gammaEnergy); 102 111 103 void SetParticle(const G4ParticleDefinition* << 112 // * fast inline functions * >> 113 inline void SetCurrentElement(G4int); 104 114 105 private: 115 private: 106 116 >> 117 void InitialiseConstants(); >> 118 >> 119 void CalcLPMFunctions(G4double gammaEnergy); >> 120 107 G4double ComputeBremLoss(G4double cutEnergy) 121 G4double ComputeBremLoss(G4double cutEnergy); 108 122 109 G4double ComputeXSectionPerAtom(G4double cut 123 G4double ComputeXSectionPerAtom(G4double cutEnergy); 110 124 111 G4double ComputeRelDXSectionPerAtom(G4double 125 G4double ComputeRelDXSectionPerAtom(G4double gammaEnergy); 112 126 113 // init special data per element i.e. per Z << 127 void SetParticle(const G4ParticleDefinition* p); 114 void InitialiseElementData(); << 115 128 116 // methods for initialisation and run-time e << 129 inline G4double Phi1(G4double,G4double); 117 void InitLPMFunctions(); << 130 inline G4double Phi1M2(G4double,G4double); >> 131 inline G4double Psi1(G4double,G4double); >> 132 inline G4double Psi1M2(G4double,G4double); 118 133 119 void ComputeLPMfunctions(G4double& funcXiS, << 134 // hide assignment operator 120 G4double& funcGS, << 135 G4eBremsstrahlungRelModel & operator=(const G4eBremsstrahlungRelModel &right) = delete; 121 G4double& funcPhiS, << 122 const G4double egam << 123 << 124 void GetLPMFunctions(G4double& lpmGs, << 125 G4double& lpmPhis, << 126 const G4double ss); << 127 << 128 void ComputeLPMGsPhis(G4double& funcGS, << 129 G4double& funcPhiS, << 130 const G4double varShat << 131 << 132 // for evaluating screening related function << 133 void ComputeScreeningFunctions(G4double& phi << 134 G4double& phi << 135 G4double& psi << 136 G4double& psi << 137 const G4doubl << 138 const G4doubl << 139 // hide assignment operator and cctr << 140 G4eBremsstrahlungRelModel& operator= << 141 (const G4eBremsstrahlungRelModel& right) = d << 142 G4eBremsstrahlungRelModel(const G4eBremsstr 136 G4eBremsstrahlungRelModel(const G4eBremsstrahlungRelModel&) = delete; 143 137 144 private: << 138 protected: 145 139 146 G4bool fIsUseCompleteSc << 140 G4NistManager* nist; 147 G4bool fIsInitializer = << 141 const G4ParticleDefinition* particle; 148 G4bool fUseLPM = true; << 142 G4ParticleDefinition* theGamma; >> 143 G4ParticleChangeForLoss* fParticleChange; 149 144 150 protected: << 145 G4double bremFactor; 151 146 152 G4bool fIsElectron = tr << 153 G4bool fIsScatOffElectr << 154 G4bool fIsLPMActive = f << 155 // << 156 G4int fCurrentIZ = 0; << 157 const G4ParticleDefinition* fPrimaryParticle << 158 G4ParticleDefinition* fGammaParticle = << 159 G4ParticleChangeForLoss* fParticleChange << 160 // cash 147 // cash 161 G4double fPrimaryParticle << 148 G4double particleMass; 162 G4double fPrimaryKinEnerg << 149 G4double kinEnergy; 163 G4double fPrimaryTotalEne << 150 G4double totalEnergy; 164 G4double fDensityFactor = << 151 G4double densityFactor; 165 G4double fDensityCorr = 0 << 152 G4double densityCorr; 166 G4double fLowestKinEnergy << 153 167 // scattering off electrons << 154 G4int currentZ; 168 G4double fNucTerm = 0.; << 155 G4bool isElectron; 169 G4double fSumTerm = 0.; << 170 156 171 private: 157 private: 172 158 173 // LPM related members << 159 static const G4double xgi[8], wgi[8]; 174 G4double fLPMEnergyThresh << 160 static const G4double Fel_light[5]; 175 G4double fLPMEnergy; << 161 static const G4double Finel_light[5]; >> 162 >> 163 // consts >> 164 G4double lowestKinEnergy; >> 165 G4double fMigdalConstant; >> 166 G4double fLPMconstant; >> 167 G4double energyThresholdLPM; >> 168 G4double facFel, facFinel; >> 169 G4double preS1,logTwo; 176 170 177 protected: << 171 // cash >> 172 G4double z13, z23, lnZ; >> 173 G4double Fel, Finel, fCoulomb, fMax; 178 174 179 static const G4double gBremFactor; << 175 // LPM effect 180 static const G4double gMigdalConstant; << 176 G4double lpmEnergy; >> 177 G4double xiLPM, phiLPM, gLPM; 181 178 182 private: << 179 // critical gamma energies >> 180 G4double klpm, kp; 183 181 184 static const G4int gMaxZet; << 182 // flags 185 // << 183 G4bool use_completescreening; 186 static const G4double gLPMconstant; << 187 // << 188 static const G4double gXGL[8]; << 189 static const G4double gWGL[8]; << 190 static const G4double gFelLowZet[8]; << 191 static const G4double gFinelLowZet[8]; << 192 // << 193 struct ElementData { << 194 /** @brief \f$ \ln(Z) \f$ */ << 195 G4double fLogZ; << 196 /** @brief \f$ \ln(Z)/3 + f_c \f$ */ << 197 G4double fFz; << 198 /** @brief \f$ ((Fel-fc)+Finel*invZ)\f$ * << 199 G4double fZFactor1; << 200 /** @brief \f$ (Fel-fc)\f$ */ << 201 G4double fZFactor11; << 202 /** @brief \f$ (1.0+invZ)/12 \f$ */ << 203 G4double fZFactor2; << 204 // LPM variables << 205 G4double fVarS1; << 206 G4double fILVarS1; << 207 G4double fILVarS1Cond; << 208 // constant factors to the screening funct << 209 G4double fGammaFactor; << 210 G4double fEpsilonFactor; << 211 }; << 212 // << 213 struct LPMFuncs { << 214 LPMFuncs() : fIsInitialized(false), fISDel << 215 G4bool fIsInitialized; << 216 G4double fISDelta; << 217 G4double fSLimit; << 218 std::vector<G4double> fLPMFuncG; << 219 std::vector<G4double> fLPMFuncPhi; << 220 }; << 221 // << 222 static std::shared_ptr<LPMFuncs> gLPMFuncs() << 223 static std::shared_ptr<std::vector<ElementDa << 224 std::shared_ptr<LPMFuncs> fLPMFuncs; << 225 std::shared_ptr<std::vector<ElementData*>> f << 226 }; 184 }; >> 185 >> 186 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... >> 187 >> 188 inline void G4eBremsstrahlungRelModel::SetCurrentElement(G4int Z) >> 189 { >> 190 if(Z != currentZ) { >> 191 currentZ = Z; >> 192 >> 193 z13 = nist->GetZ13(Z); >> 194 z23 = z13*z13; >> 195 lnZ = nist->GetLOGZ(Z); >> 196 >> 197 if (Z <= 4) { >> 198 Fel = Fel_light[Z]; >> 199 Finel = Finel_light[Z] ; >> 200 } >> 201 else { >> 202 G4double lnzt = lnZ/3.; >> 203 Fel = facFel - lnzt; >> 204 Finel = facFinel - 2*lnzt; >> 205 } >> 206 >> 207 fCoulomb = GetCurrentElement()->GetfCoulomb(); >> 208 G4double xz = 1.0/(G4double)Z; >> 209 fMax = Fel-fCoulomb + Finel*xz + (1. + xz)/12.; >> 210 } >> 211 } >> 212 >> 213 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... >> 214 >> 215 >> 216 inline G4double G4eBremsstrahlungRelModel::Phi1(G4double gg, G4double) >> 217 { >> 218 // Thomas-Fermi FF from Tsai, eq.(3.38) for Z>=5 >> 219 return 20.863 - 2.*G4Log(1. + sqr(0.55846*gg) ) >> 220 - 4.*( 1. - 0.6*G4Exp(-0.9*gg) - 0.4*G4Exp(-1.5*gg) ); >> 221 } >> 222 >> 223 inline G4double G4eBremsstrahlungRelModel::Phi1M2(G4double gg, G4double) >> 224 { >> 225 // Thomas-Fermi FF from Tsai, eq. (3.39) for Z>=5 >> 226 // return Phi1(gg,Z) - >> 227 return 2./(3.*(1. + 6.5*gg +6.*gg*gg) ); >> 228 } >> 229 >> 230 inline G4double G4eBremsstrahlungRelModel::Psi1(G4double eps, G4double) >> 231 { >> 232 // Thomas-Fermi FF from Tsai, eq.(3.40) for Z>=5 >> 233 return 28.340 - 2.*G4Log(1. + sqr(3.621*eps) ) >> 234 - 4.*( 1. - 0.7*G4Exp(-8*eps) - 0.3*G4Exp(-29.*eps) ); >> 235 } >> 236 >> 237 inline G4double G4eBremsstrahlungRelModel::Psi1M2(G4double eps, G4double) >> 238 { >> 239 // Thomas-Fermi FF from Tsai, eq. (3.41) for Z>=5 >> 240 return 2./(3.*(1. + 40.*eps +400.*eps*eps) ); >> 241 } >> 242 >> 243 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... >> 244 >> 245 inline >> 246 void G4eBremsstrahlungRelModel::SetLPMconstant(G4double val) >> 247 { >> 248 fLPMconstant = val; >> 249 } >> 250 >> 251 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... >> 252 >> 253 inline >> 254 G4double G4eBremsstrahlungRelModel::LPMconstant() const >> 255 { >> 256 return fLPMconstant; >> 257 } >> 258 >> 259 inline void G4eBremsstrahlungRelModel::SetLowestKinEnergy(G4double val) >> 260 { >> 261 lowestKinEnergy = val; >> 262 } >> 263 >> 264 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... >> 265 >> 266 inline G4double G4eBremsstrahlungRelModel::LowestKinEnergy() const >> 267 { >> 268 return lowestKinEnergy; >> 269 } >> 270 >> 271 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... >> 272 227 273 228 #endif 274 #endif 229 275