Total Variation Minimization

TV

LDPCはしばらくお休み（基本的に優先順位が低い）．今日はTV（Total Variation）の実装について書く．

神奈川大の斉藤先生が国内では精力的に研究していらっしゃる．どうにも魅力的なツールに見えるのだが，意外と名前負けしている（すなわち，ネーミングセンスが抜群）感が拭えない．曰く，画像を骨格画像とテクスチャ画像に分離する手法である．なんとか実装が出来たので，公開しておく．ただし，まだ国際的な最前線はおろか，斉藤先生にも追いついていないことは明記しておく．

まず，TVとは何か，これを使うと何が出来るかについてざっくり話す．さらに，文献を挙げておく．最後に実装の方針とソースコードを載せる．

TVとは，不連続を許す空間におけるノルムのはかり方です（多分）．不連続だから微分が出来ないところがあります．逆に，値が連続的でない場合も扱えます．実は，画像データはエッジが含まれているので，不連続を許す空間と相性が良さそうなのです．

実態としては，微分の絶対値和です．さっき微分できないと言いましたが，積分して微分は出来て，すなわち微分して積分は出来るのです（超おおざっぱ）．詳しくは参考文献とプログラムを見てください．

さて，このTVと相対を為すものは何でしょうか？TVは単調増加や単調減少の時にノルムが小さく，激しく増減をしているとノルムが大きくなります．この逆，つまり相対を為すのは，変化が少ないとノルムが大きく，振動が大きいほどノルムが小さくなれば言いわけです．これをGノルムと名付けましょう．概念としてはこれで良いのですが，数学的に定義するのは大変そうです．これも参考文献を見てください．

ここで，TVノルムとGノルムの重み付け和を最小にする，という課題を考えます．すると，画像をTVノルムの属する空間とGノルムの属する空間に分離できます．すなわち，骨格画像とテクスチャ画像に分離できちゃうのです．それも，数学に裏打ちされた結果として．

最後に歴史的なところを適当にかいつまんで列挙しておきます．そもそもの発想はWeveletから来ています．Mayerさんによる，Waveletの高周波成分がテクスチャ画像にちがいない，というアレゲなところから始まっています．それから，Rudinさん，Osherさん，Fatemiさんが頑張りまして，ROFモデルというのを提案しました．さらに，VesaさんとOsherさんが頑張りました．そのころにはWaveletはなりを潜めています．この辺の順序はちょっと未確認です（発表文献と年数を追っていけばいいですけど）．しばらくして，（Waveletを使わない）真っ当な数値計算法が提案されてきました．でもまだ，オイラー・ラグランジュ法を使っている「古典的」な手法でした．例えば，0除算を避けるイプシロンが必要でした（私はこのモデルを今年の春頃実装してました．一応完成はしたけど収束が遅かった）．彗星のごとくChambolleさんが登場し，非常にスマートにこの計算方法を確立してしまいました．PDE（Partial Differential Equation）なる理論を用います（が，私にはさっぱり理解できません）．どうも，凸空間の非線形問題で，解が収束したり，不動点定理を利用したりと，お祭り騒ぎのようです（が，これもよく分からない．私は動けばいい，使えればいい，ぐらいなんで・・・）．満を持して，Aujolさん登場．Aujolさん，他3名によるA2BCモデルが提案され，Chambolleさんの解法でスマートに計算される，というのが現在までのところです．

必須文献はこれだ！

どれも大学関係者じゃないと入手は難しいかも．しかし，ここにも（ほとんど）同じ内容がある！

というわけで，これさえ読めばあっという間に追いつけるわけです．

次に，実装に移りましょう．2番目と4番目に詳しく書いてあるので，それで良いと思います．今日公開するのはROFモデルだけです．A2BCはまた日を改めて実装したいと思います．方針ですが，C言語で純関数言語的に実装してあります．確保したメモリは呼び出し側で解放してください．まぁ，メモリも潤沢ですので，まじめに解放してませんが．

array_r1_tは画像と同じ構造です．array_r2_tはベクトル画像の構造です．なるべく構造を隠蔽するために，el関数とvec1, vec2関数を使って要素にアクセスします（どちらもdefine文ですが）．イテレーションはメイン関数の中で実行しています．そして，ここが縮小写像を計算しているところになります．回数は文献によると30回〜40回ぐらいらしいです．

画像データは0〜1に正規化しています．パラメータのτは1/4が最大にして最速で収束する値らしいです．これは固定．λは超くせ者です．0〜1に正規化しているので，これの数分の1が良いようです．とりあえず，16/256=1/16にしてあります．ざっくりいって，輝度が1/16ぐらいで振動しているのがテクスチャに分類されます．逆に言うと輝度が1/16以上の不連続な部分は，エッジとして保存されます．

プログラムを投下します．

とりあえず，大急ぎでTVとその実装を見てきました．手元で動かすといろいろ問題点が出てきますな（笑）．続きはまた時間が出来たときに（1週間ぐらいこの話に費やしてしまったよ・・・パラメータが分からなくて）．