Merge pull request #284 from GSharker/dev/mibi/ruled-surface

Dev/mibi/ruled surface

src/GShark.Test.XUnit/Geometry/ArcTests.cs

@@ -65,8 +65,7 @@ public void Initializes_An_Arc_By_Two_Points_And_A_Direction()
 
             // Act
             Arc arc = Arc.ByStartEndDirection(pt1, pt2, dir);
-            var pt = arc.EndPoint;
-            var pt5 = arc.StartPoint;
+
             // Assert
             arc.StartPoint.EpsilonEquals(pt1, 1e-6).Should().BeTrue();
             arc.EndPoint.EpsilonEquals(pt2, 1e-6).Should().BeTrue();
@@ -124,7 +123,7 @@ public void It_Is_A_Curve_Representation_Of_The_Arc_From_0_To_90_Deg()
             };
 
             // Act
-            NurbsCurve arc = new Arc(Plane.PlaneYZ, 20, new Interval(0.0, 1.8)).ToNurbsCurve();
+            NurbsCurve arc = new Arc(Plane.PlaneYZ, 20, new Interval(0.0, 1.8)).ToNurbs();
 
             // Assert
             arc.ControlPointLocations.Count.Should().Be(5);
@@ -163,7 +162,7 @@ public void It_Is_A_Curve_Representation_Of_ExampleArc3D()
             };
 
             // Act
-            NurbsCurve arc = _exampleArc3D.ToNurbsCurve();
+            NurbsCurve arc = _exampleArc3D.ToNurbs();
 
             // Assert
             arc.ControlPointLocations.Count.Should().Be(7);

src/GShark.Test.XUnit/Geometry/NurbsSurfaceTests.cs

@@ -374,5 +374,79 @@ public void Returns_True_If_Surface_Is_Close()
             // Assert
             surface.IsClosed(SurfaceDirection.V).Should().BeTrue();
         }
+
+        [Theory]
+        [InlineData(0.1, 0.1, new double[] {0.2655, 1, 2.442})]
+        [InlineData(0.5, 0.5, new double[] { 4.0625, 5, 4.0625 })]
+        [InlineData(1.0, 1.0, new double[] {10, 10, 0 })]
+        public void Returns_A_Ruled_Surface_Between_Two_Nurbs_Curve(double u, double v, double[] pt1)
+        {
+            // Arrange
+            Point3 expectedPt = new Point3(pt1[0], pt1[1], pt1[2]);
+            List<Point3> ptsA = new List<Point3>
+            {
+                new Point3(0, 0, 0),
+                new Point3(0, 0, 5),
+                new Point3(5, 0, 5),
+                new Point3(5, 0, 0),
+                new Point3(10, 0, 0)
+            };
+
+            List<Point3> ptsB = new List<Point3>
+            {
+                new Point3(0, 10, 0),
+                new Point3(0, 10, 5),
+                new Point3(5, 10, 5),
+                new Point3(5, 10, 0),
+                new Point3(10, 10, 0)
+            };
+
+            NurbsCurve curveA = new NurbsCurve(ptsA, 3);
+            NurbsCurve curveB = new NurbsCurve(ptsB, 2);
+
+            // Act
+            NurbsSurface ruledSurface = NurbsSurface.CreateRuledSurface(curveA, curveB);
+            Point3 pointAt = ruledSurface.PointAt(u, v);
+
+            // Assert
+            pointAt.EpsilonEquals(expectedPt, GSharkMath.MinTolerance).Should().BeTrue();
+        }
+
+        [Theory]
+        [InlineData(0.1, 0.1, new double[] { 0.0225, 1, 2.055 })]
+        [InlineData(0.5, 0.5, new double[] { 4.6875, 5, 4.6875 })]
+        [InlineData(1.0, 1.0, new double[] { 10, 10, 0 })]
+        public void Returns_A_Ruled_Surface_Between_A_Polyline_And_A_Nurbs_Curve(double u, double v, double[] pt)
+        {
+            // Arrange
+            Point3 expectedPt = new Point3(pt[0], pt[1], pt[2]);
+            List<Point3> ptsA = new List<Point3>
+            {
+                new Point3(0, 0, 0),
+                new Point3(0, 0, 5),
+                new Point3(5, 0, 5),
+                new Point3(5, 0, 0),
+                new Point3(10, 0, 0)
+            };
+
+            List<Point3> ptsB = new List<Point3>
+            {
+                new Point3(0, 10, 0),
+                new Point3(0, 10, 5),
+                new Point3(5, 10, 5),
+                new Point3(5, 10, 0),
+                new Point3(10, 10, 0)
+            };
+
+            Polyline poly = new Polyline(ptsA);
+            NurbsCurve curveB = new NurbsCurve(ptsB, 2);
+
+            // Act
+            NurbsSurface ruledSurface = NurbsSurface.CreateRuledSurface(poly.ToNurbs(), curveB);
+            Point3 pointAt = ruledSurface.PointAt(u, v);
+
+            // Assert
+            pointAt.EpsilonEquals(expectedPt, GSharkMath.MinTolerance).Should().BeTrue();
+        }
     }
 }

src/GShark.Test.XUnit/Geometry/PolygonTests.cs

@@ -146,7 +146,6 @@ public void It_Returns_A_Polygon_Transformed_In_NurbsCurve()
             poly2D.Degree.Should().Be(1);
             for (int i = 1; i < poly2D.Knots.Count - 1; i++)
             {
-                Point3 pt = poly2D.PointAt(knots[i]);
                 Planar2D[i - 1].Equals(poly2D.PointAt(knots[i])).Should().BeTrue();
             }
         }

src/GShark.Test.XUnit/Geometry/PolylineTests.cs

@@ -267,7 +267,6 @@ public void It_Returns_A_Polyline_Transformed_In_NurbsCurve()
             polyNurbs.Degree.Should().Be(1);
             for (int i = 1; i < polyNurbs.Knots.Count - 1; i++)
             {
-                Point3 pt = polyNurbs.PointAt(knots[i]);
                 pts[i - 1].EpsilonEquals(polyNurbs.PointAt(knots[i]), GSharkMath.MaxTolerance).Should().BeTrue();
             }
         }

src/GShark.Test.XUnit/Operation/ModifyTests.cs

@@ -271,7 +271,7 @@ public void Returns_A_Curve_Joining_Different_Types_Of_Curves()
             Point3 expectedPt1 = new Point3(3.34125, 0.005, 0.6125);
             Point3 expectedPt2 = new Point3(5, 2.363281, 4.570313);
             Point3 expectedPt3 = new Point3(5.351058, 5, -4.340474);
-            ICurve[] curves = { curve, ln.ToNurbs(), arc.ToNurbsCurve() };
+            ICurve[] curves = { curve, ln.ToNurbs(), arc.ToNurbs() };
 
             // Act
             ICurve joinedCurve = Modify.JoinCurves(curves);

src/GShark.Test.XUnit/Operation/OffsetTests.cs

@@ -95,7 +95,7 @@ public void Returns_The_Offset_Of_A_Curve()
             {
                 Point3 pt = offsetResult.PointAt(i);
                 Point3 closestPt = crv.ClosestPoint(pt);
-                pt.DistanceTo(closestPt).Should().BeApproximately(offset, GSharkMath.MinTolerance);
+                pt.DistanceTo(closestPt).Should().BeApproximately(offset, GSharkMath.MaxTolerance);
             }
         }
     }

src/GShark/ExtendedMethods/Curve.cs

@@ -123,10 +123,7 @@ public static List<Plane> PerpendicularFrames(this ICurve curve, List<double> uV
                 var sNext = Vector3.CrossProduct(pointsOnCurveTan[i + 1], rNext); //compute vector s[i+1] of next frame
 
                 //create output frame
-                var frameNext = new Plane();
-                frameNext.Origin = pointsOnCurve[i + 1];
-                frameNext.XAxis = rNext;
-                frameNext.YAxis = sNext;
+                var frameNext = new Plane {Origin = pointsOnCurve[i + 1], XAxis = rNext, YAxis = sNext};
                 perpFrames[i + 1] = frameNext; //output frame
             }
 

src/GShark/Geometry/Arc.cs

@@ -170,43 +170,38 @@ public static Arc ByStartEndDirection(Point3 ptStart, Point3 ptEnd, Vector3 dir)
         /// Implementation of Algorithm A7.1 from The NURBS Book by Piegl and Tiller.
         /// </summary>
         /// <returns>A nurbs curve shaped like this arc.</returns>
-        public NurbsCurve ToNurbsCurve()
+        public new NurbsCurve ToNurbs()
         {
             Vector3 axisX = Plane.XAxis;
             Vector3 axisY = Plane.YAxis;
             int numberOfArc;
             Point4[] ctrPts;
-            double[] weights;
 
             // Number of arcs.
             double piNum = 0.5 * Math.PI;
             if ((Angle - piNum) <= GSharkMath.Epsilon)
             {
                 numberOfArc = 1;
                 ctrPts = new Point4[3];
-                weights = new double[3];
             }
             else if ((Angle - piNum * 2) <= GSharkMath.Epsilon)
             {
                 numberOfArc = 2;
                 ctrPts = new Point4[5];
-                weights = new double[5];
             }
             else if ((Angle - piNum * 3) <= GSharkMath.Epsilon)
             {
                 numberOfArc = 3;
                 ctrPts = new Point4[7];
-                weights = new double[7];
             }
             else
             {
                 numberOfArc = 4;
                 ctrPts = new Point4[9];
-                weights = new double[9];
             }
 
             double detTheta = Angle / numberOfArc;
-            double weight1 = Math.Cos(detTheta / 2);
+            double weight = Math.Cos(detTheta / 2);
             Point3 p0 = Center + (axisX * (Radius * Math.Cos(Domain.T0)) + axisY * (Radius * Math.Sin(Domain.T0)));
             Vector3 t0 = axisY * Math.Cos(Domain.T0) - axisX * Math.Sin(Domain.T0);
 
@@ -229,7 +224,7 @@ public NurbsCurve ToNurbsCurve()
                 Intersect.LineLine(ln0, ln1, out _, out _, out double u0, out _);
                 Point3 p1 = p0 + (t0 * u0);
 
-                ctrPts[index + 1] = new Point4(p1, weight1);
+                ctrPts[index + 1] = new Point4(p1, weight);
                 index += 2;
 
                 if (i >= numberOfArc)

src/GShark/Geometry/NurbsSurface.cs

@@ -6,6 +6,7 @@
 using System.Collections.Generic;
 using System.Linq;
 using System.Text;
+using GShark.ExtendedMethods;
 
 namespace GShark.Geometry
 {
@@ -214,6 +215,63 @@ public static NurbsSurface CreateLoftedSurface(IList<NurbsCurve> curves, LoftTyp
             return new NurbsSurface(degreeU, degreeV, knotVectorU, knotVectorV, surfaceControlPoints);
         }
 
+        /// <summary>
+        /// Constructs a ruled surface between two curves.
+        /// <em>Follows the algorithm at page 337 of The NURBS Book by Piegl and Tiller.</em>
+        /// </summary>
+        /// <param name="curveA">The first curve.</param>
+        /// <param name="curveB">The second curve.</param>
+        /// <returns>A ruled surface.</returns>
+        public static NurbsSurface CreateRuledSurface(NurbsCurve curveA, NurbsCurve curveB)
+        {
+            // ToDo: this should be removed rebuilding the curve with less point. 
+            if (curveA.ControlPointLocations.Count != curveB.ControlPointLocations.Count)
+            {
+                throw new ArgumentException("Ruled surface only works if the curves have the same number of points.");
+            }
+
+            ICurve copyCurveA = curveA;
+            ICurve copyCurveB = curveB;
+
+            // Ensure that the two curves are defined on the same parameter range
+            if (Math.Abs(copyCurveA.Knots.Domain - copyCurveB.Knots.Domain) > GSharkMath.Epsilon)
+            {
+                Interval knotsIntervalB = new Interval(copyCurveB.Knots.First(), curveB.Knots.Last());
+                Interval knotsIntervalA = new Interval(copyCurveA.Knots.First(), copyCurveA.Knots.Last());
+                var knots = curveB.Knots.Select(k => GSharkMath.RemapValue(k, knotsIntervalB, knotsIntervalA)).ToKnot();
+                copyCurveB = new NurbsCurve(copyCurveB.Degree, knots, copyCurveB.ControlPoints);
+            }
+
+            // Raise the degree if the lower degree curve.
+            if (copyCurveA.Degree < copyCurveB.Degree)
+            {
+                copyCurveA = Modify.ElevateDegree(copyCurveA, copyCurveB.Degree);
+            }
+
+            if (copyCurveA.Degree > copyCurveB.Degree)
+            {
+                copyCurveB = Modify.ElevateDegree(copyCurveB, copyCurveA.Degree);
+            }
+
+            // If the knot vectors U1 and U2 are not identical, merge them to obtain the knot vector U.
+            if (!copyCurveA.Knots.SequenceEqual(copyCurveB.Knots))
+            {
+                KnotVector combinedKnots = copyCurveA.Knots
+                    .Concat(copyCurveB.Knots.Where(k => !copyCurveA.Knots.Contains(k)))
+                    .OrderBy(k => k).ToKnot();
+
+                KnotVector knotToInsertA = combinedKnots.Where(k => !copyCurveA.Knots.Contains(k)).ToKnot();
+                KnotVector knotToInsertB = combinedKnots.Where(k => !copyCurveB.Knots.Contains(k)).ToKnot();
+
+                // using U, apply knot refinement to both curves.
+                copyCurveA = Modify.CurveKnotRefine(copyCurveA, knotToInsertA);
+                copyCurveB = Modify.CurveKnotRefine(copyCurveB, knotToInsertB);
+            }
+
+            return new NurbsSurface(1, copyCurveA.Degree, new KnotVector(1, 2), copyCurveA.Knots,
+                new List<List<Point4>> {copyCurveA.ControlPoints, copyCurveB.ControlPoints});
+        }
+
         /// <summary>
         /// Evaluates a point at a given U and V parameters.
         /// </summary>

src/GShark/Geometry/Vector3.cs

@@ -486,7 +486,7 @@ public double this[int i]
         /// <returns>true if obj is a Vector3d and has the same coordinates as this; otherwise false.</returns>
         public override bool Equals(object obj)
         {
-            return (obj is Vector3 && this == (Vector3)obj);
+            return (obj is Vector3 vector3 && this == vector3);
         }
 
         /// <summary>
@@ -556,9 +556,9 @@ public int CompareTo(Vector3 other)
 
         int IComparable.CompareTo(object obj)
         {
-            if (obj is Vector3)
+            if (obj is Vector3 vector3)
             {
-                return CompareTo((Vector3)obj);
+                return CompareTo(vector3);
             }
 
             throw new ArgumentException("Input must be of type Vector3d", "obj");

src/GShark/Operation/Intersect.cs

@@ -290,7 +290,7 @@ public static bool PlaneCircle(Plane pl, Circle cl, out Point3[] pts)
                 return false;
             }
 
-            bool intersection = PlanePlane(pl, cl.Plane, out Line intersectionLine);
+            _ = PlanePlane(pl, cl.Plane, out Line intersectionLine);
             Point3 closestPt = intersectionLine.ClosestPoint(clPt);
             double distance = clPt.DistanceTo(intersectionLine);
 

src/GShark/Operation/Modify.cs

@@ -444,7 +444,7 @@ private static double BezDegreeReduce(Point4[] bpts, int degree, out Point4[] rb
             P[0] = bpts[0];
             P[P.Length - 1] = bpts.Last();
 
-            bool isDegreeOdd = (degree % 2 != 0) ? true : false;
+            bool isDegreeOdd = degree % 2 != 0;
 
             int r1 = (isDegreeOdd) ? r - 1 : r;
 
@@ -497,7 +497,7 @@ the maximums of geometric and parametric errors are not necessarily at the same
         /// <param name="curve">The object curve to elevate.</param>
         /// <param name="tolerance">Tolerance value declaring if the curve is degree reducible. Default value set to 10e-4, refer to Eq 5.30 for the meaning.</param>
         /// <returns>The curve after degree reduction, the curve will be degree - 1 from the input.</returns>
-        public static ICurve ReduceDegree(ICurve curve, double tolerance = 10e-4)
+        public static NurbsCurve ReduceDegree(NurbsCurve curve, double tolerance = 10e-4)
         {
             int n = curve.Knots.Count - curve.Degree - 2;
             int p = curve.Degree;
@@ -614,7 +614,7 @@ public static ICurve ReduceDegree(ICurve curve, double tolerance = 10e-4)
                         }
 
                         // Compute knot removal error bounds (Br).
-                        double Br = 0.0;
+                        double Br;
                         if (j - l < k)
                         {
                             Br = Pw[l - 2].DistanceTo(rbpts[kj + 1]);

src/GShark/Operation/Offset.cs

@@ -59,14 +59,14 @@ public static ICurve Curve(ICurve crv, double distance, Plane pln)
                 return crv;
             }
 
-            (List<double> tValues, List<Point3> pts) subdivision = Tessellation.CurveAdaptiveSample(crv);
+            var (tValues, pts) = Tessellation.CurveAdaptiveSample(crv);
 
             List<Point3> offsetPts = new List<Point3>();
-            for (int i = 0; i < subdivision.pts.Count; i++)
+            for (int i = 0; i < pts.Count; i++)
             {
-                Vector3 tangent = Evaluation.RationalCurveTangent(crv, subdivision.tValues[i]);
+                Vector3 tangent = Evaluation.RationalCurveTangent(crv, tValues[i]);
                 Vector3 vecOffset = Vector3.CrossProduct(tangent, pln.ZAxis).Amplify(distance);
-                offsetPts.Add(subdivision.pts[i] + vecOffset);
+                offsetPts.Add(pts[i] + vecOffset);
             }
 
             return Fitting.InterpolatedCurve(offsetPts, 2);