r497: no message

[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}%
In \ctsim, a phantom object, or a geometrical description of the object
of a CT study is constructed and an image can be created.  Then a
scanner geometry can be specified, and the projection data simulated.
Finally that projection data can be reconstructed using various user
controlled algorithms producing an image of the phantom or study 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 'phantom' and
the '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 the standard phantoms used in the CT literature (the Herman
and the Shepp-Logan) can be constructed.  In fact
\ctsim\ provides a shortcut to construct those phantoms for you.  It also
allows you to write a file in which the composition of your own phantom is
described.

The types of phantom elements and their definitions are taken from 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}
item cx cy dx dy r a
\end{verbatim}
The first entry defines the type of the element, one
of {\tt rectangle}, {\tt ellipse}, {\tt triangle}, {\tt sector}, or {\tt segment}.
{\tt cx}, {\tt cy}, {\tt dx} and {\tt dy} have different meanings depending on the element type.

{\tt r} is the rotation applied to the object in degrees counterclockwise,
and {\tt 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 dx and dy to define the semi-major and semi-minor axis lengths,
with the centre of the ellipse at cx and cy.  Of note, the commonly used
phantom described by Shepp and Logan\cite{SHEPP74} uses only ellipses.

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

\subsubsection{triangle}
Triangles are drawn with the centre of the base at cx,cy, with a base
width of 2*dx in x direction, and a height of dy.  Rotations are then
applied about the origin.

\subsubsection{sector}
It appears that dx and dy
define the end points of a radius of the sector, from which the radius and
the angle of the two arms of the sector are calculated.  But then
orientation and centering of the sector don't make much sense yet.

\subsubsection{segment}
Segments are the segments of a circle between a chord and the
perimeter of the circle.  This also isn't clear to me, but it appears that
perhaps the distance from chord to circle perimeter, and circle radius is
defined by dx and dy. Chord is always horizontal through the origin, then
translated and then rotated (???).

\subsection{Phantom Size}
Also note that the overall dimensions of the phantom are increased by 1\%
above the specified sizes to avoid clipping due to round-off errors.  
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{Sizes}
Understanding the scanning geometry is the most complicated aspect
of using \ctsim. For our 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. 
\ctsim, being a very flexible simulator,
gives tremendous options is 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 diameter}, \emph{scan diameter}, and \emph{focal length}.
These variables are all input into \ctsim\ in terms of ratios rather than
absolute values.

\subsubsection{Phantom Diameter}
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},}
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{\emph{Pd}.}
These relationships are diagrammed in figure 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}
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 diameter ratio},
\latexonly{$V_{dR}$.}\latexignore{\emph{VdR}.}
The \emph{view diameter} is then set as
\latexonly{$$V_d = P_d V_{dR}$$}\latexignore{\\$$\emph{Vd = Pd x VdR}$$}

By using a 
\latexonly{$V_{dR}$}\latexignore{\emph{VdR}}
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 diameter}
\latexonly{$S_{dR}$}\latexignore{\emph{SdR}}
is born. The scan diameter
\latexonly{$S_d$}\latexignore{\emph{Sd}}
is defined as
\latexonly{$$S_d = V_d S_{dR}$$}\latexignore{\\$$\emph{Sd = Vd x SdR}$$\\}
By default and for all ordinary scanning, the \emph{scan diameter ratio} is to \texttt{1}. If the \emph{scan diameter ratio} is less than \texttt{1}, you
can plan on significant artifacts.

\subsubsection{Focal Length}
The \emph{focal length},
\latexonly{$F_l$,}\latexignore{\emph{Fl},}
is the distance of the X-ray source to the center of
the phantom. The focal length is set as a ratio,
\latexonly{$F_{lR}$,}\latexignore{\emph{FlR},}
of the view radius. Focal length is
calculated as
\latexonly{$$F_l = F_{lR} (V_d / 2)$$}\latexignore{\\$$\emph{Fl = FlR x (Vd / 2)}$$}

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}
\begin{figure}
\image{10cm;0cm}{ctsimfig1.eps}
\caption{Geometry used for a 1st generation, parallel beam CT scanner}\label{fistgenfig}
\end{figure}

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 diameter ratio} and
the \emph{scan diameter 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.
\begin{figure}
\image{10cm;0cm}{ctsimfig2.eps}
\caption{Equilinear and equiangular geometries.}
\end{figure}

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

If this quantity is less than or equal to zero, then at least for some
projections  the source is inside the phantom.  Perhaps a figure will help at
this point. Consider first the case where $f_{vR} = f_{lR} =1 $, figure 3. The
square in the figure bounds the phantom and has sides $l_p$.  For this case
then,
\latexonly{$$f_l=\sqrt{2}l_p/2 = l_p/\sqrt{2}$$,
$$f_v = \sqrt{2}l_p$$,
and
$$d_{hs} = {l_p}/{2}$$.
Then
$$\mathrm{dFocalPastPhm} = ({l_p}/{2}) (\sqrt{2}-1)$$
}
\begin{figure}
\image{5cm;0cm}{ctsimfig3.eps}
\caption{Equilinear and equiangluar geometry when focal length ratio =
field of view ratio = 1.}
\end{figure}
The angle $\alpha$ is now defined as shown in figure 3, and the detector
length is adjusted to subtend the angle $2\alpha$ as shown.  Note that the
size of the detector array may have changed and the field of view is not
used.
For a circular array of detectors, the detectors are spaced around a
circle covering an angular distance of $2\alpha$.  The dotted circle in
figure 3 indicates the positions of the detectors in this case. Note that
detectors at the ends of the range would not be illuminated by the source.

Now, consider increasing the focal length ratio to two leaving the
field of view ratio as 1,  as in  Figure 4.  Now the detectors array is
denser, and the real field of view is closer to that specified, but note
again that the field of view is not used. Instead, the focal length is
used to give a distance from the centre of the phantom to the source, and
the detector array is adjusted to give an angular coverage to include the
whole phantom.
\begin{figure}
\image{10cm;0cm}{ctsimfig4.eps}
\caption{Equilinear and equiangluar geometry when focal length ratio = 2
and the field of view ratio = 1.}
\end{figure}
\begin{figure}
\image{10cm;0cm}{ctsimfig5.eps}
\caption{Equilinear and equiangluar geometry when focal length ratio = 4.}
\end{figure}

\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.
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.