**Detection of Abnormalities in a Biological Tissue by Diffuse Optical Tomography: A Computational Study**

H. Trabelsi, M. Gantri and E. Sediki *University of Tunis El Manar, Faculté des Sciences de Tunis Unité de Rayonnement, Département de Physique, Tunis Tunisia* 

## **1. Introduction**

26 Medical Imaging

Tu, Z. & Bai, X. (2010). Auto-Context and Its Application to High-Level Vision Tasks and 3D

Ugarriza, L.G., Saber, E., Vantaram, S.R., Amuso, V., Shaw, M. & Bhaskar, R. (2009).

van Rikxoort, E.M., Isgum, I., Arzhaeva, Y., Staring, M., Klein, S., Viergever, M.A., Pluim, J.

Weickert, J. (1997). Recursive separable schemes for nonlinear diffusion filters, In: *Scale-space* 

Withey, D.J. & Koles, Z.J. (2007). Medical Image Segmentation: Methods and Software.

*International Conference on Functional Biomedical Imaging* (2007), pp. 140-143. Witkin, A.P. (1983). Scale space filtering. *Proc. Int "l Joint Conference Artificial Intelligence*, pp.

Xie, J.; Jiang, Y. & Tsui, Y. (2005). Segmentation of Kidney From Ultrasound Images Based

Yi, Z.; Criminisi, A.; Shotton, J. & Blake, A. (2009). Discriminative, Semantic Segmentation of

Yoon, S.W.; Lee, C.; Kim, J.K. & Lee, M. (2008). Wavelet-based Multi-resolution Deformation

Yu, G.; Lin, P.; Cai, S. (2008). Robust Vessel Segmentation Based on Multi-resolution Fuzzy Clustering. *Lecture Notes in Computer Science*, vol. 5326/2008, pp. 338-345. Yuille, A.L. , Poggio, T.A. (1986). Scaling theorems for zero-crossings, *IEEE PAMI*, vol. 8

Zaidi, H. & El Naqa, I. (2010). PET-guided delineation of radiation therapy treatment

Zhang, N., Ruan, S., Lebonvallet, S., Liao, Q. & Zhu, Y. (2011). Kernel feature selection to

Zhao, Z., Zhang, H., Wahle, A., Thomas, M.T., Stolpen, A.H., Scholz, T.D. & Sonka, M.

quantitative analysis. *Medical Image Analysis*, vol. 13 (2009), pp. 483–493. Zheng,Y., Vega-Higuera,F.V., Zhou,S.K. & Comaniciu, D. (2010). Fast and Automatic Heart

*Medicine and Molecular Imaging*, vol. 37 (2010), pp. 2165–2187.

*Image Understanding*, vol. 115 (2011), pp. 256–269.

*Science*, vol. 6357/2010, pp. 84-91

Viergever, M. pp. 260-271, Springer-Verlag, Berlin Heidelberg.

1019-1022, Karslruhe, Germany, Aug. 1983.

Merging. *IEEE Trans on Imag Proc*, vol. 18, no 10 (2009), pp. 2275-2288. van Ginneken, B., Romeny, B.M.T.H. & Viergever, M.A. (2001). Computer aided diagnosis in

(2010), pp. 1744-1757.

(2010), pp. 39–49.

1241.

558-565.

(2008), pp. 207–214.

(1986), pp. 15–25.

Brain Image Segmentation. *IEEE Trans on Pat Anal and Mach Intel*, vol. 3, no 10

Automatic Image Segmentation by Dynamic Region Growth and Multiresolution

chest radiography: a survey. *IEEE Trans on Med Imag*, vol. 20, no 12 (2005), pp. 1228-

P.W. & van Ginneken, B. (2010). Adaptive local multi-atlas segmentation: Application to the heart and the caudate nucleus. *Medical Image Analysis*, vol. 14

*theory in computer vision 1252.* Romeny, B.T.H., Florack, L., Koenderink, J.,

Proceeding of Noninvasive Functional Source Imaging of the Brain and Heart.

on Texture and Shape Priors. *IEEE Trans on Med Imag*, vol. 24, no 1 (2005), pp. 45-57.

Brain Tno in MR Images. *Lecture Notes in Computer Science*, Volume 5762 (2009), pp.

for Medical Endoscopic Image Segmentation. *Journal of Medical Systems*, vol. 32

volumes: a survey of image segmentation techniques. *European Journal of Nuclear* 

fuse multi-spectral MRI images for brain tumor segmentation. *Computer Vision and* 

(2009). Congenital aortic disease: 4D magnetic resonance segmentation and

Isolation in 3D CT Volumes: Optimal Shape Initialization. *Lecture Notes in Computer* 

Diffuse optical tomography (DOT) has a real potential as a new medical diagnostic tool. This fact has stimulated considerable research interest in the last years. When it is possible, this technique has some advantages not available in classical medical imaging methods. Among the benefits of this technique, nonionizing radiation, low-cost instrumentation, mobility and possible functional imaging are the most important. Also, DOT could be used in combination with other imaging methods to provide high-resolution structural information. The transparency of tissues is related to absorption and scattering of light and it reaches its maximum in the near-infrared. Interaction of radiation from this range with biological tissues is widely studied [1-5]. In particular, some works was interested to detection of tumor-like inclusions within a biological tissue [6-8].

Most of models of light interaction with biological matter are based on radiative transfer theory [9-13]. The most used parameters in modelling Laser radiation interaction with biological tissue is absorption and scattering. However some other studies evoked a significant variation of refractive index of abnormal biological tissues especially in the near infrared range. More precisely, experimental results [14-17] showed that the tissue of malignant tumors could manifest an increase of the refractive index which could attain until 10 % of that of a normal tissue which encircles them. So, medical imaging by diffuse optical tomography should take advantage from the emergence of a third contrast parameter which is the refractive index. This drove to the appearance of a big number of numerical and fundamental works in the field of radiative transfer in a varying refractive index medium. Among them some papers are interested to varying refractive index biological tissues [18-22].

In this chapter, infrared laser radiation interaction with a biological tissue including a spatially variation of refractive index is studied through a time-dependent radiative transfer model. Our general concern is to contribute at the usability of the radiative transfer theory in a potential optical tomography setting in medical imaging. At this level, studying the effect of refractive index on the transmitted light through a biological rectangular layer is crucial. This should be a contribution to obtaining information about

Detection of Abnormalities in a Biological

so equation (3) becomes

discrete directions,

cosines( , ), 1,... 

*m P*

,

*n I*

*c t*

phase function is written as

*M*

*nx ny*

*m mm*

obtained by using interpolation formula

' 1, '

*s m mm m p*

*wp I*

*mm*

'

1 1 (

' ' ',

The angular redistribution terms will be noted

Tissue by Diffuse Optical Tomography: A Computational Study 29

*n n n I nI <sup>I</sup> <sup>I</sup>*

12 1 . (. ) [(1 ) ] [(1 ) ] <sup>2</sup>

*D I x y andD I* 2 2 [(1 ) ]

,

2 12 1 . (. ) <sup>2</sup> *<sup>x</sup> <sup>y</sup> n n n I nI D D n n xy n*

*m m m M* . An orientation depending on the incident ray direction is

*m mE mW m mN mS a sm p*

( )( ) ()

,, , , ,

*mm mm*

 ))

' '

2

*y I I xI I x y I*

 

In our numerical implementation, we use a rectangular domain which is divided into a set of I x J elementary uniform volumes ΔV with a uniform unitary depth. The angular discretization is obtained through a discrete ordinate technique. This yields a set of M

adopted for each cell [12]. Calculations are done by using integration of equation (1) over an

where *Dx m*, and *Dy*,*<sup>m</sup>* are the discrete angular derivative terms at the angular ordinate *ξ<sup>m</sup>* and *ηm* respectively and *wm* is a weighting factor. The discrete term of Henyey-Greenstein

*x m y m m P s m mm m p*

,, , ,

*n n D D x yS w p I*

*<sup>g</sup> <sup>p</sup>*

4 (1 2 (

*g g*

2 3/2

 

1

If the direction cosines are positive, the directional radiances are known on the faces W and S and they are unknown on the faces E and N of the (i,j)-cell and also in the centre P. Therefore, we need two complementary relations to eliminate *m N*, and *m E*, , this can be

> *mP mE m W mP mN m S*

,

(1 ) (1 )

 

> 

*II I II I* ,, , ,, ,

where α is an interpolation parameter. Using these relations, equation (3) becomes

)

 

2 2

 

(3)

 [(1 ) ]

2 2

*<sup>m</sup>* , 1.. *m M* giving a set of angular discrete direction

  (4)

2 2

   

*nn x y n*

2

elementary volume ΔV for each discrete direction. This gives

abnormal tissue in a typical DOT scheme. However, it is important to note that in a varying refractive index medium, the rays are not straight lines but curves. So even in a rectangular geometry, the varying refractive index radiative transfer equation displays angular derivative terms [23-25]. In this paper, these terms are treated with finite differencing scheme. More precisely, we present a computational model that could be suitable for basic optical tomography forward problem with spatially varying refractive index. We investigate cases concerning optical tomography applications. We consider a biological tissue-like domain submitted to a near infrared light source. Transmitted intensity detected on the boundary is computed. The effect of the embedded nonhomogenous objects on the transmitted signal is studied. Particularly, a possible detection of these objects is discussed. Also, the model is used to study infrared radiation in a multilayered heterogeneous medium. The effect of optical parameters on the transmitted signal is presented. In particular, the effect of refractive index variation is pointed out.

#### **2. Physical model and numerical implementation**

In this work, the radiative transfer in a human biological tissue is described by using a varying refractive index RTE [26, 27]

$$\begin{split} \frac{m(\bar{r})}{c} \frac{\partial l(\bar{r}, \Omega, t)}{\partial t} + \bar{\Omega} \nabla l(\bar{r}, \bar{\Omega}, t) + (\mu\_a(\bar{r}) + \mu\_s(\bar{r})) l(\bar{r}, \Omega, t) + \frac{1}{n(\bar{r})} \nabla n \nabla\_{\bar{\Omega}} l(\bar{r}, \bar{\Omega}, t) \\ - \frac{2}{n(\bar{r})} (\bar{\Omega} \nabla n) l(\bar{r}, \bar{\Omega}, t) = \mathbf{S}(\bar{r}, \bar{\Omega}, t) + \mu\_s(r) \int\_0^{2\pi} p(\bar{\Omega}, \bar{\Omega}') l(\bar{r}, \bar{\Omega}', t) d\bar{\Omega}' \end{split} \tag{1}$$

where *Ir t* (, ,) is the directional energetic radiance at the spatial position *<sup>r</sup>* and *n r*( ) is the refractive index distribution. Equation (1) takes into account the fact that the rays are not straight lines but curves. It involves terms that illustrate the expansion or the contraction of the cross section of the tube of light rays in the medium. On the boundary, the radiance is the sum of the external source contribution and the partly-reflected radiance due to the refractive index mismatch at the boundary

$$I(\vec{r}, \vec{\Omega}, t) = S(\vec{r}\_b, \vec{\Omega}, t) + RI(\vec{r}, \vec{\Omega}\_{ref}, t), \, \vec{n}\_b, \vec{\Omega} < 0 \text{ and } \vec{n}\_b, \vec{\Omega}\_{ref} = -\vec{n}\_b, \vec{\Omega}\_{ref}$$

where *br* is a position on the boundary and *nb* is an outer normal unit vector. The reflectivity *R* can be calculated for each direction using Fresnel's relations. For a twodimensional problem and in Cartesian coordinate system of the x-y plane, the terms due to the refractive index variation can be expressed as

$$\frac{1}{n}\nabla n.\nabla\_{\vec{\Omega}}I - \frac{2}{n}(\vec{\Omega}.\nabla n)I = -\frac{1}{n}\frac{\partial n}{\partial \mathbf{x}}(\sin\varphi\frac{\partial I}{\partial \varphi} + 2\cos\varphi I) + \frac{1}{n}\frac{\partial n}{\partial y}(\cos\varphi\frac{\partial I}{\partial \varphi} - 2\sin\varphi I) \, , \tag{2}$$

where cosφ and sinφ are the Cartesian coordinates of the unit direction vector in the x-y plane. In fact we assume that the radiance of out of plane directions is negligible. By using notations cos and sin , equation (2) displays the classical form of the angular redistribution term commonly appearing when dealing with spherical and cylindrical geometries with uniform refractive index [28, 29]

$$\frac{1}{2n}\nabla n \nabla\_{\vec{\Omega}} I - \frac{2}{n}(\vec{\Omega} \nabla n)I = \frac{1}{2n^2} \left| \frac{\partial n^2}{\partial \mathbf{x}} \frac{\partial}{\partial \xi} [(1 - \xi^2)I] + \frac{\partial n^2}{\partial y} \frac{\partial}{\partial \eta} [(1 - \eta^2)I] \right| \tag{3}$$

The angular redistribution terms will be noted

$$D\_x = \frac{\partial}{\partial \xi} [(1 - \xi^2)I] \, and \, D\_y = \frac{\partial}{\partial \eta} [(1 - \eta^2)I] \,.$$

so equation (3) becomes

28 Medical Imaging

abnormal tissue in a typical DOT scheme. However, it is important to note that in a varying refractive index medium, the rays are not straight lines but curves. So even in a rectangular geometry, the varying refractive index radiative transfer equation displays angular derivative terms [23-25]. In this paper, these terms are treated with finite differencing scheme. More precisely, we present a computational model that could be suitable for basic optical tomography forward problem with spatially varying refractive index. We investigate cases concerning optical tomography applications. We consider a biological tissue-like domain submitted to a near infrared light source. Transmitted intensity detected on the boundary is computed. The effect of the embedded nonhomogenous objects on the transmitted signal is studied. Particularly, a possible detection of these objects is discussed. Also, the model is used to study infrared radiation in a multilayered heterogeneous medium. The effect of optical parameters on the transmitted signal is presented. In particular, the effect of refractive index variation is pointed out.

In this work, the radiative transfer in a human biological tissue is described by using a

() (, ,) <sup>1</sup> . ( , , ) ( ( ) ( )) ( , , ) . ( , , ) ( )

 

*a s*

*nr Ir t Ir t r r Ir t n Ir t*

*s*

refractive index distribution. Equation (1) takes into account the fact that the rays are not straight lines but curves. It involves terms that illustrate the expansion or the contraction of the cross section of the tube of light rays in the medium. On the boundary, the radiance is the sum of the external source contribution and the partly-reflected radiance due to the

reflectivity *R* can be calculated for each direction using Fresnel's relations. For a twodimensional problem and in Cartesian coordinate system of the x-y plane, the terms due to

12 1 <sup>1</sup> . ( . ) (sin 2 cos ) (cos 2sin ) *nI nI n I nI <sup>I</sup> <sup>I</sup>*

where cosφ and sinφ are the Cartesian coordinates of the unit direction vector in the x-y plane. In fact we assume that the radiance of out of plane directions is negligible. By using

redistribution term commonly appearing when dealing with spherical and cylindrical

, (2)

 

*nIr t Sr t r p Ir td*

<sup>2</sup> ( . ) ( , , ) ( , , ) ( ) ( , ') ( , ', ) ' ( )

is the directional energetic radiance at the spatial position *<sup>r</sup>*

( , , ) ( , , ) ( , , ), . 0 . . *b ref b b ref <sup>b</sup> I r t S r t RI r t n and n n* ,

*n n nx n y*

*c t n r*

2 0

is an outer normal unit vector. The

, equation (2) displays the classical form of the angular

 

(1)

is the

and *n r*( )

**2. Physical model and numerical implementation** 

varying refractive index RTE [26, 27]

refractive index mismatch at the boundary

is a position on the boundary and *nb*

the refractive index variation can be expressed as

geometries with uniform refractive index [28, 29]

*n r*

where *Ir t* (, ,)

where *br*

notations

 cos and sin

$$\frac{1}{n}\nabla n.\nabla\_{\vec{\Omega}}I - \frac{2}{n}(\vec{\Omega}.\nabla n)I = \frac{1}{2n^2}\left\{\frac{\partial n^2}{\partial x}D\_x + \frac{\partial n^2}{\partial y}D\_y\right\}$$

In our numerical implementation, we use a rectangular domain which is divided into a set of I x J elementary uniform volumes ΔV with a uniform unitary depth. The angular discretization is obtained through a discrete ordinate technique. This yields a set of M discrete directions, *<sup>m</sup>* , 1.. *m M* giving a set of angular discrete direction cosines( , ), 1,... *m m m M* . An orientation depending on the incident ray direction is adopted for each cell [12]. Calculations are done by using integration of equation (1) over an elementary volume ΔV for each discrete direction. This gives

$$\begin{aligned} &\frac{m}{c}\frac{\partial I\_{m,P}}{\partial t} + \Delta y \xi\_m (I\_{m,E} - I\_{m,\mathcal{W}}) + \Delta x \eta\_m (I\_{m,N} - I\_{m,S}) + \Delta x \Delta y (\mu\_a + \mu\_s) I\_{m,p} \\ &+ \frac{1}{n} \frac{\partial n}{\partial x} D\_{x,m} + \frac{1}{n} \frac{\partial n}{\partial y} D\_{y,m} = \Delta x \Delta y (S\_{m,P} + \mu\_s w\_m p\_{mm} I\_{m,p} \\ &+ \mu\_s \sum\_{m'=1, m' \neq m}^M w\_{m'} p\_{mm'} I\_{m',p} \end{aligned} \tag{4}$$

where *Dx m*, and *Dy*,*<sup>m</sup>* are the discrete angular derivative terms at the angular ordinate *ξ<sup>m</sup>* and *ηm* respectively and *wm* is a weighting factor. The discrete term of Henyey-Greenstein phase function is written as

$$p\_{mm'} = \frac{1 - g^2}{4\pi (1 + g^2 - 2g(\xi\_m \xi\_{m'} + \eta\_m \eta\_{m'}))^{3/2}}$$

If the direction cosines are positive, the directional radiances are known on the faces W and S and they are unknown on the faces E and N of the (i,j)-cell and also in the centre P. Therefore, we need two complementary relations to eliminate *m N*, and *m E*, , this can be obtained by using interpolation formula

$$\begin{cases} I\_{m,P} = \alpha I\_{m,E} + (1 - \alpha) I\_{m,W} \\ I\_{m,P} = \alpha I\_{m,N} + (1 - \alpha) I\_{m,S} \end{cases}$$

where α is an interpolation parameter. Using these relations, equation (3) becomes

Detection of Abnormalities in a Biological

right member. So, this gives

*mi j*

2

 *mi j <sup>k</sup> <sup>n</sup> I* , ,

value:

1

, ,

*i j*

,

*xyS*

*n*

1 2

*ct ct*

<sup>1</sup> <sup>1</sup> ,

3 2

2

*m i*

obtained when the relative discrepancy value:

,

*k kk*

*m i <sup>j</sup> m i <sup>j</sup> m i <sup>j</sup> I II* 1/2, , 1, , , ,

<sup>1</sup> ( )

*BB w m m* 1/2 1/2 2 *m m*

*n n <sup>y</sup> <sup>x</sup> II I I*

*<sup>k</sup> <sup>n</sup> i j <sup>m</sup> <sup>m</sup> ai j si j*

*mi j mij*

*j s m mm si j m mm*

 

 

, , , ' '

The iteration process is repeated until a convergence criterion is attempted. To improve convergence speed, we use a successive over relaxation method. So the updated value

*mi j mi j mi j k kk n nn*

where ρ is a relaxation parameter whose value is usually between 1 and 2. The solution is

1

*I I*

 

*I* , , , ,

is smaller than a tolerance value. In that case the result is noted *m i <sup>j</sup> I* , , . In all our calculus, we have taken 10-8 as a tolerance value. As initial condition, we take a field of null radiances. Also, all our calculations are done in the case of interpolation diamond scheme (α=0.5).

If the direction cosines are not both positive, the precedent equations are valid provided that the orientation WESN of cells is done according to the direction of propagation. In all our

*k k n n*

*<sup>k</sup> <sup>n</sup>*

 *mi j mi j*

*mi j*

, ,

1 1

(1 )

*I II* , , , , , ,

*mi j mi j mi j mi j*

, , , , , 1, ,, 1

1 1 , , 11 1

*ij i j ij ij*

2 2 2 2 2 , 1, , ,1

*yn n D xn n D*

*k k i j nn n n i j <sup>m</sup> <sup>m</sup>*

 

The coefficients A and B are given through the recurrent relations

and an equally spaced distribution of angular ordinates (,)

*<sup>n</sup> <sup>y</sup> <sup>x</sup> <sup>I</sup> x y c t*

is a linear combination of the iterated value *mi j*

*updated*

Tissue by Diffuse Optical Tomography: A Computational Study 31

*AA w m m* 1/2 1/2 2 *m m*

and *wm* being a quadrature weight factor. In this paper, we use a constant weight quadrature

we use successive iterations in order to actualise the implicit internal source term in the

 

<sup>2</sup> with *<sup>M</sup> <sup>i</sup> <sup>j</sup> <sup>i</sup> <sup>j</sup> <sup>i</sup> <sup>j</sup> M i <sup>j</sup> I III* 1, , 1, , 0, , , , , .

,

with *A B* 1/2 1/2 0 ,

 

,, ,,

 

*k k <sup>M</sup> n n*

 

*<sup>k</sup> <sup>n</sup>*

1 1

*m mm wp I wp I* ,, 1 ', , 1

' 1, '

 

 

*xmi j ymi j*

, ,, , ,,

*k k*

. To solve the above equation

*I* , , and the previously computed

1

$$\begin{aligned} &\frac{\alpha}{c}\frac{\partial I\_{m,P}}{\partial t} + \frac{\Delta y \xi\_m^{\prime}}{\alpha}(I\_{m,P} - I\_{m,W}) + \frac{\Delta x \eta\_m}{\alpha}(I\_{m,P} - I\_{m,S}) + \Delta x \Delta y (\mu\_a + \mu\_s)I\_{m,p} \\ &+ \frac{1}{2\alpha^2} \left\{ \frac{\partial n^2}{\partial x} D\_x + \frac{\partial n^2}{\partial y} D\_y \right\} = \Delta x \Delta y (S\_{m,P} + \mu\_s w\_m p\_{mm} I\_{m,p} \\ &+ \mu\_s \sum\_{m'=1, m' \neq m}^M w\_{m'} p\_{mm'} I\_{m',p} \end{aligned} \tag{5}$$

For time discretization, an implicit three level second order time differencing scheme is used:

$$\frac{\partial I\_{m,P}}{\partial t} = \frac{\Im I\_{m,P}^{n+1} - 4I\_{m,P}^n + I\_{m,P}^{n-1}}{2\Delta t}; \qquad \qquad n = 1, \ldots \\ n\_{\text{max}}$$

where *t* being the discrete time step.

Theoretically, if we know the solution in the (i,j)-cell, we can do calculus over the cells (i+1,j) and (i,j+1) using the boundary conditions and the following relations

$$\begin{cases} I\_{m,w}(i+1,j) = I\_{m,E}(i,j); i = 1,...,I-1\\ I\_{m,S}(i,j+1) = I\_{m,N}(i,j); i = 1,...,J-1 \end{cases}$$

If the direction cosines are both positive, we get the following equation

$$\begin{split} &\frac{1}{\alpha} \frac{\Im I\_{m,i,j}^{n+1} - 4I\_{m,i,j}^{n} + I\_{m,i,j}^{n-1}}{\Delta t} + \frac{\Delta y \tilde{\varepsilon}\_{m}}{\alpha} (I\_{m,i,j}^{n} - I\_{m,i-1,j}^{n}) + \frac{\Delta x \eta\_{m}}{\alpha} (I\_{m,i,j}^{n} - I\_{m,i,j-1}^{n}) + \\ &\Delta x \Delta y (\mu\_{a} + \mu\_{s}) I\_{m,i,j}^{n} + \frac{1}{2n^{2}} \left\{ \frac{\partial n^{2}}{\partial x} D\_{x} + \frac{\partial n^{2}}{\partial y} D\_{y} \right\} \\ &= \Delta x \Delta y (S\_{m,i,j} + \mu\_{s} w\_{m} p\_{mm} I\_{m,i,j}^{n} + \mu\_{s} \sum\_{m'=1, m' \neq m}^{M} w\_{m} p\_{mm'} I\_{m',i,j}^{n}) \end{split}$$

#### **3. Numerical treatment of angular derivative terms**

The numerical treatment of the angular redistribution terms can be done by using the classical angular differencing scheme [23, 28].

$$D\_{x,m,i,j} = \frac{A\_{m+1/2}I\_{m+1/2,i,j} - A\_{m-1/2}I\_{m-1/2,i,j}}{w\_m}.$$
 
$$D\_{y,m,i,j} = \frac{B\_{m+1/2}I\_{m+1/2,i,j} - B\_{m-1/2}I\_{m-1/2,i,j}}{w\_m}$$

where

$$I\_{\\_m+1\'/2,\ell,j} = \frac{I}{2}(I\_{\\_m+1,\ell,j} + I\_{\\_m,\ell,j}) \text{ and } \ell$$

*m P <sup>m</sup> <sup>m</sup> mP mW mP mS a sm <sup>p</sup>*

, , , , ,

 

*n n*

*mi j mi j mi j mi j*

 

, , , 1, ,, ,, 1

max

(5)

( ) ( ) ()

; 1,... <sup>2</sup>

*x y mP s m mm mp*

For time discretization, an implicit three level second order time differencing scheme is used:

Theoretically, if we know the solution in the (i,j)-cell, we can do calculus over the cells (i+1,j)

*I i j I ij i I I ij I ij i J*

*mi j mi j mi j m n n m n n*

2 2

*<sup>n</sup> I II <sup>y</sup> <sup>x</sup> II II*

*m i j s m mm m i j s m mm m i j*

( )

 

, , , , ' ' ', ,

*n x y*

( 1, ) ( , ); 1,... 1 ( , 1) ( , ); 1,... 1

*<sup>M</sup> n n*

( )( ) <sup>2</sup>

*m mm*

' 1, '

The numerical treatment of the angular redistribution terms can be done by using the

*AI AI*

*BI BI*

*j,i,1m*

*m mi j m mi j*

*m*

*m mi j m mi j*

*m*

*w* 1/2 1/2, , 1/2 1/2, ,

*)II( <sup>2</sup>*

*j,i,m* and

*w* 1/2 1/2, , 1/2 1/2, , , ,, .

*n n D D x yS w p I*

2 , ,

*n nn mP mP mP mP I I II*

> *m w m E m S m N*

, , , ,

1 1

)

, , ,,

*<sup>n</sup> <sup>I</sup> <sup>y</sup> <sup>x</sup> II II x <sup>y</sup> <sup>I</sup>*

 

*M*

*c t*

where

where

,

*m mm*

' 1, '

*s m mm m p*

2 2

 

*wp I*

<sup>1</sup> ( <sup>2</sup>

' ' ',

*t t*

and (i,j+1) using the boundary conditions and the following relations

If the direction cosines are both positive, we get the following equation

*a s mi j x y*

*x yS w p I wp I*

*n n xy I D D*

, , 2

**3. Numerical treatment of angular derivative terms** 

*xmi j*

*ymi j*

, ,,

*<sup>1</sup> <sup>I</sup> j,i,2/1m*

*D*

*D*

3 4

*n x y*

*t* being the discrete time step.

*n nn*

1 1 ,, ,, ,,

classical angular differencing scheme [23, 28].

*c t*

 

3 4

*n*

<sup>1</sup> ( ) <sup>2</sup>

$$I\_{\;m-1/2,i,j}^{\;k} = \frac{1}{2}(I\_{\;m-1,i,j}^{\;k} + I\_{\;m,i,j}^{\;k}) \text{ with } I\_{\;M+1,i,j} = I\_{1,i,j} \text{ } \text{ } I\_{0,i,j} = I\_{\;M,i,j} \text{ } \text{ } \text{ }$$

The coefficients A and B are given through the recurrent relations

$$A\_{m+1/2} - A\_{m+1/2} = -2w\_m \xi\_{m'}$$

$$B\_{m+1/2} - B\_{m+1/2} = -2w\_m \eta\_m \text{ with } A\_{1/2} = B\_{1/2} = 0 \text{ and }$$

and *wm* being a quadrature weight factor. In this paper, we use a constant weight quadrature and an equally spaced distribution of angular ordinates (,) . To solve the above equation we use successive iterations in order to actualise the implicit internal source term in the right member. So, this gives

$$\begin{split} & \left( I\_{m,i,j}^{n+1} \right)^{k+1} = \left[ \frac{3n\_{i,j}}{2c\Delta t} + \frac{\Delta y\_j^\varepsilon \varepsilon\_m}{\alpha} + \frac{\Delta x \eta\_m}{\alpha} + \Delta x \Delta y \left( \mu\_{a,i,j} + \mu\_{s,i,j} \right) \right]^{-1} \times \\ & \left[ \frac{2n\_{i,j}}{c\Delta t} \overline{I\_{n,i,j}^{n}} - \frac{n\_{i,j}}{2c\Delta t} \overline{I\_{n,i,j}^{n-1}} + \frac{\Delta y\_j^\varepsilon \varepsilon\_m}{\alpha} \left( I\_{m,i-1,j}^{n+1} \right)^{k+1} + \frac{\Delta x \eta\_m}{\alpha} \left( I\_{m,i,j-1}^{n+1} \right)^{k+1} \right. \\ & \left. + \frac{1}{2n\_{i,j}^{n+1}} \left( \Delta y \left( n\_{i,j}^2 - n\_{i-1,j}^2 \right) \right) \mathbf{D}\_{x,m,i,j}^{k} + \Delta x \left( n\_{i,j}^2 - n\_{i,j-1}^2 \right) \mathbf{D}\_{y,m,i,j}^{k} \right) \\ & \left. + \Delta x \Delta y \left( S\_{m,i,j} + \mu\_s w\_m p\_{mm} \left( I\_{m,i,j-1}^{n+1} \right)^k + \mu\_{s,i,j} \sum\_{m'=1,m'\neq m}^{M} w\_{m'} p\_{mm'} \left( I\_{m',i,j-1}^{n+1} \right)^k \right) \right\} \end{split}$$

The iteration process is repeated until a convergence criterion is attempted. To improve convergence speed, we use a successive over relaxation method. So the updated value *mi j <sup>k</sup> <sup>n</sup> I* , , 1 is a linear combination of the iterated value *mi j <sup>k</sup> <sup>n</sup> I* , , and the previously computed value:

$$\left(I\_{m,i,j}^{n}\right)^{k+1}\_{\
u \text{quadted}} = \rho \left(I\_{m,i,j}^{n}\right)^{k} + (1-\rho) \left(I\_{m,i,j}^{n}\right)^{k+1}$$

where ρ is a relaxation parameter whose value is usually between 1 and 2. The solution is obtained when the relative discrepancy value:

$$\mathcal{E} = \left| \frac{\left(I\_{\
u,i,j}^n\right)^{k+1} - \left(I\_{\
u,i,j}^n\right)^k}{\left(I\_{\
u,i,j}^n\right)^k} \right|$$

is smaller than a tolerance value. In that case the result is noted *m i <sup>j</sup> I* , , . In all our calculus, we have taken 10-8 as a tolerance value. As initial condition, we take a field of null radiances. Also, all our calculations are done in the case of interpolation diamond scheme (α=0.5).

If the direction cosines are not both positive, the precedent equations are valid provided that the orientation WESN of cells is done according to the direction of propagation. In all our

Detection of Abnormalities in a Biological

**y-axis** 

Tissue by Diffuse Optical Tomography: A Computational Study 33

**air, n=1** 

Fig. 1. A test-medium with background properties and detector points

*a s*

 

*g*

*cm cm*

**y-varying RI** 

**Detector points**

0.5 25 0.8

1 1

**Source x-axis** 

**(a) (b)** 

 

Fig. 2. Response of a continuous varying refractive index medium (a) in the case of a linear

In this investigation, we consider a 3cmx3cm sized background medium containing a circular heterogeneous object. The area of the object is 0.5 cm2. We limit the analysis to the effect of heterogeneity by increase of the refractive index only. The background refractive index is taken 1.33. The other optical properties are the same as in the precedent investigation for the background and the object. The geometry of the medium and the position of the heterogeneous object are shown in figure 3. The detected signal on different

variation and (b) in the case of a parabolic variation.

**4.2 Second investigation: Detection of heterogeneous object** 

investigations, the injected power source is assumed to be equivalent to a forward collimated monochromatic intensity placed at a source point on the middle of the bottom side of the boundary. Results shown below are obtained by using a near infrared collimated source with a uniform equivalent intensity value of 50 mW.cm-1. To present our results, we use the detected fluence rate which is given in a (id,jd)-detector point on the boundary as

$$\boldsymbol{\Phi}\_{d} = \sum\_{m=1}^{M} (\mathbf{1} - \mathbf{R}\_{m, i\_{d}, j\_{d}}) w\_{m} \, \overline{I}\_{m, i\_{d}, j\_{d}}, \text{ with}$$

$$\mathbf{R}\_{m, i\_{d}, j\_{d}} = \begin{cases} \mathbf{1} \, \text{if } \boldsymbol{\varphi}\_{m} \succ \arcsin\left(\frac{n\_{\text{air}}}{n\_{i\_{d}, j\_{d}}}\right) \\\\ \frac{1}{2} \mathbf{I} [\!(\!n\_{i\_{d}, j\_{d}} \cos \boldsymbol{\varphi}\_{m} - \underset{\text{air}}{\text{air}} \cos \boldsymbol{\varphi}\_{l}) \\\ \!(\!n\_{i\_{d}, j\_{d}} \cos \boldsymbol{\varphi}\_{m} + \underset{\text{air}}{\text{air}} \cos \boldsymbol{\varphi}\_{l}) \end{cases} + \begin{cases} n\_{i\_{d}, j\_{d}} \cos \boldsymbol{\varphi}\_{l} - n\_{\text{air}} \cos \boldsymbol{\varphi}\_{m} \\\ n\_{i\_{d}, j\_{d}} \cos \boldsymbol{\varphi}\_{l} + n\_{\text{air}} \cos \boldsymbol{\varphi}\_{m} \end{cases} \} 1. \text{else}.$$

where

$$\varphi\_t = \arcsin\left(\frac{n\_{air}\sin\varphi\_m}{n\_{i\_d,j\_d}}\right) \text{ and } n\_{air} \approx 1\text{ .}$$

Also, we make use of a normalized detected fluence rate defined as

$$\Phi\_N = \frac{\Phi\_d}{1 / D \sum\_{d=1}^D w\_d^\prime \Phi\_d}$$

,

where *D* is the number of the detector points on one side of the boundary and w'd is a weighting factor from the generalized trapezoidal integration rule. In all calculations, we have used 28 detector points on each side. Also, all calculus is carried out by using 16 uniformly distributed discrete directions and a space grid of 121x121cells.

#### **4. Results and discussion**

#### **4.1 First test: A continuous varying refractive index medium**

As a first investigation, we study near infrared radiation transport in a rectangular medium exposed to a continuous collimated source which is placed on the bottom side of the boundary. Figure 1 shows the considered medium, it is assumed to be 3cmx3cm sized with varying refractive index. Within the medium, we consider a y-axis linear varying refractive index between the extreme values 1.33 and 1.58. It is an increasing- decreasing variation of refractive index RI. To show the effect of the refractive index on detected fluence rate, we have used a weakly absorbing and forward scattering background medium whose optical parameters are shown in figure1. Figure 2 shows the medium response detected on the right side of the boundary in the case of linear and parabolic variation respectively. The detected fluence rate curve presents distinguish increasing and decreasing regions with a peak of transmission when refractive index is close to the minimum value 1.33. Transmission on the boundary is relatively augmented when the refractive index is the closest to the air index.

investigations, the injected power source is assumed to be equivalent to a forward collimated monochromatic intensity placed at a source point on the middle of the bottom side of the boundary. Results shown below are obtained by using a near infrared collimated source with a uniform equivalent intensity value of 50 mW.cm-1. To present our results, we use the detected fluence rate which is given in a (id,jd)-detector point on the boundary as

*d m i j m mi j*

*n n nn*

 

1 cos cos cos cos

*n n nn*

*d d*

*d*

*d D w* 1 1/ ' 

where *D* is the number of the detector points on one side of the boundary and w'd is a weighting factor from the generalized trapezoidal integration rule. In all calculations, we have used 28 detector points on each side. Also, all calculus is carried out by using 16

As a first investigation, we study near infrared radiation transport in a rectangular medium exposed to a continuous collimated source which is placed on the bottom side of the boundary. Figure 1 shows the considered medium, it is assumed to be 3cmx3cm sized with varying refractive index. Within the medium, we consider a y-axis linear varying refractive index between the extreme values 1.33 and 1.58. It is an increasing- decreasing variation of refractive index RI. To show the effect of the refractive index on detected fluence rate, we have used a weakly absorbing and forward scattering background medium whose optical parameters are shown in figure1. Figure 2 shows the medium response detected on the right side of the boundary in the case of linear and parabolic variation respectively. The detected fluence rate curve presents distinguish increasing and decreasing regions with a peak of transmission when refractive index is close to the minimum value 1.33. Transmission on the boundary is relatively augmented when the refractive index is the closest to the air index.

*d d*

*i j*

*n* ,

*N D*

 

(1 )

*M*

*m*

*air <sup>m</sup>*

Also, we make use of a normalized detected fluence rate defined as

uniformly distributed discrete directions and a space grid of 121x121cells.

**4.1 First test: A continuous varying refractive index medium** 

*<sup>n</sup> if*

 

1 arcsin

*d d*

*mi j*

, ,

**4. Results and discussion** 

*R*

where

1

*d d*

,

*i j*

*n*

*d d d d d d d d*

, ,

*air m <sup>t</sup>*

sin arcsin

*n*

*d d d d*

*R wI* ,, ,,

*i j m air t i j t air m i j m air t i j t air m*

, , 2 2

and *nair* 1 .

,

[( ) ( ) ], . 2 cos cos cos cos

*else*

, with

Fig. 1. A test-medium with background properties and detector points

Fig. 2. Response of a continuous varying refractive index medium (a) in the case of a linear variation and (b) in the case of a parabolic variation.

#### **4.2 Second investigation: Detection of heterogeneous object**

In this investigation, we consider a 3cmx3cm sized background medium containing a circular heterogeneous object. The area of the object is 0.5 cm2. We limit the analysis to the effect of heterogeneity by increase of the refractive index only. The background refractive index is taken 1.33. The other optical properties are the same as in the precedent investigation for the background and the object. The geometry of the medium and the position of the heterogeneous object are shown in figure 3. The detected signal on different

Detection of Abnormalities in a Biological

inclusions.

(n1=1.37 and n2=1.51).

Tissue by Diffuse Optical Tomography: A Computational Study 35

In the course of this investigation, we simulate a post-mortem human prostate tissue-like medium containing two inclusions. The geometry of the medium and the targets of the inclusions are shown in figure 3. Each inclusion is a circle of area 0.5cm2. The background medium properties at 850nm wavelength are taken from [30, 31]. Inclusion 1 is double scattering only and inclusion 2 is double absorbing only. The refractive index of the two inclusions is taken the same and it is noted n2. Detected fluence rate on the top and on the lateral sides of the medium is computed and presented in figure 5 and figure 6. In particular, we confront results from refractive matching case (n1=n2) and from refractive mismatching case (n2>n1) on each side. According to figure 5, the introduction of refractive index with scattering and absorbing parameters gives a distinguish effect on top side detected signal. However this effect is more attenuated in lateral sides because of the strong forward

**air, n=1**

**Inclusion 1 2μs, n2** 

Fig. 4. Heterogeneous inclusions in human prostate-like tissue: geometry of the medium and

 

**1 1**

**100 , 0.94**

**1 2**

*cm cm g*

Fig. 5. Effect of heterogeneous inclusions on detected fluence rate on the top side (a) in the case of refractive matching (n1=n2=1.37) and (b) in the case of refractive mismatching

**4.3 Third investigation: Detection of inclusions in a biological tissue** 

*a s*

**Inclusion 2 2μa, n2**

**0.6**

**1.37**

**1**

*n*

scattering within the medium. This is shown in figure 6.

**2**

sides of the medium is shown for an increase of refractive index by 5% and 10% respectively. There is no significant effect of the object presence on the detected signal on the left side because it is far from the detectors on that side. However there is a visible effect of the object on the right side and on the top side where the detectors are relatively near by the object. In such cases, there is an obvious effect of the heterogeneity on the detected signal especially when the refractive index is increased by 10%. The response shows a distinct distortion of the fluence rate curve. It is occurring on the targets where the object centre is laying exactly. This distortion allows a principled decision on the existence of such object in the medium. These findings highlight the potential of refractive index as a possible detection parameter of a tumour in a surrounding safe tissue.

Fig. 3. Detection of heterogeneous objects by the effect of refractive index variation. (a) Geometry of the considered medium and object, (b) Detected fluence rate on the left side,(c) Detected fluence rate on the right side (d) Detected fluence rate on the top side.

sides of the medium is shown for an increase of refractive index by 5% and 10% respectively. There is no significant effect of the object presence on the detected signal on the left side because it is far from the detectors on that side. However there is a visible effect of the object on the right side and on the top side where the detectors are relatively near by the object. In such cases, there is an obvious effect of the heterogeneity on the detected signal especially when the refractive index is increased by 10%. The response shows a distinct distortion of the fluence rate curve. It is occurring on the targets where the object centre is laying exactly. This distortion allows a principled decision on the existence of such object in the medium. These findings highlight the potential of refractive index as a possible

**air, n=1** 

**2** 

x-axis

**(a)** 

**Background (RI: n)** 

**Object (RI: n +n)** 

**(c) (d)** 

**y-axis** 

**2** 

**(b)** 

Fig. 3. Detection of heterogeneous objects by the effect of refractive index variation. (a) Geometry of the considered medium and object, (b) Detected fluence rate on the left side,(c)

Detected fluence rate on the right side (d) Detected fluence rate on the top side.

detection parameter of a tumour in a surrounding safe tissue.

#### **4.3 Third investigation: Detection of inclusions in a biological tissue**

In the course of this investigation, we simulate a post-mortem human prostate tissue-like medium containing two inclusions. The geometry of the medium and the targets of the inclusions are shown in figure 3. Each inclusion is a circle of area 0.5cm2. The background medium properties at 850nm wavelength are taken from [30, 31]. Inclusion 1 is double scattering only and inclusion 2 is double absorbing only. The refractive index of the two inclusions is taken the same and it is noted n2. Detected fluence rate on the top and on the lateral sides of the medium is computed and presented in figure 5 and figure 6. In particular, we confront results from refractive matching case (n1=n2) and from refractive mismatching case (n2>n1) on each side. According to figure 5, the introduction of refractive index with scattering and absorbing parameters gives a distinguish effect on top side detected signal. However this effect is more attenuated in lateral sides because of the strong forward scattering within the medium. This is shown in figure 6.

Fig. 4. Heterogeneous inclusions in human prostate-like tissue: geometry of the medium and inclusions.

Fig. 5. Effect of heterogeneous inclusions on detected fluence rate on the top side (a) in the case of refractive matching (n1=n2=1.37) and (b) in the case of refractive mismatching (n1=1.37 and n2=1.51).

Detection of Abnormalities in a Biological

different refractive indices.

**Detector** 

<sup>1</sup> ' 0.00012 cm

<sup>1</sup> 0.0049cm

Heterogeneous

*u S*

*uQ*

RI=n1

pulse and the detected signal.

**Short pulse** 

more intensity is transmitted from the muscle layer.

**4.5 Fifth investigation: Detection with a short pulse** 

*u s*

*uQ*

RI=n2

<sup>1</sup> ' 13.6cm

<sup>1</sup> 0.35cm

Background skin

Tissue by Diffuse Optical Tomography: A Computational Study 37

thickness increases. When taking different refractive indices, a distinct distortion in the curves at every interface is observed especially at the fat-muscle interface. In all cases much

Fig. 8. Response of the three layered medium: (a) for uniform refractive index and (b) for

In this investigation, we study the effect of abnormalities laying in a homogenous background medium as in the just precedent paragraph but we use a short-pulsed source. It is a 100fs-laser pulse injected in the medium through a point source on the middle of the bottom side of the boundary. Figure 9 shows the considered medium and the geometry and the properties of the

Fig. 9. A skin-like tissue containing a heterogeneous cyst-like inclusion submitted to short

Fig. 6. Detected fluence rate on lateral sides in both cases refractive matching and mismatching cases: (a) Response on the right side (effect of inclusion 1) and (b) Response on the left side ( effect of inclusion 2)

#### **4.4 Fourth investigation: Effect of refractive mismatching in multilayered tissue**

As a fourth investigation, we consider a multilayered biological tissue-like medium. It contains three superposed layers of skin, fat and muscle. Figure 7 shows the considered medium. It is 3cmx3cm-sized medium with a fixed layer of skin (2mm) and with a varying fat-layer thickness of 3mm, 6mm and 9mm respectively. Optical properties of different tissues at 850nm-wavelength and refractive indices of the layers are taken from [30, 31]. Calculus is carried out as in the just previous section.

```
       1 1
    1
  0.12 , 87.5 , 0.8
  ,
 a s cm cm g
Skin n
       1 1
    2
  0.01 , 40 , 0.8
 ,
a s cm cm g
Fat n
       1 1
      3
  0.15 , 44.7 , 0.8
    ,
a s cm cm g
Muscle n
```
Figure 8 shows results concerning transmitted fluence rate on the right side of the boundary in the case of a continuous wave source. We can see that in most detector points and for the three fat-layer thicknesses transmitted fluence rate profiles are almost superposed in the case of refractive matching with evidence of detected intensity increase when the fat

Fig. 6. Detected fluence rate on lateral sides in both cases refractive matching and

**4.4 Fourth investigation: Effect of refractive mismatching in multilayered tissue** 

the left side ( effect of inclusion 2)

muscle at 850 nm-wavelength.

Calculus is carried out as in the just previous section.

**,**

*Skin n*

**,**

*Fat n*

mismatching cases: (a) Response on the right side (effect of inclusion 1) and (b) Response on

As a fourth investigation, we consider a multilayered biological tissue-like medium. It contains three superposed layers of skin, fat and muscle. Figure 7 shows the considered medium. It is 3cmx3cm-sized medium with a fixed layer of skin (2mm) and with a varying fat-layer thickness of 3mm, 6mm and 9mm respectively. Optical properties of different tissues at 850nm-wavelength and refractive indices of the layers are taken from [30, 31].

Fig. 7. A three-layered biological tissue-like medium with optical properties of skin, fat and

**0.12 , 87.5 , 0.8**

**0.01 , 40 , 0.8**

**0.15 , 44.7 , 0.8**

 **1 1**

 **1 1**

 **1 1**

*a s cm cm g*

**3**

*a s cm cm g*

*a s cm cm g*

**1**

**2**

**,**

*Muscle n*

Figure 8 shows results concerning transmitted fluence rate on the right side of the boundary in the case of a continuous wave source. We can see that in most detector points and for the three fat-layer thicknesses transmitted fluence rate profiles are almost superposed in the case of refractive matching with evidence of detected intensity increase when the fat

thickness increases. When taking different refractive indices, a distinct distortion in the curves at every interface is observed especially at the fat-muscle interface. In all cases much more intensity is transmitted from the muscle layer.

Fig. 8. Response of the three layered medium: (a) for uniform refractive index and (b) for different refractive indices.

## **4.5 Fifth investigation: Detection with a short pulse**

In this investigation, we study the effect of abnormalities laying in a homogenous background medium as in the just precedent paragraph but we use a short-pulsed source. It is a 100fs-laser pulse injected in the medium through a point source on the middle of the bottom side of the boundary. Figure 9 shows the considered medium and the geometry and the properties of the

Fig. 9. A skin-like tissue containing a heterogeneous cyst-like inclusion submitted to short pulse and the detected signal.

Detection of Abnormalities in a Biological

after the pulse.

Tissue by Diffuse Optical Tomography: A Computational Study 39

Fig. 11. Internal light distribution in the rat-liver tissue-like medium for different instants

Fig. 12. Typical trajectory of photons for different values of the wavelength-source as

obtained through Monte Carlo simulations.

a circular cyst-like object laying just forward the source. The detected signal on the opposite side of the source is shown in same figure. There is a distinct phase shift due to refractive mismatching showing the presence of the cyst within the background tissue.

#### **4.6 Sixth investigation: Spectroscopy on a rat-liver**

In this investigation, we consider a rat liver tissue-like medium as it is shown in figure 10. Optical properties of the tissue are taken from reference [30] from a set of several wavelengths. The response illustrate transmission through the considered medium on the top side. Only near infrared light displays a relative appreciable transmission. A maximum of detected fluence rate is noted at 800nm-light. A little delay of the maximum detection of radiation is observed when wavelength increases. Figure 11 shows fluence rate distribution into the medium at different moments after the pulse. In all moments, symmetric propagation of light is observed into the medium. Also, it can be observed that scattering prevails in 488nm while more and more absorption is observed in 2100 nm. At 488nm, the phenomenon of multiple scattering is dominant into the medium. The dispersion of the pulse is accentuated in all sides. Light disappears after 250ps. At 2100nm, more light is absorbed into the medium so a fraction of energy persists into the medium for more time.

Figure 12 shows typical movement of different photons in the medium at three different wavelengths. Calculations in our studied medium are carried out by using Monte Carlo simulations as it is described in reference. It can be observed that a typical 800nm-photon is almost snake. For 488nm-photon, the mean scattering free path is very small. Multiple scattering is dominant. The typical trajectory in this case is not ballistic and the photons are almost diffusive. For 2100nm, the mean scattering free path is longer than the 488nm-photon but the absorption mean free path is narrower. A typical photon in this wavelength is almost absorbed into the medium.

Fig. 10. Geometry of rat-liver tissue and detected signal for different wavelengths.

a circular cyst-like object laying just forward the source. The detected signal on the opposite side of the source is shown in same figure. There is a distinct phase shift due to refractive

In this investigation, we consider a rat liver tissue-like medium as it is shown in figure 10. Optical properties of the tissue are taken from reference [30] from a set of several wavelengths. The response illustrate transmission through the considered medium on the top side. Only near infrared light displays a relative appreciable transmission. A maximum of detected fluence rate is noted at 800nm-light. A little delay of the maximum detection of radiation is observed when wavelength increases. Figure 11 shows fluence rate distribution into the medium at different moments after the pulse. In all moments, symmetric propagation of light is observed into the medium. Also, it can be observed that scattering prevails in 488nm while more and more absorption is observed in 2100 nm. At 488nm, the phenomenon of multiple scattering is dominant into the medium. The dispersion of the pulse is accentuated in all sides. Light disappears after 250ps. At 2100nm, more light is absorbed into the medium so a fraction of energy persists into the medium for more time. Figure 12 shows typical movement of different photons in the medium at three different wavelengths. Calculations in our studied medium are carried out by using Monte Carlo simulations as it is described in reference. It can be observed that a typical 800nm-photon is almost snake. For 488nm-photon, the mean scattering free path is very small. Multiple scattering is dominant. The typical trajectory in this case is not ballistic and the photons are almost diffusive. For 2100nm, the mean scattering free path is longer than the 488nm-photon but the absorption mean free path is narrower. A typical photon in this wavelength is

mismatching showing the presence of the cyst within the background tissue.

Fig. 10. Geometry of rat-liver tissue and detected signal for different wavelengths.

**4.6 Sixth investigation: Spectroscopy on a rat-liver** 

almost absorbed into the medium.

A Rat liver tissue 2cm

**Short pulse** 

Fig. 11. Internal light distribution in the rat-liver tissue-like medium for different instants after the pulse.

Fig. 12. Typical trajectory of photons for different values of the wavelength-source as obtained through Monte Carlo simulations.

Detection of Abnormalities in a Biological

(2002)

478, (2001)

736, (2002)

(2007)

(2005)

Phys Med Biol; 42: 841–853, (1997)

communications; 282:4308-4314, (2009)

Lasers in Med Sci; 25:41-53, (2010)

tomography. Opt Express; 2 (6): 213–226, (1998)

Tissue by Diffuse Optical Tomography: A Computational Study 41

[3] Arridge S R, Hebden J C. Optical imaging in medicine: II. Modeling and reconstruction.

[4] Arridge S R, Schweiger M. A gradient-based optimisation scheme for optical

[5] Hielscher A H, Bluestone AY, Abdoulaev AY, Klose A D, Lasker M, Stewart M, Netz U,

[6] Mourant J R, Hielscher A H, Eick A A, Johnson T M, Freyer J P. Evidence of intrinsic

[7] Ntziachristos V, Hielscher A H, Yodh A G, Chance B. Diffuse Optical Tomography of

[8] Pu Y, Wang W.B , Das BB, Alfano RR.Time-resolved spectral wing emission kinetics and

[9] Klose, A, D, Netz,U, Beuthan,J and Hielscher, A, H, Optical tomography using the time-

[10] Cai W, Xiaohui Ni, Gayen S K, and Alfano R R. Analytical cumulant solution of the

[11] Gantri M, Trabelsi H, Bensalah R and Sediki E. Solution of a radiative transfer problem

[12] Trabelsi H, Gantri M and Sediki E. Visible and Near Infrared Laser Radiation in a

[13] Gantri M, Trabelsi H, Sediki E, Ben Salah R. Computational Laser Spectroscopy in a Biological tissue. Journal of Biophysics, article ID 253763, 9pages, (2010) [14] Li H, Xie S. Measurement method of the refractive index of biotissue by total internal

[15] Das B B, Liu F , Alfano R R. Time-resolved fluorescence and photon migration studies in biomedical and model random media. Rep. Prog. Phys; 60: 227-292, (1997) [16] Ohmi M, Ohnishi Y, Yoden K, Haruan H. In vitro simultaneous measurement of refractive

[17] Lai J, Li Z, Wang C, He A. Experimental measurement of the refractive index of

[18] Beuthan J, Minet O, Helfmann J, Herrig M, Muller G. The spatial variation of the

[19] Dehghani H, Brooksby B, Vishwanath, K, Pogue B W, Paulsen K D. The effects of internal

refractive index of biological cells. Phys Med Biol; 41:369–82 (1996)

Transactions on Biomedical Engineering; 47(9):1266-1270, (2000)

cells. Cancer (Cancer Cytopathology); 84 (6): 366-384, (1998)

light from turbid media. Phys. Rev. E (74) 056605, (2006)

reflection. Applied Optics; 35(10):1793-1795, (1996)

Beuthan J. Near-infrared diffuse optical tomography. Disease Markers; 18: 313–337,

differences in the light scattering properties of tumoregenic and nontumorigenic

Highly Heterogeneous Media. IEEE transactions on medical imaging; 20(6) :470-

optical imaging of human cancerous and normal prostate tissues. Optics

independent equation of radiative transfer. Part I: Forward model, *JQSRT*, 72: 715-

vector radiative transfer equation investigates backscattering of circularly polarized

in a biological tissue. An optical tomography model. AIP Conf. Proc; 935: 237-243,

Biological Tissue. A Forward Model for Medical Imaging by Optical Tomography.

index and thickness of biological tissue by the low coherence interferometry. IEEE

biological tissues by total internal reflection. Applied Optics; 44(10):1845-1849,

refractive index variation in near-infrared optical tomography: a finite element modelling approach, Physics in Medicine and Biology; 48 (16):2713–2727, (2003)

## **5. Conclusion**

In this chapter, we attempted to develop a computational way helping in detection of abnormalities in a biological tissue. This should enable predictions of eventual tumor existence when using a diffuse optical tomography scheme. The used model is based on radiative transfer theory including a possible variation of refractive index. The model is implemented to investigate some practical situations in DOT setting and near infrared spectroscopy. Obtained results showed that variation of refractive index can yield useful predictions about the target and the location of abnormal inclusions within the tissue.

## **6. Nomenclature**



#### **Greek symbols**


## **7. References**


In this chapter, we attempted to develop a computational way helping in detection of abnormalities in a biological tissue. This should enable predictions of eventual tumor existence when using a diffuse optical tomography scheme. The used model is based on radiative transfer theory including a possible variation of refractive index. The model is implemented to investigate some practical situations in DOT setting and near infrared spectroscopy. Obtained results showed that variation of refractive index can yield useful predictions about the target and the location of abnormal inclusions within the tissue.

[1] Selser J C, Yeh Y, Baskin R J. A light scattering characterization of membrane vesicles.

[2] Jacques S L, Prahl S A. Modeling optical and thermal distributions in tissue during laser

irradiation. Lasers in Surgery and Medicine; 6: 494-503, (1987)

**5. Conclusion** 

**6. Nomenclature** 

g anisotropy factor I intensity, [W.cm-2sr-1] n refractive index p(Ω,Ω') phase function r position vector, [cm] R reflectivity coefficient,

R I Refractive Index

**Greek symbols** 

**7. References** 

c Speed of light, (in vacuum), [cm/s] DOT Diffuse Optical Tomography Dx,y Angular derivative terms

rb position vector in the boundary

w'd weighting factor from trapezoidal integration, wm weight associated with discrete direction

S source term, [W.cm-3sr-1]

x,y cartesian coordinates, [cm]

α interpolation parameter Ф fluence rate, [W.cm-2]

ρ relaxation parameter, [-] θ scattering angle , [rad] ξ, η direction cosines, [-]

λ wavelength, [nm] Ω direction solid angle, [sr] μa absorption coefficient, [cm-1] μs scattering coefficient, [cm-1]

Фd Detected fluence rate, [W.cm-2]

μ's reduced scattering coefficient, [cm-1]

Biophysical Journal; 16 (4): 337-356, (1976)


**0**

**3**

*France*

**3D Ultrasound Image Segmentation: Interactive**

The recent breakthroughs in 3D medical imaging technologies open new promising perspectives in the health domain. Significant efforts have being carried out and the precision of acquisition systems keeps on being improved. These devices are becoming widespread in the hospitals and the constantly increasing queues for this kind of exams prove the interest in these technologies. On the contrary, it seems that the evolution of 3D image analysis tools does not generate the same interest. It may be because, for a long time, people believed that considering 3D images as a succession of 2D frames was a sufficient and efficient way to make precise analysis. In our opinion, it is not the case and huge advancements can be obtained by considering 3D images with their entire complexity. This approach needs the development of specific models which can deal with multiple knowledge and are able to process huge information quantities. As a consequence, this chapter will present two original and efficient approaches of computer science for 3D image analysis and more particularly 3D ultrasound

Ultrasound techniques present a certain number of advantages when compared to other acquisition techniques like magnetic resonance imaging (MRI), X-ray computed tomography (CT), etc. Ultrasounds are not ionizing, which means they are not invasive for the patients. Moreover this technique is inexpensive and allows real time acquisitions. Nevertheless, the interpretation of ultrasound images is very complex and specialists are usually needed throughout the examination. In the same way, the segmentation (extraction of specific regions in the image) of ultrasound images is a very difficult task as these medical images present

(Noble & Boukerroui (2006)) proposed a complete survey about native B-mode ultrasound image segmentation in which they tried to point out what makes a good ultrasound segmentation method. Among efficient techniques, the authors have identified methods dealing with image features such as gray level distribution, intensity gradient, phase, similarity measures and texture measures. As described in (Noble & Boukerroui (2006)), texture analysis is very efficient for ultrasound classification and segmentation. The classical method of Haralick's co-occurrence matrices (Haralick et al. (1973)) has often been used and has obtained good performances in several applications (Basset et al. (1993); Valckx &

characteristic artifacts like shadows, speckle, attenuations, missing boundaries etc.

**1. Introduction**

images.

<sup>1</sup>*Université François Rabelais Tours, Laboratoire Informatique (EA2101)*

**Texture-Based Approaches**

Julien Olivier2,1 and Ludovic Paulhac1

<sup>2</sup>*École Nationale d'Ingénieurs du Val de Loire*

