diff --git a/src/Spatial.Tests/Euclidean/Circle2DTests.cs b/src/Spatial.Tests/Euclidean/Circle2DTests.cs index 7d3d0df..d1a2cdd 100644 --- a/src/Spatial.Tests/Euclidean/Circle2DTests.cs +++ b/src/Spatial.Tests/Euclidean/Circle2DTests.cs @@ -1,4 +1,5 @@ -using System; +using System; +using System.Linq; using MathNet.Spatial.Euclidean; using NUnit.Framework; @@ -58,5 +59,125 @@ public void CircleFromThreePointsArgumentException() Assert.Throws(() => { Circle2D.FromPoints(p1, p2, p3); }); } + + //parallel to the X-axis + [TestCase("0,0", 1, "-10,-10", "+10,-10", null)] + [TestCase("0,0", 1, "-10,-1", "+10,-1", "0,-1")] + [TestCase("0,0", 1, "-10,0", "+10,0", "+1,0;-1,0")] + [TestCase("0,0", 1, "-10,+1", "+10,+1", "0,+1")] + [TestCase("0,0", 1, "-10,+10", "+10,+10", null)] + //parallel to the Y-axis + [TestCase("0,0", 1, "-10,-10", "-10,+10", null)] + [TestCase("0,0", 1, "-1,-10", "-1,+10", "-1,0")] + [TestCase("0,0", 1, "0,-10", "0,+10", "0,+1;0,-1")] + [TestCase("0,0", 1, "+1,-10", "+1,+10", "+1,0")] + [TestCase("0,0", 1, "+10,-10", "+10,+10", null)] + //general cases + [TestCase("0,0", 1, "-10,+10", "+10,+10", null)] + [TestCase("0,0", 1, "-1.414213562373095,0", "0,+1.414213562373095", "-0.707,0.707")] + [TestCase("0,0", 1, "-10,-10", "+10,+10", "+0.707,+0.707;-0.707,-0.707")] + [TestCase("0,0", 1, "0,-1.41421356", "+1.41421356,0", "+0.707,-0.707")] + [TestCase("0,0", 1, "0,-10", "+10,0", null)] + public void CircleIntersectWithLine2D(string sc, double radius, string sps, string spe, string intersections) + { + var circle = new Circle2D(Point2D.Parse(sc), radius); + var line = new Line2D(Point2D.Parse(sps), Point2D.Parse(spe)); + + var actual = circle.IntersectWith(line); + + var expected = parseToPointsArray(intersections); + for (int i = 0; i < Math.Min(actual.Length, expected.Length); i++) + { + var a = actual[i]; + var e = expected[i]; + AssertGeometry.AreEqual(a, e, 1e-3); //needs to fix for the default tolerance + } + } + //parallel to X-axis + ////segment contains the all intersections(same to the cases of circle and line) + [TestCase("0,0", 1, "-10,+10", "+10,+10", null)] + [TestCase("0,0", 1, "-10,+1", "+10,+1", "0,+1")] + [TestCase("0,0", 1, "-10,0", "+10,0", "+1,0;-1,0")] + [TestCase("0,0", 1, "-10,-1", "+10,-1", "0,-1")] + [TestCase("0,0", 1, "-10,-10", "+10,-10", null)] + ////segments cross the circle's contour just 1 time + [TestCase("0,0", 1, "+0,+10", "+10,+10", null)] + [TestCase("0,0", 1, "+0,+1", "+10,+1", "0,1")] + [TestCase("0,0", 1, "+0,+0", "+10,+0", "1,0")] + [TestCase("0,0", 1, "+0,-1", "+10,-1", "0,-1")] + [TestCase("0,0", 1, "+0,-10", "+10,-10", null)] + ////segment contains no intersections(px of the startingPoint is too big to intersect with the circle) + [TestCase("0,0", 1, "+10,+10", "+100,+10", null)] + [TestCase("0,0", 1, "+10,+1", "+100,+1", null)] + [TestCase("0,0", 1, "+10,+0", "+100,0", null)] + [TestCase("0,0", 1, "+10,-1", "+100,-1", null)] + [TestCase("0,0", 1, "+10,-10", "+100,-10", null)] + //parallel to Y-axis + ////segment contains the all intersections(same to the cases of circle and line) + [TestCase("0,0", 1, "-10,-10", "-10,+10", null)] + [TestCase("0,0", 1, "-1,-10", "-1,+10", "-1,0")] + [TestCase("0,0", 1, "+0,-10", "+0,+10", "0,+1;0,-1")] + [TestCase("0,0", 1, "+1,-10", "+1,+10", "+1,0")] + [TestCase("0,0", 1, "+10,-10", "+10,+10", null)] + ////segments cross the circle's contour just 1 time + [TestCase("0,0", 1, "+10,0", "+10,+10", null)] + [TestCase("0,0", 1, "+1,0", "+1,+10", "+1,0")] + [TestCase("0,0", 1, "+0,0", "+0,+10", "0,+1")] + [TestCase("0,0", 1, "-1,0", "-1,+10", "-1,0")] + [TestCase("0,0", 1, "-10,0", "-10,+10", null)] + ////segment contains no intersections(py of the startingPoint is too big to intersect with the circle) + [TestCase("0,0", 1, "+10,+10", "+10,+100", null)] + [TestCase("0,0", 1, "+1,+10", "+1,+100", null)] + [TestCase("0,0", 1, "+0,+10", "+0,+100", null)] + [TestCase("0,0", 1, "-1,+10", "-1,+100", null)] + [TestCase("0,0", 1, "-10,+10", "-10,+100", null)] + //general cases + ////segment contains the all intersections(same to the cases of circle and line) + [TestCase("0,0", 1, "-10,+10", "+10,+10", null)] + [TestCase("0,0", 1, "-1.414213562373095,0", "0,+1.414213562373095", "-0.707,0.707")] + [TestCase("0,0", 1, "-10,-10", "+10,+10", "+0.707,+0.707;-0.707,-0.707")] + [TestCase("0,0", 1, "0,-1.41421356", "+1.41421356,0", "+0.707,-0.707")] + [TestCase("0,0", 1, "0,-10", "+10,0", null)] + ////segments cross the circle's contour just 1 time + [TestCase("0,0", 1, "+10,0", "+10,+10", null)] + [TestCase("0,0", 1, "+1,0", "+1,+10", "+1,0")] + [TestCase("0,0", 1, "+0,0", "+0,+10", "0,+1")] + [TestCase("0,0", 1, "-1,0", "-1,+10", "-1,0")] + [TestCase("0,0", 1, "-10,0", "-10,+10", null)] + ////segment contains no intersections(py of the startingPoint is too big to intersect with the circle) + [TestCase("0,0", 1, "+10,+10", "+10,+100", null)] + [TestCase("0,0", 1, "+1,+10", "+1,+100", null)] + [TestCase("0,0", 1, "+0,+10", "+0,+100", null)] + [TestCase("0,0", 1, "-1,+10", "-1,+100", null)] + [TestCase("0,0", 1, "-10,+10", "-10,+100", null)] + public void CircleIntersectWithLineSegment2D(string sCenter, double radius, string sStart, string sEnd, string intersections) + { + var circle = new Circle2D(Point2D.Parse(sCenter), radius); + var segment = new LineSegment2D(Point2D.Parse(sStart), Point2D.Parse(sEnd)); + + var actual = circle.IntersectWith(segment); + + var expected = parseToPointsArray(intersections); + for (int i = 0; i < Math.Min(actual.Length, expected.Length); i++) + { + var a = actual[i]; + var e = expected[i]; + AssertGeometry.AreEqual(a, e, 1e-3); //FIXME! + } + } + + private Point2D[] parseToPointsArray(string input) + { + if (input == null) + { + return new Point2D[] { }; + } + + var result = input.Split(';') + .Select(s => Point2D.Parse(s)) + .ToArray(); + + return result; + } } } diff --git a/src/Spatial/Euclidean/Circle2D.cs b/src/Spatial/Euclidean/Circle2D.cs index 96edea3..35c9650 100644 --- a/src/Spatial/Euclidean/Circle2D.cs +++ b/src/Spatial/Euclidean/Circle2D.cs @@ -1,9 +1,11 @@ using MathNet.Spatial.Internals; using System; using System.Diagnostics.Contracts; +using System.Linq; using System.Xml; using System.Xml.Schema; using System.Xml.Serialization; +using MathNet.Numerics; using HashCode = MathNet.Spatial.Internals.HashCode; namespace MathNet.Spatial.Euclidean @@ -131,6 +133,67 @@ public static Circle2D FromPoints(Point2D pointA, Point2D pointB, Point2D pointC return new Circle2D(center, radius); } + /// + /// Returns the intersections of this circle with the given line. + /// + /// the given line + /// intersections as a Point2D Array, depending on the count. + public Point2D[] IntersectWith(Line2D line) + { + var ts = this.findParameterTs(line); + var result = ts + .Select(t => line.StartPoint + t * line.Direction) + .ToArray(); + return result; + } + + private double[] findParameterTs(Line2D line) + { + // These 2 equations in vector form can be described + // (p-cc)^2=r^2 (eq1) + // p=s+t*d (eq2) + // , where p is the point on the line and/or circle, + // cc is the center of the circle, + // r is the radius of the circle, + // s is the starting point of the line, + // t is the parameter and + // d is the line direction. + // Substituting (eq2) into (eq1) yields: + // ((s+t*d)-cc)^2=r^2 (eq3) + // (eq3) reduces to the following quadratic equation: a*t^2 + b*t + c==0 + + var cc = this.Center.ToVector2D(); //center of circle + var s = line.StartPoint.ToVector2D(); + var d = line.Direction; + var r = this.Radius; + + var a = 1.0; + var b = 2 * (s.DotProduct(d) - d.DotProduct(cc)); + var c = (s - cc).DotProduct(s - cc) - r * r; + + var solutions = FindRoots.Polynomial(new[] { c, b, a }); + var ts = solutions + .Where(z => z.IsReal()) + .Select(z => z.Real) + .Distinct() + .ToArray(); + return ts; + } + + /// + /// Returns the intersections of this circle with the given line segment, which lie within the segment. + /// + /// the given line-segment + /// intersections as a Point2D Array, depending on the count. + public Point2D[] IntersectWith(LineSegment2D segment) + { + var ts = findParameterTs(segment.ToLine2D()) + .Where(t => 0 <= t && t <= segment.Length); + var result = ts.Select(t => segment.StartPoint + t * segment.Direction).ToArray(); + return result; + } + + /// /// Returns a value to indicate if a pair of circles are equal /// diff --git a/src/Spatial/Euclidean/LineSegment2D.cs b/src/Spatial/Euclidean/LineSegment2D.cs index 2cc701c..22b4d4c 100644 --- a/src/Spatial/Euclidean/LineSegment2D.cs +++ b/src/Spatial/Euclidean/LineSegment2D.cs @@ -236,5 +236,9 @@ void IXmlSerializable.WriteXml(XmlWriter writer) writer.WriteElement("StartPoint", StartPoint); writer.WriteElement("EndPoint", EndPoint); } + + /// convert this to Line2D + /// converted Line2D object + public Line2D ToLine2D() => new Line2D(StartPoint, EndPoint); } }