Radon Transform works unexpectedly
23 views (last 30 days)
Show older comments
Koen Van den Brandt
on 21 Jun 2023
Edited: ProblemSolver
on 28 Jun 2023
I've been trying to understand the radon transform, and the documentation helped a lot.
There's one last thing I did not understand, I hope you can help me.
If I do a radon transform of a unitary matrix, at 0 degrees, I would expect the result would be the sum of the collumns since it's a parrallel beam projection.
However, as shown below the outer rows (row 4 and 7) show a different value, where a part of the total leaks into the neighbouring rows.
Could someone explain this to me?
Thanks in advance.
Koen
radon(ones(4,4), 0)
0 Comments
Accepted Answer
Matt J
on 28 Jun 2023
Edited: Matt J
on 28 Jun 2023
Instead of radon(), you could use one of the parallel beam forward projectors (the SIddon algorithm is probably best) from the TIGRE library,
1 Comment
ProblemSolver
on 28 Jun 2023
Edited: ProblemSolver
on 28 Jun 2023
You can used Siddon algorithm as suggested by @Matt J. But I have no experience using the Siddon algorithm. However, there is another tool within MATLAB that is available that you can use:
ASTRA Toolbox: This toolbox is generally for tomographic reconstruction and forward projection. It provides a wide range of algorithms for various tomographic geometries, including parallel beam. ASTRA also offers highly optimized CPU and GPU implementations, making it suitable for large-scale simulations (I assume you will be working with multiple images). You can use the astra_create_projection2d function to perform parallel beam forward projection.
If you through the link, you can look at the examples to really understand how each functions are working.
More Answers (1)
ProblemSolver
on 27 Jun 2023
The behavior you observed in the radon transform result is due to the discretization and interpolation process used in the algorithm. The Radon transform is computed using the projection of lines through the image and accumulating the values along these lines.
In the case of a unitary matrix (a matrix of ones), the Radon transform at 0 degrees represents the sum of the columns, as you correctly expected. However, the interpolation process used in the Radon transform algorithm assigns the accumulated values to discrete bins or pixels, which can result in some leakage of values into neighboring bins. In your example, the input matrix is 4x4, and the resulting Radon transform is a 9x1 vector. The values in the vector represent the accumulated values at specific angles, including the 0-degree angle. The non-zero values you observed in the outer rows (row 4 and 7) indicate that some of the accumulated values from the columns have leaked into the neighboring rows due to the interpolation process.
To mitigate this leakage, you can increase the resolution of the Radon transform by increasing the number of projection angles or by increasing the size of the input matrix. This will provide a more accurate representation of the desired parallel beam projection.
Here's an example that demonstrates the effect of increasing the matrix size:
inputMatrix = ones(8, 8); % Increase the size of the input matrix
radonTransform = radon(inputMatrix, 0);
disp(radonTransform);
By increasing the matrix size to 8x8, you will observe a more accurate parallel beam projection without significant leakage into neighboring rows.
I hope this explanation clarifies the behavior you observed.
3 Comments
ProblemSolver
on 28 Jun 2023
@Matt J: Thank you for providing feedback and sharing the example. I appreciate the opportunity to learn from it. After reviewing the example, I understand that increasing the resolution of the Radon transform by either increasing the number of projection angles or increasing the size of the input matrix does not necessarily mitigate the leakage or produce the desired results. I apologize for any misunderstanding or misinformation in my previous comment. It's important to consider the limitations and potential issues of using the Radon transform, particularly in cases where non-uniformities and jagged artifacts may be introduced. Your example highlights the challenges involved, and it serves as a reminder that the Radon transform may not always be suitable for certain applications. Thank you for bringing this to my attention, and I will keep this in mind going forward.
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!