issue-1501__FBP-for-3d-Euler

[robtovey]

This is an implementation of the 3D filtered back-projection suggested in the issue #1501.

Firstly, an apology that I'm not hugely familiar with things on github so apologies for the mistakes that I have made/will make.

The two files are an update of the FBP code and an example which is very similar to the other examples on FBP but also compares against a least squares solution computed by Landweber iteration to check it's performing as expected.

Mathematically, the FBP is as expected, and follows the argument on #1501. The only extension in the code is to compute the Fourier transform of the kernel straight from wiki, number 502 and to realise that 2D Euler angles/spherical coordinates are a non-uniform parametrisation of the sphere.

There are a couple of things of note:

1. There is a weird scaling factor needed but seems to be consistent. It is on line 411: fudge = 0.03 / ray_trafo.domain.cell_volume I see no reason for why the discretisation of the volume should have any influence on a Fourier scaling in data-space but maybe this is something weird coming from astra? 0.03 is also just a meaningless number.
2. I can't think of any quantitative tests to perform to verify accuracy. The accuracy of FBP is a continuum result and in 3D I think the cost of getting near this limit is probably too high. The example code includes a comparison with a Landweber iteration which should compute the discrete analogue. The plots show they are reasonably close and the error gets lower with more projections but there is a difference, part of this may be down to point 4 below.
3. I have updated documentation in both files. The only qualitative difference compared with other methods is that the FBP equation is only valid when you specify 2 Euler angles. The case with 3 is more complicated as you can have over-complete parametrisations of the sphere so I don't know what to do there but they are the same geometry. I have explained this as Parallel3dEulerGeometry : Exact reconstruction, only implemented for 2d angle partitions.
4. The scaling I used for the integral over the sphere is actually very natural and I wonder if it should/could be built into the operator to make the behaviour closer to the continuum model. In particular, sampling uniformly in spherical coordinates means you sample very densly at the poles. Thus, at the moment, doing a least-squares fit will heavily bias to fitting the data at the poles. I wonder if this it is the natural scaling to treat each projection equally (as is currently done) or keep the metric uniform on the sphere and have a non-uniform weighting? This is actually a difference between the FBP and Landweber iteration at the moment.

Anyway, I hope this is of interest and let me know what more I can do.

[pep8speaks]

Checking updated PR...

No PEP8 issues.

Comment last updated at 2019-07-19 09:09:32 UTC

[adler-j]

Thank you for this pull request! We're currently quite busy so review might take a while sadly. Would you be able to fix the PEP8 issues indicated by @pep8speaks?

[robtovey]

No problem, I know the feeling! Don't worry about rushing for my sake, I'll try and get pep happy by the end of the week then just get back to me whenever. Thanks for taking the time to update me!