Inverse Problems Course Notes — The X-Ray Transform

We motivated the study of the Radon transform with a tomographic problem: given the change in intensity of X-rays along all lines through a region, can we reconstruct the attenuation (think of this as density) in the region?

For our purposes in this problem, light travels along straight lines and in two dimensions, those lines are hyperplanes.  So inverting the Radon transform — which sums functions over hyperplanes — solves our problem in 2-D.  But in higher dimensions this is not what we need.  We need to integrate along lines, not hyperplanes.

So we introduce the X-ray transform.

Definition Given a function $f: \mathbb{R}^n \rightarrow \mathbb{C}$, the X-ray transform of $f$ is defined as

$Xf(x, \theta) = \int_{\mathbb{R}} f(x + t\theta) dt$

where $x\in \mathbb{R}^n, \theta \in S^{n-1}$.

The first thing to notice is that there is some redundancy here. For all $x, \theta, s$ we will have

$Xf(x + s\theta, \theta) = \int_{\mathbb{R}} f(x + s\theta + t\theta) dt = \int_{\mathbb{R}} f(x + \tilde{t}\theta) d\tilde{t} = Xf(x, \theta)$

So the natural domain of $Xf$ is actually smaller, we can pick just one representative from the line $\{x + s\theta | s \in \mathbb{R} \}$ without losing any information.   There is a distinguished point on this line: the point $E_{\theta}(x) \doteq x - \langle x, \theta\rangle\theta$ is the unique point on the line that is perpendicular to $\theta$, $\langle E_{\theta}(x), \theta \rangle = 0$.

Since $Xf(x, \theta) = Xf(E_{\theta}(x), \theta)$ for all $x \in \mathbb{R}^n, \theta \in S^{n-1}$, it is natural to think of $Xf$ as a function $\theta \in S^{n-1}$ and $x \in \theta^{\perp} = \{ x | \langle x, \theta \rangle = 0\}$, the hyperplane normal to $\theta$.  Formally,  the domain of $Xf$ is

$T = \{(x, \theta) \in \mathbb{R}^n\times S^{n-1} | x \in \theta^{\perp} \}$

Where $\theta^{\perp}$ is the hyperplane normal to $\theta$.  This is a manifold of dimension $2n - 2$.

Now $\theta^{\perp} \equiv \mathbb{R}^{n-1}$, and we will use standard Lebesgue measure on this space, denoting it by $dH_{\theta^{\perp}}$.  As with the hyperplanes in the Radon transform, we will keep using the fact that

$dH_{\theta^{\perp}}dt = dx$

is standard Lebesgue measure on $mathbb{R}^{n}$.

Function Spaces on $T$

To proceed we need to define the standard function spaces on this domain.  There is nothing surprising here. Let’s start with $C^{\infty}$

$C^{\infty}(T) = \{ f: T \rightarrow \mathbb{C} | \partial_{\theta^{\perp}}^\alpha \partial_{\theta}^\beta f \in C^0(T) \forall \alpha, \beta\}$

Here $\partial_{\theta^{\perp}}$ is a tangential derivative to $\theta^{\perp}$. For example, if

$\theta = (1, 0, \dots, 0)$

then take

$\partial_{\theta^{\perp}} = \partial_{x_j}, \quad\quad j = 2, \dots, n$

For any other hyperplane, we can rotate it into this position and use these operators to define $\partial_{\theta^{\perp}}$.

Now we will define$C_0^{\infty}$.  The sphere is compact, so we only need to worry about compact support in the $\theta^{\perp}$ directions.

$C_0^{\infty}(T) = \{ f \in C^\infty(T) | \exists R s.t. f(x,\theta) = 0 \forall |x| \geq R \}$

(Sorry — having wordpress-latex trouble. Will try to get this fixed.) The Schwartz space is similar:

$\mathcal{S}(T) = \{ f \in C^{\infty}(T) | \sup_{x \in \theta^{\perp}}|x^{\beta}\partial^{\alpha}_{\theta^{\perp}}\partial^{\gamma}_{\theta} f | \leq C_{\alpha, \beta, \gamma} < \infty \}$

As with the spaces defined for the radon transform, we could choose to ignore the derivatives in the angular variables without impacting the subsequent theory, but we will not do that.

In the same way, we can define $\mathcal{E}^\prime, \mathcal{D}^\prime, \mathcal{S}^\prime$, etc.  Topologies on these spaces are also defined in the natural way.

The Transpose of $X$

It is easy to see that $X$ is injective: the Radon transform $R$ is injective and $R$ can be computed from $X$ by integrating the value of $Xf$ over all the lines that make up the hyperplanes in the definition of $Rf$.  But we will do better and produce an explicit inversion formula.

First we need to find the formal transpose of $X$, i.e. the operator $X^t$ such that $\langle Xf, g\rangle = \langle f, X^tg\rangle$.

Let $f \in C_0^{\infty}(\mathbb{R}^n), g \in C^{\infty}(T)$

$\langle Xf, g\rangle = \int_{S^{n-1}}\int_{\theta^{\perp}} Xf(x,\theta)g(x,\theta) dH_{\theta^{\perp}} d\theta$

$= \int_{S^{n-1}}\int_{\theta^\perp}\int_{\mathbb{R}} f(x + t\theta) g(x, \theta) dt dH d\theta$

As $x$ varies over $\theta^\perp$ and t varies over $\mathbb{R}$ in this integral, the variable $y = x + t\theta$ varies over all of $\mathbb{R}^n$ and

$dy = dH_{\theta^\perp}dt$

is standard Lebesgue measure. If we perform this change of variables, we must remember that $g(\cdot, \theta)$ is only defined on $\theta^\perp$, but when $y = x + t\theta$, we have $x = y - t\theta = y - \langle y, \theta\rangle\theta = E_\theta(y)$. So

$\langle Xf, g\rangle = \int_{\mathbb{R}^n}f(y)\int_{S^{n-1}}g(E_\theta(y), \theta) d\theta dy$

So we will define

$X^tg(x) = \int_{S^{n-1}}g(E_\theta(x), \theta) d\theta$

Proposition The maps

$X: C^\infty_0(\mathbb{R}^n) \rightarrow C_0^\infty(T)$

$X^t: C^\infty(T) \rightarrow C^\infty(\mathbb{R}^n)$

are linear and continuous.

Proof Exercise.

The Normal Operator

We will use a different approach to develop the inversion formula for the X-ray transform.  We could parallel our development for the Radon transform, and we could rederive the Radon inversion formula using this approach.  It is good to see more than one way to look at a problem.

In linear algebra, it is often convenient to work with an operator $X^tX$.  It is a normal operator and is usually “much better behaved” than $X$.

$X^tXf(x) = \int_{S^{n-1}}Xf(E_\theta(x), \theta) d\theta$

$= \int_{S^{n-1}}\int_{\mathbb{R}}f(x + t\theta)dt d\theta$

$= 2\int_{S^{n-1}}\int_0^\infty f(x + t\theta) dt d\theta$

Because in the integral defining $X^tX$, we  catch the integral over $t \in \mathbb{R}^{-}$ twice: once for $\theta$, and once for $-\theta$.  This can be reinterpreted as an integral in polar coordinates. Let $y = t\theta$ and $t = |y|$.  Then

$X^tXf(x) = c \int_{\mathbb{R}^n} \frac{f(x-y)}{|y|^{n-1}} dy$

$= f * \frac{1}{|x|^{n-1}}$

up to constants.  This is a singular integral, but not badly behaved.  We do not need to regularize it.

To invert this convolution, we use the Fourier transform.

$\mathfrak{F}(X^tXf)(\xi) = c\hat{f}(\xi)\widehat{\frac{1}{|x|^{n-1}}}(\xi)$

for $f \in C_0^\infty(\mathbb{R}^n)$.

Claim $\widehat{\frac{1}{|x|^{n-1}}}(\xi) = c_n|\xi|^{-1}$

We will prove this claim (in fact, we’ll prove a more general result) in a supporting post.  This result immediately lets us invert the normal operator. If

$AX^tXf(x) = f(x)$

then we must have

$\widehat{Af}(\xi) = |\xi|\hat{f}(\xi)$

In other words,

$Af = (-\bigtriangleup)^{\frac{1}{2}}$

is a power of the negative Laplacian. Now we can state the main result

Theorem [The X-Ray Inversion Formula] For all $f \in C_0^{\infty}(\mathbb{R}^n)$

$c(-\bigtriangleup)^{\frac{1}{2}}X^tXf = f$

So we gain half a derivative by applying $X$. Notice that this formula does not depend on the dimension, because the dimension of the manifolds we are integrating over to define $X$ do not change.  Also note that the inversion operator is never local.

We will see that for $n > 2$ this has too many lines, and the formula is not very useful.

Parallel Development for the Radon Transform

We could develop the inversion formula for the Radon transform using its normal operator in a completely analogous way.

$R^tRf(x) = c_n\int\frac{f(y)}{|x-y|} dy = c_n f * \frac{1}{|x|}$

So the same manipulations show that

Theorem [Radon Inversion Formula II] For $f \in C_0^{\infty}(\mathbb{R}^n)$

$(-\bigtriangleup)^{\frac{n-1}{2}}R^tRf = f$

If $n$ is even, we get a square root and the operator is not local.  If $n$ is odd, the inversion operator is local.

At first glance, this looks different, but it is really the same as what we had before.  The proof is left as an exercise.  It relies on the observation that the Radon transform intertwines (functions of the) Laplacian on $\mathbb{R}^n$ with those on $\mathbb{R}$:

$R(-\bigtriangleup)^\alpha = (\frac{\partial^2}{\partial s^2})^\alpha R$

In the coming lectures, we will

• Generalize these results to distributions
• Rewrite the inversion formula to make it look more like the formula for $R$.
• Make stability estimates
• Characterize the range of $X$

The range characterization is perhaps the most interesting part.  $Xf$ depends on $2n - 2$ variables and is overdetermined.  We can extend it to $\mathbb{R}^n \times \mathbb{R}^n$ as follows

$Xf(x, \xi) = \int f(x + t\xi) dt$

and immediately see that we have a compatibility condition:

$\frac{\partial^2 Xf}{\partial x_i\partial \xi_j} = \frac{\partial^2 Xf}{\partial x_j\partial \xi_i}$

We will see the Fritz John’s theorem that states that this system of PDEs gives necessary and sufficient conditions for a function to be in the range of $X$.

I will grant this: I don’t think computing $R^tR$ is as easy as computing $X^tX$ and if I were asked to do it on the spot, my first reaction would be to use a Fourier Transform. So the parallel seems less than perfect at the lower level steps.