CORDIC アルゴリズム

# CORDIC アルゴリズム

## 1. CORDIC アルゴリズムとは

* "CORDIC" は **CO**ordinate **R**otation **DI**gital **C**omputer の略語.
* 1956年に Jack E. Volder が考案したため, Volderのアルゴリズムとも呼ばれる.
* 三角関数, 双曲線関数, 平方根関数, 指数関数, 対数関数など様々な関数の近似値を, 加算, 減算のみを使用して求めることができる.
    →複雑な回路をもたない計算機でも計算できる.
* このドキュメントでは, 三角関数を近似するアルゴリズムを説明する.

## 2. 前提知識

### 三角関数 (正弦, 余弦, 正接) の定義

![figure1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/191898/7cfe329f-a41f-c0b4-681b-b9c6adba5b93.png)

図のような$∠C$を直角とする直角三角形において, $∠A=\theta$ に対する正弦, 余弦, 正接は次式で表される.

**正弦関数**
$\sin \theta = \frac a h$

**余弦関数**
$\cos \theta = \frac b h$

**正接関数**
$\tan \theta = \frac a b$

### 逆三角関数

**逆正弦関数** (本記事では直接扱いません)
$\sin \theta = \frac a h \Leftrightarrow \arcsin \frac a h = \theta$

**逆余弦関数** (本記事では直接扱いません)
$\cos \theta = \frac b h \Leftrightarrow \arccos \frac b h = \theta$

**逆正接関数**
$\tan \theta = \frac a b \Leftrightarrow \arctan \frac a b = \theta$

## 3. 基本的な考え方

* 起点 $P_0$ から座標変換 $M_{\pm\theta_i}$ を繰り返し適用し, 最終的な座標 $P_N$ を求めます.
    * $P_0 \xrightarrow {M_{\pm\theta_0}} P_1 \xrightarrow {M_{\pm\theta_1}} P_2, ..., P_N$
        ![figure2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/191898/811c540f-1f8b-f12c-9d99-acfaeb373576.png)
    * $\theta_i$ の前についている符号は, $\pm \theta_0 \pm \theta_1 ... \pm \theta_{N-1}$ が求めたい角 $\theta^*$ になるべく近くなるように選びます.
        * 例えば, $30^\circ$ を求めたい場合は, $+45^\circ-26.57^\circ+14.04^\circ...$ となります. (つまり $M_{+\theta_0}$, $M_{-\theta_1}$, $M_{+\theta_2}, ...$ を適用します)

* $P_N$ の $x$ $y$ 座標から正弦, 余弦, 正接を計算できます.

ここで,
* $\theta_i = ∠P_iOP_{i+1}  = \arctan \frac 1 {2^i}$ とします. これは定数なので, 事前に計算しておきます.
    * 下図のように, 縦の長さが $1$, $\frac 1 2$, $\frac 1 4$, $\frac 1 8$, ... となる三角形をイメージしてください.
        ![figure3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/191898/27b00824-ac80-f43d-2676-5b2d94438a08.png)
* $∠OP_{i}P_{i+1}=90^\circ$ です.

また,
* $M_{+\theta_i}(P_i(x_i,y_i))$ は, 三角形$A_i$ と 三角形$B_i$ が $1:\frac 1 {2^i}$ の比で相似であることに着目すると, 次式で表されます.
    * $x_{i+1} = x_i - \frac {y_i} {2^i}$
    * $y_{i+1} = y_i + \frac {x_i} {2^i}$
* 同様の考え方で, $M_{-\theta_i}(P_i(x_i,y_i))$ は次式で表されます.
    * $x_{i+1} = x_i + \frac {y_i} {2^i}$
    * $y_{i+1} = y_i - \frac {x_i} {2^i}$

![figure4.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/191898/76e42b23-adfd-26c7-b956-4cdb2ae1b844.png)

## 4. アルゴリズム

3.の内容を疑似言語で記述します.

求めたい三角関数の角度 $\theta^*$ が与えられたとします.

1. **初期化**
    * $P_0 \leftarrow (1,0)$
    * $\theta^t \leftarrow 0$

2. **繰り返し**
    * $i=0, ..., N-1$ で以下を繰り返します.
        * $\theta^t < \theta^*$ のとき
            * $\theta^t \leftarrow \theta^t + \theta_i$
            * $P_{i+1} \leftarrow M_{+\theta_{i}}(P_{i})$
                * $x_{i+1} \leftarrow x_i - \frac {y_i} {2^i}$
                * $y_{i+1} \leftarrow y_i + \frac {x_i} {2^i}$
        * $\theta^t > \theta^*$ のとき
            * $\theta^t \leftarrow \theta^t - \theta_i$
            * $P_{i+1} \leftarrow M_{-\theta_{i}}(P_{i})$
                * $x_{i+1} \leftarrow x_i + \frac {y_i} {2^i}$
                * $y_{i+1} \leftarrow y_i - \frac {x_i} {2^i}$

3. **結果の返却**
    * $P_N$ の座標 $(x_N, y_N)$ から, 正弦, 余弦, 正接の近似値が求められます.
        * $\sin\theta^* \approx \sin\theta^t = \frac {y_N} {\sqrt {x_N^2+y_N^2}}$
        * $\cos\theta^* \approx \cos\theta^t = \frac {x_N} {\sqrt {x_N^2+y_N^2}}$
        * $\tan\theta^* \approx \tan\theta^t = \frac {y_N} {x_N}$
        * $\sqrt {x_N^2+y_N^2}$ の部分は, 常に同じ値になるので, 事前に計算しておきます.

## 5. 補足

1. 4.の計算過程において, $2^i$ での除算演算がありますが, コンピュータは2進数でデータを保持しているので, 

    * 固定小数点の場合は, シフト演算
    * 浮動小数点の場合は, (指数部のみの)減算

    の演算に置き換えて計算時間を短縮することができます.

2. $\lim_{N\rightarrow\infty}{\sum_i{\theta_i}} = 1.743...$ と一定の値に収束するため, 近似可能な $\theta$ の範囲はおよそ $[-100^\circ, +100^\circ]$ の範囲であることが分かります.
    ![figure5.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/191898/e415aa48-e77d-df9b-6e5b-5f6674c743cf.png)

## 6. 参考文献など

* 本記事は, [^1] のWebページをヒントに作成しています.
* 本記事に挿入されている図は [^2] を使用して作成しています.

[^1]: <https://teamcoil.sp.u-tokai.ac.jp/calculators/column/100224/index.html>, "サルでも分かるCORDICアルゴリズム"
[^2]: <https://www.geogebra.org/>, "GeoGebra"