Monte Carlo Simulations For Brachytherapy

Transcription

Monte Carlo simulations forbrachytherapyD.W.O. RogersCarleton Laboratory forRadiotherapy Physics.Physics Dept,Carleton UniversityOttawa, Canadahttp://www.physics.carleton.ca/ drogersICTP,Trieste, Nov 15, 20071

BrachytherapyRadiation therapy given at a shortdistance Radioactive seeds or sourcesimplanted either permanently ortemporarily inside tumours Commonly used for prostate, breast,cervical and ocular cancers.2/69

Typical Prostate ImplantProstate implant of 150125Ibrachytherapy seeds3/69

Characteristics of isotopes used inbrachytherapymissingCs-131Yb-1690.92-1.0120.63- 0.69from IAEA Review of Radn Oncology Physics, E Podgorsak editor 4/69

I-125 seeds cause hot spots125 seeds 5mm spacing -along length of seeds5/69

10 mm spacing 9 to 1 ratio6/69

90% of volume receives 140 GyDVHdose volumehistogram5125mmseeds,spacing5 mm spacingfraction ofvolume ofinterest(eg prostate)receivingat least thedose on the xaxisDetails depend on size of voxel in calculation7/69

Large variety of seeds (4-5mm overall)125IAmersham 6711103PdBest 2335103PdTheragenics 200125ISource TechMedical STM125I125IImagynisoStar 125018/69

TG-43 coordinate system for calculationsP(ro,θo) P(1 cm,π/2) xyP(r,θ)β θ2 θ1rθ1θθ2zLTG43 references at end of slides9/69

What are pure geometry effects? Consider the source in vacuum. What is thefluence of particles at (r,θ)?just 1/r2 for a point source for a line sourceor(faster)10/69

How big is the geometry factor?Pure 1/r2 would beunity in both plots11/69

TG43Air-kerma strength1U 1mGy m2h-1 1 cGy cm2h-1-the air kerma rate at distance d-in vacuum (i.e., no air attenuation or scatter)-for photons δ (no effect on dosebut has effect on K, typically 5 keV)SK is defined to be independent of d (for d L).Numerically reference air kerma rate of ICRU 38/60This is a measured quantity.12/69

TG43Dose-rate constant (DRC)ro 1 cm, θo 90ocGyh-1/U cm-2depends on radionuclide & source modelDRCs are measured or calculated by Monte Carlo13/69

TG43Radial dose functionX P or L (point or line)g(ro) 1Aside from geometry (1/r2), g(r) accounts fordose fall off on transverse axisdue to photon attenuation and scatter.For a point source:Radial dose functions are measured andcalculated by Monte Carlo14/69

TG43Dose anisotropy functionF(r,θ) describes variation in dose vs polar anglerelative to transverse plane (above GL variations)F(r,θo) 1.0Dose anisotropy functions are mostly calculated byMonte Carlo15/69

TG43i.e. the TG-43 formalism is an identity,- holds independently of what G(r,θ) is.Clinically, get measured SK for your individual seedsand use Monte Carlo or measured values forΛ, g(r), and F(r,θ)16/69

Monte Carlo for TG43 parameters Jeffrey F Williamson - PTRAN CCG– this code is responsible for most of thepublished values (MCNP and EGS4/EGSnrc havealso been used)– references at end Monte Carlo for brachytherapy is, in general aphoton transport problem since so low in energy Dose Kerma(strictly the collision kerma)17/69

Kerma estimatorsanalog Ejtracklengthtj18/69

Kerma estimators Analog estimatori-th regionj-th interaction region i Tracklength estimator-uses relationships between fluence φand total tracklength per unit volumetj tracklengthsin i-th regionand between kerma and fluenceTracklength estimator: much more efficient - it usesall paths through a volume, not just interactions in it19/69

Analog vs pathlength estimatorsanalog estimator (interaction) contributes to 1 voxelpathlength estimator contributes to 11 voxels20/69

Efficiency increase from tracklength estimator125125Iseeds,5 mm spacing2 mm voxelsCPU time sameBrachyDosecodeuse history by history scoring21/69

Kerma h Ejtjtj22/69

Exponential tracklength estimator tracklength estimator scores for everytrack which actually crosses the region exponential tracklength estimator scores onevery track which would have crossed the region if interaction in or past volume ‘thin voxel’ for interaction in a thin voxel, same astracklength estimator23/69

Next flight h Ejtjdjtjnext flight24/69

Next flight estimators (cont)Estimates the kerma at a point - previous are voxel averagesAt each interaction, determine angle, θi, to point of interest& probability per unit solid angle of a photon of energy Ei-1scattering to that direction. Energy-angle relationship E(θi)Problem: interactions for which di near 0 large fluctuations25/69

Bounded next flight (BNF) estimatorsConsider a sphere of radius R around point of interestdi Rdi Rwhere we have replaced e-µd/d2 with its mean value within R,assuming there are uniform interactions within R.This is the workhorse method used in Williamson’s code forTG43 parameters.26/69

Once more collided flux estimators (OMCFE) BNF requires uniform interactions within R– does not work close to small sources with steepgradients– does not work near interfaces OMCFE is basically a next flight estimator– resamples a step if it ends within R– changes the weight and then does next flightestimate for that collision– but tracks the history using the originallysample collision point27/69

No photo-electric absorption improve efficiency by not terminating histories viaphoto-electric events– sample distance to interaction from σt– select interaction type ignoring σpe– reduce weight probability particle has escaped aterminal event– energy deposition consider compton scatter– other estimators, multiply by Wi28/69

Calculating TG43 parameters: SKKair at 10 cm on the bisector, in vacuum issues– air attenuation or not? - definition is in vacuum– low energy (4.5 keV) Ti x-rays no effect in water large effect on Kair since (µen/ρ) large at 5 keV decision to exclude using thin Al sheet– on-axis or area of wide angle free air chamber (WAFAC) definition, TLDs, clinical ion chambers all on-axis primary standard is WAFAC29/69

Calculating TG43 parameters: g(r),F(r,θ)Kwater issues.throughout the medium D(r,θ)ro 1 cm,θo 90o– what G(r,θ) to use: line or point: consistency critical– voxel size effects if using voxels– what is the medium? water for the TG43 parameters (is this the best?) experiments in plastics/”solid-virtual water” etc– size of phantom (must be big enough)(see Melhus & Rivard Med Phys 33 (2006)1729-1737)30/69

g(r) – radial dose functionI-125Pd-103 3%,6%,16%Taylor et al, Med Phys 34(2007) 445-45731/69

Effect of cross sectionsDLC-189/200in MCNP4B/CIplant 3500 I-125µ low by 5-6%DeMarco et alPMB 47 /146/174)Λ is 2% too high for this seed with MCNP32/69

Improved fitting function for g(r) dataTaylor and Rogers, in preparation, 200733/69

Improved fitting function for g(r) dataIf fitting from 5 mm to 10 cm, then 5-th orderpolynomial is perfect.34/69

Improved fitting function for g(r) dataThe newfittingfunctionhandles theupturnUpturn is from Ti x-rays and/or from geometryeffects35/69

F(r,θ): Theragenics 200 Pd-103103PdSourceTheragenicsModel 200Taylor et al, Med Phys 34(2007) 445-45736/69

Λ – dose-rate constantDose to water at reference point divided by air-kerma strength103PdSourceTheragenicsModel 200scoring plane at 10 cm from sourceTaylor et al, Med Phys 34(2007) 445-45737/69

BrachyDose a code for calculating dose distributions in a voxel CTphantom with brachytherapy sources code is based on Yegin’s MG (multi-geometry) packageand EGSnrc detailed models of multiple seed types at multiplelocations and orientations can handle x-ray sources since based on EGSnrc code is very fastSimilar code: MCPI Chibani & WilliamsonMed Phys 32(2005)3688-369838/69

Yegin’s Multi-Geometry package39/69

Features of BrachyDose have modelled most (all?) major seeds ability to extract TG-43 parameters easily– benchmarked against previous TG-43 work– studied voxel size effects (Taylor et al, Med Phys34(2007) 445-457) ability to use a CT dataset scores spectrum in any voxel models 50 kV brachytherapy x-ray source– VRT increased efficiency by factor of 104 ability to use a phase space source models Ir-192 HDR in catheters adjust cross sections for sensitivity analysis40/69

BrachyDose: speed-up features use pathlength scoring (factor of 15 over analog) reuse photons escaping from one seed for everyseed (20% for 1mm, 28% for 3mm) fast table lookups of mass energy absorptioncoefficients (5% speed up)– use same approach as EGSnrc (log interpolation,small energy increments)– use XCOM data base for EGSnrc virtual geometry to reduce checks on seeds (largespeed up)41/69

BrachyDose: timing for average uncertainty of 2% on dosein each voxel in implant(g77 compiler - single 3.0 GHz 64 bit CPU) 3 mm voxels: 14 s 2 mm voxels: 30 s 1 mm voxels: 228 si.e.,3.5 min for 1 mm voxels and with an Intel compiler it is expected to befaster based on previous experience42/69

Voxel size effects Williamson’s code uses a point estimator when using voxels one must be carefuli.e., dose is averaged over voxelKawrakowMP 33(2006)1829For 1/r243/69

Voxel size effectsSize of error in associating the average dose in the voxelwith the dose at the mid-point (near a point source)44/69

Voxel size effects125ISource TechMedical STM125ITaylor et al,Med Phys34 (2007)445 - 45745/69

Voxel size effects125ISource TechMedical STM125I46/69

Voxel size effectsrecall slide 7 showed voxel size effects on DVHs47/69

Inter-seed effects in brachytherapyStandard planning systemsignore self-shieldingA few papers have investigatedinter-seed effects.Chibani and Williamson, (MedPhys 2005) used MCNP. Eachcalculation took 3 days.48/69

central axis dose for 5 mm spacing49/69

inter-seed effects doses vary by as much as– 27% in valley– 10% near peak can be misleading since volumes near peaks areoften small in these regular patterns close to worst case (seeds are 3.6% of volumecompared to 0.5% for 10 mm spacing)50/69

inter-seed effects (cont) total energy deposited in implant– 5 mm spacing -down 13%– 10 mm spacing -down 2.8% total energy deposited in nearby “critical organ”– 5 mm spacing -down 13%– 10 mm spacing -down 3% difference is more important than the absolute valuesince dose prescription determined by experience51/69

192Irsource in catheters can model an HDR source inside a catheter the source moves to an arbitrary number of positionsin an arbitrary number of catheters we see no significant effect of the catheters.52/69

Need for treatment planning TG-43 treatment planning assumes infinite waterpatient with one source at a time what are the magnitudes of the effects of– other tissues and calcifications– air surfaces (for breast treatments)– inter-seed attenuation (up to 5-10%) treatment planning with BrachyDose based on 3 DCT data sets handles all of these53/69

Effect of finite irradiated volume50kV X-Ray SourceTG-430cm 2cm15cmSimilar effects for125I& 103Pd seeds 54/69

Effect of realistic tissue50kV X-Ray Source50 kVsource55/69

CLRP TG-43 web resourcehttp://www.physics.carleton.ca/clrp/seed database Full set of tabulateddosimetry data foreach source studied Description ofcalculation methods Descriptions & scaledrawings of sources Plots compare valuescalculated in this studyto other results Links to relevantpapers & websites56/69

CLRP TG-43 web resource57/69

Xoft Axxent X-ray source Miniature X-Ray Source 25-50 kV e- on W target. Lifetime of 5h w/ beam on “Tunable” dose rates Reduced shielding requirements Source can be turned on/off FDA Approved – Jan 2006 EGSnrc can handle this unlikeother brachytherapy codes58/69

Model for Monte Carlo simulationCooling sheathX-Ray SourceMetal BackingHV cableCooling Sleeve50kv e-Thin W Target59/69

Radial dose function comparison192Ir125IAxxent 50kv60/69

g(r) – radial dose function61/69

Energy spectrum with EIIspectrumcalculatedin airincludingdetectorresponseexplains7 keV“peak”62/69

F(r,θ) – anisotropy function63/69

BCSEbrem cross section enhancementlow-energyelectrons producevery few bremphotons waste timetracking electronswith no effectSolution: enhancethe cross sectionfor bremAli & Rogers, Med Phys 34 (2007) 2143-215464/69

BCSEUBS uniformbrem splittingUBS “saturates”sooner becausemany photons inthe same history.Optimal approachis hybrid(500,100)1 mm3 ε 15,000Ali & Rogers, Med Phys 34 (2007) 2143-215465/69

Phase-space source electrons take much longer to transportthan photons in the patient we only care about photons score the phase space data outside the source andinside the “balloon” and re-use saves a factor of 4 in the calculation will be able to recycle as in BEAMnrc not needed for seeds since already reuse eachphoton with every seed66/69

Summary Monte Carlo has played a central role in brachytherapydosimetry, especially in calculating TG43 parameters Monte Carlo dose calculations are needed in brachytherapy– interseed attenuation-realistic materials-lack of full backscatter BrachyDose does a 2% cal’n in 1 mm voxels in 4 min.– extensive benchmarking for different seeds– EGSnrc inherently can model the x-ray sources67/69

AcknowledgmentsMembers of the Carleton Laboratory for RadiotherapyPhysics gratefully acknowledge the support of:NSERCCanada Research Chairs programCFIOITOGSSTTomoTherapy MDS-Nordion Varian NucletronElsayedAliGultekinYeginRandyTaylor68/69

Thank you69/69

Bibliography re Williamson’s work PTRAN CCG (Williamson group’s code)– Monte Carlo evaluation of kerma at a point for photontransport problems Med Phys 14(1987)567-576– Volume-based geometric modeling for radiationtransport calculations, Med Phys 19 (1992) 667-677– Comparison of calculated and measured heterogeneityfactors Med Phys 20 (1993) 209-222– Monte Carlo calculations of kerma to a point in thevicinity of media interfaces, PMB 38 (1993) 1825-1840 Quantitative dosimetry methods for brachytherapy,Williamson and Rivard, pp 233-294 in AAPM 2005 SummerSchool: Brachytherapy Physics, Medical Physics Publishing70/69

Bibliography re TG43-Dosimetry of interstitial brachytherapy sources:Recommendations of the AAPM Radiation TherapyCommittee Task Group 43Med Phys 22 (1995) 209-234-Update of AAPM TG43 Report: A revised protocol forbrachytherapy dose calculationsMed Phys 31 (2004) 633 -674-Supplement to the 2004 update of the AAPM TG43 ReportMed Phys 34 (2007) 2187-220571/69

Monte Carlo for TG43 parameters Jeffrey F Williamson - PTRAN_CCG - this code is responsible for most of the published values (MCNP and EGS4/EGSnrc have also been used) - references at end Monte Carlo for brachytherapy is, in general a photon transport problem since so low in energy Dose Kerma (strictly the collision kerma)