DHShapeMaker/Document.rtf at master · TheDwarfHorde/DHShapeMaker · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
{\rtf1\ansi\ansicpg1252\deff0\deflang1033{\fonttbl{\f0\fnil\fcharset0 MS Sans Serif;}}
{\*\generator Msftedit 5.41.21.2510;}\viewkind4\uc1\pard\ri-4320\sa200\sl276\slmult1\lang9\f0\fs22  private static void GetCurveControlPoints(IList<PointDouble> knots, \par
\pard\ri-4320\sa200\sl276\slmult1\tx1170\tab out PointDouble[] firstControlPoints, out PointDouble[] secondControlPoints)\par
\pard\ri-4320\sa200\sl276\slmult1         \{\par
            if (knots == null)\par
                throw new ArgumentNullException("knots");\par
            int n = knots.Count - 1;\par
            if (n < 1)\par
                throw new ArgumentException("At least two knot points required", "knots");\par
            if (n == 1)\par
            \{ // Special case: Bezier curve should be a straight line.\par
                firstControlPoints = new PointDouble[1];\par
                // 3P1 = 2P0 + P3\par
                firstControlPoints[0].X = (2 * knots[0].X + knots[1].X) / 3;\par
                firstControlPoints[0].Y = (2 * knots[0].Y + knots[1].Y) / 3;\par
 \par
                secondControlPoints = new PointDouble[1];\par
                // P2 = 2P1 \endash  P0\par
                secondControlPoints[0].X = 2 * firstControlPoints[0].X - knots[0].X;\par
                secondControlPoints[0].Y = 2 * firstControlPoints[0].Y - knots[0].Y;\par
                return;\par
            \}\par
 \par
            // Calculate first Bezier control points\par
            // Right hand side vector\par
            double[] rhs = new double[n];\par
 \par
            // Set right hand side X values\par
            for (int i = 1; i < n - 1; ++i)\par
                rhs[i] = 4 * knots[i].X + 2 * knots[i + 1].X;\par
            rhs[0] = knots[0].X + 2 * knots[1].X;\par
            rhs[n - 1] = (8 * knots[n - 1].X + knots[n].X) / 2.0;\par
            // Get first control points X-values\par
            double[] x = GetFirstControlPoints(rhs);\par
 \par
            // Set right hand side Y values\par
            for (int i = 1; i < n - 1; ++i)\par
                rhs[i] = 4 * knots[i].Y + 2 * knots[i + 1].Y;\par
            rhs[0] = knots[0].Y + 2 * knots[1].Y;\par
            rhs[n - 1] = (8 * knots[n - 1].Y + knots[n].Y) / 2.0;\par
            // Get first control points Y-values\par
            double[] y = GetFirstControlPoints(rhs);\par
 \par
            // Fill output arrays.\par
            firstControlPoints = new PointDouble[n];\par
            secondControlPoints = new PointDouble[n];\par
            for (int i = 0; i < n; ++i)\par
            \{\par
                // First control point\par
                firstControlPoints[i] = new PointDouble(x[i], y[i]);\par
                // Second control point\par
                if (i < n - 1)\par
                    secondControlPoints[i] = new PointDouble(2 * knots[i + 1].X - x[i + 1], 2 * knots[i + 1].Y - y[i + 1]);\par
                else\par
                    secondControlPoints[i] = new PointDouble((knots[n].X + x[n - 1]) / 2, (knots[n].Y + y[n - 1]) / 2);\par
            \}\par
        \}\par
 \par
        /// <summary>\par
        /// Solves a tridiagonal system for one of coordinates (x or y) of first Bezier control points.\par
        /// </summary>\par
        /// <param name="rhs">Right hand side vector.</param>\par
        /// <returns>Solution vector.</returns>\par
        private static double[] GetFirstControlPoints(double[] rhs)\par
        \{\par
            int n = rhs.Length;\par
            double[] x = new double[n]; // Solution vector.\par
\pard\sa200\sl276\slmult1\tx6210             double[] tmp = new double[n]; // Temp workspace.\par
\pard\ri-4320\sa200\sl276\slmult1  \par
            double b = 2.0;\par
            x[0] = rhs[0] / b;\par
            for (int i = 1; i < n; i++) // Decomposition and forward substitution.\par
            \{\par
                tmp[i] = 1 / b;\par
                b = (i < n - 1 ? 4.0 : 3.5) - tmp[i];\par
                x[i] = (rhs[i] - x[i - 1]) / b;\par
            \}\par
            for (int i = 1; i < n; i++)\par
                x[n - i - 1] -= tmp[n - i] * x[n - i]; // Backsubstitution.\par
 \par
            return x;\par
        \}\par
}