[Math]Audio开发中定点的sin函数实现和评测

来源：互联网发布：c语言不合法标识符编辑：程序博客网时间：2024/06/11 08:16

Audio开发中定点的sin函数实现和评测

-- By Water

在音频开发中经常会使用到sin这个函数，例如产生一个1KHz的纯音就会用到sin(2*pi*fc*n/fs)。但是在嵌入式系统或者linux kernel开发或者DSP开发中经常会有Fixed point开发需求，这样就不呢使用系统浮点函数库，就得自己来实现fixed point sin函数了。

最简单直白的方式就是创建一个sin table在使用的时候去查表获得相应的sin值。

产生sin表格的方式参考如下：

    // Step1. 产生需要的Sin Table    {        double phase = 0;        double delta = 2*3.1415926535897932384626433832795/sin_tab_size;        for (i=0; i<sin_tab_size; i++)        {            psintab[i] = (long)((0x7FFFFFFF * sin(phase)) + 0.5);            phase += delta;        }    }

然后在使用的时候就可以直接查表获得sin值：

long lib_sin_lite(unsigned long phase){    return tab_sine_1024_32bits_Q31[phase & 0x3FF];}

这样做的好处是效率比较高，但精度较低。

在Table Size为1024的时候，最大的误差大约在-44dB左右。

如果要获得更高的精度就需要增大Table Size，每增大一倍就能减少6dB的误差。

在增大到4096的时候最大误差减少到-56dB左右。

而且在Table Size比较小的时候每两个点之间的Phase差会比较大。

例如：1024个点的时候，每两个点之间的phase差是0.006。

但在48KHz 采样率下要达到1Hz的精度就需要Phase差小于2*pi/48000 = 0.000131。

如果是Hi-Res audio sampling rate达到384000Hz的话要分辨1Hz的audio就需要Phase差小于0.0000164。

这样1024个点的Sin Table就远远不能满足要求了，就需要524288个点这么大的Sin Table了。

虽然目前嵌入式系统资源不像以前使用8032时候那么苛刻，但memory在运算量能保证情况下还是能省尽量省。

这个时候最简单用小Table获得更高精度的办法就是线性插值。在已知两个点之间的phase采用线性计算获得。

例如还是1024的Table，但是Phase可以用22bits表达，这样高10bits的值通过查表获得，低12bits就用于线性计算。

long lib_sin(long phase){    int sflag=0;    const long *pstab;    long sa, sb, da, db;    if (phase < 0)    {        phase = -phase;        sflag = -1;    }    da = phase & 0xFFF;    db = 0x1000 - da;    phase >>= 12;    pstab = tab_sine_1024_32bits_Q19 + phase;    sa = *pstab++;    sb = *pstab;    if (sflag < 0)    {        return -(sa*db + sb*da);    }    return (sa*db + sb*da);}

这样就能把1024 Table的最大误差缩小到-104.7dB。两个点之间的Phase差为0.0000015 <0.0000164 完全满足384Khz以及一般音频的应用。