這里主要先介紹如何利用CORDIC算法計算固定角度\(\phi\)的\(cos(\phi)\)、\(sin(\phi)\)值。
一般利用MATLAB計算三角函數時，用\(cos\)舉例，只需要輸入相應的\(cos(\phi)\)便自動計算出來了。但是如果是硬件處理或者沒有那么方便的函數時，該如何計算\(cos(\phi)\)的值呢？
有一種最傻瓜的方式是用rom存儲\(0^o\)到\(90^o\)所有的余弦值，然后用查表的方法計算，但隨着精度要求的提升，需要存儲的值會越來越多，這是不合適的。那么有沒有一種用較少資源且能較快計算出高精度三角函數值的方法呢？有！這就是CORDIC算法可以做的事，算法不難，有一點三角函數和移位運算的基礎就能看懂了，核心思想是數形結合計算三角函數，計算過程把乘法運算變成移位運算。
為什么說是數形結合呢，讓我們看下如何通過旋轉求\(cos\phi\)和\(sin\phi\)。如下圖可以看出，如果讓初始坐標\((x_1, y_1)\)取\(x\)軸一個特定的位置，最后使其旋轉\(\phi^o\)到坐標\((x_n,y_n)\)處，且滿足\(\sqrt{x_n^2 + y_n^2} = 1\)，那么\(cos\phi\)不就等於\(x_n\)，\(sin\phi\)不就等於\(y_n\)了嗎，這是后話。
下面首先引入一個在圓上的坐標旋轉公式（不一定單位圓）：
注：\((1)\)的推導
其中\(A\)為旋轉圓半徑，\(\alpha\)為\((x_1, y_1)\)的角度值，圖中紅線與x軸夾角。則坐標\((x_2, y_2)\)可以寫為：
推導完畢。
我們繼續看式\((1)\)，\((1)\)也可以寫成：
由於對\((x_2, y_2)\)相當於同乘了一個常數\(cos\theta\)，不影響旋轉角度，只影響旋轉后坐標到原點的距離。我們先不看它，得到：
一：乘法變移位
之前說的我們要把乘法運算變成移位運算，所以我們找到\(tan\theta\)與\(2^{-i}\)之間的對應關系，注意由於是變成移位操作，所以對應旋轉的角度也是幾個固定的值，但是通過旋轉這幾個固定的角度，旋轉\(i\)次，最終也一定能轉到我們需要的角度\(\phi\)上(\(-99.7^o\le\phi \le 99.7^o\))。
於是把\((3)\)再改寫為：
這樣，旋轉\(\theta^o\)就變成了移位、相加的操作。注意\(d(i) = ±1\)表示旋轉的逆、順時針。
比如要旋轉\(\phi = 66^o\)，可以先轉\(+45^o\)；\(45^o < 66^o\)，再轉\(+26.565^o\)；\(45^o+26.565^o \ge 66^o\)，再轉\(-14.036^o \cdots\)，最終會逼近\(66^o\)。而整個運算僅僅進行了\(2^{-0}、2^{-1}、2^{-2} \cdots\)移位操作和加法操作。
二：cos累乘項
現在考慮把\(cos\theta\)加回去，回到\((2)\)，且考慮旋轉方向\(d_i\)和旋轉角度\(\theta_i\)，得到：
進行下一次迭代(旋轉)，得到：
可以看到每次旋轉都可以提取出\(cos\theta_i\)。\(tan\theta_i\)已經用移位替代了。接下來只用計算\(\prod_{i=1}^{N}cos\theta_i\)就行了，且\(\prod_{i=1}^{N}cos\theta_i\)只跟迭代次數有關，大概收斂到0.607。確定了迭代次數后，可以預先把\(\prod_{i=1}^{N}cos\theta_i\)算出來。
三：累計旋轉角度與旋轉方向
現在考慮最后一個問題，如何確定每次迭代的旋轉方向\(d_i\)呢？其實定義一個累計旋轉角度\(z_i\)
令\(z_1\)等於目標角度值，然后每次迭代作個判斷就好，如果\(z_i > 0\)，說明當前旋轉還沒轉到目標角度，\(d_{i+1} = 1\)；如果\(z_i < 0\)，說明當前旋轉超過了目標角度，\(d_{i+1} = -1\)。
當我們最終轉到了目標角度\(\phi\)時，比如\(\phi = 66^o\)，可以此時\(z_i\)已經很小趨近於零了。
另外，在作比較判斷時，單次旋轉角度\(\theta_i\)則還是需要通過查一次\(arctan(2^{-i})\)表得到，但這個表比起文章開頭說的傻瓜式查表要小太多了。
四：計算cos和sin值
由開始所說的數形結合可以推出，假設只旋轉一次把\((x_1, y_1)\)轉到\((x_n, y_n)\)，有：
按第一節所說的，把一次就旋轉到位用很多次旋轉替代，每次旋轉后坐標到原點的距離比上一次都大了\(\frac{1}{cos\theta_i}\)倍，所以如果我們把\((x_1, y_1)\)設為\((\prod_{i=1}^n cos\theta_i, 0)\)，那么最后\((x_n, y_n)\)到原點的距離就為1，\(cos\phi\)就等於\(x_n\)，\(sin\phi\)就等於\(y_n\)了
這樣就只用通過\((4)(6)\)的移位、相加、一次查表，再迭代十多次，就能計算\(cos\phi\)和\(\sin\phi\)值啦！
其實用CORDIC算法還能計算arctan、sinh、cosh等值，以后學習了再來補充進階版。
五：完整MATLAB代碼
clc, clear, close all;
%% 初始計算cos累乘值
N = 16; % 設置迭代次數16次
Nprod(1) = 1;
for i = 1 : N
Nprod(i + 1) = Nprod(i) * cos(atan(2^(-(i - 1))));
end
%% Cordic算法計算cos、sin值
x(1) = Nprod(N); % 橫坐標初始值賦為cos累乘值
y(1) = 0; % 縱坐標初始值賦為0
z(1) = 66 / 180 * pi; %目標旋轉角度值，66°，注意轉化成弧度值
d(1) = 1; %旋轉方向，初始肯定為1
for i = 1 : N
x(i + 1) = x(i) - d(i) * y(i) * 2^(-(i - 1)); %移位運算，公式(4)
y(i + 1) = y(i) + d(i) * x(i) * 2^(-(i - 1)); %移位運算，公式(4)
z(i + 1) = z(i) - d(i) * atan(2^(-(i - 1))); %計算累計旋轉角度，查一次表，公式(6)
if z(i + 1) >= 0 % 判斷下一次的旋轉方向
d(i + 1) = 1;
else
d(i + 1) = -1;
end
end
COS = x(N) % 輸出cos66的值
SIN = y(N) % 輸出sin66的值
六：參考文章
[1]：https://mp.weixin.qq.com/s/c4oro0bOhdDUmBt0yyLpTA
[2]：https://blog.csdn.net/qq_39210023/article/details/77456031