**32**(4), 418-425 (1993)].

## 1.

## Introduction

During functional brain activation, near-infrared spectroscopy (NIRS) is an effective tool for noninvasive investigation of cerebral oxygenation and hemodynamics.^{1, 2, 3} It has several advantages over other functional measurement methods such as functional magnetic resonance imaging (fMRI), positron emission tomography (PET), and electroencephalogram (EEG). These advantages include good temporal resolution, measurement of both oxy- and deoxyhemoglobin (HbO and HbR), high portability, low restraints of subjects, and cost-effective equipment.^{4} Several types of NIRS equipment are now available, such as time-resolved spectroscopy (TRS), phase-resolved spectroscopy (PRS), and continuous-wave-type (CW-type) instruments.

TRS-type and PRS-type instruments have some advantages over a CW-type instrument. For example, they can obtain mean optical path length of the measurement light, which can be used to realize a cross-talk-free hemoglobin concentration change estimation. However, since these instruments use a complex mechanism, they are usually more expensive than a CW-type instrument. On the other hand, the CW-type instrument is now very popular and widely used due to several practical reasons such as implementation cost. Thus, if we can estimate the mean optical path length ratio by using only CW-type measurement data, this contributes a realization of more accurate and reliable measurement of the CW-type instrument and should be very useful to many CW-type instrument users.

In a CW-type NIRS system, a simple calculation based on measuring the attenuation change of light as it propagates through a subject’s head provides information about the concentration change of HbO and HbR in the brain. This calculation requires multiple measurement lights of different wavelengths, and the dual-or triple-wavelength measurement is widely used. The relationship between attenuation and concentration change is expressed using a linear equation based on the modified Lambert-Beer law. The coefficients in this equation include molar absorption coefficients of HbO and HbR and mean optical path lengths of the measurement light. All of these coefficients depend on the wavelength of the measurement light. The molar absorption coefficients of HbO and HbR have been reported earlier.^{5} On the other hand, the mean optical path lengths are difficult to determine using a CW-type instrument.

If the mean optical path length is not known, we cannot quantify hemoglobin concentration changes, and comparisons of NIRS measurements of different subjects and/or regions within a subject are impossible since the mean optical path length depends on subjects and a measured region within a subject. The mean optical path length depends also on the wavelength of the measurement light. Thus, if we do not know path lengths of used wavelengths, and substitute, for example, a wavelength-independent constant for these values, it inevitably introduces a source of error while calculating concentration changes. These errors are referred to as cross talk, because a change in one of the chromophores may mimic a change in the other chromophore.^{6, 7}

Since it is difficult to determine the mean optical path lengths using a CW-type instrument, they were usually measured by a TRS or PRS instrument beforehand,^{8} or literature values^{9, 10, 11} were substituted for these values.^{12} In the former case, an additional TRS or PRS instrument is necessary, and these are usually more expensive than a CW-type instrument. In the latter case, the substituted values may be different from the true value since the path length depends on a subject and a measured region. Thus, if we have a method to estimate the path length only by a CW-type instrument, it could be useful to obtain the more accurate hemoglobin changes.

In this paper, we assume a triple-wavelength $({\lambda}_{1},{\lambda}_{2},{\lambda}_{3})$ measurement system and propose a new method for estimating the wavelength dependence of the mean optical path length. Here, we assume that only HbO and HbR contribute toward measured attenuation changes and that their concentrations change homogeneously in the tissue volume sampled by the detector. We also assume that optical path length does not change during a measurement.

Since two measurement lights with different wavelengths are sufficient to estimate two chromophores (HbO and HbR), two estimations of hemoglobin concentration changes can be obtained if we use a triple-wavelength measurement system—one calculated from data obtained using measurement lights ${\lambda}_{1}$ and ${\lambda}_{2}$ and the other from lights ${\lambda}_{2}$ and ${\lambda}_{3}$ . Theoretically, these estimations should be identical; however, they are generally not identical if we assume $l\left({\lambda}_{1}\right)=l\left({\lambda}_{2}\right)=l\left({\lambda}_{3}\right)$ in the calculation, where $l\left({\lambda}_{1}\right)$ , $l\left({\lambda}_{2}\right)$ , and $l\left({\lambda}_{3}\right)$ are the mean optical path lengths of three wavelengths. The proposed method estimates optical path length ratios, $l\left({\lambda}_{2}\right)\u2215l\left({\lambda}_{1}\right)$ and $l\left({\lambda}_{3}\right)\u2215l\left({\lambda}_{1}\right)$ , by minimizing squared differences of these two estimations. It cannot provide the absolute value of a mean optical path length; therefore, the absolute values of hemoglobin concentration changes are not obtained, and the method does not lead to the comparisons of NIRS measurements of different subjects and/or regions within a subject. However, accurate hemoglobin concentration changes free from cross talk with an arbitrary scaling can be obtained by the method.

There have been several studies on the wavelength-dependent optical path length estimation. Kahl reported an estimation algorithm based on the assumption that the pulselike heartbeat response in observed absorbance changes are induced by volume changes in the arterial blood, which is saturated with oxygen.^{13} Sakaguchi estimated the wavelength dependence of mean optical path length at the exposed cortex of animals based on the same principle^{14} as ours. However, their method did not take into account the measurement noise. Our simulation shown in this paper indicates that measurement noise significantly influences the estimation accuracy. We propose in this paper a new estimation method that provides an accurate estimation even with a relatively high noise level. The proposed method works without knowing a noise level. We think that this an important aspect of practical estimation, because sometimes the signal-to-noise ratio (SNR) of NIRS measurements of the human head is very low.

## 2.

## Theory

When a uniformly turbid medium is irradiated with near-infrared light, a temporal absorbance change $\Delta A$ resulting from a small homogeneous change in the absorption coefficient $\Delta {\mu}_{a}$ can be represented by using the modified Lambert-Beer law:

where $i$ is the sampling index, $l\left(\lambda \right)$ is the mean optical path length at wavelength $\lambda $ , and $n(i,\lambda )$ is a measurement noise. The absorption coefficient change $\Delta {\mu}_{a}$ is given as## Eq. 2

$\Delta {\mu}_{a}(i,\lambda )={\epsilon}_{\mathrm{Hb}\mathrm{O}}\left(\lambda \right)\Delta \mathrm{Hb}\mathrm{O}\left(i\right)+{\epsilon}_{\mathrm{HbR}}\left(\lambda \right)\Delta \mathrm{HbR}\left(i\right),$The triple-wavelength measurement employs three light sources. If we use a vector representation of temporal attenuation changes, hemoglobin concentration changes, and measurement noises as

## Eq. 3

${\mathit{a}}_{\mathit{i}}={[\Delta A(i,{\lambda}_{1}),\Delta A(i,{\lambda}_{2}),\Delta A(i,{\lambda}_{3})]}^{T},$## Eq. 4

${\mathit{x}}_{\mathit{i}}={[\Delta \mathrm{Hb}\mathrm{O}\left(i\right),\Delta \mathrm{HbR}\left(i\right)]}^{T},$## Eq. 7

$L=\mathrm{diag}[l\left({\lambda}_{1}\right),l\left({\lambda}_{2}\right),l\left({\lambda}_{3}\right)],$## Eq. 8

$E=\left[\begin{array}{cc}{\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{1}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{1}\right)\\ {\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{2}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{2}\right)\\ {\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{3}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{3}\right)\end{array}\right].$## Eq. 9

${\stackrel{\mathrm{\u0303}}{\mathit{a}}}_{i}=LE{\stackrel{\mathrm{\u0303}}{\mathit{x}}}_{\mathit{i}}+{\stackrel{\mathrm{\u0303}}{\mathit{n}}}_{\mathit{i}},$## Eq. 10

${\stackrel{\mathrm{\u0303}}{\mathit{a}}}_{\mathit{i}}={\mathit{a}}_{\mathit{i}}-\overline{\mathit{a}},$## Eq. 11

${\stackrel{\mathrm{\u0303}}{\mathit{x}}}_{\mathit{i}}={\mathit{x}}_{\mathit{i}}-\overline{\mathit{x}},$## Eq. 12

${\stackrel{\mathrm{\u0303}}{\mathit{n}}}_{\mathit{i}}={\mathit{n}}_{\mathit{i}}-\overline{\mathit{n}}.$Here, we calculate two estimations of temporal concentration changes. The first is calculated using the temporal absorbance change observed using wavelengths ${\lambda}_{1}$ and ${\lambda}_{2}$ , and the second is calculated using wavelengths ${\lambda}_{2}$ and ${\lambda}_{3}$ . These estimations should coincide with each other except for the noise, if we can correctly predict the unknown path length matrix $L$ . If we write

## Eq. 13

${\left[\begin{array}{cc}{\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{1}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{1}\right)\\ {\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{2}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{2}\right)\end{array}\right]}^{-1}=\left(\begin{array}{cc}{u}_{11}& {u}_{12}\\ {u}_{21}& {u}_{22}\end{array}\right),$## Eq. 14

${\left[\begin{array}{cc}{\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{2}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{2}\right)\\ {\epsilon}_{\mathrm{Hb}\mathrm{O}}\left({\lambda}_{3}\right)& {\epsilon}_{\mathrm{HbR}}\left({\lambda}_{3}\right)\end{array}\right]}^{-1}=\left(\begin{array}{cc}{v}_{11}& {v}_{12}\\ {v}_{21}& {v}_{22}\end{array}\right),$## Eq. 15

${\mathit{y}}_{\mathit{i}}=U{\widehat{L}}^{-1}{\stackrel{\mathrm{\u0303}}{\mathit{a}}}_{\mathit{i}},$## Eq. 16

${\mathit{z}}_{\mathit{i}}=V{\widehat{L}}^{-1}{\stackrel{\mathrm{\u0303}}{\mathit{a}}}_{\mathit{i}},$## Eq. 19

$\widehat{L}=\mathrm{diag}[\widehat{l}\left({\lambda}_{1}\right),\widehat{l}\left({\lambda}_{2}\right),\widehat{l}\left({\lambda}_{3}\right)].$If we can correctly predict the path length matrix $L$ , two estimations should be identical. Thus, we can give the following criterion $J\left(\widehat{L}\right)$ .

## Eq. 20

$J\left(\widehat{L}\right)=\frac{1}{N}\sum _{i=1}^{N}{\Vert {\mathit{y}}_{\mathit{i}}-{\mathit{z}}_{\mathit{i}}\Vert}^{2}.$Using Eq. 9, ${\mathit{y}}_{\mathit{i}}$ and ${\mathit{z}}_{\mathit{i}}$ are represented as follows:

## Eq. 21

${\mathit{y}}_{\mathit{i}}={U}_{1}{\stackrel{\mathrm{\u0303}}{\mathit{x}}}_{\mathit{i}}+{U}_{2}{\stackrel{\mathrm{\u0303}}{\mathit{n}}}_{\mathit{i}},$## Eq. 22

${\mathit{z}}_{\mathit{i}}={V}_{1}{\stackrel{\mathrm{\u0303}}{\mathit{x}}}_{\mathit{i}}+{V}_{2}{\stackrel{\mathrm{\u0303}}{\mathit{n}}}_{\mathit{i}},$## Eq. 25

$J\left(\widehat{L}\right)=\frac{1}{N}\sum _{i=1}^{N}{\Vert ({U}_{1}-{V}_{1}){\stackrel{\mathrm{\u0303}}{\mathit{x}}}_{\mathit{i}}+({U}_{2}-{V}_{2}){\stackrel{\mathrm{\u0303}}{\mathit{n}}}_{\mathit{i}}\Vert}^{2},$## Eq. 26

$=\mathrm{tr}\left\{({U}_{1}-{V}_{1}){\Sigma}_{xx}{({U}_{1}-{V}_{1})}^{T}\right\}+2\mathrm{tr}\left\{({U}_{1}-{V}_{1}){\Sigma}_{xn}{({U}_{2}-{V}_{2})}^{T}\right\}+\mathrm{tr}\left\{({U}_{2}-{V}_{2}){\Sigma}_{nn}{({U}_{2}-{V}_{2})}^{T}\right\},$## Eq. 27

$J\left(\widehat{L}\right)=\mathrm{tr}\left\{({U}_{1}-{V}_{1}){\Sigma}_{xx}{({U}_{1}-{V}_{1})}^{T}\right\}+{\sigma}_{n}^{2}{\Vert {U}_{2}-{V}_{2}\Vert}^{2}.$If $\widehat{L}=L\u2215d$ , where $d=\Vert L\Vert $ , we have ${U}_{1}={V}_{1}=dI$ , and the first term of Eq. 27 becomes zero. Thus, the criterion $J\left(\widehat{L}\right)$ achieves its minimum value when $\widehat{L}=L\u2215d$ if no noise exists $({\sigma}_{n}=0)$ . However, since the second term of Eq. 27 works as a bias of the criterion if noise does exist, $\widehat{L}$ that minimizes $J\left(\widehat{L}\right)$ is different from $L\u2215d$ . To avoid this difficulty, we propose the following criterion ${J}^{\prime}\left(\widehat{L}\right)$ :

## Eq. 28

${J}^{\prime}\left(\widehat{L}\right)=\frac{J\left(\widehat{L}\right)}{{\Vert {U}_{2}-{V}_{2}\Vert}^{2}}.$## Eq. 29

${J}^{\prime}\left(\widehat{L}\right)=\frac{\mathrm{tr}\left\{({U}_{1}-{V}_{1}){\Sigma}_{xx}{({U}_{1}-{V}_{1})}^{T}\right\}}{{\Vert {U}_{2}-{V}_{2}\Vert}^{2}}+{\sigma}_{n}^{2}.$## 3.

## Method

## 3.1.

### Simulation of Path Length Estimation

We give here a simple simulation of the proposed path length estimation method on artificial data to show the effectiveness of the algorithm. The artificial data used in the simulation was prepared as follows. Real observed data of oxy- and deoxyhemoglobin concentration changes, which were calculated from the NIRS measurement data of 3600 samples of an adult forehead by using an ordinary estimation algorithm of modified Lambert-Beer law, were used for the hemoglobin concentration change ${\mathit{x}}_{\mathit{i}}$ . Then, ${\mathit{x}}_{\mathit{i}}$ was transformed into the temporal absorbance change ${\mathit{a}}_{\mathit{i}}$ as follows:

where $E$ is given by assuming the wavelengths of three light sources are 776, 809, and $850\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$ ( ${\lambda}_{1}$ , ${\lambda}_{2}$ , ${\lambda}_{3}$ , in this order), and $L$ is assumed to be diag(1.0, 0.9, 0.8).Noises were assumed to be mutually independent Gaussian noises with a zero mean. The noise variance ${\sigma}_{n}^{2}$ was defined as follows to give several SNR levels:

## Eq. 31

$\mathrm{SNR}=10\phantom{\rule{0.2em}{0ex}}{\mathrm{log}}_{10}\phantom{\rule{0.2em}{0ex}}\frac{\text{mean}\phantom{\rule{0.3em}{0ex}}\text{variance}\phantom{\rule{0.3em}{0ex}}\text{of}\phantom{\rule{0.3em}{0ex}}LE\mathit{x}}{{\sigma}_{n}^{2}}.$The optical path length estimations were performed by two methods: the simple method using the criterion $J$ of Eq. 20, and the proposed method using the criterion ${J}^{\prime}$ of Eq. 28. For each SNR level, 50 different noise sequences were generated, and path length estimations were performed on these data. The average and the standard deviation of estimated optical path length ratios, $l\left({\widehat{\lambda}}_{2}\right)\u2215l\left({\widehat{\lambda}}_{1}\right)$ and $l\left({\widehat{\lambda}}_{3}\right)\u2215l\left({\widehat{\lambda}}_{1}\right)$ , were calculated.

## 3.2.

### Path Length Estimation by Experimental NIRS Data

To prove the validity of our theoretical considerations, we tested the proposed path length estimation method by applying it to real data of NIRS measurements. We used a near-infrared oximeter (NIRO-200, Hamamatsu Photonics) and its standard probes. The probes are attached on the right and left sides of a subject’s forehead (see Fig. 1 ). Each probe has two detectors of different source–detector distances (36 and $42\phantom{\rule{0.3em}{0ex}}\mathrm{mm}$ ) to realize a spatially resolved spectroscopy (SRS). The subject sat still on a chair, and temporal absorbance changes were recorded for $10\phantom{\rule{0.3em}{0ex}}\mathrm{min}$ . Two healthy volunteers participated in the experiment. This study was approved by the Institutional Review Board of AIST, and the participants gave a written informed consent. The measurement sampling frequency was set at $6\phantom{\rule{0.3em}{0ex}}\mathrm{Hz}$ for 3600 samples. Three light sources with wavelengths of 776, 809, and $850\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$ ( ${\lambda}_{1}$ , ${\lambda}_{2}$ , ${\lambda}_{3}$ , in this order) were used. Very low frequency components $(f<0.005\phantom{\rule{0.3em}{0ex}}\mathrm{Hz})$ of the observed data were removed by fast Fourier transform (FFT) to eliminate the drift. Data recorded in the initial and final two minutes of the analysis were not used in calculation to avoid side-effects of low-frequency elimination. As a result, only 2160 samples were used in the estimation.

We applied two pathlength estimation methods—the proposed method based on criterion ${J}^{\prime}$ and the simple method based on criterion $J$ —to the obtained NIRS data.

## 4.

## Results and Discussion

The simulation results are shown in Fig. 2 . This figure shows that the simple method can give correct path length estimations only with a very high SNR. On the other hand, on average, the proposed method gives very accurate estimations even with a low SNR. The figure also indicates, however, that standard deviations of the estimated value distribution of the proposed method are not small with a low SNR. Thus, the estimated values widely vary around the average value when SNR is very low. Thus, in this case, the proposed method may fail to give accurate path length estimation, although it still gives better results than the simple method. To overcome this problem, we may need more data samples.

The estimated mean optical path length ratios calculated using the experimental data of NIRS measurement at the right and left side of the human subject’s forehead are shown in Table 1 . Since the probe used in the experiment has two detectors with their respective separation distances being 36 and $42\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$ from the source, the estimated ratios for both detectors are shown in the table. Table 1 shows that the estimated ratios, $l\left({\lambda}_{2}\right)\u2215l\left({\lambda}_{1}\right)$ and $l\left({\lambda}_{3}\right)\u2215l\left({\lambda}_{1}\right)$ , of the proposed method were approximately 0.94 and 0.83.

## Table 1

Estimated mean optical path length ratios at the left and right forehead of a human subject.

Source–detector distance 36mm | |||||
---|---|---|---|---|---|

Proposed method (J′) | Simple method (J) | ||||

l(λ2)∕l(λ1) | l(λ3)∕l(λ1) | l(λ2)∕l(λ1) | l(λ3)∕l(λ1) | ||

Subject 1 | Left forehead | 0.93 | 0.82 | 1.20 | 1.14 |

Right forehead | 0.97 | 0.85 | 1.27 | 1.20 | |

Subject 2 | Left forehead | 0.91 | 0.78 | 1.34 | 1.23 |

Right forehead | 0.89 | 0.78 | 1.23 | 1.13 | |

Source–detector distance $42\phantom{\rule{0.3em}{0ex}}\mathrm{mm}$ | |||||

Proposed method $\left({J}^{\prime}\right)$ | Simple method $\left(J\right)$ | ||||

$l\left({\lambda}_{2}\right)\u2215l\left({\lambda}_{1}\right)$ | $l\left({\lambda}_{3}\right)\u2215l\left({\lambda}_{1}\right)$ | $l\left({\lambda}_{2}\right)\u2215l\left({\lambda}_{1}\right)$ | $l\left({\lambda}_{3}\right)\u2215l\left({\lambda}_{1}\right)$ | ||

Subject 1 | Left forehead | 0.94 | 0.84 | 1.25 | 1.19 |

Right forehead | 0.99 | 0.86 | 1.34 | 1.26 | |

Subject 2 | Left forehead | 0.91 | 0.79 | 1.15 | 1.04 |

Right forehead | 0.97 | 0.86 | 1.15 | 1.06 |

Essenpreis
^{9} gave mean optical path lengths for wavelengths ranging from
$740\phantom{\rule{0.3em}{0ex}}\text{to}\phantom{\rule{0.3em}{0ex}}840\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$
by using picosecond laser pulses and a streak camera detector. The source–detector separation was
$40\phantom{\rule{0.3em}{0ex}}\mathrm{mm}$
. Table 2 summarizes the mean optical path lengths of seven subjects for the three wavelengths (776, 809, and
$850\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$
) based on Fig. 4(a) of the Essenpreis report. Since the mean optical path length at wavelength
$850\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$
was not mentioned in this report, this value was determined by extrapolating the other data. Path length ratios were calculated using these path lengths and are given in the right columns of Table 2. They are very similar to the estimated path length ratios in Table 1. Path length ratios estimated by the simple method are also given in the right columns of Table 1. All these values are greater than the path length ratios obtained by the proposed method. The differences between the corresponding values were approximately 0.3. This result agrees well with the simulation results (Fig. 2).

## Table 2

Measured mean optical path lengths and their ratios. These values were based on the report of Essenpreis (Ref. 9). The mean optical path length at wavelength 850nm indicated by an asterisk ( *) mark was determined by extrapolating the other data because it is not mentioned in the report.

Subject | Mean optical path length | Path length ratio | |||
---|---|---|---|---|---|

l(776) | l(809) | l(850)* | l(809)∕l(776) | l(850)∕l(776) | |

1 | 237 | 232 | 210 | 0.98 | 0.89 |

2 | 244 | 238 | 216 | 0.97 | 0.89 |

3 | 252 | 246 | 224 | 0.97 | 0.89 |

4 | 254 | 247 | 222 | 0.97 | 0.87 |

5 | 262 | 251 | 231 | 0.96 | 0.88 |

6 | 264 | 260 | 237 | 0.99 | 0.90 |

7 | 298 | 282 | 258 | 0.95 | 0.86 |

Figure 3 shows time course examples of calculated oxy- and deoxyhemoglobin concentration changes. These time courses form a part of the concentration changes at the left forehead of subject 1. Figures 3a and 3b show the oxy- and deoxyhemoglobin concentration changes calculated using the equal path length assumption $[l\left({\lambda}_{1}\right)=l\left({\lambda}_{2}\right)=l\left({\lambda}_{3}\right)]$ . The dark solid lines indicate concentration changes calculated using ${\lambda}_{1}$ and ${\lambda}_{2}$ measurement, while the gray solid lines are those obtained using ${\lambda}_{2}$ and ${\lambda}_{3}$ . These time courses do not completely agree with each other, particularly in the deoxyhemoglobin change. Figures 3c and 3d show oxy- and deoxyhemoglobin concentration changes calculated by using path length ratios estimated by the proposed method. The dark and gray solid lines also indicate time courses by ${\lambda}_{1}$ and ${\lambda}_{2}$ , and ${\lambda}_{2}$ and ${\lambda}_{3}$ , respectively. These figures show that these two estimations agree very well with each other.

## 5.

## Conclusion

We propose a new method of estimating wavelength-dependent mean optical path length ratios of measurement lights. The method is based on the principle that two possible estimations of hemoglobin concentration changes calculated using triple-wavelength measurement system data should be identical. The proposed method was applied to the experimental data of human subjects’ foreheads, and the estimated path length ratios were very similar to the literature values. The proposed method simply provides the ratio of path lengths and not their absolute values. The absolute value of the optical path length is required to compare observed responses at different positions.^{8} Our method is not applicable to this problem. However, the proposed method estimates accurate hemoglobin concentration changes that are free from cross talk.

## References

*in vivo*: techniques and applications from animal to man,” J. Biomed. Opt., 12 (5), 051402 (2007). https://doi.org/10.1117/1.2789693 1083-3668 Google Scholar