Pulse pileup model for spectral resolved X-ray photon-counting detectors with dead time and retrigger capability

We developed an analytical model to evaluate the effect of signal pileup on the recorded energy spectrum in X-ray photon-counting detectors affected by dead time and equipped with retrigger capability. The retrigger function allows the system to work in a specific non-paralyzable counting mode by counting the time-over-threshold of piled-up signals in multiples of a predefined and selectable retrigger time. The model, designed for rectangle-like-shaped signals, allows for arbitrary input energy spectra and can significantly help understand and optimize the behavior of counting detectors with spectral capabilities and retrigger mechanisms in applications involving polychromatic beams, e.g., spectral X-ray imaging and computed tomography (CT), in a time-efficient way. Dedicated numerical simulations were used to validate the model under several conditions of incoming flux and threshold energy, with excellent results.

1 Introduction

In the classic configuration of a particle-counting detector, an event is detected as soon as a signal passes a certain threshold level, and for the time-over-threshold period, the system is insensitive to additional signals. At high incoming rates, events undergo pileup, and the resulting signal frequently exceeds the threshold, paralyzing the counting capability. This behavior is well-described for monochromatic signal pulses by the commonly known paralyzable counting mode, where the recorded rate m is related to the incoming rate n—assumed to follow Poisson statistics—with the relation $m=n{e}^{-n{\tau }_P}$, with the dead time τ_P being the single pulse time width [1, 2]. The instant retrigger technology is a circuital expedient introduced by DECTRIS Ltd. first in the PILATUS3 [3] counting ASIC and then adopted in the subsequent counting ASICs IBEX [4] and KITE [5] to overcome the limits of the system paralysis, thus extending the response linearity toward higher values of the incoming rate. In particular, the retrigger mechanism consists of a (time-variant) logic that counts the time-over-threshold of piled-up signals in multiples of a predefined and selectable retrigger time τ_R, making the system work in a particular type of non-paralyzable counting mode. The analytical relation between the incoming and recorded rates for this specific behavior and monochromatic signal pulses was provided in [6] as $m=\frac{n}{e^{-n{\tau }_P}+n{\tau }_R}$, and further statistical insights (spatial and temporal variance) were provided in [7].

An additional and equally fundamental problem of pulse pileup is the distortion brought to the energy spectrum recorded by the detector. This has potential implications for all those applications interested in the spectral content of the detected radiation up to high levels of incoming fluxes, e.g., X-ray imaging and CT. In the last decade, it has been demonstrated that photon-counting systems can not only bring a series of potential improvements with respect to traditional energy-integrating systems in the field of X-ray imaging and computed tomography (CT)—such as lower noise, dose efficiency, and higher spatial resolution—but also allow for novel diagnostic techniques based on their uniquely offered spectral-resolving capabilities [8–13]. Understanding the dynamics of pileup-driven spectral distortion is, therefore, crucial for guiding the design and optimization of future photon-counting systems, both from the hardware and image processing/algorithm point of view. Monte Carlo-based numerical tools are a widespread choice to address this kind of problem as they provide an extremely flexible and accurate framework but are often data-intensive and time-consuming routines. Analytical (or semi-analytical) models can significantly help in this respect, providing a deeper understanding of the underlying physics and a better knowledge of the functional dependence of the model parameters.

Examples of counting detectors operating in the classic non-paralyzable¹ mode and for triangular-like pulse shapes can be found in [14]; the extension to the case of the classic paralyzable mode is discussed in [15, 16], addressing it in an exact way and for arbitrary pulse shapes. Additionally, [17] addresses both counting modes in a semi-analytical manner, in principle extendable to arbitrary pulse shapes. The peculiar behavior introduced by the retrigger mechanism requires, however, a dedicated form of treatment, which has not yet been addressed. In this work, we present the derivation and validation of an analytical model that relates the recorded rate to the incoming rate in counting systems featuring the retrigger capability for arbitrary input energy spectra. The goal is to provide a simple and time-efficient tool to evaluate the spectral performance of this not-yet-fully explored class of detectors up to high levels of the incoming rate. To achieve this, we extended the previously described model in [6] to the case of the polychromatic input energy spectrum. For the needs of the analytical treatment, signals are assumed to be rectangle-like shaped.

The model is validated through a comparison with Monte Carlo simulations, which mimic the realistic case study of a direct X-ray tube beam impinging on a pixelated cadmium telluride (CdTe) sensor read out by a counting ASIC over a wide range of threshold energies and for values of incoming rates relevant for clinical CT applications.

The paper is organized as follows. In Section 2.1, we resume the working principle of the retrigger mechanism; in Section 2.2, we derive the analytical model, in Section 2.3, we introduce the Monte Carlo simulation framework used for the comparison; and in Section 2.4, we describe the used case study. Validation results are presented and discussed in Section 3.1 in terms of the recorded count rate as a function of the incoming rate and in Section 3.2 in terms of the recorded spectra as a function of the incoming rate. Finally, in Section 3.3, we evaluate the impact of a realistic (non-ideal) pulse shape on the recorded spectra.

2 Materials and methods

2.1 Instant retrigger technology

To ease the reading, we provide a brief description of the retrigger mechanism principles, which have already been carefully addressed in [3, 6], with the help of Figure 1. In single X-ray photon counting circuitry, an event is detected and counted as soon as the rising edge of the incoming analog signal (shaperout) exceeds a predefined threshold (vth). Without a retrigger mechanism, during the time-over-threshold, the system is “dead” or insensitive to further incoming signals (pileup), leading to paralyzation at high fluxes (count paralyzed). Conversely, with the retrigger mechanism, as soon as the first analog signal exceeds the threshold, an additional internal digital dead time signal (deadgenout) is generated and periodically evaluated—with a configurable period that we call retrigger time τ_R (t_dead). At every evaluation, if the incoming signal is still high, the counter is increased by a unit; if it has fallen below the threshold, the dead time signal is deactivated. Possible threshold crossings during the τ_R period do not give rise to counts. In this way, a non-paralyzable counting mode is achieved (count retriggered) with a limiting recorded rate of 1/τ_R in the high flux regime. The value of τ_R is customarily set slightly longer than the time-over-threshold of the input signal to prevent the system counting individual pulses multiple times.

FIGURE 1

FIGURE 1. Signal waveforms illustrating the retrigger mechanism, obtained from [3].

2.2 Pulse pileup model

The model is derived based on the following assumptions:

i. The events’ arrival follows Poisson statistics with an average true incoming rate n.

ii. The signals have the shape of a rectangle with a width τ_P.

iii. Signals’ amplitudes are mutually independent stochastic variables and are distributed according to a probability function that corresponds to the energy spectrum s_E(ϵ).

iv. The retrigger time is longer than the pulse width, i.e., τ_R > τ_P ²

As for the monochromatic case described in [6], to the purpose of an analytical treatment of the problem, it is easier to focus on the rate loss n − m. This can be expressed for convenience as the sum of four separate contributions that we call l₁, l₂, l₃, and l₄ such that

$n-m=l_1+l_2+l_3+l_4.(1)$

l₁ and l₂ deal with the rate loss occurring in the “dead time” period after the detection of a signal or after an evaluation of the retrigger signal; l₃ deals with the possibility of an incoming event being undetected in the first place; and l₄ deals with the rate loss occurring during the time periods in between detected events.

In the following list, we describe the derivation of l₁ to l₄:

$\left.l_1\right)$ After a positive detection of an incoming event during the “live period” of the detector or after an evaluation of the retrigger signal, for a time interval of length Δτ = τ_R − τ_P, all incoming events are lost since they would not last until the next evaluation of the retrigger signal. Therefore, on average, nΔτ incoming events are lost for every recorded count m, yielding

$l_1=mn\mathrm{\Delta }\tau .(2)$

$\left.l_2\right)$ During the remaining time interval of length τ_P until the next retrigger evaluation, three possibilities exist. First, if there is no incoming event, nothing is lost. Second, if the overall signal (resulting from the sum of the individual signals) lies below the threshold, all the incoming events are lost. Third, if there is at least one incoming event and the overall signal lies above the threshold, all except one event are lost.

Before making l₂ explicit, it is necessary to introduce the following notation. Given s_E(ϵ) as the probability density function of the amplitude (which is proportional to its energy) of an individual signal, by definition the corresponding cumulative distribution function, S_E(ϵ) is

$S_E\left(ϵ\right)=\mathrm{Pr}(\mathrm{E}\le ϵ)={\int }_0^{ϵ}{s}_E\left(x\right)\mathrm{d}x.(3)$

When i signals pileup, the stochastic variable $E_i={\sum }_{k=1}^i{E}_k$ is distributed according to the convolution of the individual probability density function repeated i times, which is indicated as

$s_{E_i}\left(ϵ\right)=\underset{i}{\underbrace{s_E\left(ϵ\right)\ast s_E\left(ϵ\right)\ast \dots \ast s_E\left(ϵ\right)}}\stackrel{\text{def}}{=}{s}_E{\left(ϵ\right)}^{*i},(4)$

from which the corresponding cumulative distribution function can be deduced as

$S_{E_i}\left(ϵ\right)=\mathrm{Pr}({\mathrm{E}}_{\mathrm{i}}\le \mathrm{ϵ})={\int }_0^{ϵ}{s}_E{\left(x\right)}^{*i}\mathrm{d}x.(5)$

Let us also recall the probability of having i incoming events following the Poisson distribution in the time interval τ_P.

$\mathrm{Pr}\left(i|{\tau }_P\right)=\frac{{\left(n{\tau }_P\right)}^i{e}^{-n{\tau }_P}}{i!}.(6)$

We can now write l₂ as

$l_2=m\left\{\sum _{i=1}^{+\infty }i\mathrm{Pr}\left(i|{\tau }_P\right)S_{E_i}\left(E_{th}\right)+\sum _{i=1}^{+\infty }\left(i-1\right)\mathrm{Pr}\left(i|{\tau }_P\right)\left(1-S_{E_i}\left(E_{th}\right)\right)\right\},(7)$

where E_th is the threshold energy of the counting circuitry. Rearranging the terms, we can write that

$\begin{array}{ccc}\hfill l_2& =\hfill & m\left\{\sum _{i=1}^{+\infty }i\mathrm{Pr}\left(i|{\tau }_P\right)S_{E_i}\left(E_{th}\right)+\sum _{i=1}^{+\infty }i\mathrm{Pr}\left(i|{\tau }_P\right)\left(1-S_{E_i}\left(E_{th}\right)\right)-\sum _{i=1}^{+\infty }\mathrm{Pr}\left(i|{\tau }_P\right)\left(1-S_{E_i}\left(E_{th}\right)\right)\right\}\hfill \\ \hfill & =\hfill & m\left\{\sum _{i=1}^{+\infty }i\mathrm{Pr}\left(i|{\tau }_P\right)-\sum _{i=1}^{+\infty }\mathrm{Pr}\left(i|{\tau }_P\right)\left(1-S_{E_i}\left(E_{th}\right)\right)\right\}.\hfill \\ \hfill & \hfill & \hfill \end{array}(8)$

Considering that ${\sum }_{i=1}^{+\infty }i\mathrm{Pr}(i|{\tau }_P)={\sum }_{i=0}^{+\infty }i\mathrm{Pr}(i|{\tau }_P)=n{\tau }_P$ and defining $\mathcal{A}\stackrel{\text{def}}{=}{\sum }_{i=1}^{+\infty }\mathrm{Pr}(i|{\tau }_P)S_{E_i}(E_{th})$, we can rewrite l₂ as follows:

$l_2=m\left\{n{\tau }_P+1-\mathrm{Pr}\left(0|{\tau }_P\right)-\mathcal{A}\right\},(9)$

and since $\mathrm{Pr}(0|{\tau }_P)=e^{-n{\tau }_P}$, we obtain

$l_2=mn{\tau }_P-m+m{e}^{-n{\tau }_P}+m\mathcal{A}.(10)$

$\left.l_3\right)$ An additional term has to be added to account for the possibility that an incoming event actually passes undetected; in the case when the time interval τ_P precedes each recorded count m, the overall signal lies below the threshold energy.

$\begin{array}{ccc}\hfill l_3& =\hfill & m\sum _{i=1}^{+\infty }\mathrm{Pr}\left(i|{\tau }_P\right)S_{E_i}\left(E_{th}\right)\hfill \\ \hfill & =\hfill & m\mathcal{A}.\hfill \end{array}(11)$

$\left.l_4\right)$ During the “live time” of the detector, all those piled-up signals whose overall signal lies below the threshold are also missed. Given T as the total acquisition time, the detector live time is T(1 − mτ_R). In this time interval, the number of lost events can be approximated as ${\sum }_{i=1}^{+\infty }i\mathrm{Pr}(i|{\tau }_P)S_{E_i}(E_{th})$ for each sub-interval of the live time with size τ_P. This leads to

$l_4=\frac{1}{T}\frac{T\left(1-m{\tau }_R\right)}{{\tau }_P}\sum _{i=1}^{+\infty }i\mathrm{Pr}\left(i|{\tau }_P\right)S_{E_i}\left(E_{th}\right),(12)$

with the term 1/T arising for normalization purposes. Defining $\mathcal{B}\stackrel{\text{def}}{=}{\sum }_{i=1}^{+\infty }i\mathrm{Pr}(i|{\tau }_P)S_{E_i}(E_{th})$, we can rewrite l₄ as follows:

$l_4=\frac{1-m{\tau }_R}{{\tau }_P}\mathcal{B}.(13)$

We can finally substitute the expressions of l₁ to l₄ of Eqs. 2, 10, 11, and 13, respectively, into Eq. 1, and by solving it for m, we obtain the final relation.

$m=\frac{n-\frac{\mathcal{B}}{{\tau }_P}}{n{\tau }_R+e^{-n{\tau }_P}+2\mathcal{A}-\frac{{\tau }_R}{{\tau }_P}\mathcal{B}}.(14)$

A few remarks were observed. Since the energy spectrum s_E(ϵ) is a free parameter of the model that usually lacks an analytical formulation, the terms $\mathcal{A}$ and $\mathcal{B}$ need to be evaluated numerically (typically with negligible computational burden). In addition, we left implicit in the notation the main functional dependencies of m = m(n, Eth). In the following, we will refer to the quantity m(n; Eth) (the semicolon specifies the fixed parameter) as the count rate curve, m(Eth; n) as the recorded integral energy spectrum, and $\frac{\partial m(E_{th};n)}{\partial E_{th}}$ as the recorded differential energy spectrum.

2.3 Monte Carlo validation framework

To validate the model, we used the numerical results of the custom-developed Monte Carlo model already used in [6, 7], consisting of the following. First, a temporal axis is created, spanning a statistically relevant length of several tens of ms and with a sufficiently small time step, typically on the order of ns. a number pulse-like events, computed according to the chosen value of the incoming rate, are distributed across the temporal axis, following a Poisson distribution. This resulting signal, which corresponds to the input of the detector system, is then convolved with the pulse shape of the analog front-end electronics, which is, in principle, arbitrary but is assumed to be a rectangle function in this context. To mimic the behavior of the counting electronics, a count is recorded for every rising edge crossing a specified threshold. If the retrigger capability is enabled, a count is also added for every time interval of length τ_R for which the signal lies above the threshold. The total recorded counts are then normalized to the overall temporal axis length to obtain the recorded rate. The entire operation is repeated for several values of the incoming rate and as a function of the threshold energy. The simulation error is estimated by subdividing the total time span into a set of smaller intervals of equal size, computing the variance over their ensemble, and scaling it to the overall interval time (see Appendix A).

The amplitudes of the incoming signals are randomly and independently extracted from a probability distribution function that corresponds to the specified energy spectrum. The energy spectrum was obtained as the convolution of the direct spectrum of an X-ray tube, computed with the model described in [18], with the response function of a pixelated counting detector, estimated based on the numerical results from FLUKA³—a Monte Carlo particle transport and interaction suite [19, 20]—and a custom-developed pixelization and counting routine. Specifically, the FLUKA code computes, for a large number of incoming X-ray photons, the individual energy depositions in the sensor volume at the microscopic level, including all secondary effects like fluorescence and escape. The pixelization routine translates deposited energy into electric charge, which is drift-diffused toward the collecting electrodes. Signals collected at the pixel side (comprehensive of charge sharing effects) are added to a Gaussian random signal representing the electronic noise and subsequently histogrammed to obtain the sought energy spectrum.

2.4 Case study

The predictions of the analytical model were compared with the results of Monte Carlo simulations for a realistic case study consisting of a direct beam from an X-ray tube with a W-anode, an acceleration voltage of 120 kV_p, and 6.8 mm Al filtering, detected using a CdTe sensor with a pixel size of 330 μm and a thickness of 1500 μm, and read out by a counting ASIC. The spectrum and the threshold energy were sampled with a step of 1 keV. The front-end electronics was assumed to introduce a noise of 1 keV rms. The pulse time was τ_P = 80 ns, and the retrigger time was τ_R = 100 ns. Incoming fluxes n were investigated up to 10⁸ ph/s/pix, which corresponds to roughly 10⁹ ph/s/mm² for the used pixel size, a value that can be considered the upper boundary in clinical CT applications. Figure 2 shows the energy spectrum as detected by the pixelated detector and used as input for both the model and Monte Carlo simulations. The spectrum of the original direct X-ray tube beam is also shown as a reference.

FIGURE 2

FIGURE 2. Energy spectrum detected using a CdTe sensor with a pixel size of 330 μm and a thickness of 1500 μm, used as input for both the model and the Monte Carlo simulations. The direct X-ray tube beam spectrum (W-anode, an acceleration voltage of 120 kV_p, and 6.8 mm Al filtering) is shown as a reference, and the main fluorescence peaks originated from the W-anode of the X-ray tube and the CdTe sensor.

3 Results and discussion

3.1 Count rate curves

Figure 3 shows the comparison between the count rate curves predicted by the model and simulated with the Monte Carlo tool for incoming rates up to 10⁸ ph/s/pix and several values of threshold energies ranging from 5 keV to 140 keV. The rate curve corresponding to the total counts—intended as the ideal (no pileup and no electronic noise) rate obtained with a threshold energy of 0 keV—is also shown as a reference. For increasing threshold energies, the deviation of the recorded curves from one of the total counts, even at low incoming rates, becomes larger as the obvious consequence of losing an increasing fraction of low-energy events. On the other hand, all the recorded curves saturate for increasing incoming rates to the limit 1/τ_R, set by the retrigger mechanism. As a quantitative measure of the agreement between predicted and simulated results, we chose the L² relative error norm (L²REN), which is defined as follows:

$L^{2}REN=\frac{{\sum }_{i=1}^{\mathit{NoP}}{\left(x_i^{\mathit{sim}}-x_i^{\mathit{ana}}\right)}^2}{{\sum }_{i=1}^{\mathit{NoP}}{\left(x_i^{\mathit{sim}}\right)}^2},(15)$

where NoP is the number of analyzed points. The values are reported in Table 1. An excellent match on the sub-percent level is obtained up to the last threshold energy (90 keV in this case), whose value lies below the highest energy actually contained in the input spectrum, i.e., 120 keV. For threshold energies higher than that, only pure piled-up events are detected, and the increase in L²REN seems to suggest that the accuracy of the approximations used in the derivation of the model (in particular, the term l₄) begins, in these conditions, to fade.

FIGURE 3

FIGURE 3. Count rate curves predicted by the model and simulated with the Monte Carlo tool for different values of threshold energy, including the ideal (no pileup and no electronic noise) rate at the threshold energy of 0 keV. Error bars of the simulated data are smaller than the used graphical symbol.

TABLE 1

TABLE 1. L²REN of the count rate curves.

3.2 Spectral response

Figure 4 (top) shows the comparison between the recorded integral energy spectra predicted by the model and simulated with the Monte Carlo tool as a function of the threshold energy and for increasing values of incoming rates ranging from n = 10⁵ ph/s/pix to n = 10⁸ ph/s/pix. Figure 4 (bottom) shows the corresponding deviations expressed in percentage. Looking at the integral spectra, it is possible to observe how, for increasing incoming rates, events of energy higher than the highest one contained in the input spectrum, i.e., 120 keV, and, therefore due to pure pileup effects, are detected with increasing probability. At the same time, the integral spectra tend to the limit set by the retrigger mechanism 1/τ_R, in a monotonic (i.e., non-paralyzable) way at all values of threshold energy. The relative deviation is in the order of a few percent up to the notable value of 120 keV (between −8% and +4%), a result that we can consider excellent considering the several orders of magnitudes spanned by the spectra across this threshold energy range. For threshold energies higher than this limit and high input rates (above 2⋅10⁷ ph/s/pix), the deviation slightly increases. As previously mentioned, this could be due to the approximations used in the derivation of the model starting to lose their validity. It is nevertheless worth noting that, in a realistic X-ray imaging/CT scenario, the threshold energy is not set above a certain fraction of the X-ray tube acceleration voltage. Therefore, a discrepancy in this extreme energy range is not of practical concern. For this reason, we limited the computation of the L²REN figure of merit to the threshold energy interval of 1 keV–120 keV. The values, reported in Table 2, show the agreement that lies at the (sub-)percent level.

FIGURE 4

FIGURE 4. Recorded integral energy spectra predicted by the model (solid lines) and simulated with the Monte Carlo tool (dotted lines plus error bars) for increasing values of the incoming rate (top). Corresponding relative deviation in percentage (bottom).

TABLE 2

TABLE 2. L²REN of the recorded integral spectra in the energy range of 1–120 keV.

Figure 5 shows the comparison between the recorded differential energy spectra obtained by numerical derivation of the integral spectra of Figure 4 (top) with respect to the threshold energy and the corresponding relative deviations. The L²REN metrics, computed in the threshold energy interval of 1–120 keV, are reported in Table 3. Although higher values of deviation as compared to the recorded integral spectra have to be expected due to the “amplification” caused by the derivative operation, they nevertheless lie below 10% up to incoming rates of 2⋅10⁷ ph/s/pix and at most 20% for higher incoming rates, which are anyway rather extreme as they basically lead to an almost complete count rate saturation.

FIGURE 5

FIGURE 5. Recorded differential energy spectra obtained by deriving the curves of Figure 4 (top) with respect to the threshold energy. Corresponding relative deviation in percentage (bottom).

TABLE 3

TABLE 3. L²REN of the recorded spectra in the energy range of 1–120 keV.

Based on these results, we can conclude that the predictions of the presented analytical model faithfully represent the expectations obtained by numerical tools.

3.3 Impact of non-ideal pulse shape

Real-world analog signals have shapes that significantly differ from the ideal rectangle. Apart from few other special cases (e.g., pulses with a triangular shape in [14] or with a decreasing exponential shape in [16]), a general analytical solution for arbitrary shapes remains, at the moment, more complex, if not out of sight.

To evaluate the impact of this non-ideality on the recorded spectra, we simulated a case study intentionally violating assumption ii. In particular, we used a more realistic pulse shape consisting of an asymmetric Gaussian-like function, as shown in Figure 6. The curve was created by using, as a rising edge, a Gaussian function with a standard deviation of 36 ns and, as a falling edge, a Gaussian function with a standard deviation of 134 ns. The overall FWHM is 80 ns, matching one of the ideal rectangles. The retrigger time was increased to ${\tau }_R^{\prime }$ = 150 ns in order to reduce spurious multiple counts at low threshold energies.

FIGURE 6

FIGURE 6. Asymmetric Gaussian-like function used as a realistic pulse shape. The rising edge consists of a Gaussian function with a standard deviation of 36 ns and a falling edge but with a standard deviation of 134 ns. The ideal rectangle function is shown for comparison.

Figure 7 and Figure 8 show, respectively, the resulting integral and differential energy spectra (top) and corresponding relative deviations (bottom). Analytical data were recomputed with the new value of retrigger time. What emerges from the comparison is that while the overall trends are respected, significant deviations arise for values of incoming flux n ≥ 10⁷ ph/s/pix and threshold energies E_th < 10 keV. This is due to the fact that the model does not take into account the possibility that individual signals can have a time-over-threshold longer than the chosen retrigger time due to the longer pulse tail. The values of L²REN computed in the energy range of 10–120 keV for the recorded integral spectra vary in the range of 0.6%–6.4% (see Table 4), while for the recorded differential spectra, they monotonically increase from 8% to 106% with increasing incoming fluxes (see Table 5).

FIGURE 7

FIGURE 7. Recorded integral energy spectra predicted by the model (solid lines) and simulated with the Monte Carlo tool (dotted lines plus error bars) using the asymmetric Gaussian-like pulse shape of Figure 6 (top). The retrigger time was increased to ${\tau }_R^{\prime }$ = 150 ns. Corresponding relative deviation in percentage (bottom).

FIGURE 8

FIGURE 8. Recorded differential energy spectra obtained by deriving the curves of Figure 7 (top) with respect to the threshold energy. Corresponding relative deviation in percentage (bottom).

TABLE 4

TABLE 4. L²REN of the recorded integral spectra in the energy range of 1–120 keV using the asymmetric Gaussian-like pulse shape.

TABLE 5

TABLE 5. L²REN of the recorded spectra in the energy range of 1–120 keV using the asymmetric Gaussian-like pulse shape.

It is also legit to question whether the general discrepancy could be “compensated” by the use of an effective pulse time τ_Peff, as described in [6]. Unfortunately, as shown in the reference, τ_Peff would be a function of the threshold energy (to be empirically found, moreover) and, therefore, not compatible with the model presented in this work.

4 Conclusion

We presented the derivation of an analytical model to evaluate the effect of pulse pileup on the recorded energy spectrum measured by X-ray photon-counting detectors affected by dead time and featuring the retrigger capability—a circuital mechanism that allows the detector to work in a specific non-paralyzable counting mode. The model is valid for arbitrary incoming energy spectra and derived for the case of a rectangle-like pulse shape.

We compared the analytical predictions with the results of Monte Carlo simulations obtained with a custom-developed numerical code for the realistic clinical CT case study of a direct X-ray tube beam with an acceleration voltage of 120 kV_p, impinging on a CdTe sensor with a pixel size of 330 μm and a thickness of 1500 μm, and read out by a counting ASIC. The investigated quantities (count rate curves and recorded energy spectra) showed an excellent level of agreement, quantified in terms of L²REN, over a wide range of incoming rates relevant for X-ray CT applications, i.e., up to 10⁸ ph/s/pix, which corresponds to approximately 10⁹ ph/s/mm² with our pixel size, and over a wide range of threshold energies. Slight deviations, possibly ascribable to approximations used in the model derivation, were progressively found for the concurrence of increasing threshold energies above the highest energy contained in the incoming spectrum and mid-high incoming rates ($>5\cdot$10⁷ ph/s/pix). However, considering that in realistic scenarios, threshold energies are commonly kept below that limit, this does not seem like a practical concern. For completeness, we also investigated the impact of a non-ideal pulse shape on the recorded spectra, thus intentionally violating one of the assumptions on the basis of the model. We found that despite the overall trend being respected, significant relative deviations can arise, and this constitutes one of the main limitations to the general applicability of the proposed method.

Despite the simplification of rectangle-like shaped signals needed for the analytical treatment of the problem, having an analytical expression for the recorded count rate for such a class of detectors can significantly help better understand the yet unexplored impact of the retrigger mechanism on the spectral response and optimize the system design and parameters, all in a time-efficient way, with potential benefits in particular in the field of spectral X-ray imaging and CT applications allowed by modern multi-threshold photon-counting systems.

Data availability statement

The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

Author contributions

PZ: conception, theory, data analysis, and writing. CA: conception, Monte Carlo simulations, and proof-reading. All authors contributed to the article and approved the submitted version.

Conflict of interest

Authors PZ and CA were employed by DECTRIS Ltd.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphy.2023.1205638/full#supplementary-material

Footnotes

¹The classic non-paralyzable counting mode exhibits a relation $m=\frac{n}{1+n\tau }$ [2].

²This condition prevents counting single events multiple times, which is usually an undesired effect in a retriggered system, unless it is intentionally operated in a time-over-threshold mode, as shown in [21].

³v. 4–2.1. The physics was set to multiple Coulomb scattering with a cutoff energy of 1 keV for electrons and 100 eV for photons. Fluorescence was enabled, and no biasing was used.

References

1. Kraft P, Bergamaschi A, Broennimann C, Dinapoli R, Eikenberry EF, Henrich B, et al. Performance of single-photon-counting PILATUS detector modules. J Synchrotron Radiat (2009) 16:368–75. doi:10.1107/s0909049509009911

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Knoll GF. Radiation detection and measurement. 3rd ed. New Jersey United States: John Wiley and Sons, Inc. (2000). 119–22. Chapter 3.

Google Scholar

3. Loeliger T, Broennimann C, Donath T, Schneebeli M, Schnyder R, Trueb P. The new PILATUS3 ASIC with instant retrigger capability. In: IEEE Nuclear Science Symposium and Medical Imaging Conference Record (NSS/MIC) N6-2 (2012). p. 610–5. doi:10.1109/NSSMIC.2012.6551180

CrossRef Full Text | Google Scholar

4. Bochenek M, Bottinelli S, Broennimann C, Livi P, Loeliger T, Radicci V, et al. Ibex: Versatile readout asic with spectral imaging capability and high count rate capability. IEEE Trans Nucl Sci (2018) 65:1285–91. doi:10.1109/TNS.2018.2832464

CrossRef Full Text | Google Scholar

5. Zambon P, Bottinelli S, Schnyder R, Musarra D, Boye D, Dudina A, et al. Kite: High frame rate, high count rate pixelated electron counting ASIC for 4D STEM applications featuring high-Z sensor. Nucl Inst. Methods Phys Res A (2023) 1048(2023):167888. doi:10.1016/j.nima.2022.167888

CrossRef Full Text | Google Scholar

6. Zambon P. Dead time model for X-ray photon counting detectors with retrigger capability. Nucl Inst. Methods Phys Res A (2021) 994:165087. doi:10.1016/j.nima.2021.165087

CrossRef Full Text | Google Scholar

7. Zambon P. Temporal and spatial variance in pixelated counting detectors with dead time and retrigger capability. Nucl Inst. Methods Phys Res A (2022) 1032:166625. doi:10.1016/j.nima.2022.166625

CrossRef Full Text | Google Scholar

8. Taguchi K, Iwanczyk JS. Vision 20/20: Single photon counting x-ray detectors in medical imaging. Med Phys (2013) 40:100901. doi:10.1118/1.4820371

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Symons R, Krauss B, Sahbaee P, Cork TE, Lakshmanan MN, Bluemke DA, et al. Photon-counting CT for simultaneous imaging of multiple contrast agents in the abdomen: An in vivo study. Med Phys (2017) 44(10):5120–7. doi:10.1002/mp.12301

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Willemink MJ, Persson M, Pourmorteza A, Pelc NJ, Fleischmann D. Photon-counting CT: Technical principles and clinical prospects. Radiology (2018) 289:293–312. doi:10.1148/radiol.2018172656

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Mannil M, Hickethier T, von Spiczak J, Baer M, Henning A, Hertel M, et al. Photon-counting CT high-resolution imaging of coronary stents. Invest Radiol (2018) 53:143–9. doi:10.1097/RLI.0000000000000420

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Flohr T, Petersilka M, Henning A, Ulzheimer S, Ferda J, Schmidt B. Photon-counting CT review. Physica Med (2020) 79:126–36. doi:10.1016/j.ejmp.2020.10.030

CrossRef Full Text | Google Scholar

13. Lell MM, Kachelrieß M. Recent and upcoming technological developments in computed tomography. Invest Radiol (2020) 55:8–19. doi:10.1097/RLI.0000000000000601

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Taguchi K, Frey EC, Wang X, Iwanczyk JS, Barber WC. An analytical model of the effects of pulse pileup on the energy spectrum recorded by energy resolved photon counting x-ray detectors. Med Phys (2010) 37:3957–69. doi:10.1118/1.3429056

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Taguchi K, Zhang M, Frey EC, Wang X, Iwanczyk JS, Nygard E, et al. Modeling the performance of a photon counting x-ray detector for CT: Energy response and pulse pileup effects. Med Phys (2011) 38:1089–102. doi:10.1118/1.3539602

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Roessl E, Daerr H, Proksa R. A Fourier approach to pulse pile-up in photon-counting x-ray detectors. Med Phys (2016) 43(3):1295–8. doi:10.1118/1.4941743

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Cammin J, Kappler S, Weidinger T, Taguchi K. Evaluation of models of spectral distortions in photon-counting detectors for computed tomography. J Med Imag (2016) 3(2):023503. doi:10.1117/1.JMI.3.2.023503

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Tucker DM, Barnes GT, Chakraborty DP. Semiempirical model for generating tungsten target x-ray spectra. Med Phys (1991) 18:211–8. doi:10.1118/1.596709

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Ferrari A, Sala PR, Fasso A, Ranft J. Fluka: A multi-particle transport code (2005). Program version 2005 [INFN-TC-05-11].

Google Scholar

20. Böhlen T, Cerutti F, Chin M, Fassò A, Ferrari A, Ortega P, et al. The FLUKA code: Developments and challenges for high energy and medical applications. Nucl Data Sheets (2014) 120:211–4. doi:10.1016/j.nds.2014.07.049

CrossRef Full Text | Google Scholar

21. Bergamaschi A, Dinapoli R, Greiffenberg D, Henrich B, Johnson I, Mozzanica A, et al. Time-over-threshold readout to enhance the high flux capabilities of single-photon-counting detectors. J Synchrotron Radiat (2011) 18:923–9. doi:10.1107/S0909049511034480

PubMed Abstract | CrossRef Full Text | Google Scholar

Appendix A: Computation of the simulation error

Let T be the total simulated time interval, M be the total number of recorded events in this time interval, and m be the recorded rate, then we obtain m = M/T. By dividing the total interval T into N sub-intervals of equal size (large enough to maintain a relevant statistics), we obtain M_i independent recorded counts for each sub-interval such that

$M=\sum _{i=1}^N{M}_i.(\mathrm{A1})$

Given that M_i shares the same probability space, the variance of m can be inferred using the basic properties of random variables as follows:

$\mathrm{Var}\left[m\right]=\mathrm{Var}\left[\frac{{\sum }_{i=1}^N{M}_i}{T}\right]=\frac{1}{T^2}\mathrm{Var}\left[M_i\right]N.(\mathrm{A2})$