Approximating circles with bezier cirves


So, anyway. I have for some time been intending to so some work with elliptical geometry. I have looked at hyperbolic geometry, which works great for laying out tree structures but doesn’t work well for graphs with a lot of connections. So, elliptical.

By using the stereographic projection, I can work with circular regions. I can do hierarchy by putting circles within circles. Mobius transformations for the geometry, and generalized circles for my drawing. Great. But how do I represent a line segment? A serious issue, as drawings of graphs tend to be made up of them.

I can represent a line segment as two circles. A circle Γ1 containing the segment , and a circle Γ2 cutting Γ1 at the endpoints of the segment γ1, γ2. A generalised circle has an inside and an outside, and can be left or right handed (in left handed circles, the “inside” is outside). A circle switches sense when transformed by a transform that turns it inside out, so identifying the segment as that part of the border of Γ1 that is “inside” Γ2 will identify a line segment and will behave ok under arbitrary mobius transforms.

Great! Now, there are a pencil of circles that will serve as Γ2 (all circles passing through the endpoints off the segment). Unfortunately, Γ1 is in that pencil. What’s a reliable way of picking a suitable Γ2 that is never Γ1? Well – the circle at right-angles to Γ1 passing through γ1 and γ2 will fit the bill perfectly.

How to find it easily?

Hmm. Map the two endpoints γ1, γ2 onto (0,0) and (1,0) and the center of Γ1 onto (∞.∞) with a transform T. The line segment γ becomes the straight line from 0 to 1. Circle 2 is then a constant (c=(.5,0), r=.5). Transform that circle by T-1 to get Γ2.

Excellent! How to draw the arc?

Well, I don’t want to dick around with cosines – java Arc2D uses degrees clockwise rather than radians (it’s build for people drawing UIs, not people doing math) – so approximating it with a bezier is the way to go. Had a look at the Aleksas Riškus’ paper, but it didn’t work for me the first time I tried it, so I’ll have a go myself.

Which brings me to the meat of this post: suitable control points for a cubic bezier curve approximating an arc.


Well, let’s make it a little simpler. Let’s assume our arc, and ∴ γ1, γ2, are on the unit circle. We cant assume anything about the position, because I don’t want to fiddle about with rotations. But it’s easy enough to adjust our result to translate and scale it into position

I want a bezier curve starting at γ1, ending at γ2, symmetrical, having control points on lines tangent to the circle. There is therefore only one parameter – k, the magic number, the length of the “control arms”, the velocity of the initial traectory of the curve.

Now … I *can* make this simpler. Construct a unti vector ŷ pointing at the midpoint of γ1, γ2 (this will need to be handled specially for angles > π). Rotate our circle (yes, I know I promised not to do this) so that this vector is vertical. Our magic number k alters the height of the resulting bezier curve. It’s clear that I don’t really need to wory about the x value. That is, when solving the underlying cubic, I don’t have to concern myself with the fact that there are two coordinates. I can get the solution for k for a cubic that’s centered around the y axis given the height of the control points above the x axis, and work out that by getting the projection of either control point of the curve I want onto that unit vector.

How to handle the >180 case? Well, given the unit vector, γ1 should be to the right of it. If it’s to the left, reverse the vector.

Hmm. Ok: at this point, wikipedia becomes useful, having a straight-up solution for the problem.


Actually … I will want to break up the arc if it’s too long, because this algorithm becomes very crap for large angles. So the actual arc drawing does not have to worry about handling reverse midpoints. How to break up the arc? Well: cos(theta) is simply the dot product of the unit vectors, right? So if it’s <.5, then use a single arc. Else, if it's less than acos(π/3) = sqrt(3)/2 , then break it up into two arcs with the midpoint method. Else, break it up into two arcs by taking 90 degrees off counterclockwise from γ1. This will work reliably and is easy to calculate, at the price of the result not being precisely symmetrical: a long arc will be drawn as a series of 22.5 degree chunks and a bit of scrap.


So! My method will be

    public static Shape arc(
     Point2D c, // center of the arc
     Point2D p1, // first point
     Point2D p2  // a ray xc->p2 cuts the circle 
                 // counterclockwise from point 1)

Let’s start coding!!

ok! The algorithm to break the curve into segments is a joy and a wonder. Works great. I use statics because this is explicitly single-threaded. IO know you worry about performance last, but I’m a bit oldschool.

static Point2D.Double p1 = new Point2D.Double();
static Point2D.Double p2 = new Point2D.Double();
static Point2D.Double p1hat = new Point2D.Double();
static Point2D.Double p1hat90 = new Point2D.Double();
static Point2D.Double p2hat = new Point2D.Double();

public static Shape arc(Point2D.Double c, // center of the arc
		Point2D.Double _p1, // first point
		Point2D.Double _p2 // a ray xc->p2 cuts the circle
			// counterclockwise from point 1
) {
	double rad = Math.hypot(_p1.x - c.x, _p1.y - c.y);
	p1.setLocation(_p1.x - c.x, _p1.y - c.y);
	p2.setLocation(_p2.x - c.x, _p2.y - c.y);
	GeneralPath pp = new GeneralPath();
	pp.moveTo(_p1.x, _p1.y);
	double acos;
	double acos90;
	for (;;) {
		hat(p1, p1hat);
		hat(p2, p2hat);

	// vector 90 degrees anticlockwise from p1hat
		p1hat90.x = -p1hat.y;
		p1hat90.y = p1hat.x;
		acos = p1hat.x * p2hat.x + p1hat.y * p2hat.y;
		acos90 = p1hat90.x * p2hat.x + p1hat90.y * p2hat.y;
		if (acos >= 0 && acos90 >= 0) {
			break;
		}
		// not in the first quadrant.
		// take off 90 degrees
		arc_twosegments(pp, c, rad, p1hat, p1hat90);
		p1.setLocation(p1hat90);
	}
	if (acos >= .707) { // cos(π/4)
		arc_addcurve(pp, c, rad, p1hat, p2hat);
	} else { // a bit of a fudge factor
		arc_twosegments(pp, c, rad, p1hat, p2hat);
	}
	return pp;
}

static Point2D.Double midp = new Point2D.Double();
static Point2D.Double midphat = new Point2D.Double();

private static void arc_twosegments(GeneralPath pp, Point2D.Double c, double rad, Point2D.Double p1hat,
		Point2D.Double p2hat) {
	midp.setLocation((p1hat.x + p2hat.x) / 2, (p1hat.y + p2hat.y) / 2);
	hat(midp, midphat);
	midp.x = midphat.x * rad;
	midp.y = midphat.y * rad;
	arc_addcurve(pp, c, rad, p1hat, midphat);
	arc_addcurve(pp, c, rad, midphat, p2hat);
}

private static void arc_addcurve(GeneralPath pp, Point2D.Double c, double rad, Point2D.Double p1hat,
		Point2D.Double p2hat) {
	if (BOING) {
		// test code
		pp.lineTo(c.x + p1hat.x * rad, c.y + p1hat.y * rad);
		pp.lineTo(c.x + p2hat.x * rad, c.y + p2hat.y * rad);
		return;
	}
}

private static void hat(Point2D.Double a, Point2D.Double b) {
	double zz = Math.hypot(a.x, a.y);
	b.x = a.x / zz;
	b.y = a.y / zz;
}

Right! So now we have to write the money bit – arc_addcurve. This is the bit that does the bezier.

TBC…

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: