# Novel ion imaging concept based on time-of-flight measurements with low gain avalanche detectors

F. Ulrich-Pur, T. Bergauer, A. Hirtl, C. Irmler, S. Kaser, F. Pitters, Simon Rit

► To cite this version:

F. Ulrich-Pur, T. Bergauer, A. Hirtl, C. Irmler, S. Kaser, et al. Novel ion imaging concept based on time-of-flight measurements with low gain avalanche detectors. 23rd International Workshop on Radiation Imaging Detectors, Jun 2022, Riva del Garda, Italy. pp.C02062, 10.1088/1748-0221/18/02/C02062. [https://doi.org/10.5281/zenodo.6536362](https://doi.org/10.5281/zenodo.6536362)

HAL Id: hal-03815221

https://hal.science/hal-03815221

Submitted on 9 Mar 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
23rd INTERNATIONAL WORKSHOP ON RADIATION IMAGING DETECTORS

26-30 JUNE 2022

RIVA DEL GARDA, ITALY

# Novel ion imaging concept based on time-of-flight measurements with low gain avalanche detectors

F. Ulrich-Pur,ᵃ,ᵇ,ᶜ T. Bergauer,ᵇ A. Hirtl,ᶜ C. Irmler,ᵇ S. Kaser,ᵇ F. Pittersᵇ and S. Ritᵈ

E-mail: F.Ulrich-Pur@gsi.de

**ABSTRACT:** Treatment planning in ion beam therapy requires accurate knowledge of the relative stopping power (RSP) distribution within the patient. Currently, RSP maps are obtained via conventional x-ray computed tomography (CT) by converting the measured attenuation coefficients of photons into RSP values for ions. Alternatively, to avoid conversion errors that are inherent to this method, ion computed tomography (iCT) can be used since it allows determining the RSP directly. In typical iCT systems, which usually consist of a tracking system and a separate residual energy detector, the RSP is obtained by measuring the particle trajectory and the corresponding water equivalent path length (WEPL) of single ions travelling through the patient. In this work, we explore a novel iCT approach which does not require a residual energy detector. Instead, the WEPL is estimated indirectly by determining the change in time of flight (TOF) due to the energy loss along the ion's path. For this purpose, we have created a Geant4 model of a TOF-iCT system based on low gain avalanche detectors (LGADs), which are fast 4D-tracking detectors that can measure the time of arrival and position of individual particles with high spatial and time precision. To assess the performance of this TOF-iCT concept, we determined the RSP resolution and accuracy for different system settings using the Catphan® CTP404 sensitometry phantom. Within the set of investigated system parameters, the lower limit of the RSP accuracy was found at 0.91%, demonstrating the proof-of-principle of this novel TOF-iCT concept. The main advantage

<figure><img src="image_2.png" /></figure>

<footnote>ᵃGSI Helholtzzentrum für Schwerionenforschung GmbH, Planckstraße 1, 64291 Darmstadt, Germany</footnote>

<footnote>ᵇInstitute of High Energy Physics, Austrian Academy of Sciences, Nikolsdorfer G. 18, 1050 Wien, Austria</footnote>

<footnote>ᶜAtominstitut, TU Wien, Stadionallee 2, 1020 Wien, Austria</footnote>

<footnote>ᵈUniversité de Lyon, CREATIS, CNRS UMR5220, Inserm U1294, INSA-Lyon, Université Lyon 1, Centre Léon Bérard, 69008 Lyon, France</footnote>

<footnote>*Corresponding author.</footnote>
of using this approach is that it could potentially facilitate clinical integration due to its compact
design, which, however, requires experimental verification and an improvement of the current
WEPL calibration procedure.

KEYWORDS: Accelerator applications; Image reconstruction in medical imaging; Instrumentation for hadron therapy

arXiv ePrint: 2209.13676
# Contents

<table><tr><td>1</td><td>Introduction</td><td>1</td></tr><tr><td>2</td><td>Material and methods</td><td>2</td></tr><tr><td>2.1</td><td>Time of flight of ions in matter</td><td>2</td></tr><tr><td>2.2</td><td>Slowing down power</td><td>2</td></tr><tr><td>2.3</td><td>Relation between slowing down power and stopping power</td><td>2</td></tr><tr><td>2.4</td><td>Imaging problem and calibration method</td><td>4</td></tr><tr><td>2.5</td><td>Experimental setup</td><td>5</td></tr><tr><td>3</td><td>Results</td><td>6</td></tr><tr><td>3.1</td><td>RSDP error</td><td>6</td></tr><tr><td>3.2</td><td>Ion computed tomography</td><td>6</td></tr><tr><td>4</td><td>Discussion and conclusion</td><td>7</td></tr></table>

## 1 Introduction

The rising need for fast particle detectors that are able to simultaneously measure the interaction position and time of arrival (ToA) of a particle with high timing (≤100 ps) and spatial (≤100 µm) precision (i.e. 4D-tracking) has led to the development of low gain avalanche detectors (LGADs). LGADs are silicon detectors with an intrinsic amplification operated at gains of 10–20 to optimize the signal-to-noise ratio and therefore the time precision [1]. Even though LGADs have been mainly employed in high energy physics experiments [2], their radiation hardness and excellent 4D-tracking capabilities make them also suitable candidates for ion computed tomography (iCT) [3], which is an imaging modality for the direct determination of the relative stopping power (RSP), i.e. energy loss per unit path length expressed relative to the energy loss in water. In conventional iCT systems [4], the RSP is obtained by measuring the particle trajectory and corresponding water-equivalent path length (WEPL) of individual ions passing through the patient using a tracking system and a separate calorimeter placed downstream of the tracker.

In this study, we present a novel iCT approach, which does not require a residual calorimeter for the WEPL estimation, but instead uses the increase in time-of-flight (TOF) through the patient due to the energy loss along the ion's path. First, we will describe the theoretical background of this method and define a dedicated calibration procedure to map the increase in TOF to the desired WEPL. Finally, we will test this concept by creating an iCT scan of the Catphan® CTP404 sensitometry phantom [5] via a Geant4 [6] simulation of a realistic TOF-iCT scanner based on LGADs. Parts of this study are based on the author's PhD thesis [7].
## 2 Material and methods

### 2.1 Time of flight of ions in matter

The TOF of an ion with mass $m_0$ travelling along a path $\Gamma$ with the velocity $v(\mathbf{x}(s))$ at position $\mathbf{x}(s)$ can be calculated according to

$$ \text{TOF} = \int_{0}^{L} \frac{ds}{v(\mathbf{x}(s))} = \int_{0}^{L} \frac{ds}{c \frac{E_{\text{kin}}(\mathbf{x}(s))}{E_{\text{kin}}(\mathbf{x}(s)) + m_0 c^2} \sqrt{1 + 2 \frac{m_0 c^2}{E_{\text{kin}}(\mathbf{x}(s))}}} \quad (2.1) $$

using $c$ as the speed of light, $L = \int_{\Gamma} ds$ as the total path length and

$$ E_{\text{kin}}(\mathbf{x}(s)) = m_0 c^2 (\gamma(\mathbf{x}(s)) - 1) \quad \text{with} \quad \gamma(\mathbf{x}(s)) = \frac{1}{\sqrt{1 - \frac{v^2(\mathbf{x}(s))}{c^2}}} \quad (2.2) $$

as the relativistic kinetic energy at position $\mathbf{x}(s)$. Eq. (2.1) implies that the TOF strongly depends on the energy loss along the ion's path and therefore the stopping power (SP) of the traversed material. For materials in which significant energy loss takes place, the TOF is increased when compared to the TOF in vacuum, where the energy, thus, the velocity can be assumed to be constant ($v(\mathbf{x}(s)) = v(\mathbf{x}(s=0)) = \text{const}$). In order to understand whether this material-dependent increase in TOF with respect to $\text{TOF}_{\text{vac}}$ can be exploited for the WEPL determination, this effect has to be analysed and quantified.

### 2.2 Slowing down power

Using eq. (2.1) and assuming a small path length $L = \Delta x$, the difference of the real TOF and $\text{TOF}_{\text{vac}}$ can be determined according to

$$ \Delta \text{TOF} = \text{TOF} - \text{TOF}_{\text{vac}} = \int_{0}^{\Delta x} \frac{ds}{v(\mathbf{x}(s))} - \frac{\Delta x}{v(\mathbf{x}(s=0))} \quad (2.3) $$

Similar to the stopping power, one can then define the rate of change per unit path length of $\Delta \text{TOF}$, which, from now on, will be referred to as the “slowing down power” (SDP)

$$ \text{SDP}(E_{\text{kin}}(\mathbf{x}(s))) := \frac{\Delta \text{TOF}}{\Delta x} (E_{\text{kin}}(\mathbf{x}(s))). \quad (2.4) $$

### 2.3 Relation between slowing down power and stopping power

For an ion with a primary kinetic energy $E_0$ travelling a total path length of $L = 2\Delta x$ in small path length increments $\Delta x$ (figure 1), the true TOF can be approximated by

$$ \text{TOF} = \frac{\Delta x}{v(E(x))} + \frac{\Delta x}{v(E(x + \Delta x))} \quad (2.5) $$
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Approximation of the slowing down power. The slowing down power corresponds to the increase in TOF (dashed red curve) per unit path length $\Delta x$ caused by the energy loss along the particle's trajectory.</figcaption></figure>

which, after using a first-order Taylor expansion of eq. (2.5), reads as

$$ \text{TOF} \approx \frac{2\Delta x}{v(E(x))} - \frac{\Delta x^2}{v^2(E(x))} \frac{\partial v(E(x))}{\partial x}. \qquad (2.6) $$

In vacuum, on the other hand, where the energy loss can be neglected ($E(x) = E_0$), eq. (2.5) reduces to

$$ \text{TOF}_{\text{vac}} = \frac{2\Delta x}{v(E(x))}. \qquad (2.7) $$

To estimate the difference between the true TOF and the TOF in vacuum, eq. (2.7) has to be subtracted from eq. (2.6), which results in

$$ \text{TOF} - \text{TOF}_{\text{vac}} \approx - \frac{\Delta x^2}{v^2(E(x))} \frac{\partial v(E(x))}{\partial x}. \qquad (2.8) $$

To obtain the SDP, eq. (2.8) has to be divided by the total path length $L = 2\Delta x$, yielding

$$ \text{SDP}(E(x)) = \frac{\text{TOF} - \text{TOF}_{\text{vac}}}{2\Delta x} \approx - \frac{\Delta x}{2v^2(E(x))} \frac{\partial v(E(x))}{\partial x}. \qquad (2.9) $$

Using the approximation

$$ \frac{\partial v(E(x))}{\partial x} = \frac{\partial v(E(x))}{\partial E(x)} \frac{\partial E(x)}{\partial x} \approx v'(E(x)) \cdot \text{SP}(E(x)) \qquad (2.10) $$

with $\frac{\partial E(x)}{\partial x} \approx \text{SP}(E(x))$ as the SP at energy $E(x)$ and $v'(E(x)) = \frac{\partial v(E(x))}{\partial E(x)}$ as the partial derivative of the velocity $v(E(x))$ with respect to $E(x)$, the SDP (eq. (2.9)) can be rewritten as

$$ \text{SDP}(E(x)) \approx - \frac{\Delta x}{2v^2(E(x))} v'(E(x)) \cdot \text{SP}(E(x)). \qquad (2.11) $$

Similar to the RSP, one can then define the relative slowing down power (RSDP), which is the SDP expressed relative to water

$$ \text{RSDP} := \frac{\text{SDP}_{\text{mat}}(E(x))}{\text{SDP}_{\text{H}_2\text{O}}(E(x))} \approx \frac{-\frac{\Delta x}{2v^2(E(x))}v'(E(x)) \cdot \text{SP}_{\text{mat}}(E(x))}{-\frac{\Delta x}{2v^2(E(x))}v'(E(x)) \cdot \text{SP}_{\text{H}_2\text{O}}(E(x))} = \text{RSP}. \qquad (2.12) $$
Here, the subscripts “mat” and “H₂O” denote the material for which the corresponding SP and SDP were calculated. Eq. (2.12) indicates, that the RSP = SP<sub>mat</sub>(E(x)) / SP<sub>H₂O</sub>(E(x)) can be approximated by the RSDP since v(E(x)) and v'(E(x)) in eq. (2.12) are solely energy-dependent terms, and, therefore, cancel out. To estimate the validity of this approximation, the relative difference between the RSP and RSDP

$$ \epsilon_{\text{RSDP}} = \frac{\text{RSP} - \text{RSDP}}{\text{RSP}} \qquad (2.13) $$

was calculated for different materials and beam energies using the NIST PSTAR database [8].

## 2.4 Imaging problem and calibration method

Using eqs. (2.4) and (2.12), one can, similar to conventional iCT, define a WEPL

$$ \text{WEPL} := \int_{0}^{\text{TOF}-\text{TOF}_{\text{vac}}} \frac{\text{d}\Delta\text{TOF}}{\text{SDP}_{\text{H}_2\text{O}}(\Delta\text{TOF}(E(\mathbf{x}(s))))} = \int_{0}^{L} \text{RSDP}(\mathbf{x}(s))\text{d}s, \qquad (2.14) $$

which can be interpreted as the water equivalent path length resulting in the same total change in TOF as in the target material. Here, the SDP in water was expressed as a function of ΔTOF by mapping the total increase in TOF in water for different water thicknesses and a given initial beam energy E(x(s = 0)) = E₀ to the corresponding residual energy E(x(s)) with ΔTOF(E(x(s = 0))) = 0. This relation between the ΔTOF in water and E(x(s)) can be either obtained by using eq. (2.3) or via a Monte Carlo simulation. After inserting the relation between the RSDP and RSP (eq. (2.12)) into eq. (2.14), it becomes evident that the WEPL as defined in eq. (2.14) is approximately equal to the conventional definition of the WEPL, which is given by

$$ \text{WEPL} = \int_{0}^{L} \text{RSP}(\mathbf{x}(s)) \approx \int_{0}^{L} \text{RSDP}(\mathbf{x}(s)) \text{d}s. \qquad (2.15) $$

Therefore, by estimating the particle path and determining the corresponding total change in TOF (eq. (2.3)), the RSP of the scanned object can be obtained using the standard iCT reconstruction algorithms [9]. However, since actually measuring the total change in TOF would require precise knowledge of the particle’s velocity at each position (eq. (2.3)), we first introduce a simpler approach, which approximates eq. (2.3) by the total change in TOF with respect to the TOF in air (TOF<sub>air</sub>), i.e. the measured TOF in air for the same flight distance, but without the target material (TOF - TOF<sub>vac</sub> ≈ TOF - TOF<sub>air</sub>). Thus, determining ΔTOF reduces to repeating the TOF measurement for a given setup geometry and beam energy after removing the object to be imaged. However, in order to account for the systematic error, which is naturally introduced in this approximation, it is necessary to create a calibration procedure, which maps the measured ΔTOF with respect to TOF<sub>air</sub> to the desired WEPL. Within this study, we have used a fifth-order polynomial as a function of the WEPL to obtain the corresponding calibration curve

$$ \text{TOF} - \text{TOF}_{\text{air}}(\text{WEPL}, E_0) \approx \sum_{i=0}^{5} a_i(E_0) \cdot \text{WEPL}^i \qquad (2.16) $$

with aᵢ(E₀) as the fit parameters for a given primary beam energy E₀.
## 2.5 Experimental setup

To study the feasibility of this “sandwich” TOF-iCT approach, we created a Monte Carlo model of a realistic TOF-iCT scanner based on LGADs using Geant4 version 10.05.1 with the QGSP_BIC_EMY physics list, EM_Opions 3 and a maximum steplength of 1 mm. As shown in figure 2, the scanner consists of four silicon layers with a thickness of 300 µm each, simulating simplified LGAD detectors. The distance between the individual planes (marked as D) and the clearance between the phantom and the innermost planes (marked as C) were fixed and set to 10 cm. To model the intrinsic time and spatial resolution, the true ToA and hit position measured at the entrance on each sensor were distorted using Gaussian distributions centred at zero and with standard deviations σ<sub>T</sub> and σ<sub>xy</sub> corresponding to the time and spatial resolution, respectively. While σ<sub>T</sub> was varied between 0 and 100 ps, σ<sub>xy</sub> was set to a fixed value of 100/√12 µm, corresponding to a typical sensor pitch of 100 µm. In order to obtain the TOF inside the scanner, the mean ToA in the first two LGAD planes was subtracted from the mean ToA inside the last two LGAD planes. Both the calibration and subsequent iCT scan were performed for the previously mentioned time resolutions and three different proton beam energies, namely 170, 200 and 225 MeV.

<figure><img src="image_4.png" /><figcaption>Figure 2. Experimental setup for the calibration (left) and iCT measurement (right).</figcaption></figure>

For the calibration, we measured the increase in TOF with respect to air (eq. (2.16)) for different cubic water phantoms with thicknesses ranging from 0 to 200 mm centred between the two innermost sensors (figure 2(a)). After calculating the calibration parameters according to eq. (2.16), iCT scans of the CTP404 phantom (figure 2(b)) were obtained for all investigated beam energies and time resolutions. As shown in figure 2(b), the CTP404 is a PMMA cylinder with a diameter of 15 cm featuring six cylindrical inserts with a diameter of 12.5 mm each. The materials of the inserts include polymethylpentene (PMP), Teflon, polyoxymethylene (POM, also known as Delrin), Polystyrene, polyethylene (LDPE) and Acrylic. Similar to [3], the RSP images were reconstructed via filtered back projection [9] using an image voxel size of 1 mm³ and 360 projections in 1° steps with 7.5 × 10⁵ primary protons each. Then, the RSP values were measured inside each of the six inserts using square-shaped regions of interests (ROIs) with a side length of 6 mm (indicated as red boxes in figure 2(b)). To estimate the accuracy of the reconstructed RSP, the mean absolute percentage error (MAPE) of the RSP was calculated according to

$$MAPE = \sum_{i=1}^{n=6} \left( \frac{|RSP_{ref,i} - RSP_{simu,i}|}{RSP_{ref,i}} \right) / n, \quad (2.17)$$
with RSP<sub>simu,i</sub> as the mean RSP value of insert i and RSP<sub>ref,i</sub> as the corresponding reference RSP value taken from [3]. As a measure of RSP resolution, the quartile coefficient of dispersion (QCOD) was determined for each insert using

$$QCOD_{RSP} = \frac{Q_{3,RSP} - Q_{1,RSP}}{Q_{3,RSP} + Q_{1,RSP}}, \qquad (2.18)$$

with Q<sub>1,RSP</sub> and Q<sub>3,RSP</sub> corresponding to the first and third quartiles of the obtained RSP distribution.

# 3 Results

## 3.1 RSDP error

Figure 3(a) shows the calculated RSDP for different beam energies and materials using the NIST PSTAR database. No significant energy dependence could be observed. Also, the relative difference between the RSDP and RSP, which is depicted in figure 3(b), is below 1% for beam energies greater than 100 MeV. For lower beam energies, however, $\epsilon_{RSDP}$ ranges between -2% and 3%, which depends on the investigated RSP value.

<figure><img src="image_5.png" /><figcaption>Figure 3. RSDP for different materials and beam energies (left) and RSDP error (right).</figcaption></figure>

## 3.2 Ion computed tomography

The absolute TOF difference with respect to air was measured for different water phantom thicknesses as described in section 2.5. The corresponding calibration curves are shown in figure 4(a) for $\sigma_T = 30$ ps.

After determining the calibration parameters as described in section 2, iCT scans of the CTP404 were created. The corresponding central slices of the iCT images are shown in figure 4(b) for two different system settings. The image noise seems to increase with inferior time resolution. This also becomes apparent when looking at figure 5(b) where the QCOD of the RSP is shown for the Teflon insert (right). As indicated in figure 5(b), to improve the RSP resolution, the beam energy and $\sigma_T$ should be kept as low as possible, which is consistent with the calibration curves
in figure 4(a). Also, the RSP accuracy improves with lower beam energy and shows a systematic dependence on the intrinsic time resolution. Only at the lowest of the investigated beam energies and at $\sigma_T = 50$ ps, the RSP MAPE was 0.91%, which is below the required 1% margin [4].

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Calibration curves for $\sigma_T = 30$ ps (left) and reconstructed central slices of the CTP404 phantom recorded with 175 MeV protons and $\sigma_T = 0$ ps (centre) and $\sigma_T = 30$ ps (right). The pixels outside the field of view (FOV) have been set to 0.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 5.** RSP MAPE (left) and RSP resolution (right) for different intrinsic time resolutions and primary proton beam energies.</figcaption></figure>

## 4 Discussion and conclusion

In this study, we could show that “sandwich” TOF-iCT could, in principle, be used as an alternative imaging modality to conventional iCT [4]. The compact design of this scanner would help to reduce the cost and facilitate the integration into a clinical environment. The performance of this modality, however, was poorer when compared to an iCT systems that uses an LGAD-based TOF calorimeter for the WEPL estimation [3], which could be due to the simplified calibration procedure presented in this work. Therefore, more advanced calibration models should be investigated in the future.
Acknowledgments

This project received funding from the Austrian Research Promotion Agency (FFG), Grant Number
869878.

References

[1] H.F.-W. Sadrozinski et al., 4D tracking with ultra-fast silicon detectors, Rep. Prog. Phys. **81** (2018) 026101.

[2] T. Szumlak, Silicon detectors for the LHC Phase-II upgrade and beyond RD50 status report, Nucl. Instrum. Methods Phys. Res., Sect. A **958** (2020) 162187.

[3] F. Ulrich-Pur et al., Feasibility study of a proton CT system based on 4D-tracking and residual energy determination via time-of-flight, Phys. Med. Biol. **67** (2022) 095005.

[4] R.P. Johnson, Review of medical radiography and tomography with proton beams, Rep. Prog. Phys. **81** (2018) 016701.

[5] The Phantom Laboratory, *Catphan 600® Manual*, https://www.phantomlab.com/catphan-600.

[6] S. Agostinelli et al., Geant4—a simulation toolkit, Nucl. Instrum. Methods Phys. Res., Sect. A **506** (2003) 250.

[7] F. Ulrich-Pur, Advancing ion computed tomography by incorporating time-of-flight and 4D tracking, PhD thesis, TU Wien (2022).

[8] M. Berger et al., Stopping-power and range tables for electrons, protons, and helium ions, NIST Standard Reference Database 124, http://physics.nist.gov/Star.

[9] S. Rit et al., Filtered backprojection proton CT reconstruction along most likely paths, Med. Phys. **40** (2013) 031103.