2.5 Hz sample rate time-domain near-infrared optical tomography based on SPAD-camera image tissue hemodynamics

: Time-domain near-infrared optical tomography (TD NIROT) techniques based on diffuse light were gaining performance over the last years. They are capable of imaging tissue at several centimeters depth and reveal clinically relevant information, such as tissue oxygen saturation. In this work, we present the very first in vivo results of our SPAD camera-based TD NIROT reflectance system with a temporal resolution of ฀116 ps. It provides 2800 time of flight source-detector pairs in a compact probe of only 6 cm in diameter. Additionally, we describe a 3-step reconstruction procedure that enables accurate recovery of structural information and of the optical properties. We demonstrate the system’s performance firstly in reconstructing the 3D-structure of a heterogeneous tissue phantom with tissue-like scattering and absorption properties within a volume of 9 cm diameter and 5 cm thickness. Furthermore, we performed in vivo tomography of an index finger located within a homogeneous scattering medium. We employed a fast sampling rate of 2.5 Hz to detect changes in tissue oxygenation. Tomographic reconstructions were performed in true 3D, and without prior structural information, demonstrating the powerful capabilities of the system. This shows its potential for clinical applications. Abstract: Time-domain near-infrared optical tomography (TD NIROT) techniques based on diﬀuse light were gaining performance over the last years. They are capable of imaging tissue at several centimeters depth and reveal clinically relevant information, such as tissue oxygen saturation. In this work, we present the very ﬁrst in vivo results of our SPAD camera-based TD NIROT reﬂectance system with a temporal resolution of ∼ 116 ps. It provides 2800 time of ﬂight source-detector pairs in a compact probe of only 6 cm in diameter. Additionally, we describe a 3-step reconstruction procedure that enables accurate recovery of structural information and of the optical properties. We demonstrate the system’s performance ﬁrstly in reconstructing the 3D-structure of a heterogeneous tissue phantom with tissue-like scattering and absorption properties within a volume of 9 cm diameter and 5 cm thickness. Furthermore, we performed in vivo tomography of an index ﬁnger located within a homogeneous scattering medium. We employed a fast sampling rate of 2.5 Hz to detect changes in tissue oxygenation. Tomographic reconstructions were performed in true 3D, and without prior structural information, demonstrating the powerful capabilities of the system. This shows its potential for clinical applications. oxygenation in vivo . Pioneer successfully reconstructed the shape and the location of the ﬁnger positioned in a diﬀusive medium. With the exposure time of only 0.4 s (2.5 Hz), Pioneer detected the changes in the hemoglobin concentrations in the ﬁnger during arterial and venous occlusions. These results show a great potential of Pioneer to be applied in clinical environments where the detection of rapid changes in hemodynamics is important.


Introduction
Near-infrared optical tomography (NIROT) is a powerful tool that enables imaging biomarkers in tissue, such as the tissue oxygen saturation. This parameter is clinically highly relevant, because it reflects the balance between oxygen consumption and oxygen supply. Additionally, NIROT is capable of determining the presence and the concentration of tissue components such as oxy-, deoxy-and total hemoglobin ([O 2 Hb], [HHb], [tHb] in µM), water and lipids. Being non-invasive and relatively low cost, this technology will add a new dimension to biomedical imaging in research and clinical application [1][2][3][4][5]. NIROT relies on near-infrared light penetrating deep into tissue and interacting with different tissue components. The tomogram of the inner structure is reconstructed by use of the light that is detected at the surface after it travelled through the tissue. Among the three NIROT modalities, namely continuous wave (CW), frequency domain (FD) and time domain (TD), TD provides the richest information about the tissue. It employs time-of-flight (ToF) detectors and benefits from the ability to acquire the temporal distribution of the photons that have arrived [5].
Common implementations of TD NIROT employ a limited number of sources and detectors (e.g. [6]), because these instruments are based on optical fibers and photomultiplier tubes [7,8] or silicon photomultipliers (SiPM) [9]. The low amount of information also limits the image quality. However, these systems are bulky, and it is not possible to scale the number of detectors up. Several research groups worked on new systems to solve this problem. For example, an imaging system based on single photon avalanche diodes (SPADs) achieved a large number of virtual sources and detectors by scanning an object's surface and measuring at many positions [10,11]. However, this system has to be applied in free space, which hinders the use in a clinical environment, where ambient light is present. A handheld reflection TD NIROT system based on a time-gated intensified CCD camera helped diminishing the influence of ambient light [12] and showed sensitivity to inclusions placed at various depths. Since it included only one wavelength, it was unable to determine the oxygenation and an accurate 3D image reconstruction was not attempted.
The recent progress in CMOS technology has brought a new detector type into the field: the SPAD array [13,14]. Featured with low noise and high dynamic range, it provides a previously unachievable density of detectors and thus promises to drastically increase the spatial resolution [15] and depth sensitivity [16] of TD NIROT. Although, recent publications underline the progress in NIROT, the results were achieved under specific conditions, i.e. in transmission mode or employing prior knowledge of both depth of the inhomogeneity and optical properties of the background [15]. Such conditions are unachievable in real clinical scenarios where measurements in reflection mode are commonly required and a priori information is limited.
The goal of this work is to present a new TD NIROT system called Pioneer and to report its performance in various tests. The system works in reflection mode. It is based on in-house developed CMOS SPAD 32 × 32 array with high sensitivity to NIR light [14,17]. Here, we also aim at showing the benefits of the CMOS SPAD technology in NIROT, namely the compact design, the high number of detectors and the ToF information. We discuss several complications due to the tremendous amount of information that the sensor generates and provide solutions. We performed various experiments in vitro to assess all stages of tomogram generation, starting from system validation and going all the way to the image reconstruction of the inner structure of an object. We proceeded to in vivo imaging and monitoring tissue oxygenation during an occlusion test. We demonstrated the ability to reconstruct images rapidly and accurately by employing state-of-the-art hardware with GPUs and massively parallelized computation [18]. We believe this work paves the way for clinical application of NIROT in the near future.

The novel TD NIROT system called Pioneer
The schematic of Pioneer and a photo of the system are shown in Fig. 1(a),(c). It employs a custom designed 32 × 32 pixel TD image sensor called Piccolo [17] based on SPADs with dynamic sharing of time-to-digital converters (TDCs) [19]. Piccolo is the key element of Pioneer, since it provides a large number of detectors and hence source-detector pairs. To take fully advantage of these, we implemented a novel concept of NIROT in Pioneer. In this and the following sections, we describe how we harvested the benefits.
Pioneer employs a supercontinuum laser SuperK Extreme EXR-15 with acousto-optic tunable filter (AOTF) (NKT, Denmark). The laser emits ∼5 ps pulses with a repetition rate of 80 MHz at multiple wavelengths selectable freely in the range of 650-900 nm. An optical fiber guides the light to one of 11 source positions, which is selected by a fiber optical switch (Laser Components, Germany). The use of a tunable pulsed laser enables us to achieve two goals: (i) it provides time-correlated light pulses required for TD signal acquisition, and (ii) the multispectral light enters the tissue at exactly the same location, in our case through the same optical fiber, as described later. This enables a wavelength based coupling correction as previously described in Jiang et al. [20]. Additionally, the ability to freely choose a wavelength makes the system flexible in terms of detecting different chromophores.
Optical fibers guide the light from the optical switch to the source ring. The latter is the only element of Pioneer in contact with an object or subject. It has a rigid structure enclosed in soft black biocompatible silicone and forms a ring with an inner diameter of 25 mm ( Fig. 1(b)). This aperture defines the field-of-view (FoV) of Piccolo. The source fibers are equally distributed around the FoV in a circle of 45 mm diameter. They are fixed in the source ring and terminated with transparent windows to guide the light to the tissue. Soft biocompatible silicone Silpuran 2400 (Wacker, Germany) is used for both the optical windows and the black ring enclosure to achieve gentle but tight contact with the tissue. The manufacturing process of the source ring is described in [21]. We used multi-mode optical fibers with a core diameter of 62.5 µm and 0.23 NA (Thorlabs, US) with FC/APC connectors. We employed a NIR lens SY110M with a short focal length of 1.67 mm and a relatively large aperture of F/1.8 (Theia Technologies, US) to collect light from the tissue and to project it onto the Piccolo sensor. This open space arrangement enables us to interrogate a relatively large FoV of 25 mm diameter with this tiny sensor array. Each SPAD pixel covers an area of 1.06 × 1.06 mm 2 of object's surface. Due to the short focal length and wide angle of view, the lens enables a working distance of only ∼6 cm and a distance from an object to the image plane of ∼12 cm. The space between the object and the lens is shielded by a metal tube from ambient light. It is coated inside with Vantablack (Surrey NanoSystems, UK) to eliminate internal reflections. The tube also acts as a mechanical interface between the camera and the source ring. The case is in-house designed and 3D-printed. It contains the SPAD array and all accompanying electronics. A laptop controls this hardware by an in-house developed C++ / Qt software. Data post processing and image reconstruction is performed on a separate workstation (128 GB RAM, GPU Nvidia Gtx Titan Xp).
Pioneer is applicable inside an incubator for preterm infants. We achieved this by compact probe design keeping all bulky components next to the cot on a movable rack. The case of 21 × 12.5 × 6 cm 3 contains the Piccolo camera and its electronics. The optical components of the probe, i.e., the lens, the Vantablack tube and the optical source ring form a cylinder with a diameter of 6.5cm and 11.5 cm length. The probe is connected by wires (SS-USB and power supply) and optical fibers, which are flexible enough to pass through the openings of the cot. The whole probe weighs ∼0.85 kg, which is light enough to be handled manually during data acquisition. The probe can be placed on a pillow and held to the baby's head. It does not require a fixation. Thus, Pioneer is able to perform tomography of the neonatal brain at the bedside.

Data analysis and image reconstruction
The Pioneer system generates a huge amount of data, i.e., 32 × 32 pixels of 256 time gates per wavelength per source location. For each pixel the ToFs of the arriving photons are assembled to a histogram inside the Piccolo camera. The histograms are transferred to a powerful workstation (Intel Xeon CPU E5-1620 v4, 128 GB RAM, Titan Xp GPU, Ubuntu 16.04) for image reconstruction. Dark photon counts are eliminated by subtracting the median value of the photon counts distributed over all time gates. The instrumental response function is removed by wavelength normalization (WN) [20] or calibration on a homogeneous phantom with known optical properties. The corrected TD data are used for the image reconstruction. Here, we apply the physical model based image reconstruction. It requires an accurate description of light transport in tissue, which is called the forward problem. Here, we employ the diffusion equation (DE), which is valid in biological tissue due to its high reduced scattering coefficient µ s ′ and low absorption coefficient µ a in the NIR. The time-dependent DE is expressed as follows: where r is the position vector, t the time, κ = [3( µ s ′ ) ] −1 the diffusion coefficient [22,23] and ∅(r, t) the photon fluence rate. q is the source term and c the speed of light in the medium. The boundary is described by Robin conditions [24]. Finite element methods (FEM) are commonly used to solve DEs for complex boundaries. We apply NIRFAST, an FEM based software, to model the forward and inverse problem [18]. To achieve a fast (minutes) 3D image reconstruction and to reduce memory consumption, we convert the TD data to FD at ω=100 MHz. The DE in FD is stated with fluence rate Φ(r, ω) and source term Q(r, ω) at modulated frequency of ω: [25]. In a general 3D image reconstruction, the aim is to retrieve the 3D distribution of the optical properties µ by a least-squares minimization: where the optimal values µ * are achieved by minimizing the quadratic difference of simulated Φ S (µ) and measured Φ M forward results. The Φ M is Fourier transformed fluence rate from ∅(r, t).
The modified Tikhonov method was used for regularization [18,26] and β is the regularization parameter. The optimal solution is reached by Newton's method from the initial distribution of optical properties µ 0 [27]. For the volumetric image reconstruction, a fine tetrahedral mesh of small node distance is employed to simulate accurate forward results. A second coarser mesh is created for the inversion to ease the ill-posedness of the problem. This is also computationally less expensive and therefore helps to achieve fast image reconstruction.
To recover accurate optical properties we apply a region-prior, i.e., prior information about the regions in the medium. Structural information drastically reduces the ill-posedness of the inverse problem, since only µ s ′ and µ a of a few regions instead of thousands of nodes are recovered [28][29][30]. This information is retrieved from the volumetric reconstruction. The inverse problem is solved with the MATLAB function lsqcurvefit [28] and the forward problem by NIRFAST. Given that Pioneer features a large number of source-detector pairs (∼2800 in our experiments) and temporal data of a reasonable signal-to-noise ratio (SNR), the inverse problem is less ill-posed, and this improves the accuracy of the reconstructed optical properties.

Validation of Pioneer by phantom experiments
To validate the performance of Pioneer, we performed measurements on a silicone phantom of 120 mm diameter by 60 mm thickness with a cylindrical inclusion of 10 mm diameter and 30 mm length at a depth of 15 mm. We used Silpuran 2420 two-component silicone (Wacker, Germany) to produce this phantom. We added a defined quantity of silicone-compatible white ink ELASTOSIL RAL 9010 (Wacker, Germany) to obtain the specified µ ′ s . µ a was adjusted by carbon black powder. Since carbon black has little wavelength dependence in the NIR range, we added another ink, such as ELASTOSIL Heliogruen RAL 6004 (Wacker, Germany), to obtain a wavelength dependent µ a of the inclusion. The manufacturing of heterogeneous silicone phantoms was described previously [20,[31][32][33].
We applied the Pioneer probe to the surface of the phantom. One of the 11 light sources was switched on at a time and a temporal histogram was acquired by each detector simultaneously. It took ∼25 s to obtain 2 × 10 8 photons in the FoV. The whole measurement at 11 sources and 3 wavelengths (689 nm, 802 nm and 840 nm) took ∼14 min. Prior to the measurement, the system was calibrated using a reference homogeneous phantom with known optical properties. The optical properties of the phantom were measured with a commercial device Imagent (ISS Inc., Champaign, IL, USA), as describe e.g. here [33]. The µ ′ s and µ a of the phantom and of the inclusion are shown in Table 1. The inclusion had a higher µ a compared to the bulk, whereas the µ ′ s was similar. The term 'bulk' refers to the volume of the phantom without the inclusion, i.e., its background. The cylindrical phantom with a dimension of 90 mm diameter and depth of 50 mm was modeled by a fine mesh of with an average node distance of 0.9 mm (231'977 nodes) for obtaining the forward results, i.e., to model the light propagation in the object. This mesh represented only the central superficial part and not the whole silicone phantom (Fig. 2). Limiting the simulation process to this region of interest (ROI) facilitates the computation. In a first analysis to obtain regional information, we applied a coarser mesh of 3'784 nodes to perform a no-prior reconstruction of the spatial distribution of the µ a in the phantom. The optical properties of the bulk were measured and set as initial guess for the whole volume. Since there is little variation of µ s ′ in tissue we reconstructed the µ a only. The µ s ′ was fixed to the bulk µ s ′ in the whole volume. We set the maximum number of iterations to 15 and the initial regularization parameter was 10.
The reconstructed 3D image was created and segmented. The root mean squared error (RMSE) and Dice similarity were calculated to quantify the 3D image quality. We denote the ground-truth image as the array µ a , where µ ai , i ∈ N = {1, 2, . . . , n} is the value in the i-th mesh node, and the reconstructed image asμ a andμ a i , i ∈ N. The RMSE is a measure of the error of the reconstruction: The smaller the RMSE value is, the smaller is the deviation of the reconstructed image from the ground-truth. The Dice similarity quantifies the similarity of the structure of the reconstructed 3D image to the ground-truth [34].
where S(µ a ) and S(μ a ) are the ground-truth and reconstructed segmented images, respectively. | · | denotes the cardinality of a set. The closer the Dice value is to 1, the higher is the agreement between reconstruction and ground-truth. In the second analysis, the aim is to quantitatively obtain the µ a and µ s ′ for all the three wavelengths: 689, 802, and 840 nm in the inclusion and bulk. This second analysis is based on the regional information known from the first analysis, i.e., a shape for the inclusion. The phantom provided us with a wavelength dependent µ a contrast of the inclusion, which varied from 6.7 × compared to the bulk at 689 nm to 2.5 × at 840 nm ( Table 1)

Imaging hemodynamics of a finger in a diffusive medium
We designed an experiment to prove that Pioneer is able to image internal structures and show hemodynamic information about a diffusive medium. We manufactured a silicone phantom with a cavity with optical properties similar to tissue. The outer shape of this phantom is similar to the phantom described above. The cavity has the shape of an index finger of the left-hand. This finger of a healthy volunteer (male adult, 35-year-old, 183 cm, 85 kg) was placed in the cavity and was tightly surrounded by the soft silicone. We also inserted a small amount of water-Intralipid emulsion with Indian ink into the cavity to exclude air gaps between the silicone and the skin. The optical properties of this liquid matched the ones of the silicone. In this experiment, we used the same probe as described above (Fig. 2), i.e., 11 light-sources in a circle with a diameter of 45 mm arranged around the FoV with a diameter of 25 mm, which is covered by 257 detectors. The Pioneer probe was placed in the center of the flat surface of the cylindrical phantom, with the finger located ∼ 11 mm below the surface (Fig. 3). In addition, we repeated the experiment and measured with the commercial oximeter OxyPrem 1.4 (OxyPrem AG, Switzerland). This study was cleared by the Ethical Committee of the Kanton Zurich, Switzerland.

Measurement protocol
To elicit controlled hemodynamic changes, we induced venous and arterial occlusions by inflating a cuff that was wrapped around the upper arm of the subject, who was sitting. Venous occlusion was achieved by rapidly inflating the cuff to 60 mmHg. We performed three repetitions of 2 min of venous occlusion and 2 min of rest. We also measured the baseline conditions for 20 s before and after the series of occlusions. At 60 mmHg, the venous outflow from the tissue is blocked, while the inflow of highly oxygenated arterial blood is unimpeded.   [35]. Arterial occlusion: A higher pressure of 250 mmHg was applied to stop the arterial flow. Since an arterial occlusion causes discomfort, the measurements were limited to two repetitions of occlusion for 5 min and recovery for 3 min. Again, the baseline measurements for 20 s were added before and after the occlusion measurements.
For Pioneer, a short exposure time of 400 ms was used per source position and wavelength (725 nm and 802 nm) to capture these rapid changes. These wavelengths were chosen for their sensitivity to hemoglobin and high SNR. The total time for a complete measurement (2 wavelengths, 11 source positions) was 10 s, which includes data acquisition and transfer to the PC. The OxyPrem instrument has a time resolution of 1s.

Image reconstruction for in vivo measurements by Pioneer
Here we employed a novel three-step 3D image reconstruction procedure. It does not require calibration on a homogeneous phantom and enables imaging without prior information about the structure of the object. As a first step, we measured the average optical properties with a long exposure time. As a second step, we performed a no-prior reconstruction of the internal structure of the object. As a third step, we segmented the volume based on this reconstruction and performed a reconstruction of the optical properties of the regions, using only the measured data acquired at short exposure time to obtain a high sample rate.
The first step was to determine the average optical properties. TD data were transformed into FD by a fast Fourier transform. The FD data Φ consist of the amplitude, i.e. ln(I · r 2 ), where I is the amplitude of the light intensity at 100 MHz and r the scalar source detector distance, and the phase θ [20]. The average µ a and µ s ′ were calculated by the multi-distance approach as described in detail in [36]. In short, it is based on the solution to the DE for the semi-infinite boundary condition under the assumption that the coupling from the tissue to the active area of the camera is the same for all pixels [37]. This theory shows that ln(I · r 2 ) ∝ r and θ ∝ r are both linear functions and that the µ a and µ s ′ can be calculated from the slopes of these lines. We determined these slopes in our data by a linear fit (MATLAB function fit). We calculated the average slope for 11 sources. This enabled us to determine the average µ a and µ s ′ constituting our initial guess for the 3D image reconstruction.
As second step a volumetric image was reconstructed. We accumulated photon counts for an exposure time of 25 s for each source and wavelength. We employed these data for reconstructing the internal structure. For this reconstruction, two cylindrical meshes with dimensions of 90 mm diameter and 40 mm thickness were created in NIRFAST. The fine mesh consisting of 185608 nodes was used for the forward calculation, and the coarser mesh consisting of 2274 nodes for the inverse problem. We applied the previously proposed wavelength-normalization (WN) approach to avoid a tedious calibration for the system response functions [20]. We reconstructed the image from the ratio of the FD data at two different wavelengths. The Tikhonov regularization with an initial parameter of β = 10 was applied for iterative optimization. The µ s ′ was assumed constant in the whole volume and only the distribution of µ a was recovered to obtain the shape of the finger in the phantom. The reconstructed volume was segmented by setting a threshold for µ a of twice the median µ a . The segmented mesh was utilized as region-prior for reconstructing µ a and µ s ′ of the finger, and to access the hemodynamic in vivo. As third step, the optical properties µ a and µ s ′ in the finger region were recovered using TD data with a much shorter acquisition time of only 400 ms. The image reconstruction was based on region-priors, i.e., the anatomic information recovered in the second step. This enables the detection of rapid changes in the oxygenation of the finger.

Evaluation of the results of image reconstruction
For the evaluation of the 3D image reconstruction we only evaluated the structural information with Dice similarity, since the ground truth optical properties of the finger are unknown. We converted changes in the µ a into changes in [

Results of in-phantom validation
The Pioneer system measures the ToF of photons detected in all pixels of the Piccolo camera, i.e. all pixels are active at the same time. The onboard FPGA assembles this timing information from each pixel in histograms, the distribution of ToFs, which shows the number of detected photons per time bin. An example of such a response is shown in Fig. 4(a).
The typical distributions of the intensity and phase delay (i.e. ToF) over the circular FoV are shown in Fig. 4(c),(d). As expected, the intensity exponentially decreases with distance between the source and detector, whereas the phase delay increases linearly. Note that only the central part of the sensor was illuminated, thus some pixels appear outside of the FoV. There are a few 'hot' pixels in Piccolo, which are masked and therefore visible as pixels with missing data in the 3D surface plots (Fig. 4(c),(d)). 257 detectors provide valid data based on their SNR. This results in a total of ∼2'800 source-detector pairs.
The inner structure of the heterogeneous phantom was reconstructed from the acquired data. The optimization ended after 13 iterations when it reached the criterion of no further improvement. The whole process took <10 minutes. The location and shape of the inclusion were successfully recovered (Fig. 5(a)). The RMSE was 0.0023 mm −1 (Eq. (4)). After the segmentation, we obtained a high Dice similarity of 0.65 between the reconstructed image and the ground truth ( Fig. 6(a)). As second analysis, we used the segmented volume as region-prior for the reconstruction of the optical properties of the two regions. The reconstructed and ground-truth optical properties are shown in Fig. 5(b) in comparison with the initial guesses. In addition, it is seen that a better agreement appears for µ a than for µ s ′ because the contrast between the regions is higher in µ a .

Results of the in vivo occlusion tests
The 3D image reconstruction took 74 s for one iteration and the optimal solution was reached after the 24th iteration. The segmented region shows good agreement with the true finger shape ( Fig. 6(a),(b)). To calculate the Dice similarity, we converted the reconstructed distribution of the µ a into a region of voxels and compared it to the reference model. The reference model was constructed from a 3D-scan of the mold that we employed to create the phantom. The Dice similarity between reconstructed 3D image and the reference model was 0.83.
The hemodynamic in the finger was recorded at a sample rate of 2.5 Hz. The acquisition of one complete set of data took 10 s: 0.4s x 11 sources x 2 wavelengths. The optical properties of the finger region were reconstructed for every measurement. [O 2 Hb] and [HHb] were calculated from the absorption coefficients of hemoglobin at two wavelengths [18]. The concentrations were set to zero at the first data point The correlation coefficients between the inflation status (0 for rest and 1 for inflated) of the cuff and ∆[O 2 Hb] were −0.68 (p<1%) and 0.68 (p<1%) for the arterial and venous occlusion tests, respectively. The correlation coefficients between the inflation status and ∆[HHb] were 0.66 (p<1%) and 0.50 (p<1%) for the arterial and venous occlusion tests, respectively. These high correlations that the expected hemodynamic changes were consistently detected.

Discussion
In this paper, we described the TD NIROT system Pioneer, which is based on Piccolo camera. It has been shown that an increase in the number of source-detector pairs enables a higher spatial resolution of NIROT [15]. Thus, we designed the Pioneer system to take advantage of 1024 pixels of the ToF camera Piccolo, of 11 light sources and of a free choice of the wavelength in the NIR wavelength range. In this work, we used ∼2.8k source-detector pairs, only a fraction of the total number of available source-detector pairs. While Lyons et.al. demonstrated the potential of such a system in a specific 2D imaging task with a priori knowledge of the depth of a highly absorbing inclusion in transmission mode [15], we solved the more general problem of showing the 3D inner structure of an object in reflection mode and without priors.
The SPAD camera Piccolo enabled a novel concept of a TD NIROT system, which resulted in a very compact instrument. This differs fundamentally from previously developed NIROT systems, which used separate detectors, such as photomultiplier tubes, single SPADs or silicon photomultipliers, and optical fibers for emission and detection [6,8,9]. The novel Pioneer concept implements a scalable number of detectors, i.e., a larger SPAD camera with 36'288 pixels can be used [38]. In addition, the number of light sources may be expanded (up to 24 in our case). Notably, a similar arrangement of a ToF camera and light sources was employed by Sawosz et al. [12], but the image reconstruction was not performed. Our results demonstrate that it is possible to handle the abundance of data generated by the Pioneer system and to reconstruct reasonable images. Even dynamic measurements are feasible, which has been demonstrated by Jiang et al. in phantoms [39]. In this paper, we demonstrated the dynamic in vivo measurement with Pioneer for the first time. We employed repetitive measurements with a short exposure time of 400 ms to observe hemodynamic changes in an index finger during arterial and venous occlusions. The finger was located inside a silicone phantom with the optical properties similar to biological tissue. With the proposed 3D image reconstruction method, a volumetric image was obtained, and the finger region was correctly determined.
Furthermore, we demonstrated the capability of Pioneer to 3D image hemodynamic changes. With the reconstructed regional information, the optical properties were recovered at two wavelengths and converted to [ [HHb] actually only occur in the finger, we expected the OxyPrem to underestimate these changes. Since the Pioneer system provides images of the finger, we expected it to be more sensitive to the finger region and hence to yield larger concentration changes than OxyPrem. Possible explanations could be that i) the multi-distance approach of OxyPrem leads to a higher sensitivity to the tissue of the finger than being assumed by homogeneity, ii) OxyPrem assumes the optical properties of the neonate's head, while the finger has different optical properties and iii) the two instruments did not measure simultaneously. Still the similar shape of the time series obtained by the two instruments confirms that the reconstructed images provide reasonable physiological information.
The no-prior image reconstruction is known to be inaccurate in the reconstructed values due to its ill-posed and ill-conditioned nature [18]. The tomograms tend to be blurry if no prior information is provided. As mentioned above, we have applied a strong regularization on the scattering. However, no further prior information was supplied to the reconstruction algorithm. Due to the tremendous amount of the acquired data, we were able to achieve high localization accuracy. This shows again that an increase in the amount of the obtained data enables higher image quality.
The inverse problem becomes well-posed and well-conditioned for the region-prior reconstruction [28][29][30]. Here we demonstrated that an image reconstructed without priors can be obtained first and then be employed as anatomical prior for a quantitatively correct image reconstruction based on this region prior. The challenge of this approach is to achieve an accurate segmentation result. We used a threshold of twice the median µ a without priors for the segmentation and this yielded a reasonable structural 3D image. In a clinical environment, the region-prior can also be obtained from other imaging techniques, such as ultrasound, X-ray CT or MRI. The structure of an object can also be resolved in the future by machine learning methods [40,41]. The region-prior approach enabled us to reconstruct the optical properties of the two regions in the phantom at all tree wavelengths with a high precision.
We would like to emphasize, that we intended to challenge the Pioneer system and the reconstruction algorithms, and therefore intentionally set the initial guess values far from the correct ones. In the future, a more careful choice of an initial guess will facilitate the reconstruction process and is expected to result in an even more accurate solution.

Conclusion
To conclude, we reported in this paper the first in vivo results of our TD NIROT system Pioneer, which employs several cutting-edge technologies. Its Piccolo camera simultaneously counts photons in a large number of ToF detectors, which have a temporal resolution in the order of ps and a high sensitivity to near-infrared light. This enabled the collection of a massive amount of information from the object under investigation. The state-of-the-art reconstruction software enabled accurate and fast no-prior reconstruction of the structure of the object. Furthermore, we presented here the 3-step reconstruction approach: 1) The average optical properties of the object were determined and used as an initial guess. 2) The structure of the object was recovered without any priors. 3) This structural information was used to segment the object and to accurately determine the µ a and µ s ′ in each region at all 3 wavelengths. We employed this approach to image dynamic changes of tissue oxygenation in vivo. Pioneer successfully reconstructed the shape and the location of the finger positioned in a diffusive medium. With the exposure time of only 0.4 s (2.5 Hz), Pioneer detected the changes in the hemoglobin concentrations in the finger during arterial and venous occlusions. These results show a great potential of Pioneer to be applied in clinical environments where the detection of rapid changes in hemodynamics is important.