330 likes | 374 Views
Learn Recommendation Settings and Utilities using PHITS code system with examples and input files for different calculation conditions.
E N D
PHITS Multi-Purpose Particle and Heavy Ion Transport code System Exercises using Recommendation Settings and Utilities Jun. 2013 revised title 1
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 2
What is Recommendation Settings? FAQ I cannot understand the appropriate parameter setting even if I read PHITS manual carefully We prepared several examples of PHITS input files for various calculation conditions, and uploaded to PHITS website Our Reply http://phits.jaea.go.jp/examples.html Recommendation Settings 3
List of Recommendation Settings • DetectorResponse.inp Detector response calculation (event-by-event information) • Shielding.inp Shielding calculation using [t-track], example of complex geometry • ParticleTherapy.inp Dose & LET distribution calculation for charged particle therapy • PhotonTherapy.inp Dose & particle fluence calculation for X-ray therapy • SemiConductor.inp Dose calculation inside a tiny region for SER estimation • NuclearReaction.inp Double differential cross sections by calculating fluence of secondaries • H10multiplier.inp H*(10) in water using [t-track] in combination with [multiplier] • Counter.inp Dose deposited from primary and secondary particles using [counter] Recommendation Settings 4
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 5
ParticleTherapy.inp Water Phantom 10 cm Carbon 200 MeV/u 10 cm Tally • [t-deposit] (file=dose) Depth-dose distribution in water • [t-deposit] (file=dose-equivalent.out) Depth-dose-equivalent distribution in water • [t-LET] Probability density of LET in water • [t-SED] Probability density of lineal energy in water ParticleTherapy 6
Dose and Dose Equivalent Dose.eps Dose-equivalent.eps [ T - Deposit ] dedxfnc = 1 [ T - Deposit ] mesh = r-z r-type = 2 rmin = 0.000000 rmax = 1.000000 nr = 1 z-type = 2 zmin = 0.000000 zmax = 10.00000 nz = 100 dedxfnc = 0 R Output quantities are weighted by the function written in usrdfn1.f Z Q(L) relationship is given as default r-z mesh Dose is converted into Dose equivalent! ParticleTherapy 7
Probability Density of LET & y [T-LET] [T-SED] LET-distribution.eps (page 1: front surface) y-distribution.eps (page 1: front surface) LET of C(200 MeV/u) = 16 keV/mm y*f(y) broad distribution even for mono-energetic incidence L*f(L) has Sharp peak ParticleTherapy 8
Change Incident Particle Y Z Water Phantom X 10 cm Carbon 200 MeV/u Proton (0,0,-20) 10 cm [ S o u r c e ] s-type = 1 proj = 12C e0 = 200.00 r0 = 0.0000 x0 = 0.0000 y0 = 0.0000 z0 = -20.000 z1 = -20.000 dir = 1.0000 [ T - Deposit ] … [ T - Deposit ] … [ T – L E T ] ... [ T – S E D ] ... [ T - Deposit ] … [ T - Deposit ] off … [ T – L E T ] off ... [ T – S E D ] off ... Proton Execute Dose.eps 10cm water is too short to stop 200 MeV proton! ParticleTherapy 9
Change Geometry Y Z Water Phantom 10 cm X 10 cm Carbon 200 MeV/u Proton (0,0,-20) 30 cm 10 cm [ C e l l ] 1 1 -1.0000000E+00 1 -2 -3 2 0 #1 -999 3 -1 999 [ S u r f a c e ] 1 pz 0.0000000E+00 2 pz 1.0000000E+01 3 cz 5.0000000E+00 999 so 1.0000000E+02 Execute 3.0000000E+01 Forgot to change tally region!!! Dose.eps General description 10
Change Tally Region Y Z Water Phantom 10 cm X 10 cm Carbon 200 MeV/u Proton (0,0,-20) 30 cm [ T - Deposit ] title = depth-dose d mesh = r-z x0 = 0.000000 y0 = 0.000000 r-type = 2 rmin = 0.000000 rmax = 1.000000 nr = 1 z-type = 2 zmin = 0.000000 zmax = 10.00000 nz = 100 [ C e l l ] 1 1 -1.0000000E+00 1 -2 -3 2 0 #1 -999 3 -1 999 [ S u r f a c e ] 1 pz 0.0000000E+00 2 pz 2.0000000E+01 3 cz 5.0000000E+00 999 so 1.0000000E+02 set: c1[30.0] Statistics is not enough!! c1 c1 Dose.eps ParticleTherapy 11
Restart Calculation Y Z Water Phantom 10 cm X 10 cm Carbon 200 MeV/u Proton (0,0,-20) 30 cm [ P a r a m e t e r s ] icntl = 0 maxcas = 100 maxbch = 2 Execute 10 istdev = -1 Finally, we obtained the depth-dose distribution in water irradiated by 200 MeV proton with enough statistics!! Dose.eps ParticleTherapy 12
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 13
H10multiplier.inp Air 100 cm Concrete 400 MeV Neutron 100 cm 100 cm 3 ways to calculate “Dose” in PHITS • [t-track] Calculate dose from fluences of particle multiplied with dose conversion coefficients. In this case, dose conversion coefficients for H*(10) are defined in [multiplier] section • [t-heat] Calculate dose using Kerma Approximation • [t-deposit] Calculate dose from the ionization energy of charged particles H10multiplier 14
Results of Each Tally [T-Deposit] [T-Track] [T-Heat] • [t-track] Fluence multiplied with H*(10) dose conversion coefficients → No gap between concrete and air • [t-heat] Kerma approximation for low-energy neutrons and photons →There is a gap between concrete and air • [t-deposit] Ionizing energy of charged particles → There is a gap between concrete and air. Large uncertainties in air H10multiplier 15
Calculation of H*(10) and Effective Dose Fluence calculated by [t-track] can be multiplied with several functions that are pre-defined or defined in [multiplier] section [ Multiplier ] number = -250 interpolation = log ne = 140 1.13000E-08 9.14000E+00*3600*1.0e-6 1.42000E-08 9.44000E+00*3600*1.0e-6 1.79000E-08 9.83000E+00*3600*1.0e-6 ... [ T - T r a c k ] title = H*(10) in xyz mesh in uSv/h ... multiplier = all part = neutron emax = 1000.0 mat mset1mset2 all ( 1.0 -250 ) (1.0 -102) multiplier = all part = photon emax = 1000.0 mat mset1mset2 all ( 1.0 -251) (1.0 -114) H*(10) (Page 1) Effective Dose (Page 2) Track.eps They seem to be almost the same! • H*(10) conversion coefficient for neutron(-250) and photon (-251) defined in [multiplier] • Predefined Effective dose conversion coefficients for neutron (-102) and photon (-114) H10multiplier 16
Change Axis [ T - T r a c k ] title = H*(10) in xy part = ( neutron photon ) mesh = xyz x-type = 2 xmin = -60.0 xmax = 60.0 nx = 60 y-type = 2 ymin = -60.0 ymax = 60.0 ny = 1 z-type = 2 zmin = 0.0 zmax = 120.0 nz = 60 e-type = 1 ne = 1 1.00000E-10 1.00000E+03 unit = 1 material = all axis = xz file = track.out Execute 1 H*(10) (Page 1) Effective Dose (Page 2) Avoid to output too much graphs Track.eps Tough to compare because each scale is automatically adjusted 2D-plot (contour map) to 1D-plot (histogram) z H10multiplier 17
Adjust Scale [ T - T r a c k ] title = H*(10) in xy part = ( neutron photon ) mesh = xyz x-type = 2 xmin = -60.0 xmax = 60.0 nx = 60 y-type = 2 ymin = -60.0 ymax = 60.0 ny = 1 z-type = 2 zmin = 0.0 zmax = 120.0 nz = 60 e-type = 1 ne = 1 1.00000E-10 1.00000E+03 unit = 1 material = all axis = xz file = track.out Execute 1 H*(10) (Page 1) Effective Dose (Page 2) Track.eps H*(10) < Effective dose z Add ANGEL parameter (min & max for y-axis) angel = ymin(5.e-05) ymax(5.e-4) H10multiplier 18
Execute Only ANGEL Copy “Track.out” to “Track2.out”, and Edit Right click “Track2.out” → Sendto → ANGEL #------------------------------------------------ 80行目 #newpage: # no. = 1 ie = 1 ix = 1 iy = 1 … x: z[cm] y: Flux [1/cm^2/source] p: xlin ylog afac(0.8) form(0.9) p: ymin(5.e-05) ymax(5.0e-4) h: n x y(p1-group),hh0l n # z-lower z-upper flux r.err 0.0000E+00 1.0000E+00 1.2809E-04 0.0263 … #------------------------------------------------- 232行目 newpage: # no. = 2 ie = 1 ix = 1 iy = 1 … x: z[cm] y: Flux [1/cm^2/source] p: xlin ylog afac(0.8) form(0.9) p: ymin(5.e-05) ymax(5.0e-4) h: n x y(p1-group),hh0l n # z-lower z-upper flux r.err 0.0000E+00 1.0000E+00 1.9351E-04 0.0174 ... H10 # Comment out (Same format as PHITS input file) Track2.eps Change Legend Do not use “(“ and “)” Specify line color (r: red, b: blue, g: greenetc) Do not insert space between “hh0l” and color ID r Effective H10multiplier 19
Contributions from High- and Low-Energy Particles [ T - T r a c k ] title = H*(10) in xy part = ( neutron photon ) mesh = xyz x-type = 2 xmin = -60.0 xmax = 60.0 nx = 1 y-type = 2 ymin = -60.0 ymax = 60.0 ny = 1 z-type = 2 zmin = 0.0 zmax = 120.0 nz = 60 e-type = 1 ne = 1 1.00000E-10 1.00000E+03 unit = 1 material = all axis = z file = track.out Low Energy (page1) Low Energy (page3) ≅ Execute High Energy (page2) High Energy (page4) < 2 20.01.0E+03 Divide the energy bin into high- and low-energy region angel = ymin(1.e-05) ymax(2.e-4) H*(10) Effective Dose H10multiplier 20
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 21
Contents of Utilities • Animation Create an animation of particle trajectories • Rotate3dshow Rotate the geometry depicted by [t-3dshow] • SimpleGEO Instruction for how to use SimpleGEO for PHITS input generator • Autorun Shell script for successively executing PHITS by slightly changing calculation conditions SimpleGEO Rotate3dshow Animation Utilities 22
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 23
Required Software ImageMagick (http://www.imagemagick.org/script/binary-releases.php) Software that can convert a multiple-page EPS file into a GIF animation • Procedures How to Create Animation 1. Execute PHITS Calculate the time-dependence of particle fluences or deposition energies by introducing “t-type” in PHITS input file 2. Edit EPS file generated by PHITS 1. Open EPS file “track.eps” with your text editor 2. Search the word “PageBoundingBox” twice 3. Change the first 2 numbers to 0, and save the file 3. Convert the EPS file to GIF animation using ImageMagick 1.Open “command prompt” 2. Move to the folder including the EPS file 3. Type ‘convert -dispose background -rotate 90 XXX.eps XXX.gif‘ %%PageBoundingBox: 27 36 571 803 %%PageBoundingBox: 0 0 571 803 Animation 24
Increase the Time Resolution animation.inp [ T - T r a c k ] C -- Contour figure Tally -- mesh = xyz ... t-type = 2 nt = 20 # Number of frame tmin = 0.00 # Initial time (nsec) tmax = 40 # Final time (nsec) ... angel = cmin(1.e-05) cmax(1.e+00) epsout = 1 60 20 frame Increase the frame number Execute PHITS Edit EPS file Convert EPS to GIF 60 frame Animation 25
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 26
How to use SimpleGEO for PHITS • Required Software SimpleGEO&Its plugin package (only for Windows) Download site: http://theis.web.cern.ch/theis/simplegeo/ • Procedures 1. Make geometry using SimpleGEO You have to define void (or air) region too. See manual of SimpleGEO 2. Export the geometry to PHITS readable format (File → Export → PHITS) Only [cell] & [surface] sections are generated 3. Make PHITS input file except for [cell] & [surface] sections Include the exported file using “infl:” command 4. Execute PHITS using the input file 5. Visualize the tally results obtained from “mesh=xyz” in SimpleGEO Load “DaVis3D” in SimpleGEO (Macros → Load Plugins → DaVis3D) Select the xyz-mesh tally results, and visualize in SimpleGEO SimpleGEO 27
Change Geometry in SimpleGEO doll.pht [ C e l l ] c Body 00001 26 -1 -1 c Eyes 00002 11 -7.87 -6 : -3 c Head 00003 6 -1.9 -2 #2 c Void 00004 0 -4 +1 #2 #3 c Outervoid 00005 -1 -5 +4 [ C e l l ] c Body 00001 26 -1 -1 c Eyes 00002 11 -7.87 -6 : -3 c Head 00003 6 -1.9 -2 #2 c Void 00004 0 -4 +1 #2 #3 +7 +8 c Outervoid 00005 -1 -5 +4 c leg1 00006 26 -1 -7 c leg2 00007 26 -1 -8 [ S u r f a c e ] c Body 1 RCC 0.00 0.00 0.00 0.00 0.00 10.00 5.00 c Headball 2 SPH 0.00 0.00 15.00 5.00 c LeftEye 3 SPH -3.00 4.00 16.00 1.00 c OuterSphere 4 SPH 0.00 0.00 0.00 100.00 c Outmostsphere 5 SPH 0.00 0.00 0.00 100.00 c RightEye 6 SPH 3.00 4.00 16.00 1.00 [ S u r f a c e ] c Body 1 RCC 0.00 0.00 0.00 0.00 0.00 10.00 5.00 c Headball 2 SPH 0.00 0.00 15.00 5.00 c LeftEye 3 SPH -3.00 4.00 16.00 1.00 c OuterSphere 4 SPH 0.00 0.00 0.00 100.00 c Outmostsphere 5 SPH 0.00 0.00 0.00 100.00 c RightEye 6 SPH 3.00 4.00 16.00 1.00 c leg1 7 RCC -2.50 0.00 -10.00 0.00 0.00 10.00 2.00 c leg2 8 RCC 2.50 0.00 -10.00 0.00 0.00 10.00 2.00 Export to PHITS SimpleGEO 28
Visualize Tally Result in 3D SimpleGEO.inp [ T - Deposit ] title = [t-deposit] in xyz mesh mesh = xyz # mesh type is xyz scoring mesh x-type = 2 # x-mesh is linear given by xmin, xmax and nx xmin = -5.000000 # minimum value of x-mesh points xmax = 5.000000 # maximum value of x-mesh points nx = 20 # number of x-mesh points y-type = 2 # y-mesh is linear given by ymin, ymax and ny ymin = -5.000000 # minimum value of y-mesh points ymax = 5.000000 # maximum value of y-mesh points ny = 20 # number of y-mesh points z-type = 2 # z-mesh is linear given by zmin, zmax and nz zmin = 0.000000 # minimum value of z-mesh points zmax = 20.00000 # maximum value of z-mesh points nz = 40 # number of z-mesh points -10.00000 Execute Extend the tally region SimpleGEO • Select ‘Macros -> Load Plugins -> DaVis3D’ • Press ‘DaVis3D‘ button • Select ‘deposity-xy.out’ and press ‘Load data’ SimpleGEO 29
Recommendation Settings Utilities Summary and Homework Contents • ParticleTherapy • H10mupliplier • Animation of Particle Trajectories • SimpleGEO Contents 30
It is difficult to make PHITS input file all by yourself It is better to select an appropriate recommendation setting, and edit the file for your calculation condition You can enjoy PHITS more using the utilities! Summary Dose distribution in PHITS-shaped water phantom irradiated by 1 GeV protons http://phits.jaea.go.jp 31
Homework Calculate the depth-dose distributions for various radial distances inside cylindrical water phantom irradiated by 11B 250 MeV/u beam Separate the depth-dose distributions into the contributions from neutron, proton, He, Li, Be, B Divide the phantom into 2 layers, composed of water and Aluminum (2.7 g/cm3) respectively. Obtain better statistic data by increasing the history number or by using the restart function Hints • Use 1st [t-deposit] tally in “ParticleTherapy” in “recommendation” folder • Increase the number of the radial bins in the [t-deposit] tally • Comment out unnecessary tallies using “off” command for speed up • Setup appropriate source particles and energy in [source] section • Specify particle type in the [t-deposit] tally • Add new material (Al), surface, and cells to set up geometry Summary 32
Example of Simulation Results r < 1 cm 1 cm < r < 2 cm Depth dose distributions for r < 1 cm and 1 cm < r < 2 cm Let’s think about … • What kinds of particles contribute to the dose behind Bragg peak? • Why dose suddenly increase at the depth of 5 cm? • Why we cannot see the neutron contributions in these graphs? Summary 33