2011-09-11 160 views
9

我要實現環境ASIN,ACOS和ATAN這裏我只以下的數學工具:逼近反三角函數

  • 正弦
  • 餘弦
  • 基本定點運算(浮點數不可用)

我也已經有相當不錯的平方根功能。

我可以使用它們來實現合理高效的反三角函數嗎?

我不需要太高的精度(反正浮點數的精度非常有限),基本的近似值可以。

我已經半決定去查表,但我想知道是否有一些整潔的選項(這不需要幾百行代碼來實現基本的數學)。

編輯:

要澄清一些事情:我需要以每秒35幀運行函數每秒數百幀時間。

+0

[三角函數如何工作?](http://stackoverflow.com/questions/345085/how-do-trigonometric-functions-work) –

+0

建議的副本更多關於三角函數如何工作(就像這是關於反三角函數 – Teepeemm

回答

2

你需要一個大的精度arcsin(x)功能?如果否,則可以在N個節點中計算arcsin,並將值保存在內存中。我建議使用行aproximation。如果x = A*x_(N) + (1-A)*x_(N+1)x = A*arcsin(x_(N)) + (1-A)*arcsin(x_(N+1))其中arcsin(x_(N))是已知的。

+0

是的,這是我在OP中討論的表格查詢。沒有看到爲什麼我會在運行時計算出來的原因,我會simoly將值烘焙到程序中,所以實際的計算只是兩個值之間的插值。 –

2

您可能想要使用近似值:使用infinite series,直到解決方案足夠接近您爲止。

例如:
arcsin(z) = Sigma((2n!)/((2^2n)*(n!)^2)*((z^(2n+1))/(2n+1)))其中n在[0,無窮大)

1

也許一些種像牛頓拉普遜智能蠻力。

所以解決ASIN()你最速下降下去了罪()

+0

,你可以從一個小查找表中選擇起點來加速計算 –

8

在固定點e環境(S15.16)我成功地使用了CORDIC算法(參見維基百科的一般描述)來計算atan2(y,x),然後使用知名的函數標識從asin()和acos根:

asin(x) = atan2 (x, sqrt ((1.0 + x) * (1.0 - x))) 
acos(x) = atan2 (sqrt ((1.0 + x) * (1.0 - x)), x) 

事實證明,找到一個有用的關於atan2()的迭代的CORDIC迭代比我想象的更難。下面的網站似乎包含足夠詳細的說明,同時還討論了兩種替代方法,多項式逼近與查找表:

http://ch.mathworks.com/examples/matlab-fixed-point-designer/615-calculate-fixed-point-arctangent

+0

從wikip edia,CORDIC甚至不使用trig函數(整潔!);我想象你所做的是搜索;給定sin(),cos()似乎牛頓 - 拉夫森或者其他一些會更好? (需要更少的迭代,儘管迭代的代價會有所不同。) – petrelharp

+0

我建議查看CORDIC的原因是因爲它只需要定點算術。 CORDIC最常見的用途可能是實施sin/cos,這是我第一次瞭解它(1987年)。但是也可以用CORDIC來計算其他很多功能,例如atan2。由於我沒有任何代碼可以用CORDIC來計算atan2,所以我試圖找到一個足夠詳細的網站,以便有人可以在其上實現一個實現。我上面發佈的鏈接是我可以在幾分鐘內通過搜索引擎找到的最佳網頁。 – njuffa

1

將以下代碼添加到固定點應該很容易。它使用rational approximation來計算歸一化到[0 1)區間的反正切(你可以乘以Pi/2得到真正的反正切)。然後,您可以使用well known identities從反正切中獲取弧線/弧線。

normalized_atan(x) ~ (b x + x^2)/(1 + 2 b x + x^2) 

where b = 0.596227 

的最大誤差是0.1620º

#include <stdint.h> 
#include <math.h> 

// Approximates atan(x) normalized to the [-1,1] range 
// with a maximum error of 0.1620 degrees. 

float norm_atan(float x) 
{ 
    static const uint32_t sign_mask = 0x80000000; 
    static const float b = 0.596227f; 

    // Extract the sign bit 
    uint32_t ux_s = sign_mask & (uint32_t &)x; 

    // Calculate the arctangent in the first quadrant 
    float bx_a = ::fabs(b * x); 
    float num = bx_a + x * x; 
    float atan_1q = num/(1.f + bx_a + num); 

    // Restore the sign bit 
    uint32_t atan_2q = ux_s | (uint32_t &)atan_1q; 
    return (float &)atan_2q; 
} 

// Approximates atan2(y, x) normalized to the [0,4) range 
// with a maximum error of 0.1620 degrees 

float norm_atan2(float y, float x) 
{ 
    static const uint32_t sign_mask = 0x80000000; 
    static const float b = 0.596227f; 

    // Extract the sign bits 
    uint32_t ux_s = sign_mask & (uint32_t &)x; 
    uint32_t uy_s = sign_mask & (uint32_t &)y; 

    // Determine the quadrant offset 
    float q = (float)((~ux_s & uy_s) >> 29 | ux_s >> 30); 

    // Calculate the arctangent in the first quadrant 
    float bxy_a = ::fabs(b * x * y); 
    float num = bxy_a + y * y; 
    float atan_1q = num/(x * x + bxy_a + num); 

    // Translate it to the proper quadrant 
    uint32_t uatan_2q = (ux_s^uy_s) | (uint32_t &)atan_1q; 
    return q + (float &)uatan_2q; 
} 

如果你需要更精確,有一個3階有理函數:

normalized_atan(x) ~ (c x + x^2 + x^3)/(1 + (c + 1) x + (c + 1) x^2 + x^3) 

where c = (1 + sqrt(17))/8 

其中有0.00811最大逼近誤差º

1

在此提交我的回答other similar question.

NVIDIA已經和我已經用我自己的用途,舉幾個例子一些重要的資源:acosasinatan2等等...

這些算法產生足夠精確的結果。下面是他們的代碼拷貝直線上升Python的例子在粘貼:

import math 
def nVidia_acos(x): 
    negate = float(x<0) 
    x=abs(x) 
    ret = -0.0187293 
    ret = ret * x 
    ret = ret + 0.0742610 
    ret = ret * x 
    ret = ret - 0.2121144 
    ret = ret * x 
    ret = ret + 1.5707288 
    ret = ret * math.sqrt(1.0-x) 
    ret = ret - 2 * negate * ret 
    return negate * 3.14159265358979 + ret 

而且這裏的結果進行比較:

nVidia_acos(0.5) result: 1.0471513828611643 
math.acos(0.5) result: 1.0471975511965976 

那是相當接近!乘以57.29577951即可獲得度數的結果,這也是來自他們的「度數」公式。

-1

只有連續函數可以用多項式近似。並且arcsin(x)在點x = 1.same arccos(x)中是不連續的。但是在這種情況下,範圍減少到區間1,sqrt(1/2)可以避免這種情況。我們有arcsin(x)= pi/2-arccos(x),arccos(x)= pi/2-arcsin(x)。你可以使用matlab進行最小最大逼近。僅在範圍[0,sqrt )](如果arcsin的請求角度大於sqrt(1/2)的值,則找到cos(x)。僅適用於x的反正切函數1.arctan(x)= pi/2-arctan(1/x)。