Tigre algorithm wrapping

[MargaretDuff]

Changes

• Wrapped all tigre algrorithms so they can be called as a Reconstructor for example:

from cil.plugins.tigre import tigre_algo_wrapper
algo = tigre_algo_wrapper(name='sart', initial=None, image_geometry=ig, data=absorption, niter=10)
img, qual = algo.run()

This has been tested on all geometries for most of the currently available tigre algorithms, see here https://github.com/MargaretDuff/cil_test_notebooks/blob/master/os-sart-tigre/basic_tigre_wrapping.ipynb and the unit tests.

There is a bug with those Tigre algorithms that enforce TV regularisation using the "im3ddeoniser" in tigre, if you have a 2d geometry and if you are running on multiple gpus. A warning will show for users and setting the gpuids in the kwargs works ( see the test notebook). A bug report has been put into tigre CERN/TIGRE#681

Discarded options:

Alternative discarded solution

### Option 1 - wrapped tigre's OS SART, SIRT and SART as CIL algorithms so they could be called e.g. ```python from cil.plugins.tigre import ART, SIRT, SART, OSSART OSSART_algorithm_wrap = OSSART(initial=None, image_geometry=ig, data=b, blocksize=10) OSSART_algorithm_wrap.run(5) solution= OSSART_algorithm_wrap.solution ``` To do: - [x] Test on 2D geometries and cone beam geometries - [ ] Test the callback functionality and algorithm stopping - [ ] Decide if we want to explicitly expose more options from Tigre e.g. OrderStrategy for OS SART. Any unused kwargs will be passed. - [x] Investigate memory-efficient ways of calculating the objective - [ ] Decide what the user needs to provide for the algorithm to run. At minimum, data (with acquisition geometry) and image geometry. Image geometry could be extracted from an initial, if that is passed. If both image geometry and initial is passed, do we check their geometries are the same? - [ ] Look at ways we can save/reuse memory - [ ] Write unit tests - [x] Write doc strings - including explaining that this is just for CT forward problems.

See here: https://github.com/MargaretDuff/cil_test_notebooks/blob/master/os-sart-tigre/cil_data_tigre_sart_sirt.ipynb

Related issues/links

Closes #2159

Checklist

• I have performed a self-review of my code
• I have added docstrings in line with the guidance in the developer guide
• I have updated the relevant documentation
• I have implemented unit tests that cover any new or modified functionality
• CHANGELOG.md has been updated with any functionality change
• Request review from all relevant developers
• Change pull request label to 'Waiting for review'

Contribution Notes

Please read and adhere to the developer guide and local patterns and conventions.

• The content of this Pull Request (the Contribution) is intentionally submitted for inclusion in CIL (the Work) under the terms and conditions of the Apache-2.0 License
• I confirm that the contribution does not violate any intellectual property rights of third parties

[MargaretDuff]

Discussed with developers today:

• We would like to try a general wrapping of tigre algorithms e.g. tigre.algo( name, initial, image_geometry, data, niter, **kwargs) which then converts the geometry, data and initial to tigre compatible objects and passes to tigre as name(initial, data, tigre_geom, tigre_angles, niter, **kwargs). Currently all of the tigre algorithms take this form. We would provide this basic wrapping and pass to tigre documentation for information and exploration
• There is an open question about what we then pass back to CIL. We could return just the result ( which we could pass back to a CIL object) or return the result and the tigre algorithm or something else.
• Unit tests would just check the CIL parts of this - we would make no claims on the accuracy of the algorithms.

[MargaretDuff]

Discussed with developers:

• We prefer option 2, it works with all current tigre algorithms and it is clearer that it is Tigre's algorithm and the credit goes to them
• We want to wrap it as a processor with an init, set_input and get_output. get_output could take a number of iterations and still be decided where the initial is set
• We should not restrict which tigre algorithms can be passed to this, just maybe have a list of tested algorithms in the docstring.
• We should link to the tigre documentation in the docstring

[MargaretDuff]

@gfardell - I was having another think about this - we discussed writing the wrapper as a processor but to force the code to have an init, set_input and get_output seems challenging e.g. in the init we would pass: image geometry, method, iteration number and maybe initial, set_input would pass the data and would actually initialise the algorithm and then process would do everything. To me it feels a bit backward and you wouldn't get any computational benefit from this set up

I was wondering if we could write it as a reconstructor like in the recon class with just an init and run method? It also already has the machinery to check the image and acquisition geometries.

I think it might make sense if the processors go from the same type of data e.g. image to image and acquisition to acquisition but the reconstructors go from acquisition to image

[gfardell]

Thank you @MargaretDuff. I think this is a nice solution to wrapping different algorithms. I think the comments are all fairly minor, just about useability and cleaning up some code. We can discuss if you have the resources to make the changes yourself or if you want one of us to take it over from here.

[gfardell]

This might all be to support an older version of tigre so probably not necessary here.

[gfardell]

The argument order seems strange to me. I would say that the algorithm string, data and number of iterations are necessary, and that Image geometry and initial are optional.

There might be other CIL conventions that motivated this order and would be worth discussing.

I would also go with algorithm_name and num_iterations or iterations.

[MargaretDuff]

I think the choice was the ordering that tigre uses, with the addition of 'name' at the beginning. Agree algorithm_name is better.

[gfardell]

if image_geometry is None then a default should be used.

[MargaretDuff]

Ok - we will just require passing the data

[gfardell]

I think we should only care about data here, and I suspect the parent class will handle the checking that it's AcquisitionData and it exists - it's worth checking that's a useful error message though.

[MargaretDuff]

Makes sense, I have changed this so it only raises an error if data is not passed

[gfardell]

In the case where initial is None, allocating zeros and then copying isn't necessary.

I'd also consider carefully if we need a copy at all. Does tigre change the value? Or is this array reused for the solution, in which case we can hijack this behaviour for our out.

[gfardell]

I feel like this is something we should be exposing elsewhere and have a consistent interface (i.e. take a list of IDs).

But in this case it seems like a kwarg passed directly to the algorithm, which suits out 'light touch' approach, are there other kwargs they might pass?

[MargaretDuff]

Yes, they can pass a huge number and it differs per algorithm which ones they can pass

[gfardell]

I haven't understood what this is about, is it that you don't want to inherit the method from the parent class? If it's not something we want in the API maybe we could separate it's use in the parent so we have a private method for checking the input, and remove the interface for this method?

[MargaretDuff]

Yes, it's to make it comptible for being a "reconstructor" and is called by the parent class at initialisation but I don't want the user to use it (at the moment). Agree, making it private method would be useful

[gfardell]

I think it's worth being certain we need to support both of these.

[gfardell]

is result is the data array this would still pass. I think it should be under the isinstance(result,tuple) logic to unpack the tuple.

[MargaretDuff]

Have made this clearer

[gfardell]

If out=None then we should wrap the numpy array as a CIL ImageData and remove the allocate and fill.

If out is passed and we are only using it to fill the result maybe that's something we should consider carefully. I would be tempted to not take out in this case, but maybe that's not consistent with CIL elsewhere.

[MargaretDuff]

I have addressed the simple comments. I think there are two major issues left - how out is dealt with and how it behaves as a reconstructor i.e. what of the parent we want it to inherit etc.

Happy for the core developer team to take this over.