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