Influence of the background in Compton camera images for proton therapy treatment monitoring

Objective. Background events are one of the most relevant contributions to image degradation in Compton camera imaging for hadron therapy treatment monitoring. A study of the background and its contribution to image degradation is important to define future strategies to reduce the background in the system. Approach. In this simulation study, the percentage of different kinds of events and their contribution to the reconstructed image in a two-layer Compton camera have been evaluated. To this end, GATE v8.2 simulations of a proton beam impinging on a PMMA phantom have been carried out, for different proton beam energies and at different beam intensities. Main results. For a simulated Compton camera made of Lanthanum (III) Bromide monolithic crystals, coincidences caused by neutrons arriving from the phantom are the most common type of background produced by secondary radiations in the Compton camera, causing between 13% and 33% of the detected coincidences, depending on the beam energy. Results also show that random coincidences are a significant cause of image degradation at high beam intensities, and their influence in the reconstructed images is studied for values of the time coincidence windows from 500 ps to 100 ns. Significance. Results indicate the timing capabilities required to retrieve the fall-off position with good precision. Still, the noise observed in the image when no randoms are considered make us consider further background rejection methods.


Introduction
Proton therapy is a technique that consists on the irradiation of tumors by means of proton beams. It has been deeply studied since its therapeutic interests were first described in 1946(Wilson 1946. The dose deposition by protons as they pass through the patient is highly localized near the end of the range, and consequently the damage to the surrounding tissue is minimized. Even so, many factors have an impact on the range uncertainty, which may put healthy tissue adjacent to the tumor at risk or result in tumor underdose. In order to avoid this, safety margins are set during treatment planning (Paganetti 2012). The aim of reducing these safety margins motivates research on online monitoring systems, and among them, the use of prompt gamma rays (PG) emitted during nuclear de-excitation processes (Krimmer et al 2018, Wrońska andDauvergne 2021). Compton cameras, which were first proposed for medical application in Todd et al (1974), have been studied, for proton beam range verification, employing different types of detectors . Double-sided strip silicon detectors are chosen in some Compton camera layouts dedicated to proton beam monitoring, along with scintillator crystals for the absorber layers, as the designs described in Roellinghoff et al (2011), Seo et al (2011), Krimmer et al (2015), Aldawood et al (2017) which employ LYSO, NaI, BGO and LaBr 3 crystals, respectively. Other designs make use entirely of CZT detectors (McCleskey et al 2015), LaBr 3 detectors (Llosá et al 2016), Ce:GAGG scintillators (Koide et al 2018), LYSO fibers (Kasper et al 2020) or combine CZT with BGO absorbers (Golnik et al 2016).
In a two-layer Compton camera, events suitable for Compton image reconstruction are those in which a photon undergoes a Compton scatter in the first plane and photoelectric effect, Compton scattering or pair production in the second plane. The remaining events are part of the background. Various studies have 2. Materials and methods 2.1. Simulation description Monte Carlo simulations have been carried out with GATE version 8.2, a dedicated toolkit for simulations in medical imaging and radiotherapy based on Geant4 open-source code (Agostinelli et al 2003, Jan et al 2004, Sarrut et al 2014. The simulated Compton camera has two planes with areas larger than the experimental prototype (with crystals of size 25.8 × 25.8 × 5 mm 3 ) contemplating a future prototype with enhanced efficiency, so that both planes are made up of 100 × 100 × 5 mm 3 parallel scintillator crystals of LaBr 3 , separated by 50 mm (center to center).
A PMMA phantom (ρ = 1.19 g cm −3 ) of 300 mm height and 30 mm radius is irradiated by a monoenergetic proton beam as a primary source for three proton beam energies. Primary protons of 90, 150 and 230 MeV at different beam intensities have been delivered into the PMMA phantom with total amount of 5 × 10 9 in each case, a delivered number of protons similar to Hueso-González et al (2015), Draeger et al (2018). No energy spread of the beam and an ideally small beam size (1 mm standard deviation in the transverse direction) has been employed. The first plane of the Compton camera has been set, in the z-axis, at 100 mm from the center of the phantom (figure 1) to improve the detection efficiency in the simulation. In addition to the photons arriving to the detectors, the Compton camera may also detect secondary protons, neutrons, and other particles generated in nuclear reactions within the phantom and the detector layers. The treatment room is not simulated in the study. The physical materials that actually would surround the simulated system add a background component -known as room component-due to the interaction of particles in the beam nozzle and surrounding materials (walls, ceiling, gantry) (Smeets et al 2012, Pinto et al 2014. The 90, 150 and 230 MeV proton beam energies have projected ranges of 56, 136 and 284 mm in PMMA, respectively (Berger et al 2005). In each case, the Compton camera has been shifted along the x-axis accordingly (i.e. parallel to the beam direction) to center the fall-off region in the crystal of the Compton camera.
Finally, to study the influence of different phantom sizes and the distance of the Compton camera to the beam axis, simulations with phantoms of 100 and 200 mm radius have been performed, for which the Compton camera is placed further away in each case. We have chosen to maintain the distance between the first plane and the phantom surface (70 mm) and the separation between planes (50 mm), in such a manner that in the two new scenarios, the first plane of the Compton camera is placed at 170 mm and 270 mm from the beam axis, respectively.

Physics list selection
Proton simulations entail hadronic processes and therefore require complex models. Extensive studies have been carried out evaluating and comparing the performance of GEANT4 physics lists in different versions . In our case (GEANT4 version 10.5.1), reference physics list have been implemented and compared in the proton simulations: QGSP_BIC_HP_EMZ, QGSP_BIC_AllHP_EMZ and QGSP_BERT_HP_EMZ. The physics list QGSP_BIC_HP_EMZ includes the binary cascade model and a high precision neutron package for neutrons of energies below 20 MeV. It models the discrete production of PG from interactions of low energy protons with nuclei (Verburg et al 2012), although some studies have reported that PG emission yield is overestimated by 40% (Schumann et al 2015, Pinto et al 2016. To simulate the electromagnetic physics, the selected list is the GEANT4 pre-defined physics constructor G4EmStandardPhysics_option4known as emstandard_opt4, here renamed to EMZ-which is the most accurate standard electromagnetic . For that reason, and the fact that none of the physics lists reproduce the photon emission accurately, we have based the selection of the physics, not in the literature but in our experimental results. The chosen physics list is the one matching as close as possible to the measured data. In order to choose the physics list, the spectrum obtained in a simulated LaBr 3 crystal has been compared to the experimental energy spectrum measured with one of the system detector planes of MACACO (figure 3(a)). The spectrum was measured in an experimental beam test in the AGOR cyclotron at KVI-Center for Advanced Radiation Technology (University of Groningen), where a 150 MeV proton beam was delivered into a PMMA phantom located in front of MACACO . The physics list selected to carry out the simulations, according to the results shown in section 3.1, is QGSP_BIC_HP_EMZ.

Beam time structure
In order to simulate a beam structure, the time of emission of each event in the simulation is modified to distribute the launched primary protons in pulses of particles. The chosen model is based on the IBA Isochronous cyclotron C230, in which the latest experimental tests with MACACO III have been carried out, and the method presented by Botas (2013). The beam consists of a set of spots, the duration of which can vary widely depending on the energy and dose administered. Each spot is further divided into pulses (also known as bunches). In the IBA Isochronous cyclotron, pulses are separated by 10 ns and have a width of 2 ns FWHM . Equation (1), the parameters of which are listed in table 1, shows the function employed to modify the time of emission of the proton, t prot , in the ith event. The first two terms detail the beam micro-structure. The first term distributes the event in pulses according the integer value calculated in the division of the ith event number and the number of particles in the pulse, which is straightly related to the beam intensity and randomly extracted by a Poisson distribution (Fontana et al 2019). The Box-Muller transform in the second term generates a normal distribution with a source of uniformly distributed random numbers, giving a Gaussian shape to the pulse. The third term, related to the beam macro-structure, distributes the event in spots

Data processing and background classification
Every interaction taking place in the Compton camera in the course of the simulation-referred to as hits in GATE terminology-is recorded into the simulation output. The hits provide the exact information of the particle type, position, energy deposition and process suffered in each interaction. For a given event, the hits registered in a detector volume are combined into a single, which is intended to reproduce a pulse of an experimental device. In our monolithic crystal, the energy assigned to that single is the sum of the deposited energy of every hit of the event that takes place inside a volume. Likewise, the position of the single is the combined position weighted with the energy deposited by each hit. The event is entirely simulated from beginning to end and timestamps registered in the detectors are relative to the generation of the primary particle in the event. Relative timestamps are transformed to absolute time by summing constant t prot value in the given event, described in section 2.3. At the detector level, the effect of the dark counts has been modeled by generating artificial low energy singles, distributed around 40 keV, with a frequency of 400 kcps, which corresponds to the frequency reported by Hamamatsu Photonics for the MPPC array S13360-3025CS, currently used in MACACO III prototype. A plain dead time model has also been implemented. It consists in the rejection of singles when they are detected in a plane that has previously detected another single, in a time interval smaller than 30 ns which is a value considered as a lower limit of dead time as it is of the order of the recovery time of a APD cell (Grodzicka et al 2011, Rosado andHidalgo 2015). Intrinsic activity of LaBr 3 (less than 0.3 Bq cm −3 ) has not been modeled since it leads to a negligible number of coincidences within the typical activities in hadron therapy.
Since this study is intended to evaluate the contribution of the background events to the reconstructed image, no energy and spatial resolution have been applied to the singles, which further degrades the image but have no influence in the background and random events. The data obtained in the simulations are processed offline. In each individual detector, an energy threshold of 10 keV has been imposed, below which singles have not been accepted. Singles taking place in two different planes within a time interval (TCW) are sorted into coincidences. The hits associated to each coincidence are examined to extract the information about the processes involved.
Concerning the categorization of events in the Compton camera layout, coincidences are classified by examining the interactions or hits in the scintillator planes. Signal events are defined as events where a photon suffers a Compton scattering in the first scintillator plane followed by any interaction of the Compton-scattered photon in the second scintillator plane (which may be one or more Compton scatters, photoelectric effect or pair production). These events provide the necessary information to retrieve the Compton cone containing the position of the photon emission (Wilderman et al 1998). The axis of the cone is defined by the interaction points in the two planes, the apex is defined by the interaction point in the first plane and the aperture angle is given by the Compton kinematics recovered from the energy information measured by the crystals. The remaining coincidences do not allow us to retrieve the Compton cone in which the photon emission point is located and therefore constitute the background. They are detailed hereafter and outlined in figure 2. The abbreviation that will represent each group throughout the study is indicated in italics in the following list.
(i) Back-scattered events (back): the first Compton scattering of the incoming photon takes place in the axially furthest plane of the Compton camera and the scattered photon is detected afterwards in the first plane. A time resolution of 150 ps would be needed to properly associate the order of photon interactions. The sequence of the interactions may be incorrectly assessed due to the limited time resolution in experimental systems. (ii) Coincidences where the primary gamma has gone through multiple scatterings in the first plane on which it has interacted (mult).
(iii) Events with gamma conversion into electron-positron pair in the first interaction (g-o). This process involves photons with energies greater than 1.022 MeV and it becomes more intense as the photon energy increases. The coincidence event is triggered either by the detection of one element from the electronpositron pair or by the detection of an annihilation photon.
(iv) Other photonic interactions in the detectors, such as coincidences started by photoelectric effect in the first plane or a Compton scatter where the recoil electron has acquired enough energy to escape the detector and trigger a detection in other plane (other).
(v) Coincidences where the primary gamma has been previously scattered in the phantom before arriving to the detector and therefore it cannot be correlated to the delivered dose profile (phant).
(vi) Coincidences triggered by positrons and electrons emitted by the phantom (e + /e − ).
(vii) Coincidences triggered by neutrons generated in proton beam interactions in the phantom (neut).
(viii) Coincidences originated by the detection of secondary protons produced by nuclear reactions in the phantom (prot).
(ix) Accidental coincidences of any kind triggered by any type of incoming particles generated in different events (randoms).
As stated in section 2.1, the proton beam energies chosen for the study are 90, 150 and 230 MeV, which range approximately from the lowest to the highest proton beam energies typically used in hadrontherapy (Mohan 2022). The proton beam intensity (I) has been varied from a low beam current of 16 pA or 10 8 protons per second (prot/s) to a beam intensity of 1.6 nA or 10 10 prot/s. The beam intensity value of 10 10 prot/s is of the order of those used clinically, around 1.2 × 10 10 prot/s (Pausch et al 2016). Other delivery systems as Varian Pro-Beam can use an even lower dose rate of 1.2 × 10 9 prot/s (Polf et al 2021).
To study the dependence of the random coincidences with the TCW, this parameter has been set to 500 ps, 1 ns and then varied from 10 ns to 100 ns in steps of 10 ns (section 3.2). A TCW of 1 ns has been applied to report the background composition in section 3.3, for all proton beam energies and intensities considered. Values of 500 ps, 10 ns and 50 ns TCWs have been chosen to evaluate the effect of randoms in the reconstructed images in section 3.4. The particular case of a proton beam energy of 150 MeV, I = 10 9 prot/s and TCW = 50 ns (experimental case of MACACO II and MACACO III tests) is selected in figure 9 to show the contribution of different kind of background events to image degradation. Coincidences that have a total energy deposition below 0.05 MeV and above 9.95 MeV have been removed of the analysis as they are not useful for the spectral image reconstruction.

Image reconstruction
Images are reconstructed with the list-mode maximum likelihood expectation maximization algorithm, by means of the spectral reconstruction method for two-interaction events . This method estimates both spatial and energy distributions. The image is built in a four-dimensional field of view (FoV), with three spatial coordinates of the emission point and the emission energy as the fourth coordinate. Each of the photon emission energies considered in the fourth coordinate defines a different conical surface characterized by the Compton formula. The physical model implemented in the System Matrix introduces processes for photons up to 10 MeV .
A spatial FoV of 201 × 201 × 1 and spatial voxels of 3 × 3 × 3 mm 3 are employed for all images. The FoV is parallel to the Compton camera planes, positioned at the center of the phantom (z = 0). The transverse view of the FoV (x-y coordinates) is the relevant one to the study, showing the photon emission distribution along the beam direction. Regarding the energy dimension, 99 energy bins of 0.1 MeV have been set in the energy range from 0.05 to 9.95 MeV. A Gaussian filter with a standard deviation of 2 pixels in each transverse spatial direction was applied after the final iteration. From the four-dimensional image computed by the spectral reconstruction method, the reconstructed energy spectrum can be recovered integrating out the three coordinates of the spatial image. Similarly, the spatial distribution at the whole energy range can be obtained by integrating out the spectral domain. The spatial distribution from specific energy bins is extracted as well. In particular, we examine the spatial image associated to energy bins containing energies between 4.3 and 4.5 MeV, so that, in this work, the spectral range is restricted by this method and not by an energy threshold. This interval contains the 4.4 MeV spectral line which is highly correlated to the Bragg peak (Verburg et al 2013, Kelleter et al 2017).
The contrast-to-noise ratio (CNR), calculated as CNR = |S − μ|/σ, is used as a figure of merit in the reconstructed spatial images . S is defined as the maximum intensity in the PG emission area. The mean value, μ, and standard deviation of the background, σ, are calculated by taking into account the intensity values from the voxels located in the adjacent zones to the PG emission area, leaving the borders of the image out of the calculation.
Longitudinal profiles in the spatial images are plotted after summing the ten central bins in the beam transverse direction, which cover the PG emission area. From the innermost peak along the proton beam pathway, the position at 80% (R80) and 50% (R50) of the maximum height are calculated after removing the baseline (Ortega et al 2015). The uncertainties of the estimated ranges are assessed from the linear interpolation error and are strongly limited by the voxel size. The 15th iteration was selected as from this iteration onward, the measured R80 and R50 values in the reconstructed images stabilize. Figure 3(a) shows the spectra measured experimentally and simulated for the physics list selection (see section 2.2). The QGSP_BERT_HP_EMZ list manifests an oversimplified PG emission model, as it has been stated in Jarlskog and Paganetti (2008). The 4.4 MeV peak is not visible in the spectrum, but its double escape peak (3417 keV) can be slightly distinguished. Comparing the regions around the identifiable peaks (511, 718 and 3417 keV), the spectrum profile provided by QGSP_BIC_HP_EMZ list shows a better agreement with our measured data with respect to QGSP_BIC_AllHP_EMZ, and consequently this physics list has been chosen for proton simulations. With the chosen physics lists, the energy spectra of the particles emitted by the PMMA phantom are shown in figure 3(b).

Effect of the beam current
The absolute number of total coincidences for 10 8 simulated primary protons at different energies and for different beam intensities are plotted in figure 4, where the number of random coincidences found for each case is plotted in crosses.
The rise of random coincidences with the TCW and the proton beam energy can be observed. The decrease in the number of coincidences at I = 10 10 prot/s in the simulation of 230 MeV at larger TCWs is due to the saturation effect in the detectors caused by the modeled dead time.

Background composition
The coincidence event structure, as defined in figure 2, is shown in figure 5 for a fixed TCW of 1 ns. Neutroncaused coincidences is the largest background category for the lower intensities, being responsible for nearly 20% of events for 150 MeV simulations, and this category increases with the proton beam energy. At I = 10 10 prot/s, random coincidences start to dominate, for medium and high proton energies. On the contrary, coincidences generated by protons constitute the least common category recorded. In spite of the high generation of secondary protons, most of these coincidences are withdrawn as the energy deposited in the scintillators generally falls above the upper energy cut in the event selection for the spectral reconstruction, 9.95 MeV. The energy deposition of protons in the first plane (E 1 ) versus the deposition in the second plane (E 2 ) is shown in figure 6. As seen in this figure, protons that escape from the phantom and pass through both planes are more likely to have high energies (generally greater than 20 MeV), and deposit a large amount through hadronic ionizations in the scintillator. Figure 6 is characterized by two accumulations of events. The left accumulation, which ranges approximately from the energy pairs (E 1 , E 2 ) = (13 MeV, 13 MeV) to (E 1 , E 2 ) = (20 MeV, 35 MeV), belongs to protons that completely pass through the detectors-they deposit more energy in the second plane, as more energy is deposited in the first plane-. In the right accumulation, ranging approx. from (E 1 , E 2 ) = (20 MeV, 40 MeV) to (E 1 , E 2 ) = (35 MeV, 5 MeV), protons are absorbed in the second plane.
The most predominant group generated by photons is the gamma conversion into electron-positron pair (g-o). The proportion of signal events ranges from 23.3% (90 MeV proton beam, lowest beam intensity) to 14.5% (230 MeV proton beam, highest beam intensity).

Reconstructed images
We make special emphasis on the case of the 150 MeV proton beam, as it is a medium proton beam energy and the beam energy at which we carried out past beam tests . Figure 7 shows the spatial reconstructed images in the 150 MeV proton beam simulation, for the 4.3-4.5 MeV spectral bins. Large (50 ns), medium (10 ns) and small (1 ns, 500 ps) values of TCW have been selected to sort the simulated data. In addition, the image reconstructed just with signal events is presented in the fifth column. In all cases, CNR values are indicated at the top of the image. Finally, the respective longitudinal projections are plotted in the last column. An improvement of the CNR at shorter TCW, as well as a visual enhancement of the images and profiles can be appreciated for high beam intensities, where the background contribution is greater. For I = 10 8 prot/s, the proportion of randoms is not significant, and thus using a very small TCW does not bring a significant improvement to the image and projections. A small TCW is specially necessary at high intensities, as the fall-off of the reconstructed distributions can no longer be properly distinguished with a TCW of 10 ns and a proton beam intensity of 10 10 prot/s.
In figure 8(a), we plot the R80 and R50 values extracted from the projections in figure 7, where we see that only small TCWs retrieve the fall-off position at all three intensities. The values belonging to the case of TCW = 50 ns and I = 10 10 prot/s are omitted because their position could not be determined due to the high noise. Their corresponding CNR values, plotted in figure 8(b), are specially low at I = 10 10 prot/s, for medium and high TCWs. With a proton beam intensity of 10 9 prot/s, the energy spectra recovered after the spectral reconstruction are represented in figure 8(c). This plot shows that the contribution of background coincidences to the reconstructed spectrum exists along all the energy domain. Spectra variations between different TCWs (which bring more random coincidences at larger values) are slight.
Under these same conditions and employing a TCW of 50 ns, with which we see the effects of random coincidences more clearly, figure 9(a) shows the reconstructed images taking into account different sources of background separately, where we can appreciate their influence in the images and in the fall-off region. Images integrating over the spectral domain and their analogous 4.3-4.5 MeV energy bins are shown above and below in each category. The categories only background and signal + all background have been reconstructed without and with the signal events, respectively. The effect of eliminating randoms can be seen in no randoms. The ideal case, where the image is reconstructed with just signal events, is shown in the category signal. In the signal + photonic bkg, signal events have been reconstructed along with the background generated only by photons of the same event in the detector (groups from to (i) to (v), described in section 2.4). Categories signal + neutrons, signal + protons, e − , e + and signal + randoms follow the same approach. Figure 9(b) plots the R80, R50 and CNR values of the images in figure 9(a) belonging to the transverse spatial image of the 4.3-4.5 MeV energy bins (second and fourth rows projections in figure 9(a)), proving the need of removing random events (i.e. having TCWs as small as possible) and photonic background to enhance the quality of reconstructed images.   Figure 10(a) shows the pie charts of the three simulated scenarios in which the phantom radius and the Compton camera distance are both varied, being the first pie chart already presented in figure 5.

Phantom size and Compton camera distance
As expected, the rate of coincidences produced by gammas scattered in the phantom (phant) increases in larger phantoms. In addition to differences in the solid angle, the fact that the phantom is larger has led to changes in the production of secondary gammas and neutrons, which makes oscillations in the percentages of the neut category. The enlargement of the target dimensions as a cause of the difference in the production of neutrons and secondary gammas has been verified by simulating the Compton camera at such distances from the beam axis but with a constant radius of 30 mm. Results are shown in figure 11, in which we see that the neutron percentages do not vary considerably.
In figure 10(b) the reconstructed spatial images and their projections are shown, for the three scenarios detailed in figure 10(a). We observe that, at 270 mm (phantom radius 200 mm) the PG distribution when all events are reconstructed (first row, third column) is no longer visible due to the decrease in statistics. Even so, when the image reconstruction is performed by only signal events, the PG fall-off is observed in spite of the fewer events employed in the image reconstruction. Since the amount of random events at TCW = 1 ns is low, this result encourages the investigation of strategies that reduce the background beyond the improvement of temporal resolution, specially when dealing with reduced statistics scenarios.

Discussion
Background events are known to constitute a significant source of degradation in Compton cameras for hadron therapy treatment monitoring (Krimmer et al 2018, Parodi and Polf 2018, Kohlhase et al 2019. The studies conducted in this work disentangle the contribution of different sources of background and their effect in the reconstructed images. Figures 4 and 5 show how the background composition depends on the beam energy. In our approach, background generated by the detection of non-PG particles produced by the interaction of the proton beam on the phantom has been contemplated and differentiated from the background events generated by interactions of photons in the detectors. From the former kind, and rejecting protons based on energy cuts, figure 5 shows that neutrons are the main cause of coincidences produced by secondary radiation for 150 and 230 MeV (with 19.4% of events with deposited energy below 10 MeV at 150 MeV and I = 10 9 prot/s). Due to well-known limitations of nuclear reaction models used in Monte Carlo codes and uncertainties in experiments, several projects (Farah et al 2014, Englbrecht et al 2021 find difficulties on validating the neutron models, specially when considering high energy neutrons. De Saint-Hubert et al (2022) estimates the difference of fluence of neutrons simulated with three Monte Carlo codes (FLUKA, MCNPX, Geant4) with a proton beam in a water target. The total neutron fluence in lateral position with respect to the beam has a variation of a 20% between codes. Although these results could be not directly transferable to our case because we have used different beam energies, target material and version of Geant4, we could roughly expect that, in the aforementioned case of 150 MeV and I = 10 9 prot/s, the percentage of coincidences produced by neutrons vary between 15% and 23%.
Charged particles arriving from the phantom (electrons, positrons and protons which deposit energies under 10 MeV) are the major responsible of secondary radiations in simulations at 90 MeV, the lowest simulated . The Compton camera distance to the beam axis is also shown in the left title. The first row is reconstructed with all the events (signal and all background) and the second row only with signal events. The number of events used in the image reconstruction process is indicated in the right title. Figure 11. Percentage of coincidence events in a 150 MeV proton beam at a medium intensity of 10 9 prot/s and TCW = 1 ns with a constant phantom radius, at different distances of the Compton camera from the beam axis. proton beam energy (15.8%). While protons deposit high energies in the detectors and can be discarded easily, neutrons leave signals at energies similar to photons. The increase of neutron background with the proton beam energy has a direct consequence in the rise of the rate of singles and therefore in the random coincidences.
The study presented here takes into account previous studies in this area and goes beyond with a different approach.
In Rohling et al (2017), the methodology used for event classification is similar to the one employed in this work, but limited only to events produced by PG and neglecting randoms. The general agreement is good, although the proportions obtained for the different event categories vary due to the different geometries and materials of the simulated systems. For instance, they report a 34% of pair production and 17% of back-scattered events in simulations of PG, corresponding to proton beams 100 and 120 MeV at 1.2 × 10 10 prot/s. In the case of the proton beam at 90 MeV and 10 10 prot/s simulated by us, the equivalent percentages are 34% and 14%, respectively, once we remove from the calculation all the coincidences generated by particles that are not photons nor randoms, which (Rohling et al 2017) do not consider. Simulations performed by Solevi et al (2016) were carried out to support the understanding of measured data by an early MACACO prototype from a PMMA phantom irradiated by a 150 MeV proton beam. While they also do not consider coincidences generated by secondary particles, they report 30% of coincidences due to pair production and 10% of back-scattered events. Panthi et al (2020) perform an overview of the fluence of different types of particles that reach the Compton camera and interact within the crystals, though the displayed results are not analyzed in coincidence. As in our case, neutrons are found the main type of secondary radiation reaching the Compton camera. The background study carried out by Ortega et al (2015) centers on the event categorization with a different strategy, since their purpose is to quantify the random coincidences of double and triple events at different TCWs. It shows that images reconstructed with triple events produce much larger discrepancies than those reconstructed with double events, due to their lower statistics. The results are consistent with ours when it shows that the effect of random coincidences is sensitive to TCW and beam energy, due to the increase of the detected rate. In particular, their percentage of coincidences generated from neutrons range from 10% for 100 MeV (TCW = 0.75 ns) to 40% (180 MeV and TCW = 2 ns), which is comparable to us (between 12.1% and 32.8% at TCW = 1 ns, depending on the beam energy and beam intensity, figure 5). In Fontana et al (2019), the variation of some types of coincidences as a function of the gamma ray energy is studied. The increase of electron escape events at high photon energies is evaluated. Fontana et al (2019) also quantify the appearance of random coincidences at high beam intensities. In this study, random coincidences and beam time structure effects are included. Their impact on the Compton camera background at high beam rates have been taken into account. At 160 MeV and for the beam intensity value we have used of 10 10 prot/s, they report a signal percentage of 5% with no event selection in their particular Compton camera. They state that a beam intensity reduction at the beginning of the treatment would be required to increase the signal-to-noise ratio and perform range verification.
In the context of hadron therapy monitoring, Compton cameras would aim to determine the fall-off position with the highest possible accuracy and precision, which gives information of the range of the beam on the patient. Thus, we seek to assess how the background (generated by randoms or other types) affects the longitudinal profile fall-off. Figures 7 and 8 are intended to evaluate the effect of different TCW values on the longitudinal projection and also on the determination of the R80 and R50 position. The precision on the determination of the range by the R80 and R50 determination highly depends on the statistics of the measurements (ultimately, the number of primary protons simulated) (Fontana et al 2019, Polf et al 2021. Thus, our results are only meaningful for the number of primary protons of our dataset (5 × 10 9 simulated primary protons). When a TCW of 500 ps is applied, R80 and R50 values at all intensities are closer to the values calculated from the image reconstructed just with the signal events. Even so, other types of background that are not random events affect to the determination of the range. The improvement is also quantified with the CNR values. At I = 10 10 prot/s, which is the intensity closest to that used in hadron therapy, CNRs increase from 0.7 to 12.7 (images with protons of 150 MeV). Figure 9 reinforces the point that random coincidences should be minimized with TCWs as short as possible as they contribute the most to image degradation. In addition, the background derived from the room component, mentioned in section 2.1, could further increase the number of random coincidences detected. Photonic background also contributes significantly to the image degradation. Although coincidences generated by charged particles and neutrons (with the subsequent photons generated in the neutron capture) do not affect significantly the image and range estimation by themselves, we must take into consideration that, as said before, a high rate of this type of events is also contributing to generate random coincidences. While the R80 and R50 positions in these two last cases are within the error of the analogous positions determined in the signal image, the R80 and R50 values in signal+rand are approximately 3 and 2 mm offset, respectively. Still, the noise present in the image in the category no randoms, which has all the other types of background, encourages us to consider the reduction of other types of background as well. Same conclusion arises from section 3.5, in which even with a low presence of random events, distributions with low statistics and all background events are largely degraded, but with the background subtracted (and consequently with even lower statistics), the PG fall-off is distinguishable again.
In addition to the background rejection methods mentioned in section 1, the emerging use of neural networks without requiring additional instrumentation, is also being explored by this  and other groups (Zoglauer and Boggs 2007, Polf et al 2021, Kozani and Magiera 2022. Future work will focus on the improvement of the timing capabilities of the MACACO experimental prototype and other strategies to reduce the background in this system, such as the detection of secondary charged particles or the use of neutron absorbers before the system.

Conclusion
The work presented here carries out a detailed study of the sources of background existing in hadron therapy treatment monitoring with Compton cameras and their impact on the reconstructed images. The study shows how detectors with time resolution below 1 ns are needed to retrieve the fall-off position with better precision, as while at low intensities random coincidences are negligible, at high intensities they can represent between 7.8% and 25.2% of total coincidences. The contribution of background generated by secondary particles arriving from the phantom (neutrons and charged particles) in all cases represents at least a quarter of the total events detected, and therefore, other types of background rejection would be desirable to enhance the performance of Compton cameras in such conditions.