© EDP Sciences 2019

1 Introduction

The use of voxelized volumes (CT scans, anthropomorphic phantoms) with Monte Carlo (MC) simulations for dosimetric purposes is widespread in medical physics (Abella et al., 2010; Candela-Juan et al., 2013; Figueira et al., 2013).

The dosimetric results of MC simulations can be expressed through dose maps (3D matrices), where the voxelized volume is divided following a scoring grid with a user-defined resolution. The dose scoring element of the grid is generally called dosel by analogy with the term voxel (Sarrut and Guigues, 2008). The deposited energy is tracked inside each dosel to compute the absorbed dose.

In order to perform precise comparisons among voxelized volumes, it can be useful to modify the dosel resolution of each resulting absorbed dose map to the same value regardless of the original resolution of the voxelized volumes. In addition to voxelized volume comparisons, another challenge of MC dosimetry is to reduce simulation computing time. Since the dosel statistical uncertainty is correlated with its volume, decreasing the dosel resolution obviously decreases the dose uncertainty. A reduction of dose map resolution can thus provide a notable gain in simulation computing time.

Resampling methods can be applied at three different levels to obtain a dose map with a different resolution than that of the original voxelized volume:

• pre-simulation: resampling the voxelized volume before the simulation;

• per-simulation: generating the dose map following a defined scoring grid resolution during the simulation;

• post-simulation: resampling the dose map after the simulation.

A simple pre-simulation solution consists of resampling the Hounsfield units (Schneider et al., 2000) (for a CT scan) or index numbers (for voxelized phantoms). In the case of a CT scan, averaging the Hounsfield number of a group of voxels gives interesting results (Mora et al., 2001; Ai-Dong et al., 2005). Applying these methods to the index numbers¹ of voxelized phantoms will not work at interfaces between different index numbers because they are not related to the material density of the voxel (like the Hounsfield units). It will thus skew the reconstruction of the materials in the simulation. It should be noted that more complex solution can be used to resample voxelized phantom, like assigning the majority index value between the resampled voxels to the new voxel.

Post-simulation methods can be easier to implement, as they can consist of resampling a dose map. According to the definition of the absorbed dose (i.e., deposited energy over the mass of the medium), these methods will work only with materials of similar density. If the material compositions are too different, the dose resulting from the resampling can diverge significantly from the absorbed dose in each of the materials.

Concerning per-simulation methods, superimposing a scoring grid over the geometrical voxel grid of the volume allows to generate dose maps with no limitation of dosel resolution. This method preserves all the physical interactions by using the original geometrical information of the voxelized volume, allowing it to work with any type of voxelized volumes. It will bring MC codes to compute absorbed dose inside dosels which contain more than one material. Applying the definition of the absorbed dose to such heterogeneous volumes is a complex issue and the method used to reconstruct the dose in such volumes will necessarily impact its accuracy.

In the present work, we propose a comprehensive study of two per-simulation algorithms allowing the calculation of the heterogeneous dose based on weighting the dose in the voxel by either its volume or mass fraction. The first of the two algorithms is the default dose calculation of the GATE Monte Carlo software (v8.0), a GEANT4 framework dedicated to medical physics (Sarrut et al., 2014), and the second algorithm is proposed by the authors. A large number of publications used the GATE software for dosimetric and radiological protection calculations (Perrot et al., 2014; Bouzid et al., 2015; Marcatili et al., 2015). As the choice of a given algorithm can strongly affect the computed dose values, it is interesting to compare the two algorithms in terms of bias and uncertainty on the reconstructed doses in heterogeneous voxelized volumes.

2 Materials and methods

Two methods for calculating the absorbed dose in a heterogeneous medium are described in this paper. They are based either on weighting of the absorbed dose in each of the dosel materials by its volume or mass fraction and they will be referred to as Volume-Weighting (V-W) and Mass-Weighting (M-W) algorithms respectively. The Volume-Weighting algorithm is the one used by default in the MC simulation code GATE. The Mass-Weighting has been recently added to the GATE code by the authors of the present work.

Two grid types are present in the MC simulation, the voxel geometrical grid of the volume and the dosel grid that scores the dose, each having its own resolution. The scoring grid is superimposed on the geometrical grid. GATE uses a system called DoseActor to compute the absorbed dose (Sarrut et al., 2014). This system handles the scoring grid and its resolution to store the absorbed dose in a 3D matrix called a dose map. All the simulations performed in this work used the built-in physics list emlivermore with the photon and electron cuts set to a range of 1 mm. The compositions and densities of the air (G4_AIR) and vacuum (G4_Galactic) materials came from the GEANT4 NIST material database. The compositions and densities of the tissues and organs of the voxelized phantom are taken from the publication 110 of ICRP (ICRP, 2009).

2.1 Volume-Weighting algorithm

To compute the dose inside a dosel (D_dosel) the V-W algorithm weights the absorbed dose (D_i) of each of the M voxels contained in the dosel by their respective volume fraction (V_i/V_dosel). (1) with V_i the portion of volume of the voxel i contained in the dosel. Using the definition of the absorbed dose, equation (1) can be expressed in terms of voxel deposited energy E_i and mass m_i. (2)

2.2 Mass-Weighting algorithm

The M-W algorithm takes into account the heterogeneity of the dosel by using the absorbed dose definition. At the initialization step of the simulation, the M-W algorithm computes the mass of each dosel (m_dosel) by taking into account the portion of mass (m_i) of each of the M voxels contained in it. (3) with ρ_i the density of the voxel i contained in the dosel. The computation takes less than one minute on a full voxelized anthropomorphic phantom and should only be performed once per voxelized volume. Then, during the simulation, the algorithm divides the deposited energy in each dosel (E_dosel) by their respective masses to get the dosel absorbed dose (D_dosel). (4)

Equation (4) can be written in an equivalent form where the voxel dose is weighted by its mass fraction: (5)

If the dosel and voxel resolutions are equal, a unique material is attached to each dosel (homogeneous case) and so m_dosel = m_i and V_dosel = V_i. Hence, the two algorithms follow the absorbed dose definition and compute the same dose values.

2.3 Relative statistical uncertainty

In GATE Monte Carlo simulations, the relative statistical uncertainty ϵ_k of the absorbed dose in the dosel k is computed using (Chetty et al., 2006): (6) with N the number of primary particles and d_k,j is the absorbed dose contribution of the primary particle j in the dosel k.

3 Applications

The introduced algorithms were compared in two simulation setups. First, a simple geometry with two voxels of different size was considered to study the reconstructed heterogeneous dose with respect to the absorbed dose in each material (Fig. 1a). In a second setup, the algorithms were applied to an anthropomorphic voxelized phantom to analyze the effect of dosel resolution on the algorithm’s dose reconstruction (Fig. 1b).

Fig. 1 Illustration of the two applications studied.

3.1 Simple geometry

The simple geometry consists of a cube of size 10 × 10 × 10 cm³ containing two voxels with different materials irradiated tomographically with a square beam (with dimensions greater than the dimensions of the cube) of 100 keV X-rays at eight angles separated by 45°. A beam of 10⁷ primary X-rays was simulated in each simulation. The pairs of materials considered in this study are: water-vacuum, water-bone and water-air. To study the dose variation in the entire cube reconstructed with the algorithms, several MC simulations with different proportions of the two materials were generated.

3.1.1 Results

3.1.1.1 Water-vacuum case

Figure 2a illustrates the absorbed dose in each of the materials and the reconstructed dose in the whole cube as a function of the volume ratio of each material. Knowing that the only mass contribution comes from the water volume, it is expected that the M-W dose will follow the water dose. Thus, it should remain the same regardless of the material proportions in the volume. The shape of the V-W dose is determined by the proportion of water volume in the simulation. Thus, the resulting dose corresponds to the water dose weighted by the fraction of the water volume. Figure 2b shows that the statistical uncertainties of the two algorithms follow those of the water since there is no contribution from the vacuum.

Fig. 2 (a) Absorbed dose in water and vacuum and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-vacuum case. (b) Relative statistical uncertainty on dose (statistical uncertainty of vacuum is not represented). : water; : vacuum; : Mass-Weighting; : Volume-Weighting.

3.1.1.2 Water-bone case

Figure 3a displays the absorbed dose in water and bone² and the heterogeneous dose reconstructed by the algorithms. With this combination of materials, the M-W dose is always higher than the one given by V-W. It is explained by the significant contribution of absorbed dose in bone due to its higher density (1.92 g/cm³). However, the difference between the two algorithms is significantly lower than in the water-vacuum case because of the smaller difference between the material densities. Figure 3b shows that the relative uncertainty of dose in water and bone are of the same order of magnitude. But for a bone volume ratio lower than about one-third of the total volume, the number of interactions becomes higher in the water volume than in the bone volume. As the M-W algorithm gives more weight to the bone in the dose calculation, this results in a slightly higher M-W statistical uncertainty at low bone volume ratio compared to the V-W one (up to 4.8% for a ratio of 0.3).

Fig. 3 (a) Absorbed dose in water and bone and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-bone case. (b) Relative statistical uncertainty on dose. : water; : vacuum; : Mass-Weighting; : Volume-Weighting.

3.1.1.3 Water-air case

Figure 4a shows the results for a volume containing water and air. As for the water-bone case, the dose reconstructed by the M-W algorithm is higher than the one of V-W (up to 17%). Due to the very low density of air (1.20479 mg/cm³), the M-W dose receives almost exclusively the contribution of the absorbed dose in water, while the V-W dose also receives a contribution of the absorbed dose in air that depends on its volume fraction. The difference of statistical uncertainty seen on Figure 4b reflects the contribution of air in each of the algorithms. Because of the low density of air, the small number of particles interacting leads to a high statistical uncertainty on its absorbed dose and thus on the uncertainty of the V-W dose. While having almost no contribution coming from the air, the M-W relative uncertainty follows that of the water and can be up to 90% lower than the one of the V-W algorithm.

Fig. 4 (a) Absorbed dose in water and air and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-air case. (b) Relative statistical uncertainty on dose. : water; : vacuum; : Mass-Weighting; : Volume-Weighting.

3.1.2 Discussion

The two algorithms present some significant discrepancies on the reconstructed dose depending on the materials involved at the interface, especially for materials with a large difference of density. Knowing that there is more deposited energy in heavier materials, higher dose values are computed by the M-W algorithm. The behavior of the M-W algorithm is particularly useful at air-tissue interfaces, providing a more accurate estimation of absorbed dose in tissue and a statistical uncertainty not influenced by the air. These basic setups show that the reconstructed dose for a heterogeneous dosel always falls between the absorbed dose of each separate material.

3.2 Voxelized phantom

The algorithms were studied on the female version of the ICRP 110 anthropomorphic voxelized phantom (ICRP, 2009). A 10 pixels slice, located at the thoracic level, was used in the MC simulations. Its resolution is 299 × 137 × 10 for a voxel size of 1.775 × 1.775 × 4.84 mm³. The slice was placed in air and irradiated tomographically by 80 keV X-rays at eight angles. 8 × 10⁹ X-rays were simulated in order to have an average statistical uncertainty on the dose inferior to 1% at the original voxel resolution. The dose maps were computed with the two algorithms for the following dosel resolutions: 299 × 137 × 10, 150 × 68 × 5 and 75 × 34 × 3, corresponding respectively to a dosel volume of 15.1, 121 and 817 mm³ and a resampling factor respectively equal to 1, 8 and 64. Each of the following dose and uncertainty maps were generated in the same simulation.

3.2.1 Results

3.2.1.1 299 × 137 × 10 dosel resolution

Figure 5 shows dose and relative statistical uncertainty maps at the dosel resolution 299 × 137 × 10 (same as the original voxel resolution). No matter which algorithm is used, the generated maps are exactly the same. This was expected because each dosel contains only one voxel and thus one material. The average relative statistical uncertainty on Figure 5b is presented in Table 1 and is equal to 0.97%.

Fig. 5 Dosel resolution 299 × 137 × 10 (same resolution as the ICRP phantom).

3.2.1.2 150 × 68 × 5 dosel resolution

Figures 6a and 6b show the dose maps obtained with the V-W and M-W algorithms for a dosel resolution eight times lower than the original voxel resolution. In order to show the discrepancies between the two algorithms, Figures 7a and 7c display respectively the ratio of reconstructed dose and the ratio of their relative uncertainties. Figure 7a shows that the V-W reconstructed dose is lower at heavy-light material interfaces (bone-tissue, tissue-air) compared to the M-W algorithm. One can see in Table 2 that the dose difference can be up to 49% (3.6% on average). As seen in Table 1, at the dosel resolution 150 × 68 × 5 the mean relative statistical uncertainty on the reconstructed dose is 36% lower with the M-W algorithm. Figure 7c shows that the M-W algorithm gives a smaller uncertainty than the V-W algorithm at air-tissue and lung-tissue interfaces, but the M-W uncertainty is slightly higher than the V-W one at the bone-tissue interface. Table 2 shows that the M-W relative uncertainty can be up to 88% lower at the air-tissue interface and 9% higher than the V-W one at the bone-tissue interface, which is consistent with the heterogeneous algorithms behavior described in the simple geometry application.

Fig. 6 Dose maps from Volume-Weighting (V-W) and Mass-Weighting (M-W) algorithms for scaling factors of 8 and 64. Dosel resolution of Figures (a) and (b): 150 × 68 × 5. Dosel resolution of Figures (c) and (d): 75 × 34 × 3.

Fig. 7 Ratio of dose maps and relative uncertainty maps from Volume-Weighting (V-W) and Mass-Weighting (M-W) algorithms for scaling factors of 8 (×8) and 64 (×64). Dosel resolution of Figures (a) and (c): 150 × 68 × 5. Dosel resolution of Figures (b) and (d): 75 × 34 × 3.

3.2.1.3 75 × 34 × 3 dosel resolution

Figures 6c and 6d show that at a 64 times smaller dosel resolution (75 × 34 × 3) the discrepancies between the two algorithms are more pronounced. Table 2 shows that, in this case, the dose reconstructed with the V-W algorithm can be up to 43% lower than the M-W one (5.4% on average). Table 1 shows that the mean relative statistical uncertainty on the dose map obtained with the M-W algorithm is 49% smaller compared to the V-W one. At this dosel resolution, to get a mean relative statistical uncertainty on dose of about 1%, the M-W and V-W algorithms need respectively about 35 times and 9 times less computation time compared to the MC simulation at the original resolution.

3.2.1.4 Comparison with the original resolution dose map

In order to compare the dose maps reconstructed with the algorithms with the dose map at the original resolution, each dose map of 150 × 68 × 5 and 75 × 34 × 3 dosels was resampled to the 299 × 137 × 10 resolution using the nearest neighbor interpolation. This allows to highlight the biases induced by the algorithms in comparison to the reference dose map at the voxel original resolution. Table 3 presents the mean dose ratio between the original resolution dose map and the resampled one from the V-W and M-W algorithms for different materials. Figure 8 shows the ratio between the original resolution dose map and the mass and volume weighted resampled ones. The main discrepancies between the heterogeneous algorithms and the original dose map are located at the tissue-bone interface. A comparison shows that the V-W algorithm underestimates the dose on the entire phantom slice by 0.6% at 150 × 68 × 5 and 0.4% at 75 × 34 × 3 on the average. However, when the skin region is compared, the V-W algorithm underestimates the dose by 15 and 17% on average respectively for the 150 × 68 × 5 and 75 × 34 × 3 dosel resolutions. On the other hand, the M-W algorithm overestimates the dose by a factor of 1.2 and 2.4% on average respectively for the 150 × 68 × 5 and 75 × 34 × 3 dosel resolutions. Considering only the skin region, the M-W algorithm underestimates the original dose by a mean factor of 4.1 and 7.0% for the respective dosel resolutions 150 × 68 × 5 and 75 × 34 × 3. The presence of air in the bronchi leads also the V-W algorithm to underestimate the bronchial dose at air-bronchi interfaces, while this problem does not occur with the M-W algorithm. Despite a slight overestimation of the dose in the whole phantom slice and in the majority of the structures, the M-W algorithm gives a better estimation of the dose in the skin and a constant lower standard deviation of the mean dose ratio than that obtained with the V-W algorithm.

Fig. 8 Ratio of the dose maps at original resolution (Dorig) and the resampled one with the Volume-Weighting and Mass-Weighting algorithms (DV-W or DM-W). Figures (a) and (b) are resampled from the resolution 150 × 68 × 5 to 299 × 137 × 10. Figures (c) and (d) are resampled from the resolution 75 × 34 × 3 to 299 × 137 × 10.

3.2.2 Discussion

The simulation setups described previously show that the two proposed methods to calculate dose in heterogeneous dosels produce different results depending on the densities of the materials.

These methods provide an effective solution to reduce statistical variance (up to an average factor of 2.7 with the M-W algorithm with a resampling factor of 8) and a good reconstruction of the absorbed dose in heterogeneous dosels (less than 2% underestimation on average with the M-W algorithm with a resampling factor of 8) but at the expense of a loss in the spatial accuracy of dose deposition. The V-W algorithm gives a good estimation of the mean dose in a whole slice of the phantom. However, V-W can locally underestimate the reconstructed dose, especially at materials interface. Table 3 shows that for all tissues the M-W algorithm gives a higher dose value than the V-W one. This behavior makes the M-W algorithm particularly suitable for radiation protection purposes, which always prefer the most unfavorable estimation with respect to the ALARA principle³. Nevertheless, results on the two-voxels cube have shown that the M-W algorithm can lead to neglecting the dose in the more abundant material of the dosel if its density is too low compared to the other components.

The ability of the M-W algorithm to compute skin dose with a low underestimation (about 7% on average at dosel resolution resampled by a factor of 64) compared to the V-W algorithm (which produces about 10% more underestimation than M-W) is of particular interest for skin dosimetry in radiology procedures.

In addition to the statistical variance reduction caused by the increase of the dosel size, the main dose uncertainty reduction between algorithms occurs at materials interfaces, where the M-W algorithm provides a relative uncertainty which can be up to 80% lower than computed by the V-W algorithm.

It should finally be noted that in the case where only part of a heterogeneous dosel is irradiated, for example in microbeam radiotherapy, the difference between the dose reconstructed by the algorithms and the absorbed dose in the homogeneous voxels can be much larger than those presented in this work.

4 Conclusions

The comprehensive study of the two available algorithms for dose calculation in the GATE v8.0 Monte Carlo software allows to evaluate the bias and uncertainties in the dose maps that contain heterogeneous voxels. At material interfaces like skin-air, the Volume-Weighting method can lead to an underestimation of the dose as high as 20%. The Mass-Weighting method, integrated recently into GATE, enables to limit the bias in the dose reconstruction at material interfaces, and to reduce the uncertainty up to 80% depending on the voxel size. This algorithm should therefore be favored for radiological protection studies based on dose calculations in heterogeneous volumes.

References

All Tables

All Figures

	Fig. 1 Illustration of the two applications studied.
In the text

	Fig. 2 (a) Absorbed dose in water and vacuum and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-vacuum case. (b) Relative statistical uncertainty on dose (statistical uncertainty of vacuum is not represented). : water; : vacuum; : Mass-Weighting; : Volume-Weighting.
In the text

	Fig. 3 (a) Absorbed dose in water and bone and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-bone case. (b) Relative statistical uncertainty on dose. : water; : vacuum; : Mass-Weighting; : Volume-Weighting.
In the text

	Fig. 4 (a) Absorbed dose in water and air and heterogeneous doses computed from the Mass-Weighting and Volume-Weighting algorithms for the water-air case. (b) Relative statistical uncertainty on dose. : water; : vacuum; : Mass-Weighting; : Volume-Weighting.
In the text

	Fig. 5 Dosel resolution 299 × 137 × 10 (same resolution as the ICRP phantom).
In the text

	Fig. 6 Dose maps from Volume-Weighting (V-W) and Mass-Weighting (M-W) algorithms for scaling factors of 8 and 64. Dosel resolution of Figures (a) and (b): 150 × 68 × 5. Dosel resolution of Figures (c) and (d): 75 × 34 × 3.
In the text

	Fig. 7 Ratio of dose maps and relative uncertainty maps from Volume-Weighting (V-W) and Mass-Weighting (M-W) algorithms for scaling factors of 8 (×8) and 64 (×64). Dosel resolution of Figures (a) and (c): 150 × 68 × 5. Dosel resolution of Figures (b) and (d): 75 × 34 × 3.
In the text

	Fig. 8 Ratio of the dose maps at original resolution (Dorig) and the resampled one with the Volume-Weighting and Mass-Weighting algorithms (DV-W or DM-W). Figures (a) and (b) are resampled from the resolution 150 × 68 × 5 to 299 × 137 × 10. Figures (c) and (d) are resampled from the resolution 75 × 34 × 3 to 299 × 137 × 10.
In the text