Line 66: | Line 66: | ||
&= \int_0^{\pi}\underbrace{[\int_{-\infty}^{\infty}|\rho|P_{\theta}e^{2\pi j\rho(x\cos(\theta) +y\sin(\theta))}d\rho]}_{g_{\theta}(x\cos(\theta) + y\sin(\theta))}d\theta | &= \int_0^{\pi}\underbrace{[\int_{-\infty}^{\infty}|\rho|P_{\theta}e^{2\pi j\rho(x\cos(\theta) +y\sin(\theta))}d\rho]}_{g_{\theta}(x\cos(\theta) + y\sin(\theta))}d\theta | ||
\end{align} | \end{align} | ||
− | </math><br /> | + | </math> |
− | Pulling out <math>g_{\theta}(t) | + | <br /> |
+ | Pulling out <math>g_{\theta}(t)</math>,<br /> | ||
<math>\begin{align} | <math>\begin{align} | ||
g_{\theta}(t) | g_{\theta}(t) | ||
− | &= \int_{\infty}^{\infty}|\rho | P_{\theta}(\rho )e^{2\pi j\rho t} d\rho \\ | + | &= \int_{-\infty}^{\infty}|\rho | P_{\theta}(\rho )e^{2\pi j\rho t} d\rho \ \text{, by comparison with the CTFT}\\ |
− | &= CTFT^{-1}\{|\rho |P_{\theta}(\rho)\} | + | &= CTFT^{-1}\{|\rho |P_{\theta}(\rho)\} \\ |
&= h(t) * p_{\theta}(r) \\ | &= h(t) * p_{\theta}(r) \\ | ||
\end{align}</math> | \end{align}</math> | ||
− | + | <br /> | |
+ | ---- | ||
=Projection Filter Analysis= | =Projection Filter Analysis= | ||
Now let's focus on <math>g_{\theta}(r)</math>.<br /> | Now let's focus on <math>g_{\theta}(r)</math>.<br /> | ||
Line 89: | Line 91: | ||
h(r) &= f_c^2[2sinc(2f_c t)-sinc^2(tf_c)] | h(r) &= f_c^2[2sinc(2f_c t)-sinc^2(tf_c)] | ||
\end{align}</math> | \end{align}</math> | ||
+ | ---- | ||
=Back Projection Analysis= | =Back Projection Analysis= | ||
− | + | As mentioned above, the back projection function is<br /> | |
+ | <math>f(x,y) = \int_{0}^{\pi} b | ||
+ | For each angle <math>\theta</math>, the back | ||
---- | ---- |
Revision as of 13:44, 21 December 2014
Convolution/Fourier Back Projection Algorithm
A slecture by ECE student Sahil Sanghani
Partly based on the ECE 637 material of Professor Bouman.
Contents
Introduction
Convolution Back Projection (CBP) offers a reconstruction method that is not computationally expensive. Although the method is based on the Fourier Slice Theorem, there's never a transformation to the frequency domain. Theoretically CBP involves extruding every projection back through the origin and then summing the results. This operation is called back projection. The projections in Figure 1 were all assumed to be the same regardless of $ \theta $. In Figure 1a, the extrusion is demonstrated. In Figure 1b, the summing is demonstrated.
Summary
There are 3 steps to reconstructing an object from its projections using CBP.
1. Measure the projections $ p_{\theta}(r) $. (i.e. using CT)
2. Filter the projections $ g_{\theta}(r) = h(r) * p_{\theta}(r) $
3. Back project filtered projections
$ f(x,y) = \int_0^{\pi}g_{\theta}(x\cos(\theta)+x\sin(\theta))d\theta $
An infinite number of filtered back projections will result in a perfect reconstruction of the original image of the object given it is band limited. Since it is impossible to take that many projections, a good practice is to take at least $ n $ back projections for a $ n $ by $ n $ image.
Derivation
From the Fourier Slice Theorem, we get the following relationships.
$ \begin{align} u &= \rho\cos(\theta) \\ v &= \rho\sin(\theta) \end{align} $
Now let's calculate the Jacobian of the polar coordinate transformation.
$ dudv = \left | \frac{\partial (u,v)}{\partial (\theta,\rho)} \right \vert d\theta d\rho $
$ \begin{align} \left | \frac{\partial (u,v)}{\partial (\theta,\rho)} \right \vert &= det\begin{bmatrix} \frac{\partial u}{\partial \theta} & \frac{\partial u}{\partial \rho} \\ \frac{\partial v}{\partial \theta} & \frac{\partial u}{\partial \rho} \end{bmatrix} \\ &= det\begin{bmatrix} \frac{\partial (\rho\cos(\theta))}{\partial \theta} & \frac{\partial (\rho\cos(\theta))}{\partial \rho} \\ \frac{\partial (\rho\sin(\theta))}{\partial \theta} & \frac{\partial (\rho\sin(\theta))}{\partial \rho} \end{bmatrix} \\ &= det\begin{bmatrix} -\rho\sin(\theta) & \cos(\theta) \\ \rho\cos(\theta) & \sin(\theta) \end{bmatrix} \\ &= |-\rho\sin^{2}(\theta) - \rho\cos^{2}(\theta)| \\ &= |-\rho(\sin^{2}(\theta) + \cos^{2}(\theta))| \\ &= |-\rho| \\ &= |\rho| \\ \end{align} $
$ \Rightarrow dudv = |\rho|d\theta d\rho $
Starting with the formula for the inverse CSFT and using the Fourier Slice theorem, we will end up with
$ \begin{align} f(x,y) &= \int_{-\infty}^{\infty}\int_{-\infty}^{\infty}F(u,v)e^{2\pi j(xu+yv)}dudv \\ &= \int_{-\infty}^{\infty}\int_{0}^{\pi}F(\rho\cos(\theta),\rho\sin(\theta))e^{2\pi j(x\rho\cos(\theta) +y\rho\sin(\theta))}|\rho|d\theta d\rho \\ &= \int_0^{\pi}\underbrace{[\int_{-\infty}^{\infty}|\rho|P_{\theta}e^{2\pi j\rho(x\cos(\theta) +y\sin(\theta))}d\rho]}_{g_{\theta}(x\cos(\theta) + y\sin(\theta))}d\theta \end{align} $
Pulling out $ g_{\theta}(t) $,
$ \begin{align} g_{\theta}(t) &= \int_{-\infty}^{\infty}|\rho | P_{\theta}(\rho )e^{2\pi j\rho t} d\rho \ \text{, by comparison with the CTFT}\\ &= CTFT^{-1}\{|\rho |P_{\theta}(\rho)\} \\ &= h(t) * p_{\theta}(r) \\ \end{align} $
Projection Filter Analysis
Now let's focus on $ g_{\theta}(r) $.
$ g_{\theta}(r) = h(r) * p_{\theta}(r) \ $
The frequency response of the filter is
$ H(\rho) = CTFT{(h(r))} = |\rho| \ $
After graphing the frequency response, it is apparent that the filter is a high pass filter.
This filter can be represented by a rect function minus a triangle function, so its equation is
$ \begin{align} H(\rho) &= f_c[rect(f/(2f_c))-\wedge(f/f_c)] \\ h(r) &= f_c^2[2sinc(2f_c t)-sinc^2(tf_c)] \end{align} $
Back Projection Analysis
As mentioned above, the back projection function is
$ f(x,y) = \int_{0}^{\pi} b For each angle <math>\theta $, the back
References:
[1] C. A. Bouman. ECE 637. Class Lecture. Digital Image Processing I. Faculty of Electrical Engineering, Purdue University. Spring 2013.