So I figured I might as well blog about it.
(CSG stands for Constructive Solid Geometry by the way)
I'm doing this from memory, so hopefully i'm not missing anything, please don't hesitate to tell me if something isn't clear to you.
About Numerical Stability
Before I begin, I need to stress that with these kinds of things you
really have to worry about floating point errors.
As you might know, floating point numbers have limited precision and
cannot possible hold a number with infinite precision.
If you divide 1 by 3, you won't get 0.3333333 with 3's into infinity,
instead you'll just get 0.333333343 (single precision).
Above all, consistency is more important than accuracy with geometry
on a computer. Otherwise you always end up with coordinates that are
almost the same, but never quite the same.
You'll need to use 'epsilon values', which is a small bias
(for example 0.01), when comparing coordinates or distances.
This is because you cannot assume that two calculations that are
equal from a mathematical sense, will produce the same result.
For example:
a/b does not (necessarily) give the same result as a * (1/b)
a*(b*c) does not (necessarily) give the same result as (a*b)*c
etc.
So, for example, when you determine the distance of a point to a plane you need to define for yourself that a point will be:
- on the negative side of a plane when it's distance is < -epsilon.
- on the positive side if it's distance is > epsilon.
- on the plane itself if it's distance is between >= -epsilon and <= epsilon.
Be consistent!
Half-edges
I use a half-edge data structure to represent my geometry, you can find a basic explanation of the half-edge data structure here.
For my CSG algorithm you'll need a couple of basic methods.
First of all you'll need a method to create a pair of half-edges (with a new vertex) in between two other half-edges, which essentially splits a full edge in two and inserts the new vertex in between.
Polygon cutting
Another thing you'll need is a method to cut a polygon with a plane.
It works like this:
You cut a polygon by looping trough all your vertices, and for every vertex determine if it's on the positive side (outside), the negative side (inside) or -on- the plane you're currently cutting with. (using epsilon values like I mentioned before)
- If all vertices are on the plane, then you've got a duplicate plane and you'll need to remove that plane (quite frankly, you should've detected that before you start cutting)
- If all your vertices are on the 'outside' or are on the plane then you need to treat the entire polygon as being cut.
- If all your vertices are on the 'inside' or are on the plane then your entire polygon wasn't cut by the plane.
- In all other situations you'll need to cut your polygon into two pieces.
You do this by finding the two vertices where the polygon was intersected by the plane.
If the vertices you find just happen to lie on the cutting plane you won't need to calculate any intersection, but often you'll need to calculate the intersection point of an edge with the cutting plane, and then use the half-edge insert method I mentioned earlier.
The big advantage of half-edges here is that when you split the half-edge, you automatically do it for both polygons on both sides of the edge.
And since you won't be calculating an intersection when a point is already on the plane, you won't be splitting it twice and everything will be nice and consistent!
When you have an edge to split, you may have an edge that goes from the inside to the outside, or an edge that goes from the outside to the inside.
Now what is -really- important is that you always calculate the intersection point in the same direction.
What this means is that if you do:
myNewVertex = myCuttingPlane.Intersect(myInsideVertex, myOutsideVertex);
Then you don't mix it with:
myNewVertex = myCuttingPlane.Intersect(myOutsideVertex, myInsideVertex);
When you calculate your intersection point, stick with one direction!
Otherwise you'll end up with some vertices that are biased towards one side of the plane, and other vertices that are biased at the other side of the plane.
So make sure you're, once again, consistent here.
In the part 2 I'll explain how to build geometry from
brushes (convex polytopes), which are the building blocks of my CSG algorithm.
Update:
Created a repository for the code.