Modelling of energy-dependent spectral resolution for SPECT Monte Carlo simulations using SIMIND.

PURPOSE: Monte Carlo (MC) modelling techniques have been used extensively in Nuclear Medicine (NM). The theoretical energy resolution relationship ( ∝ 1 / E ), does not accurately predict the gamma camera detector response across all energies. This study aimed to validate the accuracy of an energy resolution model for the SIMIND MC simulation code emulating the Siemens Symbia T16 dual-head gamma camera.
METHODS: Measured intrinsic energy resolution data (full width half maximum (FWHM) values), for Ba-133, Lu-177, Am-241, Ga-67, Tc-99m, I-123, I-131 and F-18 sources in air, were used to create a fitted model of the energy response of the gamma camera. Both the fitted and theoretical models were used to simulate intrinsic and extrinsic energy spectra using three different scenarios (source in air; source in simple scatter phantom and a clinical voxel-based digital patient phantom).
RESULTS: The results showed the theoretical model underestimated the FWHM values at energies above 160.0 keV up to 23.5 keV. In contrast, the fitted model better predicted the measured FWHM values with differences less than 3.3 keV. The I-131 in-scatter energy spectrum simulated with the fitted model better matched the measured energy spectrum. Higher energy photopeaks, (I-123: 528.9 keV and I-131: 636.9 keV) simulated with the fitted model, more accurately resembled the measured photopeaks. The voxel-based digital patient phantom energy spectra, simulated with the fitted and theoretical models, showed the potential impact of an incorrect energy resolution model when simulating isotopes with multiple photopeaks.
CONCLUSION: Modelling of energy resolution with the proposed fitted model enables the SIMIND user to accurately simulate NM images. A great improvement was seen for high-energy photon emitting isotopes (e.g. I-131), as well as isotopes with multiple photopeaks (e.g. Lu-177, I-131 and Ga-67) in comparison to the theoretical model. This will result in accurate evaluation of radioactivity quantification, which is vital for dosimetric purposes.

Introduction

Monte Carlo (MC) modelling techniques have been used extensively in the field of nuclear medicine (NM) over the past 50 years. These techniques are ideal in NM because of the stochastic nature of radiation emission, transport, and the detection process. Essential applications of MC modelling techniques in NM include collimator and detector design, image reconstruction and scatter correction techniques, as well as internal dosimetry. MC techniques are also useful in studying phenomena that may be impractical, or impossible to measure, such as collimator septal penetration and scatter. Various MC simulation codes with application to NM imaging, such as SIMIND [1], SimSET [2] and GATE [3], are available, illustrating the significance of MC simulations in NM.

The importance of patient-specific targeted radionuclide therapy for cancer therapy has increased over the past decade [4, 5, 6, 7, 8]. Thus, accurate image quantification has become increasingly important in NM. MC techniques have been used extensively to evaluate quantification procedures in NM imaging [1, 9, 10, 11]. When making use of MC simulation software, it is essential to ensure that the intrinsic characteristics of the imaging system are defined correctly. One such characteristic is the energy resolution, which is described as the system's ability to distinguish between two radiation energies with a small energy difference [12]. Factors such as variations in the number of light scintillation photons generated in the crystal, the number of photoelectrons produced in the photocathode, as well as the multiplication factor of the photomultiplier tubes, all contribute to variations in the signal amplitude produced by the photomultiplier tube. For every keV of absorbed energy in a NaI(Tl) crystal, approximately 38 light photons are emitted. This relates to a single light photon per approximately 30 keV of absorbed energy [13]. Depending on the quantum efficiency of the photocathode, number of photoelectrons are released from the photocathode, upon absorption of a light photon. For events where the crystal has absorbed the same amount of energy, the resultant photomultiplier tube signal amplitude will vary depending on the number of photoelectrons released. If we assume Poisson statistics, the variations can be described by the standard deviation, . Ultimately, this results in an energy spectrum, which is not a narrow line, but a Gaussian-shaped response, as shown in Figure 1. A system will thus be able to resolve two energies if they are separated by at least the value of the system's full width half maximum (FWHM), which is the width of the distribution at half the maximum value of the Gaussian photopeak, for a given photon energy.

Figure 1

Definition of FWHM.

Considering the number of photoelectrons produced depends on the absorbed energy, (, one can assume the simple energy resolution relationship [13, 14],

It has been noted that this simplified relationship does not accurately predict the gamma camera detector response across all energies. When making use of the relationship above, the FWHM values at higher energies are underestimated [15, 16, 17]. Any under- or overestimation of the energy resolution can ultimately affect the accuracy of the simulated image, which may lead to poor agreement between measured and simulated images. Inaccuracies could influence the quantification accuracy of simulated images, which is vital in the process of performing internal dosimetry [17, 18].

Hakimabad et al. [19] proposed a nonlinear response function (Eq. (2)) to model the energy resolution of a 3 × 3 in. NaI scintillation detector;with the photon energy, and , and curve-fitting parameters. The use of this function to model the energy resolution of a gamma camera has not yet been reported in literature and was used in this study.

The study aimed to validate the accuracy of the proposed energy resolution model for the SIMIND MC code emulating the Siemens Symbia T16 dual-head gamma camera. This was achieved by fitting the proposed energy resolution model to measured energy resolution values as a function of photon energy, ranging between 27.0 and 637.0 keV. The accuracy of the fitted model was then validated by comparing measured and simulated energy spectra obtained from in-air and in-scatter measurements. Finally, the energy spectra from isotopes simulated in a voxel-based digital patient phantom were compared and evaluated.

Materials and methods

A dual-head Siemens Symbia T16 (Siemens Healthcare, Erlangen, Germany) SPECT gamma camera is used for routine clinical studies at Universitas Academic Hospital, Bloemfontein, South Africa. The system is equipped with low energy high resolution (LEHR), medium energy low penetration (ME) and high energy (HE) collimators which were used for the respective extrinsic experimental measurements in this study. The Siemens Symbia T16 gamma camera was modelled using version 6.2 of the SIMIND MC code [1]. The SIMIND MC code, written in FORTRAN-90, is a program, which simulates a standard clinical SPECT camera. With simple modifications, SIMIND can be used to perform a variety of calculations or measurements typically encountered in the field of NM.

When creating the virtual Siemens Symbia T16 system using the SIMIND MC code, all detector (crystal and collimator) characteristics were defined according to the Siemens Symbia T Series specification sheet [20]. The thickness of the aluminium cover was set to 0.1 cm. SIMIND does not model the photomultiplier tubes nor any of the electronics; thus, a layer of backscatter material, mimicking the photomultiplier and electronic assembly, was simulated. This material has a density of 1.47 g.cm−3 and consists of boron (4.0%), oxygen (54.0%), sodium (2.8%), aluminium (1.2%), silicon (37.7%) and potassium (0.3%). The thickness of backscatter material which mimicked the photomultiplier and electronic assembly of the physical gamma camera was predetermined and set to 7.5 cm. The intrinsic spatial resolution of the physical gamma camera was obtained from acceptance test results and used as input into SIMIND. This ensured that the intrinsic characteristics of the virtual gamma camera are the same as that of the physical gamma camera.

The physical properties of the isotopes used in this study and incorporated in SIMIND are shown in Table 1 [21]. The backscatter photon energy refers to the energy of the photon, originating from the backscattering medium (photomultiplier and electronic assembly), which has scattered at a 180° angle, before being absorbed in the crystal [13].

Table 1

Characteristics of the isotopes used in this study, including photopeak energies, abundances and backscatter photon energies [21].

Isotope	Photopeak energy (keV)	Abundance (%)	Backscatter photon energy (keV)
Am-241†	59.5∗	35.92	48.3
Ba-133	30.8§	48.10	27.5
	35.5§	11.34	31.2
	80.9∗	33.31	61.4
	356.0	62.05	148.7
	383.8	8.94	153.4
F-18	511.0∗	193.72	170.3
Ga-67 ◊	93.3∗	38.10	68.3
	184.8∗	20.96	107.2
	300.2∗	16.60	138.0
	393.5	4.59	154.9
I-123†,‡,◊	27.3§	35.34	24.7
	159.0∗	83.25	98.0
	528.9	1.28	172.3
I-131†,‡,◊	364.5∗	81.20	150.2
	636.9∗	7.12	182.3
	722.9	1.79	188.8
Lu-177†,◊	55.2§	2.19	45.4
	112.9∗	6.20	78.3
	208.4∗	10.38	114.8
Tc-99m†,‡	140.5∗	88.50	90.7

energy photopeaks included in the fitted model.

isotopes considered for the extrinsic in-air energy spectra measurements and simulations.

isotopes considered for the extrinsic in-scatter energy spectra measurements and simulations.

isotopes considered for the voxel-based digital patient phantom simulations.

average energy and abundance for multiple emissions.

Modelling of intrinsic energy resolution

Intrinsic energy spectra: measurements and modelling

Radioactive sources listed in Table 1 were positioned at 30 cm from one of the gamma camera detectors without a collimator and energy spectra were acquired from the source in air. Gamma camera acceptance tests showed little difference between the two detectors; thus, only a single detector was used in this study. The opposing detector was extended as far back as possible to minimise the contribution of backscatter from the second detector. The experimental setup is shown in Figure 2. Intrinsic energy spectra, with a digital resolution of 1.3 keV per channel, were acquired to obtain 30 000 counts in the channel with the highest count contribution. All intrinsic energy spectra acquisitions were repeated three times. The acquired energy spectra were exported from the manufacturer's computer and analysed using the public domain software, ImageJ (version 1.52r) [22]. The main emission photopeaks were fitted with a Gaussian function, and the FWHM (keV) values were calculated for each photopeak. The measured intrinsic FWHM values (keV) were plot as a function of photon energy, and Eq. (2) was fitted (using ImageJ) to the FWHM values.

Figure 2

Experimental setup for measurement of in-air intrinsic energy spectra.

Computed tomography (CT) images were acquired of the setup, for use in the simulation. The CT acquisition was carried out in a 512 × 512 image matrix with a pixel size of 0.127 × 0.127 cm2 and a slice thickness of 0.5 cm.

Intrinsic energy spectra: Monte Carlo simulations

SIMIND allows for three methods to incorporate energy resolution, based on (i) a theoretical model defined by Eq. (1), normalised to Tc-99m, (ii) a predetermined fixed FWHM value for a specific photon energy, or (iii) a fitted model estimated from measured data, defined by Eq. (2). When making use of method (i) and (ii), SIMIND requires the measured FWHM value of Tc-99m, and that of the relevant photon energy, respectively. Method (iii) requires the user to define the fitting parameters of the function (Eq. (2)) used to model the energy resolution. In this study, energy resolution estimation methods (i) and (iii) were utilised. Neither of these two methods require the user to determine the gamma camera energy resolution for the isotope of interest.

The CT images acquired with the CT of the physical gamma camera were segmented using ITK-SNAP version 3.6.0 [23]. The segmented images were used to create a voxel-based digital phantom, as described by Ramonaheng et al. [24]. This ensured that the setup for the simulation was identical to that of the measurement. Each segmented region was assigned a unique value, as shown in Figure 3. This value, together with the radioactivity concentration to be simulated, was defined in a text file. The CT images were converted from Hounsfield units to density values using an in-house software program developed in Visual Basic. The software program incorporates a bi-linear conversion model, which is based on measurements acquired with the Gammex RMI-465 phantom [25, 26]. The density images, along with the segmented images and text file, served as input for the SIMIND MC code.

Figure 3

Transverse slice of the segmented CT image of the source, Styrofoam block and patient bed.

Intrinsic energy spectra were simulated for each of the isotopes listed in Table 1 using both the theoretical and the fitted models. The stop condition for each simulation was determined from the duration time acquiring the energy spectra on the physical gamma camera. A high number of histories (>1 billion) were simulated for each energy spectrum. The FWHM was calculated for each main emission photopeak, in the same manner as for the measured energy spectra.

All MC simulations were performed on the computer cluster of the High Performance Computing unit situated at the University of the Free States. The cluster has 36 computer nodes which perform calculations with 5560 CPU cores and 13.8 terabytes of system memory. The simulations in this study made use of 2 of the 36 computer nodes, with 100 CPU cores per node.

Validation of the fitted energy resolution model

To validate the fitted model, a series of extrinsic (detector fitted with an appropriate collimator) energy spectra, were both acquired, on the Siemens Symbia T16 gamma camera, and simulated. These energy spectra included in-air and in-scatter measurements and simulations for the isotopes indicated in Table 1 (in-air (†) and in-scatter (‡)). All in-air acquisitions were repeated three times. However, due to long acquisition times and the lack of availability of the clinically used gamma camera, in-scatter acquisitions were only performed once.

Extrinsic energy spectra: measurements

A series of in-air extrinsic energy spectra were acquired on the physical gamma camera for Am-241, I-123 and Tc-99m using the LEHR collimator, I-123 and Lu-177 with the ME collimator and I-131 using the HE collimator, at a distance of 30 cm from the detector. The experimental setup is shown in Figure 4(a). Additionally, in-scatter extrinsic energy spectra were acquired for Tc-99m and I-123 using the LEHR collimator, I-123 with the ME collimator and I-131 using the HE collimator. Similar to the in-air measurements, the individual sources were positioned 30 cm from the detector, with an added 15 cm of Perspex®, as shown in Figure 4(b). The energy spectra were processed in the same manner as for intrinsic energy spectra. FWHM values were calculated for each main emission photopeak as listed in Table 1 (in-air(†)). As before, CT images were acquired of each experimental setup for use in the simulation procedure.

Figure 4

Experimental setup for measurement of extrinsic (a) in-air and (b) in-scatter.

Extrinsic energy spectra: Monte Carlo simulations

In-air and in-scatter extrinsic energy spectra were simulated, with the appropriate collimators, for isotopes as listed in Table 1 (in-air (†) and in-scatter (‡)), using both the theoretical and fitted models. The stop condition and number of photons simulated per energy spectra were the same as that for the intrinsic energy spectra simulations. The FWHM values for the main photopeaks were calculated for the in-air extrinsic energy spectra.

The overall simulated and measured extrinsic energy spectra were visually evaluated and compared for both the in-air and in-scatter extrinsic energy spectra.

Comparison of the energy resolution models in a simulated voxel-based digital patient phantom

Energy spectra obtained from simulations of a voxel-based digital patient phantom, with both the theoretical and fitted energy resolution models, were compared. To create a voxel-based digital patient phantom, CT image data of a retrospective patient study carried out in Universitas Academic Hospital, was randomly and anonymously selected from the Symbia T16 patient database. The voxel-based digital patient phantom was created as described before (Section 2.1.2). Three different sized spheres (0.5 cm, 3.0 cm and 5.0 cm diameters), mimicking spherical tumours, were digitally added to the thorax of the voxel-based digital patient phantom (Figure 5). Since accurate image quantification plays an increasingly important role in NM, isotopes typically used in targeted radionuclide treatment procedures as well as isotopes with multiple photopeaks (listed in Table 1 (◊)) were considered.

Figure 5

Schematic showing the position of the three spherical tumours in the patient phantom.

I-123, Lu-177, Ga-67 and I-131 radioactivity concentrations were assigned respectively to the tumours, lungs, liver and soft tissue of the patient phantom (tumour to background ratio of 100:1, tumour to liver ratio of 100:7.5 and tumour to lung ratio of 100:3.3). The radioactivity concentration values were based on clinical SPECT patient data. Anterior whole-body images were simulated for, I-123 with the LEHR collimator, I-123 and Lu-177 with the ME collimator, and I-131 with the HE collimator. Energy spectra generated from the simulated data, with both the theoretical and fitted energy resolution models, were visually evaluated.

Results

The measured intrinsic energy resolution expressed as FWHM (keV) for photopeaks listed in Table 1 were plot as a function of energy in Figure 6. As a comparison, the theoretical and fitted modelled energy resolution values are also reported on the graph. The fitted model was determined aswith an r2 value of 0.993.

Figure 6

Comparison between the measured and simulated (theoretical and fitted models) intrinsic energy resolution.

It is important to note that all the isotopes with multiple photon energies that could not be resolved by the gamma camera (i.e. Ba-133: 30.6, 30.9, 35.0 and 35.9 keV; Ba-133: 356.0 and 383.4 keV; I-123: 27.2, 27.4, 31.1 and 31.7 keV; Lu-177: 54.6 and 55.7 keV) were not considered in the model.

Table 2 shows the FWHM values obtained from the measured and simulated (using both the theoretical and fitted models) energy spectra for different photopeaks. FWHM differences between the simulated and measured results are reported for each of the isotopes' photopeak energies. FWHM average and standard deviation values for the three measurements of each photopeak energy are tabulated.

Table 2

Intrinsic measured and simulated FWHM values and differences between measured and simulated data for each isotope at their photopeak energies. Simulations were performed with both the theoretical and fitted models.

Isotope	Photopeak energy (keV)	Measured	Simulated
			Theoretical model		Fitted model
		FWHM (keV)	FWHM (keV)	Difference (keV)	FWHM (keV)	Difference (keV)
I-123	27.3	6.5 ± 0.03	7.0	0.5	5.7	-0.8
Ba-133	30.8	6.4 ± 0.01	7.4	1.0	6.0	-0.4
Lu-177	55.2	8.9 ± 0.23	10.1	1.2	9.0	0.1
Am-241	59.5	7.3 ± 0.04	9.0	1.7	8.2	0.9
Ba-133	80.9	9.1 ± 0.09	10.6	1.5	10.3	1.2
Ga-67	93.3	10.9 ± 0.04	11.3	0.4	11.3	0.4
Lu-177	112.9	12.5 ± 0.19	12.7	0.2	13.2	0.7
Tc-99m	140.5	13.4 ± 0.02	13.6	0.2	14.9	1.5
I-123	159.0	15.1 ± 0.07	14.4	-0.7	16.4	1.3
Ga-67	184.8	18.4 ± 0.22	15.9	-2.5	19.1	0.7
Lu-177	208.4	19.7 ± 0.24	16.5	-3.2	20.2	0.5
Ga-67	300.2	29.0 ± 0.16	19.8	-9.2	27.4	-1.6
Ba-133	356.0	36.4 ± 0.27	22.0	-14.4	33.1	-3.3
I-131	364.5	33.9 ± 0.06	21.8	-12.1	32.5	-1.4
F-18	511.0	41.4 ± 0.14	25.8	-15.6	43.5	2.1
I-131	636.9	52.7 ± 0.74	29.2	-23.5	53.8	1.1

According to Table 2, the largest standard deviation of the FWHM values calculated from the measured data was 0.74 keV for the 636.9 keV photopeak of I-131 (relative standard deviation of 1.5%). Note that the energy spectra were acquired with a digital resolution of 1.3 keV per channel; thus, the standard deviations obtained indicate that the inherent noise in the measurements was acceptable.

The measured and simulated intrinsic energy spectra, based on both the theoretical and fitted models, for each of the isotopes listed in Table 2 are compared in Figure 7(a) – (h), respectively. All energy spectra were normalised and peaked to their respective main emission photopeaks for comparing the relative energy resolution of the spectra. The prominent photopeaks (P), backscatter peaks (B) and Compton edges (C) are indicated on the energy spectra.

Figure 7

Comparison of the measured and simulated intrinsic energy spectra for the theoretical and fitted models for (a) Ba-133, (b) Lu-177, (c) Am-241, (d) Ga-67, (e) Tc-99m, (f) I-123, (g) I-131 and (h) F-18. The photon energy photopeak (P), backscatter peak (B) and Compton edge (C) are indicated on the energy spectra.

Table 3 shows the calculated FWHM values and differences between the measured and simulated extrinsic in-air data for isotopes listed in Table 1(†). As for the intrinsic measurements, average and standard deviation values for the three measurements of each photopeak energy are reported.

Table 3

Extrinsic in-air FWHM values and differences between the measured and simulated data, for respective isotopes' photopeak energies with the highest abundance. Simulations were performed with the theoretical and the fitted models.

Isotope	Photopeak energy (keV)	Measured	Simulated
			Theoretical model		Fitted model
		FWHM (keV)	FWHM (keV)	Difference (keV)	FWHM (keV)	Difference (keV)
Lu-177 ME	55.2	8.7 ± 0.22	10.3	1.6	8.9	0.2
Am-241 LEHR	59.5	6.6 ± 0.01	9.0	2.4	8.3	1.7
Lu-177 ME	112.9	12.3 ± 0.18	12.9	0.6	13.5	1.2
Tc-99m LEHR	140.5	13.1 ± 0.03	13.6	0.5	14.9	1.8
I-123 LEHR	159.0	17.1 ± 0.05	17.4	0.3	20.2	3.1
I-123 ME	159.0	14.8 ± 0.04	14.8	0.0	16.6	1.8
Lu-177 ME	208.4	19.2 ± 0.11	16.6	-2.6	20.2	1.0
I-131 HE	364.5	38.6 ± 0.08	23.5	-15.1	36.5	-2.1

HE = high energy; LEHR = low energy high resolution; ME = medium energy low penetration.

Figures 8 and 9 show the extrinsic in-air and in-scatter energy spectra for isotopes listed in Table 1 (in-air (†) and in-scatter (‡)). The standard deviation values reported in Table 3 for the average measured in-air extrinsic FWHM values indicate a small measurement error similar to the intrinsic measured data in Table 2.

Figure 8

Comparison of the measured and simulated extrinsic energy spectra for (a) Am-241 with the LEHR collimator, (b) Tc-99m with the LEHR collimator, (c) I-123 with the LEHR collimator, (d) I-123 with the ME collimator, (e) Lu-177 with the ME collimator and (f) I-131 with the HE collimator. HE = high energy; LEHR = low energy high resolution; ME = medium energy low penetration.

Figure 9

Comparison of the measured and simulated extrinsic energy spectra with 15 cm scatter media for (a) Tc-99m with the LEHR collimator, (b) I-123 with the LEHR collimator, (c) I-123 with the ME collimator, and (d) I-131 with the HE collimator. HE = high energy; LEHR = low energy high resolution; ME = medium energy low penetration.

Figure 10 shows the simulated extrinsic energy spectra obtained with the voxel-based digital patient phantom for I-123 using both the LEHR and ME collimators, for Lu-177 and Ga-67 using the ME collimator and for I-131 with the HE collimator. These spectra were not normalised to their respective main emission photopeaks. The raw counts, obtained directly from the simulated image, are represented on the y-axis. The difference in sensitivity, as a result of the different energy resolution models, is clearly visible.

Figure 10

Comparison of the simulated extrinsic energy spectra obtained from the theoretical and fitted models of the patient phantom for (a) I-123 with the LEHR collimator, (b) I-123 with the ME collimator, (c) Lu-177 with the ME collimator, (d) Ga-67 with the ME collimator and (e) I-131 with the HE collimator. HE = high energy; LEHR = low energy high resolution; ME = medium energy low penetration.

Discussion

The measured and simulated FWHM values compare well at energies below 160.0 keV for simulations with both the theoretical and fitted models (Figure 6). For photon energies below 160.0 keV, the maximum absolute difference between the measured and simulated FWHM values using the theoretical and fitted models were 1.7 keV for Am-241 (59.5 keV) and 1.5 keV for Tc-99m (140.5 keV), respectively (Table 2). It is important to note that the simulations performed with the theoretical model made use of a predetermined intrinsic energy resolution value for Tc-99m, thus explaining the small difference (0.2 keV) between the measured and simulated intrinsic FWHM value. From Figure 6, it is evident that the fitted model accurately predicts the FWHM values above 160.0 keV, with a maximum absolute difference between the measured and simulated FWHM values of 3.3 keV for the 356.0 keV photopeak of Ba-133. In contrast, the theoretical model gradually underestimates the FWHM values as the photopeak energy increases above 160.0 keV, with a maximum difference of 23.5 keV for the 636.9 keV photopeak of I-131. This underestimation may lead to inaccurate radioactivity quantification when simulating NM images with high-energy photon emitting isotopes (e.g. I-131), as well as isotopes with multiple photopeaks (e.g. Lu-177, I-131 and Ga-67).

Energy spectra shown in Figure 7 are normalised and peaked to their respective main emission photopeaks allowing for easy comparison of energy resolution, however, this makes visual comparison of the measured and simulated energy spectra difficult. Notwithstanding this limitation, measured and simulated energy spectra, using the fitted model, show good agreement. It is important to note that the discrepancies noted on all the measured and simulated energy spectra at the low and high photon energies are due to the cut-off of photon energies at approximately 20.0 keV, respectively, for the Siemens Symbia gamma camera. Thus, photon energies below this are not reflected in any of the measured energy spectra. The cut-off of photon energies below 20.0 keV resulted in decreased amplitudes for the combined 30.8 keV and 35.5 keV photopeak for Ba-133 as well as for the 27.3 keV photopeak of I-123. The small offset between the measured and simulated photopeaks, visible on the spectra of Lu-177, F-18 and the high-energy I-123 and I-131 photopeaks is due to the nonlinear energy response of the detector to Compton and photoelectric events [13], which is not considered in the Monte Carlo simulation. This was also reported by Ramonaheng et al., for Lu-177 [24].

The difference noted between the measured and simulated energy spectra can be attributed to the addition of backscatter photons, originating from the second detector, which was included in the measurement, but not in the simulation. The simulated energy spectra only reflect backscatter originating from the detector components. Figure 7(b) shows the backscatter peak at 114.8 keV, originating from the 208.4 keV Lu-177 photopeak, which can't be distinguished from the 112.9 keV photopeak. This explains the difference between the amplitudes of the measured and simulated energy spectra at 112.9 keV.

The simulated energy spectra obtained with the theoretical model, however, show larger discrepancies when compared to the measured data. This is primarily due to the theoretical model not taking into account the variances of the photomultiplier tube gain and electronics as pointed out by Cherry et al. [13]. Significant differences were visible, on Figure 7, at the higher energy photopeaks (Ba-133: 356.0 keV; Ga-67: 300.2 keV; I-123: 528.9 keV; I-131: 364.5 keV and 636.9 keV and F-18: 511.0 keV).

The small discrepancies visible in the Compton region for Lu-177 can be attributed to minor differences in the physical and simulated source geometries for this isotope. The large differences in the Compton regions for I-131 and F-18 measured and simulated with the theoretical model, is a result of incorrect modelling of the energy resolution at high energies. The limitations of the theoretical model were also reported by Rault et al. [16].

In general, Figure 7 shows that measured energy spectra at the higher energy range (>160.0 keV) were better simulated using the fitted model.

From Table 3, it is evident that the theoretical model underestimates the FWHM values at high photon energies (Lu-177 and I-131), with the largest absolute difference of 15.1 keV for the 364.5 keV photopeak of I-131. The FWHM values obtained with the fitted model is in better agreement with the measured data for Lu-177 (208.4 keV) and I-131 (364.5 keV). The largest discrepancy obtained between the FWHM values of the measured and simulated energy spectra, using the fitted model, was 3.1 keV.

The effect of septal penetration for I-123 with the LEHR collimator is evident, as the extrinsic in-air FWHM value at 159.0 keV (Table 3) is larger than that of the intrinsic FWHM value (Table 2). I-123 with the ME collimator shows a good agreement between the intrinsic and extrinsic FWHM values at 159.0 keV. This is because the effect of septal penetration is less when using the ME collimator. As with I-123 with the LEHR collimator, the effect of septal penetration is evident in the increase in the extrinsic FWHM of I-131 using the HE collimator, in comparison to the intrinsic FWHM.

Differences noted between the measured and simulated extrinsic energy spectra (Figure 8 and 9) were similar to that of the intrinsic energy spectra (Figure 7). As mentioned previously, the discrepancy noted between the measured and simulated energy spectra at low photon energies is due to the photon energy cut-off at 20.0 keV. Backscatter peaks in the extrinsic in-air measured energy spectra (Figure 9) were less prominent than in the measured intrinsic energy spectra (Figure 7) due to collimation. Due to the energy spectra being normalised and peaked to their respective main emission photopeaks, small differences are noted between the measured and simulated high-energy photopeaks of I-123 and I-131 for both the in-air and in-scatter energy spectra. The measured energy spectra at these high-energy photopeaks were better simulated using the fitted model. As for the intrinsic energy spectra, the slight offset noted between measured and simulated photopeaks can be attributed to the nonlinear energy response of the detector.

Figures 8(f) and 9(d) show an overall underestimation of the I-131 lower energy Compton regions when simulating with the theoretical model. This can be attributed to the normalisation of the energy spectrum at the 364.5 keV photopeak. The fitted model better emulates both the in-air and in-scatter measured I-131 energy spectra.

Figure 10 (a) and (b) show minor differences between the energy spectra simulated with the theoretical and fitted energy resolution models for I-123 with both the LEHR and the ME collimator. However, the amplitude of the 159.0 keV photopeak obtained using the fitted model exceeds that of the theoretical model. This may lead to a difference in image quality and quantification. Figure 10 (c to e) shows a large discrepancy between the simulated high-energy photopeaks of Lu-177, Ga-67 and I-131. Differences are also noted in the scatter regions below and above the 208.4 keV Lu-177 photopeak, the 184.8 keV, 300.2 keV and 393.5 keV Ga-67 photopeaks and the 364.5 keV I-131 photopeak. If image quantification includes energy window-based scatter correction, inaccurate quantification may result due to these differences. These differences are more pronounced with increasing photopeak energy, as expected. It is important to note the difference at the high-energy photopeaks for I-123 (528.9 keV) and I-131 (636.9 keV). There is a strong drive to correct for collimator septal scatter and penetration as part of quantification procedures. Thus, incorrect simulation of these photons may result in inaccurate quantification compensation techniques.

Conclusion

The measured and simulated FWHM values compare well at energies below 160.0 keV for both simulations with the theoretical and fitted models. At energies above 160.0 keV, the theoretical model gradually underestimates the FWHM values as the photopeak energy increases, with differences up to 23.5 keV. In contrast, the fitted model accurately predicts the FWHM values, across all photopeak energies, with a maximum absolute difference of 3.3 keV. Intrinsic energy spectra simulated with the fitted model for Ba-133, Lu-177, Am-241, Ga-67, Tc-99m, I-123, I-131 and F-18 compared well to the measured energy spectra. Intrinsic energy spectra simulated with the theoretical model resulted in large discrepancies for I-131 and F-18 as well as for the high-energy photopeaks of Ba-133 and I-123.

Similar results were obtained for the extrinsic energy spectra as for the intrinsic data when simulated using the theoretical model, with FWHM differences up to 15.1 keV. Overall, the energy spectra simulated with the fitted model compared well with the measured energy spectra. The most notable difference between the simulated spectra using the theoretical model and the measured data was seen for the I-131 spectrum. Minor differences noted between measured and simulated energy spectra, in the Compton region, could be attributed to the normalisation at respective photopeaks.

When a scatter medium was introduced, both the theoretical and fitted models show good agreement to the measured energy spectra for I-123 with the LEHR and ME collimators. However, for I-131, the overall energy spectrum simulated with the fitted model better matched the measured energy spectrum. Both the I-123 and I-131 higher energy photopeaks (528.9 keV and 636.9 keV) simulated with the fitted model, resembled the measured photopeaks more accurately.

The energy spectra obtained from simulations in a clinical scenario were compared when using the theoretical and fitted models for I-123, Lu-177, Ga-67 and I-131. Differences observed were more pronounced as the photon energy increased. Image quantification may be affected due to these differences. Accurate image quantification becomes increasingly important with the demand for patient-specific targeted radionuclide dosimetry and therapy.

This study did not make provision for estimating the quantification error that may be introduced due to incorrect modelling of the energy resolution of the gamma camera. A follow-up study should be conducted that considers different patient geometries, with a complete quantification protocol, including scatter and attenuation correction of planar and SPECT images.

This study shows that modelling of the energy resolution across all energies is essential when simulating NM studies with high-energy photon emitting isotopes (e.g. I-131), as well as isotopes with multiple photopeaks (e.g. Lu-177, Ga-67 and I-131). The fitted energy resolution model proposed in this study showed that it was accurate, also under circumstances where scatter was introduced. The comparison between the theoretical and fitted models in the simulated clinical environment shows that possible errors can be introduced in the simulation if the energy resolution model is not selected with great care. These errors can affect the accuracy of radioactivity quantification, which is vital for dosimetry purposes in patient-specific targeted radionuclide therapy.

Accurate modelling of the Siemens Symbia T16 dual-head gamma cameras' energy resolution with the fitted model will result in accurately simulated images for all isotopes used in NM, using the SIMIND MC code. These simulated images can be used to optimise image processing software as well as radioactivity quantification and internal dosimetry.

Ethics approval

This study was performed with approval from the Health Sciences Research Ethics Committee at the University of the Free State, Ethics number UFS-HSD2019/2135/2502. The HSREC functions in compliance with, but not limited to, the following documents and guidelines: The SA National Health Act. No. 61 of 2003; Ethics in Health Research: Principles, Structures and Processes (2015); SA GCP(2006); Declaration of Helsinki; The Belmont Report; The US Office of Human Research Protections 45 CFR 461 (for non-exempt research with human participants conducted or supported by the US Department of Health and Human Services- (HHS), 21 CFR 50, 21 CFR 56; CIOMS; ICH-GCP-E6 Sections 1–4; The International Conference on Harmonization and Technical Requirements for Registration of Pharmaceuticals for Human Use (ICH Tripartite), Guidelines of the SA Medicines Control Council as well as Laws and Regulations with regard to the Control of Medicines, Constitution of the HSREC of the Faculty of Health Sciences at the University of the Free State.

Declarations

Author contribution statement

Morphis, M; van Staden, J: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

du Raan, H; Ljungberg, M: Conceived and designed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement

This work was supported by the National Research Foundation of South Africa (NRF)/Swedish Foundation for International Cooperation in Research and Higher Education (STINT) Science and Technology Research Collaboration (160811184621), as well as by the South African Medical Research Council (SAMRC-RFA-UFSP-01-2013/HARD). Partial financial support was received from the Nuclear Technologies in Medicine and the Biosciences Initiative (NTeMBI) managed by Necsa.

Data availability statement

Data included in article/supplementary material/referenced in article.

Declaration of interests statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.