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)

  • 6
    $\begingroup$ 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? $\endgroup$
    – MarcoB
    Jun 28, 2022 at 14:38
  • $\begingroup$ 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. $\endgroup$
    – NotaChoice
    Jun 28, 2022 at 15:14
  • $\begingroup$ 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. $\endgroup$
    – NotaChoice
    Jun 28, 2022 at 15:43


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.