如何查找不同频率的正弦从固定大小的查找表?(How to look up sine of diffe

2019-08-02 23:13发布

我以48kHz采样的正弦波,我的正弦波的频率范围可以与约100Hz的步骤从0到20000赫兹变化。 我使用的是查找表的方法。 于是,我产生4096个样本为4096个不同相位的正弦波。 我觉得这背后的总体思路,以递增步长和使用不同的步长为不同frequncy。 所以我做了以下(伪代码)。 但我不知道该步长将如何涉及到我要生成的正​​弦波的采样频率是多少? 例如,如果我的频率为15000赫兹这将是步长,我要穿越? 是我的样品尺寸(4096)太低呢?

 // Pseudocode
 uint16_t audio_sample[4096] = {...};
 NSTEP = freq; //???How is the step size going to be related to the freq here

 for(int i = 0; i < 4096; i = i+NSTEP)
 {
     sine_f(i) = audio_sample[i];
 }

提前致谢。

Answer 1:

你在正确的轨道上 - 首先,我们需要生成一个正弦波LUT:

const int Fs = 48000;       // sample rate (Hz)
const int LUT_SIZE = 1024;  // lookup table size

int16_t LUT[LUT_SIZE];      // our sine wave LUT

for (int i = 0; i < LUT_SIZE; ++i)
{
    LUT[i] = (int16_t)roundf(SHRT_MAX * sinf(2.0f * M_PI * (float)i / LUT_SIZE));
}                           // fill LUT with 16 bit sine wave sample values

请注意,我们只需要在初始化过程中一旦产生这种LUT,例如。

现在,我们有一个正弦波LUT,我们可以用它来产生我们希望利用相位累加器的任何频率:

const int BUFF_SIZE = 4096;  // size of output buffer (samples)
int16_t buff[BUFF_SIZE];     // output buffer

const int f = 1000;          // frequency we want to generate (Hz)

const float delta_phi = (float) f / Fs * LUT_SIZE;
                             // phase increment

float phase = 0.0f;          // phase accumulator

// generate buffer of output
for (int i = 0; i < BUFF_SIZE; ++i)
{
    int phase_i = (int)phase;        // get integer part of our phase
    buff[i] = LUT[phase_i];          // get sample value from LUT
    phase += delta_phi;              // increment phase
    if (phase >= (float)LUT_SIZE)    // handle wraparound
        phase -= (float)LUT_SIZE;
}

注:更高质量的输出,你可以在使用LUT值之间的线性插值phase_iphase_i + 1 ,但上述方法对于大多数音频应用足够好。



Answer 2:

随着查找表的方法不妨一有效地使用内存和存储仅正弦波的第一象限。

Then second quadrant = sin[180 - angle]  ; // for angles 90..180 degrees
third quadrant = -sin[angle-180];          // for 180..270
fourth quadrant = -sin[-angle+360]         // for 270..360

我会建议线性插值,也有对正弦代基于矢量旋转的方法(即产生正弦和余弦同时)

x_0 = 1, y_0 = 0;
x_n+1 = x_n * cos(alpha) - y_n * sin(alpha)
y_n+1 = x_n * sin(alpha) + y_n * cos(alpha), 

其中频率的阿尔法= 相位差 == 2PI * FHZ / Ts时,与FHZ是要生产的频率(以赫兹为单位)和Ts的采样时间(或1 / TS =采样frequenzy例如。44100赫兹)。

这导致谐振反馈滤波器方法,其传递函数F(z)具有在单位圆2个共轭极点(Z = E ^ jomegaT)。

y[n] = y[n-1] - 2*cos(alpha) * y[n-2], with
y[0] = 0,
y[1] = sin(alpha)

有趣的部分是一个可以在飞行中改变α(COS(阿尔法))。 这种IIR滤波器方法的缺点是,它的定义是不稳定的。 浮点和特别是固定点的不精确性累积并导致任一指数衰减,或大小的指数增长。 也就是说然而,可以用允许轻微的相位失真固化。

而不是作为在CORDIC旋转具有每次迭代放大率已知:

K = sqrt(1+sin(alpha)^2);

x' = x - y*sin(alpha);
y' = y + x*sin(alpha);

one can do

x' = x - y * sin(alpha);
y'' = y + x' * sin(alpha);

其不产生用于正圆(X”,Y'),但稳定的椭圆即使定点算术。 (注意,这假定阿尔法的相对小的值,也意味着相对低的频率。)



Answer 3:

如果你想精度高,可以用一个三角身份既有小LUT,干净的正弦波。

static float gTrigSinHi[256], gTrigSinLo[256];
static float gTrigCosHi[256], gTrigCosLo[256];

////////////////////////////////////////
// Sets up the fast trig functions
void FastTrigInit()
{
    unsigned i;
    for(i = 0; i < 256; ++i)
    {
        gTrigSinHi[i] = sin(2.0 * M_PI / 0xFFFF * (i << 8));
        gTrigSinLo[i] = sin(2.0 * M_PI / 0xFFFF * (i << 0));
        gTrigCosHi[i] = cos(2.0 * M_PI / 0xFFFF * (i << 8));
        gTrigCosLo[i] = cos(2.0 * M_PI / 0xFFFF * (i << 0));
    }
}

////////////////////////////////////////
// Implements sin as
//      sin(u+v) = sin(u)*cos(v) + cos(u)*sin(v)
float FastSin(unsigned short val)
{
    unsigned char hi = (val >> 8);
    unsigned char lo = (val & 0xFF);
    return gTrigSinHi[hi] * gTrigCosLo[lo]
        + gTrigCosHi[hi] * gTrigSinLo[lo];
}

////////////////////////////////////////
// Implements cos as
//      cos(u+v) = cos(u)*cos(v) - sin(u)*sin(v)
float FastCos(unsigned short val)
{
    unsigned char hi = (val >> 8);
    unsigned char lo = (val & 0xFF);
    return gTrigCosHi[hi] * gTrigCosLo[lo]
        - gTrigSinHi[hi] * gTrigSinLo[lo];
}


Answer 4:

很好的答案,这是经典的软件DDS。 面对同样的问题,这些天。 有没有必要使用花车

        UInt16 f = 400; // Hz, up to 16384 :)
        UInt16 delta = (UInt16)(((UInt32)f * LUT_SIZE ) / fmt_SamplesPerSec);
        UInt16 phase = 0;

        for (int i = 0; i < BUFF_SIZE; ++i)
        {
            buff[i] = LUT[phase];
            phase += delta;
            phase &= LUT_SIZE-1;
        }

让相位回绕LUT大小的面膜。 而不在乎使用象限,因为我的目的,我有这方面的要求已经是一个巨大的MIPS。



Answer 5:

据“http://en.wikipedia.org/wiki/X86_instruction_listings”如果你有x80387或更高版本有一个正弦指令,所以就直接调用它。 你只需要弄清楚如何在一些行内汇编语言添加到您的程序。 这样,您就不必担心,如果你的输入值不完全匹配,什么是您的表。



文章来源: How to look up sine of different frequencies from a fixed sized lookup table?