Scatter simulation #44

NikEfth · 2016-09-12T23:16:06Z

So, this is the PR for the scatter simulation.

New classes ScatterSimulation and SingleScatterSimulation which have taken over the old ScatterEstimation. This is their "minimal" implementation having only the functionality that the older class used to have. I removed, for the time being, the subsampling stuff, in order to first see why/where there is the mistake. Once this has been sorted I will start adding more things bit-by-bit.

Please don't review the new ScatterEstimation class, yet. I will create a branch, from this one, on which I will work on the ScatterEstimation and let you know when it will be ready for reviewing.

In addition, there are some minor changes in other places, like clone functions in ExamInfo.
In my laptop all tests succeeded and Travis is doing good, so far.

ATB,
Nikos

Transfer code from the other branches.

Files which were deleted when merged with UCL/master

the run_scatter_tests.sh.

New clone functions added in ExamInfo.

KrisThielemans · 2016-09-13T13:48:08Z

Hi
Most of these are small. some are pedantic...

Some of these functions should probably have PET in their name (Maybe only the ScatterSimulation bits?)

Not sure about authors of these files (and corresponding copyright). Probably need to include Harry, maybe Pablo and Nikos Dikaios. Depends on what you copied from where.
Copyright (C) 2013 - 2016 University College London, use a comma unless changes made 2014, 2015 as well :-;

ScatterEstimationByBin

Best to rename to avoid name clash with old class. Drop ByBin?
Why do we need reconstruct_iterative and reconstruct_analytic (with such differrent calling sequence)?
get_image_from_file doxygen says it’s set_image_from_file. Also some doxygen below that for functions which don’t exist, which will confuse it.
post_processing() I think we should move most of its functionality to setup(). post_processing is only called after parsing, so wouldn’t be appropriate within python/Matlab.
ProjDataInfoCylindricalNoArcCorr * proj_data_info_2d_ptr; No raw pointers as class members.
ExamInfo::clone() . no need for the static cast (just return new ExamInfo(*this))
There are multiplicative_data_3d and normalisation_data_3d (what’s the difference?) and norm_projdata_3d_sptr. The latter one should normally not be a class member as it’s presumably used to create normalisation_data_3d
mask_atten_image_filename etc. Why the “atten” in the variable name?
variable multimulti2d is a strange name :-;

ScatterSimulation

doxygen containis typo SignelScatterSimulation (important for linking). Also needs to be prefixed by stir:: to get the link to work.
Lots of inline functions. I think we should move it to .cxx. Inline is good for small functions called often, but otherwise slow down the compiler (and mean more recompilations if you change their implementation). (could leave the Compton ones inline if you prefer)
set_activity_image_sptr return Succeeded, set_activity_image doesn’t (calls error). Best to let both do the same,and I’m these days in favour of the latter. (i.e. let the user catch exceptions)
get_output_proj_data(shared_ptr& arg), change to shared_ptr get_output_proj_data() const;
set_attenuation_threshold(const float& arg) etc. ok to use references but it’s quite non-standard and actually slower for such small entities (as you’ll be passing pointers around)
post_processing creates a ProjDataInMemory for input_projdata_2d_sptr if data is 2D, just copy pointer?
reconstruction_template_par_filename usage. I see you had to do some explicit KeyParser stuff here after all. Maybe that means we had to have it elsewhere after all? I guess same for ScatterSimulation
Zoom factors currently still hard-wired
reconstruct_iterative. Cannot have paths like "./extras/recon_". They will likely fail on Windows. Shouldn’t the writing be left to the reconstruction object? (i.e. set its output filename)
apply_mask_in_place. This seems far too complicated. It was done like that in the script as we didn’t have a masking facility there. If the mask is binary, shouldn’t we just multiply?

other

scatter_detection_modelling.cxx, scatter_estimate_for_one_scatter_point.cxx, single_scatter_integrals.cxx, upsample_and_fit_scatter_estimate.cxx sadly have a load of white space differences. Not sure what to do about this, but I’d prefer to revert those…
single_scatter_estimate.cxx . has white space diffs, but we probably should rename single_scatter_estimate etc?

- [x] Choose set_up by dynamic_cast rather than name. - [x] post_processing() I think we should move most of its functionality to setup(). post_processing is only called after parsing, so wouldn’t be appropriate within python/Matlab. - [x] reconstruction_template_par_filename usage. I see you had to do some explicit KeyParser stuff here after all. Maybe that means we had to have it elsewhere after all? I guess same for ScatterSimulation

SingleScatterSimulation ScatterEstimationByBin: - [x] post_processing() I think we should move most of its functionality to setup(). post_processing is only called after parsing, so wouldn’t be appropriate within python/Matlab. - [x] ProjDataInfoCylindricalNoArcCorr * proj_data_info_2d_ptr; No raw pointers as class members. ExamInfo::clone() . no need for the static cast (just return new ExamInfo(*this)) - [x] Choose set_up by dynamic_cast rather than name. ScatterSimulation: - [x] Lots of inline functions. I think we should move it to .cxx. Inline is good for small functions called often, but otherwise slow down the compiler (and mean more recompilations if you change their implementation). (could leave the Compton ones inline if you prefer) - [x] set_activity_image_sptr return Succeeded, set_activity_image doesn’t (calls error). Best to let both do the same,and I’m these days in favour of the latter. (i.e. let the user catch exceptions) - [x] set_attenuation_threshold(const float& arg) etc. ok to use references but it’s quite non-standard and actually slower for such small entities (as you’ll be passing pointers around) - [x] post_processing creates a ProjDataInMemory for input_projdata_2d_sptr if data is 2D, just copy pointer? — I presume that this comment is for ScatterEstimationaByBin. Well not quite. If it is 2D then I cannot upsample. So (better) through an error? - [x] reconstruction_template_par_filename usage. I see you had to do some explicit KeyParser stuff here after all. Maybe that means we had to have it elsewhere after all? I guess same for ScatterSimulation - Most of this code is now in the post_processing().

Next steps are to run the estimation process.

What is left to be tested is the return to 3D. - [x] There are multiplicative_data_3d and normalisation_data_3d (what’s the difference?) and norm_projdata_3d_sptr. The latter one should normally not be a class member as it’s presumably used to create normalisation_data_3d

the results are disappointing.

KrisThielemans

hi. most comments are tiny. I think there is one bug (in attenuation coefficients).

there's still a few files with mostly white-space diffs. I'd prefer that you revert those such that I can see what actually changed (unless you're 100% sure that you current editor settings are the final ones that we're going to use, and you email me a diff ignoring white-space...)

thanks!

examples/Siemens-mMR/unlist_and_randoms.sh

recon_test_pack/run_scatter_tests.sh

src/include/stir/scatter/ScatterEstimationByBin.h

src/scatter_buildblock/ScatterEstimationByBin.cxx

src/scatter_buildblock/upsample_and_fit_scatter_estimate.cxx

src/scatter_buildblock/single_scatter_integrals.cxx

src/scatter_buildblock/scatter_registries.cxx

src/scatter_buildblock/single_scatter_estimate.cxx

In addition to the comments in the github side: - [x] mask_atten_image_filename etc. Why the “atten” in the variable name? - [x] variable multimulti2d is a strange name :-; - [x] get_output_proj_data(shared_ptr& arg), change to shared_ptr get_output_proj_data() const; - [x] Zoom factors currently still hard-wired - [x] Every member initialised in post_proccessing must have a set function - [x] Use only BinNormalisation throught out the code.

NikEfth

Needs some clean up but most comments have been implemented.

NikEfth · 2016-09-19T17:01:56Z

src/scatter_buildblock/ScatterEstimationByBin.cxx

+    this->multiplicative_data_2d_sptr.reset(
+                new ChainedBinNormalisation(attenuation_correction, normalisation_coeffs));
+
+    iterative_object->get_objective_function_sptr()->set_normalisation_sptr(multiplicative_data_2d_sptr);


I cannot fully understand the comment.

NikEfth · 2016-09-21T11:23:19Z

src/scatter_buildblock/ScatterEstimationByBin.cxx

+        int len_x = this->atten_image_sptr.get()[0][min_z][min_y].get_length();
+
+        if (len_y != len_x)
+            error(boost::format("The voxels in the x (%1%) and y (%2%)"


I have to check this more. I think I saw it somewhere

I had some mistakes mainly in the normalisations. I think the estimation it working now.

fitting to the input_data than (input_data - randoms). In addition, I removed some obsolete variables.

minor amendments and addition of the code to sabsample the scanner

Use boost::filesystem to create the path and filenames in debug mode. Now, paths are cross-platform more minor aesthetic corrections

… Correction made.

In estimation the appending to boost path was corrected in more places.

The run_scatter_test was a bit weird as in the simultion sets the zoomed image both in for the scatter points and as the original attenuation image. I am not sure wether this is bug or a feature to speed up the process. I left it as is for now.

KrisThielemans · 2018-01-23T14:37:10Z

hi nikos, I haven't checked why you need boost filesystem and boost system, but I'm somewhat reluctant to use this. We'd introduce a dependency on compiled libraries for boost, and that's creating us a ton of problems for SIRF sadly (you would think it'd be simple but it isn't).

could you tell me what you're using this for? We likely have equivalents in utilities.h.

We were using the environment variable RANDOM before, but this is a system variable that returns a random number. This caused failurs on OSX.

KrisThielemans · 2020-04-13T10:32:15Z

Success! https://travis-ci.org/github/UCL/STIR/builds/674304290.

The problem was that the test-script was using the env variable RANDOM. This turns out to be used to let you get a random number... (try echo $RANDOM a few times). Looks like you can overwrite it in Linux, but not OSX. Oh well.

Code was fine all along. I've also (re)checked that results of the scatter-simulation are identical to what we have in master.

I'll extend the tests a bit, possibly do a bit more renaming, and then it should be ready.

However, while checking everything, I discovered a problem when using the zooming of the attenuation image. I'll file a separate issue for that as it's also in master. Not sure yet if I'll fix this in this PR or not.

KrisThielemans · 2020-04-16T16:56:21Z

A few more thing to finish

fix (potential) axial shift in scatter simulations #495
remove used of static variables such as max_single_scatter_cos_angle that actually depend on energy etc (as unsafe when the Scatter Simulation is run for different energy windows)
rename some classes to include the PET name (as wouldn't work for SPECT etc)
add a simple test for ScatterEstimation

- Currently only check in Debug mode to avoid impact on computation time. - move computation of shift_detector_coordinates_to_origin to set_up (from process_data)

see UCL#495 - in downsample_density_image*, adjust zoom_z to avoid problems with UCL#495, also allow specifying default zoom_z (-1) if new_size_z is given - add functions to check consistency in z-location (currently still disabled) Also add get_*image() functions to avoid using shared_ptr too much

we had some static variables, and forget to clear the cached detection_points_vector

These are tests for a symmetric object, checking if the result is symmetric. Also enabled the consistency check mentioned in UCL#495. This means that run_scatter_tests/sh will now fail.

This is now done by default in ScatterSimulation. Removing this avoids problems with UCL#495. - changed recon_test_pack and PET_simulation script/par-file - changed scatter-output file in recon_test_pack to take the new zooming into account. Difference with respect to the previous version is about 4%. (Note that the previous version was incorrect due to UCL#495)

- some clarifications in sample scatter_simulation.par - removed obsolete file

hopefully a better name...

now consistent with the rest of STIR

Use scatter_estimate(Bin) as opposed to in terms of 2 detectors such that it will apply for SPECT at some point.

original values did not reflect this was PET specific. Also used more standard-form keywords for the simulate_scatter .par file. Also fixed some issues in example files and STIR-UsersGuide

This was not necessary earlier-on as SingleScatterSimulation::process_data() was calling set_up(), but no longer...

- use standard keywords for reconstruction and scatter simulation - allow specifying these in the parameter file itself (as opposed to requiring a separate file)

3 Travis jobs gave errors of 0.124% while our threshold was 0.1%.

KrisThielemans · 2020-04-21T07:02:16Z

Current results (not changed since a while now) on the NEMA data from examples/Siemens-mMR/process_NEMA_data.sh
(500s, OSEM 14 subsets, 42 subiterations, 8mm Gaussian post-filter)

Good, but some non-uniformity still visible. We don't have the bed attenuation currently though.

This commit finally closes UCL#495

…kip]

KrisThielemans · 2020-04-21T13:57:50Z

I'm going to stop with this. There's still improvements to be made, but time has run out. It's working ok as far as I can see, and tests succeeded.

I'll make an issue with extra things to do.

…3_05_07 Merge master to tof 2023 05 07 + Siemens support

NikEfth added 5 commits September 12, 2016 12:02

Working on ScatterSimulation.

5490bec

Transfer code from the other branches.

Restored missed files

96750e1

Files which were deleted when merged with UCL/master

Cleaning up ScatterSimulation, enchance simulate_scatter.cxx and updata

981156f

the run_scatter_tests.sh.

ScatterSimulation almost finished.

4af0b46

New clone functions added in ExamInfo.

Removed some unessesary includes and improved some doxygen strings.

0e43d5d

NikEfth added 7 commits September 13, 2016 22:11

The set_up() process is successfull.

d3b387c

Next steps are to run the estimation process.

Estimation is working,

5efe7db

the results are disappointing.

Forgot to uncomment one command

d560807

Corrected the parfiles path

2ac2e5e

KrisThielemans requested changes Sep 17, 2016

View reviewed changes

NikEfth commented Sep 21, 2016

View reviewed changes

NikEfth added 13 commits September 21, 2016 17:32

Fixed some bugs in the upsampling process.

1210468

I had some mistakes mainly in the normalisations. I think the estimation it working now.

Corrected the error during upsample_and_fit_scatter_estimate I was

2d71dbc

fitting to the input_data than (input_data - randoms). In addition, I removed some obsolete variables.

BUG corrected: I started reconstructing from subset 1.

e5facf7

Comment out the filtering of negative values when subtract two sinograms

ae80372

Merge remote-tracking branch 'upstream/master' into Scatter_simulation

31b83e5

Additions on Scatter simulation

9afe197

minor amendments and addition of the code to sabsample the scanner

Dropped the "byBin" suffix in the name ScatterEstimation name

5ac913e

Use boost::filesystem to create the path and filenames in debug mode. Now, paths are cross-platform more minor aesthetic corrections

Removed std::make_shared

6ce834b

typedef struct compilation error fixed

20d3fa5

filesytem append was introduced in boost 1.55 but we build on v.1.54.…

ee47ffd

… Correction made.

Further clean-up in simulation,

650ff8e

In estimation the appending to boost path was corrected in more places.

Travis installs suggested package boost-filesystem and boost-system

57ea910

KrisThielemans force-pushed the Scatter_simulation branch from 8fd93fd to e6a41de Compare April 13, 2020 00:40

recon_test_pack/ run_scatter_test.sh fixes

e04c3af

We were using the environment variable RANDOM before, but this is a system variable that returns a random number. This caused failurs on OSX.

KrisThielemans force-pushed the Scatter_simulation branch from e6a41de to e04c3af Compare April 13, 2020 08:03

add scatter test on symmetry of output

0c96397

KrisThielemans mentioned this pull request Apr 14, 2020

(potential) axial shift in scatter simulations #495

Open

KrisThielemans added 16 commits April 16, 2020 18:22

Scatter caching safer for old OpenMP versions

819e2dc

ScatterSimulation: add more checks if set_up

818adca

- Currently only check in Debug mode to avoid impact on computation time. - move computation of shift_detector_coordinates_to_origin to set_up (from process_data)

don't call set_up in SingleScatterSimulation::process_data()

0caf7bf

make ScatterSimulation safe for reusing the same object

70249f9

we had some static variables, and forget to clear the cached detection_points_vector

add multiple tests on symmetry for scatter simulation

9e1988e

These are tests for a symmetric object, checking if the result is symmetric. Also enabled the consistency check mentioned in UCL#495. This means that run_scatter_tests/sh will now fail.

minor clean-up of example .par files for scatter

1bdb54b

- some clarifications in sample scatter_simulation.par - removed obsolete file

scatter simulation: renamed "random" to "randomly_place_scatter_points"

d5b9b74

hopefully a better name...

renamed get_ExamInfo_sptr in scatter code to get_exam_info_sptr

40e3c02

now consistent with the rest of STIR

removed some obsolete/accidentally checked-in files

33e89de

Minor clean-up of ScatterSimulation

8125e31

Use scatter_estimate(Bin) as opposed to in terms of 2 detectors such that it will apply for SPECT at some point.

changed parsing keywords to PET Single Scatter Simulation

96baf65

original values did not reflect this was PET specific. Also used more standard-form keywords for the simulate_scatter .par file. Also fixed some issues in example files and STIR-UsersGuide

fix calling scatter_simulation::set_up() at every iteration

ea527b9

This was not necessary earlier-on as SingleScatterSimulation::process_data() was calling set_up(), but no longer...

fix ScatterEstimation parsing

0a4b2e4

- use standard keywords for reconstruction and scatter simulation - allow specifying these in the parameter file itself (as opposed to requiring a separate file)

slightly increased threshold for comparing scatter output

8146aa1

3 Travis jobs gave errors of 0.124% while our threshold was 0.1%.

KrisThielemans added 2 commits April 21, 2020 07:46

some improvements on documentation for scatter [ci skip]

ac4ec11

This commit finally closes UCL#495

moved scatter-estimation examples files to the examples/scatter [ci s…

0a54fca

…kip]

KrisThielemans merged commit d510555 into UCL:release_4.0.0 Apr 21, 2020

KrisThielemans mentioned this pull request Apr 22, 2020

updating scatter #511

Open

KrisThielemans added a commit to KrisThielemans/STIR-1 that referenced this pull request Sep 17, 2023

Merge pull request UCL#44 from KrisThielemans/merge_master_to_TOF_202…

c5ce566

…3_05_07 Merge master to tof 2023 05 07 + Siemens support

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Scatter simulation #44

Scatter simulation #44

NikEfth commented Sep 12, 2016 •

edited

Loading

KrisThielemans commented Sep 13, 2016

KrisThielemans left a comment

NikEfth left a comment

NikEfth Sep 19, 2016

NikEfth Sep 21, 2016

KrisThielemans commented Jan 23, 2018

KrisThielemans commented Apr 13, 2020

KrisThielemans commented Apr 16, 2020 •

edited

Loading

KrisThielemans commented Apr 21, 2020

KrisThielemans commented Apr 21, 2020

Scatter simulation #44

Scatter simulation #44

Conversation

NikEfth commented Sep 12, 2016 • edited Loading

KrisThielemans commented Sep 13, 2016

KrisThielemans left a comment

Choose a reason for hiding this comment

NikEfth left a comment

Choose a reason for hiding this comment

NikEfth Sep 19, 2016

Choose a reason for hiding this comment

NikEfth Sep 21, 2016

Choose a reason for hiding this comment

KrisThielemans commented Jan 23, 2018

KrisThielemans commented Apr 13, 2020

KrisThielemans commented Apr 16, 2020 • edited Loading

KrisThielemans commented Apr 21, 2020

KrisThielemans commented Apr 21, 2020

NikEfth commented Sep 12, 2016 •

edited

Loading

KrisThielemans commented Apr 16, 2020 •

edited

Loading