寻找一条直线和三次样条曲线的交点

Question

寻找一条直线和三次样条曲线的交点

3

我需要一种程序化的方法来查找一条从原点开始的直线f(x)与由4个点定义的三次样条曲线相交的点。该直线保证与样条曲线的中心线段在X1和X2之间相交。

我尝试了许多方法，但似乎无法获得预期结果。我怀疑问题出在我对复数的处理上。

有人能找到我的代码问题，或者建议其他方法吗？

    private Vector2 CubicInterpolatedIntersection(float y0, float y1, 
                            float y2, float y3, float lineSlope, lineYOffset)
    {
        // f(x) = lineSlope * x + lineYOffset
        // g(x) = (a0 * x^3) + (a1 * x^2) + (a2 * x) + a3

        // Calculate Catmull-Rom coefficients for g(x) equation as found 
        // in reference (1). These 
        double a0, a1, a2, a3;
        a0 = -0.5 * y0 + 1.5 * y1 - 1.5 * y2 + 0.5 * y3;
        a1 = y0 - 2.5 * y1 + 2 * y2 - 0.5 * y3;
        a2 = -0.5 * y0 + 0.5 * y2;
        a3 = y1;

        // To find POI, let 'g(x) - f(x) = 0'. Simplified equation is:
        // (a0 * x^3) + (a1 * x^2) + ((a2 - lineSlope) * x) 
        //                         + (a3 - lineYOffset) = 0

        // Put in standard form 'ax^3 + bx^2 + cx + d = 0'
        double a, b, c, d;
        a = a0;
        b = a1;
        c = a2 - lineSlope;
        d = a3 - lineYOffset;


        // Solve for roots using cubic equation as found in reference (2).
        //  x = {q + [q^2 + (r-p^2)^3]^(1/2)}^(1/3) 
        //         + {q - [q^2 + (r-p^2)^3]^(1/2)}^(1/3) + p
        //  Where...
        double p, q, r;
        p = -b / (3 * a);
        q = p * p * p + (b * c - 3 * a * d) / (6 * a * a);
        r = c / (3 * a);

        //Solve the equation starting from the center.
        double x, x2;
        x = r - p * p;
        x = x * x * x + q * q;
        // Here's where I'm not sure. The cubic equation contains a square 
        // root. So if x is negative at this point, then we need to proceed
        // with complex numbers.
        if (x >= 0)
        {
            x = Math.Sqrt(x);
            x = CubicRoot(q + x) + CubicRoot(q - x) + p;
        }
        else
        {
            x = Math.Sqrt(Math.Abs(x)); 
            // x now represents the imaginary component of 
            //                       a complex number 'a + b*i'
            // We need to take the cubic root of 'q + x' and 'q - x'
            // Since x is imaginary, we have two complex numbers in 
            // standard form. Therefore, we take the cubic root of 
            // the magnitude of these complex numbers
            x = CubicRoot(Math.Sqrt(q * q + x * x)) + 
                 CubicRoot(Math.Sqrt(q * q + -x * -x)) + p;
        }

        // At this point, x should hold the x-value of 
        //      the point of intersection.
        // Now use it to solve g(x).

        x2 = x * x;
        return new Vector2((float)Math.Abs(x), 
           (float)Math.Abs(a0 * x * x2 + a1 * x2 + a2 * x + a3));
    }

References:

- Madison Brown

2个回答

0

对于一个二次方程，至少存在一个实数根。使用这种方法来找到根 http://en.wikipedia.org/wiki/Cubic_function#General_formula_of_roots

- Vineeth Chitteti

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Edoot · Accepted Answer

代码

if (x >= 0)
{  // One real root and two imaginaries.
    x = Math.Sqrt(x);
    x = CubicRoot(q + x) + CubicRoot(q - x) + p;
}
else
{  // Three real roots.
    x = Math.Sqrt(Math.Abs(x)); 


    x_1 = Math.Sign(q)*2*(q*q + x*x)^(1/6)*Math.Cos(Math.Atan(x/q)/3) + p;
    x_2 = Math.Sign(q)*2*(q*q + x*x)^(1/6)*Math.Cos(Math.Atan(x/q)/3 + Math.PI*2/3) + p;
    x_3 = Math.Sign(q)*2*(q*q + x*x)^(1/6)*Math.Cos(Math.Atan(x/q)/3 + Math.PI*4/3) + p;
}

你可以使用 Math.Pow( , 1/6) 或者 Math.Sqrt(CubicRoot( )) 或者 Math.Sqrt(Cbrt( )) 来计算 ( )^(1/6)。请参考Microsoft论坛上的讨论。

当 q = 0 时要小心。 ( Atan(x/0) = pi/2 弧度。 Cos(Atan(x/0)/3) = Sqrt(3)/2 )