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: G4TransitionRadiation.cc,v 1.10 2010-10-14 18:38:21 vnivanch Exp $ >> 27 // GEANT4 tag $Name: not supported by cvs2svn $ >> 28 // 26 // G4TransitionRadiation class -- implementati 29 // G4TransitionRadiation class -- implementation file 27 30 28 // GEANT 4 class implementation file --- Copyr 31 // GEANT 4 class implementation file --- Copyright CERN 1995 >> 32 // CERN Geneva Switzerland 29 33 30 // For information related to this code, pleas 34 // For information related to this code, please, contact 31 // CERN, CN Division, ASD Group 35 // CERN, CN Division, ASD Group 32 // History: 36 // History: 33 // 1st version 11.09.97 V. Grichine (Vladimir. 37 // 1st version 11.09.97 V. Grichine (Vladimir.Grichine@cern.ch ) 34 // 2nd version 16.12.97 V. Grichine 38 // 2nd version 16.12.97 V. Grichine 35 // 3rd version 28.07.05, P.Gumplinger add G4Pr 39 // 3rd version 28.07.05, P.Gumplinger add G4ProcessType to constructor 36 40 37 //#include <cmath> << 38 41 39 #include "G4TransitionRadiation.hh" << 42 #include <cmath> 40 43 >> 44 #include "G4TransitionRadiation.hh" >> 45 #include "G4Material.hh" 41 #include "G4EmProcessSubType.hh" 46 #include "G4EmProcessSubType.hh" 42 47 >> 48 // Local constants >> 49 >> 50 const G4int G4TransitionRadiation::fSympsonNumber = 100 ; >> 51 const G4int G4TransitionRadiation::fGammaNumber = 15 ; >> 52 const G4int G4TransitionRadiation::fPointNumber = 100 ; >> 53 >> 54 43 ////////////////////////////////////////////// 55 /////////////////////////////////////////////////////////////////////// >> 56 // 44 // Constructor for selected couple of material 57 // Constructor for selected couple of materials 45 G4TransitionRadiation::G4TransitionRadiation(c << 58 // 46 G << 59 >> 60 G4TransitionRadiation:: >> 61 G4TransitionRadiation( const G4String& processName, G4ProcessType type ) 47 : G4VDiscreteProcess(processName, type) 62 : G4VDiscreteProcess(processName, type) 48 { 63 { 49 SetProcessSubType(fTransitionRadiation); 64 SetProcessSubType(fTransitionRadiation); 50 fMatIndex1 = fMatIndex2 = 0; 65 fMatIndex1 = fMatIndex2 = 0; 51 66 52 fGamma = fEnergy = fVarAngle = fMinEnergy = << 67 fGamma = fEnergy = fVarAngle = fMinEnergy = fMaxEnergy = fMaxTheta = fSigma1 = fSigma2 = 0.0; 53 fSigma1 = fSigma2 = 0.0; << 54 } 68 } 55 69 56 ////////////////////////////////////////////// 70 ////////////////////////////////////////////////////////////////////// >> 71 // 57 // Destructor 72 // Destructor 58 G4TransitionRadiation::~G4TransitionRadiation( << 73 // 59 74 60 void G4TransitionRadiation::ProcessDescription << 75 G4TransitionRadiation::~G4TransitionRadiation() 61 { << 76 {} 62 out << "Base class for simulation of x-ray t << 63 } << 64 77 65 G4bool G4TransitionRadiation::IsApplicable( << 78 G4bool 66 const G4ParticleDefinition& aParticleType) << 79 G4TransitionRadiation::IsApplicable(const G4ParticleDefinition& aParticleType) 67 { 80 { 68 return (aParticleType.GetPDGCharge() != 0.0) << 81 return ( aParticleType.GetPDGCharge() != 0.0 ); 69 } 82 } 70 83 71 G4double G4TransitionRadiation::GetMeanFreePat << 84 G4double G4TransitionRadiation::GetMeanFreePath(const G4Track&, 72 << 85 G4double, >> 86 G4ForceCondition* condition) 73 { 87 { 74 *condition = Forced; 88 *condition = Forced; 75 return DBL_MAX; // so TR doesn't limit mean << 89 return DBL_MAX; // so TR doesn't limit mean free path 76 } 90 } 77 91 78 G4VParticleChange* G4TransitionRadiation::Post 92 G4VParticleChange* G4TransitionRadiation::PostStepDoIt(const G4Track&, 79 << 93 const G4Step&) 80 { 94 { 81 ClearNumberOfInteractionLengthLeft(); 95 ClearNumberOfInteractionLengthLeft(); 82 return &aParticleChange; 96 return &aParticleChange; 83 } 97 } 84 98 85 ////////////////////////////////////////////// 99 /////////////////////////////////////////////////////////////////// >> 100 // 86 // Sympson integral of TR spectral-angle densi 101 // Sympson integral of TR spectral-angle density over energy between 87 // the limits energy 1 and energy2 at fixed va 102 // the limits energy 1 and energy2 at fixed varAngle = 1 - std::cos(Theta) 88 G4double G4TransitionRadiation::IntegralOverEn << 103 89 << 104 G4double 90 << 105 G4TransitionRadiation::IntegralOverEnergy( G4double energy1, 91 { << 106 G4double energy2, 92 G4int i; << 107 G4double varAngle ) const 93 G4double h, sumEven = 0.0, sumOdd = 0.0; << 108 { 94 h = 0.5 * (energy2 - energy1) / fSympsonNumb << 109 G4int i ; 95 for(i = 1; i < fSympsonNumber; i++) << 110 G4double h , sumEven = 0.0 , sumOdd = 0.0 ; >> 111 h = 0.5*(energy2 - energy1)/fSympsonNumber ; >> 112 for(i=1;i<fSympsonNumber;i++) 96 { 113 { 97 sumEven += SpectralAngleTRdensity(energy1 << 114 sumEven += SpectralAngleTRdensity(energy1 + 2*i*h,varAngle) ; 98 sumOdd += SpectralAngleTRdensity(energy1 + << 115 sumOdd += SpectralAngleTRdensity(energy1 + (2*i - 1)*h,varAngle) ; 99 } 116 } 100 sumOdd += << 117 sumOdd += SpectralAngleTRdensity(energy1 + (2*fSympsonNumber - 1)*h,varAngle) ; 101 SpectralAngleTRdensity(energy1 + (2 * fSym << 118 return h*( SpectralAngleTRdensity(energy1,varAngle) 102 return h * << 119 + SpectralAngleTRdensity(energy2,varAngle) 103 (SpectralAngleTRdensity(energy1, varA << 120 + 4.0*sumOdd + 2.0*sumEven )/3.0 ; 104 SpectralAngleTRdensity(energy2, varA << 105 2.0 * sumEven) / << 106 3.0; << 107 } 121 } 108 122 >> 123 >> 124 109 ////////////////////////////////////////////// 125 /////////////////////////////////////////////////////////////////// >> 126 // 110 // Sympson integral of TR spectral-angle densi 127 // Sympson integral of TR spectral-angle density over energy between 111 // the limits varAngle1 and varAngle2 at fixed 128 // the limits varAngle1 and varAngle2 at fixed energy 112 G4double G4TransitionRadiation::IntegralOverAn << 129 113 << 130 G4double 114 << 131 G4TransitionRadiation::IntegralOverAngle( G4double energy, 115 { << 132 G4double varAngle1, 116 G4int i; << 133 G4double varAngle2 ) const 117 G4double h, sumEven = 0.0, sumOdd = 0.0; << 134 { 118 h = 0.5 * (varAngle2 - varAngle1) / fSympson << 135 G4int i ; 119 for(i = 1; i < fSympsonNumber; ++i) << 136 G4double h , sumEven = 0.0 , sumOdd = 0.0 ; >> 137 h = 0.5*(varAngle2 - varAngle1)/fSympsonNumber ; >> 138 for(i=1;i<fSympsonNumber;i++) 120 { 139 { 121 sumEven += SpectralAngleTRdensity(energy, << 140 sumEven += SpectralAngleTRdensity(energy,varAngle1 + 2*i*h) ; 122 sumOdd += SpectralAngleTRdensity(energy, v << 141 sumOdd += SpectralAngleTRdensity(energy,varAngle1 + (2*i - 1)*h) ; 123 } 142 } 124 sumOdd += << 143 sumOdd += SpectralAngleTRdensity(energy,varAngle1 + (2*fSympsonNumber - 1)*h) ; 125 SpectralAngleTRdensity(energy, varAngle1 + << 126 144 127 return h * << 145 return h*( SpectralAngleTRdensity(energy,varAngle1) 128 (SpectralAngleTRdensity(energy, varAn << 146 + SpectralAngleTRdensity(energy,varAngle2) 129 SpectralAngleTRdensity(energy, varAn << 147 + 4.0*sumOdd + 2.0*sumEven )/3.0 ; 130 2.0 * sumEven) / << 131 3.0; << 132 } 148 } 133 149 134 ////////////////////////////////////////////// 150 /////////////////////////////////////////////////////////////////// >> 151 // 135 // The number of transition radiation photons 152 // The number of transition radiation photons generated in the 136 // angle interval between varAngle1 and varAng 153 // angle interval between varAngle1 and varAngle2 137 G4double G4TransitionRadiation::AngleIntegralD << 154 // 138 G4double varAngle1, G4double varAngle2) cons << 155 139 { << 156 G4double G4TransitionRadiation:: 140 G4int i; << 157 AngleIntegralDistribution( G4double varAngle1, 141 G4double h, sumEven = 0.0, sumOdd = 0.0; << 158 G4double varAngle2 ) const 142 h = 0.5 * (varAngle2 - varAngle1) / fSympson << 159 { 143 for(i = 1; i < fSympsonNumber; ++i) << 160 G4int i ; >> 161 G4double h , sumEven = 0.0 , sumOdd = 0.0 ; >> 162 h = 0.5*(varAngle2 - varAngle1)/fSympsonNumber ; >> 163 for(i=1;i<fSympsonNumber;i++) 144 { 164 { 145 sumEven += IntegralOverEnergy(fMinEnergy, << 165 sumEven += IntegralOverEnergy(fMinEnergy, 146 fMinEnergy + << 166 fMinEnergy +0.3*(fMaxEnergy-fMinEnergy), 147 varAngle1 + << 167 varAngle1 + 2*i*h) 148 IntegralOverEnergy(fMinEnergy + << 168 + IntegralOverEnergy(fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 149 fMaxEnergy, << 169 fMaxEnergy, 150 sumOdd += IntegralOverEnergy(fMinEnergy, << 170 varAngle1 + 2*i*h); 151 fMinEnergy + << 171 sumOdd += IntegralOverEnergy(fMinEnergy, 152 varAngle1 + ( << 172 fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 153 IntegralOverEnergy(fMinEnergy + << 173 varAngle1 + (2*i - 1)*h) 154 fMaxEnergy, v << 174 + IntegralOverEnergy(fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), >> 175 fMaxEnergy, >> 176 varAngle1 + (2*i - 1)*h) ; 155 } 177 } 156 sumOdd += << 178 sumOdd += IntegralOverEnergy(fMinEnergy, 157 IntegralOverEnergy(fMinEnergy, fMinEnergy << 179 fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 158 varAngle1 + (2 * fSymps << 180 varAngle1 + (2*fSympsonNumber - 1)*h) 159 IntegralOverEnergy(fMinEnergy + 0.3 * (fMa << 181 + IntegralOverEnergy(fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 160 varAngle1 + (2 * fSymps << 182 fMaxEnergy, 161 << 183 varAngle1 + (2*fSympsonNumber - 1)*h) ; 162 return h * << 184 163 (IntegralOverEnergy(fMinEnergy, << 185 return h*(IntegralOverEnergy(fMinEnergy, 164 fMinEnergy + 0.3 << 186 fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 165 varAngle1) + << 187 varAngle1) 166 IntegralOverEnergy(fMinEnergy + 0.3 << 188 + IntegralOverEnergy(fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 167 fMaxEnergy, varAn << 189 fMaxEnergy, 168 IntegralOverEnergy(fMinEnergy, << 190 varAngle1) 169 fMinEnergy + 0.3 << 191 + IntegralOverEnergy(fMinEnergy, 170 varAngle2) + << 192 fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 171 IntegralOverEnergy(fMinEnergy + 0.3 << 193 varAngle2) 172 fMaxEnergy, varAn << 194 + IntegralOverEnergy(fMinEnergy + 0.3*(fMaxEnergy - fMinEnergy), 173 4.0 * sumOdd + 2.0 * sumEven) / << 195 fMaxEnergy, 174 3.0; << 196 varAngle2) >> 197 + 4.0*sumOdd + 2.0*sumEven )/3.0 ; 175 } 198 } 176 199 177 ////////////////////////////////////////////// 200 /////////////////////////////////////////////////////////////////// >> 201 // 178 // The number of transition radiation photons, 202 // The number of transition radiation photons, generated in the 179 // energy interval between energy1 and energy2 203 // energy interval between energy1 and energy2 180 G4double G4TransitionRadiation::EnergyIntegral << 204 // 181 G4double energy1, G4double energy2) const << 205 182 { << 206 G4double G4TransitionRadiation:: 183 G4int i; << 207 EnergyIntegralDistribution( G4double energy1, 184 G4double h, sumEven = 0.0, sumOdd = 0.0; << 208 G4double energy2 ) const 185 h = 0.5 * (energy2 - energy1) / fSympsonNumb << 209 { 186 for(i = 1; i < fSympsonNumber; ++i) << 210 G4int i ; >> 211 G4double h , sumEven = 0.0 , sumOdd = 0.0 ; >> 212 h = 0.5*(energy2 - energy1)/fSympsonNumber ; >> 213 for(i=1;i<fSympsonNumber;i++) 187 { 214 { 188 sumEven += << 215 sumEven += IntegralOverAngle(energy1 + 2*i*h,0.0,0.01*fMaxTheta ) 189 IntegralOverAngle(energy1 + 2 * i * h, 0 << 216 + IntegralOverAngle(energy1 + 2*i*h,0.01*fMaxTheta,fMaxTheta); 190 IntegralOverAngle(energy1 + 2 * i * h, 0 << 217 sumOdd += IntegralOverAngle(energy1 + (2*i - 1)*h,0.0,0.01*fMaxTheta) 191 sumOdd += << 218 + IntegralOverAngle(energy1 + (2*i - 1)*h,0.01*fMaxTheta,fMaxTheta) ; 192 IntegralOverAngle(energy1 + (2 * i - 1) << 193 IntegralOverAngle(energy1 + (2 * i - 1) << 194 } 219 } 195 sumOdd += IntegralOverAngle(energy1 + (2 * f << 220 sumOdd += IntegralOverAngle(energy1 + (2*fSympsonNumber - 1)*h, 196 0.01 * fMaxTheta << 221 0.0,0.01*fMaxTheta) 197 IntegralOverAngle(energy1 + (2 * f << 222 + IntegralOverAngle(energy1 + (2*fSympsonNumber - 1)*h, 198 0.01 * fMaxTheta << 223 0.01*fMaxTheta,fMaxTheta) ; 199 << 224 200 return h * << 225 return h*(IntegralOverAngle(energy1,0.0,0.01*fMaxTheta) 201 (IntegralOverAngle(energy1, 0.0, 0.01 << 226 + IntegralOverAngle(energy1,0.01*fMaxTheta,fMaxTheta) 202 IntegralOverAngle(energy1, 0.01 * fM << 227 + IntegralOverAngle(energy2,0.0,0.01*fMaxTheta) 203 IntegralOverAngle(energy2, 0.0, 0.01 << 228 + IntegralOverAngle(energy2,0.01*fMaxTheta,fMaxTheta) 204 IntegralOverAngle(energy2, 0.01 * fM << 229 + 4.0*sumOdd + 2.0*sumEven )/3.0 ; 205 4.0 * sumOdd + 2.0 * sumEven) / << 206 3.0; << 207 } 230 } >> 231 >> 232 >> 233 >> 234 >> 235 // end of G4TransitionRadiation implementation file -------------------------- 208 236