Added a method to get detector effeciencies without creating an event

[MaxTousss]

Summary

Added a method in the CPP helpers to get efficiency from two detection bins, without creating an event. It is only quality of life, since it make no sense to create an event when we simply want to get all the efficiency coefficients.

Testing performed

Compilation, eye and logic.

Related issues

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
• The code builds and runs on my machine

Contribution Notes

This is only basic polymorphism, no need to bother everyone.

Please tick the following:

• The content of this Pull Request (the Contribution) is intentionally submitted for inclusion in the ETSI software (the Work) under the terms and conditions of the Apache-2.0 License.

[gschramm]

would be nice to have a similar change in the corresponding python helper here:
https://github.com/ETSInitiative/PETSIRD/blob/main/python/petsird/helpers/__init__.py

[MaxTousss]

Polymorphism is not a thing in python, so it would need some gymnastic which I am not sure how best to take on or defines two functions with different name. Which you believe is best?

[nkarakatsanis]

Looks useful to me to have this extra morphology of get_detector_efficiencies in the standard. Implementation looks good to me too. Thanks @MaxTousss

[KrisThielemans]

Polymorphism is not a thing in python, so it would need some gymnastic which I am not sure how best to take on or defines two functions with different name. Which you believe is best?

Agreed this creates headaches. I guess we could do a type-check, but it's going to slow everything down. Having different names might be best. @naegelejd @casperdcl any suggestions?

[gschramm]

Polymorphism is not a thing in python, so it would need some gymnastic which I am not sure how best to take on or defines two functions with different name. Which you believe is best?

Agreed this creates headaches. I guess we could do a type-check, but it's going to slow everything down. Having different names might be best. @naegelejd @casperdcl any suggestions?

Completely missed the re-definition of the function (thought the the aim was to only keep the new version).
Having two names is probably the way to go in python.

On top, having a few vectorized helpers would be nice to speed up calculations.

[MaxTousss]

I just modified the python code to have two function to get eff: one that take detector ids and a new one that take event in input.

The python analysis/generator helpers were modified following this change.

Note: I did not touch the matlab case. While I can guess how to do it, I can't test it.

Funny discovery: the generator helpers create eff. that are larger than 1.0! XD

[gschramm]

I just modified the python code to have two function to get eff: one that take detector ids and a new one that take event in input.

The python analysis/generator helpers were modified following this change.

Note: I did not touch the matlab case. While I can guess how to do it, I can't test it.

Funny discovery: the generator helpers create eff. that are larger than 1.0! XD

Hehe. I think in real scanners efficiencies can have arb. magnitudes because: (i) vendor include / can correct for global scales at multiple places in the recon pipeline, (ii) it depends on which units we use for the image that is being reconstructed (counts per voxel, activity per ml, ...), (iii) LORs that contain edge crystals (in the presence of big gaps between modules) behave very differently compared to LORs connecting center crystals. so depending on your assumed global scales, some effs. can be > 1, some <1.

Long story short, the effs. are correction factors that depend on many things which easily introduce global scaling factors.

[casperdcl]

I just pushed a Pythonic overload (which also reverts the breaking change caused by renaming); feel free to revert if you prefer.

Polymorphism is not a thing in python

quack yes it is 🦆🐍

[KrisThielemans]

Thanks @casperdcl I guess this makes Python happy, but it will also slow every call down due to the type-check. What about something like this

def get_detection_efficiencies_for_bins(...):
   bla

def get_detection_efficiencies_for_event(...):
   return get_detection_efficiencies_for_detection_bins(something)

# current code with the @typing.overload
# but let get_detection_efficiencies() call one of the above

having a few vectorized helpers would be nice to speed up calculations.

sure, feel free...

we could have "stupid" implementations now (as I did for expand_detection_bins. Even if we don't do it now, probably worth thinking if it's best to pass the 2 detection_bins as pair as opposed to 2 arguments.

[casperdcl]

it will also slow every call down due to the type-check

Do you mean the single isinstance call? That's super quick compared to all the other cruft. Is it really a bottleneck?

[KrisThielemans]

I have no idea. If you say it isn't, I'll believe you.

[casperdcl]

if isinstance(..., int) takes 20ns on my machine :)

[MaxTousss]

Wait, since when does python as polymorphism in that form? At the same time as typing hint was introduced?! Woahhh XD

[MaxTousss]

While 20 ns seems like nothing, it is called O(N^2) where N is the number of detectors. So, it will take a couple of seconds for a whole scanner, but I guess it is marginal compared to everything else in get_detection_efficiencies() XD

[KrisThielemans]

Wait, since when does python as polymorphism in that form? At the same time as typing hint was introduced?! Woahhh XD

It's all smoke and mirrors in my opinion. @typing.overload makes sure that typing is ok, but the actual run-time implementation is the same as before (i.e. do a manual type check). Of course, there is no other way for a run-time interpreter, unless you'd have JIT.

[KrisThielemans]

While 20 ns seems like nothing, it is called O(N^2) where N is the number of detectors. So, it will take a couple of seconds for a whole scanner, but I guess it is marginal compared to everything else in get_detection_efficiencies() XD

ok. Let's keep @casperdcl 's version then.

Can we still think about the vectorised version in this PR? Main reason I ask is to think about the signature:

DetectionBinPair = Tuple[petsird.DetectionBin, petsird.DetectionBin]
get_detection_efficiencies( scanner: petsird.ScannerInformation,
        type_of_module_pair: petsird.TypeOfModulePair,
        event_pair_or_detection_bin: typing.Union[Iterable[petsird.CoincidenceEvent],
                          Iterable[DetectionBinPair]]) -> List[float]

or something like that? or return a numpy array, or a generator, or whatever. This needs more thought and we leave it for the next PR. But for thisPR, the question is if specifying a list of pairs, or 2 lists of bins will be more natural, and hence implementing the current single efficiency case accordingly.

choices, choices...

[KrisThielemans]

Of course, need to fix the run-time error first

  File "/home/runner/work/PETSIRD/PETSIRD/python/petsird/helpers/generator.py", line 270, in get_events
    if get_detection_efficiency(header.scanner, type_of_module_pair,
       ~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
                                event) > 0:
                                ^^^^^^
  File "/home/runner/work/PETSIRD/PETSIRD/python/petsird/helpers/__init__.py", line 110, in get_detection_efficiency
    if isinstance(event_or_detection_bin_1, petsird.DetectionBin):
       ~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/usr/share/miniconda/envs/yardl/lib/python3.13/typing.py", line 1375, in __instancecheck__
    return self.__subclasscheck__(type(obj))
           ~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^
  File "/usr/share/miniconda/envs/yardl/lib/python3.13/typing.py", line 1378, in __subclasscheck__
    raise TypeError("Subscripted generics cannot be used with"
                    " class and instance checks")
TypeError: Subscripted generics cannot be used with class and instance checks

[gschramm]

TBH: I wouldn't worry too much about the runtime right now. I see the current python helpers as proof of concept. To get sth really efficient and "fast" on real world data sets, we need vectorized helper functions and a more efficient implementation of the yardl generated reading functions anyway.

@KrisThielemans: since we will most of the time process many events contained in an eventblock, why not having future helpers that that multiple events as input and output numpy / cupy / xp arrays?

@MaxTousss: isinstance seems to be part of python since 2.0 already ;)

[MaxTousss]

@gschramm My bad, when I said "in that form", I meant the @typing.overload thing. Like, in my heart, it is not really polymorphism if I need to have an explicit if/else in the function that deal with the different type of input.

For the signature thing, hmmm, my first though is that it should follow the format that event_pair_or_detection_bin will have once it is also "vectorised"? And since we need to convert the detection_bins into expand_detection_bin, it seems to be better to aim for numpy array for easy&lazy gain in computing time.

[gschramm]

@gschramm My bad, when I said "in that form", I meant the @typing.overload thing. Like, in my heart, it is not really polymorphism if I need to have an explicit if/else in the function that deal with the different type of input.

For the signature thing, hmmm, my first though is that it should follow the format that event_pair_or_detection_bin will have once it is also "vectorised"? And since we need to convert the detection_bins into expand_detection_bin, it seems to be better to aim for numpy array for easy&lazy gain in computing time.

I completely agree @typing.overload is no real polymorphism.
I would vote for doing all calculations (expansions, ...) on vectors in python.
The most basic low-lower helper should "just" return the raw fields that are stored in an event / a series of events.

In an ideal scenario (but I don't think that this is currently supported by the way that yardl encodes the data),
we would have "big" chuncks of events that all consists of a contigous series of event words with fixed length (as most vendors do right now). If we have that, we could use numpy's memmap to more efficiently load those chunks into arrays (after decompressing the chunk).

[casperdcl]

fixed types methinks

[KrisThielemans]

TBH: I wouldn't worry too much about the runtime right now. I see the current python helpers as proof of concept. To get sth really efficient and "fast" on real world data sets, we need vectorized helper functions and a more efficient implementation of the yardl generated reading functions anyway.

Of course. Ideally we make our helpers to future-proof the user for any changes, but that'll only work as far as it goes.

@KrisThielemans: since we will most of the time process many events contained in an eventblock, why not having future helpers that that multiple events as input and output numpy / cupy / xp arrays?

Not very easy if "events" are structures in numpy, but certainly should aim for that!

For the signature thing, hmmm, my first though is that it should follow the format that event_pair_or_detection_bin will have once it is also "vectorised"? And since we need to convert the detection_bins into expand_detection_bin, it seems to be better to aim for numpy array for easy&lazy gain in computing time.

good point, as we have

PETSIRD/python/petsird/helpers/__init__.py

Lines 33 to 36 in ed570d1

	def expand_detection_bins(
	scanner: petsird.ScannerInformation, type_of_module: petsird.TypeOfModule,
	detection_bins: typing.Iterable[petsird.DetectionBin]
	) -> list[petsird.ExpandedDetectionBin]:

I think it is best to stick to the current 2 argument signature.

I therefore suggest to merge as-is. ok?

[gschramm]

Not very easy if "events" are structures in numpy, but certainly should aim for that!

(with my limited C/Cpp knowledge), I thought if the event structure is fixed (let's say two ints and a float), we can define a custom np dtype a read like:

record_dtype = np.dtype([
    ('i1', np.int32),     
    ('i2', np.int32),  
    ('f',  np.float32), 
], align=False) 

arr = np.fromfile('records.bin', dtype=record_dtype, count=N)

or sth similar using ctypes if we have a pointer to the start of the eventblock.

Then we should be able to access the float, and two int arrays like:

floats = arr['f']
i1s    = arr['i1']
i2s    = arr['i2']

Obviously that only works if the data is not compressed.

[KrisThielemans]

@gschramm when using ctypes, everything can be done. But it's not so easy for yardl to generate all this of course. If it can be done, in my opinion, it is quite unsatisfactory, as it doesn't give you access to any methods defined on the type, but that seems just to be a python/numpy limitation. This would probably be discussed with @naegelejd in the relevant yardl issue.

Obviously that only works if the data is not compressed

The data is compressed, so the yardl API will first have to decompress it, or it needs to do it on the fly, via some kind of plug-in.

In the mean time, I'll merge this PR!