Optical tomography using the SCIRun problem solving environment : Preliminary results for three-dimensional geometries and parallel processing

We present a 3D implementation of the UCL imaging package for absorption and scatter reconstruction from time-resolved data (TOAST), embedded in the SCIRun interactive simulation and visualization package developed at the University of Utah. SCIRun is a scientific programming environment that allows the interactive construction, debugging, and steering of large-scale scientific computations. While the capabilities of SCIRun’s interactive approach are not yet fully exploited in the current TOAST implementation, an immediate benefit of the combined TOAST/SCIRun package is the availability of optimized parallel finite element forward solvers, and the use of SCIRun’s existing 3D visualisation tools. A reconstruction of a segmented 3D head model is used as an example for demonstrating the capability of TOAST/SCIRun of simulating anatomically shaped meshes. c ©1999 Optical Society of America OCIS codes: (170.3010) Image reconstruction techniques; (100.6950) Tomographic image processing; (100.6890) Three-dimensional image processing References and links 1. S. R. Arridge, M. Schweiger and D. T. Delpy, “Iterative reconstruction of near-infrared absorption images,” In Inverse Problems in Scattering and Imaging, M. A. Fiddy, ed., Proc. SPIE 1767, 372–383 (1992). 2. M. Schweiger, S. R. Arridge and D. T. Delpy, “Application of the finite-element method for the forward and inverse models in optical tomography,” J. Math. Imag. Vision 3, 263–283 (1993). 3. K. D. Paulsen and H. Jiang, “Spatially-varying optical property reconstruction using a finite element diffusion equation approximation,” Med. Phys. 22, 691–701 (1995). 4. B. W. Pogue, M. S. Patterson, H. Jiang and K. D. Paulsen, “Initial assessment of a simple system for frequency domain diffuse optical tomography,” Phys. Med. Biol. 40, 1709–1729 (1995). 5. R. Model, M. Orlt, M. Walzel and R. Hünlich, “Reconstruction algorithm for near-infrared imaging in turbid media by means of time-domain data,” Appl. Opt. 14, 313–324 (1997). 6. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Problems 15, R1–R53 (1999). 7. J. C. Hebden, S. R. Arridge and D. T. Delpy, “Optical imaging in medicine: I. Experimental techniques,” Phys. Med. Biol. 42, 825–840 (1997). 8. S. R. Arridge and J. C. Hebden, “Optical imaging in medicine: II. Modelling and reconstruction,” Phys. Med. Biol. 42, 841–853 (1997). 9. H. Jiang, K. D. Paulsen, U. L. Osterberg, B. W. Pogue and M. S. Patterson, “Simultaneous reconstruction of absorption and scattering profiles in turbid media from near-infrared frequencydomain data,” Opt. Lett. 20, 2128–2130 (1995). 10. S. R. Arridge and M. Schweiger, “Sensitivity to prior knowledge in optical tomographic reconstruction,” In Optical Tomography, Photon Migration, and Spectroscopy of Tissue and Model (C) 1999 OSA 12 April 1999 / Vol. 4, No. 8 / OPTICS EXPRESS 263 #9178 $15.00 US Received March 01, 1999; Revised April 09, 1999 Media: Theory, Human Studies, and Instrumentation, B. Chance and R. R. Alfano, eds., Proc SPIE 2389, 378–388 (1995). 11. J. C. Hebden, F. E. W. Schmidt, M. E. Fry, M. Schweiger, E. C. Hillman, D. T. Delpy and S. R. Arridge, “Simultaneous reconstruction of absorption and scattering images by multichannel measurement of purely temporal data,” Opt. Lett. 24 (1999). 12. M. Schweiger and S. R. Arridge, “The UCL optical tomography software system (TOAST)” http://www.medphys.ucl.ac.uk/toast/index.htm 13. M. Schweiger and S. R. Arridge, “Comparison of 2D and 3D reconstruction algorithms in optical tomography,” Appl. Opt. 37, 7419–7428 (1998). 14. S. G. Parker, D. M. Weinstein and C. R. Johnson, “The SCIRun computational steering software system,” In Modern Software Tools in Scientific Computing, edited by E. Arge, A. M. Bruaset and H. P. Langtangen, 1–40 (Birkhauser Press, 1997). 15. S. G. Parker and C. R. Johnson, “SCIRun: A scientific programming environment for computational steering,” Supercomputing ‘95 (1995) http://www.supercomp.org/sc95/proceedings/499 SPAR/SC95.HTM 16. J. Vetter and K. Schwan, “Models for computational steering,” in: Proceedings of the Third International Conference on Configurable Distributed Systems (1996) 17. S. G. Parker, M. Miller, C. D. Hansen and C. R. Johnson, “An integrated problem solving environment: The SCIRun computational steering system,” In Proceedings of the 31st Hawaii International Conference on System Sciences (HICSS-31) (IEEE Computer Society Press, 1998). 18. T. A. McCormick, T. A. DeFanti and M. D. Brown, “Visualisation in scientific computing,” Computer Graphics 21, v–ix, 1–14, A1–E8 (1987) 19. S. R. Arridge and M. Schweiger. “Photon measurement density functions. Part 2: Finite element calculations.” Appl. Opt. 34, 8026–8037, (1995). 20. M. Schweiger and S. R. Arridge, “Optimal data types in Optical Tomography”, In Information Processing in Medical Imaging, edited by J. Duncan and G. Gindi, 71–84 (Springer, New York,


Introduction
One of the approaches to image reconstruction in Optical Tomography is also one of the earliest proposed viz. to treat the problem as the minimization of an objective function representing the sum-squared difference between measured and modeled data [1][2][3][4][5].Under the diffusion approximation, the solution sought is μa , μ s = arg min µa,µ s ||y M − F M (µ a , µ s )|| + G(µ a , µ s ) (1) where F M is a function giving modelled data of type M, y M is measured data of this type, and G is a functional representing prior knowledge.For recent detailed reviews of Optical Tomography see [6][7][8].A suitable forward model is given by a Finite Element Method (FEM) wherein (µ a , µ s ) are represented in the basis of simple locally supported functions on a mesh.Successful reconstructions have been performed with this approach from simulated and simple phantom data [9][10][11].The UCL software system for Time-resolved Optical Absorption and Scattering Tomography (TOAST) is available online [12].The measurement types we use are normalized integral transforms of the time dependent y M , which includes for example Laplace transforms and temporal moments (variance, skew, etc.) Most reported results with FEM have been limited to two-dimensional problems, due to the high computational cost of 3D reconstructions and the limited availability of 3D mesh generation and visualization software.Recently we showed 3D reconstructions from simulated data of a cylindrical test case with small embedded absorption inhomogeneities [13].The extension to more general complex geometries does not pose any principal problems with regard to the finite element forward and inverse solvers.However the mesh generation and visualization of the resulting images becomes considerably more difficult.
To minimize the additional development effort required to extend TOAST to a general 3D reconstruction package, we are now integrating TOAST into SCIRun, a versatile and efficient 3D problem solving environment developed by the Center for Scientific Computing and Imaging at the University of Utah [14,15].The modular structure of SCIRun allows rapid incorporation of new applications, and its support for parallel platforms is attractive for delivering the performance required for 3D reconstruction in clinical optical tomography applications.In this paper we show an overview of the TOAST interface in the SCIRun framework and present preliminary 3D reconstructions of a head model obtained with the system, together with a discussion of a parallel processing implementation.

The SCIRun programming environment
SCIRun is a computational steering problem solving environment [16] that allows the interactive construction, debugging, and steering of large-scale scientific computations [14,17].A scientific application is constructed by connecting computational elements (modules) to form a program (network).This program may contain several computational elements as well as several visualization elements, all of which work together in orchestrating a solution to a scientific problem.SCIRun has been designed to allow the interactive modification of geometric inputs and computational parameters, so that the results of these changes provide immediate feedback to the investigator [18].SCIRun data flow example for loading and visualization of reconstruction results as scalar fields Among the set of tools are a visual programming language, visual interactive manipulation devices, general class libraries, domain-specific component libraries, a convenient development environment, and an optimized runtime system.SCIRun composes programs from different computational algorithms (modules) using a dataflow style "boxes and wires" approach.An example of a dataflow network used for a 2D reconstruction with simultaneous visualization is shown in Fig. 1.The SCIRun implementation of TOAST currently comprises the following modules: • TOASTReadMesh: Read a mesh from an input stream and provide it on an output port.
• TOASTReadData: Read data from an input stream and provide it on an output port.
• TOASTGenData: Given mesh and QM description, generate data on the fly.
• TOASTAddNoise: Add a specified amount of noise to a data set.
• TOASTDataArray: Combine a series of data sets to a single composite data array which can be passed to the reconstruction module.• TOAST: The main reconstruction module.This accepts a mesh, QM specification and data array, and provides absorption and scatter reconstruction images on the output ports.Images are emitted for each iteration, so the progress of the reconstruction can be controlled during the program run.• TOASTResetParam: Reset the µ a and µ s parameter values stored in the mesh to an initial guess prior to reconstruction.
For each module, a user interface can be created via a Tcl/Tk script which allows interactive input of parameters relevant to the operation of the module.As an example, Fig. 2 shows the interface for the Toast module.The data flow between modules is symbolized by connections in the network diagram.User interface for the TOAST module.
Some TOAST specific data types were added, for example for source detector (QM) file specification (yellow) or measurement data (blue).The reconstructed images emitted by the Toast module are standard SCIRun Surface data, which allows the use of existing visualization modules.
While the full potential of SCIRun's interactive approach is not yet exploited in the TOAST implementation, an immediate benefit of the incorporation of TOAST into SCIRun is the availability of optimized parallel finite element forward solvers, and the use of 3D visualization tools already present in the standard SCIRun environment.

3D reconstructions in a head model
The main motivation for implementing TOAST in the SCIRun framework is the availability of 3D mesh generation and visualization tools, together with fast parallel matrix solvers which provides a significant reconstruction speed improvement on an appropriate hardware.
As a test case we consider a head mesh composed of a regular grid of tetrahedra.The mesh is constructed from a segmented 3D MRI data set and scaled to the size of a newborn baby.Available meshing tools allow the generation of meshes of various resolution.The left image in Fig. 3 shows the a rendering of the highest-resolution mesh used in the simulations presented here.The dimensions of the bounding box of the mesh are 60 × 80 × 88 mm.Each mesh element is labelled according to the physiological region to which it is assigned.The right image in Fig. 3 shows a cross-section through a lowerresolution mesh, where the regions are identified by color-coding.Optical properties assigned are similar to those previously reported [19].
The TOAST tools for data generation and reconstruction use SCIRun's parallel conjugate gradient matrix solver with diagonal preconditioning.Further performance improvements will be achieved with a parallel implementation of an incomplete Cholesky preconditioner which is currently being developed.Left: high-resolution head mesh (98142 nodes, 452243 elements) used for forward data calculation; right: cross section of a low-resolution mesh showing the 4 regions skin and muscle , bone , grey matter and white matter .
To demonstrate the effect of the parallel solver we computed forward data of steady-state intensity (log E) and the third central moment (skew) [20] using a medium resolution head mesh (145776 tetrahedra) for a data set of 24 sources and 24 detectors, arranged in 3 horizontal rings on the mesh surface, all above eye level.The calculations were executed on 1-8 processors of a 14-processor SGI Onyx VTX platform (194MHz MIPS R10000) with 4 GB of main memory, and on a dual Pentium-II 450 with 1 GB of memory running Linux.The timing results are shown in Fig. 4.
Note that the effect of parallelization and preconditioning is more evident for the skew calculation because this requires 4 times as many CG solution runs as the intensity calculation, thus reducing the influence of the setup and system matrix collection steps which are of comparable computational cost for both cases.
The relative performance gains for multiprocessor systems shown for the forward problem directly translate to the reconstruction problem, because the repeated solution of the forward problem is the single most time-consuming operation in the inverse solver.For orientation, the head mesh is rendered as a point cloud.
Figure 6 shows a simple test case of reconstruction of absorbing and scattering perturbations inside the brain after 40 iterations of a nonlinear conjugate gradient method.The initial guess was the background tissue parameters without the perturbation.Target perturbations were spheres of 5 mm radius embedded in the grey matter region.Two of the objects were absorption perturbations at 3 and 4 times the background value, and one was a scatter perturbation at 8 times the background value.Isosurfaces of PMDF for skew data.Green: absorption PMDF (isosurface value -2), blue: diffusion PMDF (isosurface value -0.23).

Figure 1 .
Figure 1.SCIRun data flow example for loading and visualization of reconstruction results as scalar fields

Figure 2 .
Figure 2.User interface for the TOAST module.

Figure 3 .
Figure 3.Left: high-resolution head mesh (98142 nodes, 452243 elements) used for forward data calculation; right: cross section of a low-resolution mesh showing the 4 regions skin and muscle , bone , grey matter and white matter .

Figure 4 .
Figure 4.Timings for forward data generation of intensity (left) and skew data (right).Tested platforms are: SGI Onyx VTX (1-8 processors) and dual Pentium-II 450 (1 and 2 processors).Each bar is divided into setup and system matrix assembly time (bottom) and conjugate gradient solver time (top).

Figure 5
Figure5shows the Photon Measurement Density Function (PMDF)[19] of the skew of the temporal response function, for a given source-detector pair.The PMDF indicates the region of sensitivity of a measurement to a change in the optical properties.For orientation, the head mesh is rendered as a point cloud.Figure6shows a simple test case of reconstruction of absorbing and scattering perturbations inside the brain after 40 iterations of a nonlinear conjugate gradient method.The initial guess was the background tissue parameters without the perturbation.Target perturbations were spheres of 5 mm radius embedded in the grey matter region.Two of the objects were absorption perturbations at 3 and 4 times the background value, and one was a scatter perturbation at 8 times the background value.