微 处 理 机MICROPROCESSORS
・大规模集成电路设计、制造与应用・CORDIC算法在FPGA中的实现王智霞,王广生(北京工业大学电控学院,北京100022) 摘 要:CORDIC算法是在许多角度计算方面有着广泛应用的经典算法,通过考虑FPGA的结构、精度局限和速度要求,采用流水线技术(pipeline),在FPGA上用CORDIC算法实现了对于大吞吐量数据的向量倾角的计算,并对实际应用中内部步骤寄存器精度的选取给出了较为详细的方法。关键词:坐标旋转数字计算;FPGA;流水线中图分类号:TN4 文献标识码:B 文章编号:1002-2279(2007)01-0004-04
FPGABasedRealizationofCORDICAlgorithmWANGZhi-xia,WANGGuang-sheng(Beijinguniversityoftechnology,Beijing100022,China)
Abstract:CORDICalgorithmisaclassicalgorithmwithmanyapplications.Consideringthearchi2tecture,precisionandspeedofFPGA,thepipelinetechnologyisusedincomputinglargenumberofvectoranglevalues.Thisarticleprovidesthewayofconfirmthewideofinnerprocessingregister.Keywords:CORDIC;FPGA;Pipeline
1 引 言FPGA以其灵活性和使用方便在现今的数字领域已经得到了广泛的应用。但FPGA实现数字系统也有其自身的局限性,其一是器件资源的门阵列规模的限制,其二是单元延迟限制,所以,这就需要设计者充分考虑器件的实际工作能力。角度的旋转计算在数字领域尤其是数字通信领域是一种应用非常广泛的计算,如果用传统的除法器、乘法器等计算方法,需要占用大量的FPGA资源,这样就不能满足设计者的要求,需要设计者考虑其他的算法实现这种类型的计算。CORDIC算法在硬件电路的实现上只用到了加法器和移位器,这样就大大节约了FPGA的资源,从而可以满足设计者的要求。
2 CORDIC算法简介CORDIC(CoordinateRotationDigitalComputer),又名:坐标旋转数字计算,是J.Voider等人于1959
年在设计美国航空导航控制系统的过程中提出来的一种算法。下面就简要地介绍一下CORDIC算法的
基本数学思想。如图1所示,向量OA逆时针旋转θ度角得到向量OB,这个关系可以用矩阵表示如式1
[1]
:
XjYj=cosθ -sinθsinθ cosθXiYi=cosθ1 -tanθtanθ 1Xi
Yj
(1)
图1 向量旋转坐标图如果假设θ是由n个θn角度叠加而成的,那么根据式(1)得出每一步的叠加操作需要按照式(2)行Xn+1Yn+1=cosθn1 -tanθntanθn 1Xn
Yn
(2)
利用式2经过n步叠加可以表示由向量OA旋转到向量OB,如下表示:
XjYj=cosθ0・cos
θ
1...cosθn1 -tanθ0tanθ0 1...
1 -tanθntanθn 1Xi
Yi
(3)
作者简介:王智霞(1980-),女,北京人,硕士研究生,主研方向:嵌入式系统设计与实现,大规模集成电路系统设计,数字通信。收稿日期:2005-02-23
第1期2007年2月 No.1Feb.,2007王智霞等:CORDIC算法在FPGA中的实现
由于计算机进行计算采用二进制形式,所以我们选取θn=arctan(12n),这样选取θn方便了tanθn的计算,即tanθn=Sn2-n,式3前面的cosθn可以取累乘的极限即:
K=∏∞n=0cosθn=∏∞n=0cos[arctan(12n)]≈0.607253
如果我们在设计的系统中提前计算K,那么当抛开K不算时,式(2)就可以表示成式(4)[1]:
Xn+1Yn+1=1 -Sn2-nSn2-n 1Xn
Yn
θ=∑mi=0Siθi (Si={-1,1})(4)
至此,我们可以得出结论,由向量OA,在先计算K的情况下,我们可以由式4逐步的计算旋转角度后得出向量OB。计算的精度由n的大小决定,式(4)中的
S
n
由每一步的具体情况而定。
3 CORDIC算法的应用利用CORDIC算法可以进行许多运算,比如计算角度正弦、角度余弦、复数的乘法、向量的角度和模的运算等等。下面就介绍两种CORDIC算法的实际应用实例:
(1)已知:角度θ,求
:cosθ?
这个问题我们可以转换为利用CORDIC这种向量旋转的思想进行解析,建模如下(见图2):
图2 向量旋转坐标图起始向量OA为[X
i,Yi,θi]
=[1,0,0]
终止向量OB为[X
j,Yj,θj]=[cosθ,sinθ,
θ]
由OA经过n步旋转到OB,即可以得到cos
θ
设:Zn表示经过n步旋转后,得到的结果与θ
的差值,即:Z
n=θ-∑ni=0θi,通过这样的假设,
就可以
得到:
Sn=1 ZnΕ0-1 Zn<0(5)将(4)和(5)式结合,就可以得到它的逻辑表述:
IFZ(n)>=0X(n+1)=X(n)-Y(n)×2
-n
Y(n+1)=Y(n)+X(n)×2
-n
Z(n+1)=Z(n)-arctan(2
-n)
ELSEX(n+1)=X(n)+Y(n)×2
-n
Y(n+1)=Y(n)-X(n)×2
-n
Z(n+1)=Z(n)+arctan(2
-n)
(2)已知:点A(Xi,Yi),求:向量OA在坐标系
中的倾角θ?
这个问题我们可以转换为利用CORDIC这种向量旋转的思想进行解析,建模如下(见图3):
图3 坐标系旋转图如图3所示,为了方便对角度θ的计算,我们旋转了一下坐标轴,可以推断出起始向量OA为[X
i,
Yi,Zi]=[Xi,Yi,0]终止向量OB为[X
j,Yj,Zj]=[R,0,-θ]
设Z
n=∑∞n=0θn,就可以得到:
Sn=1 YnΦ0-1 Yn>0(6
)
将式(4)和式(6)结合,就可以得到它的逻辑表述:
IFY(n)<=0X(n+1)=X(n)-Y(n)×2
-n
Y(n+1)=Y(n)+X(n)×2
-n
Z(n+1)=Z(n)-arctan(2
-n)
ELSEX(n+1)=X(n)+Y(n)×2
-n
Y(n+1)=Y(n)-X(n)×2
-n
Z(n+1)=Z(n)+arctan(2
-n)
4 CORDIC的FPGA实现下面就上述的第二种应用实例,描述一下用verilog语言实现CORDIC算法的具体步骤和实现方法。(1)系统要求接口类型和精度
1)输入XY的精度为(11,10,T):10位整数、无小数位的有符号数。2)输出角度为Zπ,Z为弧度值的π前的系数,精度为(8,1,T):6位小数、1位整数的有符号数,它
・5・ 1期微 处 理 机
的范围是(-1,1)。3)系统的接口,见表1,其中clk为80MHz表1 cordic模块的接口表接口名称长度类型clk1INrst1INenin1INxin11INyin11INokout1OUTzout8OUT(2)计算中间寄存器宽度的选取分析系统要求的输入X和Y的精度和输出的Z的精度进行比较,可以分析出,输出和输入的精度不对称,即:Z的精度为(8,1,T)对应的最小可以表示的角度为2.8125°XY的精度为(11,10,T),对应的可以表示的最小的角度为:45°-arctan(10221023)=0.028017可见,输入和输出误差的允许值相差非常远,这样就得到了两种精度的分析。1)如果允许输出的角度Z的最后一位不准确,那么可以以2.8125为精度进行分析,分析如下:假设计算中用到的角度精度为(9,0,T),即:8位小数,1位符号位的值,这个值与z值相同,都表示弧度值的π前的系数,小数位的最高位表示1/4π(45度),所以它能表示的最小角度为0.3515625度。根据上述假设,arctan(2-n)可以表示为下表2:表2 arctan2-n值实际值取值计算值误差arctan1/226.5650500000000100110026.718750000000.1537arctan1/414.0362400000000010100014.062500000000.02626arctan1/87.125020000000000101007.031250000000.09377arctan1/163.576330000000000010103.515624500000.060706arctan1/321.789910000000000001011.757812500000.032097误差合计0.366533 从表2中我们不难看出,由于arctan(2-n)用精度为(9,0,T)表示,有部分的误差,误差的最大合计为0.366533度,这样可以得出如下结论:用(9,0,T)表示中间计算的角度寄存器的精度,运行5步叠加结果的最大误差为:0.366533(前者的误差)+1.78991(后者的误差)=2.156443度上述计算出的最大误差为2.156度,小于2.8125度,符合z的精度范围,计算结果只会在最后一位上与z的真实结果有出入,误差在允许范围之内。2)如果不允许输出角度Z有任何的偏差,6位小数位必须全部准确,那么需要用0.028017精度进行计算,分析方法同上,最后得出的选取中间计算角度的精度是(15,0,T),即:14位小数,1位符号位的值,这个值与z值相同,都表示弧度值的前系数,小数位的最高位表示1/4π(45度),需要进行12步计算,这个过程要比第一个过程需要增加数倍的资源,
而且速度会明显变慢。根据综合后的经验得知9位加法器的延迟时间大约为3ns,而16位加法器的延迟时间大约为5ns,
系统输入的clk为80MHz,周期为1/80MHz=12.5
ns,16位加法器(XnY
n
的计算)的延迟时间可以满足
clk的速度,所以考虑到运算速度的要求,使用第一种精度和步数选取结果,即:Zn的精度为(9,0,T),5
步迭代运算。(3)一个计算过程的处理
一个计算过程的实现步骤如图4所示。
图4 计算过程步骤图1)由于系统要求输入点A(X,Y)的范围为4个象限,而文章的二.2中讲述的方法只是针对ⅠⅣ象限的,所以我们需要在开始的计算前判断点A
(X,Y)所在的象限,并将它转换到Ⅰ象限;另外,由
于我们是判断Y(n)正负来确定S(n)的正负,角度小于45的角,角度的变化对于Y值的影响比较明显,所以如果所求向量角度在45～90度的范围内,
需要先将X和Y互换,转换到0～45度的范围内进行计算,这样可以减少误差。下面介绍一个具体的例子,说明一下具体的变换过程:如图5已知点A在第三象限,经过变换后变换到第一象限变成(-X,-Y),经过判断-X<