在不使用标准函数或C99的情况下计算atan2

Question

在不使用标准函数或C99的情况下计算atan2

7

我正在使用一款三轴加速度计计算角度，但我的编译器没有atan或atan2函数。它有一个保留的内存插槽，但它调用了一个我在任何文件中都找不到的函数。

我的编译器是Keil µVision 4，运行ARMCC编译器。编译器有math.h文件，但该函数是extern类型，不存在：

  extern _ARMABI double atan2(double /*y*/, double /*x*/);

有没有可以包含反正切函数的库或函数？又或者，有没有其他函数可以用于从加速度计计算角度？我需要完整的三轴校准角度。

编辑：我希望避免使用预先计算的值的表格。

- Jeffa

我无法验证这一点，但如果您正在使用带有uVision的GNU ARM编译器，则IDE中应该有一个选项可用于使用描述此处的数学库。您还可能需要使用标准库而不是排除某些浮点函数的微型库。 - tinman

这个函数在 ARMCC 编译器自带的库中肯定存在，并且即使使用 microlib，它也应该能够正常工作，因为数学函数是相同的（区别在于低级软件浮点函数）。 - Al Grant

4个回答

9

实现自己的arctan2并不是很困难。使用此公式将arctan2转换为arctan，然后可以使用这个无限级数计算arctan。如果您对该无限级数求和足够多的项，就会非常接近库函数arctan2所做的事情。

这里有一个类似于exp()的实现，您可以将其用作参考。

- Pavan Manjunath

5

无穷级数对于探索函数的数学属性非常有用，但不适合计算。它们收敛缓慢，需要许多项才能得到准确结果，并且在使用浮点运算时难以计算而不积累误差。通常，数学库通过使用最小化最大误差的最小二乘多项式来计算这些函数，这样可以只用少量项就可以获得精确度。 - Eric Postpischil

2

此外，函数被划分为不同的区间，有各种各样的原因。arctan 应该被划分为 |x| < 1 和 |x| > 1（|x| = 1 可以在任一分区中）。对于 |x| < 1，可以使用最小最大多项式来近似计算 arctan(x)。对于 |x| > 1，可以使用以 1/x 为输入的最小最大多项式来近似计算 arctan(x)。 - Eric Postpischil

那将是一个高效的实现吗？GCC库中使用的算法是什么？ - captain

@Jay 这不是高效的实现方式，但非常容易实现。我不确定 libc 库，但请查看 Eric 的答案，其中一个开源实现广泛使用预定义常量来加速计算。 - Pavan Manjunath

@PavanManjunath 不使用泰勒级数的原因是误差呈多项式增长。更好的近似方法是使用Remez交换算法来找到一个多项式（同阶）在所需区间内给出有界绝对误差，但在区间外没有限制。 - rwong

4

这里有一个开源的atan实现，链接地址为这里.

- Eric Postpischil

3

苹果有一个开源门户网站？我从未想过。+1 - andand

0

数学函数的实际实现（或者如果存在硬件浮点单元的存根）应该在libm中。使用GCC，可以通过向编译器传递-lm来指示这一点，但我不知道您特定工具的操作方法。

- Ignacio Vazquez-Abrams

lm选项不受Keil µVision支持。"未识别的选项"。

- Jeffa

2

据我所知，Keil是一种完全不同的编译器，与GCC没有任何关系。 - Lundin

使用GCC和类似的编译器，“-l”不表示链接器选项，而是表示库。 “-l <foo>”表示名为“lib<foo>.<suffix>” 的库，其中包括适当的后缀，例如“.a”或“.dylib”，具体取决于编译器版本和目标平台。当然，命令行上列出的库必须传递给链接器。用于传递链接器选项的GCC开关是“-X”。 - Eric Postpischil

@EricPostpischil。是的，您非常具体和正确！ :) 感谢您的纠正。我删除了我的含糊评论。 - Pavan Manjunath

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

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

- rcor · Accepted Answer

以下代码使用有理逼近来获取在[0 1)区间内归一化的反正切值（您可以将结果乘以Pi / 2以获得实际的反正切值）。

归一化反正切(x) ~ (b x + x^2) / (1 + 2 b x + x^2)

其中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 normalized_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 normalized_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;
}

如果您需要更高的精度，可以使用第三阶有理函数：

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

其中 c = (1 + sqrt(17)) / 8

这个函数的最大逼近误差为0.00811度。