Fast photon-boundary intersection computation for Monte Carlo simulation of photon migration

Xiaofen Zhao; Hongyan Liu; Bin Zhang; Fei Liu; Jianwen Luo; Jing Bai

doi:10.1117/1.OE.52.1.019001

3 January 2013 Fast photon-boundary intersection computation for Monte Carlo simulation of photon migration

Xiaofen Zhao, Hongyan Liu, Bin Zhang, Fei Liu, Jianwen Luo, Jing Bai

Author Affiliations +

Optical Engineering, Vol. 52, Issue 1, 019001 (January 2013). https://doi.org/10.1117/1.OE.52.1.019001

Abstract

Monte Carlo (MC) method is generally used as a “gold standard” technique to simulate photon transport in biomedical optics. However, it is quite time-consuming since abundant photon propagations need to be simulated in order to achieve an accurate result. In the case of complicated geometry, the computation speed is bound up with the calculation of the intersection between the photon transmission path and media boundary. The ray-triangle-based method is often used to calculate the photon-boundary intersection in the shape-based MC simulation for light propagation, but it is still relatively time-consuming. We present a fast way to determine the photon-boundary intersection. Triangle meshes are used to describe the boundary structure. A line segment instead of a ray is used to check if there exists a photon-boundary intersection, as the next location of the photon in light transports is determined by the step size. Results suggest that by simply replacing the conventional ray-triangle-based method with the proposed line segment-triangle-based method, the MC simulation for light propagation in the mouse model can be speeded up by more than 35%.

1. Introduction

Optical molecular imaging is developing rapidly as a promising technology in the molecular imaging field. The establishment of a mathematical model of light propagation in biological media is essential in optical tomography. Diffuse equation (DE) is commonly used as the mathematical model, however, solving DE needs certain boundary conditions.¹^,² A Monte Carlo (MC) method is relatively flexible, as it does not require an assumption on boundary conditions.

MC in the optical field is basically a statistical method that allows for constructing a stochastic model to simulate a large amount of light propagation trajectories. Since first introduced by Wilson and Adam into the field of laser-tissue interactions,³ it has been widely used to simulate the light propagation in media.⁴^–⁹ Much work has been done to improve the efficiency of a MC simulation for modeling light propagation in turbid media.¹⁰^–¹⁴ Based on the information generated by a single MC baseline simulation, the scaling methods were developed for a wide range of optical properties.¹⁵^,¹⁶ Liebert et al. applied a similar concept to simulate the time-resolved fluorescence imaging in layered turbid media using the MC method.¹⁷ However, a single MC process is still the most time-consuming part owing to a large number of photon trajectory calculations.

In order to improve the efficiency of a single MC simulation for light propagation, Wang et al. developed an MC modeling of photon transmission in multi-layered tissues (MCML)¹⁸ in 1995, while Tinet et al.¹⁰ proposed a semi-analytical Monte Carlo method for time-resolved light propagation. Each interaction point of a photon that transports in a medium contributes directly to a detector area, thus fewer photons are required to produce a given accuracy but with complex calculations. Zolek et al.¹¹ developed an optimized MC algorithm by approximation of the logarithmic and trigonometric functions. Boas et al. proposed a voxel-based model (tMCimg) to improve the capability of modeling turbid structures.¹² However, Binzoni et al. pointed out that the voxel-based model (tMCimg) was difficult to incorporate precise boundary information.¹³ Margallo-Balbas et al., therefore, explored a shape-based method (triMC3D)¹⁴ and used triangle meshes to define the complex boundary structures. This shape-based method gives a more accurate description of the interface than the voxel-based approaches. However, in the case with complicated geometry, the computation speed is bound up with the calculation of the photon-boundary intersection. In order to judge whether the photon hits its boundary, the ray-triangle-based method proposed by Möller et al.¹⁹ is often used to calculate the photon-boundary intersection.¹⁴^,²⁰ Even though this method can save the memory required, it needs to calculate the intersection between the ray and every triangle; therefore it is quite time-consuming. In fact, an intersection test between a line segment and every triangle can be taken before the calculation of the photon-boundary intersection, as the next location of the photon in light transports is determined by the step size. Thus, the photon-boundary intersection needs to be calculated only if an intersection exists within a triangle. When there is more than one intersection, it means that there are many boundaries along the line segment. Then only the intersection with the minimum distance to the photon’s current position is calculated. This method is more effective than the method that considers all photon-boundary intersections, and makes it possible to reduce computation time.

In this paper, a fast way to determine the photon-boundary intersection in the shape-based MC simulation is presented. Triangle meshes are used to describe the boundary structure. A line segment, instead of the ray¹⁴^,²⁰ is used to check if a photon-boundary intersection exists, as the intersection will be calculated only when it does exist. Therefore, only a part of triangles have to be taken into account when calculating the photon-boundary intersection, which greatly reduces the computation time. The conventional method and the proposed method for the intersection test in MC simulation for light propagation are detailed first. Then, the validation of the proposed algorithm is conducted in a semi-infinite medium. The fluence distribution within the medium and the remitted flux of photons obtained from MC with the proposed method are compared with those obtained from the analytic solution of DE and the conventional MC method, respectively. Finally, several sets of MC simulations under different conditions are executed to test the computational speed of the proposed method.

The paper is organized with the conventional ray-triangle-based method and the proposed line segment-triangle-based methods described in Sec. 2. In Sec. 3, three sets of simulation studies are conducted to validate the algorithm and test the computational speed. The results are shown and discussed. In the last section, the conclusion is presented.

2. Method

2.1.

The Ray-Triangle-Based Method

A ray $R (t)$ with origin $O$ and normalized direction $D$ is defined as:

Eq. (1)

R (t) = O + t D,

where

t

is the unknown distance the ray travels. And a triangle is defined by three vertices

V_{0}

,

V_{1}

, and

V_{2}

. A point,

T (u, v)

, in the plane where the triangle lies is given as:

Eq. (2)

T (u, v) = (1 - u - v) V_{0} + u V_{1} + v V_{2},

where

(u, v)

are the coordinates within the triangle which satisfy

u \geq 0

,

v \geq 0

and

u + v \leq 1

. The intersection between the ray,

R (t)

, and the triangle,

T (u, v)

, can be calculated from Eqs. (1) and (2) as follows:

Eq. (3)

O + t D = (1 - u - v) V_{0} + u V_{1} + v V_{2} .

Reordering the formula, we can get¹⁹

Eq. (4)

[\begin{matrix} - D, & V_{1} - V_{0}, & V_{2} - V_{0} \end{matrix}] [\begin{array}{l} t \\ u \\ v \end{array}] = O - V_{0},

where

D = s / μ_{t} \cdot \vec{n_{μ}}

,

s

is the free path of the photon,

\vec{n_{μ}}

is the direction vector of the photon migration and

μ_{t}

is the total interaction coefficient equal to the sum of the absorption coefficient and the scattering coefficient. It means that

t

, the distance from the ray origin to the intersection point, and

(u, v)

, the intersection coordinates, can be obtained by solving Eq. (4).¹⁹ But this approach is very time-consuming.

The ray-triangle intersection algorithm tries to determine if the ray intersects with the triangle. If the distance $t$ is in the range of (0,1), the ray-triangle intersection is within the distance of the next free path of the photon. In addition, when $(u, v)$ satisfies $u \geq 0$ , $v \geq 0$ and $u + v \leq 1$ , the ray-triangle intersection is within the triangle. Only when the solution meets these two conditions above at the same time, the next free path of the photon will intersect with the triangle defined by the three vertices $V_{0}$ , $V_{1}$ and $V_{2}$ . In the ray-triangle-based method, solving Eq. (4) is implemented as a translation and change of base of the ray origin and a unit triangle.¹⁹ The procedure of the conventional ray-triangle-based method, when determining whether the photon hits the boundary in the MC simulation, was demonstrated in a flowchart, as shown in Fig. 1.

Fig. 1

Flowchart of the ray-triangle-based method.

2.2.

Line-Segment Triangle-Based Method

However, the ray-triangle-based method has to take all triangles into account when calculating the intersection, which increases the computation time. In fact, when we refer to the path of the photon, it is a line segment rather than a ray, since every step of the photon migration is limited by the step size, as shown in Fig. 2. By taking advantage of this, the proposed method uses the line segment instead of the ray to compute the intersection. Therefore, the proposed method first checks whether a line segment-triangle intersection exists. Only when the intersection exists, can the calculation of the intersection continue.

Fig. 2

The photon migration within a limited step. (a) the photon which does not intersect with the triangle in the next step of migration. (b) photon which intersects with the triangle during the next step of migration.

The standard equation of the plane in which the triangle $V_{0} V_{1} V_{2}$ lies is given by

Eq. (5)

A \cdot x + B \cdot y + C \cdot z + D D = 0,

where

(A, B, C)

is the unit normal vector of the triangle and

D D

is the perpendicular distance from the origin to the plane. Before the start of photon transport simulation, all the standard equations of the planes in which the triangles lie are solved and stored.

By incorporating Eqs. (1) and (5), the following equation is obtained

Eq. (6)

t = - \frac{O \cdot (A, B, C) + D D}{D \cdot (A, B, C)},

when

t

is in the range of (0, 1), the line segment intersects with the plane in which the triangle lies.

If an intersection does exist, the next step is to determine whether the intersection is within the triangle. As illustrated in Fig. 3, supposing there is a triangular prism with the triangle that describes the media boundary being the base and three rectangular sides perpendicular to the triangle, the problem is transformed to determine whether the point of intersection is inside the triangular prism. If the point of intersection is within the triangular prism, the intersection is also within the triangle in the plane. As shown in Fig. 3, two line segments with the same origin $O$ intersect with the plane in which a triangle (the blue part) lies respectively. The points of intersection are $P 1$ and $P 2$ , and $P 2$ is the point where photon hits the boundary.

Fig. 3

Illustration of how to determine whether the intersection is within the triangle in the proposed method.

The equations of three rectangular sides are defined respectively as:

Eq. (7)

S_{0} \cdot x + S_{1} \cdot y + S_{2} \cdot z + S_{3} = 0 S_{4} \cdot x + S_{5} \cdot y + S_{6} \cdot z + S_{7} = 0 S_{8} \cdot x + S_{9} \cdot y + S_{10} \cdot z + S_{11} = 0,

where

(S_{0}, S_{1}, S_{2})

,

(S_{4}, S_{5}, S_{6})

,

(S_{8}, S_{9}, S_{10})

are the unit normal vector of the three rectangular sides respectively, and

S_{3}

,

S_{7}

,

S_{11}

are the perpendicular distances from the origin to each plane. Coordinates of the median point,

(x_{0}, y_{0}, z_{0})

, which are inside the triangle, can be obtained from the three vertices. The intersection point,

(x_{i}, y_{i}, z_{i})

, can be solved by substituting

t

into Eq. (1). If conditions in Eq. (8) are all met, the intersection is inside the triangular prism.

Eq. (8)

(S_{0} \cdot x_{i} + S_{1} \cdot y_{i} + S_{2} \cdot z_{i} + S_{3}) \cdot (S_{0} \cdot x_{0} + S_{1} \cdot y_{0} + S_{2} \cdot z_{0} + S_{3}) > 0 (S_{4} \cdot x_{i} + S_{5} \cdot y_{i} + S_{6} \cdot z_{i} + S_{7}) \cdot (S_{4} \cdot x_{0} + S_{5} \cdot y_{0} + S_{6} \cdot z_{0} + S_{7}) > 0 (S_{8} \cdot x_{i} + S_{9} \cdot y_{i} + S_{10} \cdot z_{i} + S_{11}) \cdot (S_{8} \cdot x_{0} + S_{9} \cdot y_{0} + S_{10} \cdot z_{0} + S_{11}) > 0 .

The procedure of the proposed line segment-triangle-based method in calculating the photon-boundary intersection in MC simulation was illustrated in a flowchart, as shown in Fig. 4.

Fig. 4

Flowchart of the line segment-triangle-based method.

To sum up, the key concept of the line segment-triangle method is based on the idea that the step size can reduce the possible triangles intersected with the photon trajectory. Besides, the proposed algorithm can work out each intersection coordinates with fewer calculations. Thus the computational consumption of the proposed method can be greatly reduced.

The MC model with the proposed method is coded in C programming language and compiled with Microsoft Visual C++ 6.0. For comparison, both the MC methods were performed on the same PC workstation (2.00 GHz Intel Xeon processor, 24 GB RAM) and used the same compiler. The analytical solution of the DE is calculated with MATLAB 7.6.0.324 (R2008a, The MathWorks). The Henyey-Greenstein (HG) phase function is employed for sampling scattering directions in the simulations.²⁰ Fresnel’s equations are used to calculate the reflections and refractions at the boundaries where the refractive index mismatch occurs. These internal functions are just the same as in the conventional MC method.

3. Results and Discussion

3.1.

Validation of the Algorithm in a Semi-infinite Medium

In this work, we simulated the photon propagation in biological medium using the well-established MC procedures.⁹^,¹⁰^,¹²^,¹⁴^,²¹^–²³ Two sets of simulation studies were conducted to validate the proposed algorithm in a semi-infinite medium in this section. The first set of simulation studies was to compare the results obtained from the proposed method in MC simulation with those obtained from the analytic solution of DE. Owing to the diffusion approximation,¹^,²⁴^–²⁶ the simulations were executed in a relatively ideal model available for DE to yield an analytic solution accurate enough to be a validation standard.¹²^,²¹^,²² The second set was to validate the proposed algorithm by comparing the results with those obtained from the conventional MC method. Compared with DE, the main advantage of the MC method is that it does not require an extra simplifying assumption. Thus this makes it capable of dealing with any situation. Considering the time consumed, the second set of simulations were conducted in a homogenous medium of the same geometry and the same source settings as that used in the first set of simulations. The optical properties were set according to that of biological tissues. There were $10^{7}$ photons launched in both MC simulation sets.

3.1.1.

Comparison to the analytical solution of DE

The structure of the medium used in the simulations was a cubic model whose surface was described by 1794 triangles generated by a commercially available software COMSOL Multiphysics 3.5 (COMSOL, Inc., Burlington, MA). The medium had a dimension of $60 \times 60 \times 40 {mm}^{3}$ with its center positioned at $(0, 0, 0) mm$ . The optical parameters of the homogeneous medium were as follows: absorption coefficient $μ_{a} = 0.005 {mm}^{- 1}$ , scattering coefficient $μ_{s} = 1.00 {mm}^{- 1}$ , relative refractive index $n = 1.33$ , anisotropy factor $g = 0.1$ . The ambient medium is the air with refractive index $n = 1$ . A collimated point source was positioned at $(0, 0, 20) mm$ with the direction vector of $(0, 0, - 1)$ .

The medium was approximately semi-infinite as the source was adequately far from the edges. The fluence distribution within the medium and the remitted flux of photons of the proposed MC method and the DE method are compared in Fig. 5. A decent agreement between both methods is observed. In Fig. 5(b), the mismatch near the source is due to the isotropic point source hypothesis within the DE. The slight discrepancy near the edge at 30 mm results from the influence of boundary conditions within the DE for the differences will narrow when enlarging the medium size.¹²^,²¹

Fig. 5

Comparisons of the photon fluence and photon flux generated by the proposed MC method and DE method. (a) The fluence field within the medium of each method, with contours drawn every 5 dB in the plane $y = 30 mm$ . (b) The flux of photons exiting the medium. A base 10 logarithmic scale is used for the $y$ -axis.

3.1.2.

Comparison to the results obtained from the conventional MC method

In this section, the comparison between the results from the MC simulation with the proposed method and with the conventional method was made. The comparison was conducted under the same settings as that used in Sec. 3.1.1, such as geometry size, pseudorandom generator, source settings and photon number. The optical parameters of the biological tissue were as follows: absorption coefficient $μ_{a} = 0.005 {mm}^{- 1}$ , scattering coefficient $μ_{s} = 1.00 {mm}^{- 1}$ , relative refractive index $n = 1.33$ , anisotropy factor $g = 0.9$ . The ambient medium was the air with refractive index $n = 1$ .

The fluence distribution within the medium and the remitted flux of photons of the proposed method and the conventional method are compared in Fig. 6. A good agreement between both methods is observed. Figure 6(a) shows the contour plot in the plane $y = 30 mm$ of the fluence distribution within the medium obtained from each method. The contour interval is 5 dB. The slight discrepancy is due to the number of photons used in the simulations since the MC method is a statistical model based on the stochastic nature of light propagation. Figure 6(b) compares the remitted flux of photons on the surface along the distance from the source. It should be noted that since both methods use surface triangles as detectors in the simulations, the remitted flux value for a target point is derived from the nearby triangles through interpolation. This helps to smooth the data to a certain degree.

Fig. 6

Comparisons of the photon fluence and photon flux generated by the proposed MC method and the conventional MC method. (a) The fluence field within the medium of each method, with contours drawn every 5 dB in the plane $y = 30 mm$ . (b) The flux of photons exiting the medium. A base 10 logarithmic scale is used for the $y$ -axis.

3.2.

Comparisons of Computation Time in Complex 3-D Geometry

In the following experiments, we compared the computation time of the proposed method and the conventional method. The three-dimensional (3-D) geometric structure of a mouse model obtained from an in-vivo experiment was used in this simulation. The torso part from the oxter to the abdomen of the mouse model with a height of 40 mm was chosen as the region to be investigated. The geometry was discretized in COMSOL Multiphysics using triangle meshes, as shown in Fig. 7. The optical parameters of the mouse model were listed in Table 1. The surrounding medium was still air with refractive index $n = 1$ . A collimated point source was positioned at $(10, 0, 10) mm$ with the direction vector of $(- 1, 0, 0)$ . In addition to the geometry and optical characteristics of the medium, the total number of the triangle meshes describing the boundary of the mouse model also had a great impact on the computational speed. In our study, we compared the computation time in four cases with different numbers of triangle meshes when simulating with $10^{4}$ photons: $TriNum = (4952, 5164, 6084, 13602)$ . The results are shown in Table 2. We also compared the computational time of the conventional and proposed methods with four different numbers of simulated photons: $PhotNum = (10^{4}, 10^{5}, 10^{6}, 5 \times 10^{6})$ , when 6084 triangle meshes were used. The results are shown in Table 3. As shown in Tables 2 and 3, the computation time can be decreased by 35% to 45% by simply replacing the conventional method with the proposed method in MC simulations for light propagation in the mouse model. The proposed line segment-triangle-based method is equivalent to the conventional ray-triangle-based method in mathematics, which guarantees accuracy but the former is faster than the latter in the shape-based MC simulation for light propagation.

Fig. 7

The mouse model used in the MC simulation.

Table 1

The optical parameters of the mouse model.

$μ_{a}$	$μ_{s}$	$n$	$G$
$0.03 {mm}^{- 1}$	$1.00 {mm}^{- 1}$	1.33	0.9

Table 2

Comparison of computation time of Monte Carlo simulations with different numbers of triangle meshes and the same number of photons (104).

Number of triangle meshes	Computation time of the conventional method (s)	Computation time of the proposed method (s)	Percentage improvement
4952	282	181	35.82%
5164	302	183	39.40%
6084	357	227	36.41%
13602	856	532	37.85%

Table 3

Comparison of computation time of Monte Carlo simulations with different number of photons and the same number of triangle meshes (6084).

Number of photons	Computation time of the conventional method (s)	Computation time of the proposed method (s)	Percentage improvement
$10^{4}$	357	227	36.41%
$10^{5}$	3913	2457	37.21%
$10^{6}$	39279	24274	38.20%
$5 \times 10^{6}$	221150	121489	45.06%

4. Conclusion

In this paper, we proposed a fast photon-boundary method in MC simulation of light propagation in complicated geometry. The proposed method uses the line segment instead of the ray in the conventional method to check if there exists a photon-boundary intersection. The MC algorithm using the proposed method was validated in a semi-infinite medium by comparing the results with those from the analytical solution of the DE and with those obtained from the MC simulation using the conventional method. Finally, the computation time of MC simulations for light propagation in complex 3-D geometry using the proposed method and the conventional method are compared, under the conditions of four different numbers of triangle meshes and four different sets of simulated photons. The results indicate that the MC simulation with the proposed method for intersection test is faster than that with the conventional method.

In our future work, the tetrahedral-based model introduced in both Refs. 22 and 23 will be considered for its higher efficiency in the intersection test. In addition, hardware acceleration based on graphics processing units (GPU)¹⁹^,²¹^,²⁷ is gradually applied in MC simulation. The proposed method may accelerate the MC simulation of photon propagation more significantly when implemented in the GPU. In the future, we will also try to apply the presented method in GPU-based acceleration for MC simulation.

Acknowledgments

This work is supported by the National Basic Research Program of China (973) under Grant No. 2011CB707701; the National Natural Science Foundation of China under Grant No. 81071191, 81271617; the National Major Scientific Instrument and Equipment Development Project under Grant No. 2011YQ030114; National Science and technology support program under Grant No. 2012BAI23B00.

References

1.

A. Ishimaru, Wave Propagation and Scattering in Random Media, Academic Press, New York (1978). Google Scholar

2.

R. C. Haskellet al., “Boundary conditions for the diffusion equation in radiative transfer,” J. Opt. Soc. Am. A, 11 (10), 2727 –2741 (1994). http://dx.doi.org/10.1364/JOSAA.11.002727 JOAOD6 0740-3232 Google Scholar

3.

B. C. WilsonG. Adam, “A Monte Carlo model for the absorption and flux distributions of light in tissue,” Med. Phys., 10 (6), 824 –830 (1983). http://dx.doi.org/10.1118/1.595361 MPHYA6 0094-2405 Google Scholar

4.

S. A. Prahlet al., “A Monte Carlo model of light propagation in tissue,” Dosimetry of Laser Radiation in Medicine and Biology, 102 –111 1989). Google Scholar

5.

M. Keijzeret al., “Light distributions in artery tissue: Monte Carlo simulations for finite-diameter laser beams,” Lasers Surg. Med., 9 (2), 148 –154 (1989). http://dx.doi.org/10.1002/(ISSN)1096-9101 LSMEDI 0196-8092 Google Scholar

6.

S. T. Flocket al., “Monte-Carlo modeling of light-propagation in highly scattering tissues—I: model predictions and comparison with diffusion theory,” IEEE Trans. Biomed. Eng., 36 (12), 1162 –1168 (1989). http://dx.doi.org/10.1109/TBME.1989.1173624 IEBEAX 0018-9294 Google Scholar

7.

S. T. FlockB. C. WilsonM. S. Patterson, “Monte Carlo modeling of light propagation in highly scattering tissues—II: comparison with measurements in phantoms,” IEEE Trans. Biomed. Eng., 36 (12), 1169 –1173 (1989). http://dx.doi.org/10.1109/10.42107 IEBEAX 0018-9294 Google Scholar

8.

M. KeijzerJ. W. PickeringM. J. C. van Gemert, “Laser beam diameter for port wine stain treatment,” Lasers Surg. Med., 11 (6), 601 –605 (1991). http://dx.doi.org/10.1002/(ISSN)1096-9101 LSMEDI 0196-8092 Google Scholar

9.

S. L. JacquesL. H. Wang, “Chapter 5: Monte Carlo modeling of light transport in tissues,” Optical Thermal Response of Laser Irradiated Tissue, Plenum Press, New York (1995). Google Scholar

10.

E. TinetS. AvrillierJ. M. Tualle, “Fast semianalytical Monte Carlo simulation for time-resolved light propagation in turbid media,” J. Opt. Soc. Am. A., 13 (9), 1903 –1915 (1996). http://dx.doi.org/10.1364/JOSAA.13.001903 1084-7529 Google Scholar

11.

N. S. ZolekA. LiebertR. Maniewski, “Optimization of the Monte Carlo code for modeling of photon migration in tissue,” Comput. Meth. Prog. Bio., 84 (1), 50 –57 (2006). http://dx.doi.org/10.1016/j.cmpb.2006.07.007 CMPBEK 0169-2607 Google Scholar

12.

D. A. Boaset al., “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express., 10 (3), 159 –170 (2002). http://dx.doi.org/10.1364/OE.10.000159 OPEXFF 1094-4087 Google Scholar

13.

T. Binzoniet al., “Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization,” Comput. Meth. Prog. Bio., 89 (1), 14 –23 (2008). http://dx.doi.org/10.1016/j.cmpb.2007.10.008 CMPBEK 0169-2607 Google Scholar

14.

E. Margallo-BalbasP. J. French, “Shape based Monte Carlo code for light transport in complex heterogeneous tissues,” Opt. Express., 15 (21), 14086 –14098 (2007). http://dx.doi.org/10.1364/OE.15.014086 OPEXFF 1094-4087 Google Scholar

15.

Q. LiuN. Ramanujam, “Scaling method for fast Monte Carlo simulation of diffuse reflectance spectra from multilayered turbid media,” J. Opt. Soc. Am. A., 24 (4), 1011 –1025 (2007). http://dx.doi.org/10.1364/JOSAA.24.001011 1084-7529 Google Scholar

16.

G. M. PalmerN. Ramanujam, “Monte Carlo-based inverse model for calculating tissue optical properties. Part I: theory and validation on synthetic phantoms,” Appl. Opt., 45 (5), 1062 –1071 (2006). http://dx.doi.org/10.1364/AO.45.001062 APOPAI 0003-6935 Google Scholar

17.

A. Liebertet al., “Monte Carlo algorithm for efficient simulation of time-resolved fluorescence in layered turbid media,” Opt. Express, 16 (17), 13188 –13202 (2008). http://dx.doi.org/10.1364/OE.16.013188 OPEXFF 1094-4087 Google Scholar

18.

L. H. WangS. L. JacquesL. Zheng, “MCML-Monte Carlo modeling of light transport in multi-layered tissues,” Comput. Meth. Prog. Bio., 47 131 –146 (1995). http://dx.doi.org/10.1016/0169-2607(95)01640-F CMPBEK 0169-2607 Google Scholar

19.

T. MöllerB. Trumbore, “Fast, minimum storage ray/triangle intersection,” J. Graph. Tool., 2 (1), 21 –28 (1997). http://dx.doi.org/10.1145/1198555.1198746 Google Scholar

20.

N. Renet al., “GPU-based Monte Carlo simulation for light propagation in complex heterogeneous tissues,” Opt. Express, 18 (7), 6811 –6823 (2010). http://dx.doi.org/10.1364/OE.18.006811 OPEXFF 1094-4087 Google Scholar

21.

Q. FangD. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing unit,” Opt. Express, 17 (22), 20178 –20190 (2009). http://dx.doi.org/10.1364/OE.17.020178 OPEXFF 1094-4087 Google Scholar

22.

Q. Fang, “Mesh-based Monte Carlo method using fast ray-tracing in Plucker coordinates,” Biomed. Opt. Express, 1 (1), 165 –175 (2010). http://dx.doi.org/10.1364/BOE.1.000165 BOEICL 2156-7085 Google Scholar

23.

H. ShenG. Wang, “A tetrahedron-based inhomogeneous Monte Carlo optical simulator,” Phy. Med. Bio., 55 (4), 947 –962 (2010). http://dx.doi.org/10.1088/0031-9155/55/4/003 PHMBA7 0031-9155 Google Scholar

24.

M. S. PattersonB. ChanceB. C. Wilson, “Time resolved reflectance and transmittance for the non-invasive measurement of tissue optical properties,” Appl. Opt., 28 (12), 2331 –2336 (1989). http://dx.doi.org/10.1364/AO.28.002331 APOPAI 0003-6935 Google Scholar

25.

T. J. FarrellM. S. PattersonB. Wilson, “A diffusion theory model of spatially resolved, steady state diffuse reflectance for the noninvasive determination of tissue optical properties in vivo,” Med. Phys., 19 (4), 879 –888 (1992). http://dx.doi.org/10.1118/1.596777 MPHYA6 0094-2405 Google Scholar

26.

A. KienleM. S. Patterson, “Improved solutions of the steady-state and the time-resolved diffusion nuations for reflectance from semi-infinite turbid medium,” J. Opt. Soc. Am., 14 (1), 246 –254 (1997). http://dx.doi.org/10.1364/JOSAA.14.000246 JOSAAH 0030-3941 Google Scholar

27.

E. AlerstamT. SvenssonS. Anderson-Engels, “Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration,” J. Biomed. Opt. Lett., 13 (6), 060504 (2008). http://dx.doi.org/10.1117/1.3041496 JBOPFO 1083-3668 Google Scholar

Biography

Xiaofen Zhao received her BS in biomedical engineering from Jilin University, Jilin, China, in 2010. She is currently a postgraduate student in the Biomedical Engineering Department at Tsinghua University, Beijing, China. Her research interest is in fluorescence molecular tomography for small animal imaging.

Hongyan Liu received her BS in biomedical engineering from Tsinghua University, Beijing, China, in 2011. Her research interest consists of fluorescence molecular tomography for small animal imaging.

Bin Zhang received his MS in mechanical and electronic engineering from University of Science and Technology of China in 2007. From 2007 to 2009, he was a research assistant in Shenzhen Institute of Advanced Integration Technology, Chinese Academy of Sciences. He is now a PhD candidate in the Department of Biomedical Engineering, Tsinghua University, Beijing, China. His research interest centers on fluorescence molecular tomography for small animals.

Fei Liu received her BS in biomedical engineering from Zhejiang University, Zhejiang, China, in 2008. She is currently a PhD candidate in the Biomedical Engineering Department, Tsinghua University, Beijing, China. Her research interest is in fluorescence molecular tomography for small animal imaging.

Jianwen Luo received his BS and PhD (with honors) in biomedical engineering from Tsinghua University, Beijing, China, in 2000 and 2005, respectively. He was a postdoctoral research scientist from 2005 to 2009, and an associate research scientist from 2009 to 2011 in the Department of Biomedical Engineering at Columbia University, New York. He became a professor at Tsinghua University in 2011, and was enrolled in the Thousand Young Talents Program in 2012. He serves as an Advisory Editorial Board member of the Journal of Ultrasound in Medicine, Associate Editor of Medical Physics, and Faculty Member of Faculty of 1000. His research interest is in biomedical imaging, including ultrasound and elasticity imaging.

Jing Bai received her MS and PhD from Drexel University, Philadelphia, Pennsylvania, in 1983 and 1985, respectively. From 1985 to 1987, she was a research associate and assistant professor at the Biomedical Engineering and Science Institute, Drexel University. In 1988, 1991, and 2000, she became an associate professor, professor, and Cheung Kong chair professor at Biomedical Engineering Department of Tsinghua University, Beijing, China. From 1997, she has been an associate editor for IEEE Transactions on Information Technology in Biomedicine. Her current research interests include mathematical modeling and simulation of the cardiovascular system, optimization of cardiac assist devices, medical ultrasound, telemedicine, home health care network and home monitoring devices, and infrared imaging. She has authored or coauthored 10 books and more than 300 journal papers.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 Unported License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation Download Citation

Xiaofen Zhao, Hongyan Liu, Bin Zhang, Fei Liu, Jianwen Luo, and Jing Bai "Fast photon-boundary intersection computation for Monte Carlo simulation of photon migration," Optical Engineering 52(1), 019001 (3 January 2013). https://doi.org/10.1117/1.OE.52.1.019001

Published: 3 January 2013

Access the abstract

JOURNAL ARTICLE
7 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

KEYWORDS

Monte Carlo methods

Photon transport

Computer simulations

Optical computing

Mouse models

Geometrical optics

Refractive index

1.

Introduction

2.

Method

2.1.

The Ray-Triangle-Based Method

Eq. (1)

Eq. (2)

Eq. (3)

Eq. (4)

Fig. 1

2.2.

Line-Segment Triangle-Based Method

Fig. 2

Eq. (5)

Eq. (6)

Fig. 3

Eq. (7)

Eq. (8)

Fig. 4

3.

Results and Discussion

3.1.

Validation of the Algorithm in a Semi-infinite Medium

3.1.1.

Comparison to the analytical solution of DE

Fig. 5

3.1.2.

Comparison to the results obtained from the conventional MC method

Fig. 6

3.2.

Comparisons of Computation Time in Complex 3-D Geometry

Fig. 7

Table 1

Table 2

Table 3

4.

Conclusion

Acknowledgments

References

Biography

Show All Keywords

Keywords/Phrases

Search In:

Publication Years