OSA's Digital Library

Biomedical Optics Express

Biomedical Optics Express

  • Editor: Joseph A. Izatt
  • Vol. 1, Iss. 1 — Aug. 2, 2010
  • pp: 165–175
« Show journal navigation

Mesh-based Monte Carlo method using fast ray-tracing in Plücker coordinates

Qianqian Fang  »View Author Affiliations


Biomedical Optics Express, Vol. 1, Issue 1, pp. 165-175 (2010)
http://dx.doi.org/10.1364/BOE.1.000165


View Full Text Article

Acrobat PDF (3747 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

We describe a fast mesh-based Monte Carlo (MC) photon migration algorithm for static and time-resolved imaging in 3D complex media. Compared with previous works using voxel-based media discretization, a mesh-based approach can be more accurate in modeling targets with curved boundaries or locally refined structures. We implement an efficient ray-tracing technique using Plücker Coordinates. The Barycentric coordinates computed from Plücker-formed ray-tracing enables us to use linear Lagrange basis functions to model both media properties and fluence distribution, leading to further improvement in accuracy. The Plücker-coordinate ray-polygon intersection test can be extended to hexahedral or high-order elements. Excellent agreement is found when comparing mesh-based MC with the analytical diffusion model and 3D voxel-based MC code in both homogeneous and heterogeneous cases. Realistic time-resolved imaging results are observed for a complex human brain anatomy using mesh-based MC. We also include multi-threading support in the software and will port it to a graphics processing unit platform in the near future.

© 2010 OSA

1. Introductions

In bio-optics applications, the Monte Carlo (MC) method is often used to model photon migration inside human tissue [1

1. B. C. Wilson and G. Adam, “A Monte Carlo model for the absorption and flux distributions of light in tissue,” Med. Phys. 10(6), 824–830 (1983). [CrossRef] [PubMed]

5

5. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

]. Because MC effectively solves the Radiative Transport Equation (RTE) via random sampling, it offers excellent accuracy when simulating photon propagation inside general complex media. In many cases, it was chosen as the gold standard when validating a new algorithm [6

6. J. Ripoll, M. Nieto-Vesperinas, S. R. Arridge, and H. Dehghani, “Boundary conditions for light propagation in diffusive media with nonscattering regions,” J. Opt. Soc. Am. A 17(9), 1671–1681 (2000). [CrossRef] [PubMed]

9

9. Y. Xu, Q. Zhang, and H. Jiang, “Optical image reconstruction of non-scattering and low scattering heterogeneities in turbid media based on the diffusion approximation model,” J. Opt. A, Pure Appl. Opt. 6(1), 29–35 (2004). [CrossRef]

] or approximation [10

10. A. H. Hielscher, R. E. Alcouffe, and R. L. Barbour, “Comparison of finite-difference transport and diffusion calculations for photon migration in homogeneous and heterogeneous tissues,” Phys. Med. Biol. 43(5), 1285–1302 (1998). [CrossRef] [PubMed]

12

12. G. Ma, J. F. Delorme, P. Gallant, and D. A. Boas, “Comparison of simplified Monte Carlo simulation and diffusion approximation for the fluorescence signal from phantoms with typical mouse tissue optical properties,” Appl. Opt. 46(10), 1686–1692 (2007). [CrossRef] [PubMed]

]. MC has several additional merits such as being easy-to-program and that it’s straightforward to parallelize [13

13. D. R. Kirkby and D. T. Delpy, “Parallel operation of Monte Carlo simulations on a diverse network of computers,” Phys. Med. Biol. 42(6), 1203–1208 (1997). [CrossRef] [PubMed]

,14

14. A. Colasanti, G. Guida, A. Kisslinger, R. Liuzzi, M. Quarto, P. Riccio, G. Roberti, and F. Villani, “Multiple processor version of a Monte Carlo code for photon transport in turbid media,” Comput. Phys. Commun. 132(1–2), 84–93 (2000). [CrossRef]

]. Compared with finite-element (FE) [15

15. S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20(2 Pt 1), 299–309 (1993). [CrossRef] [PubMed]

] diffusion solvers used in many diffuse optical imaging applications, MC produces more accurate solutions, especially when simulating low-scattering media where the diffusion approximation becomes invalid. In many other cases, MC attracts users by avoiding the complex pre-/post-processing steps (such as mesh-generation) required by FE or finite-volume (FV) [16

16. K. Ren, G. S. Abdoulaev, G. Bal, and A. H. Hielscher, “Algorithm for solving the equation of radiative transfer in the frequency domain,” Opt. Lett. 29(6), 578–580 (2004). [CrossRef] [PubMed]

] methods. In particular, the expense of generating high quality 3D meshes for complex media had been a non-trivial task due to the lack of appropriate software tools. Only in recent years, the development of image-based 3D mesh generation software [17

17. CGAL, Computational Geometry Algorithms Library, URL: http://www.cgal.org

,18

18. Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging 2009, 1142–1145 (2009).

] had started to make mesh generation from volumetric data practical and easily accessible for general data processing.

The major drawbacks of MC are the low computational efficiency and the lack of ability to accurately model curved boundaries [19

19. T. Binzoni, T. S. Leung, R. Giust, D. Rüfenacht, and A. H. Gandjbakhche, “Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization,” Comput. Methods Programs Biomed. 89(1), 14–23 (2008). [CrossRef] [PubMed]

]. A traditional MC simulation can easily take several hours or even over a day to obtain a solution with the desired accuracy [11

11. A. Custo, W. M. Wells 3rd, A. H. Barnett, E. M. Hillman, and D. A. Boas, “Effective scattering coefficient of the cerebral spinal fluid in adult head models for diffuse optical imaging,” Appl. Opt. 45(19), 4747–4755 (2006). [CrossRef] [PubMed]

,20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

]. In the last couple of years, the fast development of multi-core and many-core processors, especially the graphics processing units (GPU), has opened new possibilities for highly efficient MC simulations by exploring massively parallel computing. As a technology preview, Alerstam et al. [21

21. E. Alerstam, T. Svensson, and S. Andersson-Engels, “Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration,” J. Biomed. Opt. 13(6), 060504 (2008). [CrossRef] [PubMed]

] had shown a 1000x acceleration using GPU in a simple homogeneous domain. More realistic GPU-accelerated MC code for 3D heterogeneous media and time-resolved imaging was recently reported by Fang and Boas [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

], obtaining a better than 300x acceleration on a low-cost graphics card compared with single-threaded CPU-based MC simulations. With a modern graphics card, these GPU-based MC codes are capable of producing 3D solutions in less than a minute.

In addition to the dramatic acceleration in speed, significant effort has been made to improve the capability of modeling complex heterogeneous structures. The most widely used MC photon migration code developed by Wang and Jacques, MCML [5

5. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

], can handle layered media with a circular symmetry. Boas et al. [20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

] and Fang et al. [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

] had developed methods to use a voxelated space to represent arbitrarily complex media, showing great flexibilities in a range of applications [23

23. Y. Fukui, Y. Ajichi, and E. Okada, “Monte Carlo prediction of near-infrared light propagation in realistic adult and neonatal head models,” Appl. Opt. 42(16), 2881–2887 (2003). [CrossRef] [PubMed]

,11

11. A. Custo, W. M. Wells 3rd, A. H. Barnett, E. M. Hillman, and D. A. Boas, “Effective scattering coefficient of the cerebral spinal fluid in adult head models for diffuse optical imaging,” Appl. Opt. 45(19), 4747–4755 (2006). [CrossRef] [PubMed]

,12

12. G. Ma, J. F. Delorme, P. Gallant, and D. A. Boas, “Comparison of simplified Monte Carlo simulation and diffusion approximation for the fluorescence signal from phantoms with typical mouse tissue optical properties,” Appl. Opt. 46(10), 1686–1692 (2007). [CrossRef] [PubMed]

]. However, to model regions with curved boundaries with these tools, one has to increase grid density which requires more memory and computation. Margallo-Balbás and French [24

24. E. Margallo-Balbás and P. J. French, “Shape based Monte Carlo code for light transport in complex heterogeneous Tissues,” Opt. Express 15(21), 14086–14098 (2007). [CrossRef] [PubMed]

] as well as Ren et al. [25

25. N. Ren, J. Liang, X. Qu, J. Li, B. Lu, and J. Tian, “GPU-based Monte Carlo simulation for light propagation in complex heterogeneous tissues,” Opt. Express 18(7), 6811–6823 (2010). [CrossRef] [PubMed]

] have explored triangular-surface-based MC approaches which allow an improved approximation for domain interfaces. However, to test for ray-surface intersection, one has to scan a range of triangles based on a partitioning scheme, which introduces redundant calculations. In addition, the surface representation is unable to model continuously varying complex media. An algorithm that combines the strengths of MC and FE methods, i.e both accuracy for general media and flexibility for representing arbitrary domain shapes, is highly desired.

Inspired by the works from computer graphics [26

26. N. Platis and T. Theoharis, “Fast ray-tetrahedron intersection using Plücker coordinates,” Journal of Graphics Tools 8(4), 37–48 (2003).

], we have developed a mesh-based Monte Carlo method (referred to as MMCM hereafter) by making use of a fast ray-tracing algorithm with Plücker coordinate representation [27

27. H.S.M. Coxeter, Non-Euclidean geometry, Univ. Toronto Press pp. 88–90 (1965)

]. With this approach, one cannot only model domains with smooth boundaries, but also simulate continuously varying random media, or meshes with complex or mixed element types (tetrahedral, hexahedral, etc). We implemented MMCM in C and multi-threaded programming and developed a software package [28

28. Q. Fang, “Mesh-based Monte Carlo – the software”, URL: http://mcx.sourceforge.net/mmc/

], including the core MMCM simulation code, validation suite, mesh processing scripts and a toolbox to compute analytical diffusion solutions for a heterogeneous sphere. Combined with the open-source 3D mesh-generation toolbox developed previously [18

18. Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging 2009, 1142–1145 (2009).

], we expect this software to be an accurate and efficient option for MC simulations in biological tissues.

We also noticed an independent study and software (TIM-OS) recently published by Shen and Wang [29

29. H. Shen and G. Wang, “A tetrahedron-based inhomogeneous Monte Carlo optical simulator,” Phys. Med. Biol. 55(4), 947–962 (2010). [CrossRef] [PubMed]

] exploring a similar idea, i.e using tetrahedral meshes for MC simulations. The key differences between TIM-OS and MMCM are 1) MMCM uses Plücker coordinates while TIM-OS uses a Cartesian representation, 2) MMCM allows linear Lagrange basis functions to represent optical properties and fluence continuum while TIM-OS only supports piece-wise constant basis, 3) MMCM can handle more complicated element shapes such as high-order tetrahedral or hexahedral elements, and 4) MMCM can perform both static and time-resolved simulations. Additionally, we provide rigorous validations using analytical solutions and voxel-based MC for both homogeneous and heterogeneous cases to show the accuracy improvement using the mesh-based approach.

In the remaining portion of the paper, we first present the basic ray-polygon intersection test using Plücker coordinates and describe the algorithm work-flow. Then we discuss extensions of the algorithm such as using linear basis functions, supporting complex mesh elements and parallelization. In the Results section, we validate the algorithm using a homogeneous domain with and without a spherical inclusion, and compare the results with analytical diffusion solutions as well as those from a voxel-based MC simulation; then we show an example of time-resolved simulation from a complex human brain atlas. In the last section, we summarize our findings and discuss further expansion of the algorithm, particularly the massively parallel version running on GPU processors.

2. Methods

2.1 Ray-polygon intersection test using Plücker coordinates

A 3D ray (with a direction) can be uniquely determined by specifying two distinct 3D points, p = (x 1,y 1,z 1) and q = (x 2,y 2,z 2). As a result, the Cartesian form of a 3D ray can be typically written as
x  =p+(qp)t
(1)
where t is a scalar. In the Plücker coordinates [27

27. H.S.M. Coxeter, Non-Euclidean geometry, Univ. Toronto Press pp. 88–90 (1965)

], a 3D ray, R, can be expressed by a vector-pair R(d:m) where d represents displacement and m represents moment, defined as
d  =qpm=p×q
(2)
respectively. One of the attractive properties of the Plücker coordinates is its simplicity in a ray-polygon intersection test. In Fig. 1(a)
Fig. 1 (a) Diagram of Plücker-formed ray-polygon intersection test and (b) other commonly used polyhedral elements.
, we show a diagram where a ray enters (R 1), exits (R 2) or misses (R 3) a 3D triangle (ABC). We assume the edges of a 3D triangle (e 1(U e1:V e1), e 2(U e2:V e2) and e 3(U e3:V e3)) are oriented in counter-clock-wise order. A ray r(U r:V r) can be tested for intersection with the triangle by computing the permuted inner products w i as
wi=UrVei+VrUei
(3)
from which we can assert that

  • wi0 for all ir enters the triangle
  • wi0 for all ir exits the triangle (4)
  • wi=0 for all ir is coplanar with the triangle
  • else → r misses the triangle

p=iuipi.
(6)

For a tetrahedron that is known to intersect with a ray, the above test can quickly identify the indices of the faces where the ray enters and exits the tetrahedron; with a few additional floating-point operations, we can get the coordinates of the entry/exit points. Note that the above test is not limited to tetrahedrons. It can be easily extended to any convex polyhedral elements [31

31. C. A. Felippa, Advanced Finite Element Methods, URL: http://www.colorado.edu/engineering/cas/courses.d/AFEM.d/

]. A few examples of other commonly used polyhedral elements are shown in Fig. 1(b).

2.2 Ray-tracing based MC algorithm

In MMCM, the heterogeneous media is represented by a 3D volumetric mesh, which can be either a uniform grid (with cubic or cuboid elements) such as the one for voxel-based MC, or an unstructured mesh with polyhedral elements (tetrahedra, hexahedra etc). For a given element in the mesh, the indices of the neighboring elements through each shared face are pre-computed as the “face-neighbor list” and stored in the input files. A source position vector (s 0) and an incident directional vector (c 0) are specified by the user.

After initializing the mesh information, including nodes, elements and face-neighbor list, an initial intersection test is performed for all elements to identify the initial element (e 0) that encloses the source (e 0 can also be pre-computed). Starting from e 0, we simulate the photon propagation using the well-established MC simulation procedures [5

5. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

,20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

,22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. Briefly, we produce a random scattering length and scattering (azimuth and zenith) angles based on the exponential distribution and Henyey-Greenstein phase function [5

5. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

], respectively, and move each photon along the scattering trajectory recursively. For each element along the path, we perform the above ray-tracing test and calculate the distance (d) between the current photon position (p) to the exit point (p x) of the enclosing element along the propagation vector (c). If d is greater than the remaining scattering length, we will move the photon to the end of the scattering path, and generate a new set of scattering length and angles and repeat the above process. If d is less than the remaining scattering path, we move the photon to the intersection point, p x, and update the enclosing element ID (e) by looking up the face-neighbor list. Then we repeat the above ray-element test for the new element to further propagate the photon. Encountering an empty face-neighbor (denoted by a zero value) indicates a photon passing through a boundary face. We can either reflect the photon at the exiting face of the element based on the Snell's law, or terminate the photon and launch a new one from the source. The escaped photon energy is stored to an array associated with surface nodes. Because the domain boundary is represented by a surface mesh, MMCM can be more accurate than the boundary reflection scheme developed for the voxelated-space [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

].

2.3. Time-resolved photon migration

MMCM supports time-resolved photon migration using a similar approach as in [20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

] and [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. Briefly, we create multiple copies of photon weight array at all nodes for each time gate of the simulation. For each photon, we keep track of the elapsed time of the propagation and accumulate the photon weight loss to the array of the matched time gate. MMCM can produce static solutions when specifying a single time gate with sufficient length. To produce the volumetric fluence distribution, we apply a normalization scheme as in [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

] to get the temporal point-spread function, or TPSF.

2.4. Parallelization and optimization

The parallelization of MMCM is surprisingly straightforward. A simple multi-threaded implementation using OpenMP [32

32. OpenMP Architecture Review Board, OpenMP Application Program Interface, Version 3.0, URL: http://www.openmp.org/mp-documents/spec30.pdf (2008)

] is included in our software. Because the essential part of the calculation in MMCM involves inner and cross products of short vectors, for modern CPUs, this calculation is well suited for Streaming SIMD Extensions (SSE) instructions (inner product is supported in SSE4). Therefore in our implementation, vector operations are written in SSE assembly and can be enabled by compiling switches. These operations are expected to be significantly more efficient on graphics hardware [33

33. NVIDIA, CUDA Compute Unified Device Architecture - Programming Guide, Version 3.0 (2010)

]. Further discussion about the massively parallel version of MMCM for the GPU platform using CUDA [33

33. NVIDIA, CUDA Compute Unified Device Architecture - Programming Guide, Version 3.0 (2010)

] and OpenCL [34

34. Khronos Group, OpenCL 1.0 Specification, URL: http://www.khronos.org/opencl/ (2008)

] can be found in the last section.

3. Results

A multi-threaded implementation of MMCM algorithm for multi-core CPU was developed and released under an open-source license [28

28. Q. Fang, “Mesh-based Monte Carlo – the software”, URL: http://mcx.sourceforge.net/mmc/

]. In the current version of the software, we provide support for tetrahedral elements. The extensions to more general elements or mixed types will be added in the future. Two multi-threaded random number generators (RNG) are included: a thread-safe 48bit linear congruential RNG and a Logistic-lattice RNG [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. In this software package, we also provide a Matlab toolbox to compute the analytical solution for a sphere inside a semi-infinite medium or infinite slab based on diffusion theory.

In this section, we first assess the accuracy improvement for MMCM by comparing its output with the analytical solution from diffusion theory and the output from a voxel-based MC software, Monte Carlo eXtreme (MCX) [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. We also explore memory utilization and simulation speed with respect to different mesh configurations. Lastly, we show a realistic simulation using complex meshes generated from a human brain atlas. All the computations were performed on a computer running Ubuntu 9.10 (64bit) with an Intel Xeon E5520 (2.3G) Quad-core processor. The voxel-based MC simulator, MCX, ran on the same computer with an nVidia GTX 470 GPU. For all the cases, the FE meshes were produced using our previously published 3D mesh generation software “iso2mesh” [18

18. Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging 2009, 1142–1145 (2009).

] where part of the meshing tasks were performed by the Computational Geometry Algorithms Library (CGAL) [17

17. CGAL, Computational Geometry Algorithms Library, URL: http://www.cgal.org

] and tetgen [35

35. H. Si, and K. Gaertner, “Meshing piecewise linear complexes by constrained Delaunay tetrahedralizations,” In Proc. of the 14th Int. Meshing Roundtable, pp. 147–163, (2005)

].

In the first example, we validate the algorithm using a homogeneous medium in a 60x60x60mm cubic domain. The configuration is similar to that in [20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

] and [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]: the point source is located at (30.1,30.2,0) mm (the non-integer coordinates are used to avoid simulating photons near the edges/faces of this particular mesh) with an incident direction vector of (0,0,1); the medium has an absorption coefficient μ a = 0.005/mm, scattering coefficient μ s = 1/mm, refractive index n = 1 and anisotropy g = 0.01. We simulate 3x107 photons for a time window of 0 to 5 ns and a resolution of 0.1 ns. For MMCM, we generate a tetrahedral mesh by splitting each 2x2x2 cube in a 60x60x60 grid into 5 tetrahedra (i.e. a T5 mesh [36

36. P. Fleischmann, B. Haindl, R. Kosik, and S. Selberherr, “Investigation of a mesh criterion for three-dimensional finite element diffusion simulation,” Proc. Int. Conf. Simulation Semiconductor Processes Devices (SISPAD '99). 71–74 (1999)

]). The resulting mesh contains 29791 nodes and 135000 elements.

The temporal response (TPSF) recorded at location (30,14,10) mm is shown in Fig. 2(a)
Fig. 2 Comparisons between MMCM with voxel-based MC code, MCX, and the diffusion model for a homogeneous cubic domain: (a) the time-domain response (TPSF) at position (30,14,10) mm, and (b) CW fluence contour plots (with 10db spacing) in plane y = 30 mm.
. The continuous wave (CW) fluence profiles [20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

] along plane y = 30 mm, represented by contour lines with 10db spacing, are shown in Fig. 2(b). For comparisons, we also plot the analytical solution from the diffusion model [20

20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

] and the simulation results from MCX (running for 3x108 photons) in Fig. 2(a).

In the second example, we changed the background μ a to 0.002/mm and embed a sphere centered at (30,30,30) mm with a radius of 10 mm inside the domain. The optical properties for the sphere are μ a = 0.05/mm, μ s = 5/mm, n = 1.37 and g = 0.9. The source position is set to (30,30,0) mm with an incident vector (0,0,1). We also set n = 1.37 in the background medium to ignore reflection. Two sets of meshes are generated for this case: 1) a high-density uniform mesh and 2) a mesh with higher density around the sphere boundary and source. The node and element numbers are summarized in Table 1

Table 1. Memory utility and speed comparisons for MMCM and MCX for a heterogeneous simulation.

table-icon
View This Table
| View All Tables
. The cross-cut views of the 2 meshes are shown in Figs. 3(a) and (c)
Fig. 3 Validations of MMCM in heterogeneous media. We show the mesh cross-cut-views and the fluence contour plots (with 10db spacing) of two mesh configurations: (a)-(b) a high-density uniform mesh and (c)-(d) a mesh with higher density at the surface of the spherical inclusion and near source.
. We simulate 3x108 photons with MCX using a 60x60x60 voxel grid, and 3x107 photons for MMCM using the above meshes. By interpolating the MC solutions to a grid using a piecewise-linear basis functions, the contour plots at plane y = 30 mm are shown in Figs. 3(b) and (d) for the two mesh configurations, respectively. Meanwhile, we compute the analytical solution for a sphere inside an infinite-slab extended from the diffusion models in [37

37. D. A. Boas, M. A. O’Leary, B. Chance, and A. G. Yodh, “Scattering of diffuse photon density waves by spherical inhomogeneities within turbid media: analytic solution and applications,” Proc. Natl. Acad. Sci. U.S.A. 91(11), 4887–4891 (1994). [CrossRef] [PubMed]

] and overlap it in both figures. Because of the presence of the vertical boundaries, we only show the analytical solution near the source and inside the sphere for the comparison. The speed and memory utilities for MCX and MMCM using various number of CPU cores/threads are summarized in Table 1.

In the last example, we use MMCM to simulate photon migration in a 3D complex FE mesh generated from a segmented human brain atlas [38

38. D. L. Collins, A. P. Zijdenbos, V. Kollokian, J. G. Sled, N. J. Kabani, C. J. Holmes, and A. C. Evans, “Design and construction of a realistic digital brain phantom,” IEEE Trans. Med. Imaging 17(3), 463–468 (1998). [CrossRef] [PubMed]

,22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. Using the “iso2mesh” toolbox, we first extract four triangular surfaces for scalp, cerebro-spinal fluid (CSF), gray matter and white matter, respectively. Then we apply 4 iterations of the Laplacian + HC smoothing algorithm [39

39. J. Vollmer, R. Mencel, and H. Müller, ““Improved Laplacian smoothing of noisy surface meshes,” In Proc,” EuroGraphics 99, 131–138 (1999).

] to each surface, from which we generate a 3D volumetric FE mesh. We configure the mesh algorithm to generate a denser mesh near the top of the head or close to the complex cortex surface. This helps greatly in reducing overall memory demands of the computation. The total number of nodes and elements of the mesh are 69865 and 425224, respectively. The brain tissue optical properties are summarized in Table 2

Table 2. Optical parameters in various brain tissue types for the human brain atlas simulation. The properties for scalp/skull and cerebro-spinal fluid (CSF) are based on [11] and those for gray and white matters are based on [40] at 630 nm.

table-icon
View This Table
| View All Tables
and are chosen to be identical to those used in [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]. A point source is positioned at (75.7,67.0,168.2) mm pointing toward the center of the head.

In Figs. 4(a) and (b)
Fig. 4 Test of MMCM with a complex 3D brain atlas: (a) 3D-cut and (b) sagittal-cut views of the FE mesh; (c) the CW fluence and (d) time-resolved fluence (TPSF), both in log10-scale, extracted at plane x = 76 mm (Media 1). The tissue layers, from exterior to interior, are scalp/skull, cerebro-spinal fluid (CSF), gray matter and white matter, respectively. The boundaries of the tissue layers are overlapped in (c) and (d).
, we show the 3D-cut and the sagittal-cut views of the generated mesh. Again, a total of 3x107 photons are simulated for a time-window of 0-3 ns. The CW and the animated time-resolved solutions (in log-scale) extracted at plane x = 76 mm are shown in Figs. 4 (c) and (d), respectively. The computation for mesh generation is under 2 minutes and the MMCM simulation time is 40.5 minutes using 4 CPU cores.

4. Discussion

For simulations in the homogeneous medium, the results from MMCM match well with voxel-based MCX output in both temporal and spatial profiles. In Fig. 2(a), both MCX and MMCM solutions show lower fluence in the later temporal windows compared with the diffusion model. This artifact was caused by the domain truncation due to a finite boundary in the two simulations [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

] and disappeared when using a larger volume (not shown). In the heterogeneous case, the MMCM solutions inside the sphere and near the source for both meshes match accurately with the analytical diffusion model. In comparison, the solution from the voxel-space simulation deviates from the analytical/MMCM solutions inside the sphere, confirming that an inaccurate geometric representation of the target may negatively impact the accuracy of the solution [19

19. T. Binzoni, T. S. Leung, R. Giust, D. Rüfenacht, and A. H. Gandjbakhche, “Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization,” Comput. Methods Programs Biomed. 89(1), 14–23 (2008). [CrossRef] [PubMed]

].

Comparing Figs. 3(b) and 3(d), we find that the two MMCM solutions using different meshes are quite similar. From Table 1, their run-times are also similar. This is expected because the computational expense in MMCM is roughly proportional to the mesh density multiplied by the fluence distribution, not to the absolute numbers of elements or nodes. However, the essential advantage of the second mesh is its memory efficiency. From Table 1, we can see that the total memory required by the second mesh is only 1/3 of the first one. This may have a significant implication when running simulations on a memory-limited device, such as a GPU. Although MCX does not need extra memory to store the basic mesh information, for each additional time-gate, it needs to duplicate the entire grid. By comparison, the memory cost for each additional time-gate for MMCM is very small. For a device with 1GB memory (a typical size for a high-end graphics card), the above difference can translate to over 12000 time-gates using the coarser MMCM mesh compared to only 1200 time-gates using MCX.

Comparing speeds with different CPU cores, we find that MMCM offers excellent parallel acceleration, which is almost proportional to the number of CPU cores. It is not surprising that MCX is a few hundred times faster than MMCM, as we are comparing parallel implementations of different scales and architectures. However, it is reasonable to believe that by porting the MMCM ray-tracing code to CUDA or OpenCL, we can expect to see a comparable acceleration. In addition to running simulations with many parallel threads, the built-in short-vector operations in GPU hardware are highly efficient; this may further improve the speed of MMCM. On the other hand, each ray-tracing step in MMCM requires accessing mesh data (such as node coordinates and elements) from global memory, which may introduce overhead. Fortunately, by running a large number of threads, this latency can be efficiently “hidden” [33

33. NVIDIA, CUDA Compute Unified Device Architecture - Programming Guide, Version 3.0 (2010)

]. To further reduce memory access latency, we can also exploit the constant memory in the GPU and optimal ordering of the mesh nodes to improve data caching efficiency.

In the last example, the calculated CW and time-domain solutions generally agree with what we have observed in [22

22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

]: photons travel a longer distance in the low-scattering, low-absorption CSF layer, while attenuating rapidly toward the center of the brain. This indicates that both iso2mesh and MMCM are capable of handling real-world complex heterogeneous media.

In summary, we have described a mesh-based Monte Carlo algorithm and developed a free software package that combines the accuracy of MC and the flexibility of FE method for 3D photon migration modeling. We validate this code for both homogeneous and heterogeneous media and show that the algorithm is capable of producing more accurate solutions when simulating objects with curved boundaries.

Acknowledgement

QF would like to thank the funding support from Harvard Catalyst Pilot Grant. He also would like to thank David Boas for the helpful discussions on the analytical diffusion solutions of a sphere.

References and links

1.

B. C. Wilson and G. Adam, “A Monte Carlo model for the absorption and flux distributions of light in tissue,” Med. Phys. 10(6), 824–830 (1983). [CrossRef] [PubMed]

2.

P. van der Zee, and D. T. Delpy, “Simulation of the point spread function for light in tissue by a Monte Carlo method,” in Oxygen Transport to Tissue I. A. Silver and A. Silver, eds. (Plenum, New York,), Vol. IX. (1987)

3.

Y. Hasegawa, Y. Yamada, M. Tamura, and Y. Nomura, “Monte Carlo simulation of light transmission through living tissues,” Appl. Opt. 30(31), 4515–4520 (1991). [CrossRef]

4.

M. Hiraoka, M. Firbank, M. Essenpreis, M. Cope, S. R. Arridge, P. van der Zee, and D. T. Delpy, “A Monte Carlo investigation of optical pathlength in inhomogeneous tissue and its application to near-infrared spectroscopy,” Phys. Med. Biol. 38(12), 1859–1876 (1993). [CrossRef] [PubMed]

5.

L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

6.

J. Ripoll, M. Nieto-Vesperinas, S. R. Arridge, and H. Dehghani, “Boundary conditions for light propagation in diffusive media with nonscattering regions,” J. Opt. Soc. Am. A 17(9), 1671–1681 (2000). [CrossRef] [PubMed]

7.

S. R. Arridge, H. Dehghani, M. Schweiger, and E. Okada, “The finite element model for the propagation of light in scattering media: a direct method for domains with nonscattering regions,” Med. Phys. 27(1), 252–264 (2000). [CrossRef] [PubMed]

8.

F. Martelli, D. Contini, A. Taddeucci, and G. Zaccanti, “Photon migration through a turbid slab described by a model based on diffusion approximation. II. Comparison with Monte Carlo results,” Appl. Opt. 36(19), 4600–4612 (1997). [CrossRef] [PubMed]

9.

Y. Xu, Q. Zhang, and H. Jiang, “Optical image reconstruction of non-scattering and low scattering heterogeneities in turbid media based on the diffusion approximation model,” J. Opt. A, Pure Appl. Opt. 6(1), 29–35 (2004). [CrossRef]

10.

A. H. Hielscher, R. E. Alcouffe, and R. L. Barbour, “Comparison of finite-difference transport and diffusion calculations for photon migration in homogeneous and heterogeneous tissues,” Phys. Med. Biol. 43(5), 1285–1302 (1998). [CrossRef] [PubMed]

11.

A. Custo, W. M. Wells 3rd, A. H. Barnett, E. M. Hillman, and D. A. Boas, “Effective scattering coefficient of the cerebral spinal fluid in adult head models for diffuse optical imaging,” Appl. Opt. 45(19), 4747–4755 (2006). [CrossRef] [PubMed]

12.

G. Ma, J. F. Delorme, P. Gallant, and D. A. Boas, “Comparison of simplified Monte Carlo simulation and diffusion approximation for the fluorescence signal from phantoms with typical mouse tissue optical properties,” Appl. Opt. 46(10), 1686–1692 (2007). [CrossRef] [PubMed]

13.

D. R. Kirkby and D. T. Delpy, “Parallel operation of Monte Carlo simulations on a diverse network of computers,” Phys. Med. Biol. 42(6), 1203–1208 (1997). [CrossRef] [PubMed]

14.

A. Colasanti, G. Guida, A. Kisslinger, R. Liuzzi, M. Quarto, P. Riccio, G. Roberti, and F. Villani, “Multiple processor version of a Monte Carlo code for photon transport in turbid media,” Comput. Phys. Commun. 132(1–2), 84–93 (2000). [CrossRef]

15.

S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20(2 Pt 1), 299–309 (1993). [CrossRef] [PubMed]

16.

K. Ren, G. S. Abdoulaev, G. Bal, and A. H. Hielscher, “Algorithm for solving the equation of radiative transfer in the frequency domain,” Opt. Lett. 29(6), 578–580 (2004). [CrossRef] [PubMed]

17.

CGAL, Computational Geometry Algorithms Library, URL: http://www.cgal.org

18.

Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging 2009, 1142–1145 (2009).

19.

T. Binzoni, T. S. Leung, R. Giust, D. Rüfenacht, and A. H. Gandjbakhche, “Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization,” Comput. Methods Programs Biomed. 89(1), 14–23 (2008). [CrossRef] [PubMed]

20.

D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]

21.

E. Alerstam, T. Svensson, and S. Andersson-Engels, “Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration,” J. Biomed. Opt. 13(6), 060504 (2008). [CrossRef] [PubMed]

22.

Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]

23.

Y. Fukui, Y. Ajichi, and E. Okada, “Monte Carlo prediction of near-infrared light propagation in realistic adult and neonatal head models,” Appl. Opt. 42(16), 2881–2887 (2003). [CrossRef] [PubMed]

24.

E. Margallo-Balbás and P. J. French, “Shape based Monte Carlo code for light transport in complex heterogeneous Tissues,” Opt. Express 15(21), 14086–14098 (2007). [CrossRef] [PubMed]

25.

N. Ren, J. Liang, X. Qu, J. Li, B. Lu, and J. Tian, “GPU-based Monte Carlo simulation for light propagation in complex heterogeneous tissues,” Opt. Express 18(7), 6811–6823 (2010). [CrossRef] [PubMed]

26.

N. Platis and T. Theoharis, “Fast ray-tetrahedron intersection using Plücker coordinates,” Journal of Graphics Tools 8(4), 37–48 (2003).

27.

H.S.M. Coxeter, Non-Euclidean geometry, Univ. Toronto Press pp. 88–90 (1965)

28.

Q. Fang, “Mesh-based Monte Carlo – the software”, URL: http://mcx.sourceforge.net/mmc/

29.

H. Shen and G. Wang, “A tetrahedron-based inhomogeneous Monte Carlo optical simulator,” Phys. Med. Biol. 55(4), 947–962 (2010). [CrossRef] [PubMed]

30.

M. Meyer, H. Lee, A. Barr, and M. Desbrun, “Generalized Barycentric Coordinates on Irregular Polygons,” Journal of Graphics Tools 7(1), 13–22 (2002).

31.

C. A. Felippa, Advanced Finite Element Methods, URL: http://www.colorado.edu/engineering/cas/courses.d/AFEM.d/

32.

OpenMP Architecture Review Board, OpenMP Application Program Interface, Version 3.0, URL: http://www.openmp.org/mp-documents/spec30.pdf (2008)

33.

NVIDIA, CUDA Compute Unified Device Architecture - Programming Guide, Version 3.0 (2010)

34.

Khronos Group, OpenCL 1.0 Specification, URL: http://www.khronos.org/opencl/ (2008)

35.

H. Si, and K. Gaertner, “Meshing piecewise linear complexes by constrained Delaunay tetrahedralizations,” In Proc. of the 14th Int. Meshing Roundtable, pp. 147–163, (2005)

36.

P. Fleischmann, B. Haindl, R. Kosik, and S. Selberherr, “Investigation of a mesh criterion for three-dimensional finite element diffusion simulation,” Proc. Int. Conf. Simulation Semiconductor Processes Devices (SISPAD '99). 71–74 (1999)

37.

D. A. Boas, M. A. O’Leary, B. Chance, and A. G. Yodh, “Scattering of diffuse photon density waves by spherical inhomogeneities within turbid media: analytic solution and applications,” Proc. Natl. Acad. Sci. U.S.A. 91(11), 4887–4891 (1994). [CrossRef] [PubMed]

38.

D. L. Collins, A. P. Zijdenbos, V. Kollokian, J. G. Sled, N. J. Kabani, C. J. Holmes, and A. C. Evans, “Design and construction of a realistic digital brain phantom,” IEEE Trans. Med. Imaging 17(3), 463–468 (1998). [CrossRef] [PubMed]

39.

J. Vollmer, R. Mencel, and H. Müller, ““Improved Laplacian smoothing of noisy surface meshes,” In Proc,” EuroGraphics 99, 131–138 (1999).

40.

A. N. Yaroslavsky, P. C. Schulze, I. V. Yaroslavsky, R. Schober, F. Ulrich, and H.-J. Schwarzmaier, “Optical properties of selected native and coagulated human brain tissues in vitro in the visible and near infrared spectral range,” Phys. Med. Biol. 47(12), 2059–2073 (2002). [CrossRef] [PubMed]

OCIS Codes
(170.3660) Medical optics and biotechnology : Light propagation in tissues
(170.5280) Medical optics and biotechnology : Photon migration
(170.7050) Medical optics and biotechnology : Turbid media

ToC Category:
Optics of Tissue and Turbid Media

History
Original Manuscript: June 8, 2010
Revised Manuscript: July 11, 2010
Manuscript Accepted: July 13, 2010
Published: July 15, 2010

Citation
Qianqian Fang, "Mesh-based Monte Carlo method using fast ray-tracing in Plücker coordinates," Biomed. Opt. Express 1, 165-175 (2010)
http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-1-1-165


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. B. C. Wilson and G. Adam, “A Monte Carlo model for the absorption and flux distributions of light in tissue,” Med. Phys. 10(6), 824–830 (1983). [CrossRef] [PubMed]
  2. P. van der Zee, and D. T. Delpy, “Simulation of the point spread function for light in tissue by a Monte Carlo method,” in Oxygen Transport to Tissue I. A. Silver and A. Silver, eds. (Plenum, New York,), Vol. IX. (1987)
  3. Y. Hasegawa, Y. Yamada, M. Tamura, and Y. Nomura, “Monte Carlo simulation of light transmission through living tissues,” Appl. Opt. 30(31), 4515–4520 (1991). [CrossRef]
  4. M. Hiraoka, M. Firbank, M. Essenpreis, M. Cope, S. R. Arridge, P. van der Zee, and D. T. Delpy, “A Monte Carlo investigation of optical pathlength in inhomogeneous tissue and its application to near-infrared spectroscopy,” Phys. Med. Biol. 38(12), 1859–1876 (1993). [CrossRef] [PubMed]
  5. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “MCML - Monte Carlo modeling of light transport in multilayered tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]
  6. J. Ripoll, M. Nieto-Vesperinas, S. R. Arridge, and H. Dehghani, “Boundary conditions for light propagation in diffusive media with nonscattering regions,” J. Opt. Soc. Am. A 17(9), 1671–1681 (2000). [CrossRef] [PubMed]
  7. S. R. Arridge, H. Dehghani, M. Schweiger, and E. Okada, “The finite element model for the propagation of light in scattering media: a direct method for domains with nonscattering regions,” Med. Phys. 27(1), 252–264 (2000). [CrossRef] [PubMed]
  8. F. Martelli, D. Contini, A. Taddeucci, and G. Zaccanti, “Photon migration through a turbid slab described by a model based on diffusion approximation. II. Comparison with Monte Carlo results,” Appl. Opt. 36(19), 4600–4612 (1997). [CrossRef] [PubMed]
  9. Y. Xu, Q. Zhang, and H. Jiang, “Optical image reconstruction of non-scattering and low scattering heterogeneities in turbid media based on the diffusion approximation model,” J. Opt. A, Pure Appl. Opt. 6(1), 29–35 (2004). [CrossRef]
  10. A. H. Hielscher, R. E. Alcouffe, and R. L. Barbour, “Comparison of finite-difference transport and diffusion calculations for photon migration in homogeneous and heterogeneous tissues,” Phys. Med. Biol. 43(5), 1285–1302 (1998). [CrossRef] [PubMed]
  11. A. Custo, W. M. Wells, A. H. Barnett, E. M. Hillman, and D. A. Boas, “Effective scattering coefficient of the cerebral spinal fluid in adult head models for diffuse optical imaging,” Appl. Opt. 45(19), 4747–4755 (2006). [CrossRef] [PubMed]
  12. G. Ma, J. F. Delorme, P. Gallant, and D. A. Boas, “Comparison of simplified Monte Carlo simulation and diffusion approximation for the fluorescence signal from phantoms with typical mouse tissue optical properties,” Appl. Opt. 46(10), 1686–1692 (2007). [CrossRef] [PubMed]
  13. D. R. Kirkby and D. T. Delpy, “Parallel operation of Monte Carlo simulations on a diverse network of computers,” Phys. Med. Biol. 42(6), 1203–1208 (1997). [CrossRef] [PubMed]
  14. A. Colasanti, G. Guida, A. Kisslinger, R. Liuzzi, M. Quarto, P. Riccio, G. Roberti, and F. Villani, “Multiple processor version of a Monte Carlo code for photon transport in turbid media,” Comput. Phys. Commun. 132(1–2), 84–93 (2000). [CrossRef]
  15. S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20(2 Pt 1), 299–309 (1993). [CrossRef] [PubMed]
  16. K. Ren, G. S. Abdoulaev, G. Bal, and A. H. Hielscher, “Algorithm for solving the equation of radiative transfer in the frequency domain,” Opt. Lett. 29(6), 578–580 (2004). [CrossRef] [PubMed]
  17. CGAL, Computational Geometry Algorithms Library, URL: http://www.cgal.org
  18. Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging 2009, 1142–1145 (2009).
  19. T. Binzoni, T. S. Leung, R. Giust, D. Rüfenacht, and A. H. Gandjbakhche, “Light transport in tissue by 3D Monte Carlo: influence of boundary voxelization,” Comput. Methods Programs Biomed. 89(1), 14–23 (2008). [CrossRef] [PubMed]
  20. D. A. Boas, J. P. Culver, J. J. Stott, and A. K. Dunn, “Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head,” Opt. Express 10(3), 159–170 (2002). [PubMed]
  21. E. Alerstam, T. Svensson, and S. Andersson-Engels, “Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration,” J. Biomed. Opt. 13(6), 060504 (2008). [CrossRef] [PubMed]
  22. Q. Fang and D. A. Boas, “Monte Carlo simulation of photon migration in 3D turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef] [PubMed]
  23. Y. Fukui, Y. Ajichi, and E. Okada, “Monte Carlo prediction of near-infrared light propagation in realistic adult and neonatal head models,” Appl. Opt. 42(16), 2881–2887 (2003). [CrossRef] [PubMed]
  24. E. Margallo-Balbás and P. J. French, “Shape based Monte Carlo code for light transport in complex heterogeneous Tissues,” Opt. Express 15(21), 14086–14098 (2007). [CrossRef] [PubMed]
  25. N. Ren, J. Liang, X. Qu, J. Li, B. Lu, and J. Tian, “GPU-based Monte Carlo simulation for light propagation in complex heterogeneous tissues,” Opt. Express 18(7), 6811–6823 (2010). [CrossRef] [PubMed]
  26. N. Platis and T. Theoharis, “Fast ray-tetrahedron intersection using Plücker coordinates,” Journal of Graphics Tools 8(4), 37–48 (2003).
  27. H.S.M. Coxeter, Non-Euclidean geometry, Univ. Toronto Press pp. 88–90 (1965)
  28. Q. Fang, “Mesh-based Monte Carlo – the software”, URL: http://mcx.sourceforge.net/mmc/
  29. H. Shen and G. Wang, “A tetrahedron-based inhomogeneous Monte Carlo optical simulator,” Phys. Med. Biol. 55(4), 947–962 (2010). [CrossRef] [PubMed]
  30. M. Meyer, H. Lee, A. Barr, and M. Desbrun, “Generalized Barycentric Coordinates on Irregular Polygons,” Journal of Graphics Tools 7(1), 13–22 (2002).
  31. C. A. Felippa, Advanced Finite Element Methods, URL: http://www.colorado.edu/engineering/cas/courses.d/AFEM.d/
  32. OpenMP Architecture Review Board, OpenMP Application Program Interface, Version 3.0, URL: http://www.openmp.org/mp-documents/spec30.pdf (2008)
  33. NVIDIA, CUDA Compute Unified Device Architecture - Programming Guide, Version 3.0 (2010)
  34. Khronos Group, OpenCL 1.0 Specification, URL: http://www.khronos.org/opencl/ (2008)
  35. H. Si, and K. Gaertner, “Meshing piecewise linear complexes by constrained Delaunay tetrahedralizations,” In Proc. of the 14th Int. Meshing Roundtable, pp. 147–163, (2005)
  36. P. Fleischmann, B. Haindl, R. Kosik, and S. Selberherr, “Investigation of a mesh criterion for three-dimensional finite element diffusion simulation,” Proc. Int. Conf. Simulation Semiconductor Processes Devices (SISPAD '99). 71–74 (1999)
  37. D. A. Boas, M. A. O’Leary, B. Chance, and A. G. Yodh, “Scattering of diffuse photon density waves by spherical inhomogeneities within turbid media: analytic solution and applications,” Proc. Natl. Acad. Sci. U.S.A. 91(11), 4887–4891 (1994). [CrossRef] [PubMed]
  38. D. L. Collins, A. P. Zijdenbos, V. Kollokian, J. G. Sled, N. J. Kabani, C. J. Holmes, and A. C. Evans, “Design and construction of a realistic digital brain phantom,” IEEE Trans. Med. Imaging 17(3), 463–468 (1998). [CrossRef] [PubMed]
  39. J. Vollmer, R. Mencel, and H. Müller, ““Improved Laplacian smoothing of noisy surface meshes,” In Proc,” EuroGraphics 99, 131–138 (1999).
  40. A. N. Yaroslavsky, P. C. Schulze, I. V. Yaroslavsky, R. Schober, F. Ulrich, and H.-J. Schwarzmaier, “Optical properties of selected native and coagulated human brain tissues in vitro in the visible and near infrared spectral range,” Phys. Med. Biol. 47(12), 2059–2073 (2002). [CrossRef] [PubMed]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Figures

Fig. 1 Fig. 2 Fig. 3
 
Fig. 4
 

Multimedia

Multimedia FilesRecommended Software
» Media 1: AVI (1689 KB)      QuickTime

« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited