直接数字频率合成技术及其C++的实现

DDFS-Direct Digital Frequency Synthesizer 直接数字频率合成技术可以用来产生任意波形的周期信号。所谓的DDFS简单的说是查表法，内部维护一个Lookup Table储存一个周期的波形。那么这个查找表其实就是给出了相位到函数值的一个映射关系:

这里相位ω又是时间t的线性函数。

而调整输出频率实际上就是调整相位函数的系数kΔ。

举个具体的例子，比如我们的查找表存储的是正弦波的一个周期。查找表的长度为1024。正弦波的相位是以2π为周期的，那么查找表相邻的两项间的相位差为。也就是说查找表的相位分别为：

需要注意的是最后一个元素对应的相位并不是2π，因为相位与相位0其实是相同的，查找表中只能出现一次。

我们还需要一个相位累加器和一个频率寄存器，相位累加器用来记录当前的相位值是多少。比如说我们用一个16位无符号整数寄存器来记录这个相位值，那么0到65537(注意这里是65537)就对应0到2π这个相位周期。相位随时间变化是线性的，所以我们每次只需将相位累加器累加一个固定的相位就可以了。累加的这个值就存在频率寄存器中，比如频率寄存器中的值为5，那么相位累加器的值依次为：0,5,10,15,20,...,65530,65535,3,8,13,...。这里超过65536的值比如65537自动溢出成为0，刚好利用到了整型变量的这种周期性。

频率寄存器的值设为多少合适呢，如果我们的采样频率为 Fs，我们希望的输出频率为f。查找表的长度为M，相位累加器的长度为N，频率寄存器的值为k。那么相位累加器的最小相位变化量为。输出波形一个周期有Fs/f 个数据点。数据点间的相位差为。那么有如下等式：

实际计算出的k值不一定是个整数，四舍五入就可以了。举个具体的例子：采样频率为1kHz，我们希望输出10Hz的正弦波，相位累加器的长度为65536。那么频率寄存器的值应为。

可以看出，相位累加器的长度越长，频率分辨率就越高。因此实际使用时相位累加器通常会设计为32位或更高位数。查找表的长度通常就短很多，一般都是10位或12位。这时只需将相位累加器的最高几位作为查找表的指标就可以了。

直接数字频率合成技术并不仅仅是用于硬件生成特定函数波形，在数值计算领域，也可以采用这项技术。还是以生成正弦波为例。虽然大多数编程语言中都提供了sin()和cos()函数，但是这个函数的计算速度相对还是较慢的。如果你的程序中需要调用几百万次，那么花费的时间还是很可观的。另外，随着t的增加，sin(t)的计算精度也会下降。这时，采用DDS技术的优势就体现出来了。首先查表法保证它的速度比任何的现有的三角函数计算代码都要快的多，其次，所有的运算都是整数运算保证了不存在任何累计计算误差。

下面是实现代码，可以生成多种周期函数波形：

#ifndef DDS_H_INCLUDED
#define DDS_H_INCLUDED

class LookupTable
{
    friend class DDS;
private:
    unsigned int m_length;
    unsigned int m_shift;
    double *m_table;
public:
    LookupTable();
    LookupTable(unsigned int N);
    ~LookupTable();
    void sine(double amplitude = 1.0, double baseline = 0.0);
    void rect(double dutycircle = 0.5, double amplitude = 1.0, double baseline = 0.0);
    void sawtooth(double dutycircle = 0.5, double amplitude = 1.0, double baseline = 0.0);
    void arbitrary(double table[]);
};

class DDS
{
private:
    LookupTable *m_pLut;
    unsigned int m_iter;
    unsigned int m_stride;
    unsigned int m_phase;
public:
    DDS(LookupTable *pLut, unsigned int stride = 1, unsigned int phase = 0);
    void setPhase(double phase = 0.0);
    void setFreq(double freq);
    void setFreq(double freq, double fs);
    inline double value();
};
inline double DDS::value(void)
{
    unsigned int index = (m_iter + m_phase);
    index = index >> m_pLut->m_shift;
    double value = m_pLut->m_table[index];
    m_iter += m_stride;
    return value;
}
#endif // DDS_H_INCLUDED

#include "dds.h"
#include <cstdlib>
#include <cmath>
LookupTable::LookupTable()
{
    m_length = 0;
    m_table = NULL;
}
LookupTable::~LookupTable()
{
    delete[] m_table;
}
LookupTable::LookupTable(unsigned int N)
{
    if(N < 1)
    {
        N = 1;
    }
    if(N > 16)
    {
        N = 16;
    }
    m_length = 0x01 << N;
    m_shift = 32 - N;
    m_table = new double[m_length];
}

void LookupTable::sine(double amplitude, double baseline)
{
    if(m_length != 0)
    {
        for(unsigned int i = 0; i < m_length; i++)
        {
            m_table[i] = amplitude * sin(2 * M_PI / m_length * i) + baseline;
        }
    }
}

void LookupTable::rect(double dutycircle, double amplitude, double baseline)
{
    if(m_length != 0)
    {
        for(unsigned int i = 0; i < m_length; i++)
        {
             double value = (1.0 * i / m_length < dutycircle) ? 0 : 1;
             m_table[i] = amplitude * value + baseline;
        }
    }
}

void LookupTable::sawtooth(double dutycircle, double amplitude, double baseline)
{
    double pos, value;
    if(m_length != 0)
    {
        for(unsigned int i = 0; i < m_length; i++)
        {
            pos = 1.0 * i / m_length;
            if(pos < dutycircle)
            {
                value = pos / dutycircle;
            }
            else
            {
                value = (1.0 - pos) / (1.0 - dutycircle);
            }
            m_table[i] = amplitude * value + baseline;
        }
    }
}
void LookupTable::arbitrary(double table[])
{
    if(m_length != 0)
    {
        for(unsigned int i = 0; i < m_length; i++)
        {
            m_table[i] = table[i];
        }
    }
}
DDS::DDS(LookupTable *pLut, unsigned int stride, unsigned int phase)
{
    m_iter = 0;
    m_pLut = pLut;
    m_stride = stride;
    m_phase = phase;
}
void DDS::setPhase(double phase)
{
    m_phase = round(phase /(2 * M_PI) * 4294967296.00);
}
void DDS::setFreq(double omega)
{
    double stride = omega / (2 *M_PI) * 4294967296.00;
    m_stride = round(stride);
}
void DDS::setFreq(double freq, double fs)
{
    double omega = 2 *M_PI *freq / fs;
    setFreq(omega);
}

#include <iostream>
#include "dds.h"
using namespace std;

int main()
{
    LookupTable sincos(10);
    LookupTable sawtooth(10);
    sincos.sine(2, 0);
    sawtooth.sawtooth(0.3, 1, 0);
    DDS dds1(&sincos);
    DDS dds2(&sincos);
    DDS dds3(&sawtooth);
    dds1.setFreq(10, 1000);
    dds1.setPhase(0.0);
    dds2.setPhase(3.1415926/2);
    dds2.setFreq(10, 1000);
    dds3.setFreq(10, 1000);
    for(int i =0; i <200; i++)
    {
        cout << i << ", " << dds1.value() << ", " ;
        cout << dds2.value() << ", " << dds3.value() << endl;
    }
    return 0;
}

将输出结果绘出图形如下：

希望这个代码对大家有用。