ARTICLE DETAIL

资讯详情

深耕网站建设、视觉设计与SEO优化的一线实战洞察。

CORDIC算法在FPGA中的实现:移位加减法搞定三角函数与坐标变换

CORDIC算法在FPGA中的实现:移位加减法搞定三角函数与坐标变换 简介面向FPGA开发者的CORDIC算法学习与实现资料包覆盖坐标旋转数字计算机算法的原理、仿真与硬件实现适合数字信号处理、无线通信及图像处理方向的工程师和进阶学习者。资源共7个文件包含2份PDF原理/Xilinx应用文档、2个Verilog实现与testbench仿真代码以及3个MATLAB辅助脚本用于增益系数、反正切表及正弦/余弦计算可帮助理解迭代收敛并快速验证设计压缩包整体32.25MB结构清晰。CORDIC通过移位与加法即可逼近三角函数、坐标变换等运算在FPGA中资源占用少、实时性强非常契合无线通信、FFT/DFT等场景。目前已有1611人学习下载是从算法原理到Verilog落地实践的实用参考资料。 前阵子在调一块接收机的中频处理逻辑需要同时输出信号的幅度和相位。起初图省事直接挂了两个查表IP资源倒是没爆但延迟和精度都不太理想。折腾了一圈最后老老实实回到CORDIC算法上来——用移位和加法把反正切、求模、sin/cos全办了。这个算法老但老得有道理。这篇就当是把我自己重新走了一遍的弯路整理出来给准备在FPGA或单片机里做三角函数、坐标变换的朋友做个参考尤其是那些和我一样一开始被“迭代”“收敛”“增益补偿”这些词劝退过的人。CORDIC全称是Coordinate Rotation Digital Computer坐标旋转数字计算机。它最值钱的地方在于算三角函数、双曲函数、开方、对数、乃至线性方程都可以只靠移位和加减法完成。这对没有硬件乘法器、或者不想为乘法器付出太多资源成本的场景非常友好。整套算法的复杂度可控精度随着迭代次数线性提升属于那种“只要肯翻手册就一定能调明白”的方案。1. 为什么硬件里算三角函数绕不开CORDIC1.1 查表法和多项式逼近的困境想在不依赖处理器软件库的情况下算一个sin或者arctan很多人的第一反应是查表。小角度、低精度的场景查表完全够用比如只有0到90度、步进0.1度、输出16位也就900个点BRAM随便塞。但一旦角度范围变大、分辨率变高表的大小会迅速膨胀。就拿16位角度输入来说直接查表需要65536个深度再加上插值用的系数表资源开销和设计复杂度都上去了。多项式逼近是另一个思路。泰勒展开在角度接近0的时候收敛很快可一旦角度往两端跑需要保留的项数就明显增加。每一项都涉及乘加运算多次迭代乘法器的组合逻辑链拉得很长时序很容易崩。而且多项式系数的定点化处理也是个体力活系数量化误差和截断误差叠在一起经常是调了一下午系数最后精度还是差那么两三个LSB。1.2 CORDIC把乘法换成了旋转CORDIC的聪明之处在于绕开了直接计算函数值而是把问题转化为“旋转一个向量到目标角度”。在圆周坐标系统下每步迭代只做一次角度旋转旋转的角度序列是预先设计好的固定值——arctan(2的负i次方)。这个序列的巧妙之处在于沿坐标轴方向的位移可以表示成移位而移位在硬件里是免费的不占乘法器资源。整个迭代过程只需要三类操作比较决定旋转方向、移位完成因子2^-i的缩放、加减法更新x、y、z三个累加器。这在FPGA里意味着什么意味着LUT和寄存器就能撑起整个运算单元DSP48模块可以腾出来去做滤波器、FFT等更需要乘法器的模块。我在实际项目里经常把这个特性当作选型的第一理由不是CORDIC算得最快而是它能把宝贵的乘法器留给真正非要乘法不可的地方。2. 向量旋转的数学底子迭代逼近如何省掉乘法器2.1 旋转矩阵与角度分解如果有一个点(x, y)想把它逆时针旋转角度θ标准做法是乘一个旋转矩阵x x·cosθ - y·sinθy x·sinθ y·cosθ直接这么做需要四次乘法和两次加法而且cosθ和sinθ本身还得用别的办法算出来。CORDIC把θ拆成一串可以“偷懒”的小角度之和每次只旋转一个预先算好的角度α_i其中满足tan(α_i) 2^-i也就是α_i arctan(2^-i)。旋转矩阵除以cosθ之后会变成下面这个形式x_{i1} x_i - d_i · y_i · 2^-iy_{i1} y_i d_i · x_i · 2^-iz_{i1} z_i - d_i · α_i这里的d_i是旋转方向取值为1或-1。乘上2^-i在二进制里就是右移i位硬件上一条移位总线就搞定不需要乘法器。之所以能在没有cos项的情况下继续旋转是因为cos(α_i)这个公共因子可以放到最后统一补偿这就是后面要说的增益问题。2.2 为什么偏偏选arctan(2^-i)这个角度序列选arctan(2^-i)而不是其他序列是CORDIC算法最核心的一步设计。第一这个序列满足“每次旋转后剩余角度可以被后续更小的角度逐步逼近”。换句话说任何在收敛范围内的角度θ都能表示成若干个±α_i之和且误差随迭代次数增加而减小。第二它保证了硬件实现只需要移位正好落在“没有乘法器也能算”的红线上。这里有个容易踩的坑——这个方法依赖角度范围的收敛性。圆周系统下所有α_i加起来大约等于99.88度。也就是说输入角度如果超出了大约±99.88度直接喂给CORDIC核心无法保证收敛。实际工程里几乎不会只处理这么小的范围所以通常会在前端做角度折叠把0到360度折叠到0到90度先利用三角函数的对称性把角度换算到第一象限再根据象限把旋转方向和符号修正回来。这个折叠逻辑本身很简单但忘记做的人不在少数。2.3 迭代次数与增益补偿精度和硬件的平衡点每做一次迭代向量的模长都会被拉长一个因子sqrt(1 2^(-2i))。迭代n次之后总增益趋近于一个常数大约1.646760258。为了得到正确的幅度输出必须在最后乘上这个系数的倒数也就是0.607252935。这就是CORDIC里著名的“增益补偿”。补偿的实现方式有两种。第一种是直接乘一个定点数这是最常见也最稳定的一种第二种是把补偿因子拆成移位和加减法做成“乘法器免费”的版本但代价是逻辑层数变多。在FPGA里如果DSP48有空余我建议直接用乘法器做补偿简单可靠如果DSP48紧张再考虑移位加减法补偿。迭代次数我一般取16到20轮16轮时幅度误差大约在0.005%量级20轮会更干净一些但每多一轮就多一级流水线资源消耗也随之增加需要根据系统精度指标来权衡。3. 圆周系统下的两种工作模式旋转模式和向量模式3.1 旋转模式给定角度算sin/cos旋转模式的目标是输入一个角度z_0通过一系列迭代把z_0逐步“消耗”到0同时得到旋转后的坐标。初始时把x_0设为增益补偿因子Ky_0设为0。每次迭代判断z_i的符号z_i 0就逆时针旋转d_i 1z_i 0就顺时针旋转d_i -1。迭代结束后x_n就是cos(z_0)y_n就是sin(z_0)。实际使用中有一个细节值得注意如果直接把x_0设成K而不是1那么最后输出就直接是余弦和正弦值省掉一次独立的乘法。如果不这么做就需要在迭代结束后统一乘上K的倒数。两种方式结果一样但前者的运算开销更低延迟更短。3.2 向量模式把坐标变成模长和相位向量模式解决的是另一个问题已知一个向量的x和y求它的模长和相位角。这也是我在接收机项目里最需要的功能。初始时把z_0设为0每次迭代根据y_i的符号决定旋转方向y_i 0就顺时针转y_i 0就逆时针转。目标是把y_i逐步“压”到0。迭代结束后x_n就是模长乘以增益z_n就是向量的相位角。这里面有个隐藏的好处模长输出x_n是正的所以输出的相位角的取值范围是-99.88度到99.88度。配合上象限判断逻辑就可以还原出0到360度甚至±180度的完整相位。我习惯在向量模式后面挂一个小的象限修正状态机几行代码的事能省去很多后续软件层的麻烦。3.3 线性与双曲系统不止是三角函数很多人不知道CORDIC还有另外两套坐标系统。线性系统下迭代公式里的2^-i因子变成纯加法/减法可以计算乘除法双曲系统下把迭代序列换成双曲旋转序列可以算双曲函数、指数、对数、开方。实际上那些硬件上做开方做得特别快的实现内部往往就是双曲CORDIC。不过要提醒一句双曲系统有一个特殊问题——某些迭代需要重复一次才能保证收敛比如i4、13、40这些序号。如果直接照搬圆周系统的迭代结构精度会莫名奇妙掉一截。所以如果只是做坐标变换建议先把圆周系统吃透双曲系统等确实需要开方或对数运算时再单独拉出来调。4. 从公式到可综合电路定点化关键细节4.1 定点数格式选择CORDIC对定点格式比较敏感。我的经验是先把输入输出位宽定下来再倒推内部迭代位宽。假设输入输出都是16位x、y路径的位宽建议在内部放到20位左右因为每轮迭代虽然只有加法和移位但增益累积可能带来溢出。角度累加器z的位宽最好覆盖整个0到360度的范围用有符号定点数表示高位留出象限信息。迭代过程中x和y的右移位数随i增大而增大。移位后需要按有符号数做符号扩展否则负数会出大问题。这一步经常被忽略尤其是从MATLAB浮点原型转到Verilog/VHDL的时候最容易出bug。4.2 迭代结构串行还是流水线如果资源极其紧张可以用串行结构一个迭代单元循环使用状态机控制迭代次数。这种结构面积最小但吞吐率低适合低速控制类应用。如果需要连续出数据比如通信基带里的相位解调就得用流水线结构每一级迭代独占一组寄存器数据像流水一样往下走每个时钟周期都能吐一个结果。流水线的代价是寄存器数量线性增长。16级流水大约需要16组x、y、z寄存器每级位宽20位的话光这个单元就要接近1000个寄存器。在资源评估阶段最好提前计算不要等综合完发现时序违例再回去拆。4.3 我实测过的精度和资源数据之前在一款中端FPGA上做过一个16级流水线CORDIC16位定点输入输出角度路和向量路各一套。角度模式下的最大绝对误差大约在2到3个LSB向量模式下相位误差在0.02度以内模长误差在0.1%以内。逻辑资源消耗大约1500个LUT配合不到100个寄存器块不同厂商器件差异较大仅供参考。这个量级对于绝大多数信号处理链路来说非常友好。如果精度不够优先检查两点一是增益补偿是否做在了正确的位置二是迭代过程中是否出现过中间值溢出。溢出这个问题很隐蔽因为大多数时候结果看起来“差不多”但到了边界角度就会突然冒出一个无规律的跳变。加宽内部位宽一两比特往往就能消除。5. 实际工程里的几个容易翻车的地方5.1 角度折叠逻辑的正确做法如果输入角度是0到360度的无符号定点数先拆出最高两位作为象限标志再根据象限把角度映射到0到90度区间。第一象限0到90度保持不变第二象限90到180度用180度减去原始角度第三象限180到270度减去180度第四象限270到360度用360度减去原始角度。旋转结束后根据原始象限对cos和sin的符号做修正。这个逻辑看着简单但边界情况特别容易出错尤其是正好落在90度、180度这些整点上。我建议在仿真阶段把边界角度全部枚举一遍而不是只挑几个随机角度测。很多时候bug就藏在“315度 tiny offset”这类组合里。5.2 旋转方向判断的时序问题在流水线结构里每一级的d_i判断依赖当前级的z_i或y_i符号这是一个典型的组合逻辑路径。如果位宽比较大比如24位以上路径延迟可能会成为时序瓶颈。解决办法是在每级寄存器的输出端直接接比较器不要等进入下一级再判断或者在某些级数上插入额外的流水寄存器牺牲一个周期的延迟换取更好的时序收敛。另外如果用有符号数做判断不要用最高位直接当作符号位去接MUX因为补码的符号位不能直接决定比较结果——你还需要排除“负零”之类的情况。正确做法是老老实实做一次有符号比较或者用最高位和零标志组合判断。5.3 与CORDIC电路相关的资源复用技巧如果你的系统里同时需要多个通道的坐标变换不妨共享一套CORDIC核心用时分复用方式轮转处理各通道数据。比如四通道信号每个通道的采样率不是特别高就可以用4倍时钟频率跑一个CORDIC输出侧用FIFO把各通道结果分开。这样资源消耗直接从四套变成一套加少量缓存。另一个技巧是如果同时需要sin/cos和arctan而两者不会并发使用可以设计成模式可配置的CORDIC通过一个mode信号在旋转模式和向量模式之间切换。虽然流水线级数不能省但至少可以把两套逻辑合并成一套节省一部分控制逻辑。6. 绕开手柄的点我最初被卡住的三个细节我最早接触CORDIC的时候在三个地方卡了很久现在拿出来说说希望后来的朋友能少花点冤枉时间。第一个是迭代的起始序号。在圆周系统里迭代序号从i0开始但在双曲系统里需要从i1开始才能保证双曲旋转收敛。如果我当时能早点意识到“不是所有CORDIC都从0开始”就不会拿着一套双曲迭代公式去对圆周系统的参考实现对半天对不上。第二个是增益补偿的时机。补偿必须在所有迭代都完成之后做不能在迭代中间提前补偿否则每级旋转角度都会偏离设计的arctan(2^-i)序列导致收敛失败。这个错误我犯过一次结果输出的sin值在角度大时会明显偏离理想曲线排查起来很痛苦。第三个是仿真时一定要用“最坏情况”的角度组合。很多参考代码在测试时只测了0度、30度、45度、60度、90度这些常规角度看着精度很好。一旦把角度放到接近收敛边界比如97度或者放到0度附近的小角度问题立刻暴露。建议仿真时加一个随机角度扫描把所有角度的误差曲线都拉出来看一遍哪里尖峰高就重点排查哪里。个人体会是CORDIC这个算法属于那种“上手不快但一旦吃透用起来非常稳”的类型。它的数学原理不复杂难的是定点化实现时各种细枝末节。只要你踩过了上面这些坑后续再遇到坐标变换、相位计算、三角函数硬件化基本都能从容应对。后续有机会的话我打算再写一篇双曲CORDIC在开方和对数运算里的具体用法那个内容也不少到时候再和大家细聊。本文还有配套的精品资源点击获取
返回列表