2012-12-06 128 views
2

我想优化我的MCU的sin/cos以计算地理距离。公式的这一部分尤其是使用三角:如何改善正弦/余弦函数?

double e = (MyTan(lat2/2 + quarter_pi)/MyTan(lat1/2 + quarter_pi)); 

所以我试图建立自己的sin/cos查找表-PIPI如下:

#define PARTPERDEGREE 10 
double mysinlut[PARTPERDEGREE * 90 + 1]; 
double mycoslut[PARTPERDEGREE * 90 + 1]; 
void MySinCosCreate() 
{ 
    int i; 
    double angle, angleinc; 

    // Each degree also divided into 10 parts 
    angleinc = (M_PI/180)/PARTPERDEGREE; 
    for (i = 0, angle = 0.0; i <= (PARTPERDEGREE * 90 + 1); ++i, angle += angleinc) 
    { 
     mysinlut[i] = sin(angle); 
    } 

    angleinc = (M_PI/180)/PARTPERDEGREE; 
    for (i = 0, angle = 0.0; i <= (PARTPERDEGREE * 90 + 1); ++i, angle += angleinc) 
    { 
     mycoslut[i] = cos(angle); 
    } 
} 



double MySin(double rad) 
{ 
    int ix; 
    int sign = 1; 


    if(rad > (M_PI/2)) 
     rad = M_PI/2 - (rad - M_PI/2); 

    if(rad < -(M_PI/2)) 
     rad = -M_PI/2 - (rad + M_PI/2); 

    if(rad < 0) 
    { 
     sign = -1; 
     rad *= -1; 
    } 

    ix = (rad * 180)/M_PI * PARTPERDEGREE; 

    return sign * mysinlut[ix]; 
} 

double MyCos(double rad) 
{ 
    int ix; 
    int sign = 1; 


    if(rad > M_PI/2) 
    { 
     rad = M_PI/2 - (rad - M_PI/2); 
     sign = -1; 
    } 
    else if(rad < -(M_PI/2)) 
    { 
     rad = M_PI/2 + (rad + M_PI/2); 
     sign = -1; 
    } 
    else if(rad > -M_PI/2 && rad < M_PI/2) 
    { 
     rad = abs(rad); 
     sign = 1; 
    } 

    ix = (rad * 180)/M_PI * PARTPERDEGREE; 

    return sign * mycoslut[ix]; 
} 

double MyTan(double rad) 
{ 
    return MySin(rad)/MyCos(rad); 
} 

你可以看到表的分辨率为10份每度。我可以增加一点,但它没有多大帮助,看起来我需要一些插值。任何人都可以为我的功能提出一些实际改进,以获得更好的结以下是e的234个不同结果的图表。蓝色系列具有理想的正弦/余弦,红色来自LUT。

enter image description here

+0

你可以尝试使用导数sin(x + h)〜sin x + h * cos x'来近似,这应该给出更好的结果。 –

+0

@Daniel Fischer:你的意思是填充查找表吗?你能否用一些小的伪代码来回答,这会给我提示如何实现它? – Pablo

回答

4

看起来你的查找表太粗糙了。如果你不能使你的桌子更精细,使用衍生物来近似值应该会让你获得更好的结果。我们有

sin (x+h) ≈ sin x + h*cos x 
cos (x+h) ≈ cos x - h*sin x 

为小h。 (您可以使用更高的导数或通过使用表中的两个值来获得更好的近似值,在这两个值之间您的(计算的)角度位于其间,但这需要更长的时间,并且我收集速度是LUT的首要原因。)

所以你归一化的角度后,用

ix = (rad * 180)/M_PI * PARTPERDEGREE; 

使用

double h = rad - ix*angleinc; 
return sign*(mysinlut[ix] + h*mycoslut[ix]); 

RESP。

return sign*(mycoslut[ix] - h*mysinlut[ix]); 

这应该不会太慢,而且应该给予LUT点之间更好的近似值。

+0

不错,让我试试这个! – Pablo

+0

这使得LUT的结果与理想的结果相差甚远!我真的很感谢你的帮助! – Pablo

+0

实际上,当结果是切换标志时会有一些小问题...我会尽量简短地描述它。如果您能跟进并帮助我解决问题,我将不胜感激。 – Pablo

1

赎罪旧的,但良好的递推关系/ COS插补描述here

+0

这是我应该在LUT上应用的东西,当从表中加载值或者这可以完全替代LUT?我认为插值可以使用表格中的相邻值完成。 – Pablo