4D treatment planning for scanned ion beams
© Bert and Rietzel. 2007
Received: 05 April 2007
Accepted: 03 July 2007
Published: 03 July 2007
Skip to main content
© Bert and Rietzel. 2007
Received: 05 April 2007
Accepted: 03 July 2007
Published: 03 July 2007
At Gesellschaft für Schwerionenforschung (GSI) more than 330 patients have been treated with scanned carbon ion beams in a pilot project. To date, only stationary tumors have been treated. In the presence of motion, scanned ion beam therapy is not yet possible because of interplay effects between scanned beam and target motion which can cause severe mis-dosage. We have started a project to treat tumors that are subject to respiratory motion. A prototype beam application system for target tracking with the scanned pencil beam has been developed and commissioned.
To facilitate treatment planning for tumors that are subject to organ motion, we have extended our standard treatment planning system TRiP to full 4D functionality. The 4D version of TRiP allows to calculate dose distributions in the presence of motion. Furthermore, for motion mitigation techniques tracking, gating, rescanning, and internal margins optimization of treatment parameters has been implemented. 4D calculations are based on 4D computed tomography data, deformable registration maps, organ motion traces, and beam scanning parameters.
We describe the methods of our 4D treatment planning approach and demonstrate functionality of the system for phantom as well as patient data.
Intrafractional target motion has a relevant impact on the precision of treatment delivery in conformal radiotherapy. Even if treatment planning margins are sufficient to encompass the full extent of motion, intrafractional motion degrades dose gradients to surrounding healthy tissue [1, 2]. For intensity modulated therapy like IMRT or scanned particle therapy, the relative motion between target and multi-leaf collimator or scanned beam can have a severe impact on the delivered dose. These interplay effects between target and beam motion usually cause hot and cold spots in the delivered dose distributions [1, 3–10].
To mitigate the impact of intrafractional motion, several techniques were proposed but have not yet been used clinically with scanned particle beams: beam gating [11, 12], rescanning [13, 14], tracking [15–17], and internal margins (IM) to generate an internal target volume (ITV) . Besides tracking, all other techniques require IMs due to unmitigated or residual motion. For particle beams, the use of ITVs requires explicit consideration of possible range changes because ranges are most often influenced by organ motion [19, 20]. Adjustments of the beam range, e.g. by compensator smearing, have to be applied to ensure coverage of the distal field edge for all motion states of the target [19, 21].
Precise motion mitigation techniques require quantification of the motion, for example by time resolved computed tomography (4DCT) [22–26] or 4D magnetic resonance tomography (MRT) . Both techniques sample periodical motion in several motion phases. The motion phases correspond to quasi-static 3D volumes, e.g. standard CT volumes. Usually sampling is based on an external motion surrogate [28–30]).
Several investigators have used time-resolved volumetric imaging to extend treatment planning capabilities for tumor sites influenced by respiratory motion [21, 31–34]. The main idea was reported by Keall et al.  and Rietzel et al. : Dose calculations are performed per motion phase of 4DCT data. Deformation maps obtained by non-rigidly registering motion phases are then used to transform the resulting sub-dose distributions to a reference motion phase for effective dose calculation by time weighted summation.
At Gesellschaft für Schwerionenforschung (GSI) charged particle therapy is performed with an intensity-modulated, raster-scanned carbon ion beam in collaboration with the University Hospital Heidelberg, the German Cancer Research Center, and the Forschungszentrum Rossendorf. To date, more than 330 patients with tumors in the head, neck, spinal cord, and pelvis region have been stereotactically treated within the pilot project [35, 36]. GSI plans to treat tumors affected by respiratory motion by tracking. Consequently we have started a project to develop beam application and treatment planning capabilities for this technique [15, 37, 38]. For a short introduction, the following paragraphs summarize standard treatment delivery and treatment planning for scanned carbon ion beams as well as the current status of the tracking project.
Treatment delivery is performed with the active raster-scanner system . In beam's eye view, the planning target volume (PTV) is divided into j slices of iso-energies E j (typical slice distance: 3 mm water-equivalent). Each iso-energy slice contains a regular grid (typical grid spacing: 2–3 mm) of i beam positions (x i , y i ). Individual pencil beams with a focus of 3–9 mm (full width at half maximum) are applied per grid position. To achieve the desired dose distribution, the number of carbon ions N ij is optimized for each position. During treatments, a synchrotron accelerates carbon ion pencil beams to the required beam energy E j . Particle extraction from the synchrotron is performed in beam pulses with a length of 2.2 s followed by 3.3 s of acceleration for the next energy. For each iso-energy slice (IES) a new pulse has to be requested from the synchrotron because beam energy can – up to now – not be changed within a pulse. For each IES the pencil beam is scanned by a magnetic deflection system across all beam positions with N ij > 0. The scanning process is controlled by fluence monitors. They measure the number of particles deposited per beam position and request transition to the next position as soon as the required number of particles N ij has been reached. The beam is not turned off during transition to the next grid position. Scanning speed (up to 11m/s) is thus dependent on N ij and the synchrotron extraction profile.
d (E j , z) quantifies the energy loss distribution for a certain beam energy E j with respect to traversed amount of tissue z in water-equivalent units. During optimization, the energy levels E j and the raster grid (x i , y i ) are set to cover the extent of the PTV in beams-eye-view. The minimal beam-width is determined by the grid-spacing (uniform in x and y) according to σ > 1.27 (xi + 1- x i ). To calculate the longitudinal extent (zmin - zmax) of the target, CT numbers are converted into particle ranges based on a Hounsfield look-up table . Optimization of N ij is performed by least square minimization such that D (x, y, z) meets the prescribed dose.
The combination of scanned particle beams and target motion represents a double-dynamic system that requires a dedicated solution for 4D treatment planning. We extended our treatment planning system TRiP  to full 4D functionality to allow dose calculation and parameter optimization in the presence of motion. In principle the 4D functionality can handle all types of motion. In the following we will, however, focus on respiratory motion as this is our initial intended application.
In order to compare different motion mitigation techniques like tracking, gating, rescanning, and internal margins, the 4D version of TRiP allows to
generate particle specific ITVs for treatment plan optimization of gating and rescanning,
optimize parameters for motion compensation (tracking),
calculate physical dose distributions in the presence of motion as well as for tracking, gating, rescannning, and the use of internal margins.
In contrast to previous simulation studies of our group [15, 17], we have implemented 4D treatment planning based on multiple volumetric data sets, e.g. CT data sets. In addition, ITV generation and optimization of tracking parameters including methods to correct for target rotation and deformation were realized. Calculation of physical dose distributions can be performed for patient data with patient specific, non-rigid motion.
The purpose of this contribution is the technical description of the 4D treatment planning extensions. The functionality will be presented for phantom simulations as well as for an example patient data set. Experimental validation and application to clinical data of lung cancer patients will be reported elsewhere.
For scanned particle beams, 4D dose calculations require temporal correlation of beam motion and organ motion considering possible changes in particle range. Dose calculations are based on a reference motion phase independent of the motion mitigation technique (see fig. 1). The following sections describe the calculation of dose distributions in the presence of motion as well as the parameters that are required for these calculations.
For treatment planning we assume organ motion to be non-rigid as well as represented by time-resolved volumetric data that allows precise calculation of particle ranges, e.g. CT data. For respiratory motion measured by 4DCT, motion is assumed to be periodical. Temporally, organ motion has to be measured in correlation to beam motion (section Beam motion parameters). Spatially, organ motion is described with respect to the particle beam and with respect to a reference motion phase (fig. 1).
3D information of the motion like amplitude, trajectory, and volumetric changes are included in the 4DCT phases. Quantitatively, the motion is parametrized by B-splines which describe the non-rigid motion components between 4DCT phases. Optimization and application of these transformation maps are performed with vtkCISG . Details on calculation and validation of the transformation data have been reported previously by Rietzel et al. [33, 45].
Temporal changes from 3D data set to 3D data set were implemented via motion traces. For example 4DCT period and initial respiratory phase can be given by a motion trajectory. The 4D version of TRiP allows handling of measured motion trajectories, modeling of sinusoidal motion, or modeling motion according to Lujan et al. .
With beam motion we refer to the time dependent movement of the raster-scanned pencil beam as it traverses the target volume grid-position by grid-position and slice-by-slice (see fig. 2 and section Carbon ion therapy at GSI). Pencil beam motion is quantitatively determined by the particle intensity profile extracted from the synchrotron, the number of particles per grid-position N ij , the order of beam positions (x i , y i ) within each iso-energy slice j, and the order in which iso-energy slices are irradiated.
Particle extraction is not exactly deterministic or reproducible (see fig. 2 in ). Slight changes in the acceleration and especially the extraction process lead to changes that do not affect irradiations of stationary targets but cause changes in the temporal progress of the scanning process. For precise dose calculations in the presence of organ motion, beam intensity and irradiation time of each individual beam position (typical duration < 10 ms per position) have to be considered in temporal correlation to organ motion.
The 4D version of TRiP can handle measured intensity distributions and can model the extraction characteristics at GSI as well as the characteristics of the Heidelberg Ion-Therapy center (HIT, under construction) . In contrast to GSI's synchrotron with so called slow extraction, knock-out-extraction  will be used at HIT. This extraction method allows intermitted extraction within one pulse and thus optimal gated irradiations. Modeling of the extraction pattern for gated irradiations has to be performed for each motion trajectory and gating window combination individually because the pulse structure (~1.5s acceleration for each iso-energy slice followed by a maximal pulse length of 10 s) is fixed. Treatment planning can thus also be used to estimate realistic treatment times for gated irradiations.
For all motion mitigation techniques, treatment delivery is based on a reference motion phase. For respiratory motion, this reference motion phase typically corresponds to the end-exhale phase of the 4DCT data. Using the standard functionality of TRiP , a reference treatment plan (x, y, E, N) is optimized on the reference motion phase. This reference treatment plan is then applied to the moving target by raster-scanning. The reference treatment plan is modified by compensation parameters Δ(x, y, z, N) at time of delivery for tracking, interrupted for delivery by gating, applied multiple times for rescanning, and unchanged for mitigation by internal margins.
In contrast to other motion mitigation techniques, tracking changes beam parameters. Corresponding sub-treatment plans require consideration of compensation parameters. Adjustment of lateral beam positions (Δ x, Δ y) and numbers of particles Δ N can readily be applied (e.g. x new = x + Δ x). The longitudinal change Δ z corresponds to a shift of the depth dose distribution (d (E, z + Δ z)) and has to be considered in the summation of dose contributions (eq.1).
Sub-dose distributions represent the cumulative dose delivered within a specific motion phase. For evaluation of a treatment scenario (e.g. mitigation technique, target motion parameters, extraction rate), the effective dose distribution of the complete treatment plan has to be calculated. Because the quasi-static motion phases are not anatomically registered the sub-dose distributions can not simply be summed but need to be transformed to the reference motion phase. Transformation maps which quantitatively describe the non-rigid motion (see section Organ motion parameters) are used to transform each sub-dose distribution to the reference motion phase. The transformed sub-dose distributions can then be summed time weighted to calculate the effective dose distribution.
Changes in the optimization of treatment plans for moving targets in comparison to the algorithms used for standard irradiations (see section Carbon ion therapy at GSI and ) depend on the motion mitigation technique. The following sections describe our implementations of rescanning, gating, internal margins, and tracking.
In contrast to tracking, gating and rescanning do not require change of treatment parameters during delivery. However, particle specific ITVs have to be generated for both techniques. Particle specific ITVs account for possible range changes due to organ motion . ITVs ensure target dose coverage in each 4DCT phase for rescanning or for the subset of 4DCT phases included in the gating window.
Timing of the beam delivery sequence is modified for gating, where the beam is turned on only within the gating window, e.g. close to the end-exhale breathing phase. The gating window has to be defined for the optimization. Rescanning also changes timing of the beam delivery sequence because the same irradiation pattern is applied several times. The number of particles per grid position is therefore divided by the number of rescans. Typically the beam fluence is constant but the scanning speed increases. Besides ITV generation, the definition of the number of rescans is the main task for optimization of rescanning.
Required input: native 4DCT and ITV contour
Calculation of the maximal lateral ITV extension in beam's eye view (BEV) and setup of the raster grid (x i , y i )
Definition of iso-energy slices (IES with accelerator energies E j ): The IES are arranged from the distal to the proximal water-equivalent extension of the ITV. Because tumor motion influences particle ranges the extreme water equivalent ranges of all 4DCT phases have to be considered per grid position to ensure target dose coverage to the distal edge independent of the motion phase during irradiation.
Calculation of the water-equivalent depth for each CT voxel as required for dose calculation (d (E, z) in eq.1): Because the depth is influenced by organ motion the maximum depth from all 4DCT phases has to be used per voxel.
v. Optimization of N ij at each grid position (x i , y i ) based on the maximum water-equivalent depth data to achieve the required dose distribution.
The additional treatment parameters for gating and rescanning (gating window and number of rescans) can be determined or even optimized by performing corresponding dose calculations. Simulated organ motion trajectories and particle extraction rate are necessary for this step. By variation of the parameters for organ motion and beam application coverage of the CTV can be analyzed for different gating and rescanning parameters.
For tracking, target motion is mitigated by adaptation of the reference treatment plan parameters during treatment delivery. Calculation of the compensation parameters Δ(x, y, z) is based on a 4DCT data set and corresponding transformation maps (cf. section Organ motion parameters). Compensation parameters have to be calculated during treatment planning because calculations are too time-consuming to be performed online during treatment delivery. Because motion trajectory and temporal scanning progress (determined by the synchrotron extraction) are not known at the time of treatment planning, compensation parameters have to be calculated for all possible interplay combinations. During treatment delivery, the motion phase is continuously measured, ideally with the same system as used for 4DCT acquisition. Based on the currently irradiated grid position the corresponding pre-calculated compensation parameter set Δ(x, y, z) is used for beam adaptation. Fluctuations in synchrotron extraction have no impact on compensation parameter sets because the intensity controlled raster-scanning process determines which grid-position is irradiated. The motion detection system continuously determines the actual motion phase or interrupts the irradiation if the current motion state is not included in the pre-calculated parameter sets.
Mitigation of dose changes Δ D due to dose contributions that differ from those in the reference treatment plan can be achieved by adaptation of the number of particles Δ N at each 3D grid position. Because the target is scanned only once in a pre-determined manner, Δ D depends on the irradiation order of iso-energy slices and the scan path within each slice. Dose contribution mitigation is only possible for dose changes resulting from grid positions irradiated previously. The irradiations should therefore start with the highest beam energy at the most distal slice because dose contributions to different grid positions are significantly higher in the entrance channel than in the fragment tail.
Determination of the CT coordinate of each grid position in the reference treatment plan by conversion from the water-equivalent system to the CT system based on the reference 4DCT phase (fig. 1, grey volume)
Motion vector determination (fig. 1): the transformation maps (section Organ motion parameters) provide the geometrical transformation into all other 4DCT phases.
Lateral compensation (Δ x, Δ y): corresponds to the motion vector components
Longitudinal compensation component Δ z : corresponds to the change in particle range between the original grid position in the reference 4DCT phase and the transformed grid position (fig. 1, white circle in red volume) in the corresponding 4DCT phase.
The calculation in (i) is repeated for all possible motion phases. This dose calculation has to consider the changed 4DCT phase as well as the compensation based on the Bragg peak position, Δ(x, y, z) (fig. 7b).
The change in dose Δ D due to different contributions is the difference between the dose contributions of (i) and (ii).
During treatment, changes in deposited dose Δ D are used to determine the required change in particle deposition Δ N. Each grid position causes dose changes Δ D at several other grid-positions irradiated afterwards which are taken from the pre-calculated data depending on the motion phase valid at the time of delivery. Consequently each grid position suffers from dose changes from grid-positions irradiated previously. For determination of Δ N, the cumulative dose changes from all previously irradiated grid positions is used (∑Δ D). The adjustment Δ N is given by Δ N = - ND -1∑Δ D where N and D are parameters of the reference treatment plan. If N + Δ N is negative, i.e. if too much dose has been applied, no particles are delivered at the grid position but the over-dosage can not be corrected for. This shows the lack of optimization for the described approach, but this is unavoidable if arbitrary motion trajectories of beam and target can occur. The main goal is avoidance of under-dosage.
We have implemented 4D treatment planning for charged particle radiotherapy with scanned pencil beams. Dose calculations in the presence of motion as well as optimizations for gating, rescanning, and tracking are based on time-resolved anatomical data, motion trajectories, and extraction characteristics of the accelerator.
Theoretically, tracking of the target with the scanned pencil beam should result in the best possible sparing of surrounding, healthy tissues. Whereas for gating and rescanning, required ITV margins lead to an increased PTV that encompasses surrounding tissues. For gating, the size of the gating window and the resulting residual target motion will determine ITV expansions. Then a trade off between treatment delivery time and target conformity has to be made. Comparing the motion mitigation techniques, it should be noted that increasing conformity will elevate technical complexity. As long as tracking has not been developed for clinical routine use, we clearly favor gating in comparison to rescanning due to the increase in conformity. Currently, treatment planning studies based on 4DCT patient data are performed to explore and quantitatively evaluate the differences between motion mitigation techniques.
In general, treatment planning – independent whether in 3D or 4D – relies on the validity of the input data. In current practice, these data are most often assumed to be ground truth throughout the treatment course. Generation of the required input data is out of the scope of this paper. We will therefore only briefly discuss possible implications on 4D treatment planning.
4DCT samples moving anatomy in several discrete 3D motion phases. Usually motion is detected during data acquisition by an external monitoring system to sort CT data according to respiratory phases. As pointed out by several authors, 4DCT is not free of residual motion artifacts for example due to irregular respiration [22, 23, 25, 26, 49]. Such artifacts will have an impact on 4D planning, manifested by wrong position of the target and possible changes in particle ranges. Possible improvements in 4DCT data acquisition are currently under investigation [50–52].
Target motion trajectories could be measured directly by fluoroscopy [53, 54] or indirectly by external motion detection systems [28–30]. Techniques and limitations of motion detection systems are out of the scope of this contribution; we assume a reliable detection of motion phases. For retrospective dose calculations, motion trajectories have to be recorded during irradiations only. For optimization of motion mitigation strategies, the treatment delivery system has to react to actual motion phases online. Fluoroscopic tracking has been successfully used in Japan . Treatments are gated with millimeter precision based on trajectories of fiducial markers close to the target. Currently, fluoroscopic tracking of the target or nearby structures without fiducial markers is under investigation. If external motion detection systems have to be used, ideally, the same motion detection system would be used during treatment delivery as was used for 4DCT data acquisition.
Besides target motion beam motion during scanned beam application has to be considered to model interplay effects. Recording of the irradiation time of each beam position has already been implemented at GSI. Treatment times for individual pencil beam positions are typically below 10 ms which usually results in less than 0.1 mm of motion for typical respiratory parameters.
We extended GSI's treatment planning system TRiP to full 4D functionality. The new modules facilitate 4D dose calculation and optimization for tracking, gating, rescanning, and internal margins. Calculations and optimizations are based on 4DCT information, organ motion, and trajectory of the scanned ion pencil beam.
The authors thank Prof. Dr. Gerhard Kraft for fruitful discussions and supervision of this project and Dr. Michael Krämer for constant support regarding TRiP. The authors thank Siemens Medical Solutions, Particle Therapy, for partial funding of this project.
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.