A slice of TOAST

---

TOAST is a software package for finite element modelling and image reconstruction in optical tomography written by Martin Schweiger and Simon Arridge. It includes some other programs for creating meshes and visualisation. For more information and versions of toast which can be downloaded, see Toast's homepage.

This tutorial assumes that you have a working version of toast installed on a PC running Linux. The instructions would be similar for a Sun. Instructions for compiling and installing toast are available.

Anything written like

this

is a command which should be typed onto the command line.

---

Contents

• Getting Started
• Forward Model
• Image Reconstruction
• Difference imaging
• 3D modelling and reconstruction
• Problems?

---

Getting Started

• Setting up your environment

  If you're running windows, you need to log on to a unix computer using these instructions.

  You need to be running C shell. If you're not sure, type

  csh

  Now you need to find out where toast is running from

  which toast

  This should give you the path of the toast reconstruction program, something like /home/toast/bin/toast. You need the path for the toast directory which is everything up to, but not including /bin/. Here it would be /home/toast. We'll call this path /my/toast/directory. If I say type /my/toast/directory, replace /my/toast/directory by the path you've just found. Check you've got it right

  ls /my/toast/directory/bin

  This should list the contents of the toast home directory. It should give you a list of directories (like bin and docs) and files (including toast.env). Now you can set up your toast environment by typing

  source /my/toast/directory/toast.env

  Now you should be able to type

  echo $TOASTDIR

  and it should print the path you originally found back to your terminal.

  If you want this to be done automatically every time you log in, then add this to your .cshrc file, which should be in your home directory. It's a hidden file so to see it, you need to type "ls -a":

  setenv TOASTDIR /my/toast/directory
  source $TOASTDIR/toast.env
  

• Make a mesh

  First, we'll make a circular finite element mesh, using the program meshmod, which comes with toast.

  meshmod

  Select Mesh > New > Circle. Change the radius to 35 mm, and select 12 element rings, 10 sectors and one boundary layer. Save the new mesh as circle.msh using Mesh > Save As. Exit meshmod.

  You now have a mesh. Open the file using any text editor e.g emacs, and look at it. The file format is described here.

• Process the mesh

  This mesh probably won't have enough nodes and elements to solve the problem accurately. We'll convert it to a quadratic mesh. This means that toast will fit a quadratic function inside each element rather than a linear function. Type

  lin2quadmesh < circle.msh > circle_quad.msh

  This will have created a new quadratic mesh called circle_quad.msh from the original linear mesh.

  The reconstruction using this mesh can be made to run much faster by optimising the mesh. This just renumbers the nodes and elements without changing the geometry. Load your quadratic mesh into meshmod.

  meshmod circle_quad.msh

  A window should open, looking something like this:

  

  This screen gives a lot of important information about the mesh. The name, dimension, elements and nodes should be self-explanatory. BBox gives the size of the smallest cube-like which encompasses the mesh (in this case it's a 70 x 70 mm square, as it's holding a circle of radius 35 mm. Mesh size is the summed total of all the element areas. Element size shows the largest and smallest element areas. If the range from the smallest to the largest element size is too big, the mesh might be bad. If any elements have negative area, the mesh is obviously bad. Sysmat SBW gives the bandwidth of the system matrix.

  Now we'll optimise the mesh. Select Edit > Optimise and highlight "minimum degree". Then click "Apply" and "Done".

  Then we'll put a blob in the mesh. First add one large blob, bigger than the mesh to initialise the properties of the entire mesh. Click on Edit > Insert Blob and insert a circle centred on 0 0 with radius 50 mm, with mua = 0.01, mus=1, n=1.4, region 0. Click on "Apply" and "Done". Now we'll add a blob. With the Insert Blob window still open, insert a circle centred on 10 10 with radius 5 mm, with mua = 0.05, mus=5, n=1.4, region 1. Finish with "Apply" and "Done". Look in the command window at the number of nodes in the blob you created. This example contains 60 nodes. If the blob has less than about 10 nodes, then you need a finer mesh.

  Save the mesh "Mesh > Save As". Save it as circle_quad_blobs.opt. We usually call unoptimised meshes ".msh" and optimised meshes ".opt".

  Exit from meshmod

• Make a qm file

  A qm file, defined here, contains the positions of all the sources and detector and defines what measurements were made from them. You can make a qm file suitable for the 2D mesh which you've created using

  makeqm

  You can either get the program to read in a mesh or assume a circular outline. Because we have a circluar mesh, we can do either. Type "2" to assume a circular outline. Enter 35 for the radius. Give a filename of "circle.qm". 16 sources, with zero angular offset and point profiles for source and detector. Again, open up your qm file in a text editor and try to understand the format.

---

Forward model

• Definition file

  The options for the forward modelling and reconstruction in toast are kept in definition files. Download this definition file for the forward modelling. Save it as a file called "fwd.def". The format of the .def file is described in detail here. Look at the definition file and the different options available.

  Open up the fwd.def and change the names of the qmfile and mesh so they match those of the qmfile and mesh which you generated above.

• Forward model

  Once the qmfile, mesh and .def file are all correct, you can simulate forward data.

  femdata -o fwd.def -I &

  This line calls the program femdata and tells it to use the definition file called fwd.def. -I makes femdata save the internal fields of light intensity or meantime (the sensitivity maps). Femdata's progress is saved to a log file (named in the definition file) which can be viewed using

  tail -f test_fwd.log

  This will say whether there are any errors and after only a second or two, depending on your computer, it should say "===== BYE =====". femdata will have generated four files, two data files one each for intensity and mean time and two image files containing the forward fields.

  Look at the dat files. Open one of the files (probably .int, .mean or .dat, depending on the filename in your .def file) in a text editor. Its format is described here.

• View forward fields

  The forward fields are saved as toast images in the .nim format. They can be viewed using the toast program toastim. Load in the intensity fields.

  toastim fwdsim.nim

  You'll need to adjust the display using Image > Mapping > Logarithmic to view the images on a log scale and Image > Colour > RGB Colours from File. Then you should get something like this:

    

  This shows the toastim window, with some useful information, and an image of the log of the light intensity due to the last source. The image is clearly distorted due to the blob you inserted. To see the other sources, scroll using the arrows on the top right of the toastim window. Do the same thing with the meantime data.

  Note. It is important at this point to check the quality of the solution. Any negative values for intensity or meantime mean that there is an error in the solution, most probably caused by the mesh being too coarse. Read the two image files into a text editor and search for negative values. You will not be able to view the log of the intensity field, as described above, if there are negative values in the solution. See what happens when you use a bad mesh - repeat everything you've done so far, using a coarser mesh.

---

Image reconstruction

• Definition file

  Use this definition file for reconstruction. It has a few more options selected in the [INVERSE_SOLVER] section and the mesh is reset to being homogeneous in [INIT_PARAM]. This is equivalent to assuming that the background properties are known, but we have no information about the perturbation.

• Reconstruction

  A successful forward model suggests that the qmfile, mesh and .def file are correct and agree with each other. It is not worth starting a reconstruction unless you have already got reasonable results from a forward simulation.

  To perform a reconstruction, run

  toast -o recon.def &

  and again follow what's happening using

  tail -f test_recon.log

  This should continually update the screen, showing the current status of the reconstruction. It should take a few minutes.

• View solution

  As before, review the solution using toastim. Scrolling with the arrows shows images obtained from different iterations. The images after 30 iterations are shown below, scaled slightly using Image > Contrast. the right hand side image is of absorption, and the left is of scatter. It is normal for the scatter image to appear to have better spatial resolution than the absorption, but it can also show more image noise. The images very clearly show an increase in both absorption and scatter in the correct place. The sharp line which appears to bisect the blob is an artefact from the regularity of the finite element mesh.

    

---

Difference imaging

This tutorial so far has described modelling and reconstructing absolute data. This means that no reference data is available, so images are reconstructed from a single set of data. This method is very sensitive to errors in the measurement, modelling and image reconstruction and so, for images of test phantoms, and particularly for clinical data, we tend to use difference imaging. In this method, we have two datasets, ideally one of a homogeneous reference object and one with the perturbations in place. We then reconstruct an image of the difference between the two states. The errors then largely cancel, generally giving a more reliable image.

To generate simulated reference data, modify the forward definition file so that the mesh is set to homogeneous by modifying the options in [INIT_PARAM] like this:

[INIT_PARAM]
RESET_MUA = HOMOG 0.01
RESET_P2 = HOMOG 1 MUS
RESET_N = HOMOG 1.4
RESET_A KEIJZER

Also change the names of the log file and datafiles, so you don't overwrite the original. Run a forward simulation, as above, and you should see much more homogeneous fields.

To reconstruct difference data, add another line to the reconstruction definition file in each of the [DTYPE] sections so that a line like

FILE = datafilename

becomes

FILE = datafilename
REFFILE = referencefilename

Note the spelling - REFFILE has two Fs!

Now run the reconstruction and view the images as before. For this example, because you've created perfect simulated data, these difference images are identical to the absolute you obtained earlier. Open up the log file in a text editor and have a look at it. There should be a line saying "| Reading reference data", confirming that you did do a successful difference reconstruction.

---

3D modelling and reconstruction

• Mesh and qm file

  Once you can use toast in 2D, try a 3D reconstruction. Make a spherical finite element mesh, using mkmesh_sphere.

  mkmesh_sphere

  The message "mkmesh_sphere: Create spherical 3D mesh using tetrahedral elements" should appear.

  • Enter a name for your mesh.
  • Give it a radius of 35 (units are mm).
  • Resolution = 3 (this controls the fineness of the mesh - bigger values are finer. Try a few values between about 2 and 5).
  • The next section wants to know values for the background properties of mua (in mm-1), mus (mm-1), refractive index (n) and region number. Enter these as a string, separated by spaces: 0.01 1 1.4 0

  You can make spherical and cylindrical meshes using mkmesh_sphere and mkmesh_cyl, circular, rectangular or brick shaped meshes using meshmod (see below) and other meshes using Netgen.

  Now make a qm file, suitable for this mesh.

  makeqm3d

  Enter the mesh name where prompted. Enter "1" to generate a plane of optodes The plane is defined by ax + by + cz + d = 0 so if, for example, you want a plane of optodes at z=10, a=0, b=0, c=1 and d=-10. Enter "0 0 1 -10". Choose 16 optodes with an angular offset of zero and a sub-surface placement of zero. Generate a similar plane at z=-10. Enter "w" and save the qm file, then "q" to quit.

• Forward Model

  Use this definition file, change the mesh and qm file names and run it as before. It will take a lot longer, probably a few minutes. Look at the data in toastim. Scroll through the different slices using the 3D: Single Cross Section window.

• Reconstruction

  Get this definition file and change the filenames again. Run toast, and monitor the log file. This might take 20-30 minutes per iteration. Also type

  top

  This tells you how much of the processor and memory you are using. If the total memory goes above about 80%, the reconstruction will start to run very slow and it may crash. If this happens, you need a smaller mesh, fewer measurements, or more memory. You might also be able to change some of the options in the .def file so the reconstruction requires less memory.

  For better reconstructions, try using pixel basis (explained here): Change

  BASIS = FWDMESH
  IMAGE_FILTER = MEDIAN EACH_UPDATE 10
  

  to

  BASIS = PIXEL 12 12 12

  This will give smoother images, but it will be slower and need more memory. To reduce the memory required, use a coarser basis like 10 10 10.

---

Problems?

• Toast didn't do what you wanted it to

  Check the spelling of any modifications you made to the .def file. Look at the .log file. The .def file contains the instructions which you sent to toast, the .log file contains what it actually did.

• Femdata gets stuck in a simulation

  Could be a mismatch between the qmfile, mesh, data or .def file. For example, if the optodes don't lie near the mesh, this can cause problems. Check the .log file.

• Toast gets stuck in a reconstruction

  If a forward simulation works, and you're using measured data, there could be an error or a mismatch between the data and qm files. If this is OK, check the mesh.

• Rotten images

  Check the mesh gives positive fields everywhere. Try using difference data if you were using absolute data before. Try using MRI instead of optical tomography.

---