diff --git a/src/Triangle.Examples/Examples/Example10.cs b/src/Triangle.Examples/Examples/Example10.cs index 1bd276e..a02664b 100644 --- a/src/Triangle.Examples/Examples/Example10.cs +++ b/src/Triangle.Examples/Examples/Example10.cs @@ -1,133 +1,89 @@ -using System; -using System.Collections.Generic; -using System.Linq; -using TriangleNet.Geometry; -using TriangleNet.Meshing; -using TriangleNet.Meshing.Algorithm; -using TriangleNet.Rendering.Text; -using TriangleNet.Tools; - + namespace TriangleNet.Examples { + using System; + using System.Collections.Generic; + using System.Linq; + using TriangleNet.Geometry; + /// - /// Scattered data interpolation without USE_Z or USE_ATTRIBS. + /// Troubleshooting: finding degenerate boundary triangles. /// - internal class Example10 + public class Example9 { - // The function we are sampling. - private static readonly Func F = p => Math.Sin(p.X) * Math.Cos(p.Y); - - // The mesh size, for a structured grid (SIZE x SIZE) points. - private const int SIZE = 20; - public static bool Run(bool print = false) { - // The input domain. - var r = new Rectangle(0d, 0d, 10d, 10d); - - var mesh = GetScatteredDataMesh(r, out double[] data); - //var mesh = GetStructuredDataMesh(r, out double[] data); - - if (print) SvgImage.Save(mesh, "example-10.svg", 500); - - // The points to interpolate. - var xy = Generate.RandomPoints(50, r); - - var xyData = InterpolateData((Mesh)mesh, data, xy); - - double error = xy.Max(p => Math.Abs(xyData[p.ID] - F(p))); - - // L2 error - //double error = Math.Sqrt(xy.Sum(p => Math.Pow(xyData[p.ID] - F(p), 2))); - - // Define tolerance dependent on mesh dimensions and size. - double tolerance = 0.5 * Math.Max(r.Width, r.Height) / SIZE; - - return error < tolerance; - } - - private static IMesh GetStructuredDataMesh(Rectangle domain, out double[] data) - { - var mesh = GenericMesher.StructuredMesh(domain, SIZE, SIZE); - - mesh.Renumber(); - - // Generate function values for mesh points. - data = new double[mesh.Vertices.Count]; - - foreach (var item in mesh.Vertices) + var pts = new List { - data[item.ID] = F(item); - } + // The 4 corners of the rectangle. + new Vertex(1.5, 1.0), + new Vertex(1.5, -1.0), + new Vertex(-1.5, -1.0), + new Vertex(-1.5, 1.0), - return mesh; - } + // The edge midpoints. + new Vertex(0.0, 1.0), + new Vertex(0.0, -1.0), + new Vertex(1.5, 0.0), + new Vertex(-1.5, 0.0) + }; - private static IMesh GetScatteredDataMesh(Rectangle domain, out double[] data) - { - var r = new Rectangle(domain); + var r = new Random(78403); - double h = domain.Width / SIZE; + // The original rectangle. + var poly = Rotate(pts, 0); - // Generate a rectangle boundary point set (20 points on each side). - var input = Generate.Rectangle(r, 0.5); - - // Making sure we add some margin to the boundary. - h = -h / 2; - r.Resize(h, h); - - int n = Math.Max(1, SIZE * SIZE - input.Points.Count); - - // Add more input points (more sampling points, better interpolation). - input.Points.AddRange(Generate.RandomPoints(n, r)); - - var mesher = new GenericMesher(new Dwyer()); - - // Generate mesh. - var mesh = mesher.Triangulate(input.Points); - - mesh.Renumber(); - - // Generate function values for mesh points. - data = new double[mesh.Vertices.Count]; - - foreach (var item in mesh.Vertices) + for (int i = 0; i < 10; i++) { - data[item.ID] = F(item); - } + var mesh = poly.Triangulate(); - return mesh; - } + var list = MeshValidator.GetDegenerateBoundaryTriangles(mesh); - private static double[] InterpolateData(Mesh mesh, double[] data, IEnumerable xy) - { - // The interpolated values. - var values = new double[xy.Count()]; - - var qtree = new TriangleQuadTree(mesh); - - int i = 0; - - foreach (var p in xy) - { - var tri = qtree.Query(p.X, p.Y); - - // For easy access of the interpolated values. - p.ID = i; - - if (tri == null) + if (print && list.Any()) { - values[i] = float.NaN; - } - else - { - values[i] = Interpolation.InterpolatePoint(tri, p, data); + Console.WriteLine("Iteration {0}: found {1} degenerate triangle(s) of {2}.", + i, list.Count(), mesh.Triangles.Count); + + foreach (var t in list) + { + Console.WriteLine(" [{0} {1} {2}]", + t.GetVertexID(0), + t.GetVertexID(1), + t.GetVertexID(2)); + } } - i++; + // Random rotation. + poly = Rotate(pts, Math.PI * r.NextDouble()); } - return values; + return true; + } + + /// + /// Rotate given point set around the origin. + /// + private static IPolygon Rotate(List points, double radians) + { + var poly = new Polygon(points.Count); + + int id = 0; + + foreach (var p in points) + { + double x = p.X; + double y = p.Y; + + double s = Math.Sin(radians); + double c = Math.Cos(radians); + + double xr = c * x - s * y; + double yr = s * x + c * y; + + poly.Points.Add(new Vertex(xr, yr) { ID = id++ }); + } + + return poly; } } -} \ No newline at end of file +} diff --git a/src/Triangle.Examples/Examples/Example11.cs b/src/Triangle.Examples/Examples/Example11.cs new file mode 100644 index 0000000..1bd276e --- /dev/null +++ b/src/Triangle.Examples/Examples/Example11.cs @@ -0,0 +1,133 @@ +using System; +using System.Collections.Generic; +using System.Linq; +using TriangleNet.Geometry; +using TriangleNet.Meshing; +using TriangleNet.Meshing.Algorithm; +using TriangleNet.Rendering.Text; +using TriangleNet.Tools; + +namespace TriangleNet.Examples +{ + /// + /// Scattered data interpolation without USE_Z or USE_ATTRIBS. + /// + internal class Example10 + { + // The function we are sampling. + private static readonly Func F = p => Math.Sin(p.X) * Math.Cos(p.Y); + + // The mesh size, for a structured grid (SIZE x SIZE) points. + private const int SIZE = 20; + + public static bool Run(bool print = false) + { + // The input domain. + var r = new Rectangle(0d, 0d, 10d, 10d); + + var mesh = GetScatteredDataMesh(r, out double[] data); + //var mesh = GetStructuredDataMesh(r, out double[] data); + + if (print) SvgImage.Save(mesh, "example-10.svg", 500); + + // The points to interpolate. + var xy = Generate.RandomPoints(50, r); + + var xyData = InterpolateData((Mesh)mesh, data, xy); + + double error = xy.Max(p => Math.Abs(xyData[p.ID] - F(p))); + + // L2 error + //double error = Math.Sqrt(xy.Sum(p => Math.Pow(xyData[p.ID] - F(p), 2))); + + // Define tolerance dependent on mesh dimensions and size. + double tolerance = 0.5 * Math.Max(r.Width, r.Height) / SIZE; + + return error < tolerance; + } + + private static IMesh GetStructuredDataMesh(Rectangle domain, out double[] data) + { + var mesh = GenericMesher.StructuredMesh(domain, SIZE, SIZE); + + mesh.Renumber(); + + // Generate function values for mesh points. + data = new double[mesh.Vertices.Count]; + + foreach (var item in mesh.Vertices) + { + data[item.ID] = F(item); + } + + return mesh; + } + + private static IMesh GetScatteredDataMesh(Rectangle domain, out double[] data) + { + var r = new Rectangle(domain); + + double h = domain.Width / SIZE; + + // Generate a rectangle boundary point set (20 points on each side). + var input = Generate.Rectangle(r, 0.5); + + // Making sure we add some margin to the boundary. + h = -h / 2; + r.Resize(h, h); + + int n = Math.Max(1, SIZE * SIZE - input.Points.Count); + + // Add more input points (more sampling points, better interpolation). + input.Points.AddRange(Generate.RandomPoints(n, r)); + + var mesher = new GenericMesher(new Dwyer()); + + // Generate mesh. + var mesh = mesher.Triangulate(input.Points); + + mesh.Renumber(); + + // Generate function values for mesh points. + data = new double[mesh.Vertices.Count]; + + foreach (var item in mesh.Vertices) + { + data[item.ID] = F(item); + } + + return mesh; + } + + private static double[] InterpolateData(Mesh mesh, double[] data, IEnumerable xy) + { + // The interpolated values. + var values = new double[xy.Count()]; + + var qtree = new TriangleQuadTree(mesh); + + int i = 0; + + foreach (var p in xy) + { + var tri = qtree.Query(p.X, p.Y); + + // For easy access of the interpolated values. + p.ID = i; + + if (tri == null) + { + values[i] = float.NaN; + } + else + { + values[i] = Interpolation.InterpolatePoint(tri, p, data); + } + + i++; + } + + return values; + } + } +} \ No newline at end of file diff --git a/src/Triangle.Examples/Examples/Example2.cs b/src/Triangle.Examples/Examples/Example2.cs deleted file mode 100644 index c497ccf..0000000 --- a/src/Triangle.Examples/Examples/Example2.cs +++ /dev/null @@ -1,48 +0,0 @@ -namespace TriangleNet.Examples -{ - using TriangleNet.Geometry; - using TriangleNet.Meshing; - using TriangleNet.Rendering.Text; - - /// - /// Triangulate a polygon with hole and set minimum angle constraint. - /// - public static class Example2 - { - public static bool Run(bool print = false) - { - // Generate the input geometry. - var poly = CreatePolygon(); - - // Set minimum angle quality option. - var quality = new QualityOptions() { MinimumAngle = 30.0 }; - - // Generate mesh using the polygons Triangulate extension method. - var mesh = poly.Triangulate(quality); - - if (print) SvgImage.Save(mesh, "example-2.svg", 500); - - return mesh.Triangles.Count > 0; - } - - public static IPolygon CreatePolygon(double h = 0.2) - { - // Generate the input geometry. - var poly = new Polygon(); - - // Center point. - var center = new Point(0, 0); - - // Inner contour (hole). - poly.Add(Generate.Circle(1.0, center, h, 1), center); - - // Internal contour. - poly.Add(Generate.Circle(2.0, center, h, 2)); - - // Outer contour. - poly.Add(Generate.Circle(3.0, center, h, 3)); - - return poly; - } - } -} \ No newline at end of file diff --git a/src/Triangle.Examples/Examples/Example3.cs b/src/Triangle.Examples/Examples/Example3.cs index 7e0182d..c497ccf 100644 --- a/src/Triangle.Examples/Examples/Example3.cs +++ b/src/Triangle.Examples/Examples/Example3.cs @@ -1,54 +1,48 @@ - -namespace TriangleNet.Examples +namespace TriangleNet.Examples { - using System; using TriangleNet.Geometry; using TriangleNet.Meshing; using TriangleNet.Rendering.Text; - using TriangleNet.Smoothing; /// - /// Triangulate a polygon with hole with maximum area constraint, followed by mesh smoothing. + /// Triangulate a polygon with hole and set minimum angle constraint. /// - public class Example3 + public static class Example2 { public static bool Run(bool print = false) { - // Generate mesh. - var mesh = CreateMesh(); + // Generate the input geometry. + var poly = CreatePolygon(); - if (print) SvgImage.Save(mesh, "example-3.svg", 500); + // Set minimum angle quality option. + var quality = new QualityOptions() { MinimumAngle = 30.0 }; + + // Generate mesh using the polygons Triangulate extension method. + var mesh = poly.Triangulate(quality); + + if (print) SvgImage.Save(mesh, "example-2.svg", 500); return mesh.Triangles.Count > 0; } - public static IMesh CreateMesh() + public static IPolygon CreatePolygon(double h = 0.2) { // Generate the input geometry. - var poly = Example2.CreatePolygon(); + var poly = new Polygon(); - // Since we want to do CVT smoothing, ensure that the mesh - // is conforming Delaunay. - var options = new ConstraintOptions() { ConformingDelaunay = true }; + // Center point. + var center = new Point(0, 0); - // Set maximum area quality option (we don't need to set a minimum - // angle, since smoothing will improve the triangle shapes). - var quality = new QualityOptions() - { - // The boundary segments have a length of 0.2, so we set a - // maximum area constraint assuming equilateral triangles. - MaximumArea = (Math.Sqrt(3) / 4 * 0.2 * 0.2) * 1.45 - }; + // Inner contour (hole). + poly.Add(Generate.Circle(1.0, center, h, 1), center); - // Generate mesh using the polygons Triangulate extension method. - var mesh = poly.Triangulate(options, quality); + // Internal contour. + poly.Add(Generate.Circle(2.0, center, h, 2)); - var smoother = new SimpleSmoother(); + // Outer contour. + poly.Add(Generate.Circle(3.0, center, h, 3)); - // Smooth mesh. - smoother.Smooth(mesh, 25); - - return mesh; + return poly; } } -} +} \ No newline at end of file diff --git a/src/Triangle.Examples/Examples/Example4.cs b/src/Triangle.Examples/Examples/Example4.cs index d3c122a..7e0182d 100644 --- a/src/Triangle.Examples/Examples/Example4.cs +++ b/src/Triangle.Examples/Examples/Example4.cs @@ -1,71 +1,54 @@  namespace TriangleNet.Examples { - using TriangleNet; + using System; using TriangleNet.Geometry; using TriangleNet.Meshing; using TriangleNet.Rendering.Text; using TriangleNet.Smoothing; /// - /// Refine only a part of a polygon mesh by using region pointers and an area constraint. + /// Triangulate a polygon with hole with maximum area constraint, followed by mesh smoothing. /// - public class Example4 + public class Example3 { public static bool Run(bool print = false) { - // Generate the input geometry. - var poly = CreatePolygon(); - - // Define regions (first one defines the area constraint). - poly.Regions.Add(new RegionPointer(1.5, 0.0, 1, 0.01)); - poly.Regions.Add(new RegionPointer(2.5, 0.0, 2)); + // Generate mesh. + var mesh = CreateMesh(); - // Set quality and constraint options. - var options = new ConstraintOptions() - { - ConformingDelaunay = true - }; - - var quality = new QualityOptions() - { - MinimumAngle = 25.0, - VariableArea = true - }; - - //quality.UserTest = (t, area) => t.Label == 1 && area > 0.01; - - var mesh = poly.Triangulate(options, quality); - - var smoother = new SimpleSmoother(); - - smoother.Smooth(mesh, 5); - - if (print) SvgImage.Save(mesh, "example-4.svg", 500); + if (print) SvgImage.Save(mesh, "example-3.svg", 500); return mesh.Triangles.Count > 0; } - public static IPolygon CreatePolygon() + public static IMesh CreateMesh() { - // Generate three concentric circles. - var poly = new Polygon(); + // Generate the input geometry. + var poly = Example2.CreatePolygon(); - // Center point. - var center = new Point(0, 0); + // Since we want to do CVT smoothing, ensure that the mesh + // is conforming Delaunay. + var options = new ConstraintOptions() { ConformingDelaunay = true }; - // Inner contour (hole). - poly.Add(Generate.Circle(1.0, center, 0.1, 1), center); + // Set maximum area quality option (we don't need to set a minimum + // angle, since smoothing will improve the triangle shapes). + var quality = new QualityOptions() + { + // The boundary segments have a length of 0.2, so we set a + // maximum area constraint assuming equilateral triangles. + MaximumArea = (Math.Sqrt(3) / 4 * 0.2 * 0.2) * 1.45 + }; - // Internal contour. - poly.Add(Generate.Circle(2.0, center, 0.1, 2)); + // Generate mesh using the polygons Triangulate extension method. + var mesh = poly.Triangulate(options, quality); - // Outer contour. - poly.Add(Generate.Circle(3.0, center, 0.3, 3)); + var smoother = new SimpleSmoother(); - // Note that the outer contour has a larger segment size! + // Smooth mesh. + smoother.Smooth(mesh, 25); - return poly; + return mesh; } } } diff --git a/src/Triangle.Examples/Examples/Example5.cs b/src/Triangle.Examples/Examples/Example5.cs index 17540af..d3c122a 100644 --- a/src/Triangle.Examples/Examples/Example5.cs +++ b/src/Triangle.Examples/Examples/Example5.cs @@ -1,98 +1,71 @@  namespace TriangleNet.Examples { - using System.Collections.Generic; using TriangleNet; using TriangleNet.Geometry; using TriangleNet.Meshing; - using TriangleNet.Meshing.Iterators; using TriangleNet.Rendering.Text; + using TriangleNet.Smoothing; /// - /// Two ways finding boundary triangles. + /// Refine only a part of a polygon mesh by using region pointers and an area constraint. /// - public static class Example5 + public class Example4 { public static bool Run(bool print = false) { - var mesh = Example3.CreateMesh(); + // Generate the input geometry. + var poly = CreatePolygon(); + + // Define regions (first one defines the area constraint). + poly.Regions.Add(new RegionPointer(1.5, 0.0, 1, 0.01)); + poly.Regions.Add(new RegionPointer(2.5, 0.0, 2)); - FindBoundary1(mesh); + // Set quality and constraint options. + var options = new ConstraintOptions() + { + ConformingDelaunay = true + }; - if (print) SvgImage.Save(mesh, "example-5-1.svg", 500, true, false); + var quality = new QualityOptions() + { + MinimumAngle = 25.0, + VariableArea = true + }; - FindBoundary2(mesh); + //quality.UserTest = (t, area) => t.Label == 1 && area > 0.01; - if (print) SvgImage.Save(mesh, "example-5-2.svg", 500, true, false); + var mesh = poly.Triangulate(options, quality); + + var smoother = new SimpleSmoother(); + + smoother.Smooth(mesh, 5); + + if (print) SvgImage.Save(mesh, "example-4.svg", 500); return mesh.Triangles.Count > 0; } - /// - /// Find boundary triangles using segments. - /// - private static void FindBoundary1(IMesh mesh, bool neigbours = true) + public static IPolygon CreatePolygon() { - mesh.Renumber(); + // Generate three concentric circles. + var poly = new Polygon(); - var cache = new List(mesh.Segments.Count + 1); + // Center point. + var center = new Point(0, 0); - var circulator = new VertexCirculator((Mesh)mesh); + // Inner contour (hole). + poly.Add(Generate.Circle(1.0, center, 0.1, 1), center); - foreach (var s in mesh.Segments) - { - int label = s.Label; + // Internal contour. + poly.Add(Generate.Circle(2.0, center, 0.1, 2)); - for (int i = 0; i < 2; i++) - { - var vertex = s.GetVertex(i); + // Outer contour. + poly.Add(Generate.Circle(3.0, center, 0.3, 3)); - // Check the vertex ID to see if it was processed already. - if (vertex.ID >= 0) - { - var star = circulator.EnumerateTriangles(vertex); + // Note that the outer contour has a larger segment size! - foreach (var triangle in star) - { - triangle.Label = label; - } - - // Mark the vertex as "processed". - vertex.ID = -vertex.ID; - - cache.Add(vertex); - } - } - } - - // Undo the vertex ID changes. - foreach (var vertex in cache) - { - vertex.ID = -vertex.ID; - } - } - - /// - /// Find boundary triangles using vertices. - /// - private static void FindBoundary2(IMesh mesh) - { - var circulator = new VertexCirculator((Mesh)mesh); - - foreach (var vertex in mesh.Vertices) - { - int label = vertex.Label; - - if (label > 0) - { - var star = circulator.EnumerateTriangles(vertex); - - foreach (var triangle in star) - { - triangle.Label = label; - } - } - } + return poly; } } } diff --git a/src/Triangle.Examples/Examples/Example6.cs b/src/Triangle.Examples/Examples/Example6.cs index 3fccd89..17540af 100644 --- a/src/Triangle.Examples/Examples/Example6.cs +++ b/src/Triangle.Examples/Examples/Example6.cs @@ -1,56 +1,98 @@ -namespace TriangleNet.Examples + +namespace TriangleNet.Examples { - using System.Linq; + using System.Collections.Generic; using TriangleNet; using TriangleNet.Geometry; + using TriangleNet.Meshing; using TriangleNet.Meshing.Iterators; using TriangleNet.Rendering.Text; - using TriangleNet.Tools; - using TriangleNet.Topology; /// - /// Boolean operations on mesh regions (intersection, difference, xor). + /// Two ways finding boundary triangles. /// - public static class Example6 + public static class Example5 { public static bool Run(bool print = false) { - // Generate the input geometry. - var polygon = new Polygon(8, true); + var mesh = Example3.CreateMesh(); - // Two intersecting rectangles. - var A = Generate.Rectangle(0d, 0d, 4d, 4d, label: 1); - var B = Generate.Rectangle(1d, 1d, 4d, 4d, label: 2); + FindBoundary1(mesh); - polygon.Add(A); - polygon.Add(B); + if (print) SvgImage.Save(mesh, "example-5-1.svg", 500, true, false); - // Generate mesh. - var mesh = (Mesh)polygon.Triangulate(); + FindBoundary2(mesh); - if (print) SvgImage.Save(mesh, "example-6.svg", 500); + if (print) SvgImage.Save(mesh, "example-5-2.svg", 500, true, false); - // Find a seeding triangle (in this case, the point (2, 2) lies in - // both rectangles). - var seed = (new TriangleQuadTree(mesh)).Query(2.0, 2.0) as Triangle; + return mesh.Triangles.Count > 0; + } - var iterator = new RegionIterator(mesh); + /// + /// Find boundary triangles using segments. + /// + private static void FindBoundary1(IMesh mesh, bool neigbours = true) + { + mesh.Renumber(); - iterator.Process(seed, t => t.Label ^= 1, 1); - iterator.Process(seed, t => t.Label ^= 2, 2); + var cache = new List(mesh.Segments.Count + 1); - // At this point, all triangles will have label 1, 2 or 3 (= 1 xor 2). + var circulator = new VertexCirculator((Mesh)mesh); - // The intersection of A and B. - var intersection = mesh.Triangles.Where(t => t.Label == 3); + foreach (var s in mesh.Segments) + { + int label = s.Label; - // The difference A \ B. - var difference = mesh.Triangles.Where(t => t.Label == 1); + for (int i = 0; i < 2; i++) + { + var vertex = s.GetVertex(i); - // The xor of A and B. - var xor = mesh.Triangles.Where(t => t.Label == 1 || t.Label == 2); + // Check the vertex ID to see if it was processed already. + if (vertex.ID >= 0) + { + var star = circulator.EnumerateTriangles(vertex); - return intersection.Any() && difference.Any() && xor.Any(); + foreach (var triangle in star) + { + triangle.Label = label; + } + + // Mark the vertex as "processed". + vertex.ID = -vertex.ID; + + cache.Add(vertex); + } + } + } + + // Undo the vertex ID changes. + foreach (var vertex in cache) + { + vertex.ID = -vertex.ID; + } + } + + /// + /// Find boundary triangles using vertices. + /// + private static void FindBoundary2(IMesh mesh) + { + var circulator = new VertexCirculator((Mesh)mesh); + + foreach (var vertex in mesh.Vertices) + { + int label = vertex.Label; + + if (label > 0) + { + var star = circulator.EnumerateTriangles(vertex); + + foreach (var triangle in star) + { + triangle.Label = label; + } + } + } } } -} \ No newline at end of file +} diff --git a/src/Triangle.Examples/Examples/Example7.cs b/src/Triangle.Examples/Examples/Example7.cs index cdc2f8f..3fccd89 100644 --- a/src/Triangle.Examples/Examples/Example7.cs +++ b/src/Triangle.Examples/Examples/Example7.cs @@ -1,73 +1,56 @@ - -namespace TriangleNet.Examples +namespace TriangleNet.Examples { - using System; + using System.Linq; + using TriangleNet; using TriangleNet.Geometry; - using TriangleNet.Meshing; using TriangleNet.Meshing.Iterators; using TriangleNet.Rendering.Text; + using TriangleNet.Tools; + using TriangleNet.Topology; /// - /// Using a user test function to define a maximum edge length constraint. + /// Boolean operations on mesh regions (intersection, difference, xor). /// - public static class Example7 + public static class Example6 { - const double MAX_EDGE_LENGTH = 0.2; - public static bool Run(bool print = false) { - var poly = new Polygon(); - // Generate the input geometry. - poly.Add(Generate.Rectangle(0.0, 0.0, 1.0, 1.0)); + var polygon = new Polygon(8, true); - // Set minimum angle quality option, ignoring holes. - var quality = new QualityOptions() - { - UserTest = MaxEdgeLength - }; + // Two intersecting rectangles. + var A = Generate.Rectangle(0d, 0d, 4d, 4d, label: 1); + var B = Generate.Rectangle(1d, 1d, 4d, 4d, label: 2); - // Generate mesh using the polygons Triangulate extension method. - var mesh = (Mesh)poly.Triangulate(quality); + polygon.Add(A); + polygon.Add(B); - // Validate. - foreach (var e in EdgeIterator.EnumerateEdges(mesh)) - { - double length = Math.Sqrt(DistSqr(e.GetVertex(0), e.GetVertex(1))); + // Generate mesh. + var mesh = (Mesh)polygon.Triangulate(); - if (length > MAX_EDGE_LENGTH) - { - return false; - } - } + if (print) SvgImage.Save(mesh, "example-6.svg", 500); - if (print) SvgImage.Save(mesh, "example-7.svg", 500); + // Find a seeding triangle (in this case, the point (2, 2) lies in + // both rectangles). + var seed = (new TriangleQuadTree(mesh)).Query(2.0, 2.0) as Triangle; - return true; - } + var iterator = new RegionIterator(mesh); - static bool MaxEdgeLength(ITriangle tri, double area) - { - var p0 = tri.GetVertex(0); - var p1 = tri.GetVertex(1); - var p2 = tri.GetVertex(2); + iterator.Process(seed, t => t.Label ^= 1, 1); + iterator.Process(seed, t => t.Label ^= 2, 2); - var s1 = DistSqr(p0, p1); - var s2 = DistSqr(p1, p2); - var s3 = DistSqr(p2, p0); + // At this point, all triangles will have label 1, 2 or 3 (= 1 xor 2). - // Comparing against squared max leg length. - var maxlen = MAX_EDGE_LENGTH * MAX_EDGE_LENGTH; + // The intersection of A and B. + var intersection = mesh.Triangles.Where(t => t.Label == 3); - return s1 > maxlen || s2 > maxlen || s3 > maxlen; - } + // The difference A \ B. + var difference = mesh.Triangles.Where(t => t.Label == 1); - static double DistSqr(Vertex a, Vertex b) - { - var dx = a.X - b.X; - var dy = a.Y - b.Y; + // The xor of A and B. + var xor = mesh.Triangles.Where(t => t.Label == 1 || t.Label == 2); - return dx * dx + dy * dy; + return intersection.Any() && difference.Any() && xor.Any(); } } -} +} \ No newline at end of file diff --git a/src/Triangle.Examples/Examples/Example8.cs b/src/Triangle.Examples/Examples/Example8.cs index f0a16a5..cdc2f8f 100644 --- a/src/Triangle.Examples/Examples/Example8.cs +++ b/src/Triangle.Examples/Examples/Example8.cs @@ -2,83 +2,72 @@ namespace TriangleNet.Examples { using System; - using System.Collections.Generic; - using TriangleNet; + using TriangleNet.Geometry; + using TriangleNet.Meshing; using TriangleNet.Meshing.Iterators; - using TriangleNet.Tools; + using TriangleNet.Rendering.Text; /// - /// Compute the adjacency matrix of the mesh vertices. + /// Using a user test function to define a maximum edge length constraint. /// - public class Example8 + public static class Example7 { - public static bool Run() + const double MAX_EDGE_LENGTH = 0.2; + + public static bool Run(bool print = false) { - var mesh = (Mesh)Example3.CreateMesh(); + var poly = new Polygon(); - return FindAdjacencyMatrix(mesh); - } + // Generate the input geometry. + poly.Add(Generate.Rectangle(0.0, 0.0, 1.0, 1.0)); - private static bool FindAdjacencyMatrix(Mesh mesh) - { - mesh.Renumber(); - - var ap = new List(mesh.Vertices.Count); // Column pointers. - var ai = new List(4 * mesh.Vertices.Count); // Row indices. - - var circulator = new VertexCirculator(mesh); - - int k = 0; - - foreach (var vertex in mesh.Vertices) + // Set minimum angle quality option, ignoring holes. + var quality = new QualityOptions() { - var star = circulator.EnumerateVertices(vertex); + UserTest = MaxEdgeLength + }; - ap.Add(k); + // Generate mesh using the polygons Triangulate extension method. + var mesh = (Mesh)poly.Triangulate(quality); - // Each vertex is adjacent to itself. - ai.Add(vertex.ID); - k++; - - foreach (var item in star) - { - ai.Add(item.ID); - k++; - } - } - - ap.Add(k); - - var matrix1 = new AdjacencyMatrix(ap.ToArray(), ai.ToArray()); - var matrix2 = new AdjacencyMatrix(mesh); - - // Column pointers should be exactly the same. - if (!CompareArray(matrix1.ColumnPointers, matrix2.ColumnPointers)) + // Validate. + foreach (var e in EdgeIterator.EnumerateEdges(mesh)) { - return false; - } + double length = Math.Sqrt(DistSqr(e.GetVertex(0), e.GetVertex(1))); - return true; - } - - private static bool CompareArray(int[] a, int[] b) - { - int length = a.Length; - - if (b.Length != length) - { - return false; - } - - for (int i = 0; i < length; i++) - { - if (a[i] != b[i]) + if (length > MAX_EDGE_LENGTH) { return false; } } + if (print) SvgImage.Save(mesh, "example-7.svg", 500); + return true; } + + static bool MaxEdgeLength(ITriangle tri, double area) + { + var p0 = tri.GetVertex(0); + var p1 = tri.GetVertex(1); + var p2 = tri.GetVertex(2); + + var s1 = DistSqr(p0, p1); + var s2 = DistSqr(p1, p2); + var s3 = DistSqr(p2, p0); + + // Comparing against squared max leg length. + var maxlen = MAX_EDGE_LENGTH * MAX_EDGE_LENGTH; + + return s1 > maxlen || s2 > maxlen || s3 > maxlen; + } + + static double DistSqr(Vertex a, Vertex b) + { + var dx = a.X - b.X; + var dy = a.Y - b.Y; + + return dx * dx + dy * dy; + } } } diff --git a/src/Triangle.Examples/Examples/Example9.cs b/src/Triangle.Examples/Examples/Example9.cs index a02664b..f0a16a5 100644 --- a/src/Triangle.Examples/Examples/Example9.cs +++ b/src/Triangle.Examples/Examples/Example9.cs @@ -3,87 +3,82 @@ namespace TriangleNet.Examples { using System; using System.Collections.Generic; - using System.Linq; - using TriangleNet.Geometry; + using TriangleNet; + using TriangleNet.Meshing.Iterators; + using TriangleNet.Tools; /// - /// Troubleshooting: finding degenerate boundary triangles. + /// Compute the adjacency matrix of the mesh vertices. /// - public class Example9 + public class Example8 { - public static bool Run(bool print = false) + public static bool Run() { - var pts = new List + var mesh = (Mesh)Example3.CreateMesh(); + + return FindAdjacencyMatrix(mesh); + } + + private static bool FindAdjacencyMatrix(Mesh mesh) + { + mesh.Renumber(); + + var ap = new List(mesh.Vertices.Count); // Column pointers. + var ai = new List(4 * mesh.Vertices.Count); // Row indices. + + var circulator = new VertexCirculator(mesh); + + int k = 0; + + foreach (var vertex in mesh.Vertices) { - // The 4 corners of the rectangle. - new Vertex(1.5, 1.0), - new Vertex(1.5, -1.0), - new Vertex(-1.5, -1.0), - new Vertex(-1.5, 1.0), + var star = circulator.EnumerateVertices(vertex); - // The edge midpoints. - new Vertex(0.0, 1.0), - new Vertex(0.0, -1.0), - new Vertex(1.5, 0.0), - new Vertex(-1.5, 0.0) - }; + ap.Add(k); - var r = new Random(78403); + // Each vertex is adjacent to itself. + ai.Add(vertex.ID); + k++; - // The original rectangle. - var poly = Rotate(pts, 0); - - for (int i = 0; i < 10; i++) - { - var mesh = poly.Triangulate(); - - var list = MeshValidator.GetDegenerateBoundaryTriangles(mesh); - - if (print && list.Any()) + foreach (var item in star) { - Console.WriteLine("Iteration {0}: found {1} degenerate triangle(s) of {2}.", - i, list.Count(), mesh.Triangles.Count); - - foreach (var t in list) - { - Console.WriteLine(" [{0} {1} {2}]", - t.GetVertexID(0), - t.GetVertexID(1), - t.GetVertexID(2)); - } + ai.Add(item.ID); + k++; } + } - // Random rotation. - poly = Rotate(pts, Math.PI * r.NextDouble()); + ap.Add(k); + + var matrix1 = new AdjacencyMatrix(ap.ToArray(), ai.ToArray()); + var matrix2 = new AdjacencyMatrix(mesh); + + // Column pointers should be exactly the same. + if (!CompareArray(matrix1.ColumnPointers, matrix2.ColumnPointers)) + { + return false; } return true; } - /// - /// Rotate given point set around the origin. - /// - private static IPolygon Rotate(List points, double radians) + private static bool CompareArray(int[] a, int[] b) { - var poly = new Polygon(points.Count); + int length = a.Length; - int id = 0; - - foreach (var p in points) + if (b.Length != length) { - double x = p.X; - double y = p.Y; - - double s = Math.Sin(radians); - double c = Math.Cos(radians); - - double xr = c * x - s * y; - double yr = s * x + c * y; - - poly.Points.Add(new Vertex(xr, yr) { ID = id++ }); + return false; } - return poly; + for (int i = 0; i < length; i++) + { + if (a[i] != b[i]) + { + return false; + } + } + + return true; } } }