Geant4_10
demo.py
Go to the documentation of this file.
1 #!/usr/bin/python
2 # ==================================================================
3 # python script for Geant4Py test
4 #
5 # gtest01
6 # - check basic control flow
7 # ==================================================================
8 from Geant4 import *
9 import demo_wp
10 import g4py.MedicalBeam
11 import ROOT
12 
13 # ==================================================================
14 # ROOT PART #
15 # ==================================================================
16 
17 # ------------------------------------------------------------------
18 def init_root():
19 # ------------------------------------------------------------------
20  ROOT.gROOT.Reset()
21 
22  # plot style
23  ROOT.gStyle.SetTextFont(42)
24  ROOT.gStyle.SetTitleFont(42, "X")
25  ROOT.gStyle.SetLabelFont(42, "X")
26  ROOT.gStyle.SetTitleFont(42, "Y")
27  ROOT.gStyle.SetLabelFont(42, "Y")
28 
29  global gCanvas
30  gCanvas= ROOT.TCanvas("water_phantom_plots",
31  "Water Phantom Demo Plots",
32  620, 30, 800, 800)
33 
34 # ------------------------------------------------------------------
35 def hini():
36 # ------------------------------------------------------------------
37  global gPad1
38  gPad1= ROOT.TPad("2D", "2D", 0.02, 0.5, 0.98, 1.)
39  gPad1.Draw()
40  gPad1.cd()
41 
42  ROOT.gStyle.SetPalette(1);
43 
44  global hist_dose2d
45  hist_dose2d= ROOT.TH2D("2D Dose", "Dose Distribution",
46  200, 0., 400.,
47  81, -81., 81.)
48  hist_dose2d.SetXTitle("Z (mm)")
49  hist_dose2d.SetYTitle("X (mm)")
50  hist_dose2d.SetStats(0)
51  hist_dose2d.Draw("colz")
52 
53  gCanvas.cd()
54  global gPad2
55  gPad2= ROOT.TPad("Z", "Z", 0.02, 0., 0.98, 0.5)
56  gPad2.Draw()
57  gPad2.cd()
58 
59  global hist_dosez
60  hist_dosez= ROOT.TH1D("Z Dose", "Depth Dose", 200, 0., 400.)
61  hist_dosez.SetXTitle("(mm)")
62  hist_dosez.SetYTitle("Accumulated Dose (MeV)")
63  hist_dosez.Draw()
64 
65 # ------------------------------------------------------------------
66 def hshow():
67 # ------------------------------------------------------------------
68  gPad1.cd()
69  hist_dose2d.Draw("colz")
70  gPad2.cd()
71  hist_dosez.Draw()
72 
73 
74 # ==================================================================
75 # Geant4 PART #
76 # ==================================================================
77 
78 # ==================================================================
79 # user actions in python
80 # ==================================================================
82  "My Primary Generator Action"
83 
84  def __init__(self):
85  G4VUserPrimaryGeneratorAction.__init__(self)
87 
88  def GeneratePrimaries(self, event):
89  self.particleGun.GeneratePrimaryVertex(event)
90 
91 # ------------------------------------------------------------------
93  "My Run Action"
94 
95  def EndOfRunAction(self, run):
96  print "*** End of Run"
97  print "- Run sammary : (id= %d, #events= %d)" \
98  % (run.GetRunID(), run.GetNumberOfEventToBeProcessed())
99 
100 # ------------------------------------------------------------------
102  "My Event Action"
103 
104  def EndOfEventAction(self, event):
105  gPad1.Modified()
106  gPad1.Update()
107  gPad2.Modified()
108  gPad2.Update()
109  ROOT.gSystem.ProcessEvents()
110 
111 # ------------------------------------------------------------------
113  "My Stepping Action"
114 
115  def UserSteppingAction(self, step):
116  pass
117 
118 # ------------------------------------------------------------------
119 class ScoreSD(G4VSensitiveDetector):
120  "SD for score voxels"
121 
122  def __init__(self):
123  G4VSensitiveDetector.__init__(self, "ScoreVoxel")
124 
125  def ProcessHits(self, step, rohist):
126  preStepPoint= step.GetPreStepPoint()
127  if(preStepPoint.GetCharge() == 0):
128  return
129 
130  track= step.GetTrack()
131  touchable= track.GetTouchable()
132  voxel_id= touchable.GetReplicaNumber()
133  dedx= step.GetTotalEnergyDeposit()
134  xz= posXZ(voxel_id)
135  hist_dose2d.Fill(xz[1], xz[0], dedx/MeV)
136  if( abs(xz[0]) <= 100 ):
137  hist_dosez.Fill(xz[1], dedx/MeV)
138 
139 # ------------------------------------------------------------------
140 def posXZ(copyN):
141  dd= 2.*mm
142  nx= 81
143 
144  iz= copyN/nx
145  ix= copyN-iz*nx-nx/2
146 
147  x0= ix*dd
148  z0= (iz+0.5)*dd
149  return (x0,z0)
150 
151 
152 # ==================================================================
153 # main
154 # ==================================================================
155 # init ROOT...
156 init_root()
157 hini()
158 
159 # configure application
160 #app= demo_wp.MyApplication()
161 #app.Configure()
162 
163 myMaterials= demo_wp.MyMaterials()
164 myMaterials.Construct()
165 
166 myDC= demo_wp.MyDetectorConstruction()
167 gRunManager.SetUserInitialization(myDC)
168 
169 myPL= FTFP_BERT()
170 gRunManager.SetUserInitialization(myPL)
171 
172 # set user actions...
174 gRunManager.SetUserAction(myPGA)
175 
176 myRA= MyRunAction()
177 gRunManager.SetUserAction(myRA)
178 
179 myEA= MyEventAction()
180 gRunManager.SetUserAction(myEA)
181 
182 #mySA= MySteppingAction()
183 #gRunManager.SetUserAction(mySA)
184 
185 # set particle gun
186 #pg= myPGA.particleGun
187 #pg.SetParticleByName("proton")
188 #pg.SetParticleEnergy(230.*MeV)
189 #pg.SetParticleMomentumDirection(G4ThreeVector(0., 0., 1.))
190 #pg.SetParticlePosition(G4ThreeVector(0.,0.,-50.)*cm)
191 
192 # medical beam
193 beam= g4py.MedicalBeam.Construct()
194 beam.particle= "proton"
195 beam.kineticE= 230.*MeV
196 #beam.particle= "gamma"
197 #beam.kineticE= 1.77*MeV
198 beam.sourcePosition= G4ThreeVector(0.,0.,-100.*cm)
199 beam.SSD= 100.*cm
200 beam.fieldXY= [5.*cm, 5.*cm]
201 
202 # initialize
203 gRunManager.Initialize()
204 
205 # set SD (A SD should be set after geometry construction)
206 scoreSD= ScoreSD()
207 myDC.SetSDtoScoreVoxel(scoreSD)
208 
209 # visualization
210 gApplyUICommand("/control/execute vis.mac")
211 
212 # beamOn
213 gRunManager.BeamOn(100)
214 
215 #ROOT.gSystem.Run()
216 
def hini
Definition: demo.py:35
def __init__
Definition: demo.py:122
def hshow
Definition: demo.py:66
def posXZ
Definition: demo.py:140
if(nlines<=0)
def EndOfRunAction
Definition: demo.py:95
gApplyUICommand
Definition: __init__.py:168
def ProcessHits
Definition: demo.py:125
def init_root
Definition: demo.py:18
def EndOfEventAction
Definition: demo.py:104
def UserSteppingAction
Definition: demo.py:115