「2次元の点列にフィットする3次ベジェ曲線群を求める」の改良案

### はじめに
Graphics Gem1 の「AN ALGORITHM FOR AUTOMATICALLY FITTING DIGITIZED CURVES」
の改良案です。
とりあえず確実に効果のあった手法4つを記載します。
「AN ALGORITHM FOR AUTOMATICALLY FITTING DIGITIZED CURVES」については<a href="https://qiita.com/ML081/items/f8c6037ab1e85349a40e">こちら</a>で解説してあります。

### 点列の端点では接ベクトルを求める必要は無い
点列の端点では接ベクトルを求めずに曲線を最適化することにより、曲線の自由度が向上する。

### 点列の分割法
許容誤差を満たさない点列を分割するとき、対象とならない残りの点列はまとめる。<br>
Original algorithm
```math
\begin{aligned}
&1:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ}_{target\ \rightarrow\ \times}\\
&2:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\ \circ}_{target\ \rightarrow\ \times}\\
&\hspace{214px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{stack[0]}\\
&3:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ}_{target\ \rightarrow\ \circ}\\
&\hspace{114px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{stack[0]}\\
&\hspace{214px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{stack[1]}\\
&4:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ}_{\circ}\\
&\hspace{114px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{target\ \rightarrow\ \circ}\\
&\hspace{214px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{stack[0]}\\
&5:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ}_{\circ}\\
&\hspace{114px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{\circ}\\
&\hspace{214px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{target\ \rightarrow\ \circ}\\
\end{aligned}
```

New algorithm
```math
\begin{aligned}
&1:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ}_{target\ \rightarrow\ \times}\\
&2:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\ \circ}_{target\ \rightarrow\ \times}\\
&\hspace{214px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\ \circ}_{stack[0]}\\
&3:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ}_{target\ \rightarrow\ \circ}\\
&\hspace{114px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ}_{stack[0]}\\
&4:\ \underbrace{\circ\circ\circ\circ\circ\circ\circ}_{\circ}\\
&\hspace{114px}\underbrace{\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ\circ}_{target\ \rightarrow\ \circ}\\
\end{aligned}
```

### 最大回数まで弧長パラメータを更新し、フィッティングする
maxIterations まで計算すれば基本的に誤差は小さくなる。
しかし、iterations に対して誤差関数は単調減少ではないので、誤差が増えるようなら処理を中断する。
許容誤差を満たしたら、それを採用するのは勿体ない。

### 再フィッティング
点列を分割し、3次ベジェ曲線列を作成した後の処理である。
元の点列の端点以外の点は通る必要は無いので、隣り合う3つベジェ曲線をG1連続でフィッティングしなおす。<br>
```math
\begin{aligned}
&1:\ \underbrace{\circ\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}}_{fit\ \sharp 0\ and\ \sharp 1}\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\circ\\
&2:\ \circ\circ\circ\circ\circ\circ\underbrace{\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}}_{fit\ \sharp 1\ and\ \sharp 2}\circ\circ\circ\circ\circ\circ\\
&3:\ \circ\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\underbrace{\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\circ}_{fit\ \sharp 2\ and\ \sharp 3}\\
&4:\ \circ\circ\circ\circ\circ\circ\underbrace{\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}}_{fit\ \sharp 1\ and\ \sharp 2}\circ\circ\circ\circ\circ\circ\\
&5:\ \underbrace{\circ\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\color{cyan}{\circ}}_{fit\ \sharp 0\ and\ \sharp 1}\circ\circ\circ\circ\circ\color{cyan}{\circ}\circ\circ\circ\circ\circ\circ\\
\end{aligned}
```
"fit #0 and #1" では ${\bf p}_1,{\bf p}_2,{\bf p}_3,{\alpha}_0,{\alpha}_2$ を更新し、${\bf p}_0,{\bf q}_3$ の位置ベクトルと ${\bf q}_2-{\bf q}_3$ の方向ベクトルは固定する。
"fit #1 and #2" では $\beta_0,{\bf p}_2,{\bf p}_3,{\alpha}_0,{\alpha}_2$ を更新し、${\bf p}_0,{\bf q}_3$ の位置ベクトルと ${\bf p}_1-{\bf p}_0,{\bf q}_2-{\bf q}_3$ の方向ベクトルは固定する。
"fit #2 and #3" は点列を逆転させれば、fit #0 and #1 と同様の処理である。
"fit #0 and #1"、"fit #1 and #2" の最適化時の目的関数は凸ではないので、「準ニュートン法」などで解くのがよい。

${\alpha_0},{\alpha_2},{\beta_0}$ は $||{\bf q}_1 - {\bf q}_0||,||{\bf q}_2 - {\bf q}_3||,||{\bf p}_1 - {\bf p}_0||$ のことである。

### 追記案：誤差を指定せずに点列の曲率半径の絶対値の極小値で点列を分割する
点列にノイズがのっている場合も考えると、FFTで周波数領域に変換して高周波の振幅をカットして元に戻してノイズを削除した後に、曲率半径を計算するとよい。
ノイズのカットはFFT以外でもよい。ガウシアンフィルタとか。
ちなみにこれは試していない。備忘録。
入力する誤差に依存しないアルゴリズムなので、いいんじゃないかな？

### 追記：点列が閉ループの場合
${\bf p}_0={\bf q}_3$ とし、${\bf p}_0$ は点列の始点と一致する必要は無くなる。
つまり、${\bf p}_0$ も最適化の対象となる。

### 最後に
<a href="https://qiita.com/BlueRayi/items/7965798ba1127d269ebb">こちら</a>の記事にあるようにQiitaは変な仕様で書いていてイライラする。