Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4XTRGammaRadModel.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#include "G4XTRGammaRadModel.hh"
28
29////////////////////////////////////////////////////////////////////////////
30// Constructor, destructor
32 G4double alphaPlate, G4double alphaGas,
33 G4Material* foilMat, G4Material* gasMat,
34 G4double a, G4double b, G4int n,
35 const G4String& processName)
36 : G4VXTRenergyLoss(anEnvelope, foilMat, gasMat, a, b, n, processName)
37{
38 G4cout << "Gamma distributed X-ray TR radiator model is called" << G4endl;
39
40 // Build energy and angular integral spectra of X-ray TR photons from
41 // a radiator
42 fAlphaPlate = alphaPlate;
43 fAlphaGas = alphaGas;
44 G4cout << "fAlphaPlate = " << fAlphaPlate << " ; fAlphaGas = " << fAlphaGas
45 << G4endl;
46 fExitFlux = true;
47}
48
49///////////////////////////////////////////////////////////////////////////
51
52void G4XTRGammaRadModel::ProcessDescription(std::ostream& out) const
53{
54 out << "Rough model describing X-ray transition radiation. Thicknesses of "
55 "plates\n"
56 "and gas gaps are distributed according to gamma distributions.\n";
57}
58
59///////////////////////////////////////////////////////////////////////////
60// Rough approximation for radiator interference factor for the case of
61// fully GamDistr radiator. The plate and gas gap thicknesses are distributed
62// according to exponent. The mean values of the plate and gas gap thicknesses
63// are supposed to be about XTR formation zones but much less than
64// mean absorption length of XTR photons in coresponding material.
66 G4double varAngle)
67{
68 G4double result, Qa, Qb, Q, Za, Zb, Ma, Mb;
69
70 Za = GetPlateFormationZone(energy, gamma, varAngle);
71 Zb = GetGasFormationZone(energy, gamma, varAngle);
72
73 Ma = GetPlateLinearPhotoAbs(energy);
74 Mb = GetGasLinearPhotoAbs(energy);
75
76 Qa = (1.0 + fPlateThick * Ma / fAlphaPlate);
77 Qa = std::pow(Qa, -fAlphaPlate);
78 Qb = (1.0 + fGasThick * Mb / fAlphaGas);
79 Qb = std::pow(Qb, -fAlphaGas);
80 Q = Qa * Qb;
81
82 G4complex Ca(1.0 + 0.5 * fPlateThick * Ma / fAlphaPlate,
84 G4complex Cb(1.0 + 0.5 * fGasThick * Mb / fAlphaGas,
85 fGasThick / Zb / fAlphaGas);
86
87 G4complex Ha = std::pow(Ca, -fAlphaPlate);
88 G4complex Hb = std::pow(Cb, -fAlphaGas);
89 G4complex H = Ha * Hb;
90
91 G4complex F1 = (0.5 * (1 + Qa) * (1.0 + H) - Ha - Qa * Hb) / (1.0 - H);
92
93 G4complex F2 = (1.0 - Ha) * (Qa - Ha) * Hb / (1.0 - H) / (Q - H);
94
95 F2 *= std::pow(Q, G4double(fPlateNumber)) - std::pow(H, fPlateNumber);
96
97 result = (1. - std::pow(Q, G4double(fPlateNumber))) / (1. - Q);
98
99 G4complex stack = result * F1;
100 stack += F2;
101 stack *= 2.0 * OneInterfaceXTRdEdx(energy, gamma, varAngle);
102
103 result = std::real(stack);
104
105 return result;
106}
double G4double
Definition: G4Types.hh:83
std::complex< G4double > G4complex
Definition: G4Types.hh:88
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
G4double GetPlateLinearPhotoAbs(G4double)
G4double GetGasFormationZone(G4double, G4double, G4double)
G4complex OneInterfaceXTRdEdx(G4double energy, G4double gamma, G4double varAngle)
G4double GetPlateFormationZone(G4double, G4double, G4double)
G4double GetGasLinearPhotoAbs(G4double)
void ProcessDescription(std::ostream &) const override
G4XTRGammaRadModel(G4LogicalVolume *anEnvelope, G4double, G4double, G4Material *, G4Material *, G4double, G4double, G4int, const G4String &processName="XTRgammaRadiator")
G4double GetStackFactor(G4double energy, G4double gamma, G4double varAngle) override