​

​
A Wolfram Notebook on the Radon Transform and CT Scans

Introduction

A century ago, in 1917, the Austrian mathematician Johann Radon (1887–1956) published a paper entitled “On the Determination of Functions from Their Integral Values along Certain Manifolds.” In two dimensions, the problem studied by Radon involved the reconstruction of a function of two variables from the values of its integrals over the set of all lines in the plane. The representation of a function in terms of its integrals along lines is now called its Radon transform in his honor. The operation studied by Radon in his paper of 1917 is referred to as the inverse Radon transform.

Radon’s research was motivated by developments in pure mathematics, and yet the ideas in his paper provide the theoretical basis for computed tomography (CT), which is a vital tool for screening and diagnosis in modern medicine.

Here, I would like to explain the connection between the Radon transform and CT scans using the built-in functions for computing symbolic Radon transforms and the image processing capabilities of the Wolfram Language.

I will begin by giving examples for computing the Radon transform and its inverse in closed form using RadonTransform and InverseRadonTransform, respectively. Next, I will show how a two-dimensional image of an object can be reconstructed from its image in Radon transform space using InverseRadon. Finally, I will explain how a practical application of these ideas allows doctors to probe a diseased organ in the human body by sending x-rays at different angles through it in a CT scan machine.

The Radon Transform

Let
f(x,y)
be a function of the two variables
x
and
y
. Consider a line
L
given in normal form by the equation
pxcos(ϕ)+ysin(ϕ)
, as shown in the following figure.

Then the integral of
f
along the line is given by:

ℛ(p,ϕ)=∫
∞
-∞
f(pcos(ϕ)-ssin(ϕ),psin(ϕ)+scos(ϕ))s

As the line
L
varies, one obtains a function
ℛ(p,ϕ)
of the two variables
p
and
ϕ
. This function is the Radon transform of
f
.

For example, suppose that
f
is the following Gaussian exponential function.

In[1]:=
f[x_,y_]:=E^(-x^2-y^2)

As seen in the following plot,
f
has circular symmetry about the origin and decays rapidly to
0
as
x
and
y
increase in absolute value.

In[2]:=
Plot3D[f[x,y],{x,-2,2},{y,-2,2}]
Out[2]=

Hence, one expects that the Radon transform of
f
will have an expression that does not depend on the angle
ϕ
and decays to 0 for large values of
p
. This is indeed the case, as seen in the calculation below, which shows that the Radon transform is a Gaussian function of the variable
p
alone.