Symmetry Set, Midpoint Locus, and Medial Axis

Drag to move a control point, click to add a new control point, Shift-click to delete a control point.

The code for drawing the splines comes from Tim Lambert's web page at http://www.cse.unsw.edu.au/~lambert/splines/ and forms the core around which this applet was built.

What this applet displays and how it displays it

The medial axis of a closed curve is the set of all points inside the curve that do not have a unique closest point on the curve. Since a point, x, in the medial axis lies on the normal line to the curve at its closest point, y, a circle centered at x of radius |x - y| is tangent to the curve at y. So it is natural to extend the medial axis to the larger set consisting of all points in the plane lying at the center of a circle that is tangent to the curve at two or more points; this set is called the symmetry set of the curve and the points of tangency we will call bi-tangent points (to the curve). If we take the midpoint of all bi-tangent points in place of the centers of the bi-tangent circles, we have what is called the midpoint locus of the curve.

P. J. Giblin and S. A. Brassett, in Local Symmetry of Plane Curves December 1985 American Mathematical Monthly p. 689-707, characterize the symmetry set and midpoint locus of a smooth curve in the plane. Once two bi-tangent points are found, the equations on p. 693-694 of Giblin and Brassett allow us to move one bi-tangent point a small distance along the curve and determine where the other point must move to remain bi-tangent (to a nearby circle). They are in effect linearizing a differential equation involving the tangent vectors and curvatures at the bi-tangent points. So we can "walk along" the curve producing pairs of bi-tangent points, adding their centers to the symmetry set as we go, until an obstruction is reached. Obstructions can take several forms, but the most desirable is that we have reached the end of one piece of the symmetry set. We will call such a piece a component of the symmetry set (it is a one-manifold, but is not topologically a component, since our components can intersect each other). This is the basis of our algorithm: find pairs of bi-tangent circles and then walk along them, taking some care to avoid walking the same part of the symmetry set twice.

(The symmetry set is more complicated than we have described above; it can, for instance, contain isolated points and the components can have cusps. We have no hope of detecting such isolated points for a general curve; however, we do--usually--detect the cusps. See Giblin and Brassett's paper for a full description of the symmetry set.)

Finding the pairs of bi-tangent points is a global process--we must look everywhere along the curve if there is no detailed a priori knowledge of the geometry of the curve. Walking the components is a local process. If a point lies in the medial axis, then it is the center of a circle bi-tangent to the curve that remains wholly within the curve; determining whether or not this is true is another global process. Global processes tend to be more expensive to compute than local processes. This is why, at least for the approach we have taken, it is harder to calculate the medial axis than it is to calculate the symmetry set. The best compromise we have reached so far is to display only those components of the symmetry set containing at least one point lying in the medial axis, since this way we need only check for membership in the medial axis before beginning to walk the component, not for each point along the component. This is the portion of the symmetry set being displayed in the sample execution above.

The primary objective in creating this applet is to be able to calculate the medial axis quickly enough that one can drag a control point on the curve and watch the medial axis change almost immediately. Some loss of accuracy is acceptable, but only up to the resolution of the display. The more time spent trying to find bi-tangent points to seed the process, the longer the execution, but the lower the likelihood of missing part of the symmetry set. This is the main tradeoff in fine-tuning the algorithm. A complicating factor is that we are currently solving the differential equation used to walk along a component by very simple-mindedly linearizing the equation and using a fixed step size. This can cause the solution to diverge. We detect this divergence, but stop the walking at that point, and so partial components are often displayed. To compensate, we seed the walking with more pairs of bi-tangent points, but this slows down the execution substantially. But then this applet is a work in progress, and maybe one day we will spend enough time on it to overcome these difficulties.

Running the applet

To understand the options, it is best view the source of this page. The pertinent portion is the following:

It is a lot easier to play with the applet if the width and height, specified on the first line above, are increased, but I didn't want to prevent the viewer from seeing the applet on the first page. Don't change the second or third lines or the applet won't work. The controlPoints on the fourth line are a series of x, y values for the closed curve that is initially displayed. The values are relative to the upper left corner, unfortunately, with positive values going down and to the right.

The fifth line contains arguments controlling what type of set is calculated. By default, the symmetry set is displayed (or enter the value "symmetry"). The value "midpoint_locus" displays the midpoint locus; the X in front of "midpoint_locus" produces an unrecognized value, making the applet choose the default of the displaying the symmetry set. No error message is issued. This is a convenience, making it easy to switch back and forth between the symmetry set and the midpoint locus. The symmetry set or midpoint locus is displayed in red.

The second option in line five--inside--causes only those points that are (topologically) connected to the medial axis to be displayed. When displaying the "exact" solution, however (see next paragraph), only the media axis is displayed.

The third option in line five--Xexact_set_intersecting_normal_lines--is also disabled by an X. When enabled, the applet will calculate the symmetry set by a fairly brute force method of finding the intersection of normal lines all around the curve. The resulting set will consists only of isolated points, since calculating the symmetry set this way is very inefficient. This version of the symmetry set is displayed in blue.

The value of walkmesh on the sixth line can be small, medium, large, or huge for a walk mesh size of 50, 150, 500, or 2500 (as described below). Large seems to do the best job for a cubic with only a few control points, but the whole walk mesh should be done dynamically, removing the need for this parameter.

The seventh line contains debug options. The first option, which can be either center_func_giblin (the default) or center_func_radii determines which of two functions to use to find bi-tangent points. They correspond to the functions g and F described below. The second two options have been cleverly disabled by adding an "X" in front of them, making it easier to enable when desired. The option bitangent_circles will cause every tenth bi-tangent circle to be displayed in yellow as the algorithm "walks" along the curve. The option bitangent_normals displays the fixed-length (but non-unit) normals at the bitangent points at which walking begins.

In brief, the applet works as follows. It views the closed curve as a period-1 map from the real numbers into the plane--that is, a map from the unit interval into the plane except that arguments differing by an integer map to the same point. The cubic spline class knows how to calculate points on the curve, tangent vectors, and normal vectors, but code that calculates the symmetry set has no knowledge of these details (it is in a separate class). The applet determines for each of a fairly large set of points, equally spaced by parameter value along the curve, all the points sharing a tangent circle with it. The applet then walks each pair of bi-tangent points along the curve, first in one direction, then in the next. A bitmap to the resolution of the image is used to indicate which points have been walked; walking of a pair of bi-tangent points will not occur if that point or a neighboring point (in the bitmap) has already been walked. (If this is not done, not only are unnecessary computations performed, but the small errors in computation during a walk blur the image of one or more components.)

The bi-tangent points are located by determining the zeros of one of two functions--either the function g described on page 693 of Giblin and Brassett's paper, or the function F on page 692 of that same paper--by walking in fixed steps (the number of steps being specified by the walk_mesh parameter), then using bisection to pinpoint the zero. A seeming advantage of using g is that it is numerically stable for bi-tangent points with nearly parallel tangent vectors, while F is not. But F will still change sign near such bi-tangent points and that is all that is needed. On the other hand, F does not require that the curves near each point be oriented in any particular way, and so is easier to compute. In fact, we do not account for orientation in our computation of g, and so with our applet it is better to use F.

Arguments for expanding the above applet by a factor of three linearly are: