2# ==================================================================
3# Python module (Python3)
13# ==================================================================
15# ==================================================================
16__all__ = [ 'CalculatePhotonCrossSection', 'CalculateDEDX' ]
18# ==================================================================
20# ==================================================================
21def CalculatePhotonCrossSection(mat, elist, verbose=0,
22 plist=["compt",
"",
"phot",
"conv"]):
24 Calculate photon cross section for a given material
and
25 a list of energy, returing a list of cross sections
for
26 the components of
"Copmton scattering",
"rayleigh scattering",
27 "photoelectric effect",
"pair creation" and total one.
32 verbose: verbose level [0]
33 plist: list of process name
34 (compton/rayleigh/photoelectic/conversion) [StandardEM set]
37 "compt": Compton Scattering
38 "rayleigh": Rayleigh Scattering
39 "phot" : photoelectric effect
40 "conv" : pair Creation
45 value= xsec_list[energy_index][
"compt"]
48 print("-------------------------------------------------------------------")
49 print(
" Photon Cross Section (", mat,
")")
50 print(
"Energy Compton Raleigh Photo- Pair Total")
51 print(
" Scattering Scattering electric Creation")
52 print(
"(MeV) (cm2/g) (cm2/g) (cm2/g) (cm2/g) (cm2/g)")
53 print(
"-------------------------------------------------------------------")
59 = gEmCalculator.ComputeCrossSectionPerVolume(ekin,
"gamma", plist[0],
62 = gEmCalculator.ComputeCrossSectionPerVolume(ekin,
"gamma", plist[1],
66 = gEmCalculator.ComputeCrossSectionPerVolume(ekin,
"gamma", plist[2],
69 = gEmCalculator.ComputeCrossSectionPerVolume(ekin,
"gamma", plist[3],
72 xsec[
"tot"]= xsec[
"compt"] + xsec[
"rayleigh"] + xsec[
"phot"] + xsec[
"conv"]
74 xsection_list.append((ekin, xsec))
77 print(
" %8.3e %8.3e %8.3e %8.3e %8.3e %8.3e" \
78 % (ekin/MeV, xsec[
"compt"]/(cm2/g), xsec[
"rayleigh"]/(cm2/g),
79 xsec[
"phot"]/(cm2/g), xsec[
"conv"]/(cm2/g), xsec[
"tot"]/(cm2/g)))
88 plist=["eIoni", "eBrem", "muIoni", "muBrems", "hIoni"]):
90 Calculate stopping powers for a give particle, material
and
91 a list of energy, returing stopping power
for the components of
92 "Ionization",
"Radiation" and total one.
98 verbose: verbose level [0]
99 plist: list of process name
100 (electron ionization/electron brems/
101 muon ionization/muon brems/hadron ionization) [StandardEM set]
105 "brems": Bremsstrahlung
110 value= dedx_list[energy_index][
"ioni"]
113 print("------------------------------------------------------")
114 print(
" Stopping Power (", part,
",", mat,
")")
115 print(
" Energy Ionization Radiation Total")
116 print(
" (MeV) (MeVcm2/g) (MeVcm2/g) (MeVcm2/g)")
117 print(
"------------------------------------------------------")
121 if ( part==
"e+" or part==
"e-" ):
122 procname_ioni= plist[0]
123 procname_brems= plist[1]
124 elif ( part==
"mu+" or part==
"mu-"):
125 procname_ioni= plist[2]
126 procname_brems= plist[3]
128 procname_ioni= plist[4]
135 = gEmCalculator.ComputeDEDX(ekin, part, procname_ioni, mat) * MeV*cm2/g
137 = gEmCalculator.ComputeDEDX(ekin, part, procname_brems, mat) * MeV*cm2/g
138 dedx[
"tot"]= dedx[
"ioni"]+ dedx[
"brems"]
141 print(
" %8.3e %8.3e %8.3e %8.3e" \
142 % (ekin/MeV, dedx[
"ioni"]/(MeV*cm2/g),
143 dedx[
"brems"]/(MeV*cm2/g), dedx[
"tot"]/(MeV*cm2/g) ))
146 dedx_list.append((ekin, dedx))
void print(G4double elem)
def CalculateDEDX(part, mat, elist, verbose=0, plist=["eIoni", "eBrem", "muIoni", "muBrems", "hIoni"])
def CalculatePhotonCrossSection(mat, elist, verbose=0, plist=["compt", "", "phot", "conv"])