# Implementing the 3D Radon transform

I am wondering how to implement the Radon transform, the 3D Radon transform, that is, given a 'density' function $$f: \mathbb{R}^3\to \mathbb{R}$$ The Radon transform of $$f$$ is $$Rf(s,w)= \int_{x\cdot w= s} f(x)d\mu(x)$$ where $$\mu$$ is the Lebesgue usual measure on the plane $$x\cdot w= s$$ on which the integration is happening.
Here $$x= (x_1, x_2, x_3)$$ is a vector.

I am more interested in a symbolic implementation, let me illustrate with an example :

For a (nice) function that depends only on the radial distance $$r=\sqrt{x_1^2+ x_2^2+ x_3^2}$$ , and is piecewise polynomiale :

\begin{align*} f(r)=& \ (r-1)^2(r-3)^2= r^4- 8r^3+ 22r^2- 24r+ 9 \quad \quad \text{ if } 1\leq r\leq 3 \\ f(r)=& \ 0 \quad \quad \quad \quad \quad \text{ otherwise } \end{align*}

Then the Radon of $$f$$ depends only on $$s$$ ( that is $$R(f)(s, w_1)= R(f)(s, w_2)= R(f)(s)$$ )
We can find that for $$s\leq 1$$, by just doing some calculus, the Radon of $$f$$ is : \begin{align*} R(f)(s, w)=& \int_C f(r)d\mu= 2\pi\int_{\sqrt{1-s^2}}^{\sqrt{9-s^2}} r_1f\left(\sqrt{r_1^2+ s^2}\right)\,dr_1 \\ & \text{ (where } C \text{ is the circular sector of intersection) } \\ =& 2\pi\int_{1}^{3} rf\left(r\right)\,dr \quad \text{ by the change of variable r= \sqrt{r_1^2+ s^2} } \\ =& 2\pi\left[ \frac{r^6}{6}- \frac{8r^5}{5}+ \frac{22r^4}{4}- \frac{24r^3}{3}+ \frac{9r^2}{2} \right]_1^3 \\ =& \frac{64}{15}\pi \end{align*}

And for $$1\leq s\leq 3$$ the Radon of $$f$$ is : \begin{align*} R(f)(s, w)=& \int_D f(r)d\mu= 2\pi\int_{0}^{\sqrt{9-s^2}} r_1f\left(\sqrt{r_1^2+ s^2}\right)\,dr_1 \\ & \text{ (where } D \text{ is the disk of intersection) } \\ =& 2\pi\int_{|s|}^{3} rf\left(r\right)\,dr \quad \text{ by the change of variable r= \sqrt{r_1^2+ s^2} } \\ =& 2\pi\left[ \frac{r^6}{6}- \frac{8r^5}{5}+ \frac{22r^4}{4}- \frac{24r^3}{3}+ \frac{9r^2}{2} \right]_s^3 \\ =& \pi\left( 782,992224- \frac{s^6}{3}+ \frac{16s^5}{5}- 11s^4+ \frac{48s^3}{3}- 9s^2 \right) \end{align*} And, of course, for $$s> 3$$ the Radon of $$f$$ is zero.

You can see now, how the Radon transform is actually a piecewise polynomial of $$s$$, one can try any other polynomials $$f= f(r)$$ and check that he will always find $$R(f)(s)$$ to be a polynomial.

Is there a way to make the software do this and output the expression $$R(f)(s)$$ when given the input $$f= f(r)$$ ?

We should be able to parametrize the integration on the plane $$x\cdot w= s$$ ?

How about outputting the expression $$R(f)(s, w)$$ when given the general input $$f= f(x)= f(x_1,x_2,x_3)= f(r,\theta, \phi)$$ ?
(if you like spherical coordinates)

• What is the Mathematica question here? There is a lot of math and definitions, and no code. What are the perceived shortcomings of the built-in Radon and RadonTransform? Jun 28, 2022 at 14:38
• These are commands for the 2D radon transform as I discover them, you give it a picture or image and it compute the integral over a line. I am looking for an integral over a plane for 3D scans. Jun 28, 2022 at 15:14
• I wonder if the source code of the reference.wolfram.com/language/ref/RadonTransform.html is available or open ? It will be nice to generalize it to 3D functions. Jun 28, 2022 at 15:43