This note will describe various approaches to creating a polygonal based representation of the human cortical surface. Traditionally computer based surface modelling deals with either objects which can be represented by a number of essentially flat bounded planes (facets), for example: a table or building or objects made up of smooth curves, for example: a telephone or apple. These smooth surfaces are often formed by combinations of bounded parametric functions called patches. These patches form elegant and efficient descriptions but can be difficult to design and manipulate. Alternatively, and perhaps more commonly these smooth surfaces are approximated a large number of facets. This facet approach, while usually inefficient, is considerably easier to deal with numerically and must easier to create in an automated way from sampled data.

By comparison to the objects normally subject to computer modelling the surface of the human cortex is an incredibly complicated surface. Some idea of the geometric complexity can be seen in the following diagram. Note the folding and the sulci which typically extend inwards by one third of the radius.

Creating computer based representations of the human brain is an active area of research. A distinction should be clearly made between surface models and voxel models, the later will not be discussed here except where it is the starting point for a surface derivation as in techniques based on the tomographic methods and the marching cubes algorithm.

Note: In the computer renderings used here to represent surface models I have intentionally not used smooth shading techniques which in many cases would greatly improve the appearance of the image. The underlying geometry is illustrated better with plane shading algorithms.

Manual surface sampling

A model resulting from a manual digitisation of a plastic cast is available from ViewPoint DataLabs. The model shown below consists of over 24,000 planar facets and is the most complex they have to offer to date.

An obvious characteristic of this model is the even distribution of vertices which tend to follow the structure of the surface. In addition the facet density is function of detail, a highly desirable feature.

This approach, while totally manual, has the potential to create well designed and efficient surfaces. The primary disadvantage is the reliance on some sort of solid original. The rendering shows the relative flatness of large portions of the surface at a larger scale than the facet approximation. The sulci are also not nearly deep enough nor is the separation between the two hemispheres.

Surface reconstruction from contours

Many scanning devices, such as MRI, result in slices through an object, each each 2D slice is normally presented as a grey scale image.

Each pixel in the image corresponds to some characteristic at that point in space, for example protein density in MRI or X-ray absorption in CT scans. These images stack together to form a solid representation of the head. The contour approach attempts to trace the outline of the characteristic of interest and use these contours to form the surface.

The following example is created from contour data sourced from Brigham and Womens Hospital. There are a total of 15 contours (for one hemisphere) made up of approximately 3400 vertices.

Note that there may be multiple contours at along any slice but, as is typical of most data from tomographic techniques, the spaces between the slices is large compared to the detail in the plane of the slice..

The following surface model was created from these contours using a set of software tools (Nuages) by Bernard Geiger. These tools together form surfaces between adjacent contours based on 3D Delaunay triangulation.

This approach would seem to offer hope in forming highly convoluted surfaces, however the automatic reconstruction of the surface from disjoint contours is numerically problematic.

Voxel based reconstruction

Surface reconstruction from voxel datasets generated by MRI (Magnetic Resonance Imaging) typically use a variant of the marching cubes algorithm by Lorensen and Cline to create the surface. Before this can occur a process usually called segmentation is performed to enhance the surface boundary of interest within the volume. Segmentation will not be discussed further here except to say that it has been extensively studied with as yet no clear solutions. Unfortunately the quality of the generated surface is heavily influenced by the quality of the segmentation.

The following examples were supplied by Dennis Strelow from original data and segmentation provided by Robert Grzesczuk and David Levin. Due to the enormous number of triangular facets that result from marching cubes, the mesh is normally simplified (decimated) in such a way as to attempt to retain the geometric features.

The following image shows the surface decimated to about 22500 vertices and 45000 triangles. The most noticeable feature of such surfaces is the regular grid the polygon vertices lie on. The cubic voxels, which were the basis of this model, are obvious.

The following is the same geometry after a relaxation process has been applied which smooths out the sharp discontinuities. This relaxed surface has the same number of triangular facets but the vertices are moved so as to generate a smoother representation while still attempting to retain the features of the surface geometry.

References

Dennis W.Strelow, Clinton S. Potter, Paul C. Lauterbur
The Construction and Visualisation of Surfaces from MHRI Data
June 1996

William Erensen, Harvey E. Cline
Marching Cubes: A high Resolution 3D Surface Construction Algorithm
Computer Graphics, 21(4):163-169, 1987

William J.Schroeder, Jonathan A. Zarge, William E. Lorensen
Decimation of Triangle Meshes
Computer Graphics, 26(2):65-70, 1992

Bernard Geiger
Three-Dimensional Modelling of Human Organs and its Application to Diagnostic and Surgical Planning
PhD Thesis, 1993

J.D. Boissonnat
Shape Reconstruction from Planar Cross-Sections
Computer Vision, Graphics, and Image Processing, 44:1-29,1988

MRI dataset

The MRI dataset for this evaluation was kindly provided by Dr. F. Kruggel of the Max-Planck-Institute of Cognitive Neuroscience Inselstrasse 22-26, 04103 Leipzig, Germany. Described as:

The white matter which was used to form the basis of the implicit surface controls consisted of a clear band of intensity values, from 135 -> 215.

The control points are shown below in green for representative slices along each of the three orthogonal planes. Each control point is considered to be a sphere of 1mm radius so as to avoid overlapping.

Implicit surface modelling

In this initial experiment the implicit surface function used is usually called blobby molecules and attributed to Jim Blinn. The field D at a distance r from a single control point is given by a Gaussian

The control points were created within the white matter, an isosurface of this resulting field is an approximation to the cortex layer. In what follows the isosurface was taken to be at a level of 0.5, and A = 1.

The approximate value of B can be determined by considering two adjacent control points and the field at points between them. B should be chosen so that the isosurface from two adjacent control points is approximately cylindrical. Varying B effectively creates isosurfaces at varying distance from the white matter. High B result in isosurfaces that are close to the white matter (close to the control points), low B results in more distant isosurfaces (diffuse field and smoother isosurface).... Some example for B = 0.1 to 0.25 are illustrated below.

B = 0.1

B = 0.15

B = 0.2

B = 0.25

The second experiment was to consider each pixel in the MRI volume to be a control point, the intensity of the pixel being used to determine the value of A in the blobby molecule field equation. The mapping from intensity to A would be such that the white matter would result in high positive values, pixels clearly outside the cortex would result in negative values of A. It was hoped that this would help in the resolving of the sulci, preventing bridging between gyri.

The function used was a Gaussian with a mean of 175 and standard deviation 40, namely e-[(x - 175)/40]2 - 0.1 The zero crossings are approximately at intensities 114 and 235, in other words all intensities below 114 resulted in negative field contributions.

B = 0.05

B = 0.10

B = 0.15

B = 0.25

Credits

The implicit surfaces were rendered by polygonising the surface using a modified marching cubes algorithm written by myself....that's another story.

The resulting polygonal surfaces were rendered using Geomview on an SGI Maximum Impact.

Movies of the MRI slices in the three planes were generated using Apple QuickTime (animation) technology.

References

James F. Blinn
A Generalization of Algebraic Surface Drawing
ACM Transactions on Graphics, July 1982

Tatsumi, H., Takaoki, E., Omura, K. and Fujito, H.
A new method for 3D reconstruction from serial sections by computer graphics using meta-balls: reconstruction of "Hepatoskeletal System" formed by Ito Cells in the Cod Liver
Computers and Biomedical Research, 1990

Brian Wyvill, Geoff Wyvill
Field Functions for Implicit Surfaces
Visual Computer, Issue 5 1989

F. Kruggel, G. Lohmann
Automatic Adoption of the Stereotactical Coordinate System in MRI Datasets
Information Processing in Medical Imaging Duncan J, Gindi G (eds.), LCNS 1230, Springer 1997.

Modelling the Surface of the Human Cortex Written by Paul Bourke February 1997 A low resolution version of this model is available in the "OFF" file format as specified by the GeomView documentation: 60000.off.gz (20000 faces) or 60000.dxf.gz Higher resolution versions are available on request.

Introduction

The following informally describes a project aimed at generating a computer based surface model of a human cortex. The requirements were for a "clean" facet based model that could be used for other applications that required a surface model rather than volumetric information. In particular it was necessary to derive a model with a realistic representation of the deep sulci, a feature missing from existing surface models.

Slices

It was decided early on to use cryosection data due to its high resolution compared to MRI and other volumetric datasets available at the time.

The source for the slices used in this reconstruction is the "Human Cryotome Data" from the UCLA Laboratory of Neuro Imaging.

Examples of the slices

In total there were 317 slices.

Tracing

Every second section of the left hemisphere was traced with a third order B-Spline curve resulting in 168 slices. The software to perform this tracing was written in-house, it is a Macintosh program which reads a PICT file and optionally an existing spline file. The user may manipulate the spline tracing by moving, adding, or deleting control points. The spline curve is drawn interactively so the user has an immediate feedback as to the closeness of fit.

Needless to say the process of tracing each slice is tedious.

The overall dimensions for the model are

Surface Construction

The reconstruction of the cortical surface from the parallel B-Spline sections is still being considered. Linear interpolation between the spline slices was used here, it is straightforward but results in a surface with undesirable structure perpendicular to the slices. Another method evaluated involved stacking the spline sections together and spline curves are lofted perpendicular to the slice sections through the equivalent points on each section.

Improved reconstruction approaches include spline surface descriptions, implicit surfaces, and spherical harmonics.

The aim is to derive a surface description which is "well behaved" and can ideally be constructed at variable degrees of detail.

At this stage it is possible to estimate any position on the 2D spline surface and thus form triangular facet versions of the cortex surface of any desired resolution (up to some maximum as determined by the original sampling and slice thickness). The model seen at the start of this page was created by specifying a maximal edge length for the triangular facets of 1mm.

Editing/smoothing

Due to variations in the slice splines such as image registrations errors, the final reconstructed surface has undesirable vertical features (see left image below). A standard relaxation filter was passed over the surface to smooth it out, the effect is shown below.

The algorithm used is applied repeatedly depending on the degree of smoothing required, two iterations were used to create the final cortical model seen at the top of this page.

Images of fibres in a model of the cortex Written and created by Paul Bourke May 1997

Illustration of a computer generated connectivity between discrete areas on a human cortex. Two modelling exercises have been combined in this illustration; the formation of a realistic surface representation of the human cortex and the result of a simulation to derive pathways by which fibre bundles connect surface patches.

Generating a realistic surface representation of the human cortex involves a reconstruction from geometric data acquired from a real cortex. A common approach is to use MRI scans which come in the form of a number of 2 dimensional slices (images) through the 3 dimensional volume. While there may be sufficient resolution within a slice, the resolution perpendicular to the slices is much poorer resulting in inadequate spatial sampling in that dimension. It is thus not possible, or at least very difficult, to get a realistic reconstruction of deep sulci using MRI data. The approach taken here was to reconstruct the surface geometry from cryosection slices for which there is significantly better resolution both within a slice and between slices. The cryosection data was sourced from the UCLA Laboratory of Neuro Imaging, each slice is separated by 0.6mm and was digitally photographed at 1000x1000 pixel resolution. The upper cortical surface was traced using a mixture of automatic and manual procedures to form ribs over which 2 dimensional splines were lofted. This parametric form of the cortical surface can be sampled at any user selected resolution to give simpler surface representations approximated by triangular facets. The cortical simulation software being developed within the Brain Dynamics Laboratory requires information on how small areas (patches) on the cortex are coupled together, that is, the strength of coupling and delay between any two patches. Both of these parameters are functions of the distance of the pathway taken by the fibre bundle. The pathways were derived by modelling the cortex as a positively charged surface and a fibre as a collection of positively charged springs connected end-to-end along the pathway. The fibre path is found by iteratively relaxing the spring model to find the route of minimum potential energy. This model gives the desired properties we chose for the fibre pathways, the fibres are roughly perpendicular to the surface of the cortex, they attempt to find the shortest route between patches, and they don't intersect the surface or each other. For this type of whole brain macroscopic simulation the surface of the cortex was split into about 600 patches. Fibres from any one patch can potentially connect to any other patch resulting in a huge number of fibres. In practice however, the very long or convoluted pathways have a very small coupling value and can be safely pruned. Additional long range fibres representing the superior longitudinal fasciculus, inferior longitudinal fasciculus, and unicinate fasciculus were also added to the model. These additional longer range fibres have their pathways determined by similar techniques but their coupling strength is based upon a different function of distance.

The illustration is a view from above showing the surface of the cortex as transparent facets with black edges. The fibres are drawn as blue spline curves connecting the centers of the surface patches.