D690 Simulations will not unlist

[robbietuk]

WARNING: CListModeDataROOT: I've set the scanner from STIR settings and ignored values in the hroot header.

ERROR: the number of transaxial blocks per bucket, the number of axial crystals per block, the number of axial crystals per singles unit, the number of transaxial crystals per singles unit, 
libc++abi.dylib: terminating with uncaught exception of type std::__1::basic_string<char, std::__1::char_traits<char>, std::__1::allocator<char> >
sub_scripts/unlist.sh: line 28: 29470 Abort trap: 6           lm_to_projdata lm_to_projdata_${SGE_TASK_ID}.par

Due to mismatch between STIR and root file (unsure if .hroot or .root, probbably .root). https://github.com/UCL/STIR/blob/c11e9673cad3702aa83b288edcc27527e25fea7e/src/listmode_buildblock/CListModeDataROOT.cxx#L260-L264

mMR unlisting works fine.

[robbietuk]

@KrisThielemans This STIR error does not seem to be outputting what it should.

https://github.com/UCL/STIR/blob/c11e9673cad3702aa83b288edcc27527e25fea7e/src/listmode_buildblock/CListModeDataROOT.cxx#L239

This message is not given by STIR, actually made debugging a little harder.

[robbietuk]

The above error was resolved by changing originating system := GE Discovery 690 to originating system := User_defined_scanner in the root_header_template.hroot.

However it seems there is a discrepancy between the GATE and STIR D690 Scanners, see console output:

List mode file

WARNING: LmToProjData:
Scanner from list mode data (ROOT_defined_scanner) is different from
scanner from template projdata (GE Discovery 690).
Full definition of scanner from list mode data:
Scanner parameters:= 
Scanner type := ROOT_defined_scanner
Number of rings                          := 24
Number of detectors per ring             := 576
Inner ring diameter (cm)                 := 81.02
Average depth of interaction (cm)        := 0.94
Distance between rings (cm)              := 0.654
Default bin size (cm)                    := 0.21306
View offset (degrees)                    := 0
Maximum number of non-arc-corrected bins := 381
Default number of arc-corrected bins     := 381
Number of blocks per bucket in transaxial direction         := 1
Number of blocks per bucket in axial direction              := 4
Number of crystals per block in axial direction             := 24
Number of crystals per block in transaxial direction        := 9
Number of detector layers                                   := 1
Number of crystals per singles unit in axial direction      := 24
Number of crystals per singles unit in transaxial direction := 9
end scanner parameters:=

STIR SCANNER

Full definition of scanner from template:
Scanner parameters:= 
Scanner type := GE Discovery 690
Number of rings                          := 24
Number of detectors per ring             := 576
Inner ring diameter (cm)                 := 81.02
Average depth of interaction (cm)        := 0.94
Distance between rings (cm)              := 0.654
Default bin size (cm)                    := 0.21306
View offset (degrees)                    := -5.021
Maximum number of non-arc-corrected bins := 381
Default number of arc-corrected bins     := 331
Number of blocks per bucket in transaxial direction         := 2
Number of blocks per bucket in axial direction              := 4
Number of crystals per block in axial direction             := 6
Number of crystals per block in transaxial direction        := 9
Number of detector layers                                   := 1
Number of crystals per singles unit in axial direction      := 1
Number of crystals per singles unit in transaxial direction := 1
end scanner parameters:=

The notable difference are:

1. View offset (degrees)
2. Default number of arc-corrected bin
3. Number of blocks per bucket in transaxial direction
4. Number of crystals per block in axial direction
5. Number of crystals per singles unit in axial direction
6. Number of crystals per singles unit in transaxial direction

While 1,2 probably don't matter as much, 3,,4,5, 6 are scanner properties. Is there a mistake in the GATE D690 Geometry?

[robbietuk]

Checking my old GATE simulations and unlistings, the geometries are the same and I think Elise and I used a modified STIR scanner to match the ROOT file for all unlisting/ reconstructions.

Simulating mMR and D690 scans (same STIR generated phantom) over 10s (GATE time). mMR: Coincidences in root file = 68541, unlisted= 22984 D690: Coincidences in root file = 5000, unlisted = 3922.

[KrisThielemans]

certainly the STIR info on blocks is correct, while the .hroot isn't.

you can check this from the pettoolbox

https://github.com/UCL/GEPETToolbox/blob/master/pettoolbox/recon/matlab/D690scanner.m calls https://github.com/UCL/GEPETToolbox/blob/master/pettoolbox/recon/matlab/XRscanner.m which contains a lot of this

https://github.com/UCL/GEPETToolbox/blob/dc812cc847506b0f7d27e66fb9d5e33e70c5b24b/pettoolbox/recon/matlab/XRscanner.m#L51-L52 confirms that Default number of arc-corrected bins is correct (nUrr, obviously!).This will only matter for FBP. I'd recommend removing it from the .hroot file, although then it might fail.

View offset (degrees) is tricky. STIR ignores it (until https://github.com/UCL/STIR/pull/181 will be merged). However, the current GATE geometry will be defining it with 0 I guess (otherwise, all images would have been rotated), so this is a case where the .hroot file would have to differentate from the stir:Scanner, until we fix the GATE geometry (after merging https://github.com/UCL/STIR/pull/181)

[KrisThielemans]

differences in number of coincidences: it is in the ROOT file itself? the mMR does have long axial FOV hence higher sensitivity, but not a factor 11...

[robbietuk]

https://github.com/UCL/GEPETToolbox/blob/dc812cc847506b0f7d27e66fb9d5e33e70c5b24b/pettoolbox/recon/matlab/XRscanner.m#L51-L52 confirms that Default number of arc-corrected bins is correct (nUrr, obviously!).This will only matter for FBP. I'd recommend removing it from the .hroot file, although then it might fail.

This can and has now been removed from the .hroot and it will still unlist.

differences in number of coincidences: it is in the ROOT file itself? the mMR does have long axial FOV hence higher sensitivity, but not a factor 11...

Yes, I got those "Coincidences in root file" values were given by root tools. I wonder if the geometry.mac for the D690 is wrong, most likely explaination.

[robbietuk]

Gate D690 geometry (PETTOOLBOX) differences:

Crysrtal Length: 25mm (30mm) https://github.com/UCL/GEPETToolbox/blob/dc812cc847506b0f7d27e66fb9d5e33e70c5b24b/pettoolbox/recon/matlab/XRscanner.m#L26

Crystal sizes (computed). https://github.com/UCL/GEPETToolbox/blob/dc812cc847506b0f7d27e66fb9d5e33e70c5b24b/pettoolbox/recon/matlab/XRscanner.m#L27-L34 However, this might be to remove gaps as values are close and slightly larger (trying to close gaps).

Besides these slight differences, no major discrepancies. I am fairly sure that the D690 geometry is acceptable.

Now I will check positioning. It would be useful to have my visualisation tools working...

[KrisThielemans]

crystal sizes: yes, as we don't have the block-projector yet, I'm sure that's what Elise did, such that the STIR and GATE sims match a bit better.

D690 geo macro looks ok at first sight. Worth putting some comments in there about offset and crystal sizes, such that next time we remember.

[KrisThielemans]

positioning: I see that the mMR uses ECAT while the D690 uses Cylindrical. Maybe a difference there (although you'd hope not).

I'm quite confused though that both seem to define a scanner from z=0..L, not centred around 0. That doesn't seem to fit what we though for the phantom placement.

Not sure if vgate would allow you to visualise. Or via thinlinc.

[robbietuk]

crystal sizes: yes, as we don't have the block-projector yet, I'm sure that's what Elise did, such that the STIR and GATE sims match a bit better.

D690 geo macro looks ok at first sight. Worth putting some comments in there about offset and crystal sizes, such that next time we remember.

Done.

I managed to get my visualisation to work. The big green cubiod is suppose to be a phantom (attenuation) XCAT torso, I hope I don't look like that on the inside... Changing the phantom to a voxalised phantom cyclinder (generated with STIR) does not resolve the issue as the green cubiod is the same size. It is of a constant size, as shown by switching scanner to the mMR scanner with a smaller radius. I am still working on this.

[KrisThielemans]

it's the sideways view that's important.

could also do a set of simple geometric sphere sims, where you step the sphere along the axis. should tell you where the middle is...

[robbietuk]

Checking my old GATE simulations and unlistings, the geometries are the same and I think Elise and I used a modified STIR scanner to match the ROOT file for all unlisting/ reconstructions.

Simulating mMR and D690 scans (same STIR generated phantom) over 10s (GATE time). mMR: Coincidences in root file = 68541, unlisted= 22984 D690: Coincidences in root file = 5000, unlisted = 3922.

mMR has a significantly lower threshold than D690 in digitiser.mac

mMR https://github.com/UCL/STIR-GATE-Connection/blob/962320d202a23045a49d0c8f73555b6164f5219c/ExamplesOfScannersMacros/mMR/digitiser_mMR.mac#L24-L25

D690 https://github.com/UCL/STIR-GATE-Connection/blob/962320d202a23045a49d0c8f73555b6164f5219c/ExamplesOfScannersMacros/D690/digitiser_D690.mac#L20

This was probably put there by @ludovicabrusaferri for her lower energy unlisting. Switching back to 430. keV results in

mMR: Coincidences in root file = 22000 D690: Coincidences in root file = 5000 (unchanged) Getting closer. I dont know the mMR true lower threshold, I would assume it is ~430 keV.

[robbietuk]

it's the sideways view that's important.

Fixed with changing forceSolid -> forceWireframe and using OGLI. The example simple cylinder phantom (attenuation) is approximately at the center of both D690 and mMR scanners, there is a slight offset in STIR parameter file that I have played with.

shape type:= ellipsoidal cylinder
Ellipsoidal Cylinder Parameters:=
   radius-x (in mm):=100
   radius-y (in mm):=100
   length-z (in mm):=110
   origin (in mm):={70,10,20}
   END:=
value := 0.096

This does not seem to be a source of error here, although visualising source is a little tougher in GATE.

[KrisThielemans]

no idea what I'm looking at here... I feel we don't want to add the (confusing) STIR coordinate system for generate_image in here. I'd probably use a long cylinder such that it fails the whole image. There's then an easy test that the sinograms counts have to be symmetric.

Maybe this isn't necessary now, but we'll have to think about some tests. Might be best to do that with some others.

[robbietuk]

So the purpose of that image and those tests was to verify that the object (in the case a cylinder in white) actually existed at approximately the same position in both the scanners (yes) and could be manipulated by adjusting the offsets and sizes in STIR generate_image parameter file (yes).

I have not yet done any fine tuning tests to see where they if voxels are offset or are actually centeralised.

This was done to investigate the discrepancy between the two scanners number of counts in the root file, a factor of approximately 4 now.

mMR: Coincidences in root file = 22000 D690: Coincidences in root file = 5000 (unchanged)

I'm still not sure what is causing this larger than expecting discrepancy.

differences in number of coincidences: it is in the ROOT file itself? the mMR does have long axial FOV hence higher sensitivity, but not a factor 11...

The mMR has a smaller diameter and a larger axial FOV than the D690, maybe a factor 4 isnt unreasonable... I will check.

The big issue now is the difference between the scanner definition given by the root file and STIR's definition: https://github.com/UCL/STIR-GATE-Connection/issues/14#issuecomment-636795629

[KrisThielemans]

Sensitivity measures from literature

• mMR: 13.3 cps/kBq from https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4766138/
• D690: 7.5cps/kBq from http://dx.doi.org/10.1118/1.3635220

Looks like both were determined in the same way (long line source in the centre and 10cm off-centre). So I'm not so sure. That's of course with the energy window settings for each system etc. Another factor is that sensitivity is measured without randoms, and total ROOT counts not. Still a bit surprising of course.

[KrisThielemans]

The big issue now is the difference between the scanner definition given by the root file and STIR's definition: #14 (comment)

I thought we agreed that STIR is correct, GATE macro is correct, but hroot is wrong (with the caveat of the view offset has to be 0, which is different in the STIR file)

[robbietuk]

The first comment in this issue: https://github.com/UCL/STIR-GATE-Connection/issues/14#issue-628310168

ERROR: the number of transaxial blocks per bucket, the number of axial crystals per block, the number of axial crystals per singles unit, the number of transaxial crystals per singles unit, 
libc++abi.dylib: terminating with uncaught exception of type std::__1::basic_string<char, std::__1::char_traits<char>, std::__1::allocator<char> >
sub_scripts/unlist.sh: line 28: 29470 Abort trap: 6           lm_to_projdata lm_to_projdata_${SGE_TASK_ID}.par

was caused by the .hroot scanner definition https://github.com/UCL/STIR-GATE-Connection/blob/962320d202a23045a49d0c8f73555b6164f5219c/ExamplesOfScannersMacros/D690/unlisting/root_header.hroot#L24-L33

STIR takes these values and does some small processing https://github.com/UCL/STIR/blob/7dc238f5ec228554f594f1f2228450eb4618472b/src/include/stir/IO/InputStreamFromROOTFileForCylindricalPET.inl#L20-L107 checking them against the STIR scanner https://github.com/UCL/STIR/blob/7dc238f5ec228554f594f1f2228450eb4618472b/src/listmode_buildblock/CListModeDataROOT.cxx#L237-L297

I think this is impossible to set these .hroot values to make STIR happy with the D690 scanner. Here is an example of why in the Transaxial dimention.

• STIR scanner - Transaxial blocks per bucket = 2. STIR interprets this as = module_repeater_y, therefore number of modules_Y = 2 in the .hroot.
• STIR scanner - Transaxial crystals per block = 9 but is computed as crystal_repeater_y * module_repeater_y, therefore number of crystals_Y = 4.5? in the .hroot.
• STIR scanner - Transaxial crystals per singles unit = 1. For a Singles readout depth = 1 (what we are using), this is computed from the .hroot as module_repeater_y * submodule_repeater_y * crystal_repeater_y.

There are similar inconsitancies betweeen Axial crystals per block = 6 and Axial crystals per singles unit = 1. The number all make sense but I think it is the math in https://github.com/UCL/STIR/blob/7dc238f5ec228554f594f1f2228450eb4618472b/src/include/stir/IO/InputStreamFromROOTFileForCylindricalPET.inl#L20-L107 that is leading to this issue.

[KrisThielemans]

I thought there was a STIR issue on the interpretation of hroot vs macro vs STIR, but cannot find it. there is https://github.com/UCL/STIR/issues/96. This seems to need a new issue (then cross-reference in https://github.com/UCL/STIR/issues/96). Best to discuss there (not here).

at first glance, Transaxial crystals per block = crystal_repeater_y. The formula for Transaxial crystals per singles unit seems wrong as well (not that I understand what the singles readout depth means. would be best to add a link to a GATE page)

[KrisThielemans]

@robbietuk Regarding positioning, I found the following on ECAT at https://opengate.readthedocs.io/en/latest/defining_a_system_scanner_ct_pet_spect_optical.html#id35

For a multiple block-ring system centered axially on the base ecat, the axial position of this first block should be set to zero (see Fig. 1.20):

implying that for ECAT, the scanner is centred at 0. There's no such statement for CylindricalPET.

[robbietuk]

There was an issue with the world size being ill defined that meant that the scanner was outside of the world in places. https://github.com/UCL/STIR-GATE-Connection/blob/962320d202a23045a49d0c8f73555b6164f5219c/VoxelisedXCAT/main_muMap_job.mac#L12-L15 To fix this, the world will now be defined in the respective scanner geometry.mac files.

There is still an issue with number of counts detected that I am working on.

[robbietuk]

15 fixed a lot. Closing becase #18 describes the current state of this issue better.