Monte Carlo: Electrons and Brachytherapy

Simulating electron transport in clinical linear accelerators and computing dose distributions in brachytherapy rank among the most demanding challenges in computational medical physics. The Monte Carlo method provides the most accurate tool for solving both, overcoming fundamental limitations of analytical algorithms. In this article, we explore how Monte Carlo modelling has evolved for external electron beams and for brachytherapy source dosimetry — two applications requiring distinct levels of rigour.

Series overview: for the full roadmap and related articles, return to the complete guide on Monte Carlo in radiotherapy.

For a comprehensive overview of all Monte Carlo applications, see our complete guide to Monte Carlo techniques in radiotherapy.

Clinical Electron Beams: Why Monte Carlo Is Essential

Monte Carlo simulations are indispensable for dose calculation in electron beams — even more so than for photons. The reason is straightforward: analytical dose-calculation algorithms for electrons show significant errors in irregular fields and heterogeneous targets, even in modern treatment-planning systems. While alternative methods can rival Monte Carlo accuracy for photons, the situation is different for electrons, with numerous studies exposing flaws in conventional methods.

Clinical electron beams typically use energies ranging from 4 to 20 MeV. The LINAC configuration for electrons differs from photon mode in three key aspects: there is no photon conversion target, a scattering foil broadens the beam, and a multi-stage collimator (applicator) shapes the field close to the patient surface. This proximal collimation is necessary because electrons scatter far more in air than photons.

Early Work in Electron Beam Modelling

Berger and Seltzer, in 1978, were among the first to model the interaction of electrons with lead scattering foils — the LINAC component that most significantly influences the electron beam and where contaminant bremsstrahlung photon production occurs. An important finding: the intervening air causes significant energy degradation of the electron beam, whereas in high-energy photon beams the effect of air can be ignored. The same group that pioneered photon beam modelling was also one of the first to present a Monte Carlo model for clinical electrons.

Borrell-Carbonell et al. published simplified LINAC models in 1980, treating beam collimators as apertures with no wall interactions — an approximation that does not yield realistic particle fluence for electron beams. Rogers and Bielajew (1986) compared calculated and measured depth-dose curves for monoenergetic electrons, noting that simulations predicted a less steep dose gradient near the surface and an excessively steep dose fall-off beyond the depth of maximum dose. When the electrons were passed through the simulated exit window, scattering foils, and air, these differences diminished.

Andreo and Fransson (1989) demonstrated that stopping-power ratios are relatively insensitive to the details of the electron spectrum, but indicated that preserving the energy-angle correlation is critical — a direct implication for virtual source models. Ebert and collaborators studied simplified models of applicators and Cerrobend cutouts, identifying two main processes: electron scattering off the inner edges and bremsstrahlung production.

First Complete LINAC Models for Electron Beams

Modelling electron beams is widely recognized as more difficult than photon beams, owing to the greater sensitivity of particle fluences and absorbed-dose distributions to primary-beam details (energy, spatial, and angular distribution) and LINAC geometry — particularly the scattering foils and applicators/collimators.

The pioneering efforts to model complete electron beam geometries with the EGS4 code were made by Udale and Udale-Smith. Their models for Philips LINACs included the exit window, primary collimator, scattering foils, monitor chamber, mirror, movable photon jaws, accessory ring, and applicator.

Udale simulated five configurations of increasing complexity: from a monoenergetic pencil beam in vacuum to the full LINAC geometry. She used measured depth-dose distributions to tune the primary electron beam energy, attempting to match the 50% dose depth ($R_{50}$) and the practical range ($R_p$). Electron range rejection was employed as a variance-reduction technique to avoid transporting electrons that could not reach the LINAC exit — a practice still recommended despite the enormous increase in computer speed.

Udale recorded phase-space files at the bottom of the LINAC and, in a second step, used them for phantom dose calculations. She also extracted energy and angle distributions from the phase-space files for use in a virtual LINAC model. However, she demonstrated that some degree of correlation between position, energy, and angle must be maintained to avoid loss of simulation accuracy.

Udale-Smith compared models of several LINACs and established that some had superior designs: they produced fewer contaminant photons, lower-energy contaminant photons, fewer scattered electrons, and narrower electron angular distributions. Monte Carlo simulations are the ideal tool for such comparative equipment-design studies.

The BEAM Code and Electron Beam Applications

The advent of the BEAM code in 1995 was a landmark in Monte Carlo modelling of LINACs, including electron beams. In practice, the first results reported with BEAM were precisely for electron beams. The code provides a wide array of geometry modules, source geometries, variance-reduction techniques, scoring techniques, and tagging methods. It remains a popular and extremely useful tool.

The excellent agreement obtainable between measured and simulated dose distributions was demonstrated in the original 1995 paper by Rogers et al., using a 20 MeV beam from a research LINAC with very well-known characteristics. In clinical practice, the necessary LINAC details are often not fully available, forcing the user to “tune” the model.

In a series of studies, Ma et al. investigated multiple source models that exploit virtual point-source positions, taking into account the diffusivity of electron scattering. This makes source model derivation more complicated for electrons than for photons, which exhibit less scattering. Using realistic LINAC models, they calculated mean energy in water and stopping-power ratios.

Verhaegen et al. studied backscatter to the monitor chamber — an effect present in both electron and photon beams. The relative monitor chamber signal increase was 2% when jaws decreased from a 40 cm square field to 0 cm in a 6 MeV beam. For higher energies, the effect was smaller. A significant difference from photons: in electron beams, the spectral shape differs substantially between forward and backward electrons.

Another critical aspect is the size of the electron spot striking the scattering foils. Huang et al. (2005) proposed a slit-camera technique to derive the bremsstrahlung source size emerging from the foils, equivalent to the primary beam size. They obtained full width at half maximum (FWHM) values of 1.7–2.2 mm for 6–16 MeV beams and noted a primary beam displacement of up to 8 mm from the LINAC centre. Information about the primary electron beam remains among the hardest to estimate clinically.

Treatment Planning and Advanced Techniques

Conventional electron treatment-planning systems have well-documented errors in irregular fields and heterogeneous targets. It is widely accepted that Monte Carlo algorithms offer unmatched accuracy for electron dose calculations. Multiple treatment-planning systems with Monte Carlo modules are available and have been extensively evaluated.

Fast Monte Carlo codes with electron treatment planning as their primary application were introduced by Neuenschwander et al. (MMC, 1995) and Kawrakow et al. (VMC, 1996). Ma and collaborators studied multiple aspects: clinical beam simulation, beam characterization and modelling for dose calculations, air-gap factors for extended-distance treatments, commissioning procedures, and stopping-power ratios for dose conversion.

Modulated electron radiation therapy (MERT) was extensively studied by Ma and collaborators. In this technique, both the intensity and energy of the beams are modulated. Al-Yahya et al. designed a few-leaf electron collimator using Monte Carlo modelling during the design stage. The full range of rectangular fields deliverable by the device, combined with the available electron energies, served as input to an inverse-planning algorithm based on simulated annealing. The authors demonstrated that highly conformal treatments can be planned in this way.

A noteworthy recent development is FLASH radiotherapy. This technique exploits ultra-high dose rates to create a differential biological response between tumour and normal tissue. In this early phase, prototype LINACs are used, and Lansonneur et al. conducted one of the few studies using Monte Carlo simulations with the GATE code for FLASH electron radiotherapy. For the foundational particle transport concepts, see also our article on Monte Carlo fundamentals in radiotherapy.

Monte Carlo in Brachytherapy: Fundamentals and Evolution

In brachytherapy (BT), Monte Carlo simulation has become an essential tool, playing a central role in both clinical practice and research. The most established application is the determination of dose-rate distributions around individual sources. Modern sources contain low-energy radionuclides (mean energies < 50 keV) such as $^{103}$Pd, $^{125}$I, or $^{131}$Cs; higher-energy radionuclides such as $^{192}$Ir, $^{137}$Cs, or $^{60}$Co (mean energies of 355, 662, or 1250 keV); or miniature X-ray sources.

In low-dose-rate (LDR) brachytherapy, radioactive material and radio-opaque markers are encapsulated in permanently implantable seeds. In high-dose-rate (HDR) brachytherapy, an iridium pellet is encapsulated and welded to the tip of a remote afterloader cable. Although inverse-square-law dependence dominates dose distributions, photon attenuation and scatter buildup in the surrounding medium, combined with interactions within the source structure, produce anisotropic dose distributions.

The earliest computational efforts to obtain BT dose distributions date back to the 1960s, with Meisberger deriving tissue attenuation and buildup factors for $^{198}$Au, $^{192}$Ir, $^{137}$Cs, $^{226}$Ra, and $^{60}$Co point sources. Dale was the first to apply similar techniques to modern $^{125}$I sources in 1983. The first 3D Monte Carlo model of a BT source was performed in 1971 by Krishnaswamy for $^{252}$Cf needles.

A critical milestone: Williamson demonstrated in 1983, using 3D MC simulations, that the Sievert integral deviated by 5%–100% from MC results for monoenergetic photons with energies below 300 keV emitted by an encapsulated linear source. Burns and Raeside were the first to fully model a commercial $^{125}$I seed (model 6711), simulating the silver marker, the radioactivity distribution, and the titanium encapsulation to obtain a 2D dose distribution. Since the range of secondary electrons generated by 30 keV photons is less than 20 μm, they did not transport electrons and instead scored collision kerma using a track-length estimator.

Single-Source Dosimetry and the TG-43 Formalism

The rising popularity of LDR prostate seed implantation in the US — from 5,000 in 1995 to about 50,000 in 2002 — fuelled growth in commercially available seeds. The original TG-43 report (1995) presented consensus dosimetry parameters for one $^{103}$Pd and two $^{125}$I seeds. The TG-43U1 update (2004) covered 8 models; the 2007 supplement added 8 more; and the second supplement (2017) included the remaining commercially available low-energy sources.

The TG-43 formalism requires dosimetry data extracted from MC-computed or experimentally measured dose distributions. The parameters include the dose-rate constant $\Lambda$, the radial dose function $g_L(r)$, and the 2D anisotropy function $F(r,\theta)$. The AAPM requirement that at least one experimental and one MC determination of dosimetric parameters be published before clinical use of a source has made MC simulations a de facto standard in dosimetric practice.

An important point: MC results should not be accepted blindly. Significant differences in dose estimation result from the use of different photon cross-section databases. For dose distributions in the < 50 keV energy range, where energy deposition is dominated by photoelectric absorption, 1%–2% errors in the photoelectric cross section can generate 10%–15% dose errors at 5 cm from a seed. This led to the adoption of modern libraries derived from theoretical quantum-mechanical models.

The PSS (Primary and Scatter Separation) formalism, proposed by Russell and Ahnesjö in 1996, offers a complementary approach to TG-43. In this formalism, the primary dose distribution in water functions as the source’s “fingerprint,” independent of phantom size, allowing dose distributions in heterogeneous media to be derived from the water data.

Cross Sections and Dose Estimators

The choice of cross-section libraries is critical for MC simulations in brachytherapy. Modern libraries — EPDL97 (Lawrence Livermore), DLC-146 (RSICC), and XCOM (NIST) — are based on the same theoretical models, despite differences in format and extent of the compiled data. For BT, the most important consideration is selecting a library with accurate post-1983 photoelectric and scattering data.

In brachytherapy simulations, omitting electron transport and simulating only photon transport substantially reduces the computational burden. This simplification — approximating absorbed dose by collision kerma — is valid everywhere for low-energy sources, where electron ranges are < 0.1 mm. However, for high-energy sources ($^{192}$Ir, $^{137}$Cs, $^{60}$Co), with secondary electron ranges of 1–5 mm, the charged-particle equilibrium approximation may introduce significant errors near metal-tissue interfaces. Dose errors exceeding 15% at distances less than 1 mm from an HDR $^{192}$Ir source have been reported.

For scoring, the analogue estimator — where only collisions within the voxel of interest contribute to the dose — is the simplest approach. A far more efficient alternative is the track-length estimator, which approximates dose as collision kerma scored by $(\Delta l \cdot \Delta E) / (V) \cdot (\mu_{en}/\rho)$. Efficiency improves substantially because every voxel intercepted by a photon’s path produces a non-zero score. Relative gain factors over the analogue estimator reach 20–50 for $^{125}$I scenarios, and up to 70 ($^{103}$Pd), 90 ($^{125}$I), and 300 ($^{192}$Ir) for different BT treatments.

For a deeper understanding of photon beam modelling, which underpins many of the concepts discussed here, see our article on Monte Carlo modelling of external photon beams.

Outlook and Final Remarks

Monte Carlo simulations are powerful and increasingly essential tools in electron beam radiotherapy and brachytherapy. For electrons, the superiority over analytical algorithms is even more pronounced than for photons, given the well-known deficiencies of conventional methods. MC simulations play a decisive role in designing complex dose-delivery techniques such as modulated electron radiation therapy (MERT) and, more recently, FLASH radiotherapy.

In brachytherapy, the Monte Carlo method has been consolidated as the standard for individual source dosimetry, underpinning the TG-43 formalism that forms the basis of clinical practice. MC methodology, extensively benchmarked against experimental measurements with 1%–3% accuracy, is sufficiently mature and robust to support clinical dosimetry across the entire BT energy spectrum.

EGS/EGSnrc-based codes (including the BEAM user interface) have historically dominated this field, but other codes such as GEANT4 are seeing increasing use. An important caveat: since MC codes generally differ more in cross sections and electron transport methods than in photon transport, rigorous benchmarking of simulations is essential for electron beams.