PET and SPECT benchmarks
Two benchmarks, one for PET and one for SPECT, are included in the Gate distribution. These benchmarks check the integrity of a Gate installation or upgrade. Each benchmark consists of the macros required to run the simulation, analyze the simulation output, and generate figures. In addition, a set of baseline figures are included to compare the users' results with those from a correct run.
- 1 PET Benchmark
- 2 SPECT benchmark
- 2.1 Loading the SPECT benchmark package
- 2.2 SPECT benchmark: moving source
- 2.3 Output
- 2.4 Results: SPECT benchmark
Benchmark configuration description
The PET benchmark distribution contains the following Gate macro files:
benchmarkPET.mac camera.mac phantom.mac physics.mac sources.mac visu.mac
The main macro file is benchmarkPET.mac.
By default, the visualization is disabled. To execute the benchmark without visualization, type
To execute the benchmark with a visualization of the simulation setup, comment (with #) the following line in benchmarkPET.mac:
and uncomment the following line:
You can execute the main macro file right after you get the Gate prompt using:
Gate Idle> /control/execute benchmarkPET.mac
On a 2 GHz processor (Xeon 5130), the computing time will be around 5 hours and the size of the simulation output will be around 72 MBytes.
The PET benchmark simulates a whole-body scanner that does not correspond to any existing system. The scanner is described in the macro file camera.mac. It consists of eight detector heads arranged in a 88 cm diameter and 40 cm axial length octagonal cylinder. Each head is made of 400 detector blocks and each block is a 5x5 array of dual-layer BGO-LSO crystals.
The cylindricalPET system is used to describe the logical organization of the scanner. The relationship between the physical components of the scanner and the hierarchical levels of the system are described in Table 3.1:
|Physical component||cylindricalPET system|
Each head is shielded at both ends of the axial field-of-view with a 1 mm thick rectangular lead plate. In addition, 0.5 mm thick rectangular tungsten plates are placed in front of each head to provide a partial axial collimation (partial septa). There are three plates per head, 10 cm apart along the scanner axis (see Figure 3.1).
The LSO and BGO crystals are defined as crystal sensitive detectors. A Gaussian blurring with a 26 % Full Width at Half Maximum (FWHM) at 511 keV is applied to the energy deposited in the crystals. A detected gamma is further processed only if its energy falls between 350 and 650 keV. A coincidence is defined when two gammas are detected within 120 ns. A large coincidence time window is used in order to record a large number of random coincidences despite low activity levels.
The phantom is described in the macro file phantom.mac. It consists of a 20 cm diameter and 70 cm long cylinder made of water. The cylinder is centered in the field-of-view of the scanner.
The sources are described in the macro file sources.mac. Two 68 cm long and 1 mm in diameter line sources are inserted in the water cylinder, at a radial distance of 2 cm from the cylinder axis. The first source is filled with F-18 and the second with O-15, both with an initial activity of 100 kBq. The fast production module for positron emission is used and the half-lives of F-18 and O-15 are explicitly set to 6586.2 and 122.24 s.
Selection of the physics processes
The selection of the physics processes is described in the macro file physics.mac. The standard electromagnetic package is used for Compton scattering and for photoelectric effect. The low energy package is only used for Rayleigh scattering. To speed up the simulation, X-rays and secondary electrons are not tracked; the electron energy cut is set to an equivalent of 30 cm and the X-rays energy cut to 1 GeV.
The simulated acquisition is 4 minutes in duration, which represents about two O-15 half-lives. The source activities and the acquisition duration are set such that the PET benchmark runs in less than 5 CPU hours on a 2 GHz processor (Xeon 5130). The acquisition is divided into two 2 minute frames. After the first frame, the gantry rotates by 22.5 deg. Only coincident events are recorded and written in the benchmarkPET.root file. All other outputs are disabled. The water cylinder, the end-shieldings, and the axial collimators are defined as "phantom sensitive detectors". This allows for the selection, in the output, of the gammas that underwent a Rayleigh or Compton scattering in these materials before being detected.
Approximately 3.7 x 10^7 decays occur during the simulated acquisition and some 7.0 x 10^5 coincidences are recorded and written in the output. Based on this output, several figures and plots are calculated by the macro benchmarkPET.C.
To run the macro, type:
The plots are saved in the file benchmarkPET.gif and the figures are printed in the console.
Seven plots are produced by the macro. They consist of 1-D or 2-D histograms with one entry per true coincidence (unscattered or scattered). A reference file is provided with the benchmark (benchmarkPET_RootOutputReference.gif). It contains the plots obtained with a correct execution of the benchmark (see Figure 3.2).
The three plots in the upper row show the detection coordinates of the annihilation photons in the detector heads. The first two plots ("Transaxial detection position") are 2-D histograms of the X and Y coordinates, one per frame, and illustrate the rotation of the heads. The third plot ("Axial detection position") is a 1-D histogram of the Z coordinate for the two frames, and illustrates the effect of the axial collimators at coordinates -100, 0, and +100. In addition to the detection position, the first two histograms also contain the emission coordinates of the positrons inside the line sources.
Two decay time curves are overlayed in the left plot of the middle row ("O-15 decay"). The red curve corresponds to the O-15 source and the blue curve to the F-18 source. The result of an exponential curve fit is superimposed to the O-15 curve. These curves illustrate the decay of the two sources.
The second plot in the middle row ("Axial sensitivity") illustrates the axial sensitivity of the scanner. For all unscattered coincidences, the axial position of the corresponding line-of-response is shown as an histogram.
The third plot in the middle row ("Axial scatter fraction") shows the scatter fraction as a function of the axial position.
The plot in the lower row ("Acollinearity Angle Distribution (deg)") shows the distribution of the deviation from 180 deg. of the angle between the two annihilation gammas. Unlike the other plots, this histogram is created in Gate and is part of the output. It should be noted that the selection of the annihilation photons is not based on the physical process that produced the photon (annihilation of the emitted positron with an electron), but rather on the energy of the photon (greater than 450 keV) and on the occurrence of the photon (should be the first two photons radiated from the positron). As a result, there is a very slight contamination by Bremsstrahlung X-rays. The last plot compares the delay coincidences obtained from the coincidence sorter with the real random coincidences sorted from the prompt coincidences with the rule enventID1 != eventID2.
Six physical variables are calculated by the macro from the output:
- the total number of generated events
- the total number of random coincidences
- the total number of true unscattered coincidences
- the total number of true scattered coincidences
- the simulated O-15 half life
- the full-width-half-maximum of the annihilation gamma acollinearity angle.
Here is an example of the macro output for one correct run of the benchmark:
Acquisition start time = 0 [s] Acquisition stop time = 240 [s] O-15 decay factor = 0.546383 F-18 decay factor = 0.987477 O-15 initial activity = 100000 [Bq] F-18 initial activity = 100000 [Bq] O-15 decays = 1.31132e+07 F-18 decays = 2.36994e+07 ==> Expected total number of decays during the acquisition is 3.68e+07 +/- 6067.34 There are 3.68187e+07 recorded decays There are 312864 true unscattered coincidences There are 23639 random coincidences There are 371679 scattered coincidences ==> there is a total of 709182 coincidences (true, scattered, and random) ==> global scatter fraction = 0.542959 ==> absolute sensitivity = 0.849742 Measured O-15 half-life = 122.476 [s] Nominal O-15 half-life = 122.24 [s] ==> difference = 0.192912 Gamma acollinearity FWHM = 0.584795 degree (expected: 0.58)
The PET benchmark has been run on twelve different system configurations. Two operating systems have been tested: Linux and Mac OS. The computing time for the PET benchmark varied between 4.5 and 10 hours (due to the clock processor).
The six PET physical variables characterizing the simulation results are shown in Table 3.2 with their mean value and standard error obtained from the twelve system configuration sample.
|variable type||mean value||𝜎|
|total decays during the acquisition||3.6813 x 10^7$||± 0.01 %|
|random coincidences||23536$||± 0.20 %|
|unscattered coincidences||312725$||± 0.12 %|
|scattered coincidences||370116$||± 0.10 %|
|simulated 15-Oxygen life time||121.85$ s||± 0.07 %|
|gamma acollinearity angle||0.58$ degree||0.01 %|
The results in Table 3.2 and, especially the standard deviations show that the main physical simulation variables are stable within less than 1 %. It is strongly recommended that the user validates a new or upgraded Gate installation using this table.
The SPECT benchmark consists in a macro corresponding to the simulation of a SPECT acquisition of a moving source described analytically (i.e., not voxelized).
Loading the SPECT benchmark package
The macro of the SPECT benchmark and the associated programs are available in the benchmarkSPECT directory:
SPECT benchmark: moving source
The gamma camera which is simulated does not correspond to any real system.
Description of the scanner
The simulated gamma camera is made of four heads at textstyle 90^o each (see Figure ), each with:
- a 2 cm thick lead shielding
- a lead collimator (hole diameter: 0.15 cm, collimator thickness: 3 cm, septal thickness: 0.6 mm)
- a 1 cm thick NaI crystal
- a 2.5 cm thick backcompartment in perspex, to roughly model the set of photomultipliers and associated electronics.
The head of the camera (i.e. the scanner) called SPECThead is 7 cm thick, 21 cm wide and 30 cm long.
Part of the macro describing the camera head:
/gate/world/daughters/name SPECThead /gate/world/daughters/insert box
The gamma camera scanner called SPECThead is created.
/gate/SPECThead/geometry/setXLength 7. cm /gate/SPECThead/geometry/setYLength 21. cm /gate/SPECThead/geometry/setZLength 30. cm
The X, Y and Z dimensions of the gamma camera head are set to 7 cm, 21 cm, and 30 cm respectively.
/gate/SPECThead/repeaters/insert ring /gate/SPECThead/ring/setRepeatNumber 4
The gamma camera head is repeated four times around the Z axis.
Description of the phantom
The object under investigation is a cylindrical phantom called Phantom (5 cm in diameter and 20 cm long) filled with water and including a cylinder called movsource, which is 2 cm in diameter and 5 cm long and which is filled with Tc-99m (emission energy: 140 keV) (see Figure 4.2).
During the acquisition, the phantom moves parallel to the Z axis at 0.04 cm/s.
Part of the macro describing the phantom:
/gate/world/daughters/name Phantom /gate/world/daughters/insert cylinder
The phantom is created.
/gate/Phantom/geometry/setRmax 5. cm /gate/Phantom/geometry/setRmin 0. cm /gate/Phantom/geometry/setHeight 20. cm /gate/Phantom/placement/setTranslation 0. 0. -6. cm
The dimensions and the position of the phantom are defined.
/gate/Phantom/moves/insert translation /gate/Phantom/translation/setSpeed 0 0 0.04 cm/s
The phantom will move at 0.04 cm/s along the Z axis.
Description of the source
The source movement is achieved using confinement (see Chapter VIII of the GATE Users Guide).
The movsource volume, which has the same dimensions as the source, is required to model the source motion.
The activity volume has the same diameter as the movsourse volume and is as long as the displacement of the movsource volume.
The activity is confined to the movsource volume. The activity of the source is set to 30 kBq.
Parts of the macro describing the confinement:
/gate/Phantom/daughters/name movsource /gate/Phantom/daughters/insert cylinder /gate/movsource/geometry/setRmax 2. cm /gate/movsource/geometry/setRmin 0. cm /gate/movsource/geometry/setHeight 5. cm /gate/movsource/placement/setTranslation 0. 0. -6. cm
The movsourse volume in which events will be confined is placed in the phantom and will have the same movement as its mother, namely the Phantom volume.
/gate/source/SourceConfinement/gps/type Volume /gate/source/SourceConfinement/gps/shape Cylinder /gate/source/SourceConfinement/gps/radius 2. cm /gate/source/SourceConfinement/gps/halfz 14.5 cm /gate/source/SourceConfinement/gps/centre 0. 0. 0. cm
The shape, the dimensions and the position of the source are defined.
The activity is confined to the movsource volume.
Description of the acquisition geometry
The table (0.6 cm deep, 8 cm wide and 34 cm long) is simulated. During the acquisition, the table moves together with the phantom (which does not move with respect to the table), at 0.04 cm/s. (see Figure 4.3 ).
Simulation parameters and data processing
Detected events are recorded in the crystal volume. Compton events are detected in the phantom (phantom and movsource), in the collimator, in the backcompartment, in the shielding and in the table.
The physics processes are modeled using the low energy electromagnetic processes package. Rayleigh, photoelectric and Compton interactions are turned on while the gamma conversion interactions are turned off. To speed up the simulation, thresholds are introduced. The X-rays are tracked until their energy falls under 20 keV. Secondary electrons are not tracked.
The digitizer processes the hits recorded by the sensitive detector. The adder module transforms hits into pulses. The blurring module simulates a Gaussian energy blurring of FWHM=10% at 140 keV.
The limited spatial resolution of the photomultipliers and associated electronic was simulated using a Gaussian blurring with a standard deviation of 2 mm. A thresholder and an upholder were used to consider only the particles detected with an energy between 20 and 190 keV.
Part of the macro describing signal processing:
/gate/digitizer/Singles/insert adder /gate/digitizer/Singles/insert blurring /gate/digitizer/Singles/blurring/setResolution 0.10 /gate/digitizer/Singles/blurring/setEnergyOfReference 140. keV
The blurring module simulates a Gaussian energy blurring of FWHM=10% at 140 keV.
Description of the acquisition protocol
A SPECT acquisition is simulated, using the SPECThead system. The acquisition consists of 64 projections, i.e. 16 projections per head. The heads move along a circular orbit with a 7 cm radius and a speed of 0.15 degree per second. Sixteen runs of 37.5 seconds are performed, to simulate the 16 positions of the 4 gamma camera heads, yielding 64 projections. Figure 4.4 illustrates the simulated configuration and examples of photon trajectories.
Part of the macro describing some "acquisition" parameters:
/gate/application/setTimeSlice 37.5 s /gate/application/setTimeStart 0. s /gate/application/setTimeStop 600. s
The output files of the SPECT benchmark can be analyzed using C programs from the ASCII output file or using ROOT from the ROOT output file. The type of ouput files to be created can be selected with the following commands:
/gate/output/root/disable /gate/output/ascii/setFileName benchSPECT /gate/output/ascii/setOutFileHitsFlag 0 /gate/output/ascii/setOutFileSinglesFlag 1
The ASCII file benchSPECTSingles.dat will be created
/gate/output/ascii/disable /gate/output/root/setFileName benchSPECT /gate/output/root/setRootSinglesFlag 1 /gate/output/root/setRootNtupleFlag 1
The ROOT file benchSPECT.root will be created.
To get both ASCII and ROOT output files, the macro should contain:
/gate/output/root/setRootSinglesFlag 1 /gate/output/root/setRootNtupleFlag 1 /gate/output/ascii/setOutFileHitsFlag 0 /gate/output/ascii/setOutFileSinglesFlag 1
By default the name of the ASCII file and ROOT output files are gateSingles.dat and gate.root respectively. To change these default names, use the command setFileName.
Two C programs are provided for processing the ASCII files produced by the macro:
- benchmark_projections.c creates Interfile image files corresponding to the total, scatter and primary photons detected in a 20-190 keV energy range.
- benchmark_spectra.c creates different spectra. These two C programs have first to be compiled.
To run the programs, use:
where GateSingles.dat is the name of the file generated when running the benchmark macro.
Help regarding the programs can be obtained using: Executable_name -h
The program creates raw image files from the ASCII file GateSingles.dat. Three images can be created:
- Total: image of primary and scattered photons.
- Scatter: image of scattered photons.
- Primary: image of primary photons. The program is an interactive program and runs as follows:
The user is asked for the dimension of the image to be created:
Enter the number of rows and colums of the images to be created 128
Image size is 128 x 128.
The user is asked for the lower energy threshold of the images to be created:
Enter the lower energy threshold (MeV) 0.02
Lower energy threshold is 0.020000 MeV.
The user is asked for the upper energy threshold of the images to be created:
Enter the upper energy threshold (MeV) 0.19
Upper energy threshold is 0.190000 MeV.
The user is asked for the pixel size of the images to be created:
Enter the size of the pixels (mm) 0.904
Pixels are 0.904 mm x 0.904 mm.
The number of projections is predefined and is 64.
Beginning of processing
End of processing
The program lists the number of detected counts in each image:
- Number of events in the image of primary and scattered photons = 17606
- Number of events in the image of primary photons = 5689
- Number of events in the image of scattered photons = 11917
The program lists the names of the file corresponding to the images that have been created:
- File "Total" corresponds to the image of primary and scattered photons.
- File "Primary" corresponds to the image of primary photons.
- File "Scatter" corresponds to the image of scattered photons.
Output files are Interfile files without header, that is raw data, in short integer coded on 2 bytes.
The program reads the ASCII file GateSingles.dat and creates spectra.
The program is an interactive program and runs as follows:
The user is asked for the lower threshold of the energy spectrum to be created:
Enter the lower energy threshold (MeV) 0.02
The lower energy threshold is 0.020000 MeV.
The user is asked for the upper threshold of the energy spectrum to be created:
Enter the upper energy threshold (MeV) 0.19
The upper energy threshold is 0.190000 MeV.
The user is asked for the number of energy windows in which the energy range should be split:
Enter the number of energy windows 100
The energy range will be divided into 100 energy windows.
The user is asked for the number of scattering volumes for which the scatter spectrum should be created. The maximum number of scattering volumes is 10.
Enter the number of scattering volumes to be considered (max is 10) 6
6 scattering volumes will be considered.
The user is asked for the first four letters of each scattering volume to be considered:
Enter the first four letters of the scattering volume number 1 shie
The user chooses shie for shielding.
Enter the first four letters of the scattering volume number 2 coll
The user chooses coll for collimator.
Enter the first four letters of the scattering volume number 3 comp
The user chooses comp for compartment.
Enter the first four letters of the scattering volume number 4 tabl
The user chooses tabl for table.
Enter the first four letters of the scattering volume number 5 movs
The user chooses movs for movsource.
Enter the first four letters of the scattering volume number 6 Phan
The user chooses Phan for Phantom.
The scattering volumes that will be considered are: shie coll comp tabl movs Phan.
Beginning of processing
End of processing
The program lists the names of the files corresponding to each spectrum that has been created:
- File fich_p corresponds to the primary spectrum.
- File fich_c corresponds to the photons which have scattered in the crystal.
- File fich_t corresponds to the total spectrum.
- File fich_1 corresponds to the photons whose last scatter interaction was in the shiel volume.
- File fich_2 corresponds to the photons whose last scatter interaction was in the coll volume.
- File fich_3 corresponds to the photons whose last scatter interaction was in the comp volume.
- File fich_4 corresponds to the photons whose last scatter interaction was in the tabl volume.
- File fich_5 corresponds to the photons whose last scatter interaction was in the movs volume.
- File fich_6 corresponds to the photons whose last scatter interaction was in the Phan volume.
- File fich_o1 corresponds to 1st order scatter File fich_o2 corresponds to 2nd scatter order.
- File fich_o3 corresponds to 3rd scatter order File fich_o4 corresponds to 4th scatter order.
- File fich_o5 corresponds to scatter orders greater than 4.
The output files are ASCII files. Each file is made of one column and as many lines as the number of energy windows. A line represents the number of detected photons in an energy window. The generated files can be easily imported in Excel to be processed. As an example, the total spectrum and the scatter spectra as a function of the scatter order of the detected events are shown in the Figure below.
Using ROOT, an analysis program can be applied to the output file benchSPECT.root.
This command runs the ROOT software and loads the file in the program.
In order to execute the same analysis as the benchmark_spectra.c program, the user should type the command:
root  .x benchSPECT.C
Results: SPECT benchmark
The results presented in that section are the compilation of the results obtained from eleven different computer installations.
The number of emitted particles is 17999920±1108 and the total number of detected counts from 20 to 190 keV is 35758±57.
The percentages of primary and scatter photons, as a function of the medium in which the last scatter event occurred, are shown in Table 4.1.
|35.6 ± 0.1 %||53.3 ± 0.1 %||3.0 ± 0.1 %||0.34 ± 2.9 %||6.4 ± 0.2 %||1.2 ± 3.3 %|
The percentages of scatter photons as a function of the scatter order are shown in Table 4.2
|Scatter order||Order 1||Order 2||Order 3||Order 4||Order > 4|
|Number of detected counts (%)||48.2 ± 0.2||26.4 ± 0.3||12.8 ± 0.15||6.6 ± 0.1||6.0 ± 0.4|
- Detected counts per projection
The minimum, maximum, mean and associated standard deviation values of detected counts per projection are (see Figure ):
Min value = 2±1
Max value = 589±10
The Figure below shows a four-peak pattern. Each peak presents a bell-shape as the source moves from one side of the gamma camera field of view to another: the number of detected counts is the highest when the source is at the center of the axial field of view, and is the lowest when the source is at one edge. As the gamma camera has 4 heads at 90° each and because of the symmetry of the object under consideration, about the same number of counts is detected by the 4 heads at each position of the source, hence the four-peak pattern.
- Count profile
The 64 projections of the image of primary and scattered photons were summed and all rows of the resulting image were summed, yielding the count profile shown in the Figure below. The minimum, maximum, mean and associated standard deviation values of the count profile:
Min value = 6±2
Max value = 380±2
Mean = 120.6±1.6
To get the simulation time automatically in the output of a Gate run, the macro should contain the following command:
Table 4.3 summarizes the PC configurations on which the benchmark was run, the OS, compiler and software versions used and the corresponding computation time of the SPECT
|PC Configuration||OS Version||Compiler||Geant4||CLHEP||Computation time|
|Dual processor 1 GHz||RedHat 9.0||gcc 3.2.2||4.5.2 p02||188.8.131.52||11H06|
|Monoprocessor 2.6 GHz||RedHat 9.0||gcc 3.2.2||4.5.2||184.108.40.206||7H51|
|Dual processor 0.8 GHz||SuSE 8.1||gcc 3.2||4.5.2||220.127.116.11||6H42|
|Monoprocessor 2.6 GHz||SuSE 9.0||gcc 3.3.1||4.5.2 p02||18.104.22.168||5H46|
|Monoprocessor 2.0 GHz||RedHat 8.0||gcc 3.2||4.4.1||22.214.171.124||5H21|
|Monoprocessor 1.0 GHz||Mac OS X 10.2.8||gcc 3.1||4.5.2 p02||126.96.36.199||5H04|
|Monoprocessor 2.4 GHz||RedHat 9.0||gcc 3.2||4.5.1||188.8.131.52||3.10|
|Monoprocessor 2.4 GHz||RedHat 7.1||gcc 3.2.2||4.5.2 p02||184.108.40.206||4H16|
|Monoprocessor 1.8 GHz||RedHat 7.3||gcc 2.95||4.5.1||220.127.116.11||3.05.07|
|Monoprocessor 2.8 GHz||Fedora Core 1||gcc 3.3.2||4.5.2 p02||18.104.22.168||2H55|
|Monoprocessor xeon 5130 @ 2.0 GHz||Fedora Core 4||gcc 4.1||4 9.0||22.214.171.124||1H05|