UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Incorporating microdosimetry into radiation therapy treatment planning with multi-scale Monte Carlo simulations Lucido, Joseph 2013

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


24-ubc_2013_fall_lucido_joseph.pdf [ 7.89MB ]
JSON: 24-1.0074272.json
JSON-LD: 24-1.0074272-ld.json
RDF/XML (Pretty): 24-1.0074272-rdf.xml
RDF/JSON: 24-1.0074272-rdf.json
Turtle: 24-1.0074272-turtle.txt
N-Triples: 24-1.0074272-rdf-ntriples.txt
Original Record: 24-1.0074272-source.json
Full Text

Full Text

Incorporating Microdosimetry into Radiation TherapyTreatment Planning with Multi-Scale Monte CarloSimulationsbyJoseph LucidoB.S.E., The University of Michigan, 2005M.S., Rutgers University, 2008A THESIS SUBMITTED IN PARTIAL FULFILLMENTOF THE REQUIREMENTS FOR THE DEGREE OFDoctor of PhilosophyinTHE FACULTY OF GRADUATE AND POST DOCTORAL STUDIES(Physics)The University Of British Columbia(Vancouver)September 2013c? Joseph Lucido, 2013AbstractIn order to choose and design optimal treatment plans for radiation therapy, it isnecessary to employ models that can predict the tissue response to the ionizingradiation. The conventional models are based on the radiological absorbed dose,which does not account for the effect of radiation damage clustering on cellularresponse ? a more mechanistic model can lead to better metrics for treatment plan-ning. This dissertation presents a novel method for performing multi-scale MonteCarlo simulations to obtain microdosimetric information for patient specific treat-ments. This is done by using a track structure Monte Carlo simulation in the re-gions of interest for scoring and a condensed history algorithm for the rest of thegeometry. Since the condensed history code does not correctly follow the tracks ofparticles below a certain energy threshold, the volume in which the track structuresimulation is performed must extend beyond the volume in which scoring is done.The effect of this extended volume on simulation accuracy and performance arediscussed, and it is shown that the watch volume must extend beyond the target bya distance equal to the range of the subthreshold electrons. This simulation methodis benchmarked against experimental measurements for several radioisotopes andrun for a volumetric arc radiotherapy plan. In addition, there is a comparison ofthe microdosimetric characteristics of two widely used track structure simulations(Geant4-DNA and NOREC), and a discussion of the use of Monte Carlo in thepatient specific treatment planning for Stereotactic Body Radiotherapy and TotalBody Irradiation.iiPrefaceThe identification of the thesis project selected for this thesis was a collaborationof the author, the research supervisor Tony Popescu and supervisory committeemember Vitali Moiseenko. The design of the research project was carried out bythe author, in consultation with the supervisor and Dr. Moiseenko.In the second chapter, the commissioning of the new Monte Carlo model of thelinear accelerator treatment, and its comparison to the current model was entirelycarried out by the author, with some guidance from Tony Popescu. The VMATQA tools presented in Section 2.2 were mainly produced by me, although oth-ers (including Tony Popescu, Julio Lobo, Tony Teke) did contribute parts of it.I made numerous contributions to the research discussed in Sec. 2.3, which hasbeen submitted for publication (?Monte Carlo Calculation of Dose Distribution inOligometastatic Patients Planned for Spine Stereotactic Ablative Radiotherapy? byV. Moiseenko et al., in submission). This study involved the VMAT QA tools that Iwrote, as well as a new method to prepare the files for the Monte Carlo simulationsusing density heterogeneity corrections (but neglecting the atomic composition cor-rections). In addition, the software tools to convert the dose to water, and to run thesimulations for this required extensive modification. The analysis and writing inthat paper were predominantly carried out by the other authors, but the discussionof the results and conclusions presented in this thesis on that topic are my own. Thediscussion of dose-to-water conversion in Section 2.4 is my own viewpoint, but it isbased on numerous conversations I had with my supervisor. The RunTBI programfrom design to implementation, validation, and analysis were entirely my own, aswas the development of the Monte Carlo simulation user code and the analysis ofit.iiiChapter 3 is primarily based on the paper ?A microdosimetric comparison ofGeant4-DNA and NOREC for electrons in water? by J. Lucido, T. Popescu, and VMoiseenko (in submission). I wrote all of the software and performed the simula-tions and data analysis. The interpretation of the results was primarily my work,done in consultation with the other authors. One of the figures was put together byS. Incerti (who oversaw the GEANT4-DNA simulations for that figure), althoughI contributed the NOREC part of the data.Parts of the fourth chapter are submitted for publication (?A Method to Per-form Multi-Scale Monte Carlo Simulations in the Clinical Setting? by J. Lucido,T. Popescu and V Moiseenko). All of the simulations were designed and carriedout by me, although the design of the study and some of its analysis was a collabo-rative effort from all three authors. The part involving the clinical simulation doesnot appear in that paper, but was entirely my work.ivTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiiList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ixGlossary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiv1 Introduction To Radiation Oncology Physics . . . . . . . . . . . . . 11.1 Photon Interaction with Matter . . . . . . . . . . . . . . . . . . . 31.1.1 Coherent Scattering . . . . . . . . . . . . . . . . . . . . . 31.1.2 Incoherent Scattering . . . . . . . . . . . . . . . . . . . . 31.1.3 The Photoelectric Effect . . . . . . . . . . . . . . . . . . 61.1.4 Pair Production . . . . . . . . . . . . . . . . . . . . . . . 71.2 Charged Particle Interaction with Matter . . . . . . . . . . . . . . 81.2.1 Elastic Scattering . . . . . . . . . . . . . . . . . . . . . . 111.2.2 Excitation . . . . . . . . . . . . . . . . . . . . . . . . . . 121.2.3 Ionization . . . . . . . . . . . . . . . . . . . . . . . . . . 131.2.4 Variations in Track Structure . . . . . . . . . . . . . . . . 141.3 Radiochemistry . . . . . . . . . . . . . . . . . . . . . . . . . . . 14v1.3.1 Hydroxyl Radical . . . . . . . . . . . . . . . . . . . . . . 151.3.2 Superoxide . . . . . . . . . . . . . . . . . . . . . . . . . 151.3.3 Nitric Oxide . . . . . . . . . . . . . . . . . . . . . . . . 161.3.4 Hydrogen Peroxide . . . . . . . . . . . . . . . . . . . . . 161.4 Radiobiology . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171.5 An Overview of Contemporary Treatment Delivery Techniques . . 181.5.1 External Beam Therapy . . . . . . . . . . . . . . . . . . 181.5.2 Brachytherapy . . . . . . . . . . . . . . . . . . . . . . . 191.6 Why Dose Is Not Enough . . . . . . . . . . . . . . . . . . . . . . 201.7 Overview of Thesis . . . . . . . . . . . . . . . . . . . . . . . . . 212 Condensed History Monte Carlo for Dosimetry . . . . . . . . . . . . 252.1 Modeling a Linear Accelerator . . . . . . . . . . . . . . . . . . . 282.2 VMAT QA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 332.3 Comparison Between AAA and Monte Carlo for Spine SABR . . 362.4 Dose To Water . . . . . . . . . . . . . . . . . . . . . . . . . . . . 462.5 Total Body Irradiation . . . . . . . . . . . . . . . . . . . . . . . . 492.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 583 Track Structure Monte Carlo for Microdosimetry . . . . . . . . . . 603.1 Microdosimetry . . . . . . . . . . . . . . . . . . . . . . . . . . . 613.2 Limitations of the Quantitative Use of Microdosimetry . . . . . . 623.3 Physics Models . . . . . . . . . . . . . . . . . . . . . . . . . . . 633.3.1 Inelastic Interactions . . . . . . . . . . . . . . . . . . . . 633.3.2 Elastic Interactions . . . . . . . . . . . . . . . . . . . . . 653.3.3 Subexcitation Electrons . . . . . . . . . . . . . . . . . . 663.4 Global Track Measures . . . . . . . . . . . . . . . . . . . . . . . 673.5 Software Design for Microdosimetry Simulation . . . . . . . . . . 683.5.1 Initial Track Variables . . . . . . . . . . . . . . . . . . . 693.5.2 Scoring . . . . . . . . . . . . . . . . . . . . . . . . . . . 693.5.3 Range Rejection . . . . . . . . . . . . . . . . . . . . . . 723.5.4 Ensuring Convergence of the Simulations . . . . . . . . . 733.6 Parallel Beam . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73vi3.7 Isotropic Sources . . . . . . . . . . . . . . . . . . . . . . . . . . 783.8 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 824 Multi-Scale Monte Carlo . . . . . . . . . . . . . . . . . . . . . . . . 844.1 Watch Region . . . . . . . . . . . . . . . . . . . . . . . . . . . . 854.2 Uniform Electrons . . . . . . . . . . . . . . . . . . . . . . . . . 854.2.1 Simulation Methods . . . . . . . . . . . . . . . . . . . . 854.2.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . 864.3 Benchmarking with Photons . . . . . . . . . . . . . . . . . . . . 894.4 Clinical Use . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 914.4.1 Simulation Method . . . . . . . . . . . . . . . . . . . . . 924.4.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . 944.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 955 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103viiList of TablesTable 2.1 The percent of the planning treatment volume receiving certaindoses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38Table 2.2 Comparison of measured and calculated doses in an inhomoge-neous phantom . . . . . . . . . . . . . . . . . . . . . . . . . . 43Table 2.3 Spine SABR spinal canal doses . . . . . . . . . . . . . . . . . 45Table 2.4 Percent change in stopping power for different path lengths . . 48Table 2.5 Calculated TBI block thicknesses for a simple phantom . . . . 56Table 3.1 Excitation Energies in Geant4-DNA and NOREC . . . . . . . 64viiiList of FiguresFigure 1.1 Mass attenuation coefficients for photon processes in water . . 4Figure 1.2 Angular distribution of incoherent scattering . . . . . . . . . 5Figure 1.3 Electron stopping power in water . . . . . . . . . . . . . . . 10Figure 2.1 An Example of a DVH . . . . . . . . . . . . . . . . . . . . . 26Figure 2.2 Definition of the coordinate system for linear accelerator mod-elling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30Figure 2.3 Percent Depth Dose Calcuations for Different Beam Models . 31Figure 2.4 Beam Profiles along the in-plane axis . . . . . . . . . . . . . 32Figure 2.5 Beam profiles along the in-plane axis, showing the field region 33Figure 2.6 Beam Profiles along the cross-plane axis . . . . . . . . . . . . 34Figure 2.7 Beam Profiles along the cross-plane axis, showing the field re-gion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35Figure 2.8 Sample CT Phantoms produced by Eclipse and iCTC . . . . . 37Figure 2.9 Spine SABR PTV DVHs calculated by different inhomogene-ity corrections . . . . . . . . . . . . . . . . . . . . . . . . . . 39Figure 2.10 Spine SABR dose distributions comparison with different in-homogeneity corrections . . . . . . . . . . . . . . . . . . . . 40Figure 2.11 Spine SABR dose distributions with no heterogeneity corrections 41Figure 2.12 Spine SABR DVHs with no heterogeneity corrections . . . . 42Figure 2.13 Geometry of the inhomogeneity test phantom . . . . . . . . . 43Figure 2.14 Spine SABR DVH!s for the spinal canal calculated with differ-ent inhomogeneity corrections . . . . . . . . . . . . . . . . . 46Figure 2.15 TBI treatment set-up . . . . . . . . . . . . . . . . . . . . . . 51ixFigure 2.16 Representation of the projection problems of the TPS! for TBI 52Figure 2.17 Comparison of the DRR! images used for planning in TBI . . . 53Figure 2.18 RunTBI user interface . . . . . . . . . . . . . . . . . . . . . 54Figure 2.19 Cross-section of RunTBI validation phantom . . . . . . . . . 55Figure 2.20 Dose-Volume Histograms for TBI . . . . . . . . . . . . . . . 57Figure 3.1 Comparison of global track measures . . . . . . . . . . . . . 68Figure 3.2 Geometry for plane parallel beam simulations . . . . . . . . . 70Figure 3.3 Geometry for isotropic microdosimetric simulations . . . . . 71Figure 3.4 Microdosimetric spectra for 150, 250, and 500 keV electrons . 74Figure 3.5 Microdosimetric spectra for 80, 50, and 25 keV electrons . . . 75Figure 3.6 Microdosimetric spectra for 5 and 10 keV electrons . . . . . . 76Figure 3.7 Dose-mean lineal energies for electrons between 5 and 500 keV 77Figure 3.8 Ratio of dose-mean lineal energies at different depths . . . . . 79Figure 3.9 Microdosimetric spectra for electrons in 30 nm targets . . . . 80Figure 3.10 Dose-mean lineal energies for 30 nm targets . . . . . . . . . . 81Figure 4.1 Comparison of the microdosimetric spectra for 0.5 ?m watchvolume . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87Figure 4.2 Comparison of microdosimetric spectra for 2.5?m watch volume 88Figure 4.3 Mean lineal energies for different watch volume sizes . . . . . 89Figure 4.4 Average time per event scored for different watch volume sizes 90Figure 4.5 Microdosimetric spectra for photons in water . . . . . . . . . 91Figure 4.6 Dose-mean lineal energy for photons in water . . . . . . . . . 92Figure 4.7 Microdosimetric spectra for a VMAT treatment plan . . . . . . 95Figure 4.8 Electron energies in the clinical phase space file . . . . . . . . 96xGlossary3DCRT 3D Conformal RadiotherapyAAA Anisotropic Analytical AlgorithmAAPM American Association of Physicists in MedicineBER Base Excision RepairCT Computed TomographyCSDA Continuous Slowing Down ApproximationDICOM Digital Image Communication Protocol for MedicineDNA deoxyribonucleic acidDSB Double Strand BreakDVH Dose Volume HistorgramELF Energy Loss FunctionEPID Electronic Portal Imaging DeviceFWHM Full Width Half MaximumHDR High Dose Rate BrachytherapyICTC Improved CTCreate input makerIMRT Intensity Modulated Radiation TherapyxiKURBUC Kyoto Unversity Radiation Biology User CodeLDR Low Dose Rate BrachytherapyLEPTS Low Energy Particle Track SimulationLET Linear Energy TransferMCNP Monte Carlo for N-ParticlesMDS Multiply Damaged SiteMLC Multi Leaf CollimatorNHEJ Non-Homologous End-JoiningNOREC New Oak Ridge Electron CodeNTCP Normal Tissue Complication ProbabilityPARTRAC PARTicle TRACk codePDD Percent Depth DosePITS Positive Ion Track Structure codePTV Planning Target VolumeQA Quality AssuranceRBE Relative Biological EffectivenessRNS Reactive Nitrogen SpeciesROS Reactive Oxygen SpeciesSABR Stereotactic Ablative RadiotherapySBRT Stereotactic Body RadiotherapySSB Single Strand BreakTBI Total Body IrradiationxiiTCP Tumor Control ProbabilityTLD Thermo-Luminescent DosimeterTRION TRi-particle ION codeVMAT Volumetric Modulated Arc TherapyxiiiAcknowledgmentsI would like to thank my supervisory committee for their contributions that madethis possible. In particular, the guidance, support, and patience of Tony Popescuand Vitali Moiseenko were invaluable to me. Equally important, I am grateful tothe members of my family who have helped in a myriad of other ways during thiseffort, including my brother, my parents, and Beata and Steven Baker. I wouldbe remiss not to include my gratefulness to my friends who help ensure that thisproject was a success.xivChapter 1Introduction To RadiationOncology PhysicsRadiation therapy is the use of ionizing radiation to treat a tumor site. The radia-tion transfers energy to the tissue, which can result in the reproductive inactivationof cells in that tissue. However, the energy transfer cannot be confined only to thetumor site because of the nature of radiation tracks. This energy deposition meansthat surrounding normal tissue can potentially develop complications. Thus, thereare two competing objectives in radiation therapy: to ensure sufficient radiationdose is delivered to control the tumor while limiting the exposure of the normaltissue. This dual mandate is the fundamental challenge in radiation therapy and thesource of any unfavorable outcomes - it has driven the long history of the develop-ment of better methods of treatment. Radiation can be delivered through the use ofexternal beams as well as both sealed (brachytherapy) and unsealed radio-isotopes.Modern clinical therapy makes use of photons, electrons, and protons (although awide range of others have been considered). A major development has been the useof multi leaf collimators (MLCs) to define the field shape for external beam therapyin order to achieve better normal tissue sparing. Contemporary delivery techniquesare discussed in more detail in Section 1.5.In order to make the optimal choice from all the options available, it is nec-essary to evaluate the merits of specific treatment plans. At the most basic level,this means making prediction about patient outcomes and requires an understand-1ing of the mechanisms underlying the patient?s biological response to the therapy.It is possible to divide these mechanisms into different stages, each with differentprocesses. The first stage involves the physical interaction of the radiation withthe material it is passing through. Energy is transferred to particular moleculesalong the radiation tracks. In living tissue, this involves both the bulk water of thecell as well as direct damage to the nuclear DNA and other biomolecules. Theseinteractions are discussed further in Section 1.1 and Section 1.2. The next stageis dominated by different forms of chemistry, which involves the radiolysis of thebulk fluid and the creation of reactive species. These species then diffuse and in-teract both with other reactive species, scavengers, and biomolecules like DNA. Amore expansive treatment is given in Section 1.3. The attempt to repair the dam-age to the nuclear DNA makes up the biological response stage, which starts inthe seconds after the radiation traverses the cell. The last stage involves the tis-sue, organ, and systematic responses to the loss of cells. This is what determinesthe outcome of the treatment, and is thus of the ultimate aim for the modeling intreatment planning, these make up Section 1.4.However, to make use of the radiobiology when making these predictions, thechemistry associated with the ionizing radiation must be understood on the scaleof nanometers (which corresponds to a few hundred DNA base-pairs). The dy-namics of the chemistry is determined by the locations and magnitudes of inelasticscattering events of the ionizing radiation with the medium. The best approachto modelling the ionization tracks is radiation transport Monte Carlo simulations.However, to deal with the complex geometry and radiation spectra that exist inclinical situations and still perform the simulations in a reasonable amount of time,it was necessary to develop a sophisticated new approach to these simulations.This thesis presents the development, implementation, and application of this newmethod. This new tool allows clinicians and researchers to obtain detailed infor-mation about the structure of the radiation tracks for specific treatment plans for anindividual patient. Since this data was not previously available, there is no consen-sus on how it should be employed. But this new approach facilitates this discussionby allowing it to be placed on a firm footing. For instance, comparing proton orheavy-ion radiation treatments with existing high energy photon treatments for aparticular site in an individual patient requires a level of understanding of the ra-2diobiology which is currently based on in vitro results whose relevance in thesesituations is far from established, but now is open to investigation prior to runninglarge scale clinical trials.1.1 Photon Interaction with MatterWhen ionizing photons travel through a material, they may interact with the atomsmaking up the material and transfer energy and momentum to atomic electrons.In the energy ranges of interest in this thesis, there are four distinct processes bywhich the photons can interact with the material: coherent (or Rayleigh) scatter-ing, incoherent (or Compton) scattering, the photoelectric effect, and pair produc-tion[40]. These are discussed individually in Section 1.1.1 to Section 1.1.4. Thecontributions of each of these processes to the total mass attenuation coefficient areshown in Figure 1.1. The interaction process often results in the transfer of someof the energy to one or more charged particles, which are ultimately responsiblefor directly transferring most of the energy to the material, through the channelsdiscussed in Section Coherent ScatteringWhen the interaction of a photon and the electrons of an atom results in the ab-sorption of the photon and the energy absorbed is too low to ionize the atom, theprocess is referred to a coherent (or Rayleigh) scattering. This results in the tem-porary excitation of the atom, which then relaxes by the emission of a photon ofnearly the same energy and momentum direction as the original absorbed photon.This is an important process only for photon energies below the ionization thresh-old of the medium, and therefore is relevant for photons with energies in the ultraviolet range or below. It plays less of a role in clinical radiotherapy using x-rays[5],as shown in Figure Incoherent ScatteringThe transfer of some (but not all) of the energy of the incident photon to an elec-tron is called incoherent or Compton scattering if the electron is unbound after thetransfer of the energy. This can occur for a free electron or an atomic electron as310?3 10?2 10?1 100 101 10210?810?610?410?2100102104Photon Energy (MeV)Mass Attenuation Coefficient (cm2/g)  CoherentIncoherentPhotoelectricPair ProductionTotalFigure 1.1: Mass attenuation coefficients for different photon processes inwater, data taken from[6].long as the photon energy is greater than the ionization energy of that electron.The likelihood of this interaction occurring is greatly diminished for more tightlybound atomic electrons, and so it is generally sufficient to consider only the nearlyfree electrons in the material, which means it is reasonable to treat the electrons asunbound for a first approximation. If the photon is scattered at an angle ? , thenmomentum conservation implies that the photon?s energy after the scattering eventis given byE ? = E1? Em0c2 (1? cos?)(1.1)where E is the photons original energy, m0 is the electrons rest mass, and c is thespeed of light. From conservation of energy, the electron is emitted at an angle? = cos?1(1+ tan?(1+?)2)(1.2)4  5e?005  0.0001302106024090270120300150330180 0  5e?005  0.0001302106024090270120300150330180 0  5e?005  0.0001302106024090270120300150330180 0  5e?005  0.0001302106024090270120300150330180 028 keV 660 keV1 MeVa) b)c)6 MeVd)Figure 1.2: Angular distribution of incoherent scattering calculated from theKlein-Nishina formula for (a) 28 keV, (b) 660 keV, (c) 1 MeV, and (d) 6MeVwhere ? = E/m0c2. If the spin and magnetic moments are considered, the Klein-Nishina formula[45] gives the differential angular cross-section:d?d?=r202(1+ cos2?)FKN (1.3)where r0 is the classical radius of the electron andFKN =(11+?(1? cos?))2(1+ ?2(1? cos?)2[1+?(1? cos?)](1+ cos2?))(1.4)The angular distributions of the differential cross section for several different ener-gies are shown in Figure 1.2.By integrating Equation 1.3 over all scattering angles, the mass attenuation co-5efficient can be calculated. For the energy ranges used in radiation therapy, themass attenuation coefficient approximately inversely proportional to the square-root of the photon energy. This analysis can further be developed by incorporatingform factors that take into account the binding energy of the electrons. Since theincoherent scattering in materials predominantly involves the least bound atomicelectrons, the linear attenuation coefficient is proportional to the ratio of the atomicnumber to the atomic mass. This ratio is nearly constant for all of the most commonelements in biological material, which implies that the attenuation due to incoher-ent scattering depends on atomic density and is roughly independent of materialcomposition. This has important consequences for radiation therapy because, ascan be seen from Figure 1.1, this process is the predominant interaction in radia-tion therapy.1.1.3 The Photoelectric EffectThe photoelectric effect involves the absorption of the entire photon of energy E?by an atomic electron and the emission of an atomic electron with an energy Ee =E??Eb where Eb is the binding energy of that electron. Unlike in the incoherentscattering case, this process favors interaction with more deeply bound electrons(in particular, electrons in the K-shell), provided the photon is energetic enoughto ionize these electrons. Because this interaction is more likely to occur withinner shell electrons, the mass attenuation coefficient is proportional to Z3/A[45].However, the energy dependence is somewhat more complex. While in general thecross-section is proportional to 1/E3, there is significant enhancement when thephoton energy is slightly above the binding energy of electrons of an inner shell.For most of the materials of interest to radiation therapy, this means the K-shell,and is therefore referred to as the K-shell edge effect, and occurs at energies lessthan 100 keV. This means that the photo-electric effect is less important for mostclinical photon beams than incoherent scattering.Since the ionization generally involves a deeply bound electron, the ion is leftin an excited state and must return to its ground state. Higher energy electrons inthe ion will drop into this vacancy, leaving behind another vacancy. This processcontinues until a valence orbital is unoccupied, and a passing thermal electron can6then recombine with the ion leaving the atom in its ground state. Since the energyof the orbitals of the ion are quantized, this means that the de-excitation of the ionresults in the emission of photons at characteristic energies. The de-excitation canalso occur by the spontaneous emission of auger electrons, but energy conservationconstrains the energies of these electrons to characteristic energies as well onceone includes the binding energy of the electron. It is possible for a de-excitationcascades to include both characteristic x-rays and auger electrons.1.1.4 Pair ProductionA photon can interact with the electric field of a nucleus to create an electron/-positron pair, a process referred to as pair production. This requires that the photonhas energy higher than the combined rest mass of the particle pair (1.02 MeV). Ifthe photon energy is above this lower limit, additional energy can be transferredto the particle pair, giving them a large kinetic energy. For photons with energygreater than 1.02 MeV, the linear attenuation coefficient increases logarithmicallywith the photon energy. Figure 1.1 demonstrates that this process is the dominantphoton interaction at energies above 10 MeV in water. On the other hand, this pro-cess depends heavily on the strength of the nuclear electric field an therefore themass attenuation coefficient is proportional Z2/A, which means for most biologicalmaterials of interest, it scales roughly as Z. However, in high-Z materials, Z/A isnot constant and there is significant enhancement of the contribution of pair pro-duction to the total attenuation of the photon beam. This is one reason why lead isa favorable choice for radiation shielding, and is an important consideration in theuse of tungsten in linear accelerator MLCs, flattening filters, and x-ray targets.Pair production is further distinguished by the creation of a positron. Thepositron will eventually recombine with an electron and release a pair of photonsin roughly opposite directions. While this can occur while the positron has a largekinetic energy (known as annihilation in flight), most commonly the positron willthermalize and recombine with a lightly bound or free electron (called annihilationat rest). This means that the annihilation photons will be near the characteristicenergy of 511 keV, that is, the rest energy of a positron (or an electron).The photon can also interact with the electric field of an electron. For lower en-7ergies this is less frequent than interactions with the nuclear electric field. In somecases, enough energy is transferred to the electron that it can overcome its bindingenergy. This results in the emission of two electrons and a positron (and possiblycharacteristic x-rays or auger electrons as the ion relaxes to the ground state). Thisis known as triplet production. In fact, higher order production can occur if suf-ficient energy is transferred by the reaction to ionize multiple electrons from theatom. However, these processes require increasingly large initial photon energies,and even triplet production is only important at the highest energies typically usedin radiation therapy.1.2 Charged Particle Interaction with MatterThe scattering of charged particles from atoms and molecules in a medium is cru-cially important because it is the mechanism by which nearly all of the energy ofan incident charged or uncharged radiation source transfers energy to the medium.Although there has been a wide range of charged particles considered for radia-tion therapy, currently the most relevant ones are electrons, positrons, and protons.For the energy ranges used in clinical radiation therapy, it is sufficient to consideronly mechanical and electric interactions with the material. In this case there aretwo general channels by which the charged particles lose energy: collisions andbremsstrahlung emissions.The charged particle can exchange discrete energy packets with individualatoms or molecules in the material, and these are known as collisional (or ioniza-tional) energy losses. In general, the mean free path between these events is veryshort relative to the range of the charged particle, and most of these events involverelatively small energy transfers relative to the charged particle?s total energy. Thismeans that if we wish to consider the effects of many of these scattering events,the energy transfer per unit thickness will approach a deterministic average knownas the collisional mass stopping power, Scoll , so long as the total energy transfer issmall compared to the kinetic energy of the charged particle (or alternatively, thedistance considered is short relative to the range of the particle, but long comparedto the mean-free path between collisions). For relativistic electrons (which are usu-ally responsible for most of the direct energy transfer to the medium), it is possible8to analytically calculate Scoll for a classical electron while still taking into accountthe binding energy of the electrons as well as charge screening effects to find[40]Scoll = 2pir20NeEm0?2(lnK2(K+2Em0)2Em0I2+K2/8? (2K+Em0)Em0 ln(2)(K+Em0)2+1?? 2??)(1.5)where:? r0 is the classical electron radius? Ne is the number of electrons per mass of the material? Em0 = m0c2 is the electron rest mass energy? ? is the ratio of the electron?s velocity to the speed of light? K = 1/2mv2 is the relativistic kinetic energy of the electron? I is the mean binding energy of the electrons in the material? ? is the screening effect correction.The collisional stopping power for electrons in water is shown in Figure 1.3. Thiscan also be referred to as the unrestricted (electronic) linear energy transfer or LET.However, photons can also be discussed in terms of LET, but this actually refers tothe LET of the secondary electrons set in motion by the photons. Specific types ofcollisional interactions are discussed in Section 1.2.1 to Section 1.2.3.A charged particle can also lose energy by the emission of electromagneticbremsstrahlung radiation. This occurs whenever there is a change in the momen-tum of charged particle. For charged particles in matter, this is generally becauseof an interaction with the nuclear or (less frequently) atomic electrical fields. Sincethe acceleration is larger when the field is larger, the probability of producing abremsstrahlung photon increases with the charge of the radiation particle and thecharge of the scattering atom. In addition, the energy and number of the photonsemitted increases with the energy of the charged particle. Specifically, for elec-910?2 10?1 100 101 102 10310?310?210?1100101102Energy (MeV)Stopping Power (MeV cm2 /g)  ScollSradStotFigure 1.3: The collisional, radiative, and total stopping power for electronsin water, with data from[7].trons, the radiative mass stopping power has been calculated to be[40]Srad = 4r20NeZK137(ln2(K+Em0)Em0? 13)(1.6)where Z is the atomic number of the scattering material. As shown in Figure 1.3,this implies that the energy loss to bremsstrahlung radiation plays a more importantrole at higher particle energies, as well as higher Z materials. In many practicalapplications, it is most relative to consider the total stopping power, Stot = Scoll +Srad .An alternative way to discuss the transport of charged particles is through theuse of the continuous slowing down approximation (CSDA). In this approximation,as the particle traverses the medium the rate of energy loss is assumed to be equal10to the total stopping power for the particle?s instantaneous energy. By integratingthe stopping power over the particle?s energy, it is possible to define the particle?srange and it can be interpreted as the total path length of the particle?s track in themedium.Since the rate of energy loss greatly increases as the particle loses energy, mostof the actual energy transfer occurs near the end of the track. This enhancement isknown as the Bragg peak. For particles that have a large rest mass relative to theaverage energy transfer in the tracks, the total change in direction of the particles isvery small, and so the Bragg peak will occur at a well defined depth in the material.This includes protons and heavy ions. This phenomenon is a major motivationfor the interest in proton therapy. Electrons and positrons, meanwhile, tend to bescattered through much larger angles which smears out the penetration depth atwhich the Bragg peak occurs. Anti-particles, such as positrons, will also undergoannihilation (as discussed in Section 1.1.4), but the (non-annihilation) interactioncross-sections are roughly the same as for their conjugate.Collisional processes can be grouped based on the state of the scattering ma-terial after the interaction. If there is a negligible transfer of energy, the processis referred to as an elastic event. This channel can be important because the scat-tering angle of the charged particle can be large. The scattering site can also haveenough energy to be excited out of its ground state. In some cases, this can re-sult in the formation chemically active reactive species, which play an importantrole in radiation chemistry and biology. Finally, the scattering site can becomeionized. More details on these charged-particle impact interactions are discussedbelow, with particular focus on electrons in water.1.2.1 Elastic ScatteringIn the literature on charged particle transport, the term ?elastic? refers more broadlyto collisions in which the magnitude of the energy transfer is negligible comparedto the kinetic energy of the incident particle. Although the energy exchange is verysmall, the momentum transfer can be very large, particularly for slow, light parti-cles. This is most relevant for describing the end of electron tracks. While there hasbeen extensive theoretical treatment of these processes ? and an increasing pool of11measurements ? there is still significant room for more experimental data[15].One component of the elastic scattering potentials comes from the Coulombpotential of the atoms and molecules in the material (this can be called the staticcomponent) as well as a contribution from any molecular dipoles in the medium.There can also be an additional exchange term in the potential if the scatteringparticle is identical to any of the particles in the molecule. These exchanges havebeen shown to be important for electrons interacting with a wide variety molecules,and the process can be modelled as rotationally symmetric[33]. For electrons inwater, the static term accounts for nearly all of the scattering potential above 500eV, while it is less important than the exchange and dipole terms below 300 eV[12].1.2.2 ExcitationThe internal state of the scattering site can become excited either mechanically orelectronically. Scattering by excitation is the most relevant scattering process whenthe charged particle?s incident energy is insufficient to produce an ionization buttoo large for elastic scattering. The method of relaxation, dictated by the type ofexcitation and the energy transfer, determines the radiobiological implications ofthe event.Molecules can experience mechanical excitations in rotational and vibrationalmodes. The available modes and energy levels are determined by the symmetryproperties of the molecule as well as the atomic composition. There tend to bemany possible mechanical excitations at low level: for water the largest rotationalstates are only a few milli-eV (meV) above the ground state, while the vibrationalmodes are in the range of 0.1 eV to 1 eV[38]. It is possible for an charged particleimpact to result in a change from one vibrational mode to another at the sameenergy accompanied by a change in the rotational state of the molecule. Since therotational energy levels are much lower the vibrational modes, this is referred toa ?vibrationally elastic? scattering. These low energy excitations do not dissipatethrough the creation of radicals or ions, but contribute to the thermal energy of themedium.The energies involved are larger for electronic excitation: for water, they be-gin at an energy of about 7.4 eV. Importantly, the scattering by charged particles12includes transitions that are not available during photon interactions. There areseveral intermediate states that the excitation can achieve before relaxation. It ispossible that the atomic electrons occupy higher orbitals. In this case, the relax-ation can come by means of the emission of characteristic x-rays or auger electrons,or the energy can be transferred to other chemical species in the medium. Alter-natively, the excitation can cause the molecule to dissociate, producing an excitedradical (OH?, H?, O?)[31]. On the other hand, an incident electron can becomeattached to the molecule, giving it a net negative charge, and the eventual dissoci-ation of the molecule into a negative ion and other products. For water this meansthat there are a variety of outcomes that includes a reactive species:? H2 + O?? H + H + O?? OH + H?? O + H + H??These excitation products play an important role in the radiobiology of living tis-sue.1.2.3 IonizationThe energy transfer can also result in the ionization of the scatterer, and the ionizedelectron can potentially gain significant kinetic energy from the transfer. Further-more, if the charged particle has sufficient energy, an inner shell electron can beionized. This leaves the ion in an excited state, and it will relax in the same man-ner as a photo-ionized atom (see Section 1.1.3). For liquid water, the minimumionization potential is 8.8 eV. There are three other levels below 32.2 eV, but thenext set of states is at energies above 500 eV[12]. Ionization plays an increasinglydominant role in the scattering process of electrons in water starting at around 15eV[31].131.2.4 Variations in Track StructureThe spatial distribution of energy transfer events varies between particle types andenergies. The spacing between events is greater for lighter mass particles thanheavier ones and even greater for photons. This can be observed in the LET[37] val-ues as well as different measures of event clustering[67]. In radiation therapy, theenergy spectra can be very complex due to the spread in the initial bremsstrahlungspectra as well as from variations caused by the scatter within the patient. Both ofthese effects mean that the characteristics of the spectra change in different partsof the patient, for instance inside and outside of the radiation field[69], as well aswithin different parts of the field. This has important consequences for the radia-tion chemistry as well as radiobiology, and needs to be addressed when modelingrealistic treatments (as discussed in Section 1.7).1.3 RadiochemistryThe excitations and ionizations created by the energy transferred from the radia-tion leads to the creation of free radicals, which are molecules or atoms that haveat least one unpaired electron in a shell. There is a huge variety of recognized rad-ical species that play a role in the chemistry of radiation damage, including both?reactive oxygen species? (ROS) and ?reactive nitrogen species? (RNS)[70]. In theimmediate aftermath of the passing of the radiation, the products are distributednear to the location of the energy deposition events, forming spurs along the trackof the high energy radiation particle. The evolution of the chemistry involves dif-fusion from these initial locations, as well as a competition between a variety ofinteractions. One type of interaction involves free radical scavengers, which resultin the formation of non-radical species. Furthermore, the free radicals can interactwith biomolecules, and a variety of other chemicals which lead to the creation ofsecondary radical species. There is a myriad of specific interactions of each ofthese types[105], and so the rest of the section will focus on the most importantradical species and their interactions, specifically the hydroxyl(?OH), superoxide(?O?2 ), and nitric oxide (NO?2) radicals, plus hydrogen peroxide (H2O2).141.3.1 Hydroxyl RadicalThe hydroxyl radical has a very short lifetime in water (less than 1 ns)[99]. Thishigh reactivity, combined with the variety of ways it can be produced from otherradicals, makes it responsible for most of the direct damage to DNA bases[80][70].The hydroxyl radical can undergo addition reactions with the double bonds of withall of four of the heterocyclic DNA bases (adenine, thymine, guanine, and cyto-sine), resulting in OH-adduct radicals[18]. It can also abstract a hydrogen fromthymine or the deoxyribose backbone, creating an allyl radical or a sugar radical,respectively[18]. These radicals can be oxidized leading to deprotonation and thecreation of a glycol of the DNA base[99]. More generally, hydroxyl radicals caninteract with biomolecules through hydrogen abstraction, electron transfer, and ad-dition reactions, resulting in another (less reactive) radical. Another critical inter-action is with hydrogen peroxide, which results in a hydrogen ion and a superoxideradical. It also readily interacts with nitrite to produce nitric oxide[70]. The en-dogenous scavengers of hydroxyl radicals include vitamin E-derivatives, uric andascorbic acids, glutathione, and acetylcysteine[97].1.3.2 SuperoxideSuperoxide itself does little damage directly to DNA, but its indirect in vivo ef-fects are far reaching[99]. It can directly interact with metals in biomoleculesthrough reduction, to create a hydroxyl radical, or oxidize a biomolecule whichcan lead to damage to that biomolecule while also producing a hydrogen per-oxide molecule[99]. Furthermore, it can react with other radicals and ions toproduce an wide variety of reactive products. It can undergo dismutation withhydrogen ions to produce hydrogen peroxide, the rate of which can be greatlyenhanced by enzymes[105]. It can also reduce nitric oxide to produce perox-ynitrite (ONOO?), which can interact with a variety of other molecules to pro-duce hydroxyl radicals[70][105]. The primary cellular scavenger of superoxide isascorbate[70]. It can be produced by interactions of molecular oxygen with hydro-gen radicals or the attachment of hydrated electrons[80]. In radiation chemistry,superoxide is most commonly produced by the attachment of a hydrated electronto molecular oxygen[105].151.3.3 Nitric OxideLike superoxide, nitric oxide is a primary radical for which its biological conse-quences result from the secondary products it produces, and there is substantialcross-talk between the reactive oxygen species (the chemistry resulting from hy-droxyl radicals, superoxide radicals, and hydrogen peroxide) and the reactive ni-trogen species[105]. In particular, reactions with superoxide result in the highlyreactive peroxyl nitrite radical discussed above. Another important interaction ishigh rate of interaction with molecular oxygen to generate nitrogen dioxide radi-cals (NO?2), which is less reactive than some of the other radicals, but more diffu-sive and selective of biomolecular targets[70]. Furthermore, the nitrogen dioxideradicals can interact with another nitric oxide radical to a form dinitrogen triox-ide (N2O3), which readily causes the deaminisation of DNA molecules[70]. Thischain of events is thought to contribute the radiosensitization of hypoxic cells vianitrogen based chemicals, including the nitroimidazole-class[70]. The scaveng-ing of reactive nitrogen species has much in common with that of reactive oxygenspecies, but the importance of reactive nitrogen has only recently been recognizedand the details of the scavenging are still being worked out[34].1.3.4 Hydrogen PeroxideHydrogen peroxide reacts poorly with most biomolecules[105], but it is longerlived than the other species and can create a wide range of more reactive species,which makes it critical to understanding in vivo free radical chemistry[80]. Mostimportantly, there are a variety of interaction pathways that can result in the cre-ation of a hydroxyl radical, manyh of which may be faciliated by catalase andperoxidase enzymes[80]. With peroxidase substrates, it can also interact to withphenol groups and thiols to form a wide range of bimolecular radicals[105]. Itcan also produce hypochlorous acid (HOCl)[105], and interact with NO?2 ions toproduce nitrogen dioxide radicals[70]. Since hydrogen peroxide is also a commonproduct of endogenous cellular mechanisms[34], there are a wide variety of cel-lular mechanisms for preventing damage from hydrogen peroxide, and it can bescavenged by many thiol groups[80].161.4 RadiobiologyThe free radicals (including reactive oxygen species, ROS, and reactive nitrogenspecies, RNS) created by the radiolysis of the nucleoplasm of a cell can damagethe DNA in several ways[52]. The nucleotide bases themselves can be altered,for instance by damage to the pyrimidine or purine rings, or removed from theDNA altogether. Additionally, damage can be done to the deoxyribose backbone,causing a single strand-break (SSB). If the radiation track yields multiple radicalsor reactive molecules near to each other, they can create a cluster of damage in theDNA. A double strand break (DSB) can be formed if there are two SSBs on oppositestands within close proximity. More complex forms of damage are also possible,where there is at least one DSB plus other types of damage. These are referredto as multiply damaged sites (MDS). The damage will cause the cell to attemptto undergo the repair process. If this process fails that can lead to cell death (thefailure to reproduce) and apoptosis[8]. To prevent these outcomes, the cell mayattempt to avoid to a potentially fatal cell division by entering cell cycle arrest andgive the DNA damage repair mechanisms an opportunity to fix the damage.There are three pathways that the cell can use to repair the damage to the DNAcaused by the ionizing radiation. Base damage or loss, as well as SSBs, are cor-rected by base excision repair, BER[52]. This removes the damaged base, andreplaces it with a new one. Since of there is a complementary nucleotide at thesite, this is an error free process for simple damage. However, if there is anotherdamage site within 2-3 base pairs, this may induce a lethal MDS from previouslynon-lethal damage[52]. In the case of DSBs, there are two different repair path-ways. If there is an undamaged sister chromosome available, it can be used as atemplate for error-free homologous repair[52]. If the matching chromosome is notused, the repair process is known as non-homologous end joining (NHEJ). Certainproteins bind to ends of the strand breaks, and processing is done to remove anyflaps remaining on the breaks. During this process other nearby damage can berepaired, including SSBs and altered bases. Finally the strands are re-ligated, al-though some base pairs may have been removed and will be missing in the finalproduct, which introduces permanent errors into the chromosome. The illegitimaterejoining of ends from different DSBs is another source of errors, resulting in a17chromosome aberration[84]. There is a strong relationship between the creation ofa chromosome aberration and cell death[52]. In all of these cases, the presence ofcomplex damage greatly increases the risk of producing a fatal error in the dam-aged DNA. This is caused by the correlation of production of the free radicals, andultimately is linked to the structure of the radiation tracks ? it depends on whathappens at the nanometer scale. Since these lengths are much smaller than typicalinter-track distances in radiation therapy[66], accurate modeling of the cellular re-sponse to the radiation requires more than just the knowledge of the total energydeposition in the tissue; knowledge of the correlations in the locations and times ofthe individual events is necessary. The measurement or calculation of these quanti-ties is known as microdosimetry[37]. The microdosimetry is important in radiationoncology because we ultimately need to know the response of the tissue or tumor,which depends on the number, location and type of cells killed.1.5 An Overview of Contemporary Treatment DeliveryTechniquesThe long-term goal of radiation oncology has been the controlled localization ofthe energy deposition of the radiation in order to get the best possible treatmentoutcome and has resulted in several different delivery techniques which are cur-rently employed. External beams achieve this object by precise field shaping(Section 1.5.1). Another option is the use of sealed radio-isotopes placed internallyin the patient, known as brachytherapy (Section 1.5.2). Historically the methodswere considered complimentary - certain disease sites favored one approach overthe other. More recently, the range of overlap has expanded. The implications ofthis are discussed in Section External Beam TherapyExternal beam therapy makes use of high energy photons, electrons, or protonsto deliver the treatment. Today, this is mostly done with linear accelerators andbeam shaping is often employed to ensure that the energy deposition by the radia-tion is predominantly delivered to the target volume. In some cases a homogenousdose is desired because the target volume is poorly defined. In these situations,18placed fields from one or more directions are used, and are shaped to spare nor-mal tissue. In particular, this may be the preferred approach if the beam can bedelivered from several directions without substantially involving sensitive normaltissue. For instance, the use of three large fields is a common treatment for rec-tal cancer in order to spare the bowel. For other types of treatment, it may benecessary to use an more complexly shaped field. Most commonly, this is doneusing a set of thin metallic leaves, called a multi-leaf collimator (MLC) which is at-tached to the treatment head of the linac. If the position of the leaves is fixed duringeach field, the technique is known as 3D conformal radiation therapy (3DCRT)[24].However, if a treatment plan is designed by using inverse planning to optimizethe fluence, it is referred to as intensity modulated radiation therapy IMRT. TheIMRT approach can allow more complex treatment plans, which can improve dosedistributions[102]. More recently, this has been combined with continuous gantryrotation during treatment[71]. This approach, referred to as volumetric modulatedarc therapy (VMAT), can spare normal tissue better and reduce delivery time. Inaddition, these techniques can be used for several kinds of particle energies ? mod-ern linear accelerators have multiple energies, and can deliver either photons andelectrons or protons. There are several means of sculpting an optimal distributionof dose using external beams.1.5.2 BrachytherapyThe main alternative to external beam therapy is brachytherapy. Brachytherapy isoften used to treat cancers of the prostate, gynecological organs, skin, and breast,and many other sites on a less frequent basis. Dose is localized by positioningsources inside the patient, where the inverse square law and attenuation concentratemost of the dose to the vicinity of the target volume. The main classification of thesources is based on the dose rate. High dose rate brachytherapy (HDR) uses a highactivity 192Ir source that is placed into a cavity with an applicator to irradiate theadjoining tissue or is inserted through catheters into the tissue. In both cases, oneHDR source moves to various positions where it dwells for a short period of time (atmost a few minutes). A patient will generally receive several fractions of treatment,typically between 2 and 10.19For low dose rate (LDR) therapy, multiple sources are placed into the desiredlocations. The sources can be placed via needle or wire directly into tissue, cavitiesor lumen, or externally for superficial treatments. The sources are either removedafter a suitable time (at least 24 hours) or the seeds are left permanently in thepatient. Over the years, many different isotopes have garnered attention as potentialLDR sources, and 103Pd and 125I, 137Cs, and 192Ir are the most widespread today.In brachytherapy, the main obstacle to achieving the ideal dose distribution isthe location of the seeds themselves. There are inherent difficulties in the initialplacement of the delivery devices (whether the catheters or sources themselves),as well as seed movement and anatomical changes (such as organ enlargementin response to the insertion). This requires that treatment planning is done bothbefore, after, and often during the delivery process to ensure that the actual dosedistributions are acceptable.1.6 Why Dose Is Not EnoughAchieving an optimal balance between the two objectives of radiation therapy re-quires that potential treatment plans be evaluated. This evaluation means it is nec-essary to project the outcome of an individual patient to a specific course of treat-ment, and so predictive metrics need to be used. The current standard metrics inclinical treatment planning are based on the radiological absorbed dose. This isthe amount of energy absorbed in a macroscopic volume of the patient per unitmass. Since the dose is averaged over a region that is large compared to the scalesrelevant to the formation of DNA damage clusters, its use assumes that there isconstant relationship between the total amount of energy deposited and the micro-scopic track structures of the radiation. Perhaps the most important use of dose isin the construction of tumor control probability, TCP, and normal tissue compli-cation probability, NTCP[103][59]. These are curves fit to the fraction of patientsshowing a particular outcome for a certain prescription dose. If the patients under-went similar treatments, then the characteristics of the track structure will not varygreatly between the patients, and so this is a reliable approach. In addition, it isuseful to decide the total amount of radiation to prescribe to a patient for a particu-lar course of treatment. However, dose is not as meaningful for the comparison of20plans of different modalities and techniques because the energy and radiobiologi-cal properties of the alternative plans may be very different. This is an increasinglyimportant issue given the adoption of so many different treatment techniques andradiation sources. In IMRT and VMAT, the scatter dose can contribute substantiallyto the dose in a larger volume of the patient, and this may have radiobiologicalconsequences[92]. In order to make treatment decisions between different tech-niques and modalities, dose is not enough.One attempt to overcome this shortcoming has been the use of the relative bi-ological effectiveness[37], RBE, which is the ratio of the doses required to achievea fixed cell survival fraction between two different radiation sources. This factorcan then be multiplied by the dose to give some idea how the living tissue willrespond. There are complications in practice, since a typical external beam un-dergoes a large amount of scatter as it passes through the patient, causing largevariations in the energy fluence between different sub-volumes of the patient. Withthe increasing complexity of treatment plans, this can be a major concern. For in-stance, a 10% variation in RBE has been shown for clinical photon[69] and protonbeams[73]. This is unacceptably large for clinical practice, since this uncertaintyis greater than the tolerance levels recommended for treatment QA[49][47]. Aone-size-fits-all correction is not accurate enough for modern treatment planning.1.7 Overview of ThesisIt is necessary to base treatment planning decisions on a more mechanistic modelthat can be used for a specific treatment plan. In order to predict the damage doneto the DNA, the radiochemistry must be modeled, and this depends on the trackstructure of the charged particles. Furthermore, these track structures must capturethe complexity of the spectrum being used, and the only approach that can get thisinformation is Monte Carlo simulations. A Monte Carlo simulation uses randomnumbers to sample transitions between states in a stochastic system. For radiationtransport, a particle?s trajectory is determined by two steps. The distance to thenext interaction is determined randomly using the total interaction cross-section,and then the nature of the scattering event is determined by random sampling fromthe different types of allowed interaction processes. This process continues until21interest in that particle ceases, because its kinetic energy has fallen below a cut-offor it has moved out of the region of interest.The objective of this thesis is to develop a Monte Carlo track structure sim-ulation tool that can be applied to external beams for clinical radiation therapyand to demonstrate the use and relevance for patient treatment plans. There aretwo approaches to simulating charged particle transport. First of all, the sim-ulation can model each interaction discretely. This is known as a track struc-ture (or class I) algorithm[67][13], and examples include NOREC[87], Geant4-DNA [36], PARTRAC[27][28], KURBUC[98], MOCA8[74], PTra[32], PITS[104]and TRION[50]. While this method makes all microdosimetric information avail-able, it is computationally very expensive because charged particles undergo a verylarge number of events. Since most of these events involve small energy and mo-mentum transfers, they do not drastically alter the particle?s trajectory. This makesit possible to average over many of these small events and treat only larger eventsdiscretely, and this creates a more efficient simulation. In addition, low energycharged particles have very short ranges, so it is often convenient to neglect the mo-tion of these particles and assume that all of their energy is lost in the immediatevicinity of the particle?s position. A simulation making use of these approxima-tions is referred to as a condensed history (or class II) algorithm[42], and havebeen implemented in many different code systems: EGSnrc[43], Geant4[1][2],FLUKA[22], MCNP[11], among others. The trade-off to this approach is a lossof information about all of the events, and about the trajectories of the low en-ergy electrons (which have an important radiobiological role[64]) but is perfectlyacceptable for calculating the dose deposited in millimeter or larger volumes. Con-densed history algorithms are routinely used for radiation therapy treatment plan-ning, and their ubiquity means that there is a broad set of useful tools and insti-tutional experience available. Systems such as BEAMnrc[82], DOSXYZnrc[101]and GAMOS[4] have been developed to simulate complex geometries like the lin-ear accelerator treatment head and CT-based patient phantoms as well as many ma-terials. To take advantage of this, some groups have used the energy depositionsfrom condensed history codes to measure the lineal energy[83][69]. While someof the events are modelled discretely by this approach, many others are treated aspart of a multi-scattering transport step. This means that the energy deposition for22those events is based on their stopping power and that the energy transferred bythese events is replaced by an average value. This approximation reduces the vari-ance of the events. Other uncertainties are introduced because these histories donot carry out the transport of the subthreshold electrons, which will also effect thecontribution from the events with a large lineal energy, and this can be seen in themicrodosimetric spectra presented in those works. There is also a long history ofusing track structure studies of charged particles[62]. The initial spectra of the par-ticles was very simple, either monoenergetic electrons or the products of interac-tions of uniform photon spectra. This does not accurately reflect the fluence spectrafor clinical radiation therapy, where scatter, energy loss, and energy straggling playa role. Other track study research has attempted to model the radiation damage tothe DNA[29] for uniform monoenergetic electron sources, or by using the sec-ondary electron spectra measured for low energy photon sources[75][96]. Therehas also been research devoted to particle tracks associated with positrons[14],protons[68][72] and heavy ions[67] [35][86][93]. Another approach has been touse a track structure simulation for a very small volume, but even this requiredprohibitive resources[83].However, neither of these algorithms is suitable for track structure simulationsof clinical radiation therapy sources in patients individually. The track structure al-gorithms are too inefficient to be practical for large volumes, and they are generallycapable of simulating interactions only in water. On the other hand, the sources ofthe efficiency of condensed history algorithms are what make them incapable ofperforming track structure simulations. However, these two approaches are com-plementary: a multi-scale approach can be taken where the charged particles incertain volumes are be simulated using the track structure code, while the con-densed history algorithm handles all of the photons and the charged particles in therest of the volume. A series of projects have attempted this in certain limited cases.Most of this has focused on x-ray sources in the keV range[66]. These simulationsare practical to accomplish because the range of the photons and electrons is notvery large. This means it is fine to use a homogenous geometry (filled with wa-ter), and a manageable number of interactions per primary particle. This approachhas been applied to brachytherapy[30][78][79], radiation protection[75], and CTimaging[76]. Another body of research has been devoted to particle tracks associ-23ated with protons[68][72] and heavy ions[67] [35][86][93]. There has also been anattempt to calculate dose distributions (but not microdosimetric quantities) in theneighborhood of metallic nanoparticles for a variety of keV and MV sources[41],using simple spectra and geometry. Another project looked at the microdosime-try of megavoltage beams[54]. However, since the condensed history algorithmsneglect the important low energy electrons, all of these multi-scale simulations donot start the track structure simulations with the correct energy spectra. Section 4.1presents a novel method using a watch volume to ensure that the correct spectrum issimulated. This thesis reports on the first use of track structure Monte Carlo simu-lations to investigate the microdosimetry of clinical therapy beams in patient-basedCT-phantoms and which makes use of the watch volume.In order to perform a multi-scale simulation, each scale must be benchmarkedindependently. Chapter 2 focuses on the use of condensed history Monte Carlo indosimetry and treatment planning. Work done to ensure that the linear acceleratoris accurately modeled is reported, as well as other work relating to the use of MonteCarlo to perform useful tasks for dosimetry and patient specific quality assurance.Chapter 3 discusses the benchmarking and validation of track structure algorithmsfor use in microdosimetry, and gives verification results. Combining these toolsin a multi-scale framework is related in Chapter 4, as well as results from the ap-plication of this framework. The final chapter draws conclusions as well as futuredirections of research.24Chapter 2Condensed History Monte Carlofor DosimetryThe radiological absorbed dose is the current standard measurement for treat-ment planning in radiation therapy and the long history has lead to widespreadfamiliarity[94]. Critically, dose can be readily measured in the clinical setting:ion chambers, thermoluminescent dosimeters (TLDs), solid state devices, radio-sensitive films, and many other methods are all available and are very convenientfor routine usage[56]. These reasons have made dose the fundamental paradigmfor clinical decision making regarding treatment plans.Evaluating the dose distribution is one of the most important aspect of treat-ment planning. Prescriptions are specified in terms of a desired dose to the tumor,and the process of treatment planning is to ensure that this is delivered and thatthere is sufficient sparing of the normal tissue. The current standard is to incor-porate a model of the treatment delivery source with a voxelized, 3-D phantomcreated from a patient CT. The dose is then calculated within each voxel for thetreatment, and then it is analyzed based on several metrics. The region containingall evidence of disease is classified as the Gross Tumor Volume (GTV). Adding re-gions with suspected microscopic disease to the GTV forms the clinical treatmentvolume (CTV). The planning target volume (PTV) is created by adding margins toaccount for uncertainties in patient set up and treatment delivery, and used as theobjective for treatment planning. One important check is the existence of regions25Figure 2.1: An example of a cumulative dose-value histogram (DVH) for aclinical treatment plan. This plan involves two separate planning targetvolumes (or PTVs, in red). The other curves refer to healthy organs atrisk of complications from the radiotherapy.where the planned dose is much higher or lower than the prescription dose in thePTV. The number of these hot or cold spots and/or their size are important to deter-mining the acceptability of proposed treatment. In most cases, the plan evaluationinvolves a cumulative dose-volume histogram (DVH), which is a tabulation of thefraction of voxels within a subvolume of the phantom that receive a dose largerthan certain fixed values. A sample DVH for a clinical treatment plan is shown inFigure 2.1. These are generally evaluated by setting upper limits on what fractionof the volume can receive more than a certain dose. For instance, Dx is the largestdose that is absorbed in at least x% of the PTV, and the fractional volume that re-ceives at least y% of the prescription dose is referred to as Vy. Planning objectivesare site-specific, and are often determined by what is achievable based on previ-ous experience. Dose-based metrics are used to optimize and evaluate complextreatment plans.26For external beam therapy, it is important to ensure that the delivered radiationcorresponds to the plan. One necessity is to check the output from the linac undervarious circumstances to confirm that the planning model is accurate. Another isto ensure that specific plans can be carried out correctly through a comparison ofthe measured and predicted. While this can be achieved by using exit dosimetry(using a film or electronic portal image device, EPID) or by externally placing adosimetry device on a patient, it is more common to perform the check using asimplified geometry (e.g. a cylindrical water phantom) and a few spot checks witha dosimeter. These calculations can be carried out via Monte Carlo simulations,which offer an accurate method for modeling the stochastic processes involved inthe scattering of radiation particles. For dosimetry, condensed history algorithmsare preferred because the measurements and comparisons are for volumes that arelarge relative to the mean free pathlength of the charged particles. In addition,Monte Carlo simulations can be used as an independent check of treatment plansgenerated using a treatment planning system. For IMRT and VMAT techniques,inverse planning optimization is done using a pencil beam-based superposition-convolution method. However, these algorithms can give inaccurate results nearinterfaces between different materials[91], which may affect the calculated dosedistribution. This makes it critical to verify the dose distributions with Monte Carlowhich can properly account for the material inhomogeneity.This chapter covers several different applications of Monte Carlo simulationsfor treatment planning and dosimetric purposes. A linear accelerators is a com-plex devices, and careful evaluation must be done to ensure that a Monte Carlomodel of it is acceptable. The commissioning and testing of a clinical accelera-tor is related in Section 2.1. In addition, the use of VMAT plans requires that thedynamic motion of many parts of the linear accelerator and potentially couch andpatient motion can be included. Incorporating this motion into the simulation isdiscussed in Section 2.2. Using these tools, Section 2.3 presents a comparison be-tween the superposition-convolution algorithm AA and Monte Carlo simulationsfor spine stereotactic body radiation therapy. Section 2.4 presents a discussion ofdifferent methods for reporting doses. The application of a novel method to plantotal body irradiation, and a Monte Carlo investigation into that method is the focusof Section Modeling a Linear AcceleratorAccurate simulations require accurate models of the components involved. Whilethe general design of clinical linear accelerators is well known, the precise detailsare less well established. This is partially due to the proprietary design informa-tion which the manufacturers do not wish to make public, as well as variations inthe installation, updates, and modifications to specific machines. Thus, there isno canonical model of the linear accelerator for Monte Carlo simulations in theacademic or clinical realms. There is particular disagreement regarding the mate-rial properties in different parts of the accelerator, particularly the actual design ofthe flattening filter. Another source of uncertainty for photon beams is the energyspectrum and beam shape for the electron beam when it hits the bremsstrahlungtarget. Consequently, it is necessary to commission a Monte Carlo model againstmeasured data to demonstrate the efficacy of the simulations. A new linac modelfor BEAMnrc[82] was provided by J. Siebers for 6 MV photons. Compared tothe model currently in use at the BC Cancer Agency, this new model has differentdimensions for the jaws and flattening filter, as well as a new set of material com-positions for several components. To ensure that the most accurate of these twomodels was put into use, the author compared the beam characteristics calculatedthrough simulations wtih both models with the experimentally measured values.Given a model of a linear accelerator, it is possible to derive estimates of theuncertain or unknown parameters by a comparison of measured and simulated data.The standard approach to commissioning a Monte Carlo model of a linac for pho-ton beams is discussed by Sheikh-Bagheri and Walters[89][88]. The first step is toestablish the energy spectrum of the electron beam incident on the bremsstrahlungtarget, which can be determined by matching the shape of the percent depth dose(PDD) along the beam axis with the experimental curve. The coordinate systemis represented in Figure 2.2. For this work, simulations were performed usingDOSXYZnrc[101] which called BEAMnrc[82] shared libraries to generate parti-cles. For each set of parameters, 5?107 independent histories followed through theDOSXYZnrc geometry, which was a cubic water phantom (50 cm on each side).For the PDD studies, the field dimensions were 10 cm by 10 cm at the isocenterplane, which was placed at a depth of 1.5 cm in the phantom and an SAD of 10028cm. To increase efficiency, a column of cubic voxels was centered along the beamaxis, each voxel having a side length of 0.25 cm, while the rest of the geometrywas divided into larger voxels. A wide range of combinations of mean energy andenergy full-width half-maximum (FWHM) of the beam were tested. The best fitwas obtained for a mean energy of 6.0 MeV, with a FWHM of 0.05 MeV. This iscompared to the measured values for a 6 MV beam from a Varian Clinac21iX, aswell as the curve created for the current beam model in Figure 2.3. The currentmodel uses a monoenergetic 6.0 MeV electron source. Since the mean energiesare the same, the simulated curves agree well at depths beyond the maximum, asshould be expected. The difference in the build up region can partially be attributedto the lower energy part of the proposed model?s source electron spectrum. On thewhole, the simulations both capture the PDD everywhere outside of the build-upregion.The next stage of the commissioning of a linear accelerator model is to opti-mize the shape of the initial electron beam just before it strikes the bremsstrahlungtarget. In this approach, the initial electron beam incident on the bremsstrahlungtarget is assumed to be elliptical in shape and normally distributed about the centralaxis[46]. The effect of these parameters will be easiest to observe in the shape ofthe dose distributions off the beam axis. To evaluate these distributions, the dosewas recorded along lines passing through the phantom parallel to the x- and y-axesat different depths in the water phantom (see Figure 2.2). These curves are thennormalized relative to the maximum dose along the beam axis in the phantom andare referred to as beam profiles. The simulated profiles were obtained for differentfield sizes. The phantom was divided into cubic voxels that were 0.25 cm on a sideand larger voxels were used in the rest of the phantom. The profiles were comparedto each other and to the measured values. The optimal FWHM maximum in boththe x and y directions was found to be 0.113 cm for the proposed model, comparedto 0.075 cm in the current one. For these parameters, the profiles at a depth of 10cm for a 10 cm by 10 cm field in the x (in-plane) direction are shown in Figure 2.4and Figure 2.5, and for the y (cross-plane) direction in Figure 2.6 and Figure 2.7.The simulations for both models agree well with the data. On the other hand, it isdifficult to compare the penumbra measurements to the simulated values becausethe gradient in that region is large compared to the size of the voxels, but the two29Figure 2.2: A schematic representation presenting the definition of the coor-dinate system used in the modelling of a linear accelerator. In this case,the source represents the focal point of the photon beam being studied,and is the origin of the coordinates. The direction of the beam is towardsthe water phantom, and the z-axis is the central axis of this beam. ThePercentage Depth Dose (PDD) is measured along this axis, where thedepth is measured relative to the top surface of the water phantom. Thebeam passses through a multi-leaf collimator (MLC). The direction ofthe leaf travel defines the x-axis (alternatively, the ?in-plane? direction),and the perpendicular axis is the y-axis, which is also referred to asthe ?cross-plane? direction. The beam profiles are dose measurementstaken along a line perpendicular to the beam axis in the phantom, andare normalized to the maximum dose measured along the beam axis.300 5 10 15 20 25 300102030405060708090100Depth (cm)Percent Depth Dose (%)  MeasuredCurrentProposedFigure 2.3: A comparison of the Percent Depth Dose curves for the currentand proposed linear accelerator models and the measured values from a6MV beam from a Clinac 21iX, for a 10 cm by 10 cm field in a waterphantom.models agree with each other quite well. In addition, the presence of horns thatare more pronounced than in the measured profiles is observable in the new model(but not the current model) in these figures. This effect in increased for larger fieldsizes, becoming a 2 percentage point difference for a 40 cm by 40 cm field.In the final analysis, the two models agreed quite well in most circumstances,and both fit well within the acceptable treatment planning tolerances[47]. The useof either model would be justified. However, since the newly proposed modelproduced the anomalous horns, it was determined that the current method betterrepresented the beams in use at the Vancouver Cancer Center, and remain the basisof the current simulation of Clinac 21iX accelerators. Recently, the linac manu-facturer Varian has decided to protect their proprietary information by providing a31?30 ?20 ?10 0 10 20 30010203040506070Position (cm)Percent Maximum Dose (%)  MeasuredCurrentProposedFigure 2.4: A comparison of the dose in-plane axis beam profiles for bothlinear accelerator models as well as the values measured on VarianClinac21iX along the in-plane axis (the x-axis) at a depth of 10 cm in awater phantom for a 10 cm by 10 cm field.phase space description file captured above the jaws rather than providing detailedspecifications for the new TrueBeam design. Since there are no parameters that areadjustable during treatments and little variation between units for this part of thetreatment head, this can allow a more accurate model to be used without tuning anyparameters. On the other hand, this method makes it difficult to model variationsin the back scatter from the jaws and MLC, which contributes to the total exposurein the ion chamber, and potentially to the patient. While this means that it is nolonger possible nor necessary to adjust the electron beam parameters, commission-ing must still be performed to check the jaw and MLC models, and to ensure thatthe provided phase spaces do agree with the actual particle fluence for the unit.32?6 ?4 ?2 0 2 4 6585960616263646566676869Position (cm)Percent Maximum Dose (%)  MeasuredCurrentProposedFigure 2.5: A comparison of the in-plane beam profiles in the field region forboth linear accelerator models as well as the values measured on VarianClinac21iX along the in-plane axis (the x-axis) at a depth of 10 cm in awater phantom for a 10 cm by 10 cm field.2.2 VMAT QAOnce the linear accelerator model has been commissioned, Monte Carlo simula-tions for patient-specific QA can be done. This process requires information aboutthe settings of the linear accelerator, MLCs, and patient anatomy that come fromthe treatment planning system. This information must be parsed into a form that isunderstandable by the Monte Carlo code. For VMAT plans, the continuous rotationof the linear accelerator coupled with the dynamic motion of the jaws and MLCs re-quired a major update to parts of the BEAMnrc and DOSXYZnrc programs that arewidely used for radiation therapy simulations. These changes were implementedby Lobo and Popescu[55]. However, it was still necessary to develop a set of toolsto convert the output of the treatment planning software (DICOM files) into the33?30 ?20 ?10 0 10 20 30010203040506070Position (cm)Percent Maximum Dose (%)  MeasuredCurrentProposedFigure 2.6: A comparison of the cross-plane beam profiles for both linearaccelerator models as well as the values measured on Varian Clinac21iXalong the cross-plane axis (the y-axis) at a depth of 10 cm in a waterphantom for a 10 cm by 10 cm field.proper input for the Monte Carlo simulation, and this required a major overhaulof the existing suite of tools in use. Some of these changes simply involved pars-ing the DICOM file and transferring that data to the input file to the simulations;however, much of it involved significant additional processing and analysis.First of all, BEAMnrc and DOSXYZnrc use different coordinate systems, andif the motion of the system is taken into account, the transformation between thembecomes dynamic. The beam source models for DOSXYZnrc created by Lobo andPopescu for DOSXYZnrc allow the simulation to have rotation of the gantry by anangle ?gantry, the couch by ?couch, and the collimator angle ?coll , as well as thetranslation of the isocenter (whether due to couch or patient motion) plus changesto the source to isocenter distance. Since this data is stored for certain times during34?6 ?4 ?2 0 2 4 6585960616263646566676869Position (cm)Percent Maximum Dose (%)  MeasuredCurrentProposedFigure 2.7: A comparison of the cross-plane beam profiles showing the fieldregion for both linear accelerator models as well as the values measuredon Varian Clinac21iX along the cross-plane axis (the y-axis) at a depthof 10 cm in a water phantom for a 10 cm by 10 cm field.the delivery of the beam known as control points, it is necessary to interpolatebetween the control points for most times during the delivery. The instantaneousrotation angles are given by:? = cos?1(sin?couch sin?gantry)(2.1)? = 90.0+ tan?1(cos?couch sin?gantrycos?gantry)(2.2)?coll = 90.0+ tan?1(?sin?couch cos?gantry cos?coll ? cos?couch sin?collsin?couch cos?gantry sin?coll ? cos?couch cos?coll)(2.3)where ? is the polar angle with respect to the patient?s posterior-anterior (PA) di-35rection, ? is the azimuthal angle with respect to counterclockwise rotation in theinferior-superior direction, and ?coll is the angle of rotation about the beam axis.It was also crucial to produce a tool that can create a customized phantombased on the patient CT data. While there existed some tools to do this, theydid not account for changes to the transformations caused by the patient beingin any orientation other than head-first supine, nor could it account for differentconventions for encoding the CT numbers in the DICOM file. In addition, it isoften useful to be able to overwrite the density in the CT with other values eitherinside or outside of a specified contour. This need drove the development of acode capable of manipulating the CT files prior to the construction of the phantom.This tool, known as ICTC was implemented and is currently used for all patient-specific external beam Monte Carlo simulations at the Vancouver Cancer Centre.A illustration of the phantoms created by ICTC, as well as the original in Eclipse,are shown in Figure Comparison Between AAA and Monte Carlo forSpine SABREclipse is a treatment planning system distributed by Varian for use in conjunc-tion with their linear accelerators. For dose calculations, Eclipse makes use of aconvolution-superposition algorithm called AAA (for anisotropic analytical algo-rithm). This algorithm calculates the dose distribution by doing a weighted summa-tion over the contributions from many narrow pencil beams. The lateral scatteringis accounted by using a kernel that is scaled according to the electron density of thematerial in that direction[10]. This approximation relies on assuming that particlefluence is in equilibrium in the area of the pencil beam, which can break downin the vicinity of inhomogeneous materials. This is a particular concern for spinestereotactic ablative radiotherapy (SABR), which is the delivery of a highly con-formal treatment in a small number of high dose fractions. The target volume issurrounded by bone, soft tissue, lung, and spinal fluid, all of which have differentdensities and radiological properties due to the atomic composition, and presents achallenge for AAA. The high conformality of the treatment is required because ofthe presence of the radiosensitivity of the spinal cord and lung, which is magnified36Figure 2.8: Sample CT based voxelized phantoms produceed by (a) Eclipse,(b) iCTC and (c) iCTC where the lung has been overwritten with a den-sity of 2.0 g/cm3 for demonstration purposes.37Patient ID MC %PTV covered AAA %PTV covered100% 90% 80% 100% 90% 80%Spine 1 70.3 86.8 90.1 79.5 88.5 90.7Spine 2 70.5 90.4 93.0 79.9 91.4 93.4Spine 3 92.9 99.8 100 98.9 99.9 100.0Spine 4 90.9 97.9 98.9 95.0 98.3 100.0Spine 5 61.8 86.8 91.1 75.3 88.9 92.2Spine 6 87.3 99.1 100 95.0 100 100.0Table 2.1: The percent of the planning treatment volume receiving certaindoses calculated by the treatment planning system?s algorithm (AAA)and Monte Carlo for the Spine SABR patients.by the high biological equivalent dose associated with the large fraction sizes.Before the adoption of the stereotactic method, it was very important to un-derstand the limitations of the treatment planning process to ensure that acceptableplans are created. A retrospective planning study was carried out for treatmentsof 35 Gy in 5 fractions using VMAT[63]. Seven patients were initially selected forthe study, but one of the patients was excluded because of the use cement to re-inforce the spine after a previous course of radiation therapy. The tools discussedin Section 2.2 were used to prepare the plans for Monte Carlo simulations. Thepurpose was to see how well the plans generated by the AAA algorithm met theplanning objectives when the dose was calculated with Monte Carlo. The treatmentplan design and optimization was performed by Eclipse and the simulations werecarried out using DOSXYZnrc[101] using the sources of Lobo and Popescu[55].According to the commercial treatment planning system, the goal of 80% cover-age of the planning target volume (PTV) with the prescription dose was achievedto within 0.5 percentage points for 5 of the 6 patients, while the Monte Carlo re-sults show that all of the patients actually had less satisfactory coverage. Of those6 plans that did meet, or nearly met, the requirements according to the treatmentplanning system, the plans actually fell short of the objective by between 4 and13.5 percentage points. The results for these patients are shown in Table 2.1.The DVHs for one of these patients is shown in Figure 2.9. More details aboutthis discrepancy are presented in Figure 2.10. Comparing the AAA calculated dis-tributions in Figure 2.10(a) to the Monte Carlo ones in Figure 2.10(b), there is a38significant disagreement in the PTV coverage. In particular, the coverage on the(patient?s) right side of the PTV is much poorer than the AAA calculation. Inaddition, the dose is lower along the left bone-tissue interface in the PTV. Theseproblems both occur within bone, and the isodose lines are in better agreement inthe soft-tissue parts of the PTV.0 20 40 60 80 100 1200102030405060708090100Percent Prescription Dose [%]Ratio of Total Structure Volume [%]  AAAMCMC Density OnlyMC DwFigure 2.9: Spine SABR PTV DVHs by using different inhomogenity correc-tions. The dose distributions were calculated with AAA, Monte Carlosimulations (MC) with full inhomogeneity correction, MC simulationsthat only corrects for density heterogeneity, and MC simulations forwhich the dose to water is reported (Dw).It is important to investigate the source of these differences. The first step wasto check results when the material is homogeneous. By using a phantom that in uni-formly filled with water in both cases, it is possible to directly compare the linearaccelerator models as well as the interpretation of the beam settings to see if theyare in agreement without the added complication of the different ways of handling39Figure 2.10: Dose distributions for a SpineSABR plan calculated by differentmethods and inhomogeinity corrections by Monte Carlo simulationsusing DOSXYZnrc[101] for the slice containing the isocenter. Thecolor wash is represents the planning treatment volume, the yellow andlight blue contours outline the spinal cord and canal. The orange, darkblue, and green contours represent the isodose lines for 50%, 90%, and100% of the prescription dose. The AAA distribution is shown in (a)and compared to Monte Carlo simulations with (b) full inhomogeneity,(c) only density inhomogeneity, and (d) dose-to-water[90] corrections.the heterogeneity within the phantom. These dose calculations agreed very well,as can be seen in the sample comparison shown in Figure 2.11 and Figure 2.12.While reassuring, this is not a surprising result because of the exhaustive cross-validation that has been performed between both dose calculation algorithms forwater phantoms.The next step was to confirm that the AAA algorithm only takes electrondensity but not atomic composition into account. A simple phantom was con-structed from SolidWater (a commercial solid material that is radiologically water-40Figure 2.11: Dose distributions for a SpineSABR plan calculated by (a) AAAand (b) Monte Carlo simulations using DOSXYZnrc[101] for phan-toms filled with uniform density water, that is with no heterogeinitycorrections. The red-volume represents the PTV, while the yellow andlight blue contours delineate the spinal cord and canal. The orange,dark blue, green, and magenta contours represent the 50%, 85%, 90%,and 100% prescription doses. Note that the overall normalization hereis to 25 Gy, although the plan is otherwise identical to the other plansdiscussed.equivalent) and bone-equivalent material, as shown in Figure 2.13. Three filmswere positioned inside the phantom - at the bone-water interfaces and in the mid-dle of the bone-equivalent material, and the doses were calculated at these pointsusing the AAA algorithm and Monte Carlo with full inhomogeneity correctionsfor a 10 cm by 10 cm field from a 6 MV beam. In addition, to compare with theAAA algorithm, the Monte Carlo was performed on a model in which the materialwas water-equivalent but the density of each voxel was determined based on theCT number. The relative doses are shown in Table 2.2, where the normalization isto the dose at the upper bone-water interface for each method. The results for theMonte Carlo with only density corrections agrees exactly with the AAA correction,but not with the other measurements. This shows that the Eclipse AAA algorithmdoes not consider the atomic composition.Likewise, the Spine SABR patient Monte Carlo phantoms were remade to con-sist of only water, but the density of that water was determined by the CT number.This means that the results only account for inhomogeneity in the density but not41Figure 2.12: Dose-Volume histograms created with no heterogeinity correc-tions for a spine SABR plan calculated by AAA (squares) and MonteCarlo simulations using DOSXYZnrc[101] (triangles) for phantomsfilled with uniform density water. The magenta curves are for the PTV,while the lighter blue line is for the spinal canal, and the slightly darkerblue color represents the spinal cord.in the atomic composition. These calculations showed very strong agreement withthe AAA values. For the sample plan being considered here, Figure 2.10(c) showsa very close match with the AAA in (a). Specifically, there is no major discrepancyat right side of the PTV or the left bone-tissue interface for the 100% prescrip-tion isodose line. Overall, the DVH (in Figure 2.9) calculated with the densitybut not material inhomogeneity is in better agreement with the AAA method thanthe Monte Carlo simulation with full heterogeneity corrections. This demonstratesthat part of the disagreement is related to the handling of different material com-positions. This result may imply that there is a straightforward way to improvesuperposition-convolution algorithms by simply accounting for the atomic make-42Location Measured Full MC Density-Only MC AAAUpper Bone-Water Interface 1 1 1 1Center of Bone 0.90 0.90 0.93 0.93Lower Bone-Water Interface 0.88 0.83 0.85 0.85Table 2.2: A comparison of the relative doses measured or calculated in aphantom constructed from bone- and water-equivalent materials at sev-eral different depths in the phantom. The doses are normalized to thedose at the upper bone-water interface. Full MC refers to theMonte Carlocalculation with full material inhomogeneity considered, while Density-Only Monte Carlo refers to simulations done with only density inhomo-geneity considered (the material was considered water-equivalent). TheAAA calculation was done by the Eclipse treatment planning system.Figure 2.13: A diagram of the geometry of test phantom used to investigatethe inhomogenity corrections used by the AAA algorithm.43up of media in the simulation. Since they already account for electron density,it would not be difficult in principle to use add an atomic composition correctionto the kernel. Using the density in the CT to assign the material already a well-established practice in Monte Carlo: this is done by ctcreate, the phantom creationtool associated with DOSXYZnrc. This approach would not be strictly correct,because the energy of the particles determines the differences in material response,but it could lead to more accurate models. These functions could be calculatedfrom existing cross section data, and would be open to validation by both experi-ments and Monte Carlo simulations.It was also important to be able to present the results in terms of the ab-sorbed dose to water instead of the actual medium calculated, in order to compareto historic data and treatment planning algorithms. The approach of Siebers[90]was adopted using pre-tabulated stopping power ratios based on the medium inwhich the Monte Carlo dose was computed. It has been recommended that a largenumber of materials (up to 55) be used in the conversion[23][19]. However, theDOSXYZnrc program, and the tools to create phantoms for it assumes that eachmaterial is indexed by a single digit integer. This means that regardless of the valuechosen for the MXMED user variable, the distribution version of the programs willcause errors if more than 10 materials are employed. To overcome this limitation,substantial changes had to be carried out to the routines responsible for the read-ing and writing of the phantom files in these programs by the author. In addition,the author also had to make modifications to the Dose-To-Water code (shared by J.Siebers) that performed the conversion. These changes were made in order to makethe Dose-To-Water code compatible with the file types used at the Vancouver Can-cer Centre. This approach was applied to the patients in the study. For the sampleplan discussed above, we show this conversion in Figure 2.9 and Figure 2.10(d).For this plan, the agreement the with the AAA method is slightly improved in thebone in the right side of the PTV, but the interface troubles on the left side are stillunresolved. A more thorough discussion of this approach is given in Section 2.4.In addition, the main constraint on spinal radiosurgery is the maximum dose tothe spinal chord because the complications associated with hot spots in this organ-of-risk can have significant quality of life impacts. One of the reasons for consider-ing spine SABR is that it can potentially improve sparing of the spinal chord, which44Patient ID Monte Carlo AAAD1.2cc [Gy] Dmax [Gy] D1.2cc [Gy] Dmax [Gy]Spine 1 11.0 14.1 11.6 15.4Spine 2 14.6 20.1 15.5 21.4Spine 3 10.7 15.3 11.5 16.5Spine 4 14.8 16.8 15.5 16.9Spine 5 12.4 15.8 13.5 16.5Spine 6 13.7 15.5 14.5 16.9Table 2.3: The spinal canal doses calculated for Spine SABR treatments. Themaximum dose and the largest dose contained in a 1.2 cc volume of thespinal canal as caluclated by the treatment planning systems AAA algo-rithm and Monte Carlo with full inhomogeneity corrections for the SpineSABR patients.could be beneficial for patients with local recurrence after previous radiotherapytreatments. However, this means that tighter constraints must be placed on the max-imum dose to the spinal canal (in this study it was 18 Gy compared to 22 Gy forthe patient who had not been previously irradiated). Other dose limits were placedon the maximum dose to the skin, great vessels, and large and small bowels, whilehot spot limits were placed on the lung, kidney and liver. The Monte Carlo sim-ulations revealed no unacceptable overdosing to any of these structures, as shownin Table 2.3. In fact, for the spinal canal, the maximum doses and hot spots wereconsistently smaller when calculated by Monte Carlo methods compared to theAAA values. A comparison of the doses to the spinal canal calculated by AAAand Monte Carlo for the sample patient discussed above is shown in Figure 2.14.This shows that the AAA method may not cause unexpected over-dosing of theorgans-at-risk.Taking these results as a whole, there should be some concern about PTVunder-coverage for plans generated by the Eclipse treatment planning system forspine SABR. On the other hand, perhaps the potential for overdosing of the organs-at-risk is not as large. This treatment technique does indeed seem feasible, althoughit would appear necessary to do a Monte Carlo check of all plans before they areaccepted. Furthermore, it is demonstrated that it is possible to get as ?fair? compar-ison of the dose calculations done by AAA with those of Monte Carlo using only450 2 4 6 8 100102030405060708090100Dose (Gy)Ratio of Total Structure Volume [%]  AAAMCMC Density OnlyMC DwFigure 2.14: Spine SABR DVH!s for the spinal canal calculated by differentinhomogeinity corrections: with AAA (blue), Monte Carlo with fullinhomogeneity corrections (red), Monte Carlo simulations with onlydensity corrections (black) and Monte Carlo simulations with full in-homogeneity corrections converted to dose to water (magenta). Thedose limit for this patient was 18 Gy.the density inhomogeneity correction. Outside of the scope of this study, since theorgan-at-risk constraints did not present a major optimization challenge, it wouldseem reasonable to consider the potential for further dose escalation when usingthis technique.2.4 Dose To WaterThere is a healthy debate about how to report the doses calculated by Monte Carlomethods in patient CT-phantoms. In Monte Carlo simulations, the dose in eachvoxel is calculated based on the medium that occupies that voxel, Dm. On the46other hand, the historical development of treatment planning methods was builton calculating the dose in a water-equivalent region (Dw). This meant that onlyone conversion factor was needed to compare calculated values with ion chambermeasurements and is a legacy of the practice pre-dating the use of CT data fortreatment planning. A consequence of this approach is that the tools and metrics oftreatment planning, such as NTCP and TCP, are based on the dose-to-water values,and some believe that this is a reason that Monte Carlo doses should actually beconverted to Dw[90].In order to make the conversion from Dm to Dw, Bragg-Gray theory is applied.This theory states that the ratio of the doses at a single point is given byDwDm=Emax?0?(E)(S(E)/?)wdEEmax?0?(E)(S(E)/?)mdE(2.4)where Emax is the largest energy of the electrons passing through the voxel, ?(E)is the electron energy fluence at the point, (S(E)/?)m is the mass stopping powerfor an electron of energy E in the medium, and (S(E)/?)w is the mass stoppingpower for an electron in water. This calculation is accurate only for volumes thatare small enough that the it is acceptable to neglect changes in the electron stoppingpower while the electron traverses the volume, and also that the electron fluence isnot altered by the different material that makes up the volume. This approximationwas applied by Siebers, et al.[90], to convert the dose in a voxel to a dose-to-water. Since the voxels used in the Monte Carlo simulations almost always havedimensions greater than 1 mm, this means that the small volume assumption fails.To demonstrate this, consider a 1 MeV electron passing traveling through water.The change in the stopping power of the electron as it moves through a mediumis shown in Table 2.4. These values are calculated based on the NIST ESTARtables[7], where the change in electron energy is taken from the unrestricted stop-ping power, and the electron trajectory divided into 103 discrete steps of equalsize. These results show that the electron fluence does change somewhat non-linearly over the volumes of interest, and therefore the commonly employed ideaof using the midpoint energy as the electrons contribution to the electron fluence47Path Length (mm) Stopping Power (MeVcm2/g) Change Energy (MeV)1 0.21 3% 0.9792.5 0.20 6% 0.9375 0.18 12% 0.90110 0.16 23% 0.814Table 2.4: The percent change in the total stopping power at different pathlength. Specifically, its percent change at different path lengths for a 1MeV electron in water, based on the stopping powers given in [7], as wellas the kinetic energy at those point.is not completely justified. In the given approach, the stopping power ratio is as-sumed to be the same everywhere in the phantom, which neglects the effects ofbeam hardening, lateral scatter, and disequilibrium. However, it was demonstratedthat stopping power ratio for monoenergetic electrons between 1 MeV and 10 keVwas nearly uniform for all of the materials considered, and there was only a slightenergy dependence above that for most of the materials[90]. Furthermore, thisstudy showed that electrons in this energy range dominated the spectrum both nearthe surface and deep in the phantom for the sources of interest in radiotherapy.Further measurements and simulations confirmed that there was little variation inthe stopping power ratios for low dose regions (that is, out of the field) of waterphantoms[39]. Ultimately, these results imply that the approximations used in thisconversion method do not cause meaningul inaccuracies.There has been discussion about the number of materials to employ in the dose-to-water conversions[23][19]. One of the reasons is to prevent interfaces with dras-tic changes in the material properties and to prevent significant electronic disequi-librium at these points. On the other hand, using a large number of materials is nottoo important in practice, because these interfaces are small relative to the voxelsizes in the simulations. Since the medium is considered uniform within the voxel,this gradual change in the conversion factor is not realized in practice. As an ex-ample, the interface effects are still present in the Dw conversion in Figure 2.10(d).In the final analysis, however, the fact that this approach is valid is not moti-vation for its use. The relevant fact is the dose delivered to the patient, which isthe dose to the medium. This can be seen as an overriding reason to present all re-48sults as dose to medium, reversing the proposed paradigm. Even the dose-responsecurves, while nominally in terms of dose to water, ultimately reflect the effects ofthe dose to medium. In addition, if a comparison with other dose calculations isrequired, a more sophisticated approach should be taken. For instance, the EclipseAAA algorithm takes density but not atomic composition into consideration (seeSection 2.3), so it is better to make the same assumption in both calculations if wewant to make a comparison. Ultimately, however, the goal of this research is toimprove treatment outcomes, which can best be achieved by accurate modeling ofthe radiation exposure. Therefore, it is usually best to not artificially manipulatethe different methods to agree, but to look at the most accurate results generatedby each approach.2.5 Total Body IrradiationPatient specific planning and QA is also used for total body irradiation (TBI)[77].The goal of this treatment method is to deliver a uniform dose throughout the pa-tient, in order to sterilize the bone marrow prior to a transplant, and is most com-monly used to treat leukemia. At the Vancouver Cancer Centre, this treatment isdelivered through a 60Co sweeping radio-isotope source. The patient is placed onthe floor to create an extended source-to-patient-surface distance (SSD) close to 2m (as opposed to the unit?s SAD of 80 cm). The dose is delivered in three frac-tions via two fields: an anterior-posterior (AP) field with the patient prone, and asupine field (PA). The exposure time is determined by the patient separation at theumbilicus, the treatment SSD, the sweeping rate, and the current source activity.A divergence compensator is employed to correct for the inverse square law ef-fects on the extremities of the patient ? since the 0-degree beam axis is set to passthrough the umbilicus, the superior (head) and inferior (feet) ends of the patient arefarther from the source than the patients abdomen.Another potential source of inhomogeneity of dose is the lungs of the patient.Since the densities inside the lungs are less than that of soft tissue, the lungs areless effective at attenuating the beam than other parts of the body. If uncorrected,this would lead to overdosing in the tissue underneath the lung. Furthermore, thelateral electronic disequilibrium at the lung/tissue interfaces can contribute to in-49creased doses at the boundary of the lungs or the chest wall. To prevent this, com-pensator blocks are used to reduce the beam fluence through the lung. The shapeand thickness of these blocks must be customized for both of the fields for eachpatient. A model of the treatment set up is shown in Figure 2.15.Traditionally, the block shape and thickness was based on exit dosimetry: ra-diographic films were placed under the patient during short pre-treatment expo-sures at the 60Co unit. The oncologist then used these films to outline the projec-tion of the volume to be blocked. These outlines were then projected to the planeof the tray that holds the block and the lead blocks were cut. The average tissuedeficit of the lungs relative to the umbilicus was then computed using the opticaldensity measurements from the film. Figure 2.17(b) shows a typical film. The out-line is drawn, and the small circles inside of that outline represent the points usedto determine the average tissue deficit.The impending retirement of the film processor demanded a different approachbased on CT images. The commercial treatment planning system, Eclipse, can pro-duce digitally reconstructed radiographs similar to the films using CTs, but couldnot place the patient phantom at an extended SSD. This meant that the rays aretraced incorrectly through the origin, as the angle between the zero-degree axis andeach point is increased for a shorter radius, which is demonstrated in Figure 2.16.This error increases the chord length through the voxel, and increases the water-equivalent thickness calculated by the algorithm. In addition, this effect widensthe boundary between the lung and normal tissue which can lead to difficulty whentrying to outline the lungs. The solution was to create a custom built software toolto create the proper images and calculate the tissue deficit which was created bythe author of this thesis.The end result of this was an application called RunTBI. The first step was thecreation of a patient CT-based phantom, and was the impetus for the creation ofthe ICTC toolset which was later applied to VMAT QA (discussed in Section 2.2).The 2D images were created by calculating the total mass along rays originatingat the source location and ending on the film plane. A sample image is shownin Figure 2.17(a). The proper scaling was done using radio-opaque ball bearingsplaced on the surface of the patient during the CT simulation (which are repre-sented by crosses in the images), and these marks were used to ensure proper50Figure 2.15: The treatment set-up geometry for total body irradiation (TBI),including a patient phantom (in green and blue), the lung compen-sator blocks (in yellow), and the divergence filter (in purple). Theblock tray is not shown. The y axis represents the beam axis at 0 de-grees. This was created using the egsview application distributed withEGSnrc[43].51Figure 2.16: A representation of the projection problems of the treatmentplanning system (TPS!) for total body irradiation (TBI). A line is drawnfrom the origin to a fixed corner of the rectangle. For a larger dis-tance between the origin and the rectangle, the chord length of the linethrough the rectangle is reduced. This demonstrates the source of inac-curacy of the DRR image generated by the treatment planning system(a) compared to the actual patient set-up (b).alignment during patient set up. These images were used by the physician to out-line the blocks. To make images that were more familiar to the physicians, thecolor scale of the images was set to match the original calibration curve for theradiographic film. The outlines were then taken by the RunTBI application andthe average water equivalent thickness inside the outlines was calculated. In addi-tion, RunTBI computed the patient separation between the posterior and anteriorsurfaces at the umbilicus based on the body contours generated by Eclipse, fromwhich the lead-equivalent tissue deficit was computed. This value was rounded tothe nearest half-millimeter, and was used to select the appropriate thickness of leadsheet from which to cut the blocks. The user interface for the program is present inFigure 2.18.52Figure 2.17: A comparison of the DRR! images used in TBI produced by (a)RunTBI and (b) radiographic film for the same patient. For the ra-diographic film, the outline of the blocks is shown. The circles in (b)represent the points at which the optical density was measured in orderto obtain the tissue deficit of the lungs.The first test of this system came from a simple phantom. A set of solid wa-ter blocks was arranged to surround a central cavity filled with air, as shown inFigure 2.19. This phantom was then scanned on the CT simulator and a plan wascreated using RunTBI. Then the phantom was irradiated by the 60Co unit to pro-duce an image. Since the cavity was a rectangular prism, and the phantom wasotherwise cubic, it was possible to calculate the average tissue deficit mathemati-cally. The results for all three approaches are shown below in Table 2.5. There isstrong consistency for the mathematical and RunTBI values, while the film methoddoes not agree. Since the film processor was nearing the end of its service lifetime,it is natural to expect that it may be less reliable. To test this hypothesis, a cali-bration film was taken using a linear accelerator to deposit several fixed amountsof exposure in the dose range of interest for this study. These showed that thefilm was registering more dose than the calibration curve would predict, and the53Figure 2.18: A screenshot of the RunTBI user interface.images looked underdeveloped. After creating a new calibration curve, the tissuedeficit calculated for this phantom agreed better with the mathematical result, butwas still not consistent. Another check was performed using the RANDO anthro-pomorphic phantom. This phantom was used in the original commissioning of thetechnique, and this data agreed to within 0.5 mm (the minimum increment of thelead thicknesses) with the RunTBI method. Redoing the film-based calculation,the computed block thickness once again compared unfavorably with the other ap-proaches.The final validation came with a set of 14 patients for which blocks were de-signed using both the film method and RunTBI. The thickness of the blocks werecompared for both fields (although from 3 of the patients, the CT scans for thePA field did not follow the correct protocol and were not used). If the thicknessrounded to the nearest half-mm, then they judged to agree. By this criterion, 56%54Figure 2.19: A cross-section of the RunTBI validation phantom. A rectangu-lar air cavity was constructed around the center of a larger rectangularphantom made of water-equivalent material (SolidWater).55Method Block ThicknessExpected 0.75 mmRunTBI 0.72 mmFilm 0.97 mmTable 2.5: A comparison of the calculated TBI block thicknesses for a simplephantom, calculated by three different methods.of the thicknesses did agree. For another 20% the fields, the thickness calculatedusing the film method rounded to a block thickness that was 0.5 mm larger thanRunTBI?s value. The rest showed a larger difference, and most of the cases thatfell into this category were taken near the end of the study. It appears that thecalibration was rapidly drifting at the end of its life, and the image quality wasnoticeably degraded for the later patients. Estimates made based on the DRRs andother measures calculated by Eclipse were in agreement with the RunTBI values.In the TBI technique, the planning objective is to deliver a uniform dose to themidplane. In order to investigate the variation in the doses to the lung between thetwo methods, a Monte Carlo simulation was performed, using the geometry shownin Figure 2.15. This simulation was carried out by a custom-designed EGSnrc++user code[44], which scored dose only in the voxels along the midplane that wereunder the projection of the compensator block. The doses inside the contouredvolumes were calculated for each of the fields, and for a 4 cm2 area under the um-bilicus, which was used for the normalization for each field. Because there is nosmooth transition between the two fields, there is no reliable way to map the de-formation between the two patient positions unless better resolution is available inthe CTs, which was not done because it would require larger doses to the patientduring the planning stage. Therefore, it is not feasible to match voxels betweenthe fields, and so the dose distributions from each field must be considered sep-arately. The cumulative dose-volume histograms for a typical patient are shownin Figure 2.20 for each field. Simulations were performed for two different lungblock thicknesses: one for the thickness calculated by the original film method,and another for the thickness from the RunTBI program, as well as with no blocks.A total of 2? 108 particle histories were performed for an isotropic 60Co pointsource. The variance (compared to the dose under the umbilicus) was 2.1% for the56Figure 2.20: The dose-volume histograms for the TBI treatment compared byusing only the the voxels inside of the lung on the midplane for totalbody irradiation for both the (a) posterior-anterior and (b) anterior-posterior fields of a patient calculated using Monte Carlo. The bluecurve represents the simulation done with the thickness calculated bythe RunTBI program while the red line is for the thickness calculatedby the film method.RunTBI-derived thickness as opposed to 2.4% for the film-based thickness, and3.1% if no block was used. Using F-statistic analysis, this difference between thetwo methods was statistically significant at a 99% confidence level. This impliesthat a better dose homogeneity was achieved.The RunTBI program improved the physician?s ability to delineate the shape ofthe lung blocks and also created a more accurate and reliable tool for determiningthe thickness of the blocks. One issue of concern was the fact that the old films werecreated using the 60Co beam, while the RunTBI produced images were based on kVimaging x-rays. It was shown that the images were not significantly changed if thedensities in the ray tracing were scaled by the ratio of the attenuation coefficientsof the material. On the other hand, in addition to the advantages of moving fromequipment at the end of its lifespan, this approach has several other advantages. Itis streamlined and efficient for the user, and flexible enough to allow the additionof sophisticated secondary checks and tools after its adoption. Being automated57can also make it a more reliable method for determining the tissue deficit; the oldmethod relied on a set of optical density measurements performed by the physicistby hand, and selected in a biased manner (that is, avoiding the shadow of the ribs)while the new method uniformly samples everywhere under the lung blocks. Thisresults in tens of thousands of points contributing to the summation. Since theadoption of the RunTBI program as the primary tool for treatment planning, it hasbeen used for over 100 patients.2.6 ConclusionIn this chapter, a wide range of advances and uses of Monte Carlo for dosimetricpurposes were presented. New tools were developed that enabled faithful simu-lations of complex treatment plans like VMAT. These simulations were able tobe carried out because the author created new software tools to parse the treatmentfiles for interpretation by the Monte Carlo code and an improved tool for construct-ing digital phantoms from patient CTs. Likewise, alternative models of the linearaccelerator were investigated, and the efficacy relative to the current models wereexplored. It was demonstrated that the new model did not match the measureddata as well as the one currently in use. These tools have improved confidence inthe reliability of the treatment delivery, and have been integrated into the routinetreatment planning process for all external radiotherapy patients at the VancouverCancer Centre.In addition, these tools also give researchers the ability to investigate newparadigms in treatment techniques using individualized plans and phantoms. Oneexample discussed in this chapter was a study of the performance of the commer-cial treatment planning algorithm for spine SABR treatments. It was shown that thepencil-beam based convolution superposition algorithm AAA would overestimatethe PTV coverage in SABR plans relative to the Monte Carlo system, and some-what overestimate the dose to the spinal canal. This reinforces the idea that all ofthe treatment plans generated by that algorithm should be checked by Monte Carlosimulations before delivery to ensure adequate coverage.A new method for the treatment planning of total body irradiation was pre-sented. This technique makes use of patient CTs rather than the conventional ra-58diographic films to design the shape and thickness of blocks to compensate forthe tissue deficit of the lungs. Using a combination a measurements, analytic cal-culations, and Monte Carlo simulations, this new technique was shown to moreaccurately calculate the thickness of the blocks, and gives the physician a betterimage on which to do the necessary contouring. Implemented as the software pack-age RunTBI, this approach has been used for nearly 200 patients at the VancouverCancer Centre.Condensed history Monte Carlo simulations have a long history in treatmentplanning, and the different projects discussed in this chapter continue this tradition.This work has resulted in improvements to the clinical treatment planning processat the Vancouver Cancer Centre. Refinements of existing methods, plus new ones,have improved both the workflow process and enhanced the the reliability of thequality assurance. Ultimately, these advances ensure that patients are receiving thebest possible therapy.59Chapter 3Track Structure Monte Carlo forMicrodosimetryWhile dose can be a useful quantity for many purposes, it does not tell the wholestory about the biological response to radiation because using the dose alone doesnot directly give information about the clustering of energy transfer events in thepatient. Microdosimetry is the study and measurement of energy depositions onthe microscopic scale ? at length scales related to the volumes and structures thatare relevant in determining damage to the DNA and other critical bio-molecules.The nature of the track structure of the radiation can be described by the micro-dosimetric spectra as well as the dose-mean lineal energy, y?d[37].To get this degree of spatial resolution from Monte Carlo simulations, it is nec-essary to employ track structure algorithms. While track structure codes have beenused for radiation transport for decades, the last few years have seen a substan-tial growth in the number of such codes. Although there has been some notableadditions in the experimental data available[16], the progress has not matched thepace on the theoretical side[38][12]. This makes it crucial to compare differentalgorithms to ensure their accuracy and to explore the consequence of their modelsand assumptions. This has also produced criticism of the limits of the applica-tion of these simulations and their quantitative interpretation (which is discussed inSection 3.2), and the best way to resolve these issues would be a rigorous compar-ison between experiments and theory. Furthermore, the usefulness of a multi-scale60approach (Chapter 4) depends on the reliability of each stage of the simulation.Thus, it is imperative that comparative studies of track structure simulations arecarried out.There are a number of such comparative studies. Two of the most widely usedtrack structure simulations codes, NOREC[87] and Geant4-DNA[36], have beencompared to other codes[17][12], as well as each other, those comparisons havefocused on the end results of the DNA damage (e.g. [51]) or the cross sectionsand mean free paths between particle scattering interactions (as in [21]). Theseare clearly important issues, but there are many reasons to also consider the micro-dosimetry of the electron transport of both codes. First and foremost, the chemicalyields (that is, the number of particles of each species produced immediately afterthe ionizing radiation has passed through the medium) ultimately depend on boththe clustering and partitioning of the events, and so the nanodosimetry is dependenton the microdosimetry. By removing the extra layers introduced by the chemistry(with a large number of adjustable parameters and required data fitting), clearerinsight into the remaining layers can be achieved, and allows direct comparisonwith experiment. The rest of this chapter is devoted to this comparison, much of itcovered in Ref. [57].3.1 MicrodosimetryMicrodosimetry refers to a framework for describing the relative biological ef-fectiveness of different radiation types based on differences in the microscopicdistribution of energy deposits[37]. The basic quantity used to describe the mi-crodosimetry of a radiation track is the total energy imparted to the material in avolume by that track, ? . From ? , it is possible to define the lineal energy for thattrack by y= ?/l?, where l? is the mean chord length of the volume. If multiple tracksare considered, it is possible to calculate the probability density that a randomlychosen track will have a lineal energy y, which is also known as the lineal energydistribution, f (y). From this distribution it is possible to calculate the frequencymean lineal energy:y? f =?y f (y)dy (3.1)61Another quantity of interest is the dose-weighted lineal energy distribution, d(y),which can be calculated via the relationshipd(y) = y f (y)y? f(3.2)from which the dose-mean lineal energy can be derivedy?d =?yd(y)dy = 1y? f?y2 f (y)dy (3.3)The dose mean lineal energy has shown to be correlated with the experimentallydetermined relative biological effectiveness (RBE) for low energy electrons byMichael et al.[61][65], so the discussion on microdosimetry will be focused onthe the dose-mean lineal energy and the microdosimetric spectrum yd(y).3.2 Limitations of the Quantitative Use ofMicrodosimetryThere has been some discussion about the limits of the use of track structure codesfor quantitative purposes. In particular, this approach makes use of length scalesthat are very small, and can give numerical results that are more precise than isphysically relevant. As a lower limit, the energy deposition of an excitation eventcannot be ascribed to a volume smaller than that of a water molecule, since the ex-citation is a collective property of all of the electrons in that molecule. In addition,there is a limit to the localization of energy in quantum mechanics. In particular,for particles moving in a volume with identical, localized scattering sites, it is notrelevant to distinguish between scattering sites which are less than a certain dis-tance apart, often referred to as the adiabatic criterion[20]. While this does notdirectly apply to electrons in water (since the dipole electric field of the water islong-range), it does suggest there is a quantum mechanical basis for concern aboutassigning energy events too precisely. The NOREC algorithm actually deals withthis issue by modelling the localization of energy after an event by a stochasticprocess [81]. A related criticism is that track structure algorithms simultaneouslyspecify both the position and momentum of a particle[95]. The consequences of62this would not be noticed in the interaction cross-sections but rather in the dis-tance and direction between the molecules that receive the energy transfers. Theseare not fundamental criticisms of track structure simulations; in this context, trackstructure codes can be considered as approximation schemes with possibly inaccu-rate sampling methods. Importantly, the degree to which a track structure code isinaccurate can be explored experimentally. This would not affect scattering cross-sections, but the discrepancy could be tested by looking at the clustering of eventsand the microdosimetric spectra.3.3 Physics ModelsAny interpretation or discussion of results comparing Geant4-DNA and NORECrests on an understanding of the specific physics models utilized and their imple-mentations. Any proposed explanation for a discrepancy in the simulation resultsultimately rests in the differences in the models or the codes that manage the sim-ulations. This section will focus on the underlying models, and Section 3.5 willdiscuss the implementation of simulation management codes.3.3.1 Inelastic InteractionsBoth simulations use the interaction cross-section to determine the mean free pathof the particles being transported, and then partition the events into different chan-nels. The ?inelastic? cross-section will refer to both electronic excitations andionizations (but will exclude mechanical transformations). NOREC models theseinteractions in the same way as the original OREC code[87], while GEANT4-DNAuses different models for electrons with energies above or below 10 keV[25]. Bothcodes[81][25] use a dielectric function, ?(q,E), to calculate the doubly differentialcross section (DDCS) using the first Born approximation:d2?dqdE? Im[?1?(q,E)](3.4)where ? is the total inelastic cross section, E is the energy transfer and q is themomentum transfer and Im[?1/?(q,E)] is referred to as the dieletric dispersion re-lationship. However, a first principles calculation of the dielectric loss function for63Transition Geant4-DNA Threshold Energy (eV) NOREC Threshold Energy (eV)A1B1 8.10 8.4B1A1 10.10 10.1Rydberg A+B 12.00 11.26Rydberg C+D 13.51 11.93Diffuse Band 14.41 14.1Table 3.1: The excitation energies in Geant4-DNA[21] and NOREC[81] forscattering by electrons in water, empirically fit to optical data.a molecule with a complex electron configuration like water has not been achievedat this time, and it is also difficult to measure experimentally[26]. However, op-tical data (q = 0, and referred to as the energy loss function, ELF) is available forliquid water, and so the dielectric function is approximated from this data usingan extended Drude model to empirically fit the data[21]. Based on the data usedand the approximation scheme employed, the energy levels of the excitations aresomewhat different in the two models(see Table 3.1) but fairly consistent with eachother. The NOREC model has smaller values for the higher energy states as wellas a reduced threshold for the A1B1 state.Using these parameters, it is possible to calculate the energy loss function foreach model. For q = 0, the ELFS! are compared in Figure 2 of Ref. [21]. TheNOREC model overestimates the peak ELF by 33% compared to the model em-ployed by Geant4-DNA (and the experimental data) for electrons with an energyless than 10 eV, and somewhat reduces the ELF for energy losses greater than 30eV. The models are in good agreement below the peak near 20 eV. This meansthat a larger fraction of the inelastic events would be at lower energies transfersfor NOREC. Therefore, the average angular scattering deflection for low energyelectrons would be smaller in the NOREC simulations than in the Geant4-DNAones. On the other hand, since this increased elastic scattering probability wouldincrease the total cross-section, the stopping power will be larger, and the totaltrack path-lengths for the electrons will be smaller in NOREC.For electron energies greater than 10 keV, Geant4-DNA uses a different modelthan NOREC for the inelastic cross sections with relativistic corrections taken intoaccount based on Ref. [9]. The same energy loss functions were employed, but the64relationship between the dielectric dispersion relationship and the differential crosssections has relativistic corrections and is divided into longitudinal and transversecomponents. The study showed that the relativistic cross sections are about 30%larger than the non-relativistic calculations, which increases for higher energies.3.3.2 Elastic InteractionsThe umbrella of elastic interactions as implemented in both codes involves threeseparate contributions. The first is traditional Coulomb scattering from a staticcharge distribution. This also includes corrections for the polarizability of the wa-ter molecules and the correlation of the charge distribution of nearby moleculesthat end up screening the dipole. Lastly, the quantum mechanics of identical parti-cles demands the inclusion of an exchange correction. For Geant4-DNA, there aretwo possible data-sets for elastic scattering. In this study, the model of Championet al.[15] has been chosen, and so the rest of the discussion in this section will referto that model. For the eleastic cross-sections, both Geant4-DNA and NOREC arebased on data measured in water vapor but scaled to the density of liquid water.Since the elastic scattering has a limited effect on the total cross-sections for elec-trons with energies greater than a few eVs, the potential inaccuracy from the phaseeffects is minor.Both codes base the calculation of the elastic scattering cross section on ascreened, static Coulomb potential. In the Geant4-DNA model considered here,the charge distribution is based on Slater type orbitals centered on the oxygenmolecule. To account for the random orientation of the water molecules in theliquid phase, it is convenient to use the spherical average approximation[3]. TheNOREC simulation uses the NIST database ELAST[85]. This approach makesuse of Dirac-Hartree-Fock self-consistent charge density distributions to calculatethe phase shift of the scattering. Both approaches also include a correction basedon electron exchange. These are based on treating the molecular electrons as afree electron gas, although the nature of the two potentials are somewhat different.The Geant4-DNA code also corrects for the fact that the charge distribution is notstatic and changes due to the response of the atomic electrons to the motion of theparticles. This is included in the calculations of the cross sections by a combined65correlation-polarization term. At large distances, only the polarization potential isconsidered. For distances less than some fixed parameter, the potential is chosento be the maximum between the correlation and polarization potentials at that dis-tance. Despite these differences, the tabulated inverse mean free paths for elasticscattering agree for the two codes. A comparison of the integrated (that is, elastic)cross sections shown in Fig. 3 of Ref. [3] and the elastic inverse mean free path forthe NOREC code in Fig. 1 of Ref. [87] shows the consistency of the approaches.3.3.3 Subexcitation ElectronsThe two codes have a different approach to particles that do not have enough energyto induce electronic excitations. For Geant4-DNA, the electrons are followed untilall kinetic energy is lost. The interaction channels involved are vibrational excita-tions, elastic scattering and electron attachment[25]. NOREC follows a somewhatdifferent approach[87]. An electron with a kinetic energy below 10 eVwill undergoa very large number of elastic collisions before its energy falls below the minimumexcitation energy (7.4 eV in this code), and this makes the particle?s path appearto be like diffusion. Therefore, once an electron?s energy falls between 7.4 eV and10 eV, the code will only perform one more transport step. A random distance anddirection (chosen isotropically) for the final event is selected, and a new energy isassigned based on the energy loss from an inelastic collision. The particle is thendiscarded at this point. This difference between the codes will have some effect onthe clustering of the electron energy. However, the small energy transfers involvedhere have a limited effect on the radiolytic yields and the distances traversed bythese electrons are not very large compared to the lengths being studied in thischapter. This means that the differences here will have a small role in the resultspresented. On the other hand, these can have a huge effect on the computational ef-ficiency of the two codes. If the full Geant4-DNA simulation is performed, a largefraction of the computer time is spent following the subexcitation electrons. A sub-stantial reduction in simulation times was observed if the user applies the NORECapproximation to the Geant4-DNA simulation without a significant change in themicrodosimetric mean values calculated.663.4 Global Track MeasuresGlobal track measures are a convenient checks to compare different radiation trans-port simulation codes. Values like the range, penetration range, and total tracklength are straightforward to calculate and easy to interpret. Although the globalparameters are not strongly related to the local quantities that are of interest to mi-crodosimetry, they are relevant to the design of track structure simulations becausethey can be used to optimize the geometry to maximize efficiency without intro-ducing artifacts from edge effects and other sources of electronic disequilibrium.In order to investigate these parameters, a simulation was performed using acustom designed application using NOREC. Mono-energetic electrons were simu-lated with initial positions at the origin with an initial velocity in the +z-direction.For each track, the range was calculated based on the maximum distance from theinitial position to any interaction. In addition, the penetration range was taken to bethe maximum z-value of any interaction in the track. The track length was the sumof the distances between consecutive interactions involving the primary electron.The simulations were run for energies between 10 eV and 1 MeV, with 2? 104independent tracks generated for each energy.These results are shown in Figure 3.1 along with results from Geant4-DNAsimulations (provided by S. Incerti from Ref. [26]) and the theoretical values ofMeesungnoen, et al.[60]. For NOREC, any electron with an energy less than 10eV was immediately discarded. The values for each of these approaches agreesvery well above 1 keV. Below this range, the theoretical result for the total tracklength begins to diverge from the simulated values, and becomes nearly an order ofmagnitude greater at 20 eV, which can be attributed to the fact that the simulationsdo not follow the electrons to the end of the track, which corresponds to an extradistance of around 10-20 nm. In addition, NOREC begins to show an increase inthe penetration depth and range compared to Geant4-DNA at these low energies,but somewhat underestimates the total path length. Since the penetration depth andrange are based on the direction of motion, while the track length is not, this showsthat the particles simulated by NOREC do not undergo as much angular deflectiondue to the scattering. On the other hand, the energy transfer for low energy elec-trons per unit track length is larger for NOREC. This is consistent with the larger67Figure 3.1: A comparison of the global track measures. The range (trian-gles), penetration range (circles), and track length (squares) for elec-trons in water for different initial energies. These calculations are basedon Monte Carlo simulations performed with NOREC (filled shapes),Geant4-DNA[26] (empty shapes), as well as the theory of track lengthcalculated by Meesungnoen[60] (dashed curve). Figure courtesy of Se-bastien Incerti.inverse mean free path discussed in Section 3.3.2. On the whole, however, there isreasonable consistency between all three models.3.5 Software Design for Microdosimetry SimulationThe goal of this study was to compare the microdosimetry and dose-mean linealenergies for the two algorithms. Both provide the user access to the information68related to each interaction, but the user must create the software to oversee andmanage the simulations. Since NOREC was distributed as a pre-compiled library,a C++ application was created in order to manage the runs. The Geant4-DNA is apackage within the open-source Geant4 library, and so the simulation were man-aged by a custom Geant4 application made by implementing concrete extensionsof abstract classes (that is, interfaces) in C++. There were three main tasks thatthe simulation management code was responsible for: determining the initial trackvariables (position, energy, direction, and particle type), performing range rejectionon the particle tracks, and scoring the quantities of interest. Further details of thespecific implementations of these tasks in each code are discussed in the followingsubsections.3.5.1 Initial Track VariablesThere were two types of sources that were required for this study. In both cases,the initial energies were monoenergetic and only electrons were considered. InSection 3.7, particles were uniformly distributed within a sphere with an isotrop-ically uniform initial velocity direction, and was used for simulations where theelectron energy was less than or equal to 5 keV. The other source, discussed inSection 3.6, was used for all initial energies greater than or equal to 5 keV. The z-component of the initial position was the set equal to one-half the penetration depth(from Section 3.4) for all particles. The x and y values were independently selectedrandomly with uniform density on the interval [-0.5,0.5] and the multiplied by thesource side length (which was 10 cm in this study). The initial trajectory was inthe +z? direction. This type of arrangement is referred to as a planar parallel beam.Both set-ups are demonstrated in Figure 3.2 and Figure ScoringThe scoring was performed in spherical sensitive volumes. These volumes arevirtual in the sense that they do not affect the radiation transport of the electrons inany way, and are only used for calculating the dose distributions. The lineal energywas calculated in the scoring spheres by summing up the total energy depositionwithin the sensitive sphere for each particle history, and then dividing by the mean69Figure 3.2: A schematic representation of the geometry used for the planeparallel beam simulations. The electrons were given an initial posi-tion randomly distributed in the source square and a uniform directiontowards the plane containing the scoring volumes, which are packedin a rectangular grid with spacing equal to the twice the radius of thespheres (that is, a grid spacing of 1 ?m). The array of spheres has a sidelength that is smaller than the source plane to assure particle equilibriumthroughout the scoring region.70Figure 3.3: A schematic representation of the geometry used for the isotropicparticle source simulations. The electrons were given an initial positionrandomly distributed within the source sphere and a isotropically dis-tributed initial direction. The sensitive sphere is placed at the center ofthe source sphere, with a radius of 0.5 ?m.71chord length of the sphere (l? = 4r/3 where r is the radius of the sphere). For thesimulations with the parallel beam source, an array of (non-overlapping) sensitivespheres was placed in the z=0 plane. Each of the spheres had a diameter of 1 ?m,and the spheres were arranged with their centers in rows and columns parallel tothe x and y axes. In order to avoid edge effects, the array side-length wassa = ss? (R(E)?1.1+d) (3.5)where ss is the side-length of the source, R(E) is the range of electrons with anenergy E, and d is the diameter of the sensitive spheres. The 1.1 is used as anextra margin to account for electrons whose range is larger than the mean, and theamount is standard in the literature[66]. In order to maintain statistical indepen-dence, scoring was only performed in the sensitive sphere that first recorded anenergy deposition for that primary particle. In the other simulation geometry (withthe uniform, isotropic source used for the lower energy electrons), there was onlyone sensitive sphere, placed at the origin, with a diameter of 30 nm. The geometriesin both situations are shown in Figure 3.2 and Figure Range RejectionIn Monte Carlo track structure simulations, range rejection is a technique to im-prove computational efficiency by discarding particles once their range is smallerthan the distance to sensitive volumes for scoring. For the simulations in this sec-tion, the distance calculated for range rejection purposes depended on the statusof the scoring array. If energy had already been deposited in one sphere in thescoring array, then that ?active? sphere?s center was used to calculate the distance.However, if no energy had been deposited in one of the spheres yet, the centerof the nearest sensitive sphere was selected. If this distance was larger than thecurrent range of that particle plus the radius of the sensitive spheres and a ten per-cent margin, the particle was discarded. For the Geant4-DNA simulations, rangerejection was performed after any particle step (such as an interaction), during theProcessHits subroutine call of the scoring object.For the NOREC simulations, an alternative range rejection approach was em-ployed. Since NOREC does not directly give the user access to the energies of the72electrons being followed, these were tabulated outside of the NOREC track part ofthe simulation. In particular, it is not possible to distinguish between the secondaryelectrons with the public functions of the library, and so the user can only be sure ofthe energy of the primary electron. This is accomplished by subtracting the energytransfer of each particle step from the initial kinetic energy of the primary particle.This means that after each interaction involving the primary particle, it was possi-ble for the user to know both the primary particle?s position and energy, becauseNOREC will only transport the primary particle once there are no more secondaryparticles to transport. Thus, the range rejection checking was done only after aninteraction involving the primary electron.3.5.4 Ensuring Convergence of the SimulationsFor each simulation, 3.5? 104 independent tracks were scored in a sensitive vol-ume. Simulations were then run in batches of 1000 particles, and after each batchthe dose-mean lineal energy (y?d) for all of the tracks was calculated. These batcheswere run until the dose-mean lineal energy varied by less than 1% for five consec-utive batches in order to ensure convergence had been achieved.3.6 Parallel BeamSimulations were performed using the parallel beam source (as shown in Figure 3.2)and an array of sensitive volumes for many energies between 5 keV and 500 keV,using sensitive spheres with a diameter of 1 ?m for both Geant4-DNA and NORECat two different depths for each energy. In the first case, the centers of the sensitivespheres were 0.5 ?m from the source plane, and in the other case the depth was setto be equal to one-half of the penetration range from Figure 3.1 (except for the 5and 10 keV sources, for which the distance was 1 ?m). The microdosimetric spec-tra for these simulations are shown in Figure 3.4 to Figure 3.6. The dose-meanlineal energies at both of these depths are shown in Figure 3.7, and the ratio be-tween the value at these depths is presented in Figure 3.8. All of the simulations inthis chapter were carried out on a Linux Cluster with 24 nodes (Intel i5 quad-coreprocessors with 2.4 GHz clock speed and 6 GB of memory).There is very strong agreement between the two simulations. The dose-mean7310?2 10?1 100 101 10200. 500 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(a)10?2 10?1 100 101 10200. 500 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(b)10?2 10?1 100 101 10200. 250 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(c)10?2 10?1 100 101 10200. 250 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(d)10?2 10?1 100 101 10200. keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(e)10?2 10?1 100 101 10200. 150 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(f)Figure 3.4: Microdosimetric spectra for electrons at a depth of 0.5 ?m (leftside) and at one-half of the penetration range (right side) for 500 keV, (a)and (b), 250 keV, (c) and (d), and 150 keV,(e) and (f), simulated usingGeant4-DNA and NOREC using 1 ?m-diameter sensitive spheres.7410?2 10?1 100 101 10200. keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(a)10?2 10?1 100 101 10200. keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(b)10?2 10?1 100 101 10200. 50 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(c)10?2 10?1 100 101 10200. keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(d)10?2 10?1 100 101 10200. 25 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(e)10?2 10?1 100 101 10200. 25 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(f)Figure 3.5: Microdosimetric spectra for for 80 keV, (a) and (b), 50 keV, (c)and (d), and 25 keV, (e) and (f), electrons simulated using Geant4-DNAand NOREC using 1 ?m-diameter sensitive spheres, at a depth of 0.5?m (left side) and at one-half of the penetration range (right side)7510?2 10?1 100 101 10200. 10 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(a)10?2 10?1 100 101 10200. 10 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(b)10?2 10?1 100 101 102012345 5 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(c)10?2 10?1 100 101 10200. 5 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(d)Figure 3.6: Microdosimetric spectra for electrons at a depth of 0.5 ?m (leftside) and 1 ?m (right side) for 10 keV, (a) and (b), 5 keV, (c) and (d),simulated using Geant4-DNA and NOREC using 1 ?m-diameter sensi-tive spheres.lineal energies agree to within 5% for all energies expect for 10 keV at a depthof 0.5 ?m. For the larger distances, all of the dose-mean lineal energies agree towithin 13%, although the difference is consistently larger than 10% for energiesless than or equal to 50 keV, and much smaller for higher energies. The most no-ticeable discrepancy occurs at 10 keV, where the NOREC results show a significantenhancement of the peak for the tracks that are entirely contained within the sen-sitive sphere. This can partly be explained by the fact observed in Section 3.4 thatNOREC particles have a reduced track length. Another difference in the spectracan be observed by the fact that the global peak in the spectra is much higher for760 100 200 300 400 5000123456Energy (keV)Dose?Mean Lineal Energy (keV/?m)  Geant4?DNANOREC(a)0 100 200 300 400 5000123456Electron Energy (keV)Dose?Mean Lineal Energy (keV/?m)  Geant4?DNANOREC(b)Figure 3.7: Dose-mean lineal energies simulated for electrons between 5 and500 keV in water using Geant4-DNA and NOREC for 1 ?m-diameterspheres at a depth of (a) 0.5 ?m and (b) one half of the penetrationrange.77the NOREC simulations, while the location and widths of those peaks generallyagree well. This difference is driven by the fact that NOREC recorded a slightlysmaller number of larger events. This is consistent with the observations above thatthe angular scattering in NOREC is smaller than for Geant4-DNA, since straightertrajectories would mean fewer scattering events within the sensitive sphere. In ad-dition, at the higher energies, there are several well-defined peaks in the sub-100eV event range for Geant4-DNA that are not present in the NOREC spectra. Theseare related to single scattering events which are higher energy excitations whichare included in Geant4-DNA but not in NOREC.In general, Geant4-DNA also showed a smaller dose-mean lineal energy at thedeeper set of sensitive spheres. Since the primary particles should still be consid-ered ?fast? at one-half of their penetration range in these simulations, there shouldnot be very significant spread in the energies of the primary electrons when theyreach the spheres. However, since the stopping power is smaller for the Geant4-DNA model, we should expect that these particles in the Geant4-DNA simulationhave higher energy at this distance. For the energy ranges studied here, this wouldcause a reduced dose-mean lineal energy as observed in this study. Another wayto consider this effect is to look at the relationship between the dose-mean linealenergy at the two different depths. The ratio of the two different depths is shownin Figure 3.8. For energies less than 200 keV, the ratios agree very well for the twodifferent simulations, but above this level the NOREC results show an increaseddose-mean lineal energy. Again, this is consistent with lower energy primary elec-trons reaching the sensitive spheres in the NOREC simulations.3.7 Isotropic SourcesSimulations were also performed for a uniform isotropic source for electrons (seeFigure 3.3) with energies from 100 eV to 5 keV for sensitive spheres with a di-ameter of 30 ?m. The electrons were randomly distributed within a sphere with adiameter equal to the twice the range of the electrons plus the diameter of the sen-sitive spheres plus a ten percent margin. The microdosimetric spectra are presentedin Figure 3.9.At lower energies, the NOREC spectra show substantially more pronounced780 100 200 300 400 5000.511.52Energy (keV)Ratio of Dose?Mean Lineal Energies  Geant4?DNANORECFigure 3.8: The ratio of the dose-mean lineal energy at a depth equal to onehalf the penetration range to the dose-mean lineal energy at a depth of0.5 ?m.central peaks than in the Geant4-DNA results. This is probably due to the differ-ences in the two approaches to handling the transport for electrons with energiesless than 10 eV. For low energies, this accounts for a sizeable fraction of the totalenergy deposit, and reduces the track length of the electrons in the NOREC sim-ulation by around 10 nm (which is nearly the radius of the spheres). In addition,the Geant4-DNA simulation includes vibrational scattering events, which are mostprevalent for low energy electrons. While these events do not transfer much energy,they can result in large scattering angles, which would tend to lead to more eventswith smaller deposits. While the difference in the peak height and width is veryprominent, it does not have a major effect on the dose-mean lineal energy below 1keV because most events are still close to the total kinetic energy of the electrons.Above 1 keV, the dose-mean lineal energies for Geant4-DNA and NOREC dif-7910?2 100 102 10400. keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(a)10?2 100 102 10400. 2 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(b)10?2 100 102 10400. 1 keVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(c)10?2 100 102 10400. eVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(d)10?2 100 102 1040123456 500 eVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(e)10?2 100 102 104024681012 300 eVLineal Energy (keV/?m)yd(y)  Geant4?DNANOREC(f)Figure 3.9: Microdosimetric spectra in a sensitive sphere with a diameter of30 nm simulated using Geant4-DNA and NOREC. The particles wereuniformly distributed in the volume with random direction. Initial ener-gies were (a) 5 keV, (b) 2 keV, (c) 1 keV, (d) 700 eV, (e) 500 eV, and (f)300 eV.800 1 2 3 4 5051015202530Energy (keV)Dose?Mean Lineal Energy (keV/?m)  Geant4?DNANORECLiamsuwan, et al.Figure 3.10: The dose-mean lineal energy in 30 nm targets for uniformlydistributed monoenergetic electrons simulated with Geant4-DNA,NOREC, and KURBUC for liquid water by Liamsuwan, et al.[53].fer by between 1% and 5%. This could be caused by the smaller inelastic inversemean free path in NOREC caused by the larger electron loss function employed,and is consistent with the results for the larger sensitive spheres for the parallelbeam discussed above. The microdosimetric spectra in this range are mostly simi-lar, although there are differences. For 1 keV and 2 keV, the NOREC peak is stillmore pronounced and shifted slightly towards larger events, as one would expect.Both the spectra and the dose-mean lineal energies are consistent for the Geant4-DNA and NOREC simulations.The dose-mean lineal energies were also compared to the results published byLiamsuwan et al.[53] for liquid water. The data shown in Figure 3.10 for Ref.[53], were measured in cylinders with a height and diameter of 30 nm, which havethe same mean chord length as the 30 nm diameter spheres, and the distance and81direction from the starting point of the track to the center of the sensitive volumeshave the same distribution. At the lower energies the agreement is very good, butat these energies the electron?s ranges are about the same size as the dimensionsof the sensitive volumes. This means the agreement there is more a reflection ofthe initial conditions and ranges rather than the properties of the transport. Atlarger energies the agreement is not as good. In Liamsuwan, et al., it was shownthat the stopping power for NOREC (and likewise Geant4-DNA, which can beinferred from Figure 3.1) is much greater than for KURBUC at energies greaterthan a few tens of eVs to 1 keV, and slightly smaller for higher energies (see Fig.1 in [53]). Since the electrons in these studies are randomly distributed aroundthe sensitive volumes, the primary electron energy spectrum in the sensitive sphereshould be broadly distributed. This could mean that the primary electrons losemore energy outside of the sensitive volume in the KURBUC simulation, and soenter the scoring region with less energy. That would lead to an shift towards largerevents in the spectrum, and consequently a larger mean lineal energy.3.8 ConclusionA thorough comparison between the microdosimetry of electrons simulated byGeant4-DNA and NOREC was carried out over most of their ranges of validity.For the nearly monoenergetic electrons crossing through the spheres at a depth of0.5 ?m, the agreement was very strong. Somewhat larger differences were foundat deeper depths. This reflects differences in the inverse mean free path lengths.At these energies, this is mostly because of the differences in the inelastic inversemean free path resulting from the extrapolation of the dielectric functions fromthe energy loss function. Despite this fact, the agreement is still fairly good. Forenergies 5 keV and below, there is also a high level of agreement. For these ener-gies, the elastic scattering channels play a bigger role, and in this case, the largerelastic scattering inverse mean free path helps create a small difference betweenthe models. In the end, it is possible to relate these discrepancies directly back todifferences in the underlying models. This insight can contribute to understandingand interpreting the results of track structure simulations. More importantly, thesegive us results that can be experimentally measured using standard techniques with82tissu-equivalent proportional counters, and allows us to empirically validate thetrack structure codes, and to explore the veracity of the underlying models andassumptions.83Chapter 4Multi-Scale Monte CarloThe use of microdosimetric information with megavoltage external beams wouldpermit the use of more radiobiologically based analysis of treatment plans thanusing dose alone. Track structure codes can give this information, but cannot berun in a reasonable amount of time for large numbers of particles. To speed up thisprocess, one can use a condensed history code to simulate the radiation transportto generate a charged particle spectra for use by the track structure codes in asmall subvolume of the simulation. This means that most of the particles willbe simulated quickly using the condensed history code, and only a portion of thetracks of the rest of the particles will be simulated with the track structure code.Of course, paying attention to the assumptions made at each scale is crucial. Forradiation transport, the most important inconsistency is that the condensed historycodes discard the low energy electrons, which make an important contribution tothe microdosimetric spectra. Amethod to ensure that the subthreshold electrons areproperly taken into account is discussed in Section 4.1, and results justifying thisapproach are presented in Section 4.2 and Section 4.3. These sections are based onRef. [58]. In addition, these tools are then applied to a clinical radiation therapyplan for a patient CT-based phantom in Section 4.4.844.1 Watch RegionAccurate track structure simulations require an accurate charged particle spectrum.The fluence through the sensitive volume where the scoring occurs needs to reflectthe contribution of low energy electrons, but the condensed history algorithm doesnot accurately track them. The solution is to perform the track structure simulationin a volume that is larger than the sensitive volume. Since the low energy electronshave a short range, if this volume (which will be referred to as the watch volume)is large enough, then the interior of the watch volume will have the correct flu-ence. Thus, the scoring should only be done within this interior portion. Ideally,the watch volume should extend beyond the scoring volume just far enough so thatthe fluence is correct ? but any additional distance causes the simulation to be-come inefficient. This extension is expected to be equal to the range of the largestsubthreshold electron in the condensed history code.4.2 Uniform ElectronsThe role of the watch volume in a multi-scale code needed to be investigated inorder to explore several questions. The first is to prove that this approach increasesthe accuracy of the simulation. Equally importantly, it is necessary to establish howto determine the optimal size of the watch volume and to determine if this approachcan feasibly be applied to actual treatment plans and geometries. In order to answerthese questions, it was best to begin with a monoenergetic electron source as thisallow the most direct interpretations of the results.4.2.1 Simulation MethodsThe condensed history part of the simulations was simulated on a Linux clusterwith an Intel quad-core 2.8 GHz with 8GB of memory, using a custom designedEGSnrc++ user code[44]. A cube of water (side length 50 cm) was centered on theorigin. A watch volume was placed at the origin. If a charged particle entered (orwas created) inside the watch volume, its trajectory information (position, direc-tion, energy, and charge) was written to a file and it was immediately discarded bythe EGSnrc++ user code. The electrons? initial position was uniformly distributedin a cube centered on the origin, with a sidelength equal to twice the range of the85electrons plus the diameter of the sensitive sphere that was used for scoring thelineal energy. The initial kinetic energy was set to 50 keV, and all particles hadthe same initial direction. Range rejection was performed based on the distance ofthe particle from the watch volume, and the global electron total energy cut-off (orthreshold) parameter, ECUT, was set 0.521 MeV, and the material data for waterwas taken from the pegs4 data file 521icru.The phase space description file was then read and used to generate the tracks inan application based on NOREC. A 1 ?m diameter sphere was placed at the originas the sensitive volume for scoring, using the methods discussed in Section 3.5.2.In this case, the total energy deposition was accumulated for all electrons arisingfrom the same original electron in the condensed history part of the code. Rangerejection was done by NOREC in the same manner as in Section 3.5.3, exceptthat the ?primary? particle refers to the electron that generated the current trackin NOREC, not the original 50 keV electron from the condensed history part ofthe code. The track structure part of the simulation was carried out on a singleWindows laptop with a quad-core Intel i5 processor (2.4 GHz) with 6.00 GB ofmemory.Simulations were run for a watch volume radius that varied from 0.5 ?m to30.5 ?m. At a radius of 0.5 ?m, the watch volume did not extend beyond thesensitive sphere, and at 30.5 ?m was larger than the range of the electrons, and theentire simulation was carried out with the track structure code. Convergence wasassured by the means elaborated in Section ResultsA comparison of the spectrum for the full track structure simulation (30.5 ?m ra-dius) and the case where the watch volume is the same as the scoring volume (0.5?m radius) is shown in Figure 4.1. For the smaller watch volume (which does notextend past the sensitive volume), the spectrum shows fewer events with a linealenergy greater than 2.3 keV/?m. For electrons that have a range much larger thanthe dimensions of the scoring volume, the lineal energy should approximately beequal to the particle?s LET. Using this approximation, these missing events cor-respond to electrons with an energy of less than 10 keV, that is, the subthreshold8610?2 10?1 100 101 10200. (keV?m)y d(y)  0.5 ?m Radius30.5 ?m RadiusFigure 4.1: A comparison of the microdosimetric spectra for simulations of50 keV electrons with a watch volume that has a radius of 0.5 ?m (equalto the radius of the sensitive sphere, blue line) and 30.5 ?m (containingall electrons, red line).electrons.Next, consider the situation with a watch volume with a radius of 2.5 ?m. Inthis situation, the watch volume extends beyond the sensitive volume by a distanceequal to 2 ?m which is about equal to the range of the largest subthreshold elec-trons. Figure 4.2 compares the 30.5 ?m spectrum to the spectrum generated fromthe electrons with a 2.5 ?m watch volume. These spectra have a strong agreementand demonstrate that the watch volume does not have to be very large in order togive accurate results.This is further demonstrated by the looking at the value of the frequency-meanand dose-mean lineal energies calculated for different watch volume sizes. Theseare shown in Figure 4.3. From this figure, it can be seen that the dose-mean linealenergy is reduced from 2.7751 keV/?m to 1.850 keV/?m when the watch volumeis reduced from 30.5 ?m to 2.5 ?m. This shows that the calculation for the 0.5?m watch volume radius underestimates the value by 33%. On the other hand, themean lineal energy is correctly predicted to within 1% whenever the watch volume8710?2 10?1 100 101 10200. (keV/?m)  2.5 ?m Radius30.5 ?m RadiusFigure 4.2: A comparison of the microdosimetric spectra for simulations of50 keV electrons with a watch volume that has a radius of 2.5 ?m (equalto the radius of the sensitive sphere plus the range of the subthresholdelectrons, dotted blue line) and 30.5 ?m (containing all electrons, solidred line).extends beyond the sensitive volume by at least 2 ?m.This accuracy is only useful if a sufficient number of events can be accumulatedto ensure statistical convergence is acceptable within a reasonable time. The totalsimulation time (including both the condensed history and track structure parts ofthe simulation) was calculated for all of the monoenergetic electron simulationsand this time was divided by the number of events scored for each track. The av-erage time per scored track is shown in Figure 4.4. The time per scored event wasreduced by a factor of 109 by reducing the watch volume radius from 30.5 ?m to0.5 ?m. Since the size of the watch volume is determined by the cut-off thresholdin the condensed history part of the simulation, this method can easily be appliedto any multi-scale radiation transport simulation. The total time for the calculationwith a 2.5 ?m watch volume was 3 hours, with the overwhelming majority of thesimulation devoted to the track structure simulation, which was carried out on lessthan state-of-the-art and somewhat outdated personal laptop. Naturally, a dedicated880 5 10 15 20 25 30 350.511.522.53Watch Volume Radius (?m)Mean Lineal Energy (keV/?m)  ydyfFigure 4.3: The frequency-mean (y? f , blue squares) and dose-mean (y?d , redcircles) lineal energies for simulations with different watch volume sizesfor 50 keV electrons.high end computer or cluster would result in at least an order of magnitude reduc-tion in computer time. Since the mean free path of higher energy charged particlesis larger than for lower energy particles, the efficiency of the track structure codewould actually increase for a more realistic simulation. In addition, the fact thatthe most computationally demanding part of ths simulation was carried out by asingle laptop shows that incorporating track structure simulations in clinical prac-tice is feasible, particularly if the track structure simulations are performed on adedicated high-end computer or cluster.4.3 Benchmarking with PhotonsThe most fundamental validation of this approach must be carried out with pho-tons, which allows a more complex comparison with experimental results. A par-allel beam of photons was incident on the surface of a water phantom, and thetracks were simulated using the EGSnrc++ code discussed above. The spectrum ofthe particles was taken to correspond to 241Am (photons only, with a mean energyof 60 keV), 137Cs (with an energy of 662 keV), 51Cr (with a mean energy of 3208910?1 100 101 10210?210?1100101102Time Per Track Scored (s)Watch Volume Radius (?m)Figure 4.4: The average time per event scored in the sensitive volume for 50keV electrons simulations with different watch volume sizes.keV), and 99mTc (with an energy of 141 keV) radioisotope sources. These isotopeswere chosen because of their relevance as therapeutic sources and/or calibrationsources. The watch volume in this case was a parallelepiped centered at a depth of2 cm in the phantom. In order to avoid edge effects and disequilibrium, the halfside-lengths of the watch volume were shorter than the beam?s half side-length by5 cm. The thickness of the watch volume in the direction of the beam was setto be 5 ?m, so that it extended 2 ?m beyond the sensitive spheres in the trackstructure simulation. For the track structure component, the spheres were arrangedin an array (as in Section 3.5.2), but no spheres were placed within 2 ?m of theedges of the watch volume. The range rejection, scoring, and convergence check-ing were performed as for the monoenergetic electrons. The ECUT parameter inthe track structure simulation was 0.521 MeV, and the global photon cut-off thresh-old, PCUT, was set to 0.01 MeV, and the material data was taken from the 521icrupegs4 data file.The spectra for these simulations are shown in Figure 4.5. As the average en-ergy of the initial spectrum increased, the width of the distribution as increased.Furthermore, the higher energy photons showed both a shift to lower energy events9010?2 10?1 100 10100. (keV/?m)y d(y)  AmCr TcCsFigure 4.5: The microdosimetric spectra for photons in water calculated fromMonte Carlo simulations of 137Cs (black), 51Cr (magenta), 99mTc (blue),and 241Am (red, and only using the photons part of the spectrum) ra-dioisotope sources incident on the surface of a water phantom. An arrayof sensitive spheres with a diameter of 1 ?m was placed at a depth of 2cm in the phantom.and a lower energy peak value. This is in agreement with the trends reportedfor the experimental results of Kliauga and Dvorak[48], Varma et al.[100], andfrom ICRU Report 36[37]. A comparison with experimental data is presented inFigure 4.6. The trend here to decreasing dose-mean lineal energies for larger ener-gies is captured by the simulations and the magnitude of the values is consistent.4.4 Clinical UseThis method can ultimately be applied to clinical treatment plans. To demonstratethis, a lung SBRT case was selected for simulation, which was performed using atwo arc VMAT treatment plan using a 6 MV photon beam on a Varian TrueBeamlinear accelerator. The microdosimetric spectra and the dose-mean lineal energywas computed in the voxel containing the treatment isocenter.910 100 200 300 400 500 600 70000.511.522.533.544.5Photon Energy (keV)Dose?Mean Lineal Energy (keV/?m)  SimulationVarma, et al.Figure 4.6: The dose-mean lineal energy for photons in water calculated fromMonte Carlo simulations of 137Cs, 51Cr, 99mTc, 241Am (using only thephotons part of the spectrum) radioisotope sources incident on the sur-face of a water phantom. An array of sensitive spheres with a diameterof 1 ?m was placed at a depth of 2 cm in the phantom. The simu-lated results are represented by triangles, and the experimental resultsof Varma, et al.[100] are shown with circles. The results are reported asa function of the average energy of the photon sources.4.4.1 Simulation MethodThis simulation was done using a modified version of the EGSnrc++ and NORECapplications discussed earlier in the chapter. The watch volume was set to be thevoxel containing the isocenter (which had a side length of 2.5 mm). In order tosimulate the dynamic motion of the linear accelerator during beam delivery, it wasnecessary to update the BeamSource library in EGSnrc++. This source calls aBEAMnrc shared library accelerator, but since this was a VMAT plan, the simu-lation required that the monitor unit index is synchronized with the gantry, colli-mator, and collimator angles as well as any change in the SSD and translation ofthe isocenter. In effect, this called for the development of an EGSnrc++ version ofthe DOSXYZnrc source 21. The track structure part of the simulation was carriedout in the same manner as discussed above, except that a 3-D array was employed92for the sensitive spheres (with a diameter of 1 ?m). No sensitive spheres wereplaced within 2 ?m of the edges of the watch volume to ensure a proper built upof the subthreshold electrons. The voxel containing the isocenter was soft-tissueequivalent which ensures that the use of water-equivalent material in NOREC isappropriate.Another crucial difference in this situation was the energy of the particles in-volved. The photons in this situation had energies greater than 1 MeV (up to 6MeV), which lead to the creation of positrons as well as electrons. In this case0.9% of the particles written to the phase space description file were positrons.Following conventional practice, the NOREC portion of the code treated thesepositrons as if they were electrons. This neglects the fact that the cross sectionsare slightly different. The main difference is the potential for annihilation in flightof the positron, which was not modelled in the simulation, but this process is muchless likely than annihilation at thermal energies (which was much smaller than the7.4 eV cutoff). This also neglects the possible reabsorption of the annihilation pho-tons within the watch volume. This is a reasonable approximation so long as thewatch volume dimensions are small relative to the mean free path of the photons,which is true in this case (for a 0.511 MeV photon in water, the linear attenuationcoefficient is 0.1 per cm[6]).The scattering of the high energy photons can also produce electrons with anenergy greater than 1 MeV. However, NOREC is not equipped with cross sectionsfor these electrons. On the other hand, the LET value for electrons in the 1-6MeV range is lower than for lower energy electrons, so these particles will makevery small contribution to the dose-mean lineal energy and microdosimetric spec-trum. Furthermore, since the mean free path for an electron in this energy range isroughly the same size as the sensitive volumes, the variance in the energy depositedby these electrons is reduced. This means that it is appropriate to approximate theircontribution by using the CSDA. In the simulation, the particles trajectory was as-sumed to be a straight line. If a scoring event had already been simulated in asensitive sphere for the current history, the intersection of the high energy electronand that sphere was calculated and this length was multiplied by the stopping powerof this electron (from Ref. [7]) to obtain the energy deposited. If no event had yetbeen recorded in a sensitive sphere for this track, the first sphere that the particle?s93trajectory intersected with was used for scoring. The use of a straight line for theparticle?s trajectory is justified because the packing of the spheres meant that mostof the electrons started inside of a sphere, and for the rest, the mean distance to thenearest sphere was smaller than the mean free path of the electrons. In addition,for tracks in which more than one charged particle was written to the phase spacedescription file for the same history, the distance between the particles was aboutthe same size as the mean free path, which meant that the ?active? sphere would beclose to the starting point of both electrons. Histories with multiple tracks writtento the phase space were a minority in this study (less than 1% of the total), and30% of the particles had a kinetic energy greater than 1 MeV.A VMAT treatment plan for a lung tumor was selected. The central mass of thetumor had soft tissue-equivalent density, which both contained and extended sev-eral voxels beyond the voxel containing the isocenter. This makes it acceptable touse NOREC to simulated the tracks through the voxel without altering the particlefluence. The treatment consisted of two 225? arcs with the MLC and gantry trajec-tories determined by the Eclipse treatment planning system using a 6 MV VarianTrue Beam linac. The total simulation was run in three hours. The condensed-history portion of the simulation was carried out on the 48 node dedicated cluster(see Section 3.6) and the track structure simulation was performed on a 2.0 GHzwindows laptop with 4 GB of memory. Despite the fact that the most computation-ally intensive part of the simulation was carried out on the slowest computer, it stillaccounted for less than 5% of the total calculation time.4.4.2 ResultsThe microdosimetric spectrum for the isocenter of the treatment plan is shownin Figure 4.7. The general shape of this spectrum is typical for a megavoltagephoton beam, but the distribution is much wider than for a monoenergetic photonsource. From this spectrum, the dose-mean lineal energy is calculated to be 2.56keV/?m. This is larger than for 137Cs in the water phantom from Figure 4.6, whichcan be related to the fact that this beam is a continuous bremsstrahlung spectrumas opposed to a monoenergetic source. In addition, the former simulation wasperformed at a depth of 2 cm, but for the VMAT plan much of the beam reached9410?2 10?1 100 10100. Energy (keV/?m)yd(y)Figure 4.7: The simulated microdosimetric spectra for a 6 MV VMAT treat-ment plan in the voxel containing the isocenter in a voxelized CT-basedpatient phantom.the watch volume after being transported through a larger distance in the phantomand the charged particles are set in motion by more highly scattered photons.To check that the handling of the high energy charged particles does not alterthe results, the track structure part of the simulation was re-run by neglecting allparticles with an energy greater than 1 MeV. For this simulation, it was shown thatdose-mean lineal energy was changed by less than 0.1% by these particles. The fre-quency of these particles in the phase space description file is shown in Figure 4.8.These particles did not substantially effect the outcomes of the simulation.4.5 ConclusionThe application of track structure Monte Carlo plans has been demonstrated andshown to be practically achievable for clinical treatment planning purposes. Thiscan accurately and efficiently be done by making use of a watch volume surround-ing any sensitive volumes were microdosimetric information is scored. The watchregion need not be very large to ensure accurate results: it must extend beyond thesensitive volume by a distance equal to the maximum range of a subthreshold elec-950 1 2 3 4 500.0050.010.0150.020.0250.030.035Energy (MeV)FrequencyFigure 4.8: The electron energies in the phase space file presented used forthe track structure simulations for the clinical VMAT plan simulation,presented as a histogram.tron in the condensed history component of the simulation. While neglecting to usea watch volume is shown to be inaccurate, the computational cost of including anacceptably large watch volume is not prohibitive. The approach is validated withexperimental data for photons from a variety of radio-isotope sources.The ultimate objective was to run simulations for actual patient plans. Thiswas carried out effectively for a 6 MV photon beam in a patient CT phantom.This simulation included modeling of the entire linear accelerator treatment head,as well as the full geometry, complete with material and density heterogeneity.This is the first reported results from a multi-scale Monte Carlo simulation for amega-voltage treatment beam with a patient phantom. The total simulation timenecessary to get a statistically significant number of tracks scored increased thetotal simulation time by 0.75%, which demonstrates that these tools can be usedregularly in a clinical setting.96Chapter 5ConclusionsThe objective of this thesis was to create a multi-scale Monte Carlo system forsimulating external beam radiotherapy using specific patient treatment plans andCT-phantoms to calculate microdosimetric quantities. The intent was to developa tool that did not rely on assumptions that negate the validity of the results andstill could be performed without requiring impractical computer resources. Themicrodosimetric characteristics of the radiotherapy treatments ultimately enablethe incorporation of more mechanistic (and therefore predictive) radiobiologicalmodelling. This modelling can improve the optimization and evaluation of treat-ment planning in order to ensure that a patient receives the most effective treatmentpossible while minimizing unnecessary complications.The multi-scale Monte Carlo system implemented in this thesis was based onthe existing condensed history algorithm EGSnrc and the track structure algorithmNOREC. The author developed a software program to oversee these simulations.To make use of the transport simulation codes available, it was necessary to cre-ate software that could score and calculate the quantities of interest, pass particletrajectory information from one code to the other, and to add range rejection tech-niques to improve the efficiency of system. Due to the novelty of this application,each of these parts of the simulation required a new and unconventional approach.In addition, new software tools had to be created that parsed the treatment plansfor the modelling of the linear accelerator. Furthermore, a new type of linear ac-celerator model for EGSnrc++ had to be developed in order to synchronize the97dynamic motion of the linear accelerator, patient, and other accessories with thesimulation. Great effort had to go into producing the software that manages theradiation transport simulations.On the other hand, this effort paid off in a variety of other ways. Many of thesoftware tools developed for this work were also useful in a variety of other clin-ical projects. For instance, the author?s experience with the EGSnrc software andMonte Carlo simulations proved valuable in order to attempt to improve the linearaccelerator model used in the Monte Carlo simulations at the Vancouver CancerCentre. In particular, the currently used model was compared to an alternativemodel proposed by a colleague. The author carried out a series of simulations us-ing both models in order to optimize both with regard to the measured data for thelinear accelerator in question. After this optimization, the results of both modelswere compared with the data to determine which agreed better with the accelera-tors in use at the clinic. In the end, it was determined that the original model was abetter match for the measured data.Another clinical research projected that benefitted from these software toolswas a study of treatment planning algorithms for spinal SBRT. For this project, theauthor had to extensively modify the code developed to simulate the dynamic mo-tion of the linear accelerator in order to apply to the BeamNRC software package.In addition, several additional tools were created that enabled the user to gather amuch wider range of information about the radiation tracks in the simulations. Thiseffort made it possible to make a detailed study of the effects of different treatmentplanning algorithms on the quality of treatment plans generated for spinal SBRT.These results showed that it is possible to use a pencil beam-convolution algorithmto generate treatment plans for this situation, but these plans should ultimately beverified with Monte Carlo simulations to ensure that the delivered treatment meetsthe planning objectives.Other spin-offs from this project were less research based but were importantadvances in the clinical practice. The author developed a set of tools that wereuseful in the planning of total body irradiation. This toolset streamlined and im-proved the process of creating lead blocks to compensate for the lower densityof the lungs. Crucially, this projected was based on patient CTs and not radio-graphic films, which allowed the department to retire the film processor. In ad-98dition, much of the software developed for the TBI project was well suited forthe patient-specific quality assurance for IMRT delivery. The author incorporatedthese tools into the program that oversaw the simulations, as well as extensivelyupdating the rest of the program to adapt to new clinical requirements.Most importantly, the author proposed the concept of the watch volume in orderto connect the two different simulation algorithms. Ensuring an accurate particlespectrum in the scoring targets had been one of the major roadblocks preventing thedevelopment of a rigorous and accurate multi-scale algorithm. By incorporating avariable sized watch volume into the simulation, the author was able to carry out astudy of the effects of the watch volume size on the particle fluence and microdosi-metric spectra in the scoring volumes. From these investigations, it was shown thatthe energy fluence spectrum of particles in the sensitive scoring volumes wouldbe incorrect if the watch volume was not large. However, the watch volume onlyneeded to extend beyond the scoring volumes by a distance equal to the maximumrange of the subthreshold electrons in order to ensure the correct fluence in thevolume of interest. Since the energy threshold can be set to be quite small (of thesame order of magnitude as the size of the sensitive volumes), the inclusion of thewatch volume does not inhibit the potential to use this approach within a reason-able amount of time with a practical amount of computer resources. By using awatch volume, the work presented in this thesis is the first method reported that ac-curately can extract microdosimetric information for mega-voltage photon beamsin inhomogeneous geometries.Once the tools were developed, bench-marking and validation of the modelcould be carried out. The author performed a series of studies of system for differ-ent photon spectra in water phantoms. The aim of this work was to verify that thesimulation results were in agreement with the existing experimental data. It wasdemonstrated that the simulation was capable of reproducing the experimental datafor a variety of different situations and energies, which gives confidence that theresults would be accurate in other more challenging cases. In particular, it was alsodemonstrated that a multi-arc VMAT treatment plan could be simulated in a patientCT-phantom. Importantly, the computational demands of incorporating the trackstructure algorithm were much less than was needed to carry out the condensedhistory simulation alone if there were only a small number of volumes of interest99in which the track structure simulation was performed.These results show that this thesis presents the first tool capable of obtainingmicrodosimetric information directly for mega-voltage radiotherapy treatments,and can do so in a practical and efficient manner. This is a major advance forthe field because there has been a long and heated debate about the relative meritsof different treatment modalities. Much of the disagreement is caused by the lackof concrete data ? all of which has been extrapolated from situations in which therelevance to actual clinical practice is not established. This approach can createknowledge that puts this discussion on a firm foundation which will facilitate animprovement in the understanding of the radiobiology of therapy. By includingthat understanding into the creation and evaluation of treatment plans, the qualityof patient care can improve.On the other hand, the novelty of this approach means that there is as yet noconsensus on the use of microdosimetric information. There are many steps con-necting the microdosimetric track structure with the final outcome to the treatment,few of them adequately understood. However, this only highlights the contributionthe work of this thesis can provide to the field: it is now possible to develop and testhypotheses and models of these stages in vivo by comparing with experiments andclinical data. It will be possible to move beyond simplistic heuristics and embracethe complexities inherent in realistic clinical situations. For instance, there are ac-tually several possibilities for utilizing the microdosimetric information. Perhapsthe most straightforward would be to use the dose-mean lineal energy, y?d , in orderto calculate the RBE using one of several existing models. Then the RBE could beapplied to the absorbed dose to compute the equivalent radiation dose which couldbe used instead of the radiological absorbed dose for treatment planning purposes.This has the virtue of being straightforward and can be quickly calculated becauseit relies on mean values which can be computed precisely without running an ex-travagant number of histories. On the other hand, at this time these models are allempirically based, and so this approach falls short of being fully mechanistic. Fur-thermore, it is unknown whether using the mean values tells enough of the story tomake the most use of the microdosimetric data. It may be that fluctuations in themicrodosimetric event size cause different chemical or biological responses to theradiation that are relevant.100Another alternative is to use a fully mechanistic Monte Carlo simulation ofthe radiation chemistry and biology. Instead of using a quantity like the linealenergy, the code presented in this thesis could instead report the location and energytransfer sizes for each ionization and excitation event in the volumes of interest, anduse these to generate the initial distribution of the free radical and ion species of theradiochemistry stage. This approach has been taken before, but never for a realisticmega-voltage photon beam in a patient CT-phantom. Currently, determining thedistribution of the different species depends on tuning several parameters to matchempirical data. With more attention devoted to this process, however, it is likelythat a canonical model can emerge. Modelling the dynamics of the chemistry,including interaction with embedded DNA molecules, has previously been appliedin simple cases to create the DNA damage spectra. This allows the computationof such information as the rates of double strand breaks and multiple damagedsites as well as strand fragment size spectra. Furthermore, it may be possible toincorporate a stochastic damage repair modelling onto these damage spectra toderive cell survival fractions. This step is more speculative than the others, butis finally amenable to realistic simulations which can contribute to establishingaccurate models. Creating such a simulation tool would allow a fully mechanisticapproach to be used for both research and clinical purposes, and lead to betterindividualization of radiation therapy. It may also be that the full simulation isnot necessary, and that stopping at one of the intermediate stages would providesufficient predictive power.Adopting this full Monte Carlo approach is a very large undertaking. It wouldinvolve collaboration across many fields and most likely large scale studies involv-ing a variety of institutions and a lengthy debate across multiple fields. The con-tribution of this thesis is to begin that process. There are two steps, however, thatare immediately apparent and feasible to undertake immediately. First, incorporat-ing the radiochemical stage modelling is straightforward and mostly understood.The existing code has been implementing specifically with that in mind. The otherfuture direction is to make use of previously used treatment plans and patient CTsfor retrospective studies to correlate a variety of microdosimetric and nanodosi-metric measures with treatment outcomes. This would enable a prospective look atthe possibilities of incorporating microdosimetry into clinical treatment planning.101With the completion of the simulation system created in this thesis, it is possible tobegin work on both of these paths.This thesis presented the implementation of a multi-scale Monte Carlo systemcapable of simulating the track structure of mega-voltage radiation therapy beamsin patient-CT phantoms. The feasibility of this approach was demonstrated, andthe results were validated against empirical data. This is the first such result thathas been reported that does this is a rigorous and accurate manner. In addition, thiswork opens up a realm of possible investigations into improving the radiobiologi-cal modelling of radiation therapy treatment planning. Ultimately, these potentialadvances can help guide clinical practice to ensure that patients receive the optimaltreatments maximizing the opportunity for tumor control and minimizing the riskof complications.102Bibliography[1] S. Agostinelli, J. Allison, J. Amako, J. Apostolakis, and et al. Geant4?asimulation toolkit. Nuclear Intruments in Physics Research A, 506:250?303, 2003. ? pages 22[2] K. Amako, J. Apostolakis, H. Araujo, P. A. Dubois, and M. Asai, et al.Geant4 developments and applications. IEEE Transactions on NuclearScience, 53:270?278, 2006. ? pages 22[3] H. Aouchiche, C. Champion, and D. Oubaziz. Electron and positron elasticscattering in gaseous and liquid water: A comparative study. Radiat. Phys.Chem., 77:107?114, 2008. ? pages 65, 66[4] P. Arce, P. Rato, M. Caadas, and J. I. Lagares. Gamos: An easy and flexibleframework for geant4 simulations. In 2008 Nuclear Science Symposium,Medical Imaging Conference and 16th Room Temperature SemiconductorDetector Workshop, 2008. ? pages 22[5] F. H. Attix. Introduction to Radiological Physics and Radiation Dosimetry.John Wiley and Sons, Inc, Hoboken (NJ), 1986. ? pages 3[6] M. Berger, J. Hubbell, S. Seltzer, J. Chang, J. Coursey, R. Sukumar,D. Zucker, and K. Olsen. Xcom: Photon cross section database (version1.5). Online at http://www.nist.gov/pml/data/xcom/index.cfm, March 26th,2013. ? pages 4, 93[7] M. J. Berger, J. S. Coursey, M. A. Zucker, and J. Chang. Estar, pstar, andastar: Computer programs for calculating stopping-power and range tablesfor electrons, protons, and helium ions (version 1.2.3). Online athttp://physics.nist.gov/Star, March 13th, 2013. ? pages 10, 47, 48, 93[8] N. E. Bolus. Basic review of radiation biology and terminology. J. Nucl.Med. Tech., 29:67?73, 2001. ? pages 17103[9] C. Bousis, D. Emfietzoglou, P. Hadjidoukas, and H. Nikjoo. A monte carlostudy of absorbed dose distributions in both the vapor and liquid phases ofwater by intermediate energy electrons based on differentcondensed-history transport schemes. Phys. Med. Biol., 53:3739, 2008. ?pages 64[10] C. M. Bragg, K. Wingate, and J. Conway. Clinical implications of theanisotropic analytical algorithm for imrt treatment planning andverification. Radiotherapy and Oncology, 86:276?284, 2008. ? pages 36[11] J. F. Briesmeister. MCNP?a general monte carlo n-particle transport code.Technical Report LA-12625-M Version 4B, Los Alamos NationalLaboratory, 1997. ? pages 22[12] C. Champion. Theoretical cross sections for electron collisions in water:structure of electron tracks. Phys. Med. Biol., 48:2147?2168, 2003. ?pages 12, 13, 60, 61[13] C. Champion. Moving from organ dose to microdosimetry: Contribution ofthe monte carlo simulations. Brazilian Archives of Biology and Technology,48:191?199, 2005. ? pages 22[14] C. Champion and C. LeLoirec. Positron follow-up in liquid water: I. a newmonte carlo track-structure code. Phys. Med. Biol., 51:1707?1723, 2006.? pages 23[15] C. Champion, S. Incerti, H. N. Tran, Z. E. Bitar, and theGeant4 DNA Collaboration. Electron and proton elastic scattering in watervapour. Nucl. Instr. and Methods in Phys. Res. B, 273:98?101, 2012. ?pages 12, 65[16] M. Dingfelder. Track-structure simulations for charged particles. HealthPhysics, 103:590?595, 2012. ? pages 60[17] M. Dingfelder, R. H. Ritchie, J. E. Turner, W. Friedland, H. G. Paretzke,and R. N. Hamm. Comparison of calculations with partrac and norec:Transport of electrons in liquid water. Radiat. Res., 169:584 ? 594, 2008.? pages 61[18] M. Dizdaroglu, P. Jaruga, M. Birincioglu, and H. Rodriguez. Freeradical-induced damage to dna: Mechanisms and measurement. FreeRadical Biology and Medicine, 32:1102?1115, 2002. ? pages 15104[19] N. Dogan, P. J. Keall, and J. V. Siebers. Clinical comparison of head andneck and prostate IMRT plans using absorbed dose to medium and absorbeddose to water. Phys. Med. Biol., 51:4967?4980, 2006. ? pages 44, 48[20] D. Emefietzoglou, G. Papamichael, K. Kosterelos, and M. Moscovitch. Amonte carlo tack structure code for electrons ( 10 ev-10kev) and protons( 0.3-10 mev) in water: partitioning of energy and collision events. Phys,Med. Biol., 45:3171 ? 3194, 2000. ? pages 62[21] D. Emfietzoglou and H. Nikjoo. The effect of model approximations onsingle-collision distributions of low-energy electrons in liquid water.Radiat. Res., 163:98?111, 2005. ? pages 61, 64[22] A. Ferrari, P. R. Sala, A. Fasso`, and J. Ranft. Fluka. CERN-library, page(http://fluka. web. cern. ch/fluka), 2005. ? pages 22[23] M. Fippel and F. Nusslin. Comments on ?converting absorbed dose tomedium to absorbed dose to water for monte carlo based photon beam dosecalculations?. Phys. Med. Biol., 45:L17?L18, 2000. ? pages 44, 48[24] B. A. Fraass. The development of conformal radiation therapy. MedicalPhysics, 22:1911?1921, 1995. ? pages 19[25] Z. Francis, S. Incerti, R. Capra, B. Mascialino, G. Montarou, V. Stepan, andC. Villagrasa. Molecular scale track structure simulations in liquid waterusing the geant4-dna monte-carlo processes. Appl. Rad. and Isotopes, 69:220?226, 2011. ? pages 63, 66[26] Z. Francis, S. Incerti, M. Karamitros, H. Tran, and C. Villagrasa. Stoppingpower and ranges of electrons, protons and alpha particles in liquid waterusing the geant4-dna package. Nucl. Instr. and Meth. in Phys. Res. B: BeamInteractions with Materials and Atoms, 269:2307 ? 2311, 2011. 12thInternational Conference on Nuclear Microprobe Technology andApplications. ? pages 64, 67, 68[27] W. Friedland, P. Jacob, P. Bernhardt, H. G. Paretzke, and M. Dingfelder.Simulation of DNA damage after proton irradiation. Radiat Res, 159:401?410, 2003. ? pages 22[28] W. Friedland, P. Jacob, H. F. Paretzke, A. Ottolenghi, F. Ballarini, andM. Liotta. Simulation of light ion induced DNA damage patterns. Radiat.Prot. Dosim., 122:116?120, 2006. ? pages 22105[29] W. Friedland, M. Dingfelder, P. Kundra?t, and P. Jacob. Track structures,dna targets and radiation effects in the biophysical monte carlo simulationcode PARTRAC. Mutation Research, 711:28?40, 2011. ? pages 23[30] M. C. Fuss, A. Mun?oz, J. C. Oller, and F. Blanco, et al. Energy depositionby a 106ru/106rh eye applicator simulated using LEPTS, a low energy particletrack simulation. Applied Radiation and Isotopes, 69:1198?1204, 2011. ?pages 23[31] J. D. Gorfinkiel, L. A. Morgan, and J. Tennyson. Electron impactdissociative excitation of water within the adiabatic nuclei approximation.J. Phys. B: At. Mol. Opt. Phys, 35:543?555, 2002. ? pages 13[32] B. Grosswendt. Formation of ionization clusters in nanometric structuresof propane based tissue-equivalent gas or liquid water by electrons andalpha-particles. Radiat. Environ. Biophys., 41:103?112, 2002. ? pages 22[33] S. Hara. The scattering of slow electrons by hydrogen molecules. J. Phys.Soc. Japan B, 22:710?718, 1967. ? pages 12[34] K. Hensley, K. A. Robinson, S. Gabbita, S. Salsman, and R. A. Floyd.Reactive oxygen species, cell signaling, and cell injury. Free RadicalBiology and Medicine, 28:1456 ? 1462, 2000. ? pages 16[35] S. Incerti, N. Gault, C. Habchi, J. Lefaix, P. Moretto, J. Poncy, T. Pouthier,and H. Seznec. A comparison of cellular irradiation techniques with alphaparticles using the geant4 monte carlo simulation toolkit. Radiat. Prot.Dosim, 122:327?329, 2006. ? pages 23, 24[36] S. Incerti, G. Baldacchino, M. Bernal, R. C. R, C. Champion, and et al.The geant4-dna project. Int. J. Model. Simul. Sci. Comput., 1:157?178,2010. ? pages 22, 61[37] Microdosimetry ICRU Report 36. International Commission on RadiationUnits and Measurements (ICRU), 1983. ? pages 14, 18, 21, 60, 61, 91[38] Y. Itikawa and N. Mason. Cross sections for electrons with watermolecules. J. Phys. Chem. Ref. Data, 34:1?22, 2005. ? pages 12, 60[39] S. Y. Jang, H. H. Liu, R. Mohan, and J. V. Siebers. Variations in energyspectra and water-to-material stopping-power ratios in three dimensionalconformal and intensity modulated photon fields. Med. Phys., 34:1388?1397, 2007. ? pages 48106[40] H. E. Johns and J. T. Cunningham. The Physics of Radiology. Charles C.Thomas, Springfield (Illinois), 4th edition, 1983. ? pages 3, 9, 10[41] B. L. Jones, S. Krishnan, and S. H. Cho. Estimation of microscopic doseenhancement factor around gold nanoparticles by monte carlo calculations.Med. Phys., 37:3809, 2010. ? pages 24[42] I. Kawrakow and A. F. Bielajew. On the condensed history technique forelectron transport. Nuclear Instruments and Methods in Physics ResearchSection B: Beam Interactions with Materials and Atoms, 142:253?280,1998. ? pages 22[43] I. Kawrakow, E. Mainegra-Hing, D. W. O. Rogers, F. Tessier, and B. R. B.Walters. The egsnrc code system: Monte carlo simulation of electron andphoton transport. Technical Report PIRS-701, National Research Councilof Canada (NRCC), 2011.(http://irs.inms.nrc.ca/software/egsnrc/documentation/pirs701.pdf). ?pages 22, 51[44] I. Kawrakow, E. Mainegra-Hing, F. Tessier, and B. R. B. Walters. Theegsnrc c++ class library. Technical Report PIRS-898, National ResearchCouncil of Canada (NRCC), 2013.(http://irs.inms.nrc.ca/software/egsnrc/documentation/pirs898/index.html).? pages 56, 85[45] F. M. Khan. The Physics of Radiation Therapy. Lippincott Williams andWilkins, Philadephia (Pa), 3rd edition, 2003. ? pages 5, 6[46] S. Kim. Characteristics of eliiptical sources in beamnrc monte carlosystem: Implementation and application. Med. Phys., 36:1046?1052, 2009.? pages 29[47] E. E. Klein, J. Hanley, J. Bayouth, and F. Yin, et al. Task group 142 report:Quality assurance of medical accelerators. Med. Phys., 36:4197?4212,2009. ? pages 21, 31[48] P. Kliauga and R. Dvorak. Microdosimetric measurements of ionization bymonoenergetic photons. Radiat. Res., 73:1?20, 1978. ? pages 91[49] G. J. Kutcher, L. Coia, M. Gillin, and W. F. Hanson, et al. ComprehensiveQA for radiation oncology: Report of AAPM radiation therapy committetask group 40. Med. Phys., 21:581?618, 1994. ? pages 21107[50] A. V. Lappa, E. A. Bigildeev, D. S. Burmistrov, and O. N. Vasilyev. TRIONcode for radiation action calculations and its application in microdosimetryand radiobiology. Radiat. Environ. Biophyis., 31:1?19, 1993. ? pages 22[51] P. Lazarakis, M. Bug, E. Gargioni, S. Guatelli, H. Rabus, and A. Rosenfeld.Comparison of nanodosimetric parameters of track structure calculated bythe monte carlo codes geant4-dna and ptra. Phys. Med. Biol., 57:1231 ?1242, 2012. ? pages 61[52] S. Lehnert, (Ed.). Biomolecular Action of Ionizing Radiation. Taylor andFrancis, 1st edition, 2008. ? pages 17, 18[53] T. Liamsuwan, D. Emfietzoglou, S. Uehara, and H. Nikjoo.Microdosimetry of low-energy electrons. Int. J. of Radiat. Biol., 88:899 ?907, 2012. ? pages 81, 82[54] H. H. Liu and F. Verhaegen. An investigation of energy spectrum and linealenergy variations in mega-voltage photon beams used for radiotherapy.Radiat. Prot. Dosim., 99:425?427, 2002. ? pages 24[55] J. Lobo and I. A. Popescu. Two new dosxyznrc sources for 4d monte carlosimulations of continuously variable beam configurations, withapplications to rapidarc, vmat, tomotherapy and cyberknife. Phys. Med.Biol., 55:4431?4443, 2010. ? pages 33, 38[56] D. A. Low, J. M. Moran, J. F. Dempsey, L. Dong, and M. Oldham.Dosimetry tools and techniques for IMRT (AAPM task group report 120).Med. Phys., 38:1313?1339, 2011. ? pages 25[57] J. Lucido, I. A. Popescu, and V. Moiseenko. A microdosimetriccomparison of geant4-dna and norec for electrons in water. Phys. Med.Biol., 1987 (submitted). ? pages 61[58] J. Lucido, I. A. Popescu, and V. Moiseenko. A method to performmulti-scale monte carlo simulations in the clinical setting. Phys. Med.Biol., 2013 (submitted). ? pages 84[59] L. B. Marks, E. D. Yorke, A. Jackson, R. K. Ten Haken, L. S. Constine,A. Eisbruch, S. M. Bentzen, J. Nam, and J. O. Deasy. Use of normal tissuecomplication probability models in the clinic. Int. J. of Radiat. Onco. Biol.Phys., 76:S10?S19, 2010. ? pages 20108[60] J. Meesungneon, J.-P. Jay-Gerin, A. Filali-Mouhim, and S. Mankhetkorn.Low energy electron penetration range in liquid water. Radiat. Res., 158:657 ? 660, 2002. ? pages 67, 68[61] B. Michael, K. Prise, M. Folkard, B. Vojnovic, B. Brocklehurst, I. Munro,and A. Hopkirk. Action spectra for single-and double-strand breakinduction in plasmid dna: Studies using synchrotron radiation. Int. J. Rad.Biol., 66:569?572, 1994. ? pages 62[62] V. Michalik. Model of DNA damage induced by radiations of variousqualities. Int. J. of Rad. Biol., 62:9?20, 1992. ? pages 23[63] V. Moiseenko, M. Liu, S. Louwen, R. Kosztyla, E. Vollans, J. Lucido,M. Fong, R. Vellani, and I. A. Popescu. Monte carlo calculation of dosedistribution in oligometastatic patients planned for spine stereotacticablative radiotherapy. Phys. Med. Biol., 2013. ? pages 38[64] H. Nikjoo and D. Goodhead. Track structure analysis illustrating theprominent role of low-energy electrons in radiobiological effects of low-letradiations. Phys. Med. Biol., 36:229, 1991. ? pages 22[65] H. Nikjoo and L. Lindborg. Rbe of low energy electrons and photons.Phys. Med. Biol., 55:R65, 2010. ? pages 62[66] H. Nikjoo, D. T. Goodhead, D. E. Charlton, and H. G. Paretzke. Energydeposition in small cylindrical targets by ultrasoft x-rays. Phys. Med. Biol.,34:691?705, 1989. ? pages 18, 23, 72[67] H. Nikjoo, S. Uehara, W. E. Wilson, M. Hoshi, and D. T. Goodhead. Trackstructure in radiation biology: theory and applications. Int. J. Radiat. Biol.,73:355?364, 1998. ? pages 14, 22, 23, 24[68] H. Nikjoo, P. O?Neill, M. Terrissol, and D. Goodhead. Quantitativemodelling of dna damage using monte carlo track structure method.Radiat. Environ. Biophys., 38:31?38, 1999. ? pages 23, 24[69] H. Okamoto, T. Kanai, Y. Kase, and Y. Matsumoto, et al. Relation betweenlineal energy distribution and relative biological effectiveness for photonbeams according to the microdosimetric kinetic model. J. Radiat. Res., 52:75?81, 2011. ? pages 14, 21, 22[70] P. O?Neill and P. Wardman. Radiation chemistry comes before radiationbiology. Int. J. Radiat. Biol., 85:9?25, 2009. ? pages 14, 15, 16109[71] K. Otto. Volumetric modulated arc therapy: IMRT in a single gantry arc.Med. Phys, 35:310, 2008. ? pages 19[72] A. Ottolenghi, M. Merzagora, L. Tallone, M. Durante, H. Paretzke, andW. Wilson. The quality of dna double-strand breaks: a monte carlosimulation of the end-structure of strand breaks produced by protons andalpha particles. Radiat. Environ. Biophys., 34:239?244, 1995. ? pages 23,24[73] H. Paganetti, A. Niemierko, M. Ancukiewicz, L. E. Gerweck, M. Goitein,J. S. Loeffler, and H. D. Suit. Relative biological effectiveness (rbe) valuesfor proton beam therapy. International Journal of Radiation Oncology*Biology* Physics, 53(2):407?421, 2002. ? pages 21[74] H. G. Paretzke. Radiation track structure theory. In G. R. Freeman, editor,In Kinetics of NonHomogeneous Processes, pages 89?170. Wiley, 1987. ?pages 22[75] E. Pomplun and M. Terrissol. Low-energy electrons inside active DNAmodels: a tool to elucidate the radiation action mechanisms. Radiat. Prot.Biophys., 33:279?292, 1994. ? pages 23[76] Y. Qiu, I. A. Popescu, C. Duzenli, and V. Moiseenko. Mega-voltage versuskilo-voltage cone bone ct used in image guided radiation therapy:comparative study of microdosimetric properties. Radiat. Prot. Dosimetry,143:477?480, 2011. ? pages 23[77] U. Quast. Total body irradiation ? review of treatment techniques ineurope. Radiotherapy and Oncology, 9:91 ? 106, 1987. ? pages 49[78] G. Raisali, L. Mirzakhanian, S. F. Masoudi, and F. Semsarha. Calculationof dna strand breaks due to direct and indirect effects of auger electronsfrom incorporated 123i and 125i radionuclides using the geant4 computercode. Int. J. Rad. Biol., 89:57?64, 2013. ? pages 23[79] B. Reniers, S. Vynckier, and F. Verhaegen. Theoretical analysis ofmicrodosimetric spectra and cluster formation for 103pd and 125i photonemiitters. Phys. Med. Biol., 49:3781?3795, 2004. ? pages 23[80] P. A. Riley. Free radicals in biology: oxidative stress and the effects ofionizing radiation. Int. J. Radiat Biol., 65:27?33, 1994. ? pages 15, 16110[81] R. Ritchie, R. Hamm, J. Turner, H. Wright, and W. Bolch. Radiationinteractions and energy transport in the condensed phase. In Physical andChemical Mechanisms in Molecular Radiation Biology, pages 99?135.Springer, 1992. ? pages 62, 63, 64[82] D. W. O. Rogers, B. Walters, and I. Kawrakow. Beamnrc users manual.Technical Report PIRS-0509(A)revL, National Research Council ofCanada (NRCC), 2011.(http://irs.inms.nrc.ca/software/beamnrc/documentation/pirs0509/). ?pages 22, 28[83] S. Rollet, P. Colautti, B. Grosswendt, D. Moro, E. Gargioni, V. Conte, andL. DeNardo. Monte carlo simulation of mini tepc microdosimetry spectra:influence of low energy electrons. Radiation Measurements, 45:1330?1333, 2010. ? pages 22, 23[84] R. K. Sachs, A. M. Chen, and D. J. Brenner. Review: Proximity effects inthe production of chromosome abberations by ionizing radiation. Int. J.Radiat. Biol, 71:1?19, 2001. ? pages 18[85] F. Salvat. NIST Electron Elastic-Scattering Cross-Section Database. PhDthesis, Universitat de Barcelona, 2010. ? pages 65[86] M. Scholz and G. Kraft. Calculation of heavy ion inactivation probabilitiesbased on track structure, x ray sensitivity and target size. Radiat. Prot.Dosim., 52:29?33, 1994. ? pages 23, 24[87] V. A. Semenenko, J. E. Turner, and T. B. Borak. NOREC, a monte carlocode for simulating electron tracks in liquid water. Radiat Environ.Biophys., 42:213?217, 2003. ? pages 22, 61, 63, 66[88] D. Sheikh-Bagheri and D. W. O. Rogers. Monte carlo calculation of ninemegavoltage photon beam spectra using the beam code. Med. Phys., 29:391?402, 2002. ? pages 28[89] D. Sheikh-Bagheri and D. W. O. Rogers. Sensitivity of megavoltagephoton beam monte carlo simulations to electron beam and otherparameters. Med. Phys., 29:379?390, 2002. ? pages 28[90] J. V. Siebers, P. J. Keall, A. E. Nahum, and R. Mohan. Converting absorbeddose to medium to absorbed dose to water for monte carlo based photonbeam dose calculations. Phys. Med. Biol., 45:983?995, 2000. ? pages 40,44, 47, 48111[91] E. Sterpin, M. Tomsej, B. De Smedt, N. Reynaert, and S. Vynckier. Montecarlo evaluation of the AAA treatment planning algorithm in aheterogeneous multilayer phantom and IMRT clinical treatments for anelekta sl25 linear accelerator. Med. Phys., 34:1665?1677, 2007. ? pages27[92] A. Syme, C. Kirkby, R. Mirzayans, M. MacKenzie, C. Field, and B. G.Fallone. Relative biological damage and electron fluence in and out of a 6mv photon field. Phys. Med. Biol., 54:6623?6633, 2009. ? pages 21[93] M. Taguchi, A. Kimura, R. Watanabe, and K. Hirota. Estimation of yieldsof hydroxyl radicals in water under various energy heavy ions. Radiat.Res., 171:254?263, 2009. ? pages 23, 24[94] J. Terrien. News from the bureau international des poids et mesures.Metrologia, 11:179?183, 1975. ? pages 25[95] R. M. Thomson and I. Kawrakow. On the monte carlo simulation ofelectron transport in the sub-1 kev energy range. Medical physics, 38:4531?4536, 2011. ? pages 62[96] E. Torres-Garcia, H. M. Garnica-Garza, and G. Ferro-Flores. Monte carlomicrodosimetry of 188RE- and 131I-labelled anti-CD20. Phys. Med. Biol.,51:N349?N356, 2006. ? pages 23[97] J.-i. Ueda, N. Saito, Y. Shimazu, and T. Ozawa. A comparison ofscavenging abilities of antioxidants against hydroxyl radicals. Archives ofbiochemistry and biophysics, 333:377?384, 1996. ? pages 15[98] S. Uehara, H. Nikjoo, and D. T. Goodhead. Cross-sections for water vapourfor monte carlo track structure code from 10 ev to 10 mev region. Phys.Med Biol., 38:1841?1858, 1993. ? pages 22[99] M. Valko, C. Rhodes, J. Moncol, M. Izakovic, M. Mazur, et al. Freeradicals, metals and antioxidants in oxidative stress-induced cancer.Chemico-biological interactions, 160:1?40, 2006. ? pages 15[100] M. Varma, J. Baum, P. Kliauga, and V. Bond. Microdosimetric parametersfor photons as a function of depth in water using wall-less and walledcounters. Radiat. Res., 88:466?475, 1981. ? pages 91, 92[101] B. Walters, I. Kawrakow, and D. W. O. Rogers. Dosxyznrc users manual.Technical Report PIRS-794revB, National Research Council of Canada112(NRCC), 2011.(http://irs.inms.nrc.ca/software/beamnrc/documentation/pirs794/). ?pages 22, 28, 38, 40, 41, 42[102] S. Webb. The physical basis of imrt and inverse planning. British Journalof Radiology, 76:679?689, 2003. ? pages 19[103] S. Webb and A. E. Nahum. A model for calculating tumor controlprobability in radiotherapy including the effects of inhomogeneousdistributions of dose and clonogenic cell density. Phys. Med. Biol., 38:653?666, 1993. ? pages 20[104] W. E. Wilson and H. Nikjoo. A monte carlo code for positive ions tracksimulation. Radiat. Environ. Biophyis., 38:97?104, 2000. ? pages 22[105] C. C. Winterbourn. Reconciling the chemistry and biology of reactiveoxygen species. Nature chem. biol., 4:278?286, 2008. ? pages 14, 15, 16113


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items