doc/ctsim-concepts.tex

   1 \chapter{Concepts}\index{Concepts}%
   2 \setheader{{\it CHAPTER \thechapter}}{}{}{}{}{{\it CHAPTER \thechapter}}%
   3 \setfooter{\thepage}{}{}{}{}{\thepage}%
   4
   5 \section{Overview}\label{conceptoverview}\index{Concepts,Overview}%
   6 In \ctsim, a phantom object, or a geometrical description of the object
   7 of a CT study is constructed and an image can be created.  Then a
   8 scanner geometry can be specified, and the projection data simulated.
   9 Finally that projection data can be reconstructed using various user
  10 controlled algorithms producing an image of the phantom or study object.
  11
  12 In order to use \ctsim\ effectively, some knowledge of how \ctsim\ works
  13 and the approach taken is required. \ctsim\ deals with a variety of
  14 object, but the two objects we need to be concerned with are the 'phantom' and
  15 the 'scanner'.
  16
  17 \section{Phantoms}\label{conceptphantom}\index{Concepts,Phantoms}%
  18 \subsection{Overview}\label{phantomoverview}\index{Concepts,Phantoms,Overview}%
  19
  20 \ctsim\ uses geometrical objects to
  21 describe the object being scanned. A phantom is composed a one or more
  22 phantom elements. These elements are simple geometric shapes,
  23 specifically, rectangles, triangles, ellipses, sectors and segments.
  24 With these the standard phantoms used in the CT literature (the Herman
  25 and the Shepp-Logan) can be constructed.  In fact
  26 \ctsim\ provides a shortcut to construct those phantoms for you.  It also
  27 allows you to write a file in which the composition of your own phantom is
  28 described.
  29
  30 The types of phantom elements and their definitions are taken from Herman's 1980
  31 book\cite{HERMAN80}.
  32
  33 \subsection{Phantom File}\label{phantomfile}\index{Concepts,Phantoms,File}
  34 Each line in the text file describes an element of the
  35 phantom.  Each line contains seven entries, in the following form:
  36 \begin{verbatim}
  37 item cx cy dx dy r a
  38 \end{verbatim}
  39 The first entry defines the type of the element, one
  40 of {\tt rectangle}, {\tt ellipse}, {\tt triangle}, {\tt sector}, or {\tt segment}.
  41 {\tt cx}, {\tt cy}, {\tt dx} and {\tt dy} have different meanings depending on the element type.
  42
  43 {\tt r} is the rotation applied to the object in degrees counterclockwise,
  44 and {\tt a} is the X-ray attenuation coefficient of the object.
  45 Where objects overlap, the attenuations of the overlapped objects are summed.
  46
  47
  48
  49 \subsection{Phantom Elements}\label{phantomelements}\index{Concepts,Phantoms,Elements}
  50
  51 \subsubsection{ellipse}
  52 Ellipses use dx and dy to define the semi-major and semi-minor axis lengths,
  53 with the centre of the ellipse at cx and cy.  Of note, the commonly used
  54 phantom described by Shepp and Logan\cite{SHEPP74} uses only ellipses.
  55
  56 \subsubsection{rectangle}
  57 Rectangles use
  58 cx  and cy to define the position of the centre of the rectangle with respect
  59 to the origin.  dx and dy  are the half-width and half-height of the
  60 rectangle.
  61
  62 \subsubsection{triangle}
  63 Triangles are drawn with the centre of the base at cx,cy, with a base
  64 width of 2*dx in x direction, and a height of dy.  Rotations are then
  65 applied about the origin.
  66
  67 \subsubsection{sector}
  68 It appears that dx and dy
  69 define the end points of a radius of the sector, from which the radius and
  70 the angle of the two arms of the sector are calculated.  But then
  71 orientation and centering of the sector don't make much sense yet.
  72
  73 \subsubsection{segment}
  74 Segments are the segments of a circle between a chord and the
  75 perimeter of the circle.  This also isn't clear to me, but it appears that
  76 perhaps the distance from chord to circle perimeter, and circle radius is
  77 defined by dx and dy. Chord is always horizontal through the origin, then
  78 translated and then rotated (???).
  79
  80 \subsection{Phantom Size}
  81 Also note that the overall dimensions of the phantom are increased by 1\%
  82 above the specified sizes to avoid clipping due to round-off errors.
  83 So, if the phantom is defined as
  84 a rectangle of size 0.1 by 0.1, the actual phantom has extent $\pm$0.101 in
  85 each direction.
  86
  87 \section{Scanner}\label{conceptscanner}\index{Concepts,Scanner}%
  88 \subsection{Geometries}
  89 This is where things get tricky.  There are two possible approaches.  The
  90 simple approach would be to define the size of a phantom which is put at
  91 the centre of the scanner.  The scanner would have it's bore size defined,
  92 or perhaps better, the field of view defined. Here, field of view would be
  93 the radius or diameter of the circular area from which data is collected
  94 and an image reconstructed.  In a real CT scanner, if the object being
  95 scanned is larger than the field of view, you get image artifacts.  And of
  96 course you can't stuff an object into a scanner if the object is larger
  97 than the bore!  In this model, the scanner size or field of view would
  98 be used as the standard length scale.
  99
 100 However, \ctsim\ takes another approach.  I believe this approach arose
 101 because the "image" of the phantom produced from the phantom description
 102 was being matched to the reconstruction image of the phantom.  That is,
 103 the dimensions of the 'before' and 'after' images  were being matched.
 104 The code has a Phantom object and a Scanner object.  The geometry of the
 105 Scanner is defined in part by the properties of the Phantom.  In fact,
 106 all dimensions are determined in terms of the phantom size, which is used
 107 as the standard length scale. Remember, as mentioned above, the
 108 phantom dimensions are also padded by 1\%.
 109
 110 The maximum of the phantom length and height is used to define the square
 111 that completely surrounds the phantom. Let $p_l$ be the width (also height)
 112 of this square. The diameter of this boundary box, $p_d$ is then
 113 \latexonly{$$p_d = \sqrt{2}(p_l)$$}
 114 \latexignore{sqrt(2) * $p_l$.}
 115 This relationship can be seen in figure 1 with the parallel geometry.
 116
 117 \subsubsection{Focal Length \& Field of View}
 118 The two important variables is the focal-length-ratio $f_lr$.
 119 This is used along with $p_d$ to
 120 define the focal length according to
 121 \latexonly{$$f_l = f_lr p_d$$}
 122 \latexignore{\\$f_l$ = $f_lr$ x $p_d$\\}
 123 where
 124 we consider the case of a first generation parallel beam CT scanner.
 125
 126 \subsubsection{Parallel Geometry}\label{geometryparallel}\index{Concepts,Scanner,Geometries,Parallel}
 127 \begin{figure}
 128 \image{10cm;0cm}{ctsimfig1.eps}
 129 \caption{Geometry used for a 1st generation, parallel beam CT scanner}\label{fistgenfig}
 130 \end{figure}
 131
 132 In figure 1A, the excursion of the source and detector need only be $l_p$,
 133 the height (or width) of the phantom's bounding square. However, if the
 134 field of view were only $l_p$, then the projection shown in figure 1B
 135 would clip the corners of the phantom.  By increasing the field of view by
 136 $\sqrt{2}$ the whole phantom is included in every projection.  Of course,
 137 if the field-of-view ratio $f_{vR}$ is larger than 1, there is no problem.
 138 However, if $f_{vR}$ is less than one and thus the scanner is smaller than
 139 the phantom, then distortions will occur without warning from the program.
 140
 141 The code also sets the detector length equal to the field of view in this
 142 case.  The focal length is chosen to be $\sqrt{2}l_p$ so the phantom will
 143 fit between the source and detector at all rotation angles, when the focal
 144 length ratio is specified as 1.  Again, what happens if the focal length
 145 ratio is chosen less than 1?
 146
 147 The other thing to note is that in this code the detector array is set to
 148 be the same size as the field-of-view $f_v$, equation (2).  So, one has to
 149 know the size of the phantom to specify a given scanner geometry with a
 150 given source-detector distance (or $f_l$ here) and a given range of
 151 excursion ($f_v$ here).
 152
 153 \subsubsection{Divergent Geometries}\label{geometrydivergent}\index{Concepts,Scanner,Geometries,Divergent}
 154 Next consider the case of equilinear (second generation) and equiangular
 155 (third, fourth, and fifth generation) geometries.
 156 The parts of the code  relevant to this
 157 discussion are the same for both modes.  In the equilinear mode, a single
 158 source produces a fan beam which is read by a linear array of detectors.  If
 159 the detectors occupy an arc of a circle, then the geometry is equiangular.
 160 See figure 2.
 161 \begin{figure}
 162 \image{10cm;0cm}{ctsimfig2.eps}
 163 \caption{Equilinear and equiangular geometries.}
 164 \end{figure}
 165
 166 For these geometries, the following logic is executed:  A variable dHalfSquare
 167 $d_{hs}$ is defined as
 168 \latexonly{$$d_{hs} = (f_v)/(2\sqrt{2}) = (l_p/2) f_{vR}$$}
 169 This is then subtracted from the focal length $f_l$ as calculated above, and
 170 assigned to a new variable
 171 \latexonly{$\mathrm{dFocalPastPhm} = f_l - d_{hs}$}.  Since $f_l$ and
 172 $d_{hs}$ are derived from the phantom dimension and the input focal length and field of view ratios, one can write,
 173 \latexonly{$$\mathrm{dFocalPastPhm} = f_l -d_{hs}
 174     = \sqrt{2}(l_p/2) f_{lR} - (l_p/2) f_{vR} = l_p(\sqrt{2}f_{lR} - f_{vR})$$}
 175 If this quantity is less than or equal to zero, then at least for some
 176 projections  the source is inside the phantom.  Perhaps a figure will help at
 177 this point. Consider first the case where $f_{vR} = f_{lR} =1 $, figure 3. The
 178 square in the figure bounds the phantom and has sides $l_p$.  For this case
 179 then,
 180 \latexonly{$$f_l=\sqrt{2}l_p/2 = l_p/\sqrt{2}$$,
 181 $$f_v = \sqrt{2}l_p$$,
 182 and
 183 $$d_{hs} = {l_p}/{2}$$.
 184 Then
 185 $$\mathrm{dFocalPastPhm} = ({l_p}/{2}) (\sqrt{2}-1)$$
 186 }
 187 \begin{figure}
 188 \image{5cm;0cm}{ctsimfig3.eps}
 189 \caption{Equilinear and equiangluar geometry when focal length ratio =
 190 field of view ratio = 1.}
 191 \end{figure}
 192 The angle $\alpha$ is now defined as shown in figure 3, and the detector
 193 length is adjusted to subtend the angle $2\alpha$ as shown.  Note that the
 194 size of the detector array may have changed and the field of view is not
 195 used.
 196 For a circular array of detectors, the detectors are spaced around a
 197 circle covering an angular distance of $2\alpha$.  The dotted circle in
 198 figure 3 indicates the positions of the detectors in this case. Note that
 199 detectors at the ends of the range would not be illuminated by the source.
 200
 201 Now, consider increasing the focal length ratio to two leaving the
 202 field of view ratio as 1,  as in  Figure 4.  Now the detectors array is
 203 denser, and the real field of view is closer to that specified, but note
 204 again that the field of view is not used. Instead, the focal length is
 205 used to give a distance from the centre of the phantom to the source, and
 206 the detector array is adjusted to give an angular coverage to include the
 207 whole phantom.
 208 \begin{figure}
 209 \image{10cm;0cm}{ctsimfig4.eps}
 210 \caption{Equilinear and equiangluar geometry when focal length ratio = 2
 211 and the field of view ratio = 1.}
 212 \end{figure}
 213 \begin{figure}
 214 \image{10cm;0cm}{ctsimfig5.eps}
 215 \caption{Equilinear and equiangluar geometry when focal length ratio = 4.}
 216 \end{figure}
 217
 218 \section{Reconstruction}\label{conceptreconstruction}\index{Concepts,Reconstruction}%
 219 \subsection{Overview}
 220 \subsection{Filtered Backprojection}
 221 \subsection{Direct Inverse Fourier}
 222 This method is not currently implemented in \ctsim, however it is planned for a
 223 future release. This method does not give as accurate result as filtered
 224 backprojection mostly due to interpolation occuring in the frequency domain rather
 225 than the spatial domain. The technique is comprised of two sequential steps:
 226 filtering projections and then backprojecting the filtered projections. Though
 227 these two steps are sequential, each view position can be processed individually.
 228 This parallelism is exploited in the MPI versions of \ctsim\ where the data from
 229 all the views are spread about amongst all of the processors. This has been testing
 230 in a 16-CPU cluster with good results.
 231
 232 \subsubsection{Filter projections}
 233 The projections for a single view have their frequency data multipled by
 234 a filter of $|w|$. \ctsim\ permits four different ways to accomplish this
 235 filtering. Two of the methods use convolution of the projection data with the
 236 inverse fourier transform of $|w|$. The other two methods perform an fourier
 237 transform of the projection data and multiply that by the $|w|$ filter and
 238 then perform an inverse fourier transform.
 239
 240 Though multiplying by $|w|$ gives the sharpest reconstructions, in practice, superior results are obtained by mutiplying the $|w|$ filter by
 241 another filter that attenuates the higher frequencies. \ctsim\ has multiple
 242 filters for this purpose.
 243
 244 \subsubsection{Backprojection of filtered projections}
 245 Backprojection is the process of ``smearing'' the filtered projections over
 246 the reconstructing image.
 247