Prompt histogram calculation

[markus-jehl]

Changes in this pull request

Added a function to translate a projdata input into a histogram of events each crystal contributed to.

Testing performed

Related issues

fixes #1666

Checklist before requesting a review

• I have performed a self-review of my code
• I have added docstrings/doxygen in line with the guidance in the developer guide
• [] I have implemented unit tests that cover any new or modified functionality (if applicable)
• The code builds and runs on my machine
• [] documentation/release_XXX.md has been updated with any functionality change (if applicable)

Contribution Notes

Please tick the following:

• The content of this Pull Request (the Contribution) is intentionally submitted for inclusion in STIR (the Work) under the terms and conditions of the Apache-2.0 License.
• I (or my institution) have signed the STIR Contribution License Agreement (not required for small changes).

[KrisThielemans]

How does this actually differ from make_fan_sum_data https://github.com/UCL/STIR/blob/master/src/buildblock/ML_norm.cxx#L1442?

There are implementation differences (use of get_all_dot_poss, which I think is good), and OpenMP (not sure why we forgot to do that for make_fan_sum_data ), but I think conceptually the same?

[markus-jehl]

How does this actually differ from make_fan_sum_data https://github.com/UCL/STIR/blob/master/src/buildblock/ML_norm.cxx#L1442?

There are implementation differences (use of get_all_dot_poss, which I think is good), and OpenMP (not sure why we forgot to do that for make_fan_sum_data ), but I think conceptually the same?

Ah, that's why it would have been good to discuss the approach first! Didn't know there was already that functionality somewhere. So I might just merge my implementation into this, i.e. add OpenMP support.

[KrisThielemans]

sorry!

I might just merge my implementation into this, i.e. add OpenMP support.

yes please, and use get_all_det_bins.... It'll come in handy later.

[markus-jehl]

Currently make_fan_sum_data says it does not work with view_mashing and TOF. Do you want to keep that behaviour? The create_crystal_histogram function can cope with both.

[KrisThielemans]

exactly, using get_all_det_bins and collapsing over TOF will make that function more useful. Essentially, your function can replace make_fan_data I believe (obviously, keeping the original name). Might be good to first do a quick check they indeed give the same result :-)

[markus-jehl]

Funnily, they don't give the same result. My function gives exactly the double the counts of the projdata, while make_fan_sum_data returns a little bit less:

But in this case I might really just replace the function completely with the new implementation.

[KrisThielemans]

hmmm. that's worrying. Hopefully this is a case of sum() being sensitive to numerical errors (although for nD arrays in STIR, it should be ok, as it computes it in stages).

you could display the "fan" and "histogram" results i suppose to see what/where the difference is.

Obviously, replacing the old function with a new one that doesn't give the same result looks rather dangerous....

[markus-jehl]

Hah! There is indeed a very small variability in my function - I will look into this a bit more. But I would expect that the sum is exactly double the bin counts (unless there is view mashing), since each bin count should give two counts in the histogram. Or am I missing something?

[KrisThielemans]

I think you're right. I can only see numerical issues. What happens if you use Array<2,double> for data_fan_sums (or your equivalent)

[KrisThielemans]

could do that in the calculation-loop only, and keep the output argument as-is.

[markus-jehl]

We're essentially dealing with integers in my example, so floating point issues don't come into it. The reason for the difference was that I'm summing from min_tangential_pos to max_tangential_pos, and data_fan_sums goes from -half_fan_size (which is min_tangential_pos + 1) to half_fan_size. This is also why in my output the histogram is not completely homogeneous:

Note, the asymmetry is just due to my test, where I fill the projdata with homogeneous values. In realistic cases the min_tangential_pos bin would be essentially zero anyways, so then the data will be more even. Therefore I think that my approach is actually more accurate, no?

[KrisThielemans]

We're essentially dealing with integers in my example, so floating point issues don't come into it.

that's not true, really. Integers get converted to floats, and then added, causing precision loss. Anyway, that wasn't the issue, luckily.

The reason for the difference was that I'm summing from min_tangential_pos to max_tangential_pos, and data_fan_sums goes from -half_fan_size (which is min_tangential_pos + 1) to half_fan_size.

This needs some careful checking. This code is used in the ML norm estimation. It is essentially the "backprojection" step, and needs to be consistent with the forward model.

One reason that the half_fan_size was used was that for the list-mode equivalent lm_fansums, at the time, we didn't have the min/max_tangential_pos information, so passed the half_fan_size as a parameter. If we can get rid if that, so much the better.

As usual, this will therefore need some careful checking. At least for the ML_norm stuff, there are some self-consistency checks that will help.

sorry that this is getting more involved.

[markus-jehl]

As discussed, I opted to use the slightly sub-optimal half_fan_data implementation for my work, since it had only a small impact on the final results. The PR is now ready from my side.

[markus-jehl]

But should we maybe create another issue to track that ML_norm should be modified at some point?