[ctsim.git] / doc / ctsim-concepts.tex

\chapter{Concepts}\index{Concepts}%
\setheader{{\it CHAPTER \thechapter}}{}{}{}{}{{\it CHAPTER \thechapter}}%
\setfooter{\thepage}{}{}{}{}{\thepage}%

\section{Overview}\label{conceptoverview}\index{Concepts,Overview}%
The operation of \ctsim\ begins with the phantom object.  A phantom
object consists of geometric elements.  A scanner is specified and the
projection data simulated.  Finally that projection data can be
reconstructed using various user controlled algorithms producing an
image of the phantom object. This reconstruction can then be
statistically compared to the original phantom object.

In order to use \ctsim\ effectively, some knowledge of how \ctsim\ works
and the approach taken is required. \ctsim\ deals with a variety of
object, but the two objects we need to be concerned with are the
\emph{phantom} and the \emph{scanner}.

\section{Phantoms}\label{conceptphantom}\index{Concepts,Phantoms}%
\subsection{Overview}\label{phantomoverview}\index{Concepts,Phantoms,Overview}%

\ctsim\ uses geometrical objects to
describe the object being scanned. A phantom is composed a one or more
phantom elements. These elements are simple geometric shapes,
specifically, rectangles, triangles, ellipses, sectors and segments.
With these elements, standard phantoms used in the CT literature can
be constructed.  In fact, \ctsim\ provides a shortcut to load the
published phantoms of Herman and Shepp-Logan.  \ctsim\ also reads text
files of user-defined phantoms.

The types of phantom elements and their definitions are taken with
permission from G.T. Herman's 1980 book\cite{HERMAN80}.

\subsection{Phantom File}\label{phantomfile}\index{Concepts,Phantoms,File}
Each line in the text file describes an element of the
phantom.  Each line contains seven entries, in the following form:
\begin{verbatim}
element-type cx cy dx dy r a
\end{verbatim}
The first entry defines the type of the element, one of
\rtfsp\texttt{rectangle}, \texttt{}, \texttt{triangle},
\rtfsp\texttt{sector}, or \texttt{segment}. \texttt{cx},
\rtfsp\texttt{cy}, \texttt{dx} and \texttt{dy} have different
meanings depending on the element type.

\rtfsp\texttt{r} is the rotation applied to the object in degrees
counterclockwise, and \texttt{a} is the X-ray attenuation
coefficient of the object. Where objects overlap, the attenuations
of the overlapped objects are summed.


\subsection{Phantom Elements}\label{phantomelements}\index{Concepts,Phantoms,Elements}

\subsubsection{ellipse}
Ellipses use \texttt{dx} and \texttt{dy} to define the semi-major and
semi-minor axis lengths, with the center of the ellipse at \texttt{cx}
and \texttt{cy}.  Of note, the commonly used phantom described by
Shepp and Logan\cite{SHEPP74} uses only ellipses.

\subsubsection{rectangle}
Rectangles use \texttt{cx} and \texttt{cy} to define the position of
the center of the rectangle with respect to the origin.  \texttt{dx}
and \texttt{dy} are the half-width and half-height of the rectangle.

\subsubsection{triangle}
Triangles are drawn with the center of the base at \texttt{(cx,cy)},
with a base half-width of \texttt{dx} and a height of \texttt{dy}.
Rotations are then applied about the center of the base.

\subsubsection{segment}
Segments are complex. They are the portion of an circle between a
chord and the perimeter of the circle.  \texttt{dy} sets the radius of
the circle. Segments start with the center of the chord located at
\texttt{(0,0)} and the chord horizontal. The half-width of the chord
is set by \texttt{dx}.  The portion of an circle lying below the chord
is then added. The imaginary center of this circle is located at
\texttt{(0,-dy)}. The segment is then rotated by \texttt{r} and then
translated by \texttt{cx,cy}.

\subsubsection{sector}
Sectors are the like a ``pie slice'' from a circle. The radius of the
circle is set by \texttt{dy}. Sectors are
defined similarly to segments. In this case, though, a chord is not
drawn.  Instead, the lines are drawn from the origin of the circle
\texttt{(0,-dy)} to the points \texttt{(-dx,0)} and \texttt{(dx,0)}.
The perimeter of the circle is then draw between those two points
below the x-axis. The sector is then rotated and translated the same
as a segment.

\subsection{Phantom Size}
The overall dimensions of the phantom are increased by 1\% above the
specified sizes to avoid clipping due to round-off errors from
sampling the polygons of the phantom elements.  So, if the phantom is
defined as a rectangle of size 0.1 by 0.1, the actual phantom has
extent 0.101 in each direction.

\section{Scanner}\label{conceptscanner}\index{Concepts,Scanner}%
\subsection{Dimensions}
Understanding the scanning geometry is the most complicated aspect of
using \ctsim. For real-world CT simulators, this is actually quite
simple. The geometry is fixed by the manufacturer during the
construction of the scanner and can not be changed. Conversely,
real-world CT scanners can only take objects up to a fixed size.

\ctsim, being a very flexible simulator,
gives tremendous options in setting up the geometry for a scan.

In general, the geometry for a scan all starts from the size of the
phantom being scanned. This is because \ctsim\ allows for statistical
comparisons between the original phantom image and it's reconstructions.
Since CT scanners scan a circular area, the first important
variable is the diameter of the circle surround the phantom, or the
\emph{phantom diameter}. Remember, as mentioned above, the
phantom dimensions are also padded by 1\%.

The other important geometry variables for scanning objects are the
\emph{view ratio}, \emph{scan ratio}, and \emph{focal length ratio}.
These variables are all input into \ctsim\ in terms of ratios rather
than absolute values.

\subsubsection{Phantom Diameter}
\begin{figure}
$$\image{5cm;0cm}{scangeometry.eps}$$
\caption{Phantom Geometry}
\end{figure}
The phantom diameter is automatically calculated by \ctsim\ from the
phantom definition. The maximum of the phantom length and height is
used to define the square that completely surrounds the phantom. Let
\latexonly{$p_l$}\latexignore{\emph{Pl}}
be the width and height of this square. The diameter of this boundary box,
\latexonly{$p_d$,}\latexignore{\emph{Pd},}
\rtfsp is then
\latexignore{\\$$\emph{Pl x sqrt(2)}$$\\}
\latexonly{$$p_d = p_l \sqrt{2}$$}
CT scanners actually collect projections around a circle rather than a
square. The diameter of this circle is also the diameter of the boundary
square
\latexonly{$p_d$.}\latexignore{\rtfsp\emph{Pd}.}
These relationships are diagrammed in figure 2.1.

\subsubsection{View Diameter}
The \emph{view diameter} is the area that is being processed during scanning of phantoms as
well as during rasterization of phantoms. By default, the \emph{view diameter}
\rtfsp is set equal to the \emph{phantom diameter}. It may be useful, especially for
experimental reasons, to process an area larger (and maybe even smaller) than
the phantom. Thus, during rasterization or during projections, \ctsim\ will
ask for a \emph{view ratio},
\latexonly{$v_r$.}\latexignore{\rtfsp \emph{VR}.}
The \emph{view diameter} is then set as
\latexonly{$$v_d = p_d v_r$$}\latexignore{\\$$\emph{Vd = Pd x VR}$$}

By using a
\latexonly{$v_r$}\latexignore{\emph{VR}}
less than 1, \ctsim\ will allow
for a \emph{view diameter} less than
\emph{phantom diameter}.
This will lead to significant artifacts. Physically, this would
be impossible and is analagous to inserting an object into the CT
scanner that is larger than the scanner itself!

\subsubsection{Scan Diameter}
By default, the entire \emph{view diameter} is scanned. For experimental
purposes, it may be desirable to scan an area either larger or smaller than
the \emph{view diameter}. Thus, the concept of \emph{scan ratio}
\latexonly{$s_r$}\latexignore{\emph{SR}}
is born. The scan diameter
\latexonly{$s_d$}\latexignore{\emph{Sd}}
is the diameter over which x-rays are collected and is defined as
\latexonly{$$s_d = v_d s_r$$}\latexignore{\\$$\emph{Sd = Vd x SR}$$\\}
By default and for all ordinary scanning, the \emph{scan ratio} is to
\texttt{1}. If the \emph{scan ratio} is less than \texttt{1},
you can expect significant artifacts.

\subsubsection{Focal Length}
The \emph{focal length},
\latexonly{$f$,}\latexignore{\emph{F},}
is the distance of the X-ray source to the center of
the phantom. The focal length is set as a ratio,
\latexonly{$f_r$,}\latexignore{\emph{FR},}
of the view radius. Focal length is
calculated as
\latexonly{$$f = (v_d / 2) f_r$$}\latexignore{\\$$\emph{F = (Vd / 2) x FR}$$}

For parallel geometry scanning, the focal length doesn't matter. However,
divergent geometry scanning (equilinear and equiangular), the \emph{focal
length ratio} should be set at \texttt{2} or more to avoid artifacts.


\subsection{Parallel Geometry}\label{geometryparallel}\index{Concepts,Scanner,Geometries,Parallel}

As mentioned above, the focal length is not used in this simple
geometry. The detector array is set to be the same size as the
\emph{scan diameter}.  For optimal scanning in this geometry, the
\emph{scan diameter} should be equal to the \emph{phantom
diameter}. This is accomplished by using the default values of
\texttt{1} for the \emph{view ratio} and the \emph{scan ratio}. If
values of less than \texttt{1} are used for these two variables,
significant distortions will occur.


\subsection{Divergent Geometries}\label{geometrydivergent}\index{Concepts,Scanner,Geometries,Divergent}
\subsubsection{Overview}
Next consider the case of equilinear (second generation) and equiangular
(third, fourth, and fifth generation) geometries. In these cases,
the x-ray beams diverge from a single source to a detector array.
In the equilinear mode, a single
source produces a fan beam which is read by a linear array of detectors.  If
the detectors occupy an arc of a circle, then the geometry is equiangular.
See figure 2.2.
\begin{figure}
\image{10cm;0cm}{divergent.eps}
\caption{Equilinear and equiangular geometries.}
\end{figure}


\subsubsection{Fan Beam Angle}
For these divergent beam geometries, the \emph{fan beam angle} needs
to be calculated. For real-world CT scanners, this is fixed at the
time of manufacture. \ctsim, however, calculates the \emph{fan beam angle},
$\alpha$ from the \emph{scan diameter} and the \emph{focal length}
\latexignore{\\$$\emph{alpha = 2 x asin ( (Sd / 2) / f)}$$\\}
\latexonly{\begin{equation}\label{alphacalc}\alpha = 2 \sin^{-1} ((s_d / 2) / f)\end{equation}}
This is illustrated in figure 2.3.
\begin{figure}
\image{10cm;0cm}{alphacalc.eps}
\caption{Calculation of $\alpha$}
\end{figure}


Empiric testing with \ctsim\ shows that for very large \emph{fan beam angles},
greater than approximately
\latexonly{$120^\circ$,}\latexignore{120 degrees,}
there are significant artifacts. The primary way to manage the
\emph{fan beam angle} is by varying the \emph{focal length} since the
\emph{scan diameter} by the size of the phantom.

To illustrate, the \emph{scan diameter} can be defined as
\latexonly{$$s_d = s_r v_r p_d$$}\latexignore{\\$$Sd = Sr x Vr x Pd$$\\}

Further, $f$ can be defined as
\latexonly{$$f = f_r (v_r p_d / 2)$$}
Plugging these equations into
\latexignore{the above equation,}\latexonly{equation~\ref{alphacalc},}
We have,
\latexonly{
\begin{eqnarray}
\alpha &= 2\,\sin^{-1} \frac{s_r v_r p_d / 2}{f_r v_r (p_d / 2)} \nonumber \\
&= 2\,\sin^{-1} (s_r / f_r)
\end{eqnarray}
}

Since in normal scanning $s_r = 1$, $\alpha$ depends only upon the \emph{focal length ratio}.

\subsubsection{Detector Array Size}
In general, you do not need to be concerned with the detector array
size. It is automatically calculated by \ctsim.

For parallel geometry, the detector length is equal to the scan
diameter.

For divergent beam geometrys, the size of the
detector array also depends upon the \emph{focal length}.
Increasing the \emph{focal length}
decreases the size of the detector array while increasing the \emph{scan
diameter} increases the detector array size.

For equiangular geometry, the detectors are spaced around a
circle covering an angular distance of
\latexonly{$\alpha$.}\latexignore{\emph{alpha}.}
The dotted circle in
\begin{figure}
\image{10cm;0cm}{equiangular.eps}
\caption{Equiangluar geometry}
\end{figure}
figure 2.4 indicates the positions of the detectors in this case.

For equilinear geometry, the detectors are space along a straight
line. The length of the line depends upon
\latexonly{$\alpha$}\latexignore{\emph{alpha}}
and the \emph{focal length}. It is calculated as
\latexonly{$$\mathrm{detLengh} = 4\,f \tan (\alpha / 2)$$}
\latexignore{\\$$\emph{detLength} = 4 x F x tan(alpha/2)$$\\}
\begin{figure}
\image{10cm;0cm}{equilinear.eps}
\caption{Equilinear geometry}
\end{figure}
An example of the this geometry is in figure 2.5.


\subsubsection{Examples of Geometry Settings}


\section{Reconstruction}\label{conceptreconstruction}\index{Concepts,Reconstruction}%
\subsection{Overview}
\subsection{Direct Inverse Fourier}
This method is not currently implemented in \ctsim, however it is
planned for a future release. This method does not give results as
accurate as filtered backprojection. The difference is due primarily
because interpolation occurs in the frequency domain rather than the
spatial domain.

\subsection{Filtered Backprojection}
The technique is comprised of two sequential steps:
filtering projections and then backprojecting the filtered projections. Though
these two steps are sequential, each view position can be processed individually.

\subsubsection{Multiple Computer Processing}
This parallelism is exploited in the MPI versions of \ctsim\ where the
data from all the views are spread about amongst all of the
processors. This has been testing in a 16-CPU cluster with good
results.

\subsubsection{Filter projections}
The projections for a single view have their frequency data multipled by
a filter of $|w|$. \ctsim\ permits four different ways to accomplish this
filtering.

Two of the methods use convolution of the projection data with the
inverse Fourier transform of $|w|$. The other two methods perform an Fourier
transform of the projection data and multiply that by the $|w|$ filter and
then perform an inverse fourier transform.

Though multiplying by $|w|$ gives the sharpest reconstructions, in
practice, superior results are obtained by mutiplying the $|w|$ filter
by another filter that attenuates the higher frequencies. \ctsim\ has
multiple filters for this purpose.

\subsubsection{Backprojection of filtered projections}
Backprojection is the process of ``smearing'' the filtered projections
over the reconstructing image. Various levels of interpolation can be
specified.  In general, the trade-off is between quality and execution
time.