Skip to content
Snippets Groups Projects
Commit 68e5a601 authored by Erik Strand's avatar Erik Strand
Browse files

Add filtered back projection results

parent 3ad32bbd
Branches
No related tags found
No related merge requests found
assets/img/project_fbp_reconstruction.png

259 KiB

...@@ -32,6 +32,7 @@ projections at different angles. ...@@ -32,6 +32,7 @@ projections at different angles.
## Fourier Reconstruction ## Fourier Reconstruction
### Fourier Slice Theorem ### Fourier Slice Theorem
Let $$f : \mathbb{R}^2 \rightarrow \mathbb{R}$$ be a density function. If we map density to Let $$f : \mathbb{R}^2 \rightarrow \mathbb{R}$$ be a density function. If we map density to
...@@ -111,7 +112,7 @@ algorithms out there that come closer to the theoretically ideal ...@@ -111,7 +112,7 @@ algorithms out there that come closer to the theoretically ideal
The popular [gridrec](https://www.ncbi.nlm.nih.gov/pubmed/23093766) method is one. The popular [gridrec](https://www.ncbi.nlm.nih.gov/pubmed/23093766) method is one.
### Implementation ### Results
I used [FFTW](http://www.fftw.org/) to compute Fourier transforms. It's written in C and is very I used [FFTW](http://www.fftw.org/) to compute Fourier transforms. It's written in C and is very
[fast](http://www.fftw.org/benchfft/). I implemented my own polar resampling routine. It uses a [fast](http://www.fftw.org/benchfft/). I implemented my own polar resampling routine. It uses a
...@@ -148,6 +149,7 @@ It would be nice if we could avoid the interpolation required for Fourier recons ...@@ -148,6 +149,7 @@ It would be nice if we could avoid the interpolation required for Fourier recons
simplest way of doing so, and still the most popular method of performing image reconstruction, is simplest way of doing so, and still the most popular method of performing image reconstruction, is
called filtered back projection. called filtered back projection.
### Theory ### Theory
Let's hop back to the continuous theory for a moment. The Fourier reconstruction technique is based Let's hop back to the continuous theory for a moment. The Fourier reconstruction technique is based
...@@ -189,6 +191,11 @@ $$ ...@@ -189,6 +191,11 @@ $$
f(x, y) = \int_0^\pi \mathcal{F}^{-1}[q_\theta](x \cos \theta + y \sin \theta) \mathrm{d} \theta f(x, y) = \int_0^\pi \mathcal{F}^{-1}[q_\theta](x \cos \theta + y \sin \theta) \mathrm{d} \theta
$$ $$
Since multiplication in the frequency domain is equivalent to convolution in the spatial domain,
$$\mathcal{F}^{-1}[q_\theta]$$ is simply a filtered version of $$p_\theta$$. Hence the name filtered
back projection.
### Discretization ### Discretization
We'll replace the continuous Fourier transforms with DFTs as before. We'll replace the continuous Fourier transforms with DFTs as before.
...@@ -210,6 +217,13 @@ them directly for quadrature. For $$r$$, on the other hand, the regular samples ...@@ -210,6 +217,13 @@ them directly for quadrature. For $$r$$, on the other hand, the regular samples
line up with the values $$x \cos \theta + y \sin \theta$$ that we want. So we'll still have to do line up with the values $$x \cos \theta + y \sin \theta$$ that we want. So we'll still have to do
some interpolation. But now it's a simple 1D interpolation problem that's much easier to do without some interpolation. But now it's a simple 1D interpolation problem that's much easier to do without
introducing as much error. In particular, we only have to do interpolation in the spatial domain, so introducing as much error. In particular, we only have to do interpolation in the spatial domain, so
errors will only accumulate locally. [Lanczos errors will only accumulate locally.
resampling](https://en.wikipedia.org/wiki/Lanczos_resampling) would probably be the best, but I'll
just use Catmull-Rom again since I already have it implemented.
### Results
I again use [FFTW](http://www.fftw.org/) for Fourier transforms. For interpolation I just take the
sample to the left of the desired location, and for integration I just use a left Riemann sum. Even
with these lazy techniques, the reconstruction looks quite good.
![brain](../assets/img/project_fbp_reconstruction.png)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment