2.1.

Optical Fiber-Free Silicon Photomultiplier

An SiPM is basically an array of hundreds of thousands fully independent microcells with a very high fill factor (typically $>50\%$), each one consisting of a silicon SPAD, together with an integrated quenching resistor. Cells are connected in parallel, thus sharing a single common anode and a single common cathode. As a photon hits the detector active area, it generates an electron–hole pair inside the high electric field region (above the breakdown level). The pair triggers an avalanche current pulse, thus giving rise to a voltage pulse across the integrated quenching resistor and marking the photon arrival time. The block scheme of an SiPM-based detection channel is shown in Fig. 1(a), where the photodetector is directly coupled to an RF amplifier to provide a macroscopic voltage pulse with low additional time jitter, which can be easily used as a synchronism to a time-correlated single-photon counting (TCSPC) board. An external power supply is needed to provide the detector biasing and the amplifier operating voltage.

Fig. 1

(a) Block schematic of a detection channel, (b) designed SiPM-based detection channel (bottom in the picture) and the same detector, hosted into a plastic screw holder (top), and (c) EEG standard cap.

We designed two miniaturized printed circuit boards (PCBs) each one hosting a $1\times 1{\mathrm{mm}}^2$ SiPM (C30742-11-050, Excelitas Technologies, Canada) and all the required front-end and amplification electronic circuits, still keeping a very small PCB dimension [$30\mathrm{mm}\times 6\mathrm{mm}$, Fig. 1(b) bottom]. The probe was designed to be placed directly in contact with the head avoiding the use of optical fibers and to be hosted into a plastic screw holder [see Fig. 1(b) top] to be compatible with the EEG standard caps (g.EEGcap, g.tec medical engineering GmbH, Austria), shown in Fig. 1(c). In this way, the probe is very easy to handle and place. The large hole in the cap allows for easy removal of hairs.

2.2.

System Architecture

To test the performances of the SiPMs detectors, we integrated them into a two-channel TD fNIRS system, presented in Ref. 14. It employs two pulsed diode laser heads, operating at two wavelengths around the isosbestic point (690 and 829 nm), with 80-MHz repetition rate (PDL, Picoquant GmbH, Germany). Multimode GI fibers ($100/140\mu \mathrm{m}$) and a $2\times 2$ fiber optical switch are used to inject light pulses into the sample. We replaced the step-index fiber bundles and the PMTs of the original system with the two SiPMs in direct contact with the sample. Two boards for TCSPC (SPC130, Becker&Hickl GmbH, Germany), hosted in a peripheral component interconnect box (PCI-Cardbox 2F, IBP Instruments Gmbh, Germany) connected to a laptop are used for the parallel acquisition of the two photon DTOF in the two channels. The overall system is controlled by a home-made software, written in the C language, interfaced to a microcontroller unit (dsPIC30F2010, Microchip Technology Inc., Arizona), which provides a trigger signal to synchronize the optical switch and the data acquisition of the TCSPC boards.

The laser power delivered to the sample under investigation was 1.5 and 0.4 mW, respectively, at 690 and 829 nm, which gives IRF with an FWHM of less than 230 ps in both cases, as shown in Fig. 2.

Fig. 2

Example of IRF at the two wavelengths.

3.1.

Light Harvesting Capability

To test the light harvesting efficiency of our single-photon detector, we measured the responsivity parameter for diffuse optics instruments as described in the basic instrument performance (BIP) assessment protocol.¹⁵ This figure is meant to quantify the overall photon harvesting efficiency when dealing with a distributed quasi-Lambertian source and accounts for both the electronic photon detection efficiency and the optical chain collection efficiency. Being designed to be used in contact with the sample under investigation, we placed the detector directly on the reference phantom surface (described in Ref. 15), thus exploiting the whole detector NA. In Table 1, the results for one SiPM detector at the two wavelengths are presented.

Table 1

Responsivity values for different detection systems.

The values obtained are consistent with the ones already measured by Dalla Mora et al.¹¹ We then compare these results with a standard PMT tube (PMT, PMC-100, Becker&Hickl GmbH, Germany), which is often employed in TD fNIRS instruments.¹⁴ In this case, light from the sample was coupled with the detector’s active area (8-mm cathode diameter) with a glass optical fiber bundle (Loptek Glasfasertechnick GmbH, Germany) with an inner diameter of 3 mm and an NA of 0.57. The PMT’s responsivity is of the same order of magnitude compared to the SiPM one, but the PMT’s active area is much bigger than the SiPM one. We have also compared these findings with another detector, which was used in a TD fNIRS system:¹⁶ a hybrid PMT (HPM-100-50, Becker&Hickl GmbH, Germany) with 3-mm cathode’s diameter and which employ the same fiber bundle described for the PMT. In this case, the calculated responsivity was one order of magnitude smaller. The increase in the responsivity values provided by SiPM with respect to hybrid-PMT is mainly due to a theoretical NA of 1, together with the absence of optical fibers and lenses that usually attenuate the signal. We can also compare our results with the ones shown in Ref. 15, where the responsivity values, for the brain imagers present in some research’s laboratories until 2014, are reported. Here, we can find responsivity of the same order of magnitude or smaller than the ones found with the SiPM detector. Since the active areas of the detectors listed are larger, we can affirm that we have a clear improvement in this parameter.

3.2.

Linearity and Accuracy

The MEDPHOT protocol¹⁷ allows one to assess the system performances on homogeneous diffusive media and to compare different photon migration instruments. In the present study, we specifically apply two MEDPHOT tests to ascertain the “linearity” and “accuracy” of the system to retrieve the absorption (${\mu }_{\mathrm{a}}$) and reduced scattering (${\mu }_{\mathrm{s}}^{\prime }$) coefficient in homogeneous media. The tests are performed on homogeneous solid phantoms with different optical properties, labeled with a letter (A to C) and a number (2 to 8) representing, respectively, the nominal ${\mu }_{\mathrm{s}}^{\prime }$ (5 to $15{\mathrm{cm}}^{-1}$, $5\text{-}{\mathrm{cm}}^{-1}$ steps) and the nominal ${\mu }_{\mathrm{a}}$ (0.07 to $0.49{\mathrm{cm}}^{-1}$, $0.07\text{-}{\mathrm{cm}}^{-1}$ steps) at 660 nm. The ${\mu }_{\mathrm{a}}$ and ${\mu }_{\mathrm{s}}^{\prime }$ values are obtained mixing the proper quantity of ${\mathrm{TiO}}_2$ and toner powder inside an epoxy resin mould. The measurements were performed with both the detection channels, in a reflectance geometry (3-cm source–detector distance) with an acquisition time of 1 s for each wavelength, with five repetitions for each wavelength and channel. To calculate the ${\mu }_{\mathrm{a}}$ and ${\mu }_{\mathrm{s}}^{\prime }$, we performed a fitting procedure between the time-resolved reflectance curve and the analytical solution of the diffusion equation in an infinite homogeneous slab¹⁸ geometry previously convoluted with the IRF.

In Fig. 3, the results at 690 nm for one channel are shown. In Figs. 3(a) and 3(d), the measured ${\mu }_{\mathrm{a}}$ and ${\mu }_{\mathrm{s}}^{\prime }$ versus their conventional true values are plotted, respectively. The conventional values have been measured by means of a TD NIRS instrument,¹⁹ which guarantees a good estimate of the optical parameters thanks to a narrow ($<100\mathrm{ps}$), clean, and stable IRF.²⁰ The points represent the average values over the different scattering phantoms for each absorption value (a) and over the different absorption phantoms for each scattering value (d). The error bars are the corresponding standard deviations. The solid lines are the linear fit of the whole measurement points. We can notice how the experimental points are well aligned along these lines (goodness of the fit: $R^2=0.996$ for ${\mu }_{\mathrm{a}}$ and $R^2=0.989$ for ${\mu }_{\mathrm{s}}^{\prime }$) showing a good linearity.²¹ In Figs. 3(b) and 3(c), the measured ${\mu }_{\mathrm{a}}$ and ${\mu }_{\mathrm{s}}^{\prime }$ versus the conventional true ${\mu }_{\mathrm{s}}^{\prime }$ and ${\mu }_{\mathrm{a}}$, for the different phantoms are plotted, respectively. The solid lines are the linear fit of the whole measurement points. For ${\mu }_{\mathrm{a}}$ [Fig. 3(b)], we have a good fit until phantom 3 ($R^2<0.972$), while with the increasing of the absorption value, we have a decrease in the goodness of the fit ($0.143<R^2<0.494$), with the exception of phantom 7, where $R^2=0.978$. For ${\mu }_{\mathrm{s}}^{\prime }$ [Figs. 3(c)], we can see that the points are always well aligned with the fit lines ($R^2>0.979$) except for phantom 3, where $R^2=0.067$. Similar results were obtained for the 829-nm wavelength and for the other detection channel (data not shown).

Fig. 3

(a and d) Linearity and (b and c) accuracy for the absorption (${\mu }_{\mathrm{a}}$) and reduced scattering (${\mu }_{\mathrm{s}}^{\prime }$) coefficients for one channel at 690 nm. Points in panel (a) and (d) are the average of different phantoms, the error bars, and the standard deviations.

3.3.

Detection of Deep Absorption Changes

To validate the application of the instrument in clinical problems related to a localized change in absorption properties, such as functional brain imaging or muscle metabolism monitoring, we performed measurements applying the nEUROPT protocol.²² In particular, we assessed the “sensitivity,” expressed by the contrast as a function of depth, and the “spatial resolution,” based on a lateral scan. The nEUROPT tests were taken on an inhomogeneous solid phantom with switchable black polyvinyl chloride inclusions of different dimensions, which mimics a realistic absorption perturbation.²³ The inclusions (black cylinders with diameter and height of 3, 4, 5, and 7 mm, corresponding to an equivalent absorption variation of 0.05, 0.10, 0.17, and $0.40{\mathrm{cm}}^{-1}$, for a volume of $1{\mathrm{cm}}^2$, respectively) were placed inside a movable cylinder (solid epoxy resin with same optical properties of the host phantom). Both lateral (Y-scan, measuring surface parallel to the inclusion movement) and depth scans (Z-scan, measuring surface perpendicular to the inclusion movement) were performed in 1-mm steps. The time-dependent contrast as a function of the inclusion position $x_{\mathrm{inc}}$, at 690 nm, was calculated as

The contrast for the Y- and Z-scan for the 5-mm black inclusion is shown in Figs. 4(a) and 4(b), respectively. The contrast was calculated in 10 different time delays (0 to 3.50 ns, 0.35-ns step), but in the figure in the first one (0 to 0.35 ns), it is not shown since the noise was much bigger than the signal in this temporal window. During the Y-scan, we can see that at the zero point position, where the inclusion is exactly under the source–detector pair’s plane, the contrast has a maximum, while it decreases with the progressive estrangement of the inclusion. Upon increasing the time delay, an enhancement of the contrast is observed. No contrast is observed in the early gates (0.35 to 1.05 ns), which are related to the region close to the phantom surface. The last temporal delay (3.15 to 3.50 ns) is noisier, indicating that the signal, in this part of the curve, is lower. In the Z-scan (b), the starting position corresponds to the inclusion border coincident with the entrance surface, just outside the bulk phantom. The peak in the contrast for both the late and early gates is now reached at different positions, which are farther from the injection surface if the time delay is higher. These results show the excellent capability of this measurement setup to detect local changes of absorption, thus validating the SiPM detector as a very great candidate for TD fNIRS in vivo measurement. For the other inclusions with different diameters, we found similar results with an increased contrast for the larger inclusions.

Fig. 4

Contrasts as function of the position obtained during (a) the $Y$ and (b) scan, for the 5-mm inclusion, in function of the different temporal windows. The first temporal window (0 to 0.35 ns) is not shown.

3.4.

Long- and Short-Term Stabilities

The IRF stability test (foreseen in the BIP protocol¹⁵) is an important factor for determining how quickly the whole system reaches the thermal equilibrium, to quantify the system warm-up time and for understanding how long this stability is maintained. We recorded the system IRFs for the two channels and the two wavelengths, starting from the system power-on for 5 h, collecting 10 measurements of 1 s every 1 min. After the constant background subtraction (dark count rate), the total photon counts and the IRF time position (as identified by its barycenter) were calculated. For the photon counts at 690 nm, obtained integrating the IRF waveform, to enter into the $\pm 3\%$ ($\pm 1\%$) range of the average of the last 30 min (long-term stability), it is necessary to wait about 120 min (210 min). For the time position of the IRF, it is necessary to wait only 10 min (50 min) to enter into the $\pm 3\%$ ($\pm 1\%$) range of the average of the last 30 min.

Since during in vivo recordings, the measurements performed are often no longer than 30 min, we have therefore evaluated the short-term stability by dividing the whole 5-h long stability, in 10 different intervals of 30 min each, and by calculating the variation of photon counts and time position of the IRF within these short periods. Only in the first time interval (soon after turning on the system), the percentage of intensity variation is around 2.5% (calculated as the standard deviation in the period). For the following intervals and for the IRF barycenter, it is always lower than the 1%. This fact suggests that it is possible to use the system after 30 min of warm-up.

To disentangle the contribution of the detector and of the laser during the warm-up period, we performed another stability measurement, in which the lasers were already kept on for 3 h, which allows us to consider them stable within 1%. In this case, to enter into the $\pm 3\%$ ($\pm 1\%$) range of the average of the last 30 min (total measurement length: 4 h), it is necessary to wait about 70 min (170 min) for the data counts and 5 min (30 min) for the barycenter. Similar results were found for the other detector and wavelength. These results show that some thermal instability can be effectively ascribed to the SiPM front-end, probably due to an ineffective heat dissipation strategy.

4.1.

Arterial Arm Occlusion

We performed an arterial (300 mmHg) cuff occlusion of the right arm muscle of a male adult volunteer. We employed two different measurement channels (i.e., source–detectors couples), one for the occluded arm and one for the other, which was kept in a basal condition for the whole experiment. The detectors were placed directly in contact with the arms by means of a black rubber pad probe, with a source–detector separation of 3 cm. The protocol consisted of 60 s of baseline, where the basal hemodynamic condition was registered, 150 s of occlusion, and 90 s of final recovery. The sampling rate was set at 1 s. The thicknesses of the superficial layer for the left and right arm were measured with a skin folder caliper and were set at 2 and 3 mm, respectively.

The baseline average values for the hemodynamic parameters for left (right) arm were: ${\mathrm{O}}_2\mathrm{Hb}=119(108)\mu \mathrm{M}$, $\mathrm{HHb}=39(50)\mu \mathrm{M}$, $\mathrm{tHb}=157(158)\mu \mathrm{M}$, and ${\mathrm{SO}}_2=75$ (68)%. In Fig. 5, the variations of the time course of ${\mathrm{O}}_2\mathrm{Hb}$, HHb and tHb concentrations, and of ${\mathrm{SO}}_2$ percentage with respect to the baseline during the arterial occlusion are shown. A median filter of the third order was applied, since it represents a good means to remove the noisy spikes without smoothing too much the signal causing a loss of information. In Fig. 5, data for the left (light line), not occluded, arm and right (thick line), occluded one, are presented, respectively. During the arterial occlusion, since both the arteries and the veins are closed, we can observe, for the occluded arm, in Fig. 5(b), an increase of HHb, and in Fig. 5(a), a decrease of ${\mathrm{O}}_2\mathrm{Hb}$. Since no new oxygenated blood can reach the arm and no deoxygenated can go outside it, but the muscular metabolism continues, ${\mathrm{SO}}_2$ decreases and tHb increases during the occlusion period, as shown in Figs. 5(d) and 5(c), respectively. When the cuff is released (around 210 s), the typical hyperemia peak is present. No significant changes in the left arm are observed.

Fig. 5

Arterial occlusion of the right arm. Variations of the hemodynamic parameters with respect to the baseline for left (light line), not occluded, and the right (thick line), occluded arm. [(a) oxy-hemoglobin, (b) deoxy-hemoglobin, (c) total-hemoglobin, and (d) saturation].

4.2.

Brain Motor Cortex Activation

We performed two motor task experiments consisting of a finger-tapping exercise with the right and left hand, respectively, involving one adult male volunteer. The two measurement channels were placed around the C3 and C4 positions according to the $10/20$ EEG international system, in order to monitor the hemodynamic changes in both the ipsi- and contralateral brain motor cortex. The protocol was based on four repeated blocks consisting of: 20 s of rest, 20 s of task, and 20 s of recovery, for 240 s of total experiment length. The acquisition rate was set at 1 s and the source–detector separation was fixed at 3 cm. The detector and the injection fibers were put in contact with the scalp by means of a standard EEG cap (EEGcap, g-TEC, Austria), as the one shown in Fig. 2(c), according to the $10/20$ EEG international standard.

The baseline average values for the hemodynamic parameters for left (right) hemisphere were: ${\mathrm{O}}_2\mathrm{Hb}=53(56)\mu \mathrm{M}$, $\mathrm{HHb}=25(27)\mu \mathrm{M}$, $\mathrm{tHb}=78(83)\mu \mathrm{M}$, and ${\mathrm{SO}}_2=68$ (68) %. In Figs. 6(a)–6(d), the time course of $\mathrm{\Delta }{\mathrm{O}}_2\mathrm{Hb}$, $\mathrm{\Delta }\mathrm{HHb}$, $\mathrm{\Delta }\mathrm{tHb}$, and $\mathrm{\Delta }{\mathrm{SO}}_2$, calculated with respect to the baseline values, respectively, for the left and right hemispheres during a right-hand finger-tapping task, are shown. The average over the four repetitions is shown and the error bars represent the standard deviation over the repetitions. A median filter of the third order was applied. For this subject, we considered a thickness of the UP layer of 1.25 cm that was manually estimated by looking at the anatomical MRI images. In the contralateral hemisphere, we can notice an evident cortical activation, i.e., an increase of ${\mathrm{O}}_2\mathrm{Hb}$ and a contextual decrease of HHb [Fig. 6(a)] that yield to an increase of tHb and of ${\mathrm{SO}}_2$ [Fig. 6(b)]. In the ipsilateral hemisphere [Figs. 6(c) and 6(d)], no activation is shown ($\mathrm{\Delta }{\mathrm{O}}_2\mathrm{Hb}$ and $\mathrm{\Delta }\mathrm{HHb}$ have non task-related fluctuations around zero), as expected during a motor cortex activity.

Fig. 6

Right-hand finger tapping. Variations of the hemodynamic parameters with respect to the baseline for (a and b) the left and (c and d) right hemispheres. The error bars are the standard deviations over the four task repetitions.

In Fig. 7, the left-hand finger-tapping task is shown. The graphs were obtained in the same way as for the previous ones. Again, we can clearly observe the activation in the contralateral hemisphere, with a slight increase of ${\mathrm{O}}_2\mathrm{Hb}$ and a corresponding decrease of HHb [Fig. 7(c)], with an increase of ${\mathrm{SO}}_2$ and tHb [Fig. 7(d)]. In the ipsilateral hemisphere [Figs. 7(a) and 7(b)], we notice instead an increase in ${\mathrm{O}}_2\mathrm{Hb}$ and an almost constant HHb. These results can be ascribed to a task related physiological response not involving a cortical activation.

Fig. 7

Left-hand finger tapping. Variations of the hemodynamic parameters with respect to the baseline for (a and b) the left and (c and d) right hemispheres. The error bars are the standard deviations over the four task repetitions.