Sub-3 mm, near-200 ps TOF/DOI-PET imaging with monolithic scintillator detectors in a 70 cm diameter tomographic setup

Recently, a monolithic scintillator detector for time-of-flight (TOF)/depth-of-interaction (DOI) positron emission tomography (PET) was developed. It has a detector spatial resolution of ~1.7 mm full-width-at-half-maximum (FWHM), a coincidence resolving time (CRT) of ~215 ps FWHM, and ~4.7 mm FWHM DOI resolution. Here, we demonstrate, for the first time, the imaging performance of this detector in a 70 cm diameter PET geometry. We built a tomographic setup representative of a whole-body clinical scanner, comprising two coaxially rotating arms, each carrying a detector module, and a central, rotating phantom table. The fully automated setup sequentially acquires all possible lines of response (LORs) of a complete detector ring, using a step-and-shoot acquisition approach. The modules contained 2  ×  2 detectors, each detector consisting of a 32 mm  ×  32 mm  ×  22 mm LYSO crystal and a digital silicon photomultiplier (dSiPM) array. The system spatial resolution was assessed using a Na-22 point source at different radial distances in the field-of-view (FOV). Using 2D filtered back projection (2D FBP, non-TOF), tangential and radial spatial resolutions of ~2.9 mm FWHM were obtained at the center of the FOV. The use of DOI information resulted in almost uniform spatial resolution throughout the FOV up to a radial distance of 25 cm, where the radial and tangential resolution are ~3.3 mm FWHM and ~4.7 mm FWHM, respectively, whereas without DOI the resolution deteriorates to ~9 mm FWHM. Additional measurements were performed with a Na-22 filled Derenzo-like phantom at different locations within the FOV. Images reconstructed with a TOF maximum-likelihood expectation-maximization (TOF ML-EM) algorithm show that the system is able to clearly resolve 3 mm diameter hot rods up to 25 cm radial distance. The excellent and uniform spatial resolution, combined with an energy resolution of 10.2% FWHM and a CRT of ~212 ps FWHM, indicates a great potential for monolithic scintillators as practical high-performance detectors in TOF/DOI-PET systems.

2 G Borghi et al spatial resolution throughout the field of view (FOV) (Surti et al 2013, Thoen et al 2013. At present, benefits of DOI estimation have been experimentally demonstrated only in preclinical PET scanners or demonstrators (González et al 2016, Lee et al 2017, whereas there are no results available for PET scanners or demonstrators having ring dimensions compatible with clinical applications. Current clinical whole-body scanners are usually built with detectors based on segmented scintillator-crystal matrices whose pixels typically have edge dimensions of 4-5 mm and a thickness of 20-30 mm. These crystal arrays are read out using either light sharing techniques or one-to-one coupling, by means of arrays of photomultiplier tubes (PMTs), position sensitive PMTs (PS-PMTs), or arrays of silicon photomultipliers (SiPMs). In state-of-the-art clinical scanners these detectors typically achieve a CRT of 300-350 ps full-width-at-halfmaximum (FWHM) and an energy resolution of about 10.5%-11% FWHM, but none of them provides DOI estimation (Vandenberghe et al 2016).
Prototype pixelated scintillation detectors with improved performance compared to the clinical state-of-theart have been presented in literature. These detectors are typically based on smaller crystals (~2 mm pixels) (Marcinkowski et al 2014, Lee 2015, Wong et al 2015), while numerous designs have been proposed to obtain DOI estimation. For instance, 3D crystal arrays (Yamaya et al 2011), stacked matrices made of different materials with different decay times (phoswich) (Schmall et al 2015, Chen-Ming 2017, phosphor-coated crystals with pulse-shape classification (Kwon et al 2016), light sharing techniques that encode the DOI in the width of the registered light distribution , and crystals with double-sided-readout (DSR) Schaart 2014, Kang et al 2015). Often, such improvements come at the expense of other key performance parameters, such as time resolution and energy resolution, and/or result in increased cost and complexity.
An interesting alternative to pixelated detectors is represented by monolithic scintillator detectors, in which a continuous scintillator (typically a cuboidal crystal) is coupled to an array of photosensors. Monolithic scintillator detectors have shown an excellent combination of spatial resolution, coincidence resolving time, and detection efficiency (Seifert et al 2012. For example, our group has recently presented a detector based on a 32 mm × 32 mm × 22 mm LYSO:Ce crystal and a digital photon counter (DPC) array which achieved a x-y spatial resolution ~1.7 mm FWHM/~1.6 mm mean absolute error (MAE), a DOI resolution ~3.7 mm FWHM/~2.2 mm MAE, and a CRT of ~215 ps FWHM (Borghi et al 2016). Also, new methods have been proposed to increase the efficiency of the calibration procedures and of the time and position estimation algorithms (Borghi et al 2015), which had previously been considered a limitation for the practical application of this technology.
Thus, monolithic detectors are becoming a serious option for clinical TOF/DOI-PET scanners. Nevertheless, monolithic and pixelated scintillator detectors have inherently different responses and it is not possible to simply compare e.g. the FWHM spatial resolution of a monolithic detector with the crystal pitch of a pixelated detector (Tabacchini et al 2017). Hence, realistic tomographic image acquisition measurements are warranted to assess the performance of monolithic scintillator detectors at the full-scanner level.
In this work we aim to experimentally predict the imaging performance of a 70 cm diameter whole-body TOF/DOI-PET scanner based on monolithic scintillator crystals coupled to DPC arrays using a new tomographic setup, which has been developed to perform full PET acquisitions using only two detector modules. The setup comprises two rotating arms, each one carrying a detector module, as well as a central rotating phantom table. The setup can be used to sequentially acquire all possible combinations of detector modules of a complete PET ring, using a step-and-shoot acquisition procedure. A detailed description of the detector modules and of the experimental setup is presented in section 2, while the time, energy, and spatial resolution obtained in the tomographic experiments are reported in section 3.

Module design and sensor settings
The two PET modules used in this work are based on DPC3200-22-44-M22 sensor modules produced by Philips Digital Photon Counting (PDPC) (figure 1). Each module hosts 2 × 2 digital silicon photomultiplier (dSiPM) arrays (model DPC3200-22-44) at a distance of ~0.5 mm one from each other and dedicated readout electronics. The DPC arrays measure 32.6 mm × 32.6 mm and are composed of 4 × 4 independent sensors (dies) at a pitch of 8 mm. Each die is subdivided into four pixels and is equipped with an on-silicon integrated TDC system. When a die acquires a light signal, it registers the intensity of the light on each pixel, i.e. the exact number of fired single photon avalanche photodiodes (SPADs), and a single timestamp. Therefore, DPC arrays can provide an 8 × 8 pixel light distribution and 16 timestamps when read out fully. More detailed descriptions of the DPC array working principle can be found in Frach et al (2009), Frach et al (2010, Schaart et al (2016) and Schulze (2014b). The modules were controlled and read out by a field-programmable gate array (FPGA) based electronic board and a computer software provided by PDPC.
Eight monolithic scintillator detectors based on 32 mm × 32 mm × 22 mm LYSO:Ce crystals (Crystal Photonics, Sanford, USA) and DPC3200-22-44 arrays were assembled to fully populate the modules. The lateral faces of the crystals were covered with a specular reflector foil (Vikuiti ESR, 3M), while the top faces were covered with Teflon tape. The crystals were permanently glued to the photosensors using UV-curing glue (DELO Photobond 4436).
The sensor settings (Schulze 2014a) were optimized to maximize the timing performance and the sensitivity for monolithic scintillator detectors operated at reduced temperatures (−20°/−15°). On all sensors, the noisiest 10% of the SPADs were disabled. The trigger threshold was set at the first photon interacting on a die (DPC notation: MT = 1), the validation interval was set to 40 ns, the validation threshold setting required at least one photon per pixel to validate an event on a die (DPC notation: 0x7F:AND) and the integration interval equaled 165 ns. Considering also the readout and reset time (680 ns), DPC sensors have a total dead time per acquired event <1 µs. However, the bandwidth of current DPC sensor tile is artificially limited by the event storage memory available on the tile FPGA: the present sensor tile can handle a maximum of about 120 kcps per chip, i.e. each monolithic crystal detector used in this work has a maximum theoretical count rate of ~120 kcps. The neighbor logic was activated so that all dies on a sensor were acquired every time that one of the dies registered a validated event.
A hardware gating signal was defined between the modules using the FPGA electronic board to register only the events occurring within a coincidence window of ~100 ns. This broad coincidence pre-selection requirement was imposed with the sole purpose of reducing the single-event rate registered on the sensors and to minimize the amount of data to be transferred to the data acquisition computer through a USB 2.0 connection (Schulze 2014a). A more accurate coincidence selection was performed during the offline analysis (see section 2.3).

Time and position estimators for monolithic scintillator detectors
Monolithic scintillator detectors require individual calibration and optimized estimation algorithms to achieve their best spatial and timing performance. In Borghi et al (2015Borghi et al ( , 2016, new calibration methods were presented that made it possible to fully characterize the response of a detector in few hours using fan-beam and flood irradiations. In the same papers, several modifications to the statistical estimators for the position-, energy-, and time-of-interaction were introduced to accelerate these estimators and make them practically applicable. The same methods and estimators were used in the present work to calibrate the detectors on the PET modules and to process the acquired data. In particular, for x, y position estimation an accelerated version of the k-nearest neighbor (k-NN) method was used, which requires the acquisition of a reference dataset for each single detector using a fan-beam irradiation. The same datasets were used to calculate the look-up-tables (LUTs) needed to correct for the position-dependent energy response of the detectors and to estimate the DOI. The spatial sampling used to estimate the position of interaction inside the crystal was 0.25 mm × 0.25 mm × 1 mm in the x, y and DOI direction, respectively.
A separate calibration dataset, acquired by means of a simple flood irradiation, was used to calibrate the maximum likelihood interaction estimation method (MLITE), the statistical method used to estimate the time of interaction (van Dam et al 2013). A short description of the calibration procedures and estimators is given in appendix. In particular, the small adaptations that were necessary to adjust the methods to the measurement conditions used in this work are highlighted there. For a complete description of the methods, the reader is referred to the papers cited.

Tomographic setup
A new setup was designed and built to perform complete tomographic acquisitions using only two detector modules (figure 2). This setup is constructed such that PET scans with different axial extents and with the modules at different radial positions can be acquired. The setup achieves a module positioning accuracy ⩽0.2-0.3 mm even at the largest radii, thus preventing image artifacts or degradations in the imaging performance that might otherwise occur from inaccurate mechanical precision.
The setup is based on a precision goniometer originally developed for an x-ray diffractometer (Bruker) with two coaxially rotating stages capable of a full rotation. The rotational accuracy is <0.01° and each rotating stage is equipped with a mechanical arm. On each arm, two perpendicular, custom-made linear stages driven by stepper motors are used to support and move a PET module, such that its radial position r and height z can be adjusted. The range of the r and z linear stages are 450 mm and 250 mm, respectively. The positioning accuracy of both stages is <0.05 mm. At the center of the goniometer, an additional rotating stage (accuracy ~0.05°) is mounted, so that the central platform onto which the sources are placed can be rotated. All stages make use of closed-loop control. A dedicated control unit based on a CNC protocol was developed that interfaces with a PC, allowing the experimenter to easily control the setup (Berkelaar MRT, Delft).
Precise alignment of the modules was obtained by means of fine-adjustment mechanisms that allow correction of the unavoidable small inaccuracies that occur when large machineries are assembled. Alignment was performed using a 3D measuring arm (Romer Absolute Arm, 6 axis) that could determine the spatial position of the various components with a precision of <0.02 mm. The detector modules are mounted using a high-precision docking system that defines their position accurately even if they are removed and re-mounted.
The setup also includes a two-stage cooling system designed to control the temperature of the DPC sensors and to dissipate the heat produced by the detectors during operation. The first cooling stage is a liquid cooling machine (Integral XT 150, LAUDA) used to stabilize the temperature of a cooling plate at the back of each module at ~15 °C. The second cooling stage consists of high-power Peltier elements mounted in between the cooling plates and the modules, which can cool the modules down to about −20 °C. The modules and Peltier elements are enclosed by insulating foam boxes purged with dry nitrogen to avoid condensation of moisture.
Complete tomographic acquisition are performed acquiring all module pair positions (or views) that are measured in a full 3D PET ring, i.e. all module pair positions that define LORs that intersect the FOV. For a given module pair position, possible LORs include all lines connecting all the points on the surface of the detectors of one module with all the points on the surface of the detectors of the other module. A given view can be obtained by moving either the arms, the central plate, or both. For each view a measurement of the same time duration has to be registered. A LabVIEW program was developed to perform such acquisition sequences in a completely automated way. No angular constraints are applied on the possible LORs, so tomographic measurements are performed in full 3D PET mode. Coincidences are saved in list mode.
For the tomographic acquisitions performed in this work, the inner diameter (module-to-module maximum distance) of the system was set to 70 cm, thus representing a whole-body clinical scanner made of 32 modules. The FOV has a diameter of 50 cm and an axial extent of 6.5 cm, corresponding to the size of a single detector module. A representation of the PET ring emulated with the tomographic setup is shown in figure 3.
Coincidence events were selected offline, using a 4 ns wide software coincidence window and an energy window corresponding with the full-width-at-ten-maximum (FWTM) of the uncorrected 511 keV photopeaks of the different detectors, which on average corresponds to a ~100 keV wide energy window (see appendix, section A.1).

System characterization 2.4.1. Energy resolution and coincidence resolving time
The energy resolution and CRT of the tomographic setup were determined from a complete tomographic acquisition of a 22 Na point source (~3 MBq) placed at the center of the FOV. The energy resolution was calculated by creating, for each of the 8 detectors, a calibrated spectrum corrected for the position-dependent detector response (see appendix) and summing the obtained spectra. Detector calibration was performed as a simple linear calibration using only the 511 keV photopeak as a calibration point. The single-detector and overall energy resolutions were determined with Gaussian fits.
A system TOF-difference spectrum was obtained by calculating the differences between the detector timestamps for all coincidences acquired during the measurement. Detector timestamps were corrected for die and tile electronic skews, which were determined as described in appendix. No further correction to align the spectra obtained at the different combinations of module positions was performed. The system CRT was determined with a Gaussian fit of the spectrum.

Spatial resolution measurements with point sources
Complete tomographic acquisitions of a single 22 Na point source (Ø 0.5 mm, ~3 MBq) were performed at different positions along the x-axis of the scanner, from the center of the FOV up to a radial offset of 25 cm, at a pitch of 5 cm. A scan time of ~3 min per module position was used at each radial offset, which resulted in a total number of 4.5 to 5.5 million coincidences per scan, depending on the source position.
The images were reconstructed projecting all the LORs on a single plane and using 2D filtered back projection without considering TOF information (2D FBP, non-TOF), with pixel dimensions of 0.5 mm × 0.5 mm, with no smoothing or apodization. The LORs were rebinned with an angular sampling of 0.5° and a radial sampling Each module has an active area of ~65 mm × 65 mm. The inner diameter of the scanner is 70 cm (module-to-module distance), while the angular distance between modules is 11.25°. The FOV has a diameter of 50 cm and an axial extent of 6.5 cm. of 0.25 mm. The radial and tangential resolutions were defined as the FWHM of a line profile across the reconstructed point sources in the two directions, respectively. The FWHMs were calculated using spline interpolation of the 1D line profiles through the point sources. In order to investigate the accuracy and the impact of DOI estimation, images were reconstructed both with and without DOI information incorporated in the reconstruction. In case no DOI estimation was used, the center of the crystal was used as a fixed value for the DOI.

Spatial resolution measurements with Derenzo-like phantom
A custom, Derenzo-like resolution phantom was built using PMMA (polymethyl methacrylate) plastic to qualitatively assess the spatial resolution and the imaging performance of the PET scanner (figure 4). The hotrod insert of the phantom has a diameter of 10 cm, a height of 7 cm and is subdivided into six sectors, with rod diameters of 2.5 mm, 3.0 mm, 3.5 mm, 4.0 mm, 5.0 mm, and 7.0 mm. The distance between the centers of adjacent rods within each sector always equals twice the rod diameter. A total activity of ~20 MBq of 22 Na was used to prepare the homogeneous aqueous solution with which the phantom was filled. Considering the total volume of solution, which is partly outside the hot rods, the estimated activity concentration was ~150-200 kBq cm −3 .
The phantom was scanned at three radial distances from the center of the FOV, viz. 0 cm, 15 cm and 20 cm. At each radial distance, the scan time at each module position equaled ~40 min, which resulted in a total number of ~50 million coincidence events for each complete tomographic scan.
Image reconstruction was performed using a list-mode 3D TOF maximum-likelihood expectation-maximization (TOF ML-EM) algorithm (Vandenberghe et al 2000) with a voxel size of 1 mm × 1 mm × 1 mm. Siddon ray tracing was used to calculate the elements of the system matrix (Siddon 1985). No resolution modelling was used for image reconstruction. Normalization, scatter and random correction were not applied, since normalization and random correction are not expected to have an influence on spatial resolution, whereas scatter correction has only a minor effect on the tails of the measured activity distribution. Attenuation correction was applied using an analytical model of the phantom, of which the geometry, composition and position inside the scanner are accurately known. Attenuation correction was part of the image reconstruction, no prior correction was performed on the measured data. A sensitivity map (in the image space) was calculated by backprojecting LORs on the image matrix. In order to calculate the sensitivity map in the image space, a number of 10 000 LORs were randomly sampled for each monolithic detector pair and backprojected on the image matrix, i.e. the end-points of each LOR were randomly sampled in the (continuous) volume of each crystal. All images were obtained with 10 EM iterations, since additional iterations did not provide any significant visual improvement. As only the transaxial resolution is of interest for this measurement, the 50 central image slices were averaged in order to increase count statistics.
All images were reconstructed twice, once using the estimated DOI and once using the center of the crystal as a fixed DOI value.

Energy resolution and coincidence resolving time
The energy resolutions of the single detectors are reported in table 1. All of them lie between 10.6% FWHM and 10.1% FWHM, except for the resolution of detector 1 on module 1, which is slightly worse. An overall system resolution of 10.2% FWHM is found, which is competitive with state-of-the-art detectors used in commercial scanners based on pixelated scintillators.
For comparison, the recent GE Healthcare Signa TOF-PET/MRI scanner, which is based on analog SiPMs, has an energy resolution of 10.3% . For the Philips Vereos TOF-PET/CT scanner, which is based on arrays of 4 mm × 4 mm LYSO:Ce crystals coupled to the same DPC sensors as used in this work, an energy resolution of 11.1% FWHM has been reported (Miller et al 2015). Similar values have been reported for other prototype detectors and scanners based on pixelated crystals coupled to DPC sensors (Degenhardt et al 2012, Marcinkowski et al 2014, Schug et al 2015 or PMTs (Wong et al 2015), which achieved energy resolutions between 10.7% FWHM and 11.4% FWHM. The slightly improved energy resolution of monolithic scintillator detectors compared to pixelated detector, when the same type of photosensor is used, can be explained by the more efficient light collection process due to the favorable aspect ratio of monolithic crystals.
In previous work on a prototype monolithic detector (similar to the ones used in this work) we reported an even better value for the energy resolution, viz. less than 10% FWHM (Borghi et al 2016). The difference can probably be explained by the higher percentage of DPC cells disabled in the present work (10% instead of 5%) and by the slightly larger sensor temperature fluctuations observed during operation (~1 °C instead of 0.1 °C-0.2 °C).
A system CRT of ~212 ps FWHM was obtained using the MLITE method discussed in section 2.2. Figure 5 shows the corresponding time difference histogram. This result is consistent with the CRT values that were measured independently with 8 different couples of monolithic crystals (see section 4 in appendix), which were all between 210 ps FWHM and 220 ps FWHM.
For comparison, CRTs of ~315 ps FWHM and of 385 ps FWHM have been reported recently for the Vereos TOF-PET/CT scanner (Miller et al 2015) and for the Signa TOF-PET/MRI scanner (GE) Levin et al (2016), respectively. For research scanners, a CRT of ~266 ps FWHM was measured by Degenhardt et al (2012) on a prototype PET scanner based on DPC sensors and one-to-one coupled 4 mm × 4 mm × 22 mm LYSO:Ce crystals. In Schug et al (2015), a CRT of ~213 ps FWHM was reported for a prototype small-diameter scanner based on DPC sensors and 4 mm × 4 mm crystals with a thickness of only 10 mm. This comparison indicates that monolithic scintillator detectors can achieve better timing performance than segmented crystals of equal thickness, due to the favorable light collection conditions and the possibility to correct for the influence of the optical transport of the scintillation photons within the crystal, using techniques such as MLITE .
Finally, the same data used to calculate the CRT with MLITE were used to determine the CRT that could be achieved using a simple, analytical method to estimate the time of interaction, viz. using the average of the first two valid timestamps . This approach provides a CRT of ~233 ps FWHM, indicating that it is possible to achieve excellent timing performance even without MLITE. It is noted that this result furthermore implies that the skew correction based on a simple flood irradiation of the detectors (see section A.4) works very well, as this analytical approach is expected to be rather sensitive to any remaining skews.

Spatial resolution 3.2.1. Spatial resolution with point sources
The radial and tangential spatial resolutions measured with point sources and 2D FBP at different radial offsets are reported in table 2 and figure 6 (red squares). A spatial resolution of ~2.9 mm FWHM is measured at the center of the FOV, while very little degradation is observed at radial distances of up to 25 cm, where the radial and tangential resolutions become ~3.3 mm FWHM and ~4.7 mm FWHM, respectively. Such a homogeneous spatial resolution across the entire FOV shows that monolithic scintillator detectors can accurately estimate the DOI.
For comparison, when the same datasets are reconstructed using no DOI information, the spatial resolution rapidly deteriorates as the point source is moved away from the center of the FOV (figure 6, black circles). At a radial distance of 25 cm, the radial and tangential resolutions become ~8.6 mm FWHM and ~9.0 mm FWHM, respectively.
The spatial resolution of our system outperforms the resolution reported for all currently available wholebody clinical scanners. For instance, radial and tangential spatial resolutions of 4.44 mm FWHM and 4.10 mm FWHM at 1 cm and of 8.44 mm FWHM and 5.23 mm FWHM at 20 cm were recently reported for the Signa TOF-PET/MRI scanner (Grant et al 2016). Similar results were reported for the Vereos TOF-PET/CT scanner, which achieves a transverse spatial resolution of 4.11 mm FWHM at 1 cm and of 5.79 mm FWHM at 20 cm, respectively (Miller et al 2015). To our knowledge, a spatial resolution of about 3 mm FWHM in a whole-body clinical PET system has so far only been achieved by Wong et al (2015), who built a prototype scanner based on small LYSO:Ce crystals (2.35 mm × 2.35 mm × 15.2 mm) coupled to PMTs. With this scanner an impressive resolution of ~2.9 mm FWHM and ~3.7 mm FWHM was measured at the center of the FOV and at 24 cm radial distance, respectively. The relatively small influence of DOI in this system can be attributed to its relatively thin (15.2 mm) crystals and relatively large scanner diameter (~87 cm). Unfortunately, both of these factors tend to reduce system sensitivity.

Derenzo-like phantom images
The TOF ML-EM reconstructed image of the Derenzo-like phantom acquired with the phantom positioned at the center of the FOV is shown in figure 7. All rods with a diameter of ⩾3 mm are clearly distinguishable, while most of the rods in the 2.5 mm sector can be recognized.   Because of the accurate DOI information provided by the monolithic scintillator detectors, essentially the same results are obtained when the phantom is imaged at different positions throughout the FOV. Specifically, when the phantom is positioned with its center at 15 cm radial offset (figure 8), the image looks quite similar to figure 7. Even when the phantom is positioned at 20 cm radial offset (figure 9), in which case the smallest rods are located at radial distances of up to ~25 cm, little degradation of image quality is observed. For a visual impression of the benefit of including DOI information in the reconstruction process, figures 7, 8 and 9 show also images of the Derenzo-like phantom reconstructed without (right) DOI correction. As expected, a comparison of the images with and without DOI information demonstrates that DOI information is more important to improve image quality at larger radial distances in the FOV. When the phantom is positioned at the center of the FOV there are no noticeable differences between the DOI and non-DOI images. However, when the phantom is positioned at 15 or 20 cm radial offset only the sectors with the rods having a diameter > 4 mm can be resolved.
The excellent quality of the images obtained with the Derenzo-like phantom confirms the outstanding spatial resolution of the tomographic system measured with the point sources and the capability of the system to maintain an almost uniform spatial resolution throughout the whole FOV thanks to the DOI information. However, a direct comparison between the two results cannot be made, since point-source spatial resolution cannot be reliably estimated using ML-EM reconstruction (Gong et al 2016).

Conclusions
Two complete PET modules, each containing four 32 mm × 32 mm × 22 mm LYSO:Ce crystals read out by digital photon counter arrays, were assembled and fully calibrated. The modules were installed on the two rotating arms of an experimental setup capable of acquiring all LORs of a complete tomographic PET acquisition in a step-and-shoot approach. The system was used to emulate a 70 cm diameter PET ring based on monolithic scintillator detectors and to experimentally estimate its energy resolution, timing performance and spatial resolution.
An average energy resolution of ~10.2% FWHM and a CRT of ~212 ps FWHM were obtained. A spatial resolution of ~2.9 mm FWHM was measured at the center of the FOV, using a single 22 Na point source and 2D FBP reconstruction. The resolution remained almost constant throughout the FOV because of the accurate DOI estimation of the detectors. At the largest radial distance tested (25 cm), the resolution increased to only 3.3 mm FWHM and 4.7 mm FWHM in the radial and tangential direction, respectively.
Furthermore, PET acquisitions of a Derenzo-like phantom were performed and reconstructed using a TOF ML-EM algorithm. The results show that the tomographic setup is able to clearly resolve 3 mm hot rods up to a radial distance of ~25 cm.
We emphasize that no resolution modelling was used in the reconstruction of any of the images presented in this work. Future research will include investigating the potential benefit of modelling the detector response in the image reconstruction. Future work should also include a detailed study of the count-rate capabilities of the detector, which could not be performed with the current experimental setup because of the limitations of the readout electronics. Nevertheless, due to the short dead time of DPC sensors (<1 µs per event) and the bandwidth of the DPC-sensor array (>100 kcps), we do not expect count-rate limitations for monolithic detectors in realistic whole-body clinical PET acquisition conditions. In conclusion, we built a proof-of-concept tomographic setup to investigate, for the first time, the imaging performance of 32 mm × 32 mm × 22 mm LYSO:Ce monolithic scintillator detectors developed at our lab in a 70 cm diameter TOF/DOI-PET geometry. These results are particularly meaningful because no group has yet shown tomographic images experimentally obtained with thick monolithic detectors in a tomographic setup having a diameter compatible with clinical applications. Moreover, to the best of our knowledge, these are the first experimental results (obtained either with monolithic or with pixelated detectors) that show how DOI information can improve the spatial resolution of a clinical PET scanner at large radial distances. The results obtained demonstrate the excellent potential of monolithic scintillator detectors as a practical high-performance detector for clinical TOF/DOI-PET.

Acknowledgments
This work was part of the EU FP7 project SUBLIMA, Grant Agreement 241711; see also www.sublima-petmr.eu.

Appendix. Time and position estimation methods and detector calibration procedures
All calibration measurements were performed in a dedicated setup built in a climate chamber (Weiss WT 450/70), which was cooled to about −20 °C, resulting in a sensor temperature of about −17°/−16° C during operation. All measurements were acquired selecting one detector on one module for calibration and using one detector on Figure 8. TOF ML-EM reconstructed images of the Derenzo-like phantom with its center positioned at 15 cm radial distance, such that the smallest rods are located at radial distances of up to ~20 cm. DOI information is used in the reconstruction of the left-hand image and discarded in the right-hand image. The diameter of the rods in the different sectors is 2.5 mm, 3.0 mm, 3.5 mm, 4.0 mm, 5.0 mm, and 7.0 mm. The pixel intensity values obtained from image reconstruction were normalized from 0 to 1. The image was obtained with 10 EM iterations. The 50 central image slices were averaged in order to increase count statistics. Figure 9. TOF ML-EM reconstructed images of the Derenzo-like phantom with its center positioned at 20 cm radial distance, such that the smallest rods are located at radial distances of up to ~25 cm. DOI information is used in the reconstruction of the left-hand image and discarded in the right-hand image. The diameter of the rods in the different sectors is 2.5 mm, 3.0 mm, 3.5 mm, 4.0 mm, 5.0 mm, and 7.0 mm. The pixel intensity values obtained from image reconstruction were normalized from 0 to 1. The image was obtained with 10 EM iterations. The 50 central image slices were averaged in order to increase count statistics. the other module as a reference. The module containing the detector to be calibrated was assembled on a pair of perpendicular stages driven by stepper motors (Physics Instruments, M-403.42S stages with C-663 controllers) which were used to correctly position it and move it when necessary.

A.1. Data pre-processing
A few pre-processing operations were performed for all measurements done in this work, both the calibration measurements reported in this section and the characterization measurements described in section 2.4.
The first step consisted of selecting only the events in which no pixel values or time stamps were missing and for which the total (uncorrected) deposited energy fell within the full-width-at-ten-maximum (FWTM) of the 511 keV photopeak. Before applying any energy correction, the mean energy resolution of the detectors (which have practically Gaussian photopeaks) is 10.6% (minimum is 10.3%, maximum is 12%), therefore on average the FWTM is ~19.3% (min ~18.8%, max 21.9%), which corresponds to an energy window ~100 keV wide (min ~96 keV, max ~112 keV).
The second step consisted of correcting the values of all light distributions for possible non-uniformities in the DPC array response which could be due to small defects in crystal wrapping or in the optical coupling between the crystals and the sensor. The homogeneity correction maps used for this operation were calculated with the method reported in Borghi et al (2016) using the fan-beam irradiation measurement described in section A.2.

A.2. Spatial response calibration and position estimation algorithms
An accelerated version of the k-nearest neighbor (k-NN) classification algorithm (Borghi et al 2016) was used for estimating the x-y position of interaction of each gamma ray inside the PET detectors. The basic idea of this method is to first perform a rapid but coarse pre-estimation of the 3D position of interaction of an incident gamma ray by using analytical parameters of the light distributions and calibrated LUTs. For this purpose we used (i) the center of gravity (COG) of the light distribution for the x and y position and (ii) the sum of the squared pixel intensities for the DOI. Once a pre-estimated position is obtained, a more accurate estimation is obtained by running the k-NN 1D algorithm on a subset of the reference dataset composed of only reference events whose pre-estimated positions of interactions are close to the pre-estimated positon of the unknown event. Since the k-NN algorithm is run on a significantly smaller subset of the training dataset, the overall computational time is greatly reduced. The same method used for DOI pre-estimation was also used for the final estimation of the DOI, except that the x-y coordinates estimated with the k-NN algorithm were used to build the LUT and to calculate the DOI (van Dam et al 2011).
The acquisition of training (or reference) events for position estimation was performed by irradiating the crystals of the detectors with a narrow fan beam of annihilation photons at precise x or y positions, while uniformly irradiating the crystals along the other direction (Borghi et al 2015). The fan beam was obtained using a 22 Na source (Ø 0.5 mm, ~3 MBq, IDB Holland BV) and a 80 mm tungsten slit collimator with a slit width of about 0.5 mm. For each detector, two datasets of reference events were acquired, one collected along a series of equally spaced positions (0.25 mm pitch) on the x axis (x-subset) and the other one on the y axis (y-subset). At each position, 12 800 full-energy events were collected. These datasets were used to calibrate the pre-positioning LUTs and as x-or y-reference datasets for the accelerated 1D k-NN algorithm. An estimation of the distribution of the x-y-z positions of interaction of the calibration events inside the crystals, which is needed to calculate the position pre-estimation LUTs and the LUT for the final DOI estimation, was obtained with Monte Carlo simulations performed using GATE (Jan et al 2004).

A.3. Energy response calibration
Monolithic scintillator detectors can have small variations in their energy response depending on the position of interaction of the gamma rays inside the crystal. An energy correction LUT was therefore calculated for each detector (Borghi et al 2016) by subdividing the crystal into 2 mm × 2 mm × 5.5 mm 'voxels' and determining a correction factor for each voxel. This correction factor was equal to the ratio between the 511 keV photopeak position of the total energy spectrum and the photopeak position of the spectrum obtained for that voxel. For this calibration procedure, the same datasets acquired with fan-beam irradiation for x, y position estimation were used.

A.4. Electronic skew estimation and MLITE calibration
Detector calibration for the MLITE algorithm , Borghi et al 2016 consisted of two steps. First, the electronic time skews between the different dies on each sensor were estimated. Then, the probability distributions of the first photon detection delays (i.e. the delays between the time of interaction of the gammaray and the time of detection of the first visible photon interacting on a die) for a grid of positions inside the crystal and for each die were determined. For this purpose, a flood irradiation was acquired for each crystal using a 22 Na point source (~3 MBq). The monolithic crystal to be calibrated was placed at a distance of about 390 mm from the source while the distance of the reference monolithic crystal was about 60 mm, so that the coincidence event would be localized in a region of the reference detector smaller than a die. About 20 million coincidence events were acquired for each detector to be calibrated and their position of interaction inside the detectors was estimated using the methods described in section A.2.
The die skews were estimated for each single detector using all events registered during the flood irradiation. To estimate the die skew between two vertically, horizontally or diagonally adjacent dies, all the events interacting in the 8 × 8 mm 2 region in between the dies and at a DOI between 0 mm and 10 mm (from the crystal top face) were selected. The differences of the timestamps acquired by the two dies were then calculated and used to estimate their probability distribution function with kernel density estimation (KDE). The mode of the distribution was then considered an estimation of the skew. Using the relative time distances between adjacent dies, time skews of all dies with respect to a single die were then calculated.
To determine the accuracy of this skew-estimation method, it was tested on a detector similar to the detectors used in this work for which the values of the die skews measured with a laser were available (Borghi et al 2016). The test showed that the average difference between the values obtained with the new method and the values measured with a laser was smaller than 20 ps (i.e. the width of the TDC bins of the DPC sensors), while the maximum differences were smaller than 40-50 ps. This result shows that the new method provides excellent estimation and its accuracy can be considered comparable with the accuracy achievable performing a measurement with a laser.
To estimate the detector FPDDs, four million reference events were selected for each detector among the events acquired during the flood irradiations. The events were selected so that their estimated position of interaction in the reference detector was above the die on which the flood irradiation was electronically collimated and their estimated DOI was comprised between 0 mm and 6 mm (from the crystal top face). The reference timestamps were calculated as the average of the first two timestamps registered on the reference detector and were used to estimate the FPDDs of the sensor being calibrated, following the procedure described in Borghi et al (2016). Using this procedure it was possible to avoid the use of small, fast reference detectors for the MLITE calibration.
To estimate the electronic skews in between tiles and modules, an additional measurement was necessary. First, all detectors of one module were measured in coincidence with a single detector of the other module using an un-collimated 22 Na point source placed exactly in the center between them (~20 cm distance from each detector). The same measurement was repeated after inverting the modules. One million events per detector combination were acquired and the resulting time spectra were calculated. The position of the timing spectra acquired for all eight combinations were used to calculate the skews between detectors.