| Geant4 Cross Reference |
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 #include "G4XTRRegularRadModel.hh"
28
29 #include "G4PhysicalConstants.hh"
30
31 ////////////////////////////////////////////////////////////////////////////
32 // Constructor, destructor
33 G4XTRRegularRadModel::G4XTRRegularRadModel(G4LogicalVolume* anEnvelope,
34 G4Material* foilMat,
35 G4Material* gasMat, G4double a,
36 G4double b, G4int n,
37 const G4String& processName)
38 : G4VXTRenergyLoss(anEnvelope, foilMat, gasMat, a, b, n, processName)
39 {
40 G4cout << " XTR Regular discrete radiator model is called" << G4endl;
41
42 fExitFlux = true;
43 }
44
45 ///////////////////////////////////////////////////////////////////////////
46 G4XTRRegularRadModel::~G4XTRRegularRadModel() = default;
47
48 ///////////////////////////////////////////////////////////////////////////
49 void G4XTRRegularRadModel::ProcessDescription(std::ostream& out) const
50 {
51 out << "Describes X-ray transition radiation with thickness of gaps and "
52 "plates\n"
53 "fixed.\n";
54 }
55
56 ///////////////////////////////////////////////////////////////////////////
57 G4double G4XTRRegularRadModel::SpectralXTRdEdx(G4double energy)
58 {
59 static constexpr G4double cofPHC = 4. * pi * hbarc;
60 G4double result, sum = 0., tmp, cof1, cof2, cofMin, theta2, theta2k;
61 G4double aMa, bMb, sigma, dump;
62 G4int k, kMax, kMin;
63
64 aMa = fPlateThick * GetPlateLinearPhotoAbs(energy);
65 bMb = fGasThick * GetGasLinearPhotoAbs(energy);
66 sigma = 0.5 * (aMa + bMb);
67 dump = std::exp(-fPlateNumber * sigma);
68 if(verboseLevel > 2)
69 G4cout << " dump = " << dump << G4endl;
70 tmp = (fSigma1 - fSigma2) / cofPHC / energy;
71 cof1 = fPlateThick * tmp;
72 cof2 = fGasThick * tmp;
73
74 cofMin = energy * (fPlateThick + fGasThick) / fGamma / fGamma;
75 cofMin += (fPlateThick * fSigma1 + fGasThick * fSigma2) / energy;
76 cofMin /= cofPHC;
77
78 theta2 = cofPHC / (energy * (fPlateThick + fGasThick));
79
80 kMin = G4int(cofMin);
81 if(cofMin > kMin)
82 kMin++;
83
84 kMax = kMin + 49;
85
86 if(verboseLevel > 2)
87 {
88 G4cout << cof1 << " " << cof2 << " " << cofMin << G4endl;
89 G4cout << "kMin = " << kMin << "; kMax = " << kMax << G4endl;
90 }
91 for(k = kMin; k <= kMax; ++k)
92 {
93 tmp = pi * fPlateThick * (k + cof2) / (fPlateThick + fGasThick);
94 result = (k - cof1) * (k - cof1) * (k + cof2) * (k + cof2);
95 if(k == kMin && kMin == G4int(cofMin))
96 {
97 sum +=
98 0.5 * std::sin(tmp) * std::sin(tmp) * std::abs(k - cofMin) / result;
99 }
100 else
101 {
102 sum += std::sin(tmp) * std::sin(tmp) * std::abs(k - cofMin) / result;
103 }
104 theta2k = std::sqrt(theta2 * std::abs(k - cofMin));
105
106 if(verboseLevel > 2)
107 {
108 G4cout << k << " " << theta2k << " "
109 << std::sin(tmp) * std::sin(tmp) * std::abs(k - cofMin) / result
110 << " " << sum << G4endl;
111 }
112 }
113 result = 2 * (cof1 + cof2) * (cof1 + cof2) * sum / energy;
114 result *= dump * (-1 + dump + 2 * fPlateNumber);
115
116 return result;
117 }
118
119 ///////////////////////////////////////////////////////////////////////////
120 // Approximation for radiator interference factor for the case of
121 // fully Regular radiator. The plate and gas gap thicknesses are fixed.
122 // The mean values of the plate and gas gap thicknesses
123 // are supposed to be about XTR formation zones but much less than
124 // mean absorption length of XTR photons in corresponding material.
125 G4double G4XTRRegularRadModel::GetStackFactor(G4double energy, G4double gamma,
126 G4double varAngle)
127 {
128 G4double aZa = fPlateThick / GetPlateFormationZone(energy, gamma, varAngle);
129 G4double bZb = fGasThick / GetGasFormationZone(energy, gamma, varAngle);
130
131 G4double aMa = fPlateThick * GetPlateLinearPhotoAbs(energy);
132 G4double bMb = fGasThick * GetGasLinearPhotoAbs(energy);
133
134 G4double Qa = std::exp(-aMa);
135 G4double Qb = std::exp(-bMb);
136 G4double Q = Qa * Qb;
137
138 G4complex Ha(std::exp(-0.5 * aMa) * std::cos(aZa),
139 -std::exp(-0.5 * aMa) * std::sin(aZa));
140
141 G4complex Hb(std::exp(-0.5 * bMb) * std::cos(bZb),
142 -std::exp(-0.5 * bMb) * std::sin(bZb));
143
144 G4complex H = Ha * Hb;
145 G4complex Hs = std::conj(H);
146
147 G4complex F2 = (1.0 - Ha) * (Qa - Ha) * Hb * (1.0 - Hs) * (Q - Hs);
148 F2 *= std::pow(Q, G4double(fPlateNumber)) - std::pow(H, fPlateNumber);
149
150 G4double result = (1. - std::pow(Q, G4double(fPlateNumber))) / (1. - Q);
151 result *= (1. - Qa) * (1. + Qa - 2. * std::sqrt(Qa) * std::cos(aZa));
152 result /= (1. - std::sqrt(Q)) * (1. - std::sqrt(Q)) +
153 4. * std::sqrt(Q) * std::sin(0.5 * (aZa + bZb)) *
154 std::sin(0.5 * (aZa + bZb));
155
156 G4double I2 = 1.;
157 I2 /= (1. - std::sqrt(Q)) * (1. - std::sqrt(Q)) +
158 4. * std::sqrt(Q) * std::sin(0.5 * (aZa + bZb)) *
159 std::sin(0.5 * (aZa + bZb));
160
161 I2 /= Q * ((std::sqrt(Q) - std::cos(aZa + bZb)) *
162 (std::sqrt(Q) - std::cos(aZa + bZb)) +
163 std::sin(aZa + bZb) * std::sin(aZa + bZb));
164
165 G4complex stack = 2. * I2 * F2;
166 stack += result;
167 stack *= OneInterfaceXTRdEdx(energy, gamma, varAngle);
168
169 return std::real(stack);
170 }
171