Conformal Dose Modulator for Proton Beam Therapy Part 1: A Simulation Study.

Purpose: Conformal dose deliveries in proton therapy utilize either a passive scattering system with a modulator or a pencil beam scanning (PBS) system. Efforts have been made to achieve conformal dose delivery by scanning a single energy layer of pencil beams through a 3D conformal modulator (3DCM), which combines a spread-out Bragg peak (SOBP) modulator consisting of a micro-pyramid array and a range compensator. The current published approach of designing such 3DCM relies on forward calculation methods to determine the geometry of the modulator. This study presents an alternative designing algorithm that inversely generates the geometry of a 3DCM paired with a corresponding fluence map, customized to patient-specific clinical indications.
Methods: Critical spacing governing the size and separation between neighboring micro-pyramids was first determined, under which the dose homogeneity at desired depths could be achieved. We designed an adaptive ring optimization method using a modified gradient descent algorithm to inversely calculate the geometry of the 3DCM. This method includes several stages that progressively optimize both target coverage and dose conformity. The output contains the geometry of the 3DCM and its corresponding proton fluence map. Monte Carlo (MC) simulation was used to validate the results.
Results: The critical size and spacing of Lucite pyramids was determined to be 0.5 cm for a 184-MeV pristine proton beam. Using MATLAB (R2020a), the inverse designing algorithm generated an optimized 3DCM geometry and a fluence distribution achieving 100% target coverage with the 90% isodose surface and a corresponding conformity index of 1.057 on a spherical target. The resulting geometry was pruned to accommodate the MC simulation software and a currently accessible 3D printing service. The pruned geometry gave 95% target coverage by 90% isodose surface with a conformity index of 1.09 by ray-tracing dose computation. The MC simulation validated the 3DCM with 95% target coverage by 87% isodose surface and a conformity index of 1.12.
Conclusion: We have demonstrated the feasibility of using a novel inverse optimization algorithm to generate 3DCM geometry and its corresponding proton beam fluence/intensity map, which could deliver highly conformal dose distribution with pencil beam scanning system using a single energy layer.

1 Introduction

In charged particle beam therapy, a treatment target volume is covered by transversely and longitudinally spreading out Bragg peaks (SOBP). In passive scattering broad-beam systems, this is achieved by combining a SOBP-modulator and a range compensator with beam collimating features. The SOBP produced by the modulator using either a modulating wheel or ridge filter has fixed width while the range-modulating compensator conforms the shape of the prescribed dose only to the distal contour. This results in unnecessary exposure to the normal tissues in the proximal region of the target. In state-of-the-art charged particle systems, pencil beam scanning techniques (PBS) have become the primary modality for treatment delivery. SOBP is achieved by scanning pencil beams with a specific spot size through a consecutive number of iso-energy layers or slices, which can effectively conform the dose to the target contour in all three dimensions.

Despite the trend of PBS becoming more prominent than scattering systems, there are some advantages of passive systems that are still clinically relevant; examples include shorter treatment time and the robustness of dose distribution for targets under the effects of respiratory motion. The efforts of combining the advantages of passive technique and PBS system have been attempted to create 3D conformal dose distributions (1, 2). Sakae et al. (2–4) and Ishizaki et al. (5) first reported the approach of using 3D modulators consisting of stacked arrays of pre-fabricated cone-shape mini ridge filters and a compensator. Yuri Simeonov et al. (1) reported the same concept but using a 3D-printed customized array of thin pins in the construction of ridge filter. These earlier efforts demonstrated the feasibility of combining PBS with a physical modulator to achieve conformal dose coverage using a single iso-energy scanning layer, which can potentially save treatment time and be more effective in motion management during treatment.

In this study, we furthered the methodology by developing an inverse algorithm and procedures that would allow robust fabrication of 3D Conformal Modulators (3DCM) customizable to general clinical indications.

The main approaches include two stages:

Stage 1: The use of mini pencil beamlets of very small spot size to inversely generate optimized geometry of 3DCM with corresponding 2D proton fluence map with a single energy that can be used to treat tumors of arbitral depth and shape.

Stage 2: The generated geometry of the 3DCM can be 3D printed. The quasi-continuous 2D fluence map will be subsequently converted to a clinically deliverable pencil beam spot sequence (positions and fluences) of a single energy layer.

The objective of this study is to carry out a proof-of-concept simulation on the first stage, i.e., the inverse design/production of the 3DCM and the 2D fluence map. A Monte Carlo (MC) simulation was used to validate the method. The 3D printing of 3DCM and the conversion of 2D fluence map to deliverable spot sequence, as well as the subsequent dosimetric verification by phantom measurement are currently under investigation and will be reported in a separate publication.

We currently limit the application to proton beam therapy, although the methodology can be extended to other charged-particle beams with appropriately integrated biological models.

2 Methods and Materials

The physical 3DCM has two basic components: a non-uniform SOBP modulator customized to the various thickness of the target volume along the beam axis, and a range-modulating compensator customized to the target contour at the distal end ( ).

Figure 1

(A) Uniform SOBP modulator—poor overall conformity; (B) adding range-modulating compensator (shifter)—improved conformity, but only in the distal target margin; (C) 3DCM: combining non-uniform SOBP modulator and range-modulating compensator—good conformity in both distal and proximal target margins.

To demonstrate the principle and feasibility, a theoretical mono-energetic proton beam with 184 MeV (no energy spread, Bragg peak at a water equivalent depth of 22.5 cm) was used with Lucite being the material for the modulator, and water media for dose computation.

To generate the optimal geometry of 3DCM and the proton fluence map, mini pencil beamlets with a 0.01 cm radius, rather than actual clinical pencil beams, are used. The dose kernel of the mini pencil beamlet was generated by the MC simulation program (FLUKA version 4-0.0) (6–9) using the PRECISIO mode, which defaults particle transport threshold at 100 keV. The resultant dose kernel was stored in a matrix with 1 × 1 × 1 mm3 resolution.

2.1 Non-Uniform SOBP Modulator

The non-uniform SOBP modulator is an assembly of pyramid-like ridge filters. Each pyramid is formed by stacking cuboids with decreasing side length. A single-energy broad beam (to be formed by array of pencil beamlets), with the beam axis perpendicular to the pyramid base, would penetrate through different thicknesses of the material. The dose distribution along the beam axis is the integration of all sub-beamlets passing through various heights of the pyramids, forming a SOBP.

The geometry of the pyramid is derived through multiple iterations. Each pyramid with its unique shape is expected to produce a homogeneous SOBP region that would not only cover the thickness of the target along the beam axis but also merge smoothly with the neighboring SOBPs produced by adjacent pyramids. The critical elements in the geometry of the pyramid array include the shape of the pyramids (base, height, and lateral graduation) and the spacing between them.

2.1.1 Critical Size of the Pyramids and Their Spacing

To obtain the critical size and spacing of the pyramids, the mini pencil beamlets were used in hexagon formation, in which the mini pencil beamlets were located at the center and at the vertices of each hexagon, where the central axes of the pyramids are expected to situate. The critical spacing is defined and obtained by gradually changing the size of the hexagon until the composite lateral dose distribution at the depth of the Bragg peak achieves desired uniformity, with individual pencil beamlets dose pattern vanished. The obtained critical spacing was then used as a guide to determine the distance between the centers of each pyramid, which would also be the maximum side length of the pyramid base. Note that no pyramid was yet used in this calculation. The purpose of this procedure is to assure all beamlets going through various parts of a pyramid and through the neighboring pyramids will merge smoothly, when the side length of the pyramid base is equal to or less than the critical spacing.

2.1.2 Geometry of the SOBP Pyramid

The shape of a pyramid is formed by stacking cuboids of various thicknesses (levels), which leads to a spread out energy distribution of outgoing beams. For each level’s net area, a circumference of beam cross-section corresponds to the weighting factor of a specific outgoing energy associated with the thickness of that level. The design criterion of pyramid levels is to produce a uniform plateau of SOBP. Accordingly, we generated another set of 40 dose kernels using the FLUKA software, each corresponding to a 184-MeV mini pencil beamlet (0.01 cm radius) passing through Lucite of thickness 0.2 cm to 8 cm with a 0.2-cm increment. These dose kernels would be integrated with certain derived weighting assignments to obtain the final depth dose curve in the direction of the beam axis. The weightings of all levels (cuboids) of a pyramid were inversely calculated using a modified version of the basic gradient descent algorithm (10):

where θ is the array of weighting factors corresponding to each beamlet, which is forced to be greater than or equal to zero to ensure physical validity; α is a user-defined step size to ensure fast convergence to a local minimum of the cost function; x is an array that contains all the depth dose curves; y is the optimization goal, a set of the anchor points of an ideal SOBP curve. In this study, the anchor points were assigned to be 1 for the depths from 15 cm to 20.2 cm and 0 for the depths from 21.5 cm onward, and the initial value of θ was assigned to be zero. These assignments ensure that the calculated weighting factors would yield a smooth plateau between depths of 15 cm and 20 cm with a sharp distal dose fall-off. The setting would keep the entrance dose as low as possible. J is the cost function that tracks the difference between optimization results and the preset goal.

2.1.3 SOBP Modulator and Its Design Validation

To demonstrate the validity of the described method, we constructed a uniform SOBP modulator with a cross-sectional size of 5 × 5 cm2, formed by arranging individual pyramids abutting each other side by side with each row offset by half of a pyramid. The side length of the pyramid base was set to the critical spacing described in Section 2.1.1. We then conducted MC simulation of the dose distribution in water produced by a 184-MeV broad proton beam with a field size of 4 × 4 cm2 penetrating through the modulator. The MC simulation was performed in the PRECISIO mode, with a resolution of 1 × 1 × 1 mm3. The resultant percentage depth dose along the central axis was compared with the intended shape of SOBP to validate the inverse weighting factors calculated by the modified version of the basic gradient descent algorithm. The lateral dose homogeneity of the SOBP plateau region was analyzed to validate the determination of the critical spacing and the side length of the pyramid base.

2.2 Generation of 3DCM With Optimized 3D Conformal Dose Distribution

With the inverse algorithm for the optimal structure of SOBP modulator (pyramid-based) at hand, we proceed to combine energy and range modulations to produce a 3D conformal modulator using a 184-MeV mono-energy proton beam.

The test target volume was a sphere of 4 cm in diameter. We first designed a virtual “mesh” perpendicular to the beam axis. The size of each cell of the mesh is 2 mm × 2 mm (less than the aforementioned critical pyramid spacing), and would serve as a placeholder for 2 mm × 2 mm pyramids. The rows of the mesh offset with each other by 1 mm. We created a new sphere by expanding the target sphere radius by one cell width. The entire mesh area was defined by the maximum cross-section of the expanded sphere.

Using the inverse optimization method described in Equation (1), the set of 40 beamlet dose kernels passing through Lucite of thickness from 0.2 cm to 8 cm with a 0.2-cm increment were placed to the center of each cell. The dose kernels with their Bragg peaks located outside the expanded target sphere were removed. The remaining dose kernels constitute the parameter x. The array of the weighting factors θ, corresponding to each of the 40 dose kernels, would be used to define the graduating shape of each individual pyramid. The optimization goal y would be the anchor points representing the desired dose at specific locations, and would be adaptively adjusted based on the optimization staging. The dose calculation uses a 2 mm × 2 mm × 2 mm grid spacing.

We developed an adaptive ring optimization technique that consists of three stages. The first stage focuses on the target coverage. Target coverage in this stage was defined as the percentage of the target volume covered by 90% isodose surface. The initial inputs for θ would be zero. The first optimization goal y would be set to one within the spherical target volume, without any constraint outside the target volume. During the optimization iteration, the weighting factor array θ resulting in the highest target coverage would be cached, and the iteration would stop once the target coverage converged.

The second stage aims to improve the conformity index by introducing adaptive ring constraints:

where V denotes the target volume. V denotes the treated volume, which is the volume enveloped by the prescribed isodose surface (90% in this study). The conformity index nCI is the product of two ratios, one being between the target volume and the intersectional volume of target and treated volume, and the other ratio being between the treated volume and the intersectional volume of target and treated volume. The conformity index defined in Equation (2) takes into consideration the intersectional volume of target and treated volume. When the target volume and the treated volume are completely overlapped, nCI = 1 is the optimal index value. In order to regulate the optimization, we created three sets of “rings” perpendicular to the beam axis at various depths surrounding the target sphere. Each set of rings contains ten rings with different diameters at ten depths. The first set of rings has diameter 10 mm larger than the target sphere’s cross-sectional diameter at the corresponding depths; the second and third set of rings have diameters 16 and 24 mm larger than the target cross-section, respectively. In addition to the original optimization goal of setting 100% target volume coverage, constraints on rings would be implemented into y. In the first iteration, the weighting factor θ computed from stage one was used to calculate the dose distribution. Both nCI and target coverage were assessed. The relative dose on each ring was evaluated. The dose goal for the first set of optimization rings was set to be 10% less than the prescribed dose, the second set to be 6% less, and the third set to be 3% less. The gradient descent algorithm would then be used to optimize θ, the array of the weighting factors. During the first iteration, the θ corresponding to the highest nCI with the coverage greater than 95% was saved, and the iteration would stop once the cost function stabilized. After the first iteration, the optimization goals of all rings would then be adaptively adjusted based on the difference between the current dose distribution from the updated θ and that of the previous iteration. The program would subsequently proceed to optimize θ using gradient descent algorithm and would cache the optimal θ that yield the highest nCI while keeping the target coverage greater than 95%. The adaptive ring optimization would continue until the target coverage dropped below 95%. The cached θ resulting in the highest nCI would be the final optimal weighting factors.

The optimal θ values contain the weightings of all dose kernels that would result in the desired width of SOBP covering the target along that beam axis. The weightings of those dose kernels passing through the same cell on the virtual mesh would be used to reconstruct the pyramid located on that cell. Essentially, the shape (height and graduation) of an individual pyramid was determined by the values of θ on its corresponding cell. The sum of the weightings belonging to each cell was assigned to a new array, which represents proton fluence weight going through the pyramid located on that cell.

presents a flow chart of the aforementioned optimization process.

Figure 2

Flow chart of generating 3DCM geometry and fluence map through the inverse optimization process.

The 3DCM geometry and the proton fluence map would be imported into FLUKA to compute the dose distribution in water using MC simulation (in the PRECISIO mode, with resolution of 1 × 1 × 1 mm3). The nCI and the target coverage would then be compared with those calculated by the inverse optimization process.

3 Results

3.1 SOBP Modulator

3.1.1 Critical Size of the Pyramids and Their Spacing

The lateral dose distribution at the depth of Bragg peak produced by a single 184-MeV mini pencil beamlet with a 0.01-cm radius in water obtained by MC simulation is shown in . The radiuses of multiple isodose lines (5%, 10%–90%, separating in 10%, and 95%) were measured (denoted as r , r , etc.). These radiuses will be used as the side lengths of the hexagons in the subsequent operation.

Figure 3

Lateral dose distribution of pencil beams in hexagon formation (units on axes are in cm). The individual dose pattern vanishes as the distance between neighboring pencil beams reduces to less than 0.5 cm. (A) Dose pattern resulted from one pencil beam; (B–L) dose pattern resulting from seven pencil beams in hexagon formation; the distance between them reduced from r 5% to r 95%.

Dose distributions at the Bragg peak depth of seven pencil beamlets located at the center and the vertices of a hexagon with various side lengths d set according to the values of r 5%, r 10%, etc. were computed and plotted ( ). The critical spacing was identified when d = r 40% (measured to be 0.5 cm, ), as the superimposed dose distribution resembles a dilated dose distribution from a single pencil beam, with individual pencil beamlet’s dose pattern totally vanished. This result was found similar to other studies (1, 2, 11, 12), and will be applied as the maximum spacing between the centers of neighboring pyramids, as well as the maximum side length of the pyramid base, to ensure the homogeneity of lateral dose distribution. Naturally, spacing and side length smaller than the critical value would result in higher capacity of fine-tuning the dose distribution, conditioned practicality, and technical feasibility.

3.1.2 Validation of SOBP Pyramid Design

The set of 40 mini pencil beamlet dose kernels used for inverse optimization in the design of pyramid shape were obtained using FLUKA MC simulation by directing a pristine 184-MeV proton mini pencil beamlet with a radius of 0.01 cm into a water medium through Lucite blocks of various thickness placed at 20 cm above the water surface. By superimposing this set of dose kernels using MATLAB, a typical SOBP curve emerged ( ). The dose uniformity at the SOBP plateau region presented with a maximum deviation of 0.2% within the plateau region with 2 mm inner margin from both sides (13), and the sharpness of distal dose fall-off, characterized by the depth interval between distal 80% and 20% dose, L (13), being 0.8 cm.

Figure 4

(A) SOBP curve formed by the integration of weighted individual Bragg peaks. (B) A pyramid is created so that the net area of each step corresponds to the weighting of each individual Bragg peaks. The base of the pyramid serves as a compensator to constrain the range the pristine proton beam penetrates.

Using such set of dose kernels and the SOBP characteristics of , the pyramid shape was obtained by the modified version of the basic gradient descent algorithm as described in Equation (1). The resulting pyramid has a base size of 0.5 cm × 0.5 cm, equal to the critical spacing, and consists of 26 layers or levels, with each layer having a height of 0.2 cm, except for the base having a height of 1.2 cm ( ). The 26 levels correspond to the number of individual modulated energies that would be superimposed to form the desired SOBP. The heightened base served as the compensator to shift the range distally.

A uniform SOBP modulator of size 5×5 cm2 was digitally constructed with one hundred such individually designed pyramids, placed side by side ( ). The dose distribution in water resulting from a 184-MeV mono-energetic proton beam of field size 4×4 cm2 was computed using FLUKA MC simulation. Eighty (80) batches of 50,000 primary particles (4,000,000 total particles) were simulated with different random seeds for each batch. The averaged resultant uncertainty was 3.9% in the range of 10% to 100% dose levels. Comparing the percentage depth dose curves taken along the central beam axis, the average difference between MATLAB calculation and the FLUKA MC simulation was 0.5% ( ), with both results normalized to the maximum dose. The MC simulation showed excellent dose uniformity at the plateau region, with a maximum deviation of 0.4%; the sharpness, L of 0.65 cm, well agreed with the results from MATLAB calculation. The transverse homogeneity at different depths within the SOBP plateau region in a 2 × 2 cm2 area centered on the beam axis was evaluated by homogeneity index, defined as the ratio of the maximum dose to the mean dose ( ). The averaged value of the homogeneity index is 1.1, with the standard deviation in the magnitude of 10-5. All MC simulation data met SOBP modulator design objectives.

Figure 5

3D rendering of proton beam modulator with pyramid as its fundamental structure.

Figure 6

Comparison between the MatLab-calculated SOBP curve and the FLUKA-simulated SOBP curve. These two curves closely coincide with each other.

Table 1

Homogeneity at different depths within the SOBP simulated by FLUKA.

Depth (cm)	15.5	16.5	17.5	18.5	19.5
Homogeneity Index	1.088	1.105	1.099	1.096	1.126
Standard Deviation	3.14e-5	3.55e-5	3.75e-5	3.77e-5	5.20e-5

3.2 Optimized 3D Conformal Dose Distribution From the 3DCM

The 3DCM for the test sphere was developed digitally using the optimization algorithm described in the previous section. The outputs of the optimization procedure included (1) the geometry of the 3DCM that can be used by a 3D printer to fabricate the physical 3DCM, and (2) the corresponding proton fluence map to be eventually transferred to the treatment planning system to construct the clinical pencil beam spot sequence (positions and fluences). The first stage of the optimization resulted in 100% target coverage, and the corresponding nCI of 1.514. Following the second stage of optimization, the conformal index nCI was further improved to 1.057 ( ). The geometry of the 3DCM has two parts: the compensator serving as the bases of all pyramids and to modulate the range to conform the dose to the distal end of the target, and the SOBP modulator formed by 489 pyramids with various shapes (heights and graduation), spreading the dose along the beam axis and conforming the dose to the proximal end of the target ( , pruned).

Figure 7

Resultant dose distribution from the optimizations (units on axes are in cm). (A) Dose distribution of the first stage of optimization with conformity 1.514; (B) dose distribution after the adaptive ring optimization; the dose conforms better to the target with conformity 1.057.

Figure 8

Sagittal cross-sectional view of the modulator (the base colored in blue serves as the compensator, and the pyramid bodies colored in orange modulate the SOBP region).

The design was validated by FLUKA MC simulation. Due to the limitation of the number of geometric bodies that can be imported into FLUKA, the minimal weighting factors (minuscule geometric bodies) were discarded. For each pyramid, the steps with weighting factors less than 5% of that pyramid’s total weighing factors were removed. With the pruned geometry, the dose distribution was re-computed using ray-tracing algorithm by MATLAB, showing 95% target coverage by 90% of the isodose surface with a dose conformity index nCI of 1.09. The pruned geometry ( ) and the proton fluence map were then imported into FLUKA for the MC simulation. A total number of 24.46 million primary particles were used for simulation. The simulation was repeated five times, resulting in an averaged uncertainty of 5.2%. The analysis of the resulted dose distribution ( ) showed that the target volume covered by the 87% isodose surface was 95.0%, with nCI = 1.12, closely consistent with the MATLAB results.

Figure 9

Axial and sagittal view of the dose distribution generated by FLUKA MC Simulation (units on axes are in cm). The thickened line represents 86% isodose line (maximum dose is normalized to 100%), covering 95.5% of the target volume.

4 Discussion

The key to developing a 3D conformal modulator (3DCM) with single proton energy relies on the spatial resolution of SOBP modulations. Small pyramid-like cones can serve as a basic unit to construct such non-uniform SOBP modulator customized to tumor shapes. Although proof of concept of such methodology has been reported by prior research groups (1–5), applying such technique in clinical settings requires a robust mechanism to fabricate patient-specific 3DCM. With forward designing of SOBP pyramids, the inter-contribution of lateral dose distribution from adjacent pyramids cannot be effectively accounted for; in addition, the interplay between the shape of pyramids and the proton fluence distribution, crucially important to improve target conformity and critical organ avoidance, cannot be easily optimized. The inverse algorithm and procedure presented in this study offer an effective approach to take into account relevant constraints to generate the optimized geometry of 3DCM with the corresponding proton fluence map.

In this paper, we have shown that the size and the spacing between pyramids are essential in assuring the quality (conformity and homogeneity) of resultant dose distribution. Critical spacing varies with proton energies. With 184 MeV as our example, r 40% or 5 mm was shown to be the critical pacing. The side length of the pyramid base should not exceed critical spacing in order to assure the lateral dose equilibrium and the smoothness of dose distribution. Naturally, a smaller side length of the pyramid base would increase the fine-tuning capability of dose distribution and dose homogeneity. However, other factors, such as the quality of 3D printer and the rigidity of modulator material, should be taken into consideration in determining the size of the pyramids.

In principle, the algorithm and methodology presented in the paper are not limited to the shapes and size of treatment targets, as long as the target volume is not fragmented along the beam path. In addition, since the setting of dose objectives to specific points or structure is done through a programmable module, constraints related to dose-limiting critical organs can be taken into account in the optimization process similar to the ways the “rings” were treated.

As the proposed strategy, the mini beamlets, rather than clinical pencil beams are used in the first stage of obtaining the optimal geometry of the 3DCM and the corresponding quasi-continuous 2D fluence map. This is to be followed by the second stage of creating a spot-scanning layer of clinical pencil beams that would closely reproduce the same fluence map. This can be achieved by inverse optimization on most of commercially available treatment planning systems. The geometric configuration of the 3DCM is not critically sensitive to the spot size of clinical pencil beams and their arrangement. As long as the fluence map can be accurately reproduced, the spot sizes and their positional arrangement are in principle inconsequential.

For proof of concept, the proton beam of 184 MeV with no energy spread was used in this study. In actual clinical application, the dose kernel of mini beamlets with energy spread consistent with the clinical beams would be generated and used for 3DCM design. Again, the actual clinical pencil beams with specific spot sizes are not to be used in the first stage of optimization process to produce the 3DCM geometry and the 2D fluence map.

In this study, all beamlets were arranged in parallel. The divergence correction based on the characteristics of actual clinical proton systems will be introduced in the next phase of the project development.

The 3D-printed 3DCM is to be placed between the beam nozzle and the patient. Like any other traditional beam modulators, a 3DCM should be aligned with the beam axis with high precision. Due to the fact that the 3DCM is in essence to be aligned with the optimized fluence map rather than with individual pencil beam, it is therefore inherently more robust.

Clinically, using a single energy layer to cover the whole target volume might incur longer beam-on time. This challenge could potentially be overcome by the higher dose rates brought forth by newer accelerator technologies. With proton systems capable of delivering non-transmission ultra-high-dose rates or FLASH treatment (14), we envision that such 3DCM would be very valuable for SRS/SBRT treatment, with the patient motion-related error eliminated or minimized.

By incorporating biological models, the application of the current method could be extended to heavier particles, such as carbon ions (15).

5 Conclusion

This proof-of-concept study demonstrated the feasibility of using a novel inverse optimization algorithm to design the geometry of 3DCM, paired with a proton fluence map to achieve optimal 3D conformal dose distribution using a single energy layer with a PBS system. The method has the flexibility of handling tumors with complex shapes. With the successful validation by MC simulation, the physical fabrication of 3DCM and dosimetric verification with phantom measurement will follow. With further development, especially the integration with commercial planning systems, the robust clinical implementation of 3DCM can be anticipated.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author Contributions

ZZ, WZ, and XW were responsible for the study design and contributed to the writing of the manuscript. XW was responsible for developing the initial concept. ZZ was additionally responsible for MATLAB and MC simulation. MB was responsible for providing proton beam parameters, validating data analysis, and manuscript editing. All authors contributed to the article and approved the submitted version.

Conflict of Interest

XW is a managing director of the company Executive Medical Physics Associates.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.