Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4GammaXTRadiator.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 //
27 // $Id: G4GammaXTRadiator.cc 97385 2016-06-02 09:59:53Z gcosmo $
28 //
29 
30 #include <complex>
31 
32 #include "G4GammaXTRadiator.hh"
33 #include "Randomize.hh"
34 
35 #include "G4Gamma.hh"
36 
38 //
39 // Constructor, destructor
40 
42  G4double alphaPlate,
43  G4double alphaGas,
44  G4Material* foilMat,G4Material* gasMat,
46  const G4String& processName) :
47  G4VXTRenergyLoss(anEnvelope,foilMat,gasMat,a,b,n,processName)
48 {
49  G4cout<<"Gamma distributed X-ray TR radiator model is called"<<G4endl ;
50 
51  // Build energy and angular integral spectra of X-ray TR photons from
52  // a radiator
53 
54  fAlphaPlate = alphaPlate ;
55  fAlphaGas = alphaGas ;
56  G4cout<<"fAlphaPlate = "<<fAlphaPlate<<" ; fAlphaGas = "<<fAlphaGas<<G4endl ;
57 
58  // BuildTable() ;
59 }
60 
62 
64 {
65  ;
66 }
67 
68 
69 
71 //
72 // Rough approximation for radiator interference factor for the case of
73 // fully GamDistr radiator. The plate and gas gap thicknesses are distributed
74 // according to exponent. The mean values of the plate and gas gap thicknesses
75 // are supposed to be about XTR formation zones but much less than
76 // mean absorption length of XTR photons in coresponding material.
77 
78 G4double
80  G4double gamma, G4double varAngle )
81 {
82  G4double result, Za, Zb, Ma, Mb ;
83 
84  Za = GetPlateFormationZone(energy,gamma,varAngle) ;
85  Zb = GetGasFormationZone(energy,gamma,varAngle) ;
86 
87  Ma = GetPlateLinearPhotoAbs(energy) ;
88  Mb = GetGasLinearPhotoAbs(energy) ;
89 
90 
92  G4complex Cb(1.0+0.5*fGasThick*Mb/fAlphaGas,fGasThick/Zb/fAlphaGas) ;
93 
94  G4complex Ha = std::pow(Ca,-fAlphaPlate) ;
95  G4complex Hb = std::pow(Cb,-fAlphaGas) ;
96  G4complex H = Ha*Hb ;
97 
98  G4complex F1 = (1.0 - Ha)*(1.0 - Hb )/(1.0 - H)
100 
101  G4complex F2 = (1.0-Ha)*(1.0-Ha)*Hb/(1.0-H)/(1.0-H)
102  * (1.0 - std::pow(H,fPlateNumber)) ;
103 
104  G4complex R = (F1 + F2)*OneInterfaceXTRdEdx(energy,gamma,varAngle) ;
105 
106  result = 2.0*std::real(R) ;
107 
108  return result ;
109 }
110 
111 
112 //
113 //
115 
116 
117 
118 
119 
120 
121 
122 
G4double G4ParticleHPJENDLHEData::G4double result
G4GammaXTRadiator(G4LogicalVolume *anEnvelope, G4double, G4double, G4Material *, G4Material *, G4double, G4double, G4int, const G4String &processName="XTRgammaRadiator")
G4double GetGasFormationZone(G4double, G4double, G4double)
G4double GetPlateLinearPhotoAbs(G4double)
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4double GetGasLinearPhotoAbs(G4double)
G4double GetStackFactor(G4double energy, G4double gamma, G4double varAngle) override
int G4int
Definition: G4Types.hh:78
G4double GetPlateFormationZone(G4double, G4double, G4double)
tuple b
Definition: test.py:12
std::complex< G4double > G4complex
Definition: G4Types.hh:81
G4GLOB_DLL std::ostream G4cout
const G4int n
G4double energy(const ThreeVector &p, const G4double m)
G4complex OneInterfaceXTRdEdx(G4double energy, G4double gamma, G4double varAngle)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76