2013-04-03 112 views
3

Diagram查找相交的線和三次樣條之間的點

我需要一種方法以編程找到的線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)); 
    } 

參考文獻:

  1. http://paulbourke.net/miscellaneous/interpolation/
  2. http://www.math.vanderbilt.edu/~schectex/courses/cubic/

回答

1

代碼

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; 
} 

你可以計算()^(1/6)Math.Pow( , 1/6)或您的Math.Sqrt(CubicRoot())Math.Sqrt(Cbrt())。請參閱following thread on Microsoft forum

請注意q = 0(Atan(x/0) = pi/2 radians. Cos(Atan(x/0)/3) = Sqrt(3)/2)

+0

嗯......似乎並沒有工作。我使用均勻間隔的4個點進行測試,距離原點的距離相等,以定義g(x)。然後我在每個pi/8間隔檢查函數的結果。結果應該類似於正弦波,但它不是不穩定的,不連續的,並且總是正的。也許還有其他問題? – 2013-04-04 02:41:48

+0

另外,我可以肯定地說,問題出在這個函數的某個地方,而不是程序中的其他地方,因爲當我用正則線性交叉函數替換它時,結果如預期的那樣。 – 2013-04-04 02:54:03

+0

你可以給我發送第3階多項式來解決嗎?該函數解決了我測試的所有多項式。 – Edoot 2013-04-04 07:15:22