I think of this as a very applied paper, because this mathematics is not new-- its all been done in computer vision, years ago. I'm just trying to present it in a way that's easy to understand, so that its very clear how to map it to actual working code.
I think the contribution of this paper, then, is that description of a practical method for curvature measurement, combined with these 3 new ways of leveraging that information in volume rendering. ----------------------------------------- Previous work png/Slide03.png Among the many people who have measured curvature in volumes, our work is partly inspired by Monga et al., who analyzed curvature with 3rd derivatives to find feature lines of extremal curvature,
Others have used principal curvature directions to indicate surface shape with textures and kinetic visualization, but our paper is most directly inspired by Hladuvka, who first showed that a 2-D space of curvature magnitudes provides an intuitive and effective way to highlight volume features.
As far as non-photorealistic volume rendering, Rheingans and Ebert mapped many NPR methods to volume rendering in their Vis 2000 paper. Oher people have worked on pen-and-ink [Trevett and Chen], or contour-based [Csébfalvi], or stippling-based [Lu et al.] methods.
More recently, NPR volume rendering has become interactive with graphics hardware. This paper [Nagy et al.], for example, is interesting for its interactive use of curvature-directed strokes and hatches. ----------------------------------------- What is curvature png/Slide04.png Curvature describes how small movements along a surface result in changes to the surface normal, roughly, the shape of the surface.
This two-dimensional space (on the right) captures all the possibilities: planar, and cap, ridge, saddle, valley, and cup. This is parameterized by kappa1 and kappa2, the principal curvatures magnitudes, which are two of the degrees of freedom in curvature.
The other degrees of freedom in curvature are the orientation of these lines, the principal curvature directions, which are the directions along which the cross-sections have maximal and minimal curvature.
Curvature is both these magnitudes and these directions. We feel it is simplest to have a general framework for measuring all curvature information, from which we can then extract those components needed for a given problem, rather than rely on the various special-case formula which exist for particular scalar curvature metrics, as is often done with level sets, for example. ----------------------------------------- Curvature measures png/Slide05.png I made these for fun, to try to build my intuition for the relationship between the two principal curvatures, and the various curvature metrics that are defined in terms of them.
These are all volume renderings. ----------------------------------------- Application area png/Slide06.png We are concerned with medical datasets, like CT or MRI, comprised of scalar values on a regular grid.
We make the simplifying assumption that the surfaces we wish to visualize, are isosurfaces of the field. This implies a few things:
But we care about changes in the gradient vector, so when we differentiate this, we see that the change in the gradient vector around x0 is found by multiplication by the Hessian matrix.
So the Hessian is the basis of our curvature computations.
To see exactly how we get curvature from the Hessian, ... ----------------------------------------- Projections, Hessian, Invariant: n png/Slide08.png I'm going to give a pictorial version of in Section 3 of our paper, which is a more detailed and rigorous proof.
If n is the surface normal ... ----------------------------------------- Projections, Hessian, Invariant: nnT png/Slide09.png ... n times n transpose is a linear operator which projects onto this line, spanned by n. ----------------------------------------- Projections, Hessian, Invariant: I - nnT png/Slide10.png The identity matrix minus that, which we'll call P, projects onto the complement of that line, namely, the tangent plane. ----------------------------------------- Projections, Hessian, Invariant: H png/Slide11.png Remember from the last slide that the Hessian tells how the gradient vector changes as a function of displacements away from the point at which its measured. Those displacements, and the resulting gradient changes, are both vectors in three dimensions. ----------------------------------------- Projections, Hessian, Invariant: HP png/Slide12.png But if we preceed the Hessian by multiplication by P, we are only looking at changes in the gradient vector due to displacements *within* the tangent plane. ----------------------------------------- Projections, Hessian, Invariant: PHP png/Slide13.png And if, in addition, we also follow by multiplication by P, we are only looking the component within the tangent plane, of changes in the gradient, due to displacements within the same tangent plane. ----------------------------------------- Projections, Hessian, Invariant: plane is the surface png/Slide14.png But in a differential sense, the tangent plane is the surface, so this quantity basically is curvature. ----------------------------------------- Projections, Hessian, Invariant: normalization png/Slide15.png Except for that normalization which converts the gradient vector to a surface normal.
Once we build in that normalization, we have G, the geometry tensor. G tells you everything you might want to know about implicit surface curvature, so it is G that we measure for our curvature-based transfer functions.
[Note: This was sneaky. The normalization which converts the gradient vector to the surface normal includes a minus sign, but that is NOT implied by my diagrams, because they would be too confusing otherwise. I figured that those in the audience who caught this sloppiness would sympathize, and those who didn't would be none the wiser....] ----------------------------------------- Projections, Hessian, Invariant: finding kappa1, kappa2 png/Slide16.png Lots of the time, we want to know the principal curvature magnitudes kappa1 and kappa2 When expressed in the coordinate frame of the principal curvature directions, G is just all zeros except for kappa1 and kappa2 along the diagonal.
But we don't know the principal curvature directions, so we can't just pluck out those two elements and be done with it. So we use two matrix invariants, the trace and the frobenious norm, which are the same regardless of the basis in which G is expressed. The trace is the sum of diagional elements, the frobenius norm is the l2 norm of the matrix considered as a big 9-vector.
So these two invariants will show up in our ... ----------------------------------------- Recipe for Curvature Measurement png/Slide17.png ... recipe for curvature measurement. Here it is.
Note that this started with measuring partial derivatives. And we do that based on ... ----------------------------------------- Derivative Measurement png/Slide18.png ... convolution, which is how we reconstruct a continuous scalar field from our discretely sampled data.
We don't have time for details, but I'd like to note that the analytical derivative of a continously reconstructed signal is the result of convolving the same data, with the derivative of the filter.
I should also mention that its helpful for the filter to have some smoothing built into it, as these do, in order to help deal with real-world noisy data. ----------------------------------------- Derivative Measurement in 3-D png/Slide19.png Our three-dimensional filters are just the usual seperable product of three 1-D filters, so the 3-D partial derivative filters are some combination of reconstruction and derivative filters along X, Y, and Z.
This is sort of filter you'll need for the gradient, these are thes sorts of things you'll need for the Hessian.
Remember that we're using different filters, convolved with the same data, to get different derivatives, so there is no pre-computing or storage overhead. Everything is computed on the fly.
So the last issue is, how do you pick these filters? ----------------------------------------- Spatial Filter Design png/Slide20.png We use the spatial filter design framework of Moller et al. This is actually the first time the framework has been used to create second derivative filters.
Instead of controlling frequency domain properties like aliasing, leaking, and so forth, this works directly in the same spatial domain in which the filter is applied, analysing the Taylor series expansion of the convolution sum.
The main quality metric you control, called accuracy, is the degree of the polynomial that can be reconstructed exactly. You can also control the derivative, continuity and support characteristics
What you get out are piece-wise polynomial filters with symmetric support, which we like for their efficiency of evaluation.
Hoping to provide practical advice on choosing a set of filters for curvature measurement, we did a little experiment to find the optimal combination of filters, with 4x4x4 sample support, which is really the bare minimum for any kind of second derivative measurement. We tested filter combinations by rendering the mean curvature on the surface of the Marschner-Lobb test dataset. ----------------------------------------- Spatial Filter Design Experiment png/Slide21.png For comparison, here's the analyical result, and here's what you can do if you can afford an 8x8x8 sample support, with degree six polynomials. These are actually the best filters available in the current implementation of the framework.
This image, though, is really the result of our experiment, showing the best possible filter combination with a 4x4x4 sample support. These have polynomial degrees 7, 4, and 1.
But we were curious about how more familiar cubic filters and their derivatives would compare. These are the results of using Catmull-Rom for reconstruction, and the first and second derivative of the cubic B-spline. Its pretty close!
So a secondary result of this experiment is that these familiar cubic polynomials, and their derivatives, do a very respectable job of measuring curvature.
If want some smoothing built into your reconstruction, replace Catmull-Rom ... ----------------------------------------- Spatial Filter Design Experiment (2) png/Slide22.png ... with the uniform cubic b-spline.
That's the combination of filters used throughout the paper. ----------------------------------------- Half-way Point png/Slide23.png Okay, half-way done.
We know how to measure curvature.
Now we're going to put it into volume rendering. ----------------------------------------- Direct Volume Rendering png/Slide24.png Transfer functions play the important role of mapping the field properties that we sample along the ray to things that we can composite, like colors and opacities.
Notice that among these properties is the geometry tensor G.
So in each of these three applications, different aspects of implicit surface curvature will be mapped through the transfer function to enrich the volume rendering in different ways. ----------------------------------------- Level Sets png/Slide25.png First up are level sets, which, in general, iteratively deform an implicit surface according to a partial different equation of the surface properties, usually involving curvature.
We looked at a particular level set method, which generalizes Perona+Malik's anisotropic image diffusion by working in the space of surface normals. This work was presented at last year's Vis conference by Tolga Tasdizen. The curvature property that matters here is total curvature, the Frobenious norm of the geometry tensor.
We have adopted direct volume rendering as a level set visualization tool:
These pictures confirm the properties intended by the algorithm, that contiguous, high curvature edges are preserved, while noisy curvatures spots are removed.
The renderings also show this interesting network of edges which changes in topological connectivity throughout the processing. That could give us ideas about new directions for mathematical analysis of this method. We view this as an example of volume rendering providing useful feedback in the design of image processing tools. ----------------------------------------- Isosurfaced Uncertainty png/Slide27.png Here is the second of the three applications of curvature-based transfer functions. I should say that this one is probably the most experimental in nature.
Let's consider these two kinds of boundaries.
Here are two flat regions, with parallel isocontours. If we place water here, it will flow down in a straight path, viewed from above. But here, the high part is not level, the isocontours aren't really parallel. and water flows in a curved path, when viewed from above.
I think we can also say that the individual iso-contours here are a poor model of the real boundary between the two regions, which is something like this.
What we have here is called "flow-line curvature": a property unique to implicit surfaces, which measures how much the surface normal tilts as you move off of, or further into the surface. We measure flow-line curvature with the same framework described before, see the paper for details.
We propose that flow-line curvature can be a qualitative tool for visualizing the uncertainty of surface models generated by isosurface extraction, in the following sense.
Now these holes are the problem: they are not a material boundary. They're from the bone being very thin here, compared to here.
And notice that those holes are lined by very high flowline curvature, and in fact their formation is actually foreshadowed by flowline curvature here. Note that there are high surface curvature spots with low flowline curvature: like this edge here, or the surface of this round bone.
These pictures end up being a visualization of the fact that its quite hard to find a CT isosurface that reliably shows the bone surface. ----------------------------------------- Non-photorealistic volume rendering png/Slide29.png The third and final use of curvature-based transfer functions, non-photorealistic volume rendering, comes in two parts: controlling the thickness of contours, and doing ridge and valley emphasis. ----------------------------------------- Contours png/Slide30.png Contours (or silhouettes) are a basic part of non-photorealistic volume rendering. They're easy to make, and they help emphasize the shapes of objects.
They're usually make with a view-dependent environment map, where you find the dot product between the view vector and the surface normal, and if it falls within some fixed range, then you color the surface differently.
Here's the result of doing that. The problem here is that the thickness of the contour varies, depending on the shape of the object. ----------------------------------------- Better contours: Idea png/Slide31.png But now we can do that with curvature based transfer functions.
We vary the world-space contour thickness, that range of surface normals that we color differently, to create the appearance, in image-space, of a constant thickness contour.
This is possible with two pieces of information, the curvature along the view direction, which is easy to extract from the geometry tensor, and the desired image-space contour thickness, T, which is the one parameter the user controls.
And you do a little bit of trigonometry, and you find that you still express the contour in terms of v-dot-n, but now this side varies with the view direction curvature. We implement these contous with a 2-D lookup table, indexed by v dot n, and the curvature. ----------------------------------------- Better contours: Results png/Slide32.png And here are results, for two different settings of T.
T really is the only parameter to vary, so these are easy to create. ----------------------------------------- Ridge and Valley emphasis: Idea png/Slide33.png Okay, last section.
Creases, or sharp edges on a surface, are either ridges or valleys. Emphasizing them is a basic part of standard technical illustration.
We do this in volume rendering by working in the two-dimensional space of the principal curvature magnitudes, shown again here.
Following convention, we'll select these regions, the ridges, and lighten them, and these regions, the valleys, will be darkened. ----------------------------------------- Ridge and Valley emphasis: Results png/Slide34.png If we use the rainbow colormap, we can first see that we're correctly detecting valleys, in yellow, along here, and ridges, in cyan, here, and all along here.
And here's the result of using the ridge and valley emaphasis, as implemented by this transfer function. ----------------------------------------- Volume NPR: results (skin) png/Slide35.png Now I'll show the results of the two volume illustration methods combined.
Here's a rendering with no methods, except for the cool-to-warm shading. ----------------------------------------- Volume NPR: results (skin, old contours) png/Slide36.png Now with standard contours ... ----------------------------------------- Volume NPR: results (skin, new contours) png/Slide37.png ... and now with new contours. Note that we're not exactly getting this contour right (on the top and back of the head), because this is a very low curvature region. Accurately measuring these low curvatures is really hard, in the presence of any noise. But I think the effect still has an appropriate illustrative effect. ----------------------------------------- Volume NPR: results (skin, ridge and valley) png/Slide38.png And here's the ridge and valley emphasis. The effect here is to really highlight all the folds of the skin, in the eyes, on the face, on the neck. ----------------------------------------- Volume NPR: results (skin, both effects combined) png/Slide39.png And now combined with the contours. ----------------------------------------- Volume NPR: results (bone) png/Slide40.png Okay, same sequence here, starting with cool-to-warm shading. The goal here will be to properly highlight the structure of the vertebrae here. ----------------------------------------- Volume NPR: results (bone, old contours) png/Slide41.png Old contours. ----------------------------------------- Volume NPR: results (bone, new contours) png/Slide42.png New contours.
The new contours are thinner in the low curvature regions, but they're actually thicker around the nose, and the spine, where the details are smaller. ----------------------------------------- Volume NPR: results (bone, ridge and valleys) png/Slide43.png Now here's the ridge and valley emphasis. This looks like kind of a mess right now, but notice how the teeth highlighted as alternating ridges and valleys, and we've picked up the sutures on the skull surface, and on the nose.
As is, the overall shape of features is a little hard to see, but this is addressed by ... ----------------------------------------- Volume NPR: results (bone, both effects combined) png/Slide44.png ... putting back in the contours. ----------------------------------------- Volume NPR: Animation png/Slide45.png Here's an animation for which the frames only finished rendering five minutes ago: vmhead0.mov. Now I was going to stand up here and say that because there is nothing stochastic in the rendering process, I'm guaranteed frame-to-frame coherence. The problem here is that these images are really undersampled, so we're getting a quite a bit of aliasing.
[Note: I've since re-rendered the animation with larger image resolution, and a softer opacity function (which increases the image-space size of features, and decreases the minimum image size to avoid undersampling): vmhead.mov] ----------------------------------------- Discussion png/Slide46.png I think its useful to know that 4x4x4 filters based on cubic polynomials are actually all you need to get decent curvature measurements,
and I think the bigggest challenge now is trying to measure the curvature of large features (like the skull in the images before), in the presence of noise, either from measurement, or from quantization.
I think the NPR work can be expanded, by playing with different styles of contours, different kinds of effects such as halos.
Finally, my implementation is really slow, but I'd love to learn about the extent to which its possible on current graphics hardware. ----------------------------------------- Final note png/Slide47.png Finally, I'd like to draw your attention to the last sentence of the paper: "Every rendered image in this paper can be regenerated exactly with open-source software and public datasets...", and here are step by step instructions for recreating all the datasets, all the rendering parameters, and all the images in the paper.
This was an experiment, a time-consuming one, based on the idea that perhaps, visualizations should be reproducable, especially if they're supposed to be scientific visualizations.
Hopefully if the logic of this isn't too far-fetched, then others will do the same.
Anyway, thanks to the NIH for funding, and for the reviewers and other members of my department for their help in preparing the paper and the talk.