计算定点正弦和余弦- MATLAB & Simulink示例金宝app - 金宝app,下载188bet金宝搏,金宝搏官方网站

计算定点正弦和余弦

打开脚本

此示例展示了如何使用定点设计器™提供的基于cordic和基于查找表的算法来近似MATLAB®正弦(罪)和cos (因为)功能。高效的定点正弦和余弦算法对许多嵌入式应用至关重要，包括电机控制，导航，信号处理和无线通信。

用CORDIC算法计算正弦和余弦

介绍

的cordiccexp，cordicsincos，cordicsin,cordiccos函数近似MATLAB罪和因为函数使用基于cordic的算法。CORDIC是坐标旋转数字计算机的缩写。Givens基于旋转的CORDIC算法(参见[1,2])是硬件效率最高的算法之一，因为它只需要迭代的移位加操作。CORDIC算法消除了显式乘法器的需要，适用于计算各种函数，如正弦、余弦、反正弦、反余弦、反正切、矢量大小、除法、平方根、双曲和对数函数。

您可以使用CORDIC旋转计算模式来计算正弦和余弦，以及极坐标到直角坐标的转换操作。在这种模式下，矢量大小和旋转角度是已知的，旋转后计算坐标(X-Y)分量。

CORDIC旋转计算模式

CORDIC旋转模式算法首先初始化具有所需旋转角度的角度累加器。接下来，每次CORDIC迭代中的旋转决策以减小剩余角度累加器的大小的方式完成。旋转决策是基于每次迭代后角度累加器中剩余角度的符号。

旋转模式下，CORDIC方程为:

$z_{i+1} = z_{i} - d_{i}*atan(2^{-i})$

$x_{i+1} = x_{i} - y_{i}*d_{i}*2^{-i}$

$y_{i+1} = y_{i} + x_{i}*d_{i}*2^{-i}$

在哪里 $d_{i} = -1$ 如果 $z_{i} < 0$ , $$ +1 $$ 否则;

$$ I = 0,1，…， n-1 $$ , $$ n $$ 是迭代的总次数。

这将提供如下结果为 $$ n $$ 方法 $+\infty$ ：

$z_{N} = 0$

$$间{N} =现代{N}(间的{0}\因为{z_ {0}} - y_{0} \罪{z_ {0}}) $$

$$ y_ {N} =现代{N} (y_{0} \因为{z_{0}} +间的{0}\罪{z_ {0}}) $$

地点:

$A_{N} = \prod_{i=0}^{N-1}{\sqrt{1+2^{-2i}}}$ ．

在旋转模式下，CORDIC算法受限于旋转角度之间 $-\pi/2$ 和 $\pi/2$ ．为了支金宝app持超出这个范围的角度cordiccexp，cordicsincos，cordicsin,cordiccos函数在CORDIC迭代完成后使用象限校正(包括可能的额外否定)。

理解`CORDICSINCOS`正弦和余弦代码

介绍

的cordicsincos函数使用CORDIC算法计算[-2*pi 2*pi)范围内输入角的正弦和余弦。这个函数取一个角度 $\theta$ (弧度)和迭代次数作为输入参数。函数返回正弦和余弦的近似值。

CORDIC计算输出由旋转器增益缩放。这个增益是通过预缩放初始值来计算的 $1 / a_ {n}$ 常数值。

初始化

的cordicsincos函数执行以下初始化步骤:

角度输入查询表inpLUT设为atan(2 .^ (0:N-1))．
$z_{0}$ 设置为 $\theta$ 输入参数值。
$x_{0}$ 设为 $1 / a_ {n}$ ．
$y_{0}$ 设为0。

初始值的明智选择使算法可以同时直接计算正弦和余弦。后 $$ n $$ 迭代，这些初始值导致以下输出为 $$ n $$ 方法 $+\infty$ ：

$x_{N} \近似cos(\theta)$

$y_{N} \近似sin(\theta)$

共享的定点和浮点CORDIC内核代码

CORDIC算法(旋转模式)内核部分的MATLAB代码如下(对于标量的情况)x，y,z）.同样的代码用于定点和浮点操作:

函数[x, y, z] = cordic_rotation_kernel(x, y, z, inpLUT, n)%对N个内核迭代执行CORDIC旋转内核算法。XTMP = x;Ytmp = y;为Idx = 1:n如果z < 0 z(:) = z + inpLUT(idx);X (:) = X + ytmp;Y (:) = Y - xtmp;其他的z(:) = z - inpLUT(idx);X (:) = X - ytmp;Y (:) = Y + xtmp;结束XTMP = bitsra(x, idx);乘以2^(-idx)的%右移位Ytmp = bitsra(y, idx);乘以2^(-idx)的%右移位结束

可视化正弦余弦旋转模式的CORDIC迭代

CORDIC算法通常通过指定(恒定)次数的迭代来运行，因为提前结束CORDIC迭代会破坏管道代码，并导致CORDIC增益 $A_{n}$ 不是常数，因为 $$ n $$ 会有所不同。

对于非常大的值 $$ n $$ ， CORDIC算法保证收敛，但并不总是单调的。如下面的示例所示，中间迭代偶尔会产生比后期迭代更准确的结果。通常可以通过增加迭代的总次数来获得更高的准确性。

例子

在下面的示例中，迭代5提供了比迭代6更好的结果估计，并且CORDIC算法在以后的迭代中收敛。

= /5;%以弧度为单位的输入角Niters = 10;迭代次数%sinTh = sin()参考结果%cost = cos(θ);参考结果%Y_sin = 0(整数，1);Sin_err = 0 (niters, 1);X_cos = 0(整数，1);Cos_err = 0 (niters, 1);流(“恐怖\ n \ n \ nNITERS \”）;流(”——\ t ----------\ n '）;为N = 1:niters [y_sin(N)， x_cos(N)] = cordicsincos(theta, N);sin_err(n) = abs(y_sin(n) - sinTh);cos_err(n) = abs(x_cos(n) - cost);如果N < 10 fprintf(' %d \t %1.8f\n'， n, cos_err(n));其他的流(' %d \t %1.8f\n'， n, cos_err(n));结束结束流(' \ n '）;

错误码------ ---------- 1 0.10191021 2 0.13966630 3 0.03464449 4 0.03846157 5 0.00020393 6 0.01776952 7 0.00888037 8 0.00436052 9 0.00208192 10 0.00093798

在柱状图上绘制CORDIC近似误差

图(1);clf;栏(1:硝石,cos_err(1:硝酸钠));包含(“迭代次数”，“字形大小”12“fontweight”，“b”）;ylabel (“错误”，“字形大小”12“fontweight”，“b”）;标题(cos(pi/5)计算的CORDIC近似误差，…“字形大小”12“fontweight”，“b”）;轴([0位0.14]);

绘制5次迭代的X-Y结果

Niter2Draw = 5;图(2)，clf, hold在情节(cos(0:0.1:π/ 2),sin(0:0.1:π/ 2),“b——”）;%半圆为i=1:Niter2Draw plot([0 x_cos(i)]，[0 y_sin(i)]，“线宽”2);% CORDIC迭代结果文本(x_cos(我),y_sin(我),int2str(我),“字形大小”12“fontweight”，“b”）;结束情节(cos(θ),sin(θ),的r *，“MarkerSize”, 20);%理想结果包含(“X (COS)”，“字形大小”12“fontweight”，“b”) ylabel (“Y”(罪)，“字形大小”12“fontweight”，“b”)标题(“cos(pi/5)计算的CORDIC迭代”，…“字形大小”12“fontweight”，“b”)轴平等的;轴广场;

计算定点正弦值`cordicsin`

在[-2*pi, 2*pi]之间创建1024个点

stepSize = pi/256;thRadDbl = (-2*pi):stepSize:(2*pi - stepSize);thRadFxp = sfi(thradfbl, 12);%带符号，12位定点值sinThRef = sin(double(thRadFxp));参考结果%

比较定点CORDIC与双精度三角函数的结果

使用12位量化输入，迭代次数从4到10不等。

为cdcSinTh = cordicsin(thRadFxp, niters);errCdcRef = sinthf - double(cdcSinTh);图;持有在;轴([-2*pi 2*pi -1.25 1.25]);情节(thRadFxp sinThRef,“b”）;情节(thRadFxp cdcSinTh,‘g’）;情节(thRadFxp errCdcRef,“r”）;ylabel (“sinθ(\)”，“字形大小”12“fontweight”，“b”）;集(gca),“XTick”2 *π:π/ 2:2 *π);集(gca),“XTickLabel”，…｛“2 *π”，“3 *π/ 2”，“-π”，“-π/ 2”，…' 0 '，π/ 2的，“π”，“3 *π/ 2”，“2 *π”});集(gca),“YTick”, 1:0.5:1);集(gca),“YTickLabel”, {“-1.0”，“-0.5”，' 0 '，“0.5”，“1.0”});ref_str =参考:罪(双(\θ));Cdc_str = sprintf(的12位CORDICSIN;N = %d'、硝石);Err_str = sprintf('错误(max = %f)'马克斯(abs (errCdcRef)));图例(ref_str, cdc_str, err_str);标题(cdc_str,“字形大小”12“fontweight”，“b”）;结束

计算N = 10时的LSB误差

图;fracLen = cdcSinTh.FractionLength;plot(thRadFxp, abs(errCdcRef) * pow2(fracLen));集(gca),“XTick”2 *π:π/ 2:2 *π);集(gca),“XTickLabel”，…｛“2 *π”，“3 *π/ 2”，“-π”，“-π/ 2”，…' 0 '，π/ 2的，“π”，“3 *π/ 2”，“2 *π”});ylabel (sprintf ('LSB错误:1 LSB = 2^{-%d}'fracLen),“字形大小”12“fontweight”，“b”）;标题(LSB错误:12位CORDICSIN;N = 10’，“字形大小”12“fontweight”，“b”）;轴([-2*pi 2*pi 0 6]);

计算噪声底限

fft_mag = abs(fft(double(cdcSinTh)));Max_mag = max(fft_mag);Mag_db = 20*log10(fft_mag/max_mag);图;持有在;情节(0:1023 mag_db);情节(0:1023 0 (1024),“r——”）;%归一化峰值(0 dB)情节(0:1023 -62。* (1024),“r——”）;%噪声底电平ylabel (数据库级的，“字形大小”12“fontweight”，“b”）;标题(62 dB噪声底层:12位CORDICSIN;N = 10’，…“字形大小”12“fontweight”，“b”）;%轴([0 1023 -120 0]);完整的FFT轴([0 round(1024*(pi/8)) -100 10]);%放大集(gca),“XTick”，[0 round(1024*pi/16) round(1024*pi/8)];集(gca),“XTickLabel”, {' 0 '，“π/ 16”，“π/ 8”});

加速固定点`CORDICSINCOS`函数与`FIACCEL`

您可以使用MATLAB®从MATLAB代码生成MEX函数fiaccel函数。通常，运行生成的MEX函数可以提高仿真速度，尽管实际速度的提高取决于所使用的仿真平台。下面的示例展示了如何加速定点cordicsincos函数使用fiaccel．

的fiaccel函数将MATLAB代码编译成一个MEX函数。此步骤需要创建临时目录并在该目录下拥有写权限。

tempdirObj = fidemo.fiTempdir(“fi_sin_cos_demo”）;

当您将迭代次数声明为常量时(例如:10)使用coder.newtype(“常数”,10)，编译后的角度查找表也将是常量，因此不会在每次迭代时计算。还有，当你打电话的时候cordicsincos_mex，则不需要为它提供迭代次数的输入参数。如果传入迭代次数，则mex函数将出错。

输入参数的数据类型决定是否cordicsincos函数执行定点或浮点计算。当MATLAB为该文件生成代码时，仅为特定的数据类型生成代码。例如，如果THETA输入参数为定点，则只生成定点代码。

inp = {thRadFxp, coder.newtype()“不变”10)};%函数的示例输入fiaccel (“cordicsincos”，“o”，“cordicsincos_mex”，“参数”,输入)

首先，通过调用计算正弦和余弦cordicsincos．

Tstart = tic;cordicsincos (thRadFxp 10);telapsed_mccordicsincos = toc(tstart);

接下来，通过调用mex函数计算正弦和余弦cordicsincos_mex．

cordicsincos_mex (thRadFxp);%加载MEX文件Tstart = tic;cordicsincos_mex (thRadFxp);telapsed_MEXcordicsincos = toc(tstart);

现在，比较一下速度。在MATLAB命令行中输入以下命令以查看平台上的速度改进:

fiaccel_speedup = telapsed_mccordicsincos /telapsed_MEXcordicsincos;

清空临时目录，使用如下命令:

清晰的cordicsincos_mex;status = tempdirObj.cleanUp;

计算`罪`和`因为`使用查阅表

有许多基于查找表的方法可用于实现定点正弦和余弦近似。下面是基于单个实值查找表和简单的最近邻线性插值的低成本方法。

基于单个查找表的方法

的罪和因为方法fi对象中的定点设计器近似于MATLAB®内置浮点数罪和因为函数，使用基于查找表的方法与简单的最近邻值之间的线性插值。这种方法允许使用一个小的实值查找表，并使用简单的算术。

使用单个实值查找表简化了索引计算和获得非常准确的结果所需的总体算法。这些简化产生了相对较高的速度性能和相对较低的内存需求。

理解基于查找表的`罪`和`因为`实现

查找表的大小和精度

查找表的两个重要设计考虑因素是它的大小和准确性。不可能为每个可能的输入值创建一个表 $$ u $$ ．由于的量子化，也不可能完全准确 sin(u) 或 $ cos(u) $$ 查找表值。

作为妥协，定点设计器罪和因为的方法FI使用8位查找表作为其实现的一部分。一个8位的表只有256个元素长，所以它是小而高效的。在许多平台上，8位也相当于一个字节或一个字的大小。与线性插值和16位输出(查找表值)精度结合使用，8位可寻址查找表提供了非常好的精度和性能。

初始化常量罪查找表值

为了实现简单，表值一致性和速度，使用全正弦波表。首先是四分之一浪潮罪函数在[0,pi/2)弧度范围内以64个均匀间隔采样。为表值选择一个带符号的16位小数定点数据类型，即tblValsNT = numerictype(1,16,15)，可产生最佳的精度结果罪输出范围[-1.0,1.0]。这些值在设置之前进行预量化，以避免溢出警告。

tblValsNT = numerictype(1,16,15);quarterSinDblFltPtVals = (sin(2*pi*((0:63) ./ 256)))';endpointQuantized_Plus1 = 1.0 - double(eps(fi(0,tblValsNT));halfSinWaveDblFltPtVals =…[quarterSinDblFltPtVals;…endpointQuantized_Plus1;…flipud (quarterSinDblFltPtVals (2)));fullSinWaveDblFltPtVals =…[halfSinWaveDblFltPtVals;-halfSinWaveDblFltPtVals];FI_SIN_LUT = fi(fullSinWaveDblFltPtVals, tblValsNT);

算法实现概述

定点设计器的实现罪和因为的方法fi对象涉及首先投射定点角度输入 $$ u $$ (以弧度为单位)转换为[0,2 pi]范围内的预定义数据类型。为此，执行模2运算以获得定点输入值inpValInRange在[0,2 pi]范围内，并转换为最佳精度的二进制点缩放的无符号16位定点类型numerictype(0, 16日,13)：

%实际值范围[0,2 *pi]的最佳UNSIGNED类型，%，映射到定点存储的整型值[0,51472]。inpInRangeNT = numerictype(0,16,13);

接下来，我们从这个范围内的定点FI角度值中获得16位存储的无符号整数值:

idxUFIX16 = fi(storedInteger(inpValInRange)， numerictype(0,16,0));

我们将存储的整数值乘以一个归一化常数65536/51472。生成的整数值将在uint16的完整索引范围内:

normConst_NT = numerictype(0,32,31);normConstant = fi(65536/51472, normConst_NT);fullScaleIdx = normConstant * idxUFIX16;idxUFIX16(:) = fullScaleIdx;

这个全尺寸无符号16位索引的前8位最高有效位(msb)idxUFIX16用于直接索引到8位正弦查找表。执行两个表查找，一个在计算表索引位置lutValBelow，另一个在下一个索引位置lutValAbove：

idxUint8MSBs = storedInteger(bitsliceget(idxUFIX16, 16,9));zeroBasedIdx = int16(idxUint8MSBs);lutValBelow = FI_SIN_LUT(zeroBasedIdx + 1);lutValAbove = FI_SIN_LUT(zeroBasedIdx + 2);

的剩余8个最低有效位(LSBs)idxUFIX16用于在这两个表值之间进行插值。LSB值被视为具有8位分数数据类型的规范化比例因子rFracNT：

rfrnt = numerictype(0,8,8);%小数余数数据类型idxFrac8LSBs = reinterpretcast(bitsliceget(idxUFIX16,8,1)， rFracNT);rFraction = idxFrac8LSBs;

实乘法用于确定两点之间的加权差。这将导致一个简单的计算(相当于一个乘积和两个和)，以获得插值的定点正弦结果:

temp = rFraction * (lutValAbove - lutValBelow);rslt = lutValBelow + temp;

例子

使用上述算法，这里是一个例子，查找表和线性插值过程中用来计算的值罪对于一个定点输入inpValInRange = 0.425弧度:

%使用任意输入值(例如，0.425弧度)inpInRangeNT = numerictype(0,16,13);%最佳精度，[0,2 *pi]弧度inpValInRange = fi(0.425, inpInRangeNT);%任意定点输入角度将其存储的整数规范化以获得完整的无符号16位整数索引idxUFIX16 = fi(storedInteger(inpValInRange)， numerictype(0,16,0));normConst_NT = numerictype(0,32,31);normConstant = fi(65536/51472, normConst_NT);fullScaleIdx = normConstant * idxUFIX16;idxUFIX16(:) = fullScaleIdx;%使用无符号8位整数索引(即8个msb)执行两次表查找idxUint8MSBs = storedInteger(bitsliceget(idxUFIX16, 16,9));zeroBasedIdx = int16(idxUint8MSBs);%从零开始的表索引值lutValBelow = FI_SIN_LUT(zeroBasedIdx + 1);%第一个表查找值lutValAbove = FI_SIN_LUT(zeroBasedIdx + 2);%第二个表查找值%使用8个lsb进行最近邻插值(作为分数余数处理)rfrnt = numerictype(0,8,8);%小数余数数据类型idxFrac8LSBs = reinterpretcast(bitsliceget(idxUFIX16,8,1)， rFracNT);rFraction = idxFrac8LSBs;%分数值用于线性插值temp = rFraction * (lutValAbove - lutValBelow);rslt = lutValBelow + temp;

这是算法结果的示意图:

X_vals = 0:(pi/128):(pi/4);xIdxLo = 0;xIdxHi = zeroBasedIdx + 4;图;持有在;axis([x_vals(xIdxLo) x_vals(xIdxHi) 0.25 0.65]);情节(x_vals (xIdxLo: xIdxHi)、双(FI_SIN_LUT (xIdxLo: xIdxHi)),“b ^——”）;情节([x_vals (zeroBasedIdx + 1) x_vals (zeroBasedIdx + 2)),…[lutValBelow lutValAbove),“k”。）;最接近的值%图(0.425,双(rslt),的r *）;%插值定点结果罪情节(0.425,(0.425),“gs”）;%双精度参考结果包含(“X”）;ylabel (“SIN (X)”）;lut_val_str =“定点查找表值”;near_str =“两个最接近的定点LUT值”;interp_str =“插值不动点结果”;ref_str =“双精度参考值”;图例(lut_val_str, near_str, interp_str, ref_str);标题(基于线性插值的SIN定点设计查找表，…“字形大小”12“fontweight”，“b”）;

计算定点正弦使用`罪`

在[-2*pi, 2*pi]之间创建1024个点

stepSize = pi/256;thRadDbl = (-2*pi):stepSize:(2*pi - stepSize);%双精度浮点数thRadFxp = sfi(thradfbl, 12);%符号，12位定点输入

比较定点SIN和双精度SIN的结果

fxpSinTh = sin(thRadFxp);%定点结果sinThRef = sin(double(thRadFxp));参考结果%errSinRef = sinthif - double(fxpSinTh);图;持有在;轴([-2*pi 2*pi -1.25 1.25]);情节(thRadFxp sinThRef,“b”）;情节(thRadFxp fxpSinTh,‘g’）;情节(thRadFxp errSinRef,“r”）;ylabel (“sinθ(\)”，“字形大小”12“fontweight”，“b”）;集(gca),“XTick”2 *π:π/ 2:2 *π);集(gca),“XTickLabel”，…｛“2 *π”，“3 *π/ 2”，“-π”，“-π/ 2”，…' 0 '，π/ 2的，“π”，“3 *π/ 2”，“2 *π”});集(gca),“YTick”, 1:0.5:1);集(gca),“YTickLabel”, {“-1.0”，“-0.5”，' 0 '，“0.5”，“1.0”});ref_str =参考:罪(双(\θ));Fxp_str = sprintf(“带12位输入的16位定点SIN”）;Err_str = sprintf('错误(max = %f)'马克斯(abs (errSinRef)));图例(ref_str, fxp_str, err_str);标题(fxp_str,“字形大小”12“fontweight”，“b”）;

计算LSB错误

图;fracLen = fxpSinTh.FractionLength;plot(thRadFxp, abs(errSinRef) * pow2(fracLen));集(gca),“XTick”2 *π:π/ 2:2 *π);集(gca),“XTickLabel”，…｛“2 *π”，“3 *π/ 2”，“-π”，“-π/ 2”，…' 0 '，π/ 2的，“π”，“3 *π/ 2”，“2 *π”});ylabel (sprintf ('LSB错误:1 LSB = 2^{-%d}'fracLen),“字形大小”12“fontweight”，“b”）;标题(LSB错误:16位定点sin12位输入，“字形大小”12“fontweight”，“b”）;轴([-2*pi 2*pi 0 8]);

计算噪声底限

fft_mag = abs(fft(double(fxpSinTh)));Max_mag = max(fft_mag);Mag_db = 20*log10(fft_mag/max_mag);图;持有在;情节(0:1023 mag_db);情节(0:1023 0 (1024),“r——”）;%归一化峰值(0 dB)情节(0:1023 -64。* (1024),“r——”）;%噪声底电平(dB)ylabel (数据库级的，“字形大小”12“fontweight”，“b”）;标题(64 dB本底噪声:16位定点SIN, 12位输入，…“字形大小”12“fontweight”，“b”）;%轴([0 1023 -120 0]);完整的FFT轴([0 round(1024*(pi/8)) -100 10]);%放大集(gca),“XTick”，[0 round(1024*pi/16) round(1024*pi/8)];集(gca),“XTickLabel”, {' 0 '，“π/ 16”，“π/ 8”});

比较不动点逼近算法的代价

定点CORDIC算法需要以下操作:

1表查找每个迭代
2班每个迭代
3增加每个迭代

采用最近邻线性插值的简化单查找表算法需要进行以下操作:

2个表查找
1乘法
2添加

在实际应用中，为定点三角函数计算选择算法通常取决于所需的精度、成本和硬件约束。

关闭所有;%关闭所有图形窗口

参考文献

Jack E. Volder, CORDIC三角计算技术，电子计算机学报，vol . 8, 1959年9月，pp330-334。
Ray Andraka，基于FPGA的计算机CORDIC算法的研究，1998年第六届现场可编程门阵列国际学术研讨会论文集，2月22-24日，1998,pp191-200

这个话题有帮助吗?

文档