10

我正在以 48 kHz 采样正弦波，我的正弦波的频率范围可以从 0 到 20000 Hz 变化，步长约为 100 Hz。我正在使用查找表方法。因此，我为 4096 个不同相位的正弦波生成 4096 个样本。我认为这背后的一般想法是增加步长并为不同的频率使用不同的步长。所以我做了以下（伪代码）。但我不确定步长与我想要生成正弦波样本的频率有什么关系？例如，如果我的频率是 15000 Hz，我必须遍历的步长是多少？我的样本量（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];
 }

提前致谢。

4

5 回答 5

16

你在正确的轨道上 - 首先我们需要生成一个正弦波 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;
}

phase_i注意：为了获得更高质量的输出，您可以在和处的 LUT 值之间使用线性插值phase_i + 1，但上述方法对于大多数音频应用来说已经足够了。

于 2012-11-20T05:51:17.577 回答
5

通过查找表方法，人们可能希望有效地使用内存并仅存储正弦波的第一象限。

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),

其中 alpha=频率的相位差== 2pi*fHz/Ts，fHz 是要产生的频率（以赫兹为单位），Ts 是采样时间（或 1/Ts = 采样频率，例如 44100 Hz）。

这导致了谐振反馈滤波器方法，其传递函数 f(z) 在单位圆 (z=e^jomegaT) 处具有两个共轭极点。

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

有趣的部分是可以随时更改 alpha (cos(alpha))。这种 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'') 产生完美的圆，但即使使用定点算法也会产生稳定的椭圆。（请注意，这假设 alpha 值相对较小，也意味着频率相对较低。）

于 2012-11-20T06:50:17.890 回答
3

If you want high precision, you can use a trig identities to have both a small LUT, and clean sine waves.

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];
}
于 2012-11-20T17:11:04.770 回答
0

很好的答案，这是经典软件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 来满足这个要求。

于 2013-07-24T21:10:28.380 回答
-1

根据“http://en.wikipedia.org/wiki/X86_instruction_listings”，如果您有 x80387 或更新版本，则有一条正弦指令，所以直接调用它即可。您只需要弄清楚如何将一些内联汇编语言添加到您的程序中。这样，如果您的输入值与表中的值不完全匹配，您就不必担心。

于 2012-11-20T05:23:38.913 回答