Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4StrawTubeXTRadiator.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
28
29#include "G4Gamma.hh"
31#include "G4SystemOfUnits.hh"
32
33////////////////////////////////////////////////////////////////////////////
34// Constructor, destructor
36 G4Material* foilMat,
37 G4Material* gasMat, G4double a,
38 G4double b, G4Material* mediumMat,
39 G4bool unishut,
40 const G4String& processName)
41 : G4VXTRenergyLoss(anEnvelope, foilMat, gasMat, a, b, 1, processName)
42{
43 if(verboseLevel > 0)
44 G4cout << "Straw tube X-ray TR radiator EM process is called" << G4endl;
45
46 if(unishut)
47 {
48 fAlphaPlate = 1. / 3.;
49 fAlphaGas = 12.4;
50 if(verboseLevel > 0)
51 G4cout << "straw uniform shooting: "
52 << "fAlphaPlate = " << fAlphaPlate
53 << " ; fAlphaGas = " << fAlphaGas << G4endl;
54 }
55 else
56 {
57 fAlphaPlate = 0.5;
58 fAlphaGas = 5.;
59 if(verboseLevel > 0)
60 G4cout << "straw isotropical shooting: "
61 << "fAlphaPlate = " << fAlphaPlate
62 << " ; fAlphaGas = " << fAlphaGas << G4endl;
63 }
64
65 // index of medium material
66 fMatIndex3 = (G4int)mediumMat->GetIndex();
67 if(verboseLevel > 0)
68 G4cout << "medium material = " << mediumMat->GetName() << G4endl;
69
70 // plasma energy squared for plate material
71 fSigma3 = fPlasmaCof * mediumMat->GetElectronDensity();
72 if(verboseLevel > 0)
73 G4cout << "medium plasma energy = " << std::sqrt(fSigma3) / eV << " eV"
74 << G4endl;
75
76 // Compute cofs for preparation of linear photo absorption in external medium
78}
79
80///////////////////////////////////////////////////////////////////////////
82
83void G4StrawTubeXTRadiator::ProcessDescription(std::ostream& out) const
84{
85 out << "Simulation of forward X-ray transition radiation for the case of\n"
86 "a straw tube radiator.\n";
87}
88
89///////////////////////////////////////////////////////////////////////////
90// Approximation for radiator interference factor for the case of
91// straw tube radiator. The plate (window, straw wall) and gas (inside straw)
92// gap thicknesses are gamma distributed.
93// The mean values of the plate and gas gap thicknesses
94// are supposed to be about XTR formation zone.
96 G4double varAngle)
97{
98 G4double result, L2, L3, M2, M3;
99
100 L2 = GetPlateFormationZone(energy, gamma, varAngle);
101 L3 = GetGasFormationZone(energy, gamma, varAngle);
102
103 M2 = GetPlateLinearPhotoAbs(energy);
104 M3 = GetGasLinearPhotoAbs(energy);
105
106 G4complex C2(1.0 + 0.5 * fPlateThick * M2 / fAlphaPlate,
107 fPlateThick / L2 / fAlphaPlate);
108 G4complex C3(1.0 + 0.5 * fGasThick * M3 / fAlphaGas,
109 fGasThick / L3 / fAlphaGas);
110
111 G4complex H2 = std::pow(C2, -fAlphaPlate);
112 G4complex H3 = std::pow(C3, -fAlphaGas);
113 G4complex H = H2 * H3;
114
115 G4complex Z1 = GetMediumComplexFZ(energy, gamma, varAngle);
116 G4complex Z2 = GetPlateComplexFZ(energy, gamma, varAngle);
117 G4complex Z3 = GetGasComplexFZ(energy, gamma, varAngle);
118
119 G4complex R = (Z1 - Z2) * (Z1 - Z2) * (1. - H2 * H) +
120 (Z2 - Z3) * (Z2 - Z3) * (1. - H3) +
121 2. * (Z1 - Z2) * (Z2 - Z3) * H2 * (1. - H3);
122
123 result = 2.0 * std::real(R) * (varAngle * energy / hbarc / hbarc);
124
125 return result;
126}
127
128////////////////////////////////////////////////////////////////////////
129// Calculates formation zone for external medium. Omega is energy !!!
131 G4double gamma,
132 G4double varAngle)
133{
134 G4double cof, lambda;
135 lambda = 1.0 / gamma / gamma + varAngle + fSigma3 / omega / omega;
136 cof = 2.0 * hbarc / omega / lambda;
137 return cof;
138}
139
140////////////////////////////////////////////////////////////////////////
141// Calculates complex formation zone for external medium. Omega is energy !!!
143 G4double gamma,
144 G4double varAngle)
145{
146 G4double cof, length, delta, real_v, image_v;
147
148 length = 0.5 * GetMediumFormationZone(omega, gamma, varAngle);
149 delta = length * GetMediumLinearPhotoAbs(omega);
150 cof = 1.0 / (1.0 + delta * delta);
151
152 real_v = length * cof;
153 image_v = real_v * delta;
154
155 G4complex zone(real_v, image_v);
156 return zone;
157}
158
159////////////////////////////////////////////////////////////////////////
160// Computes matrix of Sandia photo absorption cross section coefficients for
161// medium material
163{
164 const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
165 const G4Material* mat = (*theMaterialTable)[fMatIndex3];
167}
168
169//////////////////////////////////////////////////////////////////////
170// Returns the value of linear photo absorption coefficient (in reciprocal
171// length) for medium for given energy of X-ray photon omega
173{
174 G4double omega2, omega3, omega4;
175
176 omega2 = omega * omega;
177 omega3 = omega2 * omega;
178 omega4 = omega2 * omega2;
179
180 const G4double* SandiaCof =
182
183 G4double cross = SandiaCof[0] / omega + SandiaCof[1] / omega2 +
184 SandiaCof[2] / omega3 + SandiaCof[3] / omega4;
185 return cross;
186}
std::vector< G4Material * > G4MaterialTable
double G4double
Definition G4Types.hh:83
bool G4bool
Definition G4Types.hh:86
std::complex< G4double > G4complex
Definition G4Types.hh:88
int G4int
Definition G4Types.hh:85
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
const double C2
#define C3
G4SandiaTable * GetSandiaTable() const
G4double GetElectronDensity() const
std::size_t GetIndex() const
static G4MaterialTable * GetMaterialTable()
const G4String & GetName() const
G4double GetSandiaCofForMaterial(G4int, G4int) const
G4double GetStackFactor(G4double energy, G4double gamma, G4double varAngle) override
void ProcessDescription(std::ostream &) const override
G4StrawTubeXTRadiator(G4LogicalVolume *anEnvelope, G4Material *, G4Material *, G4double, G4double, G4Material *, G4bool unishut=false, const G4String &processName="StrawTubeXTRadiator")
G4double GetMediumFormationZone(G4double, G4double, G4double)
G4double GetMediumLinearPhotoAbs(G4double)
G4complex GetMediumComplexFZ(G4double, G4double, G4double)
G4int verboseLevel
G4double GetPlateLinearPhotoAbs(G4double)
G4double GetGasFormationZone(G4double, G4double, G4double)
G4complex GetPlateComplexFZ(G4double, G4double, G4double)
static constexpr G4double fPlasmaCof
G4double GetPlateFormationZone(G4double, G4double, G4double)
G4double GetGasLinearPhotoAbs(G4double)
G4complex GetGasComplexFZ(G4double, G4double, G4double)