Circle2D intersection with Line(Segment)2D

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<ArgumentException>(() => { 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;
+        }
     }
 }

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);            
         }
 
+        /// <summary>
+        /// Returns the intersections of this circle with the given line.
+        /// </summary>
+        /// <param name="line">the given line</param>
+        /// <returns>intersections as a Point2D Array, depending on the count.</returns>
+        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;
+        }
+
+        /// <summary>
+        /// Returns the intersections of this circle with the given line segment, which lie within the segment.
+        /// </summary>
+        /// <param name="segment">the given line-segment</param>
+        /// <returns>intersections as a Point2D Array, depending on the count.</returns>
+        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;
+        }
+
+
         /// <summary>
         /// Returns a value to indicate if a pair of circles are equal
         /// </summary>

src/Spatial/Euclidean/LineSegment2D.cs

@@ -236,5 +236,9 @@ void IXmlSerializable.WriteXml(XmlWriter writer)
             writer.WriteElement("StartPoint", StartPoint);
             writer.WriteElement("EndPoint", EndPoint);
         }
+
+        /// <summary>convert this to Line2D </summary>
+        /// <returns>converted Line2D object</returns>
+        public Line2D ToLine2D() => new Line2D(StartPoint, EndPoint);
     }
 }