Multi-wavelength spatial frequency domain diffuse optical tomography using single-pixel imaging based on lock-in photon counting

Tongxin Li; Zhuanping Qin; Xi Hou; Mai Dan; Jiao Li; Limin Zhang; Zhongxing Zhou; Feng Gao

doi:10.1364/OE.27.023138

1. Introduction

In the last few years, spatial frequency domain (SFD) imaging technique has increasingly attracted interests since it allows rapid and non-contact wide-field imaging of the tissue optical properties [1–7]. This technique relies on generating structured light with two-dimensional (2D) sinusoidal patterns and projecting them on a tissue surface. The tissue acts as a spatial filter and blurs the sinusoidal patterns. The spatially modulated reflectance images from the tissue surface are captured and demodulated for extracting direct current (DC) and alternating current (AC) amplitudes at all locations on the tissue surface. Through a model-based analysis, the extracted DC and AC amplitudes can be used for quantitatively characterizing the tissue optical properties (mainly the absorption and scattering coefficients). From the tissue optical properties at two or more wavelengths, the concentrations of biologically relevant chromophores, such as oxy- and deoxy-hemoglobin, lipid, and water can be deduced [1–3]. Knowledge of these chromophores provides a wealth of information about the composition and metabolic activity of tissue, which has been proven useful for characterizing tumor cells in vivo and potentially provide a contrast for the improvement of the intraoperative guidance or image-guided surgery in clinical practice [8–13].

Conventional SFD imaging technique has mainly used expensive electron-multiplying charge-coupled device (EMCCD) cameras for acquisition of spatially modulated reflectance images [1–7]. Collection of spatially modulated reflectance images for multiple wavelengths, however, requires switching one wavelength to other wavelengths. Recently, a compressed single-pixel SFD imaging technique has been developed to quantitatively characterize the tissue optical properties for multiple wavelengths simultaneously [14]. This technique employs frequency-encoded sources and spatially modulated sinusoidal patterns for illuminating tissue. The light reflected from the tissue surface is modulated again by a digital micromirror device (DMD) according to many sampling patterns. The reflectance corresponding to each sampling pattern is then integrated, fed into a photo-diode (PD) detector, and frequency-decoded to extract the signals for each wavelength. Together with compressed sensing-based single-pixel imaging (CS-SPI) scheme, the spatially modulated reflectance images for the multiple wavelengths are obtained, and subsequently the tissue optical properties for the multiple wavelengths are recovered [15]. Although this technique lowers the SFD imaging system costs and enables characterizing the tissue optical properties for multiple wavelengths simultaneously, frequency-decoding for wavelength-discrimination relies on application of filtering to separate different frequency components in Fourier space. Generally, filtering processes are prone to be associated with frequency cross-talks or overlapping among multi-frequency components, especially when the multiple frequencies are close to one another [16]. Moreover, selection of appropriate filtering windows requires manual intervention, which may result in information loss associated with the division of the frequency spectrum into multiple components.

This paper proposes a cost-effective single-pixel SFD imaging system based on highly-sensitive lock-in photon counting detection for simultaneously acquiring multi-wavelength spatially modulated reflectance images of turbid media. To avoid the non-linear effect in frequency-encoding light sources, the square waves with different frequencies were utilized in the proposed system for frequency-encoding light-emitting diode (LED) sources [17,18]. Two DMDs were incorporated, one for generating wide-field sinusoidal illumination patterns, the other for generating sampling patterns of spatially modulated reflectance images before the reflected light being spatially integrated. Spatially integrated light signals were frequency-decoded by a novel highly-sensitive lock-in photon counting detection, instead of employing the conventional filtering process in the Fourier space. The lock-in photon-counting detection features cost-effective and ultra-high sensitivity that suitable for detecting faint optical signals from biological tissues, and enabled not only frequency-decoding the faint multi-wavelength signals without information mixing, but also effectively rejecting the ambient light during data acquisition [19–23]. Spatially modulated reflectance images for the three wavelengths were recovered with the single-pixel imaging (SPI) method using the frequency-decoded multi-wavelength data sets. To dramatically reduce sampling patterns, the 2D discrete cosine transform (DCT) was incorporated in the SPI method [24,25]. Since the spatially modulated reflectance images were sparse in DCT domain, only a few number of sampling patterns were required for the measurement, and, thereby, the proposed DCT-SPI scheme had a potential to provide multi-wavelength spatially modulated reflectance images in real-time, that is desired for applications of dynamic SFD imaging. The recovered spatially modulated reflectance images for the three wavelengths were demodulated by a single snapshot multiple frequency demodulation (SSMD) method to extract the DC and AC amplitudes at all locations on the media surface [26,27]. We also put forward a spatial frequency domain diffuse optical tomography (SFD-DOT), using both the extracted DC and AC amplitudes with an inversion algorithm based on the first-order Rytov approximation of the diffusion equation (DE), to acquire tomographic images of turbid media. The proposed DCT-SPI-based SFD-DOT approach was experimentally verified by a series of tissue-simulating phantom experiments, and the results were compared with those of the conventional camera-based SFD imaging system.

2. Methods

2.1. Multi-wavelength single-pixel SFD imaging system based on lock-in photon counting detection

A schematic of the system is shown in Fig. 1(a). Illumination sources were three LEDs (M455F1, M532F2, M660F1, Thorlabs) emitting light at the wavelengths of 455 nm, 532 nm and 660 nm, which were controlled by three independent LED drivers (LEDD1B, Thorlabs), with the maximum optical powers of the three LEDs of 11 mW, 9.9 mW, and 14.5 mW, respectively. For the wavelength discrimination, the three LEDs were intensity-modulated by square waves for frequency-encoding with the modulation frequencies of $f^{(1)} = 803 Hz$ , $f^{(2)} = 654 Hz$ and $f^{(3)} = 505 Hz$ , respectively, and with the same modulation amplitude of 3.3 V, using a three-channel square wave generator implemented in a field-programmable gate array (FPGA) (Spartan3E, Xilinx). Light from the three LEDs were combined together and put into a multimode optical fiber, and then focused onto DMD I (DlpVisionFly4500, Texas Instruments) for generating a wide-field spatially modulated sinusoidal illumination pattern coded over 256 gray levels in the x-direction with a spatial frequency of $f_{x} = 0.1 {mm}^{- 1}$ at a medium surface. Reflected light from the medium surface were again spatially modulated by DMD II (DlpVisionFly4500, Texas Instruments) with specified sampling patterns also coded over 256 gray levels. Each sampling pattern consisted of 64 × 64 pixels for a $40 \times 40 {mm}^{2}$ field of view (FOV), thus the pixel size of the captured spatially modulated reflectance images was $0.625 \times 0.625 {mm}^{2}$ . For each sampling pattern, the spatially integrated light from DMD II was fed into a time-resolved multispectral photomultiplier tube (PMT) (H7155-01, Hamamatsu Photonics, Japan), and the signals from the PMT were frequency-decoded for wavelength discrimination by a specially designed lock-in photon-counting module, also implemented in the FPGA. Specular reflection was avoided by illuminating the medium at a small angle of ~3 degrees to a normal incidence while the DMD II received light normally from the medium. The whole experimental setup was placed in a dark environment to shield the stray light.

Fig. 1 (a) schematic of the multi-wavelength single-pixel SFD imaging system based on lock-in photon counting detection. D1-D3: LED drivers; (b) configuration of the lock-in photon counting module.

Download Full Size | PDF

2.2. Frequency-decoding for discriminating multiple wavelengths by lock-in photon counting detection

The lock-in photon counting module for frequency-decoding consisted of three independent digital phase-sensitive detection (PSD) blocks, as shown in Fig. 1(b). The spatially integrated reflectance detected by the PMT was output as discrete single-electron response (SER) pulses. For the wavelength of l, the occurrence probability of each SER pulse is proportional to the spatially integrated reflectance [19]

S {}^{(l)}{(t)} \propto \frac{4 R^{(l)}}{π} \sum_{n = 1}^{\infty} \cos [(2 n - 1) 2 π f^{(l)} t + θ^{(l)}],

where

R^{(l)}

is the spatially integrated reflectance to be measured and

θ^{(l)}

is the initial phase. To extract

R^{(l)}

and to eliminate the effect of

θ^{(l)}

, a pair of in-phase and quadrature reference signals were introduced into each digital PSD block, i.e.,

H^{(l)} (t) = \frac{4}{π} \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n - 1}}{2 n - 1} \cos [(2 n - 1) 2 π f^{(l)} t]

and

Q^{(l)} (t) = \frac{4}{π} \sum_{n = 1}^{\infty} \frac{{(- 1)}^{n - 1}}{2 n - 1} \sin [(2 n - 1) 2 π f^{(l)} t]

, with the amplitudes of both the reference signals being “±1”.

Each digital PSD block operated using a multiple-period reference-weighted counting (RWC) strategy which simply accumulated the amplitudes of the reference signals at the occurrence of each SER pulse, i.e., “1” was added or subtracted according to the reference signals of “1” or “-1”, respectively. Accumulation over an integration time of $T$ also completed a low pass filtering process for harmonic components suppression [19,20]. The process of in-phase and quadrature multiple-period accumulations multiple-period accumulations based on the RWC strategy is illustrated in Fig. 2, the up and down arrows indicate addition and subtraction of “1”, respectively, for accumulation of the weights over the integration time of T, which is expressed mathematically by

Fig. 2 Diagram of in-phase and quadrature multiple-period accumulations multiple-period accumulations based on the RWC strategy with an integration time of T.

Download Full Size | PDF

X^{(l)} = \int_{0}^{T} S^{(l)} (t) H^{(l)} (t) d t \propto {\begin{cases} R^{(l)} (1 - 2 θ^{(l)} / π), 0 < θ^{(l)} < π \\ R^{(l)} (2 θ^{(l)} / π - 3), π < θ^{(l)} < 2 π \end{cases},

Y^{(l)} = \int_{0}^{T} S^{(l)} (t) Q^{(l)} (t) d t \propto {\begin{cases} R^{(l)} (- 2 θ^{(l)} / π), 0 < θ^{(l)} < π / 2 \\ R^{(l)} (2 θ^{(l)} / π - 2), π / 2 < θ^{(l)} < 3 π / 2 \\ R^{(l)} (- 2 θ^{(l)} / π - 4), 3 π < θ^{(l)} < 2 π \end{cases},

Equations (2) and (3) indicate $X^{(l)}$ or $Y^{(l)}$ is proportional to $R^{(l)}$ with the proportionality constant as a function of $θ^{(l)}$ . A quantity proportional to the spatially integrated reflectance, $R^{(l)}$ , is extracted by

R^{(l)} \propto | X^{(l)} | + | Y^{(l)} |,

Figure 3 shows the schematic of the digital PSD block based on the RWC strategy using a pair of parallel accumulators for performing the in-phase and quadrature multiple-period accumulations of the weights over the integration time.

Fig. 3 Schematic of the digital PSD block based on the RWC strategy

Download Full Size | PDF

2.3. Recovery of spatially modulated reflectance images with the SPI method based on 2D DCT

The SPI method based on 2D DCT, i.e., the DCT-SPI scheme, is applicable to measurements in the DCT domain, where most images are sparse. DCT coefficients of a sparse image are acquired from the measurements, and the image is reconstructed by applying an inverse 2D DCT to the acquired DCT coefficients [24,28]. The frequency-decoded results, $R^{(l)}$ , obtained in the previous subsection are the DCT coefficients of the spatially modulated reflectance image as functions of the generalized frequency pair of $(u, v)$ of the sampling patterns, being expressed as

R^{(l)} (u, v) = \iint_{Ω} I^{(l)} (x, y) P (x, y; u, v) d x d y,

where

(x, y)

denotes the location on the medium surface,

Ω

is the domain in the x-y plane,

I^{(l)} (x, y)

is the spatially modulated reflectance image at the wavelength

l

,

P (x, y; u, v)

is the sampling pattern with

64 \times 64

pixels for the

40 \times 40 {mm}^{2}

FOV given by

P (x, y; u, v) \frac{2}{\sqrt{M N}} C_{u} C_{v} \cos \frac{(2 x / 0.625 + 1) u π}{2 M} \cos \frac{(2 y / 0.625 + 1) v π}{2 N},

where

M = N = 64

,

C_{0} = 1 / \sqrt{2}

and

C_{σ} = 1

(

σ = u

or

v

),

u = 0, 1, ..., M - 1

and

v = 0, 1, ..., N - 1

. Note that the number of DCT coefficients is the same as that of the pixels in the spatially modulated reflectance image. The sampling patterns were generated by varying

u

or/and

v

, and the corresponding spatial frequency pair of each sampling pattern,

(f_{x, s a m}, f_{y, s a m})

, were calculated by

f_{x, s a m} = u / 1.25 M {mm}^{- 1}

and

f_{y, s a m} = v / 1.25 N {mm}^{- 1}

from the two cosine functions in Eq. (6). The magnitudes of

P (x, y; u, v)

given by Eq. (6) take both negative and positive values that cannot be physically implemented into the sampling patterns generated by the DMD II. Thereby, we split

P (x, y; u, v)

into two parts,

P (x, y; u, v) = P_{+} (x, y; u, v) - P_{-} (x, y; u, v)

with positive and negative parts,

P_{+} (x, y; u, v) = [| P (x, y; u, v) | + P (x, y; u, v)] / 2

and

P_{-} (x, y; u, v) = [| P (x, y; u, v) | - P (x, y; u, v)] / 2

, respectively. Then Eq. (5) is modified to

R^{(l)} (u, v) = R_{+}^{(l)} (u, v) - R_{-}^{(l)} (u, v),

where

R_{+}^{(l)} (u, v) = \iint_{Ω} P_{+} (x, y; u, v) I^{(l)} (x, y) d x d y

and

R_{-}^{(l)} (u, v) = \iint_{Ω} P_{-} (x, y; u, v) I^{(l)} (x, y) d x d y

are obtained from the sampling pattern,

P (x, y; u, v)

, and the measured image,

I^{(l)} (x, y)

, after frequency-decoding.

Equation (6) indicates that fully sampling of one spatially modulated reflectance image with 64 × 64 pixels needs 2 × 64 × 64 measurements of the spatially integrated reflectances of $R_{+}^{(l)} (u, v)$ and $R_{-}^{(l)} (u, v)$ . However, according to the 2D DCT theory applied to images consisting of DCT coefficients obtained from the sampling patterns with different frequencies, the number of measurements can be reduced if the DCT coefficient at each frequency is known in advance. As the sinusoidal illumination with the AC (spatial modulation of $f_{x} = 0.1 {mm}^{- 1}$ in the x-direction) and DC components gives rise to the fluence-rate in the medium with the spatial frequency of $f_{x} = 0.1 {mm}^{- 1}$ [6], the major (or non-zero) DCT coefficients of the spatially modulated reflectance images were primarily concentrated in the frequency regions around $(f_{x, s a m}, f_{y, s a m}) = (0.1 {mm}^{- 1}, 0 {mm}^{- 1})$ and $(f_{x, s a m}, f_{y, s a m}) = (0 {mm}^{- 1}, 0 {mm}^{- 1})$ . Then we employed a total of 306 sampling patterns (3.74% of 2 × 64 × 64 sampling patterns, 153 for $P_{+} (x, y; u, v)$ and 153 for $P_{-} (x, y; u, v)$ ), which were generated from Eq. (6) by varying $u$ from 0 to 8 and $v$ from 0 to 16 corresponding to the ranges of $f_{x, s a m}$ and $f_{y, s a m}$ from $0 {mm}^{- 1}$ to $0.2 {mm}^{- 1}$ and from $0 {mm}^{- 1}$ to $0.1 {mm}^{- 1}$ , respectively, as shown in Fig. 4.

Fig. 4 Distribution of the spatial frequency pairs of the 306 sampling patterns.

Download Full Size | PDF

The spatially modulated reflectance image for the wavelength of l was recovered by applying the inverse 2D DCT operation to the measurements of $R_{+}^{(l)} (u, v)$ and $R_{-}^{(l)} (u, v)$ for the 306 sampling patters

I^{(l)} (x, y) = \sum_{u = 0}^{8} \sum_{v = 0}^{16} P (x, y; u, v) [R_{+}^{(l)} (u, v) - R_{-}^{(l)} (u, v)],

2.4. SFD-DOT reconstruction

The recovered spatially modulated reflectance image at a surface position $ρ = (x, y)$ for the wavelength of l, $I^{(l)} (ρ)$ , is a sum of the DC and AC components for the sinusoidal illumination pattern in the x-direction only as

I^{(l)} (ρ) = I_{D C}^{(l)} (ρ) + I_{A C}^{(l)} (ρ) \cos (2 π f_{x} x),

where

I_{D C}^{(l)} (ρ)

and

I_{A C}^{(l)} (ρ)

are the DC and AC amplitudes, respectively. Both

I_{D C}^{(l)} (ρ)

and

I_{A C}^{(l)} (ρ)

were extracted by the SSMD method making use of the orthogonality of harmonic functions with an excellent noise suppression capability [26,27]. The extracted

I_{D C}^{(l)} (ρ)

and

I_{A C}^{(l)} (ρ)

for the target medium as well as

I_{D C, r e f}^{(l)} (ρ)

and

I_{A C, r e f}^{(l)} (ρ)

for the homogeneous reference medium were the input data for SFD-DOT reconstruction.

For image reconstruction in the SFD-DOT, the DE of the fluence-rate, $Φ_{0}^{(l)} (r)$ , for a homogeneous medium at a position in the medium $r = (x, y, z)$ with $z$ of the depth from the medium surface, is often employed as a forward model. When the DE is combined with the modulated illumination source, the fluence-rate is also modulated as $Φ^{(l)} (r) = Φ_{0}^{(l)} (r) \cos (2 π f_{x} x)$ , where $Φ_{0}^{(l)} (r)$ is the AC amplitude of fluence-rate. Then the equation for $Φ_{0}^{(l)} (r)$ becomes a 1-D second-order Helmholtz equation with respect to the z-coordinate as the following equation [6,29]:

\frac{d^{2}}{d z^{2}} Φ_{0}^{(l)} (z) - {μ^{'}}_{e f f}^{(l)} Φ_{0}^{(l)} (z) = - \frac{q_{0}^{(l)} (z)}{κ^{(l) 2}},

where

{μ^{'}}_{e f f}^{(l)} = \sqrt{3 μ_{a}^{(l)} (μ_{a}^{(l)} + {μ^{'}}_{s}^{(l)}) + {(2 π f_{x})}^{2}}

is the effective attenuation coefficient for the SFD,

μ_{a}^{(l)}

and

{μ^{'}}_{s}^{(l)}

are the absorption and reduced scattering coefficients, respectively,

κ^{(l)} = {[3 (μ_{a}^{(l)} + {μ^{'}}_{s}^{(l)})]}^{- 1}

is the diffusion coefficient,

q_{0}^{(l)} (z) = P_{0}^{(l)} {μ^{'}}_{s}^{(l)} \exp (- κ^{(l)} z / 3)

is the amplitude of the source illumination, and

P_{0}^{(l)}

is the incident optical power.

Using the partial current boundary condition [30], $Φ_{0}^{(l)} (r)$ is calculated, and the flux of $Φ_{0}^{(l)} (r)$ at the surface, $Γ_{0}^{(l)} (ρ)$ , is given proportional to $Φ_{0}^{(l)} (r)$ at the surface as the following equation expressing the dependence on $f_{x}$ explicitly

Γ_{0}^{(l)} (f_{x}, ρ) = - A Φ_{0}^{(l)} (f_{x}, r) |_{z \to 0^{+}},

where

A = (1 - R_{e f f}) / (2 + 2 R_{e f f})

is the proportionality constant and

R_{e f f}

is the effective internal reflectance.

Our goal was to determine the absorption perturbation of a target medium with respect to the homogeneous reference medium, from difference measurements of the spatially varying DC and AC amplitudes. For SFD-DOT reconstruction, a linear integral equation is derived from the first-order Rytov approximation to the DE by considering the difference measurements [31]

Γ^{(l)} (f_{x}, ρ) = Γ_{0}^{(l)} (f_{x}, ρ) \exp [\frac{- \int_{Ω_{r}} δ μ_{a}^{(l)} (r) G_{0}^{(l)} (f_{x}, r, ρ) Φ_{0}^{(l)} (f_{x}, r) d r}{Φ_{0}^{(l)} (f_{x}, ρ)}],

where

δ μ_{a}^{(l)} (r)

is the absorption perturbation at position

r

in the medium,

Ω_{r}

is the domain in the x-y-z space;

G_{0}^{(l)} (f_{x}, r, ρ)

is the frequency-dependent Green’s function of the DE for a homogeneous reference medium [30]

G_{0}^{(l)} (f_{x}, r, ρ) = \frac{1}{4 π κ^{(l)}} [\frac{\exp (- {μ^{'}}_{e f f}^{(l)}) r_{1}^{(l)}}{r_{1}^{(l)}} - \frac{\exp (- {μ^{'}}_{e f f}^{(l)}) r_{2}^{(l)}}{r_{2}^{(l)}}],

where

r_{1}^{(l)} = \sqrt{{(z - z_{0}^{(l)})}^{2} + ρ^{2}}

and

r_{2}^{(l)} = \sqrt{{(z + z_{0}^{(l)} + 2 z_{b}^{(l)})}^{2} + ρ^{2}}

,

z_{0}^{(l)} = 1 / {μ^{'}}_{s}^{(l)}

,

z_{b}^{(l)} = 2 κ^{(l)}

, and

ρ

is the horizontal distance between the medium position

r

and the surface position

ρ

. The sensitivity of

\ln [Γ^{(l)} (f_{x}, ρ) / Γ_{0}^{(l)} (f_{x}, ρ)]

at

ρ

to the absorption perturbation of

δ μ_{a}^{(l)} (r)

at

r

is given as the term in the exponential in Eq. (13). The demodulated quantities of

I_{D C}^{(l)} (ρ) / I_{D C, r e f}^{(l)} (ρ)

and

I_{A C}^{(l)} (ρ) / I_{A C, r e f}^{(l)} (ρ)

are equivalent to the calculated quantities of

Γ^{(l)} (f_{0}, ρ) / Γ_{0}^{(l)} (f_{0}, ρ)

and

Γ^{(l)} (f_{x}, ρ) / Γ_{0}^{(l)} (f_{x}, ρ)

, respectively, where the superscript 0 refers to the case of

f_{x} = 0 {mm}^{- 1}

.

Therefore, a matrix equation relating the measured quantities and the unknown absorption perturbations is formulated as

B = W δ μ_{a}

where

B = \ln {[I_{D C}^{(l)} (ρ_{1}) / I_{D C, r e f}^{(l)} (ρ_{1}), ..., I_{D C}^{(l)} (ρ_{K}) / I_{D C, r e f}^{(l)} (ρ_{K}), I_{A C}^{(l)} (ρ_{1}) / I_{A C, r e f}^{(l)} (ρ_{1}), ..., I_{A C}^{(l)} (ρ_{K}) / I_{A C, r e f}^{(l)} (ρ_{K})]}^{T}

is a column vector consisting of the data obtained by the difference measurements between

I_{D C}^{(l)} (ρ)

and

I_{D C, r e f}^{(l)} (ρ)

, and between

I_{A C}^{(l)} (ρ)

and

I_{A C, r e f}^{(l)} (ρ)

at the

K = M \times N

pixels,

δ μ_{a} = {[δ μ_{a}^{(l)} (r_{1}), δ μ_{a}^{(l)} (r_{2}), ..., δ μ_{a}^{(l)} (r_{E})]}^{T}

is a column vector consisting of the values of

δ μ_{a}^{(l)} (r_{e})

of the e-th voxel (e = 1, 2, ..., E) in the mesh generated in the medium for numerical calculation of the DE,

W

is a sensitivity matrix with

2 K

rows and

E

columns written as

W = - [\begin{array}{l} G_{0}^{(l)} (f_{0}, r_{1}, ρ_{1}) V_{1} Φ_{0}^{(l)} (f_{0}, r_{1}) / Γ_{0}^{(l)} (f_{0}, ρ_{1}), .., G_{0}^{(l)} (f_{0}, r_{E}, ρ_{1}) V_{E} Φ_{0}^{(l)} (f_{0}, r_{E}) / Γ_{0}^{(l)} (f_{0}, ρ_{1}) \\ ⋮ \\ G_{0}^{(l)} (f_{0}, r_{1}, ρ_{K}) V_{1} Φ_{0}^{(l)} (f_{0}, r_{1}) / Γ_{0}^{(l)} (f_{0}, ρ_{K}), .., G_{0}^{(l)} (f_{0}, r_{E}, ρ_{K}) V_{E} Φ_{0}^{(l)} (f_{0}, r_{E}) / Γ_{0}^{(l)} (f_{0}, ρ_{K}) \\ G_{0}^{(l)} (f_{x}, r_{1}, ρ_{1}) V_{1} Φ_{0}^{(l)} (f_{x}, r_{1}) / Γ_{0}^{(l)} (f_{x}, ρ_{1}), .., G_{0}^{(l)} (f_{x}, r_{E}, ρ_{1}) V_{E} Φ_{0}^{(l)} (f_{x}, r_{E}) / Γ_{0}^{(l)} (f_{x}, ρ_{1}) \\ ⋮ \\ G_{0}^{(l)} (f_{x}, r_{1}, ρ_{K}) V_{1} Φ_{0}^{(l)} (f_{x}, r_{1}) / Γ_{0}^{(l)} (f_{x}, ρ_{K}), .., G_{0}^{(l)} (f_{x}, r_{E}, ρ_{K}) V_{E} Φ_{0}^{(l)} (f_{x}, r_{E}) / Γ_{0}^{(l)} (f_{x}, ρ_{K}) \end{array}],

where

V_{e}

is the volume of the e-th voxel. Equation (14) was reasonably solved with an algebraic reconstruction technology (ART), whose row-fashioned performance was nearly memory-independent and robust for the ill-posed linear inversion [32,33]. Finally, the multi-wavelength absorption tomographic images of the turbid media are reconstructed.

3. Experimental verification

3.1. System assessment

To verify the proposed lock-in photon counting detection and assess the effectiveness of the single-pixel SFD imaging system, a series of assessment was performed, including the detection linearity, the detection stability and the anti-crosstalk. A homogeneous slab phantom with a size of 70 mm (length) × 70 mm (width) × 30 mm (thickness) was used. The phantom was made of polyformaldehyde, and its absorption and reduced scattering coefficients, $μ_{a}$ and ${μ^{'}}_{s}$ , at the three wavelengths were determined to be $μ_{a} (455 / 532 / 660 nm) = 0.006 / 0.0049 / 0.004 {mm}^{- 1}$ and ${μ^{'}}_{s} (455 / 532 / 660 nm) = 1.132 / 0.99 / 0.8 {mm}^{- 1}$ , using time-resolved spectroscopy (TRS) [34]. The sampling pattern, $P_{+} (x, y; u, v)$ with $u = v = 0$ , was employed constantly in the system assessment process.

The detection linearity was assessed by changing the driver current of one LED from 400 mA to 1000 mA, at a step of 100 mA, and with its integration time of T = 1 s in lock-in photon counting detection, while the other two LEDs were powered off. Before the detection linearity assessment, the light intensities of the three LEDs with different driver currents were measured by a PD detector (PDA36A-EC, Thorlabs) with voltage as the unit. The frequency-decoding results were recorded and maximum normalized, and then evaluated for the detection linearity with the linear regression analyse. Figure 5 shows the normalized frequency-decoding results as functions of the light density (PD output voltage) for the three wavelengths, in contrast to their respective linearly regressed curves. The frequency-decoding results exhibit high linearity with the correlation coefficients of 0.999, 0.998, and 0.999 for the three wavelengths, respectively.

Fig. 5 Frequency-decoding results for the detection linearity assessment for the wavelengths of (a) 455 nm, (b) 532 nm and (c) 660 nm.

Download Full Size | PDF

To assess the detection stability of the single-pixel SFD imaging system, the driver currents of the three LEDs were set to 800 mA. The lock-in photon counting detection with the integration time of 2 s was repeated 120 times. Figure 6 shows the results of the lock-in frequency-decoding as a function of time with normalization by the mean values for highlighting their fluctuations. The results show the fluctuations of less than 0.5%, 0.5%, and 0.6% for the three wavelengths, respectively, demonstrating the detection stability of the system.

Fig. 6 Frequency-decoding results for the stability assessment of the lock-in photon counting detection.

Download Full Size | PDF

Although the chosen frequencies for frequency-encoding light sources avoided the odd harmonics of the two DMDs, harmonic crosstalks were inevitably occur in the frequency-decoding process, of which the effective solution is the choice of frequencies for frequency-encoding the three light sources. To assess the harmonic crosstalks among the three PSD blocks, a measure referred to as Crosstalk-Index (CI) was calculated, which is defined as the difference between the results of the lock-in frequency-decoding when only one LED was ON and those when the three LEDs were ON simultaneously. The integration time of the lock-in detection was 1 s, and the driver currents of the LEDs were varied from 400 mA to 1000 mA. The results of the CI are presented in Table 1 indicating that the crosstalks among these three digital PSD blocks were less than 4.72% and decreasing with the increase in the driver current, demonstrating the harmonic crosstalk was effective solved by the chosen three frequencies.

Table 1. Calculated CIs for the three wavelengths

View Table | View all tables in this article

In summary, the performances of our single-pixel SFD imaging system were confirmed.

3.2. Phantom experiments

We conducted a series of tissue phantom experiments to verify the proposed approach. A slab-shaped tissue phantom with a size of 70 mm (length) × 70 mm (width) × 30 mm (thickness) was fabricated using polyformaldehyde, as shown in Fig. 7(a). Inside the phantom, there was a cylindrical hole with a diameter of 10 mm and a height of 5 mm located at the depth from 1.5 mm to 6.5 mm, the hole was filled with the target liquid which was a mixture of India ink and Intralipid-10% diluted with deionized water in order to match background scattering. The optical properties of target liquid were determined to be $μ_{a}^{(T)}$ (455/532/660 nm) = $0.0125 / 0.0098 / 0.008 {mm}^{- 1}$ and ${μ^{'}}_{s}^{(T)} (455 / 532 / 660 nm) = 1.132 / 0.99 / 0.8 {mm}^{- 1}$ , also measured by the TRS [34,35]. In the SFD-DOT reconstruction process, the imaging region with 40 mm (length) × 40 mm (width) × 20 mm (thickness) was discretized into 163840 voxels with an equal size of $Δ x \times Δ y \times Δ z = 0.625 \times 0.625 \times 0.5 m m^{3}$ .

Fig. 7 Workflow of SFD-DOT reconstruction based on the recovered spatially modulated reflectance images: (a) sketch of the tissue phantom, (b) spatially modulated reflectance images for the three wavelengths recovered with the DCT-SPI scheme, (c) demodulated images of the DC (left) and AC (right) amplitudes for the three wavelengths, (d) horizontal cross-sections at z = 4 mm (left) and vertical cross-sections at y = 20 mm (right), from the SFD-DOT reconstructed tomographic images of the absorption coefficients for the three wavelengths. Dashed circles and rectangles indicate the true locations and sizes of the target.

Download Full Size | PDF

To achieve a difference measurement that is required for the system calibration, the homogeneous phantom described in Section 3.1 was used as the reference phantom after the system assessment. The sinusoidal illuminating pattern at the spatial frequency of 0.1 mm⁻¹ in the x-direction was generated, with all the drive currents of the three LEDs set to 950 mA. The 306 sampling patterns described in Section 2.3 were successively uploaded to DMD II for the lock-in photon counting detections with an integration time of T = 1 s. Figures 7(b)–7(d) present the workflow of the SFD-DOT reconstruction based on the recovered spatially modulated reflectance images. The spatially modulated reflectance images for the three wavelengths as shown in Fig. 7(b) were first obtained by frequency-decoding using the DCT-SPI scheme. These spatially modulated reflectance images were then demodulated with the SSMD method to obtain the DC and AC amplitudes of all the pixels, as shown in Fig. 7(c). Figure 7(d) depicts the multi-wavelength tomographic images of the absorption coefficients of the tissue phantom reconstructed by the SFD-DOT showing, the horizontal cross-sections at z = 4 mm and vertical cross-sections at y = 20 mm including the target center, and the corresponding x and z line-profiles are presented in Fig. 8. We employed a relaxation parameter of 0.5 in the ART linear inversion in which the iteration was terminated at 20 because after 20 iterations no evident improvement in the image quality was observed. In the tomographic images and line-profiles, the target location and the size can be reasonably identified. Comparing with the true position and line-profiles of the target, the peak positions of the reconstructed targets are slightly shifted toward the surface ( $z = 0 mm$ ), probably because the distributions of the sensitivity functions have maximal near the surface, decreasing with the increase in the depth, and because the geometry served as a spatial filter for the AC components in superficial region [1–3]. In this case, it is shown that the SFD imaging technique is well suited to biomedical applications where the ranges of lesions are limited up to several millimeters in thickness, such as for the skin tissue. It is also noted that there is a deviation between the true and reconstructed values. The reason mainly lies in the intrinsic ill-posedness and non-linearity of the inverse problem, and in the insufficiency of the information content in the featured data sets, as explained in Ref [33].

Fig. 8 x line-profiles of horizontal cross-sections at z = 4 mm (left) and z line-profiles of vertical cross-sections at y = 20 mm (right), from the SFD-DOT reconstructed tomographic images of the absorption coefficients of the tissue phantom for the wavelengths of (a) 455 nm, (b) 532 nm and (c) 660 nm.

Download Full Size | PDF

We compared the tomographic images obtained by our single-pixel SFD imaging system with those obtained by a conventional camera-based SFD imaging system. In the latter system, we replaced the single-pixel detector by a cooled EMCCD camera (Versarray512, Rolera-Mgi Plus, Canada) with detection by switching the wavelengths in turn. For each detection, the exposure time of EMCCD camera for capturing spatially modulated reflectance images is set to 50 ms. For a better comparation, we varied the absorption contrast of the target to the background from 1.5:1, 2:1, 2.5:1 to 4:1. The quality of the tomographic images obtained by the two SFD imaging systems were quantitatively evaluated using two metrics, referred to as the Quantitativeness Contrast (QC) and the Quantitativeness Ratio (QR). The $Q C$ is defined as the ratio of the mean value of the reconstructed absorption coefficient in the target region, $Ω_{T}$ , to that in the background region, $Ω_{B}$ , i.e., $Q C = \underset{r \in Ω_{T}}{m e a n} [μ_{a}^{(R)} (r)] / \underset{r \in Ω_{B}}{m e a n} [μ_{a}^{(R)} (r)]$ , and the $Q R$ is defined as the ratio of the mean value of the reconstructed absorption coefficient in the target region, $Ω_{T}$ , to the true target absorption coefficient, i.e., $Q R = \underset{r \in Ω_{T}}{m e a n} [μ_{a}^{(R)} (r)] / μ_{a}^{(T)}$ . The results are given in Tables 2 and 3 for the $Q C$ and $Q R$ values, respectively.

Table 2. Results of the QC for comparing the performances of the proposed single-pixel SFD imaging system and the conventional camera-based SFD imaging system.

View Table | View all tables in this article

Table 3. Results of the QR for comparing the performances of the proposed single-pixel SFD imaging system and the conventional camera-based SFD imaging system.

View Table | View all tables in this article

As expected, the two metrics were in good agreements between the two SFD imaging systems, the less-dense sampling strategy of the proposed DCT-SPI approach was proven to be effective in the single-pixel SFD imaging technique. The average differences in the QC and QR values between the two SFD imaging systems for the three wavelengths were 4.83%, 7.11% and 8.81% for the QC and 6.49%, 7.86% and 3.72% for the QR, respectively. The increase in the absorption contrast decreases the quantitativeness, mainly resulting from the saturation effect of the tomographic reconstruction, as explained before. Other factors that might adversely affect the imaging performance include the deviation of the optical properties measured by the TRS from the true values, and the instability of the LED power. An in-depth investigation is required for solving these issues.

4. Discussions

In the conventional SFD imaging, a pixel-wise analysis is performed with at least two spatial frequencies to simultaneously extract topographical images of the absorption and scattering coefficients [1–6]. This strategy might offer better spatial resolution than DOT but fails to obtain depth-resolved information. Also it is in principle less quantitative due to the application of the optically-homogeneous semi-infinite geometry in the photon-migration modeling. We herein aim at reconstructing tomographic images of the absorption coefficient by performing the model-based inversion from the whole data set, assuming that the scattering coefficient is known. The practicability of the assumption has been justified in some of the clinical applications such as intraoperative monitoring of tissue oxygenation, where tracking the dynamic variations in oxygen metabolism is a primary concern and requires fast and convenient determination of the absorption coefficient, while the scattering background is reasonably regarded as a constant and easy to be pre-determinated using the time resolved spectroscopy [36,37]. It is an open issue whether or not the SFD-DOT approach with multiple spatial frequencies could simultaneously reconstruct the tomographic images of both the absorption and scattering coefficients, regardless of the way how the reflectance images is acquired.

Thus far, 306 sampling patterns have been used for the recovery of the spatially modulated reflectance images as shown in Fig. 4, taking about an acquisition time of ~306 s for an integration time of T = 1 s. With consideration of the spatially low-passing effect of turbid medium on the light propagation, it might be feasible to accelerate the sampling process by reducing the sampling patterns as well as the integration time to some extent. Figure 9 illustrates the DOT reconstructions of the same phantom as in Fig. 7(a) at 455-nm wavelength with 130 sampling patterns under different integration times of T = 1 s, 0.5 s, and 0.25 s, respectively. The used 130 sampling patterns were generated from Eq. (6) by varying u from 0 to 5 and v from 0 to 13 that corresponds to $f_{x, s a m}$ and $f_{y, s a m}$ ranging from $0 {mm}^{-1}$ to $0.15 {mm}^{-1}$ and $0 {mm}^{-1}$ to $0.05 {mm}^{-1}$ , respectively. Accordingly, the acquisition time were determined to be 130 s, 65 s and 32.5 s. Figure 9 indicates that the target location and size can still be identified even if only 130 sampling patterns were employed for recovery of the reflectance image, whereas the decreasing integration time results in the substantial increase in reconstruction artifacts and degradation in the target quantitativeness, probably due to the decrease in the signal-to-noise ratio (SNR) of the measured data. We have also calculated the QC and QR values of the reconstructed images. The results show 3.7%, 11.2% and 12.7% reductions in the QC, and 4.3%, 12.7% and 14.4% reduction in the QR for the three acquisition times of 130 s, 65 s and 32.5 s, respectively, as compared to those for acquisition time of 306 s. Three alternative strategies can be potentially adopted in further investigations to accelerate the data acquisition: 1) the power of the LED sources can be potentially increased to improve the SNR of the photon counting detection and thereby reduce the integration time under an acceptable SNR; 2) the number of sampling patterns for the reflectance recovery could be further reduced on the basis of a dedicate analysis of the DOT resolution limit; 3) a compressed sensing method with more effective sparse representations of the reflectance images than DCT can be devised to decrease number of the sampling patterns. In terms of the above images, an SFD-DOT reconstruction can be acceptable with the proposed lock-in photon counting based single-pixel SFD imaging system using a less-dense sampling of 130 patterns (about 1.59% compression ratio) and an integration time of 0.25 s, achieving a minimum acquisition time of 32.5 s.

Fig. 9 SFD-DOT reconstructions with 130 sampling patterns with different acquisition times of 130 s, 65 s, and 32.5 s: (a) horizontal cross-sections at z = 4 mm (top) and vertical cross-sections at y = 20 mm (bottom); (b) corresponding x line-profiles (left) and z line-profiles (right). Dashed circles and rectangles indicate the true locations and sizes of the target.

Download Full Size | PDF

To demonstrate the benefits of the DCT-SPI scheme, we have compared, with reference to the EMCCD result, the spatially modulated reflectance images recovered with the DCT-SPI and the universal CS-SPI described in Ref [14], for the aforementioned homogeneous phantom under the same experimental conditions (130 sampling patterns, an integration time of T = 1 s and a spatial frequency of $f_{x} = 0.1 {mm}^{- 1}$ ). Figure 10 illustrates the spatially modulated reflectance images recovered by the DCT-SPI and CS-SPI at 455-nm wavelength, as well as their x line-profile comparison with the EMCCD result. Both the schemes decomposed the reflectance image into the DCT orthogonal basis for the sparse representation. In the DCT-SPI, the 130 known spatial frequencies were assumed in the DCT sampling patterns as described above. While the universal CS-SPI scheme generated the 130 sampling patterns at random, according to the Binary-Bernoulli distribution, and employed the Basis Pursuit for the l₁-regularized optimization solution [15]. The significant superiority of the DCT-SPI to the CS-SPI can be observed in terms of the image quality. A quantitative assessment has been performed by calculating the structural similarity indexes of the DCT-SPI and CS-SPI recovered images with regard to the EMCCD reference [38], which are 0.927 and 0.238, respectively. The results further render the proposed DCT-SPI scheme as a real-time potential scheme to capture high-fidelity spatially modulated reflectance images for DOT applications.

Fig. 10 Comparison of DCT- and CS-SPIs under the same experimental conditions: (a) spatially modulated reflectance images recovered with DCT-SPI (left) and CS-SPI (right) schemes; (b) corresponding x line-profiles.

Download Full Size | PDF

During the experiment, we have utilized the phantom and target liquids for the validation of the proposed SFD-DOT approach. Both the phantom and target liquids feature high scattering and low absorption, which fulfills the criteria of the DE forward model, i.e., the reduced scattering coefficient must be greater (≥10 × ) than the absorption coefficient [6]. However, the optical properties of biologically relevant chromophores are often strongly wavelength-dependent. At near-infrared wavelengths (~650-950 nm), biological tissue is highly scattering and lowly absorptive making scattering dominates absorption. At non near-infrared wavelengths, however, the absorption coefficients of biologically relevant chromophores are one or several orders of magnitude larger than those at near-infrared wavelengths, then the DE forward model may not appropriate for modeling photon transport [30,31]. To further apply the proposed SFD-DOT approach to biological tissue, it is desirable to model photon transport using the radiative transform equation or monte-carlo simulations, which will be valid for a greater range of the optical properties.

The maximum spatial frequency for the DE forward model to be valid is limited to $0.25 μ_{t r}$ or $0.33 μ_{t r}$ , where $μ_{t r} = μ_{a} + {μ^{'}}_{s}$ [6]. With such a frequency constraint, we have ever performed simulative and experimental reconstructions to investigate the performance of tomographic images obtained from different frequency components [39]. The investigation has concluded that using only one single frequency component for image reconstruction leaded to obvious shift phenomena, where the reconstructed region of interest moved to the surface. Multiple frequency components involved in reconstruction, however, significantly improved the image quality, as well as the spatial resolution and quantitativeness, since the increasing number of the frequency components provided more frequency information content of the sample and effectively reduced the ill-posedness of the inverse problem. Nonetheless, the image quality was determined by the choice of multiple frequency components. A larger distance among multiple frequency components in the Fourier space performed better in spatial resolution for image reconstruction. As for the target location and quantitativeness, lower spatial frequency components involved in image reconstruction performed better than higher frequency ones, because the formers were less attenuated and propagated deeper in the media. This paper has presented the tomographic images obtained from dual-frequency components of DC and AC (spatial modulation of $f_{x} = 0.1 {mm}^{- 1}$ ), our future work will demodulate more spatial frequency components simultaneously from a single spatially modulated reflectance image, and choose the most suited frequency components for SFD-DOT reconstruction to achieve the best results. The SSMD method, in principle, can be used for extracting multi-frequency DC and AC amplitudes from a single spatially modulated reflectance image containing more than two spatial frequencies with arbitrary illumination patterns. However, the performance of the SSMD method would be deteriorated because multiple frequency components in a single image may result in frequency crosstalk, and the precise demodulation of the blurred image is probably not achievable. To a great extent, the accuracy of demodulation also depends on the number of the pixels in a recovered spatially modulated reflectance image. In this paper, the SSMD method was adopted to extract dual-frequency components from a single image with 64 × 64 pixels. Future studies are required to evaluate the fidelity and performance of the SSMD method for extracting more multiple frequency components from a single image with more pixels.

It is worth mentioning that the presented single-pixel SFD imaging system has utilized a slight inclination angle of ~3 degrees for illumination to avoid specular reflections. Such an oblique illumination might phase-shifts the spatially modulated light during propagation in the medium as compared to that with the normal incidence, and leading to errors in the amplitude of fluence-rate calculated by Eq. (10), as demonstrated in Ref [40]. To investigate the effect of the oblique illumination on the accuracy of the calculated fluence-rate amplitude, a comparison was made by numerical simulations for a semi-infinite geometry with the optical properties of $μ_{a} = 0.004 {mm}^{- 1}$ and ${μ^{'}}_{s} = 0.8 {mm}^{- 1}$ , illuminated with a sinusoidal pattern at the spatial frequency of $f_{x} = 0.1 {mm}^{- 1}$ . Figure 11 shows the calculated amplitudes of fluence-rate varying with depth under the oblique and normal illuminations, respectively. The average difference in the amplitude of fluence-rate calculated between the two illuminations is 0.57%. The observation indicates that the slight inclination angle negligibly affects the fluence-rate calculation.

Fig. 11 Amplitudes of fluence-rate v.s. depth calculated for a semi-infinite geometry under oblique illumination and normal illumination of a sinusoidal pattern, respectively.

Download Full Size | PDF

In this study, we focused on the verification of the DCT-SPI scheme based on the lock-in photon counting detection and also focused on the evaluation of the relevant SFD-DOT reconstruction, for which the image reconstruction algorithm based on the first-order Rytov formulation with the ART linear inversion was uniquely employed in the phantom experiments. This algorithm has been proved to be robust to noises and accurate in determining displacements and certain gradients in the iteration process toward the solution [32,33], but the calculation result depends on the selection of relaxation parameters. An adaptive relaxation parameter adjustment method can be introduced for further enhancing the tomographic image quality that is beyond the scope of the study.

5. Conclusions

We have presented the SFD-DOT using the highly-sensitive single-pixel SFD imaging system for simultaneously acquiring multi-wavelength tomographic images of the turbid media. We utilized the square wave with different frequencies for frequency-encoding light sources, and implemented the novel lock-in photon counting detection for decoding the faint multi-wavelength light signals. The lock-in photon counting detection represents a significant advancement in simultaneously improving the detection sensitivity and the decoding accuracy. The number of the sampling patterns was significantly reduced from that necessary for the universal CS-SPI scheme by employing the 2D DCT theory. A series of the tissue-simulating phantom experiments were performed to verify the proposed SFD-DOT approach. The results showed good agreements between the tomographic images obtained by the proposed approach and those by using the conventional camera-based SFD imaging system. The proposed approach exhibits an advantage of simultaneous multi-wavelength imaging in SFD, which is useful for fast acquisition of multi- and hyperspectral images of biological tissues. Future studies will provide a more detailed insight into the factors that limit the acquisition speed. In addition, much effort will be made to facilitate the applications of this technique to real time monitor fast phenomena in biological tissue.

Funding

National Natural Science Foundation of China (81871393, 61475115, 61575140, 61475116, 81571723, 81671728, 81771880); Tianjin Municipal Government of China (16JCZDJC3120, 17JCZDJC32700, 17JCQNJC12700, 18JCYBJC29400).

Acknowledgments

The authors acknowledge the insightful suggestions and revisions given by Professor Yukio Yamada from University of Electro-Communications, Japan.

Disclosures

The author declares that there are no conflicts of interest related to this article.

References

1. D. J. Cuccia, F. Bevilacqua, A. J. Durkin, and B. J. Tromberg, “Modulated imaging: quantitative analysis and tomography of turbid media in the spatial-frequency domain,” Opt. Lett. 30(11), 1354–1356 (2005). [CrossRef] [PubMed]

2. S. D. Konecky, A. Mazhar, D. Cuccia, A. J. Durkin, J. C. Schotland, and B. J. Tromberg, “Quantitative optical tomography of sub-surface heterogeneities using spatially modulated structured light,” Opt. Express 17(17), 14780–14790 (2009). [CrossRef] [PubMed]

3. S. D. Konecky, C. M. Owen, T. Rice, P. A. Valdés, K. Kolste, B. C. Wilson, F. Leblond, D. W. Roberts, K. D. Paulsen, and B. J. Tromberg, “Spatial frequency domain tomography of protoporphyrin IX fluorescence in preclinical glioma models,” J. Biomed. Opt. 17(5), 056008 (2012). [CrossRef] [PubMed]

4. S. Bélanger, M. Abran, X. Intes, C. Casanova, and F. Lesage, “Real-time diffuse optical tomography based on structured illumination,” J. Biomed. Opt. 15(1), 016006 (2010). [CrossRef] [PubMed]

5. C. D’Andrea, N. Ducros, A. Bassi, S. Arridge, and G. Valentini, “Fast 3D optical reconstruction in turbid media using spatially modulated light,” Biomed. Opt. Express 1(2), 471–481 (2010). [CrossRef] [PubMed]

6. D. J. Cuccia, F. Bevilacqua, A. J. Durkin, F. R. Ayers, and B. J. Tromberg, “Quantitation and mapping of tissue optical properties using modulated imaging,” J. Biomed. Opt. 14(2), 024012 (2009). [CrossRef] [PubMed]

7. R. B. Saager, D. J. Cuccia, and A. J. Durkin, “Determination of optical properties of turbid media spanning visible and near-infrared regimes via spatially modulated quantitative spectroscopy,” J. Biomed. Opt. 15(1), 017012 (2010). [CrossRef] [PubMed]

8. S. V. Panasyuk, S. Yang, D. V. Faller, D. Ngo, R. A. Lew, J. E. Freeman, and A. E. Rogers, “Medical hyperspectral imaging to facilitate residual tumor identification during surgery,” Cancer Biol. Ther. 6(3), 439–446 (2007). [CrossRef] [PubMed]

9. H. Akbari, Y. Kosugi, K. Kojima, and N. Tanaka, “Detection and analysis of the intestinal ischemia using visible and invisible hyperspectral imaging,” IEEE Trans. Biomed. Eng. 57(8), 2011–2017 (2010). [CrossRef] [PubMed]

10. S. S. Welsh, M. P. Edgar, R. Bowman, P. Jonathan, B. Sun, and M. J. Padgett, “Fast full-color computational imaging with single-pixel detectors,” Opt. Express 21(20), 23068–23074 (2013). [CrossRef] [PubMed]

11. Q. Pian, R. Yao, L. Zhao, and X. Intes, “Hyperspectral time-resolved wide-field fluorescence molecular tomography based on structured light and single-pixel detection,” Opt. Lett. 40(3), 431–434 (2015). [CrossRef] [PubMed]

12. D. G. Winters, R. A. Bartels, and S. R. Domingue, “Hyperspectral imaging via labeled excitation light and background-free absorption spectroscopy,” Optica 2(11), 929–932 (2015). [CrossRef]

13. S. Jin, W. Hui, Y. Wang, K. Huang, Q. Shi, C. Ying, D. Liu, Q. Ye, W. Zhou, and J. Tian, “Hyperspectral imaging using the single-pixel Fourier transform technique,” Sci. Rep. 7(1), 45209 (2017). [CrossRef] [PubMed]

14. M. Torabzadeh, I. Y. Park, R. A. Bartels, A. J. Durkin, and B. J. Tromberg, “Compressed single pixel imaging in the spatial frequency domain,” J. Biomed. Opt. 22(3), 030501 (2017). [CrossRef] [PubMed]

15. E. J. Candès, J. Romberg, and T. Tao, “Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information,” IEEE Trans. Inf. Theory 52(2), 489–509 (2006). [CrossRef]

16. Y.-Z. Lu and R.-F. Lu, “Using composite sinusoidal patterns in structured-illumination reflectance imaging (SIRI) for enhanced detection of apple bruise,” J. Food Eng. 199, 54–64 (2017). [CrossRef]

17. Y. Wang and S. Zhang, “Comparison of the squared binary, sinusoidal pulse width modulation, and optimal pulse width modulation methods for three-dimensional shape measurement with projector defocusing,” Appl. Opt. 51(7), 861–872 (2012). [CrossRef] [PubMed]

18. L. Wang, Y. Yan, H. Ma, and Z. Jin, “Resonant fiber optic gyro based on a sinusoidal wave modulation and square wave demodulation technique,” Appl. Opt. 55(12), 3274–3278 (2016). [CrossRef] [PubMed]

19. W. Chen, X. Wang, B. Wang, Y. Wang, Y. Zhang, H. Zhao, and F. Gao, “Lock-in-photon-counting-based highly-sensitive and large-dynamic imaging system for continuous-wave diffuse optical tomography,” Biomed. Opt. Express 7(2), 499–511 (2016). [CrossRef] [PubMed]

20. J. M. Masciotti, J. M. Lasker, and A. H. Hielscher, “Digital lock-in detection for discriminating multiple modulation frequencies with high accuracy and computational efficiency,” IEEE Trans. Instrum. Meas. 57(1), 182–189 (2008). [CrossRef]

21. X.-M. Ding, B.-Y. Wang, D.-Y. Liu, Y. Zhang, J. He, H. J. Zhao, and F. Gao, “A three-wavelength multi-channel brain functional imager based on digital lock-in photon counting technique,” Proc. SPIE 10480, 104800S (2018).

22. A. Restelli, R. Abbiati, and A. Geraci, “Digital field programmable gate array-based lock-in amplifier for high performance photon counting applications,” Rev. Sci. Instrum. 76(9), 093112 (2005). [CrossRef]

23. D. Braun and A. Libchaber, “Computer-based photon-counting lock-in for phase detection at the shot-noise limit,” Opt. Lett. 27(16), 1418–1420 (2002). [CrossRef] [PubMed]

24. B.-L. Liu, Z.-H. Yang, and L.-A. Wu, “Coloured computational imaging with single-pixel detectors based on a 2D discrete cosine transform,” J. Mod. Opt. 64(3), 259–264 (2017). [CrossRef]

25. Z. Zhang, X. Ma, and J. Zhong, “Single-pixel imaging by means of Fourier spectrum acquisition,” Nat. Commun. 6(1), 6225 (2015). [CrossRef] [PubMed]

26. M. Xu, Z.-L. Cao, W.-H. Lin, X.-L. Chen, L.-F. Zheng, and B.-X. Zeng, “Single snapshot multiple frequency modulated imaging of subsurface optical properties of turbid media with structured light,” AIP Adv. 6(12), 125208 (2016). [CrossRef]

27. X. Chen, W. Lin, C. Wang, S. Chen, J. Sheng, B. Zeng, and M. Xu, “In vivo real-time imaging of cutaneous hemoglobin concentration, oxygen saturation, scattering properties, melanin content, and epidermal thickness with visible spatially modulated light,” Biomed. Opt. Express 8(12), 5468–5482 (2017). [CrossRef] [PubMed]

28. S. Dokouzyannis and T. Tziortzios, “High throughput and energy efficient two-dimensional inverse discrete cosine transform architecture,” IET Image Processing 7(5), 533–541 (2013). [CrossRef] [PubMed]

29. A. Mazhar, D. J. Cuccia, S. Gioux, A. J. Durkin, J. V. Frangioni, and B. J. Tromberg, “Structured illumination enhances resolution and contrast in thick tissue fluorescence imaging,” J. Biomed. Opt. 15(1), 010506 (2010). [CrossRef] [PubMed]

30. R. C. Haskell, L. O. Svaasand, T. T. Tsay, T. C. Feng, M. S. McAdams, and B. J. Tromberg, “Boundary conditions for the diffusion equation in radiative transfer,” J. Opt. Soc. Am. A 11(10), 2727–2741 (1994). [CrossRef] [PubMed]

31. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. 15(2), R41–R93 (1999). [CrossRef]

32. X. Intes, V. Ntziachristos, J. P. Culver, A. Yodh, and B. Chance, “Projection access order in algebraic reconstruction technique for diffuse optical tomography,” Phys. Med. Biol. 47(1), N1–N10 (2002). [CrossRef] [PubMed]

33. D. A. Boas, D. H. Brooks, E. L. Miller, C. A. DiMarzio, M. Kilmer, R. J. Gaudette, and Quan Zhang, “Imaging the body with diffuse optical tomography,” IEEE Signal Process. Mag. 18(6), 57–75 (2001). [CrossRef]

34. D.-L. Qin, H.-J. Zhao, Y. Tanikawa, and F. Gao, “Experimental determination of optical properties in turbid media by TCSPC technique,” Proc. SPIE 6434, 64342E (2007). [CrossRef]

35. H. J. Van Staveren, C. J. M. Moes, J. van Marie, S. A. Prahl, and M. J. van Gemert, “Light scattering in Intralipid-10% in the wavelength range of 400-1100 nm,” Appl. Opt. 30(31), 4507–4514 (1991). [CrossRef] [PubMed]

36. M. A. Yücel, C. M. Aasted, M. P. Petkov, D. Borsook, D. A. Boas, and L. Becerra, “Specificity of hemodynamic brain responses to painful stimuli: a functional near-infrared spectroscopy study,” Sci. Rep. 5(1), 9469 (2015). [CrossRef] [PubMed]

37. A. C. Ehlis, S. Schneider, T. Dresler, and A. J. Fallgatter, “Application of functional near-infrared spectroscopy in psychiatry,” Neuroimage 85(Pt 1), 478–488 (2014). [CrossRef] [PubMed]

38. Y.-B. Tong, Q.-S. Zhang, and Y.-P. Qi, “Image quality assessing by combining PSNR with SSIM,” J. Image Graphics 11(12), 1758–1763 (2006).

39. T. Li, Z. Qin, W. Chen, H. Zhao, P. Yan, K. Zhao, and F. Gao, “Wide-field fluorescence tomography with composited epi-illumination of multi-frequency sinusoidal patterns,” Appl. Opt. 56(29), 8283–8290 (2017). [CrossRef] [PubMed]

40. A. Bassi, D. J. Cuccia, A. J. Durkin, and B. J. Tromberg, “Spatial shift of spatially modulated light projected on turbid media,” J. Opt. Soc. Am. A 25(11), 2833–2839 (2008). [CrossRef] [PubMed]

Wavelength	LED driver current (mA)
Wavelength	400	500	600	700	800	900	1000
455 nm	2.13%	1.17%	2.15%	1.22%	0.55%	0.27%	0.8%
532 nm	3.34%	2.26%	1.39%	0.58%	0.97%	0.54%	0.54%
660 nm	4.72%	1.05%	1.60%	0.66%	0.75%	0.24%	0.51%

Absorption contrast	455nm		532nm		660nm
Absorption contrast	Single-pixel	Camera-based	Single-pixel	Camera-based	Single-pixel	Camera-based
1.5:1	1.277	1.178	1.202	1.152	1.162	1.143
2:1	1.324	1.206	1.257	1.243	1.312	1.227
2.5:1	1.759	1.793	1.618	1.517	1.616	1.698
4:1	2.199	2.218	1.955	2.290	1.775	2.236

Absorption contrast	455nm		532nm		660nm
Absorption contrast	Single-pixel	Camera-based	Single-pixel	Camera-based	Single-pixel	Camera-based
1.5:1	78.2%	72.7%	77.9%	72.9%	73.6%	71.7%
2:1	69.1%	66.0%	65.0%	63.8%	68.4%	72.2%
2.5:1	54.6%	59.9%	58.2%	60.9%	63.7%	62.9%
4:1	43.5%	45.6%	43.6%	51.7%	50.9%	53.7%

Wavelength	LED driver current (mA)
Wavelength	400	500	600	700	800	900	1000
455 nm	2.13%	1.17%	2.15%	1.22%	0.55%	0.27%	0.8%
532 nm	3.34%	2.26%	1.39%	0.58%	0.97%	0.54%	0.54%
660 nm	4.72%	1.05%	1.60%	0.66%	0.75%	0.24%	0.51%

Absorption contrast	455nm		532nm		660nm
Absorption contrast	Single-pixel	Camera-based	Single-pixel	Camera-based	Single-pixel	Camera-based
1.5:1	1.277	1.178	1.202	1.152	1.162	1.143
2:1	1.324	1.206	1.257	1.243	1.312	1.227
2.5:1	1.759	1.793	1.618	1.517	1.616	1.698
4:1	2.199	2.218	1.955	2.290	1.775	2.236

Multi-wavelength spatial frequency domain diffuse optical tomography using single-pixel imaging based on lock-in photon counting

Abstract

1. Introduction

2. Methods

2.1. Multi-wavelength single-pixel SFD imaging system based on lock-in photon counting detection

2.2. Frequency-decoding for discriminating multiple wavelengths by lock-in photon counting detection

2.3. Recovery of spatially modulated reflectance images with the SPI method based on 2D DCT

2.4. SFD-DOT reconstruction

3. Experimental verification

3.1. System assessment

3.2. Phantom experiments

4. Discussions

5. Conclusions

Funding

Acknowledgments

Disclosures

References

Cited By

Figures (11)

Tables (3)

Equations (15)

Optics Express