坐标旋转数字计算机 (coordinate rotation digital computer, CORDIC)的概念起源于笛卡尔坐标系中的角度旋转，算法采用一组基本角度的加减组合来逼近目标角度，单个基本角度的旋转可以通过移位和加法运算完成^[1]，这使得其具有精度可以调节、硬件实现简便、速度快的优点. CORDIC算法本身作为一个强大的数学运算工具，可以用来作通用的科学计算^[2-3]，因而在快速傅里叶变换、离散余弦变换和数控振荡器、直接数字频率综合器(direct digital frequency synthesizer, DDFS)等数字信号处理和通信领域有着极其广泛的应用^[4-7].

传统的直接数字频率综合器采用只读内存(read only memory, ROM)查找表完成从累加器输出相位到幅度的转换，然后经数模转换器(digital to analog converter, DAC)和低通滤波器(low pass filter, LPF)完成数字信号到模拟信号的转换^[8]. 将CORDIC算法应用于DDFS中的相幅转换模块，与ROM查找表相比，可以有效降低高精度DDFS所需占用的存储资源，提高电路速度.

传统的流水线迭代式CORDIC算法依赖于迭代次数保证算法精度，本次迭代的旋转方向需要等待上次迭代完成才能确定^[9]. 对于精度要求高的情况，需要迭代很多次，导致输出时延过长，硬件消耗大. 随着数字通信系统向实时性和高集成度趋势的发展，研究者们提出多种改进方案：史方显等^[10]提出选择迭代式的CORDIC算法，有效减少了迭代次数，但z通路的复杂运算使得其存在硬件消耗较大的问题. 张佳等^[11]设计的三角函数运算单元提高了精度，但所耗费的资源很大. 徐成等^[12]提出免除缩放因子的双步旋转机制，将两次基本旋转合并成一次迭代，减少了一半的迭代次数，但输出时延很大；祁艳杰等^[13]采用并行CORDIC算法模块代替查找表，对进位保存加法器进行改进，提高了运算速度，但所消耗的资源增加，精度较低.

针对目前CORDIC算法实现时存在的时延大和资源消耗大的问题，结合三步旋转机制、角度二极化重编码技术和合并迭代技术，提出改进的低时延低消耗CORDIC算法. 结合加法树技术，实现了20位定点小数的高精度低时延DDFS，与其他方案相比，在资源、时延和精度上都有优化.

CORDIC的基本思想如下：通过基本角度集 $\left\{ {\left. {{\theta _i}} \right\}} \right.$包含的角度的线性组合来逼近输入角度. 如图1所示，向量 $\left( {{{ x}_n},{{ y}_n}} \right)$逆时针旋转角度 $\theta $，得到与其模值相同的向量 $\left( {{{{ x_m}}},{{{ y_m}}}} \right)$，这可以通过依次旋转一组基础角度 ${\theta _i}$来完成，即通过依次旋转 ${\theta _1}$， ${\theta _2}$, ···, ${\theta _N}$来完成.

以流水线迭代方式实现CORDIC算法时，取 ${\theta _i} = \arctan \,\,({2^{ - i}})$. 初始迭代值取(1,0)乘上模校正因子后的值，即 $({x_0},{y_0},{z_0}) = (K,0,\theta )$，其中 $K$为模校正因子，

(1) $K = \prod\limits_{i = 0}^{N - 1} {\cos \;{\theta _i}} .$

根据如下迭代公式，以流水线结构完成迭代^[9]：

(2) $\left. {\begin{array}{l} {{x_{i + 1}} = {x_i} - {d_i}{2^{ - i}}{y_i},} \\ {{y_{i + 1}} = {y_i} + {d_i}{2^{ - i}}{x_i},} \\ {{z_{i + 1}} = {z_i} - {d_i}{\theta _i},} \\ {{d_i} = {\rm sgn}\; ({z_{ i}}).} \end{array}} \right\}$

引入的 ${z_i}$表示剩余未旋转角度，用 ${d_i}$来指示旋转的方向， ${d_i}$取−1时进行逆时针旋转， ${d_i}$取1时进行顺时针旋转. 由于该次迭代的方向由 ${d_i}$上次迭代结果给出，输出时延随着迭代次数的增加而增加. 流水线结构提高了数据的吞吐率，但每步迭代都占用相应的硬件资源，因此硬件消耗随着迭代次数的增加而增加.

考虑到 $\left[ {0,{{\text{π}} / 8}} \right]$内的角度 $\theta $可以用N位的二进制小数形式表示为

(3) $\theta = \sum\limits_{i = 2}^N {{b_i}{2^{ - i}}} .$

式中： ${b_i}$取0或1，对 ${b_i}$的取值进行极化处理，即令 ${r_{i + 1}} = 2{b_i} - 1$，则 ${r_i}$取−1或1. 相应地， $\theta $可以表示为

(4) $\theta = {\varphi _0} + \sum\limits_{i = 3}^{N + 1} {{r_i}{2^{ - i}}} {\text{，}}$

(5) ${\varphi _0} = \sum\limits_{i = 3}^{N + 1} {{2^{ - i}}} = \frac{1}{4} - {\left( {\frac{1}{2}} \right)^{N + 1}}.$

这种对二进制小数位值进行极化处理的编码方式，称为角度二极化重编码^[14-15].

与取 $\left\{ {\left. {\arctan \,\,({2^{ - i}})} \right\}} \right.$角度集作为基本角度集的传统CORDIC算法不同，角度二极化重编码将输入角度 $\theta $分解为一个较大的角度 ${\varphi _0}$和若干个 ${2^{ - i}}$弧度的组合^[16]. 在坐标旋转时，先旋转一个较大的角度 ${\varphi _0}$，然后每次旋转 ${2^{ - i}}$弧度，加快了收敛到目标角度的进程. 此外，对于 ${b_i}$的极化处理不需要消耗额外的硬件资源，当 ${b_i} = 1$时，对应 ${r_{i + 1}} = 1$，逆时针旋转 ${2^{ - i}}$弧度；当 ${b_i} = 0$时，对应 ${r_{i + 1}} = - 1$，顺时针旋转 ${2^{ - i}}$弧度，每次的旋转方向在得到 $\theta $时已确定，即由 $\theta $的二进制位值决定，这既消除了流水线迭代式算法中剩余旋转角度运算通路所消耗的资源，也提高了算法的并行性.

经过角度二极化重编码处理后，单次的旋转迭代公式变为

(6) $\left. \begin{array}{l} {x_{i + 1}} = {x_i} - {r_i}\tan \left( {{2^{ - i}}} \right){y_i}, \\ {y_{i + 1}} = {y_i} + {r_i}\tan \left( {{2^{ - i}}} \right){x_i}. \end{array} \right\}$

对 $\tan \,\,({2^{ - i}})$进行泰勒展开可知， $ i \geqslant m =$ $ \left\lfloor {{{(N - {{\log }_2}3)} / 3}} \right\rfloor $时， $\tan \,\,({2^{ - i}}) \approx {2^{ - i}}$，因而式(6)可以化为

(7) $\left. \begin{array}{l} {x_{i + 1}} = {x_i} - {r_i}{2^{ - i}}{y_i}, \\ {y_{i + 1}} = {y_i} + {r_i}{2^{ - i}}{x_i}. \\ \end{array} \right\}$

对于 $m \leqslant i < {N / 2}$阶段的旋转迭代，采用三步旋转的方式来压缩迭代次数. 对于这一阶段的迭代来说， ${x_i} \times {2^{ - (3i + 3)}}$和 ${y_i} \times{2^{ - (3i + 3)}}$等于机器零，按式(7)连续迭代3次并除去机器零项，可得

(8) $\left. \begin{aligned} &{x_{i + 3}} = {x_i}\left( {1 - \frac{{{r_i}{r_{i + 1}}}}{{{2^{(2i + 1)}}}} - \frac{{{r_i}{r_{i + 2}}}}{{{2^{(2i + 2)}}}} - \frac{{{r_{i + 1}}{r_{i + 2}}}}{{{2^{(2i + 3)}}}}} \right) - \sum\limits_{j = i}^{i + 2} {{r_j}{2^{ - j}}{y_i}} , \\ & {y_{i + 3}} = {y_i}\left( {1 - \frac{{{r_i}{r_{i + 1}}}}{{{2^{(2i + 1)}}}} - \frac{{{r_i}{r_{i + 2}}}}{{{2^{(2i + 2)}}}} - \frac{{{r_{i + 1}}{r_{i + 2}}}}{{{2^{(2i + 3)}}}}} \right) + \sum\limits_{j = i}^{i + 2} {{r_j}{2^{ - j}}{x_i}} . \end{aligned} \right\}$

式(8)省去了中间迭代的过程，从 $({x_i},{y_i})$经三步旋转直接得到 $({x_{i + 3}},{y_{i + 3}})$，即在1个时钟周期内完成了3次迭代. 对应的结构图如图2所示.

三步旋转可以在1个时钟周期内完成3次迭代，减少了迭代次数. 通过图2可以看出，相比于式(6)的迭代，三步旋转可以节省硬件资源.

对于 $i < m$的迭代，即 $({x_m},{y_m})$的求解，可以通过多路选择器完成. 对于 $\left[ {0,{{\text{π}} / 8}} \right]$的角度 $\theta $来说， $i$从3开始取值，则初始迭代值 $({x_0},{y_0}) = (K,0)$在逆时针旋转 ${\varphi _0}$后可以得到 $({x_3},{y_3})$，计算公式为

(9) $\left. \begin{array}{l} {x_3} = {x_0} - \tan \left( {{\varphi _0}} \right){y_0}, \\ {y_3} = {y_0} + \tan \left( {{\varphi _0}} \right){x_0}. \\ \end{array} \right\}$

根据式(6)和r₃, ···, r_m−2, r_m−1的所有可能取值线下计算出 $({x_m},{y_m})$的所有可能取值，将 $({x_m},{y_m})$的所有可能取值存在一个小的ROM查找表或寄存器组中（对于20位的定点算法的实现来说，仅需要地址位宽为3 bit的查找表或寄存器组）. 在第一阶段时，直接通过一个多路选择器，从这个小的ROM查找表或寄存器组中选出正确的初始迭代值，然后进行三步旋转迭代. 由于m的取值小，该多路选择器的地址位 ${r_3} \cdot \cdot \cdot {r_{m - 2}}{r_{m - 1}}$位宽不会太大，与多次迭代式(6)相比，大大减少了从输入到输出的时延和耗费的硬件资源.

对于 $i \geqslant m$的迭代，由式(7)推导出两次迭代合并后的公式^[17]：

(10) $\left. \begin{array}{l} {x_{i + 2}} = {x_i} - {y_i}{r_i}{2^{ - i}} - {y_i}{r_{i + 1}}{2^{ - i - 1}} - {x_i}{r_i}{r_{i + 1}}{2^{ - 2i - 1}}, \\ {y_{i + 2}} = {y_i} + {x_i}{r_i}{2^{ - i}} + {x_i}{r_{i + 1}}{2^{ - i - 1}} - {y_i}{r_i}{r_{i + 1}}{2^{ - 2i - 1}}. \end{array} \right\}$

当 $i > {{(N - 1)} / 2}$时， ${x_i}{r_i}{r_{i + 1}}{2^{ - 2i - 1}}$和 ${y_i}{r_i}{r_{i + 1}}{2^{ - 2i - 1}}$等于机器零，因此对于这一阶段的n次迭代，可以合并成如下形式：

(11) $\left. \begin{array}{l} {x_{i + n}} = {x_i} - {y_i}\displaystyle\sum\limits_{j = i}^{k + n - 1} {{r_j}{2^{ - j}}} , \\ {y_{i + n}} = {y_i} + {x_i}\displaystyle\sum\limits_{j = i}^{k + n - 1} {{r_j}{2^{ - j}}} . \end{array} \right\}$

由式(11)可知，这一阶段由 $({x_i},{y_i})$可以一步计算出 $({x_{i + n}},{y_{i + n}})$，硬件实现时，用移位相加结构一步完成这一阶段的迭代，因而输出时延极大减小.

对于 $\left[ {0,{{\text{π}} / 8}} \right]$区间角度的正余弦值的求解，可以通过上面这种多路选择器、三步旋转和合并迭代的三阶段运算完成. 对于 $\left( {\left. {{{\text{π}} / 8},{{\text{π}} / 4}} \right]} \right.$区间的角度，可以用如下的两角差公式进行求得：

(12) $\left. \begin{aligned} &\sin\,\, \left(\frac{\text{π}}{4} - \theta \right) = \frac{{\sqrt 2 }}{2}(\cos \theta - \sin \theta ), \\ &\cos \,\,\left(\frac{\text{π}}{4} - \theta \right) = \frac{{\sqrt 2 }}{2}(\cos \theta + \sin \theta ). \end{aligned} \right\}$

可以看出，对 $\left[ {0,{{\text{π}} / 8}} \right]$区间角度的正余弦值作一次加减法，再乘上固定缩放因子 ${{\sqrt 2 } / 2}$可得 $\left( {\left. {{{\text{π}} / 8},{{\text{π}} / 4}} \right]} \right.$区间角度的正余弦值. 利用三角函数的对称性，可以求得 $\left[ {0,2{\text{π}} } \right]$区间角度的正余弦值，即整个 $\left[ {0,2{\text{π}} } \right]$区间的角度可以折叠到 $\left[ {0,{{\text{π}} / 8}} \right]$区间.

根据前面描述的算法原理，算法实现的整体流程如下：根据输入角度的二进制位值，从可能的初始值中选择正确的初始值完成第1阶段的迭代；然后通过三步旋转得到第2阶段的迭代结果，再通过合并迭代完成第3阶段的迭代；最后根据两角差公式和三角函数对称性，进行角度还原. 角度还原包括一次加法运算和缩放 ${{\sqrt 2 } / 2}$运算，其中缩放 ${{\sqrt 2 } / 2}$可以通过移位相加结构实现，因为 ${{\sqrt 2 } / 2}$可以用二进制小数的形式表示为(0.101 101 010 000 010 1)₂. 利用硬件描述语言Verilog，按照如图3所示的算法框图，实现了输出位宽为16位(包括1位符号位)的三步旋转CORDIC算法，输入数据位宽取18位(高4位指示角度所属区间范围)，中间计算结果采用18位.

在相同的输出位宽(N=16)的情况下，比较流水线迭代式算法和三步旋转算法的误差. 利用Matlab实现了流水线迭代式CORDIC算法和提出的三步旋转算法，开展定点化误差仿真，仿真时输入角度的分辨率取2⁻¹⁵，角度取 $\left[ {0,{{\text{π}} / 2}} \right]$，误差分析时均采用绝对误差. 用Verilog实现了输出位宽相同的2种算法，利用Vivado软件进行综合来评估资源消耗. 2种算法的精度、资源和时延的对比如表1所示. 表中，T_clk为时针周期.

从表1 可以看出，三步旋转算法的查找表(look-up table, LUT)和寄存器(flip-flop, FF)开销总共为1 146个，与传统流水线算法相比，硬件开销下降了29.2%；三步旋转CORDIC算法的输出时延为4个时钟周期(T_clk)，与传统流水线算法相比下降了75%，这是由于三步旋转算法采用三步旋转和合并迭代的方法对迭代次数进行了压缩. 三步旋转CORDIC算法的最大幅值误差为18.7×10⁻⁵，比传统算法的误差小. 传统算法采用θ_i=arctan (2⁻ⁱ)的线性组合来逼近输入角度 $\theta $，经过有限次的旋转迭代后，总旋转的角度与输入角度的差值较大；三步旋转CORDIC算法采用角度二极化重编码的方式（见式（4））对输入角度进行分解，对于角度 ${\varphi _0}$的旋转，根据式（9）采用线下计算的方式取得更高的精度，而后直接根据输入角度的二进制位值决定旋转 ${2^{ - i}}$弧度的方向，基于这两点，三步旋转CORDIC算法与传统算法相比，能够达到更高的精度.

利用改进的三步旋转CORDIC算法，实现DDFS中的相幅转换模块，并将结果送到DAC，结构框图如图4所示. 累加器负责对频率控制字进行累加，将累加得到的角度进行区间折叠；然后将角度转换到 $\left[ {0,{{\text{π}} / 8}} \right]$，根据转换后的角度位值选择正确的值进行三步旋转，经合并迭代后得到正余弦值；再根据三角函数公式和对称性得到原角度的正余弦值，送到DAC.

相位累加器对给定的频率控制字进行累加得到相位，当累加到 $2{\text{π}}$后，将累加结果减去 $2{\text{π}}$；接着累加，可得多个周期的正余弦波形.

为了提高电路的系统速度，对算法中多输入的加法器进行优化. 采用加法树结构将加法运算流水化，可以使电路满足更严格的时序约束，提高电路的系统时钟频率. 3级加法树结构的示意图如图5所示，即将操作数两两相加；然后用流水线寄存器进行寄存，对寄存后的结果进行两两相加然后寄存，直到最后一级输出. 算法实现时考虑速度和资源的折中，可以调整单级加法树中加法器的输入个数.

基于图4，用Verilog 硬件描述语言实现了20位定点小数的DDFS，即相位累加器的位宽取24，其中包括1位符号位；然后将相位折叠到 $\left[ {0,{{\text{π}} / 8}} \right]$，即用20位小数表示折叠后的相位. 通过Modelsim仿真得到的正余弦波形如图6所示. angle_init信号是0到2π线性增加的相位，溢出后再从0开始增加. angle信号为折叠到 $\left[ {0,{{\text{π}} / 8}} \right]$后的相位，sin和cos信号是生成的两路正交信号.

在Matlab中对DDFS进行20位小数的定点化仿真，仿真所得的正弦值和余弦值绝对误差如图7所示. 图中，E_s为正弦值误差，E_c为余弦值误差， $\theta $为输入角度. 图7中黑点标示了误差的最大值，可以看出算法的幅度误差不超过8.24×10⁻⁶.

将算法实现后输出的数据导入Matlab中，对输出波形进行频谱分析，计算出无杂散动态范围(spurious free dynamic range, SFDR). 输出信号频率为120 MHz时的频谱图（500 MHz采样频率）如图8所示. 图中，f为信号频率，A_mp为信号幅度，SFDR为105.03 dB.

对于实现的Verilog代码进行综合，选择Xilinx的XC5VLX30型号的现场可编程门阵列(field programmable gate array, FPGA)芯片. 基于该算法和文献[10]、[13]、[18]算法实现的DDFS的最大幅度误差、输出时延及硬件消耗的对比，如表2所示. 可以看出，本文的硬件开销比文献[10]、[13]都小，这是由于采用角度二极化重编码省去了剩余旋转角度的计算，即消除了z通路的运算. 角度二极化重编码后旋转的角度是直接由角度的位值决定的，这使得本文的精度更高. Xin等^[18]采用循环迭代和流水线相结合的方法，重复使用移位寄存器和加法器减少了硬件开销，但是这使得输出时延很大. 虽然Xin等^[18]的硬件开销为本文的37.48%，但是输出时延是本文的4.27倍. 本文采用三步旋转和合并迭代技术极大地压缩了迭代次数，利用加法树对多输入加法器进行适当的优化，提高了电路的工作频率，达到了更低的输出时延，即38.5 ns.

针对流水线迭代式CORDIC算法硬件开销大和输出时延大的缺点，提出基于三步旋转机制的低时延低消耗CORDIC算法. 从数学上分析三步旋转算法降低时延和资源消耗的原理，即利用三步旋转和合并迭代压缩迭代次数，结合区间折叠和角度二极化重编码技术消除剩余旋转角度运算通路；用Verilog语言实现了算法，通过仿真对比三步旋转算法和流水线迭代式算法的性能. Matlab仿真结果和Vivado综合结果表明：与流水线迭代式CORDIC算法相比，在精度稍有改善的情况下，三步旋转算法的硬件开销降低了约29.2%，从输入角度到输出正余弦值的时延减少了75%. 利用该算法，结合加法树技术实现了DDFS. 仿真结果表明，与以往的实现方案相比，在精度、硬件开销和输出时延方面均有一定的优化.