Performance evaluation of a staggered three-layer DOI PET detector using a 1 mm LYSO pitch with PETsys TOFPET2 ASIC: comparison of HAMAMATSU and KETEK SiPMs

In this study, we propose a staggered three-layer depth-of-interaction (DOI) detector with a 1 mm crystal pitch and 19.8 mm total crystal thickness for a high-resolution and high-sensitivity small animal in-beam PET scanner. A three-layered stacked LYSO scintillation array (0.9 × 0.9 × 6.6 mm3 crystals, 23 × 22 mm2 surface area) read out by a SiPM array (8 × 8 channels, 3 × 3 mm2 active area/channel and 50 μm microcell size) with data acquisition, signal processing and digitization performed using the PETsys Electronics Evaluations kit (based on the TOFPET v2c ASIC) builds a DOI LYSO detector block. The performance of the DOI detector was evaluated in terms of crystal resolvability, energy resolution, and coincidence resolving time (CRT). A comparative performance evaluation of the staggered three-layer LYSO block was conducted with two different SiPM arrays from KETEK and HAMAMATSU. 100% (KETEK) and 99.8% (HAMAMATSU) of the crystals were identified, by using a flood irradiation the front- and back-side. The average energy resolutions for the 1st, 2nd, and 3rd layers were 16.5 (±2.3)%, 20.9(±4.0)%, and 32.7 (±21.0)% (KETEK) and 19.3 (±3.5)%, 21.2 (±4.1)%, and 26.6 (±10.3)% (HAMAMATSU) for the used SiPM arrays. The measured CRTs (FWHM) for the 1st, 2nd, and 3rd layers were 532 (±111) ps, 463 (±108) ps, and 447 (±111) ps (KETEK) and 402 (±46) ps, 392 (±54) ps, and 408 (±196) ps (HAMAMATSU). In conclusion, the performance of the staggered three-layer DOI detector with 1 mm LYSO pitch and 19.8 mm total crystal thickness was fully characterized. The feasibility of a highly performing readout of a high resolution DOI PET detector via SiPM arrays from KETEK and HAMAMATSU employing the PETsys TOFPET v2c ASIC could be demonstrated.


Introduction
In recent years the effort to improve the spatial resolution and quality of PET imaging systems has steadily increased. The most commonly pursued approach is using time-of-flight (TOF) PET, which was already suggested by Brownell et al (1969). The first commercially used TOF PET scintillators (CsF and BaF 2, 1-to-1 coupled to PMTs) could achieve TOF capabilities between 470 and 750 ps. However, they were still suffering from low light yield and low material densities, resulting in limited detection efficiency and spatial resolution, the latter due to the large crystal sizes required (Melcher 2000, Surti and Karp 2016, Conti and Bendriem 2019, Saint-Gobain Crystals 2020a. Nowadays modern scintillation crystals such as LSO and LYSO, in combination with silicon photomultipliers (SiPMs) and dedicated signal processing electronics, allow for coincidence resolving times (CRT) of commercial systems to reach close to 200 ps (Surti 2015, Reddin et al 2018, Conti and Bendriem 2019, van Sluis et al 2019, while providing high detection efficiencies due to their high Z and density (Saint-Gobain Crystals 2020b) and high spatial resolution (e.g. by using an Anger-type detector readout). For non-commercial laboratory setups CRT values around 100 ps could be demonstrated by using novel scintillation crystals such as LaBr 3 :Ce and CeBr 3 (Glodo et al 2006, Schaart et al 2010. The continuous improvement of the CRT results in an improved image contrast due to an improved signal-to-noise ratio (SNR) (Surti et al 2007, Surti 2015. However, the ultimate goal targets a CRT of about 10 ps (Grundacker et al 2019), which would allow to directly obtain the β + -annihilation position by measuring the detection time difference with an accuracy that is similar to the current reconstruction methods that rely on the intersection of the linesof-response (LOR) of multiple detected positron annihilation photon pairs (Lecoq et al 2020).
While large radius PET systems, such as human body scanners, strongly benefit from an improved TOF capability and its improvements on the SNR, for smaller scanner radii (such as small-animal-PET scanners) a high detector spatial resolution becomes more important. The spatial resolution of the reconstructed PET images relies on the cross section of the scintillation crystal. For the commercially available preclinical PET scanners, the scintillation cross sections range from to 1.5×1.5 mm 2 to 1.0×1.0 mm 2 depending on the manufacturers (Kunter et al 2014). The spatial resolution of PET scanners can be further enhanced by using a fine pitch of scintillation crystals, however, the crystal pitch has a tradeoff with the geometric efficiency given by the solid angle coverage that reduces with smaller crystal sizes, due to the optical insulation layers between individual crystals. The spatial resolution of small animal PET scanners is significantly deteriorated at the periphery of the PET field-of-view (FOV) caused by parallax errors introduced by the crystal thickness due to which the true LOR deviates from the reconstructed one if one or both γ rays interact between the crystal surface and the photosensor surface (Hoffman et al 1989, Pomper and Lee 2005, Ito et al 2010. However, a certain crystal thickness is required to provide a sufficient detection efficiency (system sensitivity) for 511 keV γ rays. The optimum detector thickness depends on the crystal material and is typically between 20 and 30 mm. To minimize the parallax error and to preserve the spatial resolution in the periphery of the FOV various types of depth-of-interaction (DOI) detectors have been proposed (Mohammadi et al 2019).
Consequently, by using a detector which consists of stacked crystal layers with a defined offset of the crystals in the respective layers and therefore providing DOI information (limited to the thickness of the individual crystal layers) the image resolution can be improved, while the efficiency of the detector block remains the same as for conventional single-layer detectors with identical crystal pitch and overall thickness of the full DOI detector block (Ito et al 2010). Such multi-layer scintillation crystal configurations to be used in PET were previously evaluated in various studies. E.g. Prout el al (2020) presented a phoswich detector using LYSO and BGO crystals which were read out with HAMAMATSU MPPC arrays and the TOFPET ASIC by PETsys for signal digitization. Another concept for DOI encoding that also uses the PETsys TOFPET ASIC for signal integration and digitization was given by Yoshida et al (2020). Kang et al (2019) and Takyu et al (2018), presented studies with multi-layer LYSO crystal blocks, with crystals of the respective layers shifted with respect to each other, which is the same type of DOI crystal block also evaluated in this study.
However, the aim for smaller crystal front areas and potentially stacked crystal layers also sets high demands on the performance of the scintillator readout, i.e. the photosensors and the signal processing electronics, which need to allow to resolve all individual scintillation crystals of the detector while providing sufficiently good energy and timing resolutions as well as the ability to process high count rates. Furthermore, the system should be scalable without having to sacrifice in performance (e.g. caused by potentially increased detector noise etc).
In this work, the detailed comparative performance characterization of a three-layer LYSO PET detector with a DOI resolution of 6.6 mm and a crystal pitch of 1 mm and a highly integrated readout based on different SiPMs and an ASIC-based signal processing unit (PETsys Electronics TOFPET v2c ASIC) is presented. The staggered three-layer LYSO crystal block is conceptually based on the principle described in Ito et al (2010) and was designed using GATE optical simulations at the National Institute for Radiological Science (NIRS) in Chiba, Japan (Kang et al 2018). The focus of this work is set to (i) a systematic investigation of the possibility to resolve the individual crystals of the three detector layers with two different photosensor arrays, (ii) the achievable energy resolution as well as (iii) the CRT for an electronic readout system based on SiPM arrays and the PETsys Electronics TOFPET v2c ASIC and (iv) a comparison between SiPM arrays of two manufacturers (KETEK and HAMAMATSU). The used ASIC was initially developed for TOF-PET applications (Di Francesco et al 2016, PETsys Electronics S.A. 2019. In contrast to the previous studies performed at the NIRS (Kang et al 2019, Prout et al 2020), this readout system does not rely on a front-end resistor-chain board, which generates four analog signals (from the initially 64 SiPM channels of the (8×8) SiPM array) used to determine the photon interaction position by an Anger-type calculation (Kang et al 2019). The PETsys TOFPET ASIC v2c triggers on and digitizes all 64 detector channel signals individually. All further signal evaluation is performed after the digitization (PETsys Electronics S.A. 2019), which improves the system performance already for 64 channel SiPM arrays, due to the absence of electronic tradeoffs such as the summation of the SiPM capacitance or noise. Other than used by the the authors of Prout et al (2020), the individual crystals of the DOI LYSO block evaluated in the present study have a crystal pitch of only 1.0 mm and the microcell size of the studied SiPM array is 50 μm in contrast to the 75 μm used for the MPPCs in Prout et al (2020).

The setup
The detector used for this investigation is a trapezoidal three-layer stacked LYSO crystal block consisting of 0.9×0.9×6.6 mm 3 individual crystals, optically isolated by a layer of 0.1 mm thick BaSO 4 in between each crystal. The first layer consists of 23×22 crystals, the second layer comprises 23×23 crystals and the third layer is composed of 24×24 crystals. The crystals in the second layer are shifted by half a crystal pitch with respect to the first layer along the x-axis and the crystals in the third layer have an offset of half a crystal pitch along the y-axis with respect to the first and second layers. In order to couple the DOI scintillator block to the photosensor, optical grease (St. Gobain BC-631) was applied. A KETEK PA3350WB-0808 prototype SiPM array with a breakdown voltage of 25.2 V (5.0 V recommended overvoltage/30.2 V applied voltage) and a HAMAMATSU S14161-3050HS-08 MPPC array with a breakdown voltage of 38.6 V (2.7 V recommended overvoltage/41.3 V applied voltage) (HAMAMATSU Photonics K.K 2020) were used. Both arrays have 8×8 channels with an active area of 3×3 mm 2 and 50 μm single-photon avalanche photodiodes. The SiPM channel pitch is 3.36 mm (KETEK) and 3.2 mm (HAMAMATSU), respectively (HAMAMATSU Photonics K.K 2020, KETEK GmbH 2020b). The PETsys Electronics Evaluation kit based on the TOFPET2 v2c ASIC was used as signal processing and data acquisition system. Two of these 64 channel ASICs are mounted on a front-endmodule (FEM128) to which the detector is connected. Each ASIC channel is operated independently, i.e. has its individual amplifier, discriminator and digitizer (TDC/ADC) (PETsys Electronics S.A. 2019). A simplified schematic of one ASIC input channel is shown in figure 1 (figure provided by PETsys Electronics). The input current is replicated into three independent branches (T, E and Q). These input currents are converted into a voltage signal by a transimpedance amplifier (TIA) of selectable gain G T and G E for the T and E branches signals, respectively. While G T can be set to 3000 Ω, 1500 Ω, 750 Ω or 375 Ω, the gain of the E branch is reduced by a factor of 10 with the respective values of 300 Ω, 150 Ω, 75 Ω and 38 Ω. The voltage signals in the T and E branch are fed into the discriminators by which the trigger logic is controlled. The T branch contains two comparators, which set the timing thresholds V th_T1 and V th_T2 , thus providing the logic signals do_T1 (do_T1' if delayed) and Figure 1. Simplified schematic of one ASIC channel. The input current signal is replicated into three input currents (T, E, Q branches). The input current signal is then amplified and converted into a voltage signal via a transimpedance amplifier. Comparators in the T and E branch set thresholds to the signal, which must be exceeded for an event to be triggered and digitized (PETsys Electronics S.A. 2019) (figure provided by PETsys Electronics). Reproduced with permission from PETsys Electronics.
do_T2, respectively. The E branch comprises only one comparator, using V th_E as threshold and providing do_E as logic output. For the nominal trigger logic mode, a rising edge of do_T1 will trigger Trigger_T, a rising edge on do_T2 triggers trigger_Q and a rising edge on do_E triggers trigger_E. Trigger_Q is generated if in a logical OR condition do_T1, do_T2 or do_E is active. The event, however, will only be digitized if all three triggers (trigger_T, trigger_Q and trigger_E) have a rising edge.
The input to the Q branch can be digitized if the system is used in the QDC mode to perform charge integration, which was the operation mode for all measurements performed in this work.
For the event triggering the nominal trigger mode was used. In this mode all three logic signals provided by the three comparators (do_T1, do_T2 and do_E) must be present for an event to be validated and digitized. The timestamp corresponding to a digitized signal in one ASIC channel is set at the time at which V th_t1 is exceeded. If V th_t2 is not exceeded and therefore do_T2 is not triggered, the event will be discarded with no dead time caused in the TDC/QDC. If V th_t2 is exceeded but V th_E is not, the event will be rejected with 5 clock cycles (25 ns) dead time.
The setup is pictured in figure 2 (top panel). On the left-hand side the DOI detector (wrapped with black tape) is connected to a FEM-128 front-end-module with the TOFPET v2c ASIC. On the right, a reference detector array (3×3×5 mm 3 LYSO crystals coupled one-by-one to a KETEK PA3325WB-0808 SiPM array) is Figure 2. Top: the setup used for coincidence measurements is shown. The DOI detector is visible on the left connected to one FEM-128 frontend module containing the TOFPET v2c ASIC, while the reference detector used to determine the CRT is placed on the right side of the radioactive calibration source which is centered inbetween both detectors. The FEM-128 board are connected via the blue ribbon cables to the motherboard which contains the FPGA and is responsible for the communication to the acquisition computer via a gigabit ethernet link. The inset shows a picture of the DOI LYSO crystal placed on the KETEK PA3350WB-0808 array. Bottom: schematic of a coincidence setup for the DOI detector and a small reference detector (left) and a reference detector of larger active area in the same position (right). The orange squares represent the projection of the solid angle that can be covered by the reference detector for coincidences with the DOI detector.
shown, which was used to measure the CRT. The radioactive source (here: 22 Na with an activity of 180 kBq resulting in an accepted count rate in the low kHz range depending on the exact distance between the detector front surface and the radioactive source) is placed in the center between the two detectors. The blue ribbon cables transfer the digital signals to the FPGA. For the flood map acquisition and energy resolution measurements, the data were taken in single-mode in which the reference detector is not connected to the FPGA board and only the DOI detector is used. For measurements of the CRT, the data were taken in coincidence mode in which the reference detector is configured with the identical settings (e.g. thresholds, integration windows, etc) as the DOI detector.
2.2. The reference detector array for CRT measurements When using a small reference detector compared to the total size of the DOI LYSO detector, geometrical constraints would allow for calculating the CRT only for a part of the DOI detector, if the setup is as compact as the used one with a distance between the DOI detector and the reference detector of only 60 mm. The reason is that the reference detector can only find opposing coincident channels within a solid angle spanned by all possible LORs and therefore depends on the detector's surface area and the distance to the opposing coincident DOI detector crystal. Figure 2 (bottom panel) shows a schematic view of the coincidence setup consisting of a DOI detector and a small reference detector (left) and a larger reference detector (right). The orange lines represent the envelope of all possible LORs and the orange squares represent the projection of the solid angle spanned by all possible LORs onto the DOI detector. Figure 2 (bottom panel) shows that for a small reference detector compared to the DOI detector only for a subset of all DOI channels coincidences can be found due to geometrical constraints.
To overcome this restriction either the distance between the two detectors can be increased, which reduces the coincidence rate and increases the measurement time. Alternatively, instead of the small single-channel reference detector, a full reference detector crystal array can be used to enlarge the solid angle coverage of the reference detector. The latter method was chosen in this work. The reference detector array was an 8×8 channel array (3×3×5 mm 3 LYSO crystals) one-to-one coupled to a KETEK PA3325WB-0808 SiPM array, resulting in a reference detector active area of 27 mm×27 mm. By using this reference detector array for the CRT measurements an opposite reference channel could be found for most of the 64 possible channels of the DOI detector.
The bias voltage for the reference detector was lower than actually recommended by KETEK (set to 30.1 V). The reason is the large difference in signal amplitudes between a one-by-one coupling (reference detector array) and a light sharing approach (DOI detector). As consequence, the thresholds for optimum results are different. The focus was set to maximize the performance of the DOI detector. By reducing the bias voltage and consequently the gain of the reference detector, the performance with the low thresholds used for the DOI detector could be improved. Even though, in general, a higher bias voltage and the resulting steeper rising edge of the signal results in an improved CRT (until the deteriorating influence of noise such as dark counts become dominant), no deterioration on the measured CRT of the DOI detector is expected due to the low bias voltage of the reference detector array, since its time resolution, and therefore its contribution, is deconvolved from the DOI detector's CRT.

The configuration
All measurements were performed with the signal processing system set to the QDC mode, i.e. the charge integration stage (QDC) was used. As mentioned in the previous section, the full trigger logic (nominal mode) was used, including the timing threshold V th_T1 for dark count rejection and providing an accurate time stamp (PETsys Electronics S.A. 2019). The thresholds V th_T1 and V th_E were used for event validation and to start the actual signal integration (PETsys Electronics S.A. 2019). The integration window length is set dynamically between 20 and 900 ns and closes-on an individual channel basis-after the signal has dropped below the three respective threshold levels plus a delay of 6 ns (falling edge of trigger_B). The gain of the TIAs was set to the maximum value of 3000 Ω (G T ) for the time branch and 300 Ω (G E ) for the charge branch. The integrator gain was also set to the maximum value of 3.65. While the threshold starting the integration (v th_T2 ) was set differently for the KETEK SiPM and HAMAMASTSU MPPC arrays, V th_T1 -defining the timestamps-was set to the minimum value of 2.8 mV (PETsys Electronics S.A. 2019) for the SiPM arrays of both manufacturers. For the flood maps the applied threshold settings were V th_T2 =12.6 mV, V th_E =9.8 mV and V th_T2 =35.0 mV, V th_E =32.2 mV for the KETEK and the HAMAMATSU array, respectively.
The system was actively cooled by a 12 V DC fan, which removes the warm air around the ASICs. While the temperature measured at the ASIC acquiring and processing the signals of the DOI detector was measured to be 28.4°C ± 0.3°C during the acquisition of the flood maps, the temperature was decreased to approx. 25.0°C ± 0.3°C by using a more powerful fan when measuring the CRT of the detector. In accordance with information from the manufacturer, the ASIC's gain increases with decreasing temperature. Also the SiPM breakdown voltage decreases with rising temperature (22 mV K −1 (KETEK) KETEK GmbH 2020c) and 34 mV K −1 (HAMAMATSU) (HAMAMATSU Photonics K.K 2020), which naturally leads to a higher applied overvoltage and therefore to an increase of the SiPM gain. The gain change of the system can only be given as entangled effect arising from a gain increase of both the SiPM and the ASIC and is approximated from previous experiments (in the temperature region between approx. 23°C and 40°C) to approx. 2.5% per°C. As a consequence, the thresholds V th_T2 and V th_E could be increased for the KETEK SiPM to determine the CRT and therefore for the SiPM arrays of both manufacturers identical thresholds could be used (V th_T1 =2.8 mV, V th_T2 =39.2 mV ,V th_E =28.8 mV). This assured identical threshold settings for the reference detector during the CRT measurements of both SiPM/MPPCs.

The error analysis
The errors given to the measured values presented in the next chapter are determined as follows: the statistical error given for all values of the energy resolution arises from the error on the Gaussian fit's sigma value that is fitted to the photo peak. Since the fit error of the Gaussian's centroid only contributes with about 1% to the statistical uncertainty it is neglected for the statistical error given in the following analysis.
For the energy resolution of the individual layers, 12 crystals in a central region and 12 crystals in the edge region were evaluated for each layer. The obtained 12 values for each set of energy resolutions form a distribution. The 1σ standard deviation of this distribution is considered as the 'inter-crystal variability' uncertainty (δ var ) of the energy resolution (represented by the length of the antennas in figure 6).
The dependency, however, of the energy resolution's mean value on the layer causes the distribution for the total energy resolution to be asymmetric and the concept of using the mean value and 1σ standard deviation for the 'inter-crystal variability' is not practical. Therefore, the resulting total energy resolution will be given as the median value and the 'inter-crystal variability' uncertainty will be given as the median's difference to the 10% and 90% levels of the range, representing the energy resolution intervals that contain 80% of the values.
For the CRT three main sources contribute to the uncertainties. The statistical uncertainty originates from the fit uncertainty of the gaussian fit to the central peak in the time-difference spectra that were used to derive the CRT.
Like for the error analysis of the energy resolution, there is also a contribution by the inter-crystal variability to the uncertainty of the CRT. However, for the CRT there on the individual layers, that is derived from the CRT of individual DOI crystal's CRT, also a second contribution was evaluated and is included in the d var component of the error. This second contribution is the variability of the arrival time measurements between the channel combinations to the overall CRT, when measured versus all combinations of reference detector channels found in coincidence. This contribution to the error includes a geometrical factor. It could be shown that in some events not only one but two coincident reference detector array channels directly opposing the DOI detector channel of interest (relative to the 22 Na source) exhibit a minimum of the CRT. In contrast, the CRT increases when being measured against more distant reference channels until finally no coincidences are found any more.
For the given values of the total CRT (by the first and brightest firing channel, respectively), the second component is not present, since only one opposing reference detector channel was used to derive the CRT of one specific SiPM/ASIC channel.

Arrival time distribution
As mentioned in the previous section the PETsys TOFPET v2b ASIC is designed for one-by-one crystal-to-SiPM coupling and therefore triggers each channel individually. In order to obtain data from a detector that involves light sharing, the individual SiPM/MPPC signals need to be grouped together in order to assign the multiplicity of SiPM/MPPC channels to one initial gamma event detected in one scintillation crystal. PETsys Electronics provides a post-processing script that clusters triggered ASIC channels such that all channels firing within a (selectable) coincidence window (grouping window) around one triggered channel are assigned to a common event. For a channel to be considered inside the grouping window the arrival time stamp, set when V th_T1 is exceeded, is used. Furthermore, a geometrical constraint can be set such that only SiPM/MPPC channels within a given radius around the initially triggered channel are considered. Also, multiple hits of one channel are excluded. The grouping time window needs to be set such that all firing SiPM/MPPC channels of an initial gamma are grouped together, but ideally no noise and background events are added to the actual gamma event distribution. Figure 3 shows the arrival time distribution with respect to the first firing channel and a zoom into the region from t=0 to t=200 ns (bottom row) for the KETEK PA3350WB-0808 array (left column) and the HAMAMATSU S14161-3050HS-08 (right column). The energy axis is normalized to the total detected charge. Figure 3 shows that after approx. 70 ns no further SiPM pixels are caused to fire by the scintillation light of an initial photon.
Using these results, the grouping window was set to 75 ns in order to ensure that all firing channels are grouped together.

Spatial information
The capability to spatially resolve each individual scintillation crystal of the DOI detector was shown by a flood irradiation of the detector with a 22 Na calibration source of 180 kBq activity. The radioactive source was placed in a distance of 2 cm (central position) in front of the detector and data were taken for 600 s. The raw data were consecutively clustered into events belonging to one initial γhit in the detector by running a grouping routine provided by PETsys Electronics. The time window of this routine that defines the search range for quasicoincident ASIC channels triggering after the first registered hit and thus determines what is considered to belong to the same γ event was set to 75 ns. The 2D interaction position within a detector layer was calculated via an Anger-type logic and an energy window from approx.−5 σ to +3 σ around the photo peak of the inclusive energy spectrum was applied. For the Anger calculation the 64 channels of the SiPM array were clustered into 8 rows and columns (ROW i and COL i ), respectively, with ROW i (COL i ) being the total energy detected by the respective 8 SiPM of the ith row (column). The interaction position in x (y) was calculated according to The resulting (x, y) flood maps reveal that the individual crystals can clearly be resolved (figures 4 and 5 center rows). Most of the interactions are detected in the first layer (seen as brightest crystal response marked by a white triangle in the zoomed right-hand top panel of figure 4) due to the highest probability of 511 keV γ rays to be detected within the first 6.6 mm. Vertically shifted by half a crystal size, hits in the second layer appear with medium brightness in the flood map (marked with a white circle). The faintest spots, shifted horizontally by half a crystal size with respect to the first layer (highlighted by a white rectangle), belong to interactions in the third layer, being the least probable for interactions to occur. For both SiPM arrays under study from different manufacturers the flood maps exhibit clearly resolved crystal images with only minor distortions (see figures 4 (KETEK) and 5 (HAMAMATSU)).
In order to confirm the assignment of photon interactions to their respective layers, also data were taken from a back-side irradiation of the detector block (i.e. placing the source in front of the SiPM array/frontend readout board). The bottom row of figure 4 shows a flood map obtained from such a back-side irradiation of the DOI detector read out by a KETEK PA3350WB-0808. The ordering of brightness of hits in the respective layers is reversed compared to the front-side irradiation. Due to geometrical constraints of this arrangement, imposed by the inhomogeneous matter distribution from the FEM-128 frontend boards carrying the ASIC chip and being plugged to the SiPM/MPPC arrays from the back-side, the 22 Na radiation source could not be placed in a central position for this irradiation scenario. This asymmetric irradiation explains why in the flood maps taken by . Flood map acquired with the three-layer DOI LYSO crystal block read out by a KETEK PA3350WB-0808 SiPM array operated at 5 V overvoltage (middle and bottom panels of left column), a zoom into the regions marked by the red rectangles (middle and bottom panels of right column) and the count rate profile of crystal column 8 and row 13 (top row). The middle panels show a flood map acquired in a 600 s front-side irradiation with a 22 Na source, while the data shown in the bottom panels were acquired during a 600 s back-side irradiation (see text). Crystals belonging to the first layer are indicated by a white triangle, crystals of the second layer by a circle and the squares indicate crystals of the third layer. The red triangle exemplarily highlights a region which only contains background and inter-crystal scattering events as will be described in section 'inter crystal scattering (ICS)'.
back-side irradiations (figures 4 and 5 bottom rows) the crystal response appears brighter in the lower part of the flood map than in the top part.
By evaluating the flood maps obtained by a front-side and a back-side irradiation 1611 out of 1611 (KETEK) and 1608 out of 1611 (HAMAMATSU) crystals could be identified.

Energy resolution
The relative energy resolution ( DE E ) was determined using a 180 kBq 22 Na calibration source placed in a distance of about 5 cm centrally in front of the detector. The signal charge of any detected γ interaction within a selected crystal was filled in a histogram (energy spectrum). A calibration that takes into account nonlinearities, such as SiPM saturation, was applied by using 152 Eu (121, 244 and 344 keV), 22 Na (511 and) and 137 Cs sources (662 keV) to determine the photo-peak positions of their respective γ transitions and by fitting a quadratic calibration function Figure 5. Flood map acquired with the three-layer DOI LYSO crystal block read out by a HAMAMATSU S14161-3050HS-08 MPPC array operated at 2.7 V overvoltage, irradiated with a 180 kBq 22 Na source (middle and bottom panels of left column), a zoom to the region marked by the red rectangle (middle and bottom panels of right column) and the count rate profile of crystal column 8 and row 13 (top row). The middle row shows a flood map acquired in a 600 s front-side irradiation with a 22 Na source, while the data shown in the bottom row were acquired during a 600 s back-side irradiation (see text). Crystals belonging to the first layer are indicated by a white triangle, crystals of the second layer by a circle and the square indicates crystals of the third layer. The white square in the middle left panel marks the two regions which were considered to contain 'edge' and 'central' crystals as discussed in the section 'energy resolution'.
with E keV and E a.u. being the measured energy in keV and arbitrary units, respectively, the parameters A, B and C being determined for each crystal separately by fitting and which was subsequently applied to the raw energy spectra. For a PET application the energy resolution of the two individual crystals that detect the two coincident 511 keV annihilation photons matters. Therefore, the energy resolution was exemplarily calculated for 24 individual crystals (12 at edge positions, 12 at a central position, see figure 5 center left panel) for all three layers (see table 1). However, for the egde region of layer 3 measured with the KETEK SIPM array, only seven of the selected 12 crystals showed a spectrum that could be used to derive a spectrum.The Crystal regions were defined manually and regions with less than approx. 8%-12% of the maximal amplitude were considered as background by best effort.
We compared the measured energy resolution layer-wise and for the central and edge regions separately. For a confidence level 0.05, the energy resolution of the first layer (central and edge region) was superior when measured with the KETEK SiPM array compared to the measurements performed with the HAMAMATSU MPPC array (p-value=0.045 (central) and p-value=0.045 (edge)). For the energy resolution of the second layer (central and edge region) with p-value=0.994 (central) and p-value=0.620 (edge) and the central region of the third layer (p-value=0.072) no difference could be measured with statistical significance (p-value: 0.994). However, for the edge region of the third layer, the energy resolution measured with the HAMAMATSU MPPC array was superior compared to that measured with the KETEK SIPM array (p-value=0.019).
The error of the given values increases with deeper layers due to low statistic because of the very low probability for the low energy gamma rays (121, 244 and 344 keV from 152 Eu) to be detected in the deeper layers. Furthermore, an influence of the light loss at the edges can be observed. For the third layer, where the deterioration between the edges and the center is expected to be the strongest due to the coupling between the SiPM array and the crystal, a difference of around 6.3% ( The summarized energy resolution of the individual layers as well as the energy resolution of all layers calculated from the individual layers' energy resolution is listed in table 1, while being partially illustrated in figure 6. While for the first layer a superior energy resolution could be obtained with the KETEK array, the energy resolution of the third layer was observed to be superior when obtained by the Hamamatsu MPPC array. As an explanation for the difference of the energy resolution of the first layer (where basically no influence of Inter-Crystal_Scattering is present) for the two SiPM and MPPC arrays, respectively, can be given by (i) a higher geometrical array fill factor (82% versus 74%) of the KETEK PA3350WB-0808 that allows to detect a larger fraction of the scintillation light and (ii) while a comparable photon detection efficiency (PDE) of both photosensor arrays of about 50% at the peak sensitivity wavelength is reported (HAMAMATSU Photonics K.K 2020, Wiest 2020, KETEK GmbH 2020c), the KETEK SiPM array with a photon detection peak sensitivity at 430 nm matches almost exactly the peak emission wavelength of the LYSO crystal (Saint-Gobain Crystals 2020b, Wiest 2020, KETEK GmbH 2020a, 2020c), while the sensitivity of the HAMAMATSU MPPC array reaches its peak at 450 nm (Saint-Gobain Crystals 2020b, HAMAMATSU Photonics K.K 2020), resulting in a reduced PDE at the peak emission wavelength of the LYSO crystals. For the strong deterioration of the energy resolution of third layer crystal when measured with the KETEK SiPM array, no straightforward explanation can be given. However, it seems that there is a contribution of light leakage. Due to the slightly larger area of the KETEK SiPM, it exceeds the crystal block. This causes the Teflon wrapping at the contact region between crystal and SiPM to be less efficient for the KETEK array than for the HAMAMATSU MPPC array. This might cause more scintillation photons to be lost at the edges.

Inter-crystal scattering (ICS)
Intuitively, the best energy resolution is expected to be obtained from the third layer, since these crystals are coupled closely to the photosensor (i.e. with only a thin layer of optical grease inbetween) and experience the least scintillation light loss due to absorption when passing through the crystal. However, the energy resolution degrades from the first to the third layer, while also the inter-crystal variability increases. This phenomenon can be explained with ICS. The energy spectra of the second, and even more prominent the ones of the third layer, do not show a Gaussian shaped photopeak anymore. This phenomenon is dominant in the third layer (due to the poorest signal-to-background ratio), therefore in the following the focus is set on the third layer (ICS is considered to be one source of background). The arguments given in this section, however, are also valid for the second layer.
A peak in the energy spectra of third layer crystals may have three potential origins. (i) The trivial origin is a full photo absorption of an impinging γ-ray resulting in the photopeak. (ii) A second origin could be related to photo absorption occurring in one of the above crystal layers, while appearing as well in the energy spectrum of a 3rd layer crystal. However, this can be in practice excluded, due to the defined localization of such events (clear crystal response in the flood map) and the ability to assign photo absorption clearly to the corresponding crystal. (iii) A third effect, however, results in generating a background contribution to the actual crystal energy spectra. ICS followed by full absorption of the scattered gamma in a different crystal results in a full absorption peak in the energy spectrum. Due to the light distribution measured at the photosensor, which corresponds to the response to the interactions in two crystals, the calculated interaction position cannot be assigned clearly to one specific crystal and thus forms a uniform background throughout the flood map. A region within the flood map that contains only ICS and background in the corresponding energy spectrum is exemplarily shown by the red triangle in the central right panel of figure 4. Such an energy spectrum exhibits a full absorption peak, as displayed by the ICS spectrum in figure 7 (top row, left panel). The energy spectrum of the third layer consists of the actual gamma source spectrum plus the ICS spectrum. The photopeak and the full ICS absorption peak in these spectra are not found in the same ADC channel region, since the effectively registered amount of scintillation light originating from the full ICS absorption peak and detected by the photosensor is reduced compared to a primary photo absorption event due to absorption of the scattered photons along their trajectory through the crystals. Figure 7 (central row) shows the peak positions of three neighboring crystals in the first, second and third layer of the DOI detector when read out by the KETEK PA3350WB-0808 array (left) and HAMAMATSU Figure 6. Relative energy resolution measured at 511 keV with the KETEK PA3350WB-0808 (left) and HAMAMATSU S14161-3050HS-08 (right) shown for the three crystal layers and the total energy resolution of central and edge crystals given as the median value obtained by all evaluated individual crystals (light gray shaded column). The measured energy values are shown for crystals at the edges and the central region of the detector block. A deterioration of the relative energy resolution from the first to the third layer is observed, while no deterioration towards the edges is observed within the measurement uncertainties. S14161-3050HS-08 array (right), respectively. The scintillation photon absorption can clearly be seen in the shift of the photopeak positions. Using a back-side irradiation, a clean energy spectrum with suppressed ICS can be observed in the third layer energy spectra, because in this scenario no preceding layers have to be traversed by the photon. However, a comparison of energy spectra measured in the first layer from a front-side irradiation and the third layer from back-side irradiation still shows a shift of the photopeak due to scintillation light loss by absorption of photons originating in the first layer ( figure 7 (bottom row, left panel)).
In order to consolidate this explanation, in a back-side irradiation the ICS full absorption peaks should be visible in the energy spectra of first layer crystals. However, in this scenario the photopeak corresponds to the leftmost peak (due to absorption) and the peak structure on the right if this peak corresponds to the full absorption of inter-crystal scattered photons. Figure 7 (bottom row, right panel) shows a 22 Na energy spectrum detected in a first-layer crystal obtained from a front-side irradiation (solid red) superimposed to a 22 Na energy spectrum of the same crystal, but from a back-side irradiation including the ICS peaks (red dashed histogram). Top row: ICS energy spectrum (from a 22 Na irradiation) with a non-Gaussian-shaped full absorption peak around 400 keV in a background region (however, this effect is also present as a background contribution of the spectra obtained from any other region) of the flood map (left panel) and 22 Na energy spectrum measured by a crystal in the first layer without ICS distortions, exhibiting and a Gaussian-shaped 511 keV photopeak (right panel). The 176 Lu-background is present in both spectra. Central row: energy spectra of a 22 Na source placed in front of layer 1 and measured by a crystal in the first (red), second (blue) and third (yellow) layer of the DOI LYSO detector using the KETEK PA3350WB-0808 (left panel) and HAMAMATSU S14161-08 (right panel) SiPM arrays. Bottom row left: energy spectra of a 22 Na source measured with the DOI detector read out by the KETEK PA3350WB-0808. The spectrum drawn in red was detected by one LYSO crystal in the first layer of the DOI detector when irradiated from the front side, while the yellow spectrum was registered by the crystal at the same position, but in the third layer, irradiated from the back side. Right: energy spectra measured with one crystal in the first layer of the DOI LYSO detector and read out by the KETEK PA3350WB-0808 array from a front irradiation (solid color) and back irradiation (dashed histogram) with a 22 Na source.

Analysis of time spectra
Due to the individual SiPM readout and signal processing provided by the PETsys ASIC, an arrival timestamp is assigned to each of the N firing SiPM channels, which is triggered by an initial γ-ray hit in one of the DOI LYSO detector crystals. In general, any of the N provided arrival timestamps belonging to a photon interaction can be used to determine the CRT. A reasonable choice which timestamp to use is either using the timestamp of the first registered firing SiPM channel or the timestamp of the brightest firing channel, respectively. Using the timestamp of the first firing channel is an obvious choice. However, using the timestamp of the brightest firing SiPM may also be beneficial, due to the steeper slope of the rising edge of the scintillator pulse for high amplitude signals and thus a smaller time walk for signals with non-identical amplitudes. Both methods have been investigated. The arrival time difference was calculated with respect to any channel of an 8×8 (64 channels) LYSO crystal array coupled to a PA3325WB-0808 SiPM array (on-by-one coupling) as described in section 'materials and methods' and illustrated in figure 2. Figure 8 (top row) shows the arrival time difference spectra measured with the KETEK PA3350WB-0808 (left) and HAMAMATSU S14161-3050HS-08 (right) with the timestamp of the first detected channel of the DOI LYSO detector used for the calculation. Displayed are inclusive spectra, which contain all arrival time differences from the different combinations between DOI LYSO detector channels and the channels of the reference detector. Three observations can be made: (i) a central peak corresponding to the actual arriving time difference distribution and two satellite peaks at about ±5000 ps are visible. (ii) The satellite peak at the positive side of the time axis has a tail reaching at least down to the central peak. (iii) The central peak is asymmetric, i.e. non-Gaussian shaped. The deviation from a Gaussian distribution is caused by small deviations between the signal transit times and data processing times of the contributing individual crystals. As a result, the centroids of the Figure 8. Top row: measured arrival time difference spectra of the DOI LYSO detector read out by the KETEK PA3350WB-0808 array (left) and the HAMAMATSU S14161-3050HS -08 MPPC array (right). The CRT was determined with respect to the first firing channel of the DOI detector. Bottom row: measured arrival time difference spectra of the DOI LYSO detector read out by the KETEK PA3350WB-0808 array (left) and the HAMAMATSU S14161-3050HS -08 MPPC array (right). The CRT was determined with respect to the brightest firing channel of the DOI detector. arrival time difference distributions of the individual channels are shifted with respect to each other and cause a non-Gaussian distribution in the inclusive spectrum.
In figure 8 (bottom row) the inclusive arrival time difference spectra for the KETEK PA3350WB-0808 (left) and the HAMAMATSU S14161-3050HS-08 (right) are plotted with the timestamp of the DOI LYSO detector taken from the brightest firing channel. While the satellite peaks and the non-Gaussian shape of the central peak can still be observed, the tail of the satellite peak at the positive time axis disappears. The same difference between the respective arrival time difference spectra for the two timestamping methods can also be observed for individual channel pairs, i.e. one specific channel of the DOI LYSO detector versus a specific channel of the reference detector.
The presence of the satellite peaks in the time-difference spectra was initially evaluated by Schug et al (2018) using the TOFPET2 ASIC and a one-to-one coupling between a KETEK PM3325WB and a LYSO crystal. They could show that in particular for low values of the first timing threshold (V th_T1 ), a dark count could trigger the corresponding discriminator. In cases where the trigger logic is configured such that a delay is applied to the logic signal do_T1, the dark count would be responsible for the generation of a timestamp, even though it would not exceed the other two thresholds V th_T2 and V th_E. If this ASIC channel subsequently detects a real γ hit that triggers V th_T2 and V th_E within the delay time of do_T1, the time stamp of the previously dark count would be assigned to this γ hit. These events would cause the assigned time stamp to be wrong by Δt=t delay_T1 -(t T2 −t T1 ), with the constant term t delay_T1 (Schug et al 2018), resulting in satellite peaks at roughly ±5.5 ns around the central peak for t delay_T1 =5.8 ns (which is the default setting and was also used for the studies presented in this manuscript). However, both satellite peaks studied in Schug et al (2018), had the same amplitude, whereas in the time-difference spectra presented in this study the satellite peak on the negative time axis exhibits only about roughly one tenth of the amplitude of the satellite peak on the positive time axis. This is a consequence of the light-sharing readout of the DOI LYSO crystal block. While for a one-to-one coupling configuration (as used for the reference detector) the dark count has to be detected in the same ASIC channel that also detects the subsequent hit, in the light-sharing approach the light emitted in one crystal triggers on average 10-11 ASIC channels. As a consequence, only in one of these ASIC channels a SiPM dark count needs to be detected in order to cause the satellite peak, resulting in about 10 times more events in the satellite peak corresponding to the assignment of a wrong timestamp in the DOI detector.
Furthermore, this observation also explains why the time-difference spectra obtained by the brightest firing channel show less entries in the satellite peaks compared to the one obtained by the first firing channel. While naturally for the first firing channel all of these wrongly assigned timestamps are included in the time-difference spectra, in the case of the brightest firing channel only channels which are the first and brightest firing channel can cause such a wrong timestamp to be used.

Coincidence resolving time
The CRT of the DOI detector depends on the detection time of individual gamma hits in the detector and therefore on the timestamps assigned to these hits. As already described in the previous section, the applied signal processing electronics operates all input channels individually, i.e. in case of an operation mode involving light sharing, one timestamp per electronic channel is provided. The first ASIC-generated timestamp was taken for the CRT calculation, but also a comparison to the obtained CRT using the timestamp of the brightest channel will be given.
In order to determine the CRT of the DOI LYSO detector the procedure was as follows: The CRT of two identical reference detector arrays (3×3 × 5 mm 3 LYSO crystals one-by-one coupled to the (KETEK PA3325WB-0808 SiPM array, see figure 2 was measured and the time resolution (ΔT ref ) was calculated (under the assumption that the time resolution of the identical reference detectors (ΔT ref ) is also identical) according to Since any of the 64 SiPM/MPPC channels can be the first one detecting photons, consequently providing the trigger timestamp, non-uniformities in the matching of the timestamps between different ASIC and/or photosensor channels lead to a deterioration of the CRT, due to a time jitter of the mean value of the Gaussian arrival-time-difference distributions of individual channels, resulting in a broadening of the overall Gaussian.
Instead of using the timestamp provided by the first firing SiPM of the DOI LYSO detector, the timestamp of the brightest channel was also used to determine the CRT of the detector. A comparison of the obtained results of the two methods is given in table 2. No statistically significant difference of the CRT can be seen between these two methods when the DOI crystal block is readout by the HAMAMATSU array (p-value=0.299 at a significance level of 0.05). For the readout with the KETEK SiPM array, the two methods show a different result (p-value<0.000 02), with the method of using the first firing channel being the superior method. However, for the HAMAMATSU MPPC array the statistical uncertainty decreases by about 50% when using the brightest channel as timestamp instead of the first channel. Also, the inter-crystal variability in this case decreases by approx. 25%.
In addition to the CRT with respect to the individual readout channels of the full detector block, the CRT was also determined for individual crystals of the three respective detector layers. In this case it was additionally required that for a specific crystal the first firing channel (providing the triggering timestamp) had also to be the brightest one within this photon event. This condition is fulfilled for 57.7% (KETEK) and 49.6% (HAMAMATSU) of all events. This requirement was applied to eliminate the connecting tail in the time difference distribution when using the first firing channel as timestamp (see figure 8). This tail may otherwise introduce a significant error to the Gaussian fit in cases where the statistics is low, which is especially the case when evaluating the CRT of crystals in the 3rd layer.
Using this method, the overall CRT measured with the KETEK PA3350WB-0808 array was found to be 490 (±44 stat ±82 var ) ps. However, when disentangled into the CRT values for the individual crystal layers, on average the CRT slightly improves from the first over the second to the third layer for the KETEK SiPM array. For the HAMAMATSU MPPC array no influence of the layer position on the CRT can be observed within the measurement uncertainties. However, the statistical uncertainty decreases with larger distance to the SiPM array, due to a higher detection probability in the first and second layers compared to the third layer. The CRT of crystals belonging to the three layers is CRT 1st layer =532 (±30 stat ± 81 var ) ps, CRT 2nd layer =463 (±31 stat ±77 var ) ps and CRT 3rd layer =447 (±42 stat ±69 var ) ps, respectively. The measured CRT with the HAMAMATSU S14161-3050HS-08 is CRT 1st layer =402 (±18 stat ± 28 var ) ps, CRT 2nd layer =392 (±22 stat ± 32 var ) ps and CRT 3rd layer =408 (±36 stat ± 160 var ) ps. These found results are summarized and compared to the values of the overall detector block without layer separation in figure 9 and listed in table 3. For all three types of CRT determination (via first firing channel, brightest firing channel and individual layers), the CRT values determined for the two photosensors are in good agreement with each other within the 1σ uncertainties for the total CRT calculated by using the first firing channel, and for the CRT of the second and third layer crystals, respectively. Figure 9. Left: measured CRT of the DOI LYSO detector with KETEK PA3350WB-0808 and HAMAMATSU S14161-3050HS-08 SiPM readout, respectively. The box charts belonging to the CRT distributions of the (up to) 64 channels indicate the interquartile range (25th percentile to 75th percentile) by the box size, the mean value (black square), the median (horizontal line), and the 1σ standard deviation (whisker) of the mean value. Right: comparison of the CRT between KETEK and HAMAMATSU SiPM arrays. The values given for comparison are the combined CRT for all layers with the timestamps provided by the first firing channel (dark shaded gray)(as given in the left panel) and the brightest firing channel (light shaded gray shaded) and for the three individual DOI layers. The error bars correspond to 1σ standard deviation in the measured CRT values for all investigated crystals in the respective layer.
Furthermore, the evaluation of the time difference spectra of individual crystals shows a shift of the mean value of the individual layers between 50 and 70 ps (varies for the LYSO detector from crystal to crystal) as shown in figure 10. For a refractive index of n=1.81 (at the emission wavelength of 420 nm) (Saint-Gobain Crystals 2020b) this corresponds to a crystal thickness of 10±2 mm, which is in reasonable agreement with the actual crystal thickness of 6.6 mm, taking into account that the scintillation light path is not straight from the interaction point to the photosensor, but is prolonged by various reflections.

Discussion
The presented studies give a performance evaluation of a SiPM/MPPC array readout for a three-layer high resolution PET scintillation crystal block. For the studies presented in this work SiPM/MPPC arrays from two manufacturers (KETEK and HAMAMATSU) have been investigated and compared. For both photosensors a flood map with clearly resolved and distinguishable individual crystals could be acquired, also proving a spatial resolution of 1 mm in the x-y plane given by the crystal pitch and a DOI resolution of 6.6 mm (corresponding to the crystal layer thickness) of the DOI detector.
The overall relative energy resolution measured with the KETEK PA3350WB-0808 at 511 keV was 19.1% and 21.2% with the HAMAMATSU S14161-3050HS-08. By determining the relative energy resolution of individual crystals, the average energy resolution of the three layers of 16.5 (±0.5 stat ± 1.8 var )%, 20.9 Table 2. Summary of the CRT values measured throughout this paper for the DOI LYSO detector read out by the two investigated SiPM arrays. The CRT was calculated by using the timestamp of the first firing DOI detector channel (first row) and by using the brightest channel (bottom row).

Timestamps provided by:
Coincidence resolving time (CRT) Coincidence resolving time (CRT) KETEK PA3350WB-0808 HAMAMATSU S14161-3050HS -08 First firing SiPM 678 (±13 stat ± 124 var ) ps 613 (±12 stat ± 148 var ) ps Brightest SiPM 847 (±16 stat ± 211 var ) ps 637 (±10 stat ± 97 var ) ps  Figure 10. Time difference spectra of neighboring single crystals from the first layer (red), second layer (blue) and third layer (yellow) of the DOI LYSO detector measured with the KETEK PA3350WB-0808 (left) and HAMAMATSU S14161-3050HS -08 (right) arrays. The shift of the mean value of the respective distributions is clearly visible. Furthermore, it is evident that the interaction probability for 511 keV γ rays is the highest in the first layer and decreases with further layers.
All measurements and evaluations performed within the scope of this work were performed with SiPMs biased at the overvoltage recommended by the manufacturer. A performance evaluation of the DOI LYSO detector with respect to the applied SiPM overvoltage has not been done and is not in the scope of this paper. However, it could very well affect the relative comparison between the two different SiPM array types under study.

Conclusion
Throughout the studies presented in this work it was shown that the DOI LYSO detector module consisting of a staggered three-layer LYSO detector block with 1 mm crystal pitch, SiPM readout and highly integrated ASISbased signal processing fulfills all requirements demanded by modern high-resolution PET scanners.
The tested SiPM arrays from KETEK and HAMAMATSU, respectively, generally show a comparable performance for the detector module. However, while the energy resolution measured with the KETEK SiPM array is slightly superior for the first layer to the one measured with the HAMAMATSU MPPC array, the latter shows a superior energy resolution for the third layer. Also, the CRT that was obtained with the HAMAMATSU MPPC array was superior to the one obtained with the KETEK SiPM array.
The PETsys TOFPET v2c ASIC and its accompanying front-end electronics boards provide a suitable signal processing and data acquisition for the SiPM arrays of both manufacturers.