Paul's profileWelcome to my pagePhotosBlogListsMore ![]() | Help |
Welcome to my pageFeel free to browse around |
|||||||||||
|
January 24 2-d Triangulation algorithm in C#This post builds on the code presented in the previous post for testing whether a point is within a triangle. The following code triangulates a 2-d polygon. The input to the code is a list of points representing the perimeter of the polygon. Its output is a list of triples of integers. The integers are indexes into the list of points which, if taken in order, represent a collection of triangles that completely fills the polygon. The algorithm is a simple “ear-clipping” technique, which isn’t the fastest technique, but works well enough for small polygons. Wikipedia has a concise entry on polygon triangulation here, but there’s good information elsewhere on the web. As you can see from the namespace, this code comes from the 3-d Maze program. The 3-d geometry for the walls of the maze is first built in two dimensions by triangulating the 2-d floor plan of connected wall segments and then extruding that floor plan into a 3-d solid. Due to size constraints, I can’t post the complete source code in one blog posting, but you can get the full code from my SkyDrive here.
/// <summary> /// Implementation of a greedy ear-clipping algorithm for /// triangulating simple (no holes) polygons /// </summary> public class Triangulator { public event Action<double> ProgressChanged;
private IList<Point> _vertices; private LinkedList<Vertex> _vtxIndex; private LinkedList<LinkedListNode<Vertex>> _ears; private LinkedList<LinkedListNode<Vertex>> _reflex;
/// <summary> /// Constructor. Pass in the list of points to triangulate. /// </summary> /// <param name="vertices"></param> public Triangulator(IList<Point> vertices) { this._vertices = vertices; }
/// <summary> /// Triangulate the polygon /// </summary> /// <returns></returns> public Int32Collection Triangulate() { if (_vertices.Count < 3) { throw new ArgumentException("Cannot triangulate with less than 3 vertices"); }
// Set up phase InitializeVertexList(); InitializeReflexList(); InitializeEarList();
// Run the algorithm return InnerTriangulate(); }
/// <summary> /// Determine if the vertex is an 'ear'. That is a convex vertex /// where no other points lie in the triangle. /// </summary> /// <param name="v"></param> /// <returns></returns> private bool IsEar(LinkedListNode<Vertex> v) { // Only convex vertices can be ears if (!v.Value.IsConvex) { return false; }
int vtxIndexNext = v.CircularNext().Value.idx; int vtxIndexPrevious = v.CircularPrevious().Value.idx;
Triangle t = new Triangle(); t.p1 = _vertices[vtxIndexNext]; t.p2 = _vertices[v.Value.idx]; t.p3 = _vertices[vtxIndexPrevious];
TriangleTester tt = new TriangleTester(t);
// Test to see if any of the reflex vertices is in the triangle foreach (var r in _reflex) { int vtxIndex = r.Value.idx; if (vtxIndex != vtxIndexPrevious && vtxIndex != vtxIndexNext) { // Test reflex vertices only Point test = _vertices[vtxIndex]; if( tt.IsPointInTriangle(test) != TriangleTester.TestResult.Outside ) { // Point is inside or on an edge return false; // Not an ear } } }
return true; }
/// <summary> /// Initialize the reflex list /// </summary> private void InitializeReflexList() { _reflex = new LinkedList<LinkedListNode<Vertex>>();
foreach( var v in _vtxIndex.GetNodeEnumerator()) { switch(Convexity(v)) { case Triangle.VertexType.Convex: // Do nothing break;
case Triangle.VertexType.Reflex: v.Value.reflexRef = _reflex.AddLast(v); break;
case Triangle.VertexType.Collinear: // Remove vertex // v.RemoveSelf(); // Do nothing. Let the algorithm remove Collinear ears break; } } }
/// <summary> /// If the vertex was reflex, see if it has become convex /// If it is now convex, test for an ear. /// If the vertex has become flat (collinear), then it /// should be removed. Return true in this case. /// If the vertex was convex, it will stay convex, but /// we should test to see if it has become an ear. /// </summary> /// <param name="v">The vertex to adjust</param> /// <returns>Nothing</returns> private void FixUpVertex(LinkedListNode<Vertex> v) { if (!v.Value.IsConvex) { var type = Convexity(v);
if (type == Triangle.VertexType.Reflex) { // Reflex stayed reflex return; }
// No longer reflex v.Value.reflexRef.RemoveSelf(); v.Value.reflexRef = null;
if( type == Triangle.VertexType.Collinear ) { // Need to remove this vertex entirely // Don't try to remove it here because doing so might // invalidate 'ears' on either side. Instead, add // it as an ear to the head of the list and it will // be removed in due course. v.Value.earRef = _ears.AddFirst(v); return; }
// fall through and test for "earness" }
// Test for 'earness' if (IsEar(v)) { if (!v.Value.IsEar) { // Newly created ear. v.Value.earRef = _ears.AddFirst(v); } } else { if (v.Value.IsEar) { // No longer an ear, so remove it from the ears list. v.Value.earRef.RemoveSelf(); v.Value.earRef = null; } }
return; }
/// <summary> /// Run the ear-clipping algorithm /// </summary> /// <returns></returns> private Int32Collection InnerTriangulate() { // For a polygon with N vertices, there will be N-2 triangles (3 ints per triangle) Int32Collection triangleIndices = new Int32Collection(3 * (_vertices.Count - 2));
// Keep clipping ears until only 3 vertices (single triangle) remain while (_vtxIndex.Count > 3) { if (ProgressChanged != null && _vtxIndex.Count % 64 == 0) { FireProgressChangedEvent(); }
// Take the first ear in the list var vEar = _ears.First.Value;
// Add the ear to our output triangle list var vNext = vEar.CircularNext(); var vPrevious = vEar.CircularPrevious();
triangleIndices.Add(vPrevious.Value.idx); triangleIndices.Add(vEar.Value.idx); triangleIndices.Add(vNext.Value.idx);
// Remove the ear from the ear list and the vertex from the vertex list vEar.Value.earRef.RemoveSelf(); vEar.RemoveSelf();
if (_vtxIndex.Count == 3) { // Normal early exit break; }
// Exmine the two remaining vertices. FixUpVertex(vPrevious); FixUpVertex(vNext); }
// Add the remaining triangle triangleIndices.Add(_vtxIndex.Last.Value.idx); triangleIndices.Add(_vtxIndex.First.Value.idx); triangleIndices.Add(_vtxIndex.First.Next.Value.idx);
return triangleIndices; }
private void FireProgressChangedEvent() { double dProgress = (double)(_vertices.Count - _vtxIndex.Count) / _vertices.Count; ProgressChanged(dProgress); } } January 18 Point in Triangle Test in C#Here’s some C# code for testing whether a 2-d point lies within, outside or “on an edge” of a triangle. The last case, “on the edge”, is frequently overlooked in simple implementations. As you’ll see in the next post, this implementation is a foundational part of other geometric algorithms. The tester makes use of classes like Point, Vector and helper functions like “CrossProduct” and “IsCloseToZero”. These can be found in the full code on SkyDrive here. /// <summary> /// A triangle helper class /// </summary> internal struct Triangle { public Point p1; public Point p2; public Point p3;
/// <summary> /// Describes the type of the vertex p2. /// </summary> public enum VertexType { Collinear, Convex, Reflex }
/// <summary> /// If you take the path from p1 to p2 to p3, is the turn at /// p2 a left turn or a right turn? /// If p1, p2, p3 are adjacent points on the polygon, then /// this determines whether the angle at p2 is convex, reflex /// or collienear. /// </summary> /// <returns></returns> public VertexType Convexity { get { // Compute the cross product of two vectors in the triangle // Compute the cross product of two vectors in the triangle Vector v12 = Point.Subtract(p2, p1); Vector v23 = Point.Subtract(p3, p2); double dCross = Vector.CrossProduct(v12, v23);
if (dCross.IsCloseToZero()) { return VertexType.Collinear; }
// Cross product is -ve if the interior angle is < 180 degrees return (dCross < 0) ? VertexType.Convex : VertexType.Reflex; } }
/// <summary> /// The bounding box for the triangle. Useful for hit testing. /// </summary> public Rect Bounds { get { Rect rc = new Rect(p1, p2); rc.Union(p3); return rc; } } }
/// <summary> /// Used for point in triangle tests. Pre-computes some vectors /// and a bounding box /// </summary> internal class TriangleTester { private Triangle t; private Vector v12; private Vector v23; private Vector v31; private Rect rcBounds;
public TriangleTester(Triangle t) { this.t = t; this.v12 = Point.Subtract(t.p2, t.p1); this.v23 = Point.Subtract(t.p3, t.p2); this.v31 = Point.Subtract(t.p1, t.p3); this.rcBounds = t.Bounds; }
public enum TestResult { Inside, Outside, OnEdge }
private bool InBoundingRect(Point test) { // Surprisingly, this is (slightly) faster than calling rcBounds.Contains if (test.X < rcBounds.Left || test.X > rcBounds.Right || test.Y < rcBounds.Top || test.Y > rcBounds.Bottom) { return false; }
return true; //return rcBounds.Contains(test); }
public TestResult IsPointInTriangle(Point test) { // First check if the point is even in the bounding rectangle if (!InBoundingRect(test)) { return TestResult.Outside; }
Vector p1Test = Point.Subtract(test, t.p1); Vector p2Test = Point.Subtract(test, t.p2); Vector p3Test = Point.Subtract(test, t.p3);
// The cross products will include the sin of the angle // zero indicates 0 or 180 degrees // +ve indicates "to the right of" // -ve indicates "to the left of" double cross12 = Vector.CrossProduct(p1Test, v12); double cross23 = Vector.CrossProduct(p2Test, v23); double cross31 = Vector.CrossProduct(p3Test, v31);
// Check for collinear points if (cross12.IsCloseToZero()) { // Point is collinear with the 1.2 edge return cross23 * cross31 > 0 ? TestResult.OnEdge : TestResult.Outside; }
if (cross23.IsCloseToZero()) { // Point is collinear with the 2.3 edge return cross12 * cross31 > 0 ? TestResult.OnEdge : TestResult.Outside; }
if (cross31.IsCloseToZero()) { // Point is collinear with the 3.1 edge return cross12 * cross23 > 0 ? TestResult.OnEdge : TestResult.Outside; }
// To be in the triangle, the test points must all // be +ve or all -ve. if (cross12 * cross23 > 0 && cross23 * cross31 > 0) { return TestResult.Inside; }
return TestResult.Outside; } }
August 21 Gaussian Random Number Generator in C#The other day, on one of the discussion boards I read, someone asked for a random number generator which returned numbers distributed according to the Gaussian distribution. The Gaussian distribution (sometimes called the "Normal distribution") has a distinctive "bell curve" shape with a clear peak of values clustered around a center point (the mean), with the frequency of surrounding values trailing off to either side of the peak. The width of the peak (how spread out it is) is controlled by a second number, the "standard deviation" usually denoted by the lower-case greek symbol, sigma.
Gaussian distributions occur all over the place in nature, so they're a very good fit for many natural populations. I used a Gaussian random number generator in my Hosepipe particle system to give a more natural feel to the randomness of the particles.
The regular random number generator in .Net (System.Random) gives you a 'rectangular' distribution which means that all values are equally likely. To get a Gaussian from this rectangular distribution, we need to apply a little math. The full formula for the Gaussian distribution is complex, but we can get very close with an approximation. A little research found the formula I need: http://www.pitt.edu/~wpilib/statfaq/gaussfaq.html
It's fast, has no loops and has a very good approximation to the curve.
Without further ado, here's the code:
using System; namespace Utils { public class GaussianRandom : Random { /// <summary> /// Choose the next random number. The distribution of randomness is "normal" /// or "Gaussian" with a mean value dMu and a standard deviation of dSigma /// </summary> /// <param name="dMu">The center point (mean) of the normal distribution</param> /// <param name="dSigma">The standard deviation of the distribution</param> /// <returns>Random value</returns> public double NextGaussian(double dMu, double dSigma) { return dMu + CumulativeGaussian(base.NextDouble()) * dSigma; }
private static double CumulativeGaussian(double p) { // p is a rectangular probability between 0 and 1 // convert that into a gaussian. // Apply the inverse cumulative gaussian distribution function // This is an approximation by Abramowitz and Stegun; Press, et al. // See http://www.pitt.edu/~wpilib/statfaq/gaussfaq.html // Because of the symmetry of the normal // distribution, we need only consider 0 < p < 0.5. If you have p > 0.5, // then apply the algorithm below to q = 1-p, and then negate the value // for X obtained. bool fNegate = false;
if (p > 0.5) { p = 1.0 - p; fNegate = true; }
double t = Math.Sqrt(Math.Log(1.0 / (p * p))); double tt = t * t; double ttt = tt * t; double X = t - ((c_0 + c_1 * t + c_2 * tt) / (1 + d_1 * t + d_2 * tt + d_3 * ttt));
return fNegate ? -X : X; }
private const double c_0 = 2.515517; private const double c_1 = 0.802853; private const double c_2 = 0.010328; private const double d_1 = 1.432788; private const double d_2 = 0.189269; private const double d_3 = 0.001308; } } February 10 Hosepipe Particle SystemI couldn't resist updating another old Windows Forms application of mine to WPF. This particle system models the forces of gravity and air resistance for particles launched from a 'hosepipe' at the bottom of the screen. The direction of the hose is controlled by the mouse. Here's are a few screen shots: While I'm reasonably happy with the performance, I'm sure it's possible to do better. I ran the profiler on it once and fixed a bottleneck creating brushes, but still the majority of the time is spent drawing all those particles. The physics simulation is a breeze. Here's the XBAP: http://coad.net/users/pharring/Hose/HoseXBAP.xbap Once again, thanks to Noah for hosting. January 31 Hosting of the WPF maze GeneratorNoah Coad wrote with an offer I couldn't refuse. Noah wrote:
So, I spent a few minutes refactoring my main window so it could be hosted as a control in either a desktop WPF application or a web-based XBAP. Once I'd figured out the minimal set of ClickOnce security tags, I published it to Noah's FTP server and tried it out. Here's the fun part: I did all of this - including the upload - from my laptop on the bus during my commute home tonight. It was an unusually bad commute tonight because of heavy rain, and even though I was on a bus in the express lanes, it still took about twice as long as normal. The wireless access on the bus is not 100% reliable, but was working like a champ tonight. So, here's the link: http://coad.net/users/pharring/Maze/Maze.xbap Note: You need to have .Net 3.5 installed first (although, I think the ClickOnce settings will prompt you to install it if necessary). Noah also had a couple of suggestions to improve the application:
I like both of these suggestions, so these will probably be my next improvements. After that, though, I think it's time to give mazes a rest for a while. I'm looking for my next hobby project, but it might not be quite as graphically interesting :-) The source code: Public folders
|
Books I recommed
|
|||||||||
|
|