| 1 | #!/usr/bin/python
|
|---|
| 2 | # ==================================================================
|
|---|
| 3 | # An example of ploting by EmCalculator
|
|---|
| 4 | #
|
|---|
| 5 | # Plotting photon cross sections and stopping power with ROOT
|
|---|
| 6 | # ==================================================================
|
|---|
| 7 | from Geant4 import *
|
|---|
| 8 | import g4py.NISTmaterials
|
|---|
| 9 | import g4py.ezgeom
|
|---|
| 10 |
|
|---|
| 11 | # ==================================================================
|
|---|
| 12 | # geometry setup
|
|---|
| 13 | # ==================================================================
|
|---|
| 14 |
|
|---|
| 15 | # ------------------------------------------------------------------
|
|---|
| 16 | # setup
|
|---|
| 17 | # ------------------------------------------------------------------
|
|---|
| 18 | def Configure():
|
|---|
| 19 | g4py.NISTmaterials.Construct()
|
|---|
| 20 | g4py.ezgeom.Construct()
|
|---|
| 21 |
|
|---|
| 22 | # ------------------------------------------------------------------
|
|---|
| 23 | # constructing geometry
|
|---|
| 24 | # ------------------------------------------------------------------
|
|---|
| 25 | def SetMaterial(material_name):
|
|---|
| 26 | material= gNistManager.FindOrBuildMaterial(material_name)
|
|---|
| 27 | g4py.ezgeom.SetWorldMaterial(material)
|
|---|
| 28 |
|
|---|
| 29 |
|
|---|
| 30 | # ==================================================================
|
|---|
| 31 | # plot by ROOT
|
|---|
| 32 | # ==================================================================
|
|---|
| 33 | import ROOT
|
|---|
| 34 | from math import log, log10, sqrt, ceil, floor
|
|---|
| 35 | from array import array
|
|---|
| 36 |
|
|---|
| 37 | # ------------------------------------------------------------------
|
|---|
| 38 | # caclculate plot range
|
|---|
| 39 | # ------------------------------------------------------------------
|
|---|
| 40 | def plot_range(xmin, xmax, xmargin=0.):
|
|---|
| 41 | xmaxlog= 10
|
|---|
| 42 | xminlog= -10
|
|---|
| 43 |
|
|---|
| 44 | if(xmax!=0.):
|
|---|
| 45 | xmaxlog= log10(xmax)
|
|---|
| 46 |
|
|---|
| 47 | if(xmin!=0):
|
|---|
| 48 | xminlog= log10(xmin)
|
|---|
| 49 |
|
|---|
| 50 | ixmaxlog= xmaxlog+0.5
|
|---|
| 51 | ixminlog= xminlog-0.5-xmargin
|
|---|
| 52 |
|
|---|
| 53 | return [10**ixminlog, 10**ixmaxlog]
|
|---|
| 54 |
|
|---|
| 55 | # ------------------------------------------------------------------
|
|---|
| 56 | # ROOT init
|
|---|
| 57 | # ------------------------------------------------------------------
|
|---|
| 58 | def init_root():
|
|---|
| 59 | ROOT.gROOT.Reset()
|
|---|
| 60 |
|
|---|
| 61 | # plot style
|
|---|
| 62 | ROOT.gStyle.SetTextFont(82)
|
|---|
| 63 |
|
|---|
| 64 | ROOT.gStyle.SetTitleFont(82, "X")
|
|---|
| 65 | ROOT.gStyle.SetTitleFontSize(0.04)
|
|---|
| 66 | ROOT.gStyle.SetLabelFont(82, "X")
|
|---|
| 67 | ROOT.gStyle.SetTitleFont(82, "Y")
|
|---|
| 68 | ROOT.gStyle.SetLabelFont(82, "Y")
|
|---|
| 69 |
|
|---|
| 70 | #ROOT.gStyle.SetOptTitle(0)
|
|---|
| 71 | ROOT.gStyle.SetErrorX(0)
|
|---|
| 72 |
|
|---|
| 73 | canvas= ROOT.TCanvas("g4plot", "g4plot", 620, 30, 600, 600)
|
|---|
| 74 |
|
|---|
| 75 | canvas.SetLogy()
|
|---|
| 76 | canvas.SetLogx()
|
|---|
| 77 | canvas.SetGrid()
|
|---|
| 78 |
|
|---|
| 79 | return canvas
|
|---|
| 80 |
|
|---|
| 81 | # ------------------------------------------------------------------
|
|---|
| 82 | # do a plot
|
|---|
| 83 | # ------------------------------------------------------------------
|
|---|
| 84 | def make_plot(xlist, user_title, axis_titile, q_super_impose=0):
|
|---|
| 85 |
|
|---|
| 86 | ekin_array, y_array = array('d'), array('d')
|
|---|
| 87 |
|
|---|
| 88 | for x in xlist:
|
|---|
| 89 | ekin_array.append(x[0])
|
|---|
| 90 | y_array.append(x[1])
|
|---|
| 91 |
|
|---|
| 92 | # plot range
|
|---|
| 93 | xmin= min(ekin_array)
|
|---|
| 94 | xmax= max(ekin_array)
|
|---|
| 95 | xrange= plot_range(xmin, xmax)
|
|---|
| 96 |
|
|---|
| 97 | ymin= min(y_array)
|
|---|
| 98 | ymax= max(y_array)
|
|---|
| 99 | yrange= plot_range(ymin, ymax, 2)
|
|---|
| 100 |
|
|---|
| 101 | if(q_super_impose==0):
|
|---|
| 102 | htit= user_title
|
|---|
| 103 | global frame
|
|---|
| 104 | frame= ROOT.TH1F("dumy", htit, 1, xrange[0], xrange[1]);
|
|---|
| 105 | frame.SetMinimum(yrange[0]);
|
|---|
| 106 | frame.SetMaximum(yrange[1]);
|
|---|
| 107 | frame.SetXTitle("Kinetic Energy (MeV)")
|
|---|
| 108 | frame.GetXaxis().SetLabelSize(0.025)
|
|---|
| 109 | frame.GetXaxis().SetTitleSize(0.03)
|
|---|
| 110 | frame.SetYTitle(axis_titile)
|
|---|
| 111 | frame.GetYaxis().SetLabelSize(0.025)
|
|---|
| 112 | frame.GetYaxis().SetTitleSize(0.03)
|
|---|
| 113 | frame.SetStats(0)
|
|---|
| 114 | frame.Draw()
|
|---|
| 115 |
|
|---|
| 116 | plot= ROOT.TGraph(len(ekin_array), ekin_array, y_array)
|
|---|
| 117 | plot.Draw("L")
|
|---|
| 118 | plot.SetLineColor(q_super_impose+1)
|
|---|
| 119 |
|
|---|
| 120 | return plot
|
|---|
| 121 |
|
|---|