これは楽しい！簡単にマーブリングが作れるアルゴリズム

　ちょっとタイトルを盛ってみた。というわけで、前々から一度実装してみたかったマーブリングのシミュレーションアルゴリズムを昨日一日楽しんだのでここで紹介します。

　マーブリングってのは、皆さん（たぶん）ご存じの、水面にマーブリング用の絵具を垂らして、それを棒でかき混ぜたりして模様を作って、紙に転写するあれです。誰でも簡単に美しい模様ができる表現技法マーブリング

　マーブリングの画像を見て貰うと分かりますが、マーブリングでできる模様ってのは、色と色の境界がほとんどぼけてない、シャープな模様ができるんですよね。寺町美術教室 >> マーブリングに掲載されてる画像などが良い例だと思います。

　これを流体計算で再現するのは大変なんです。あ、私はこの一年流体計算について研究してきたので、そのうちどがっと流体計算についての記事を書きます。まぁ、それはよしとして、流体計算ってのは多くの場合1種類の流体の計算で、やっても2種類までなんですよ。それ以上は非常に大変と。なので、流体計算でマーブリングみたいなのをやろうとしたら、液体自体は1種類で、色の混合計算をするか、レイヤー事に分けて計算するかのどっちかになってしまうんです。（後者のはYouTubeでComputer-Generated Marbling Textures: A GPU-Based Design SystemというIEEE Computer Graphics and Applicationsに投稿された動画が見られます。）

　んじゃ、流体計算やらなくて、もっと別な手段でやればいいじゃない、というのがこの記事で紹介するMathematical MarblingというShufang Luらによる手法です。昨年IEEE Computer Graphics and Applicationsに投稿された論文で、読んで字のごとく数学的にマーブリング画像を生成します。この手法は流体計算なんかと違って、とにかく実装が簡単です。すぐに作って遊べます。そして数式が簡潔で美しいです。高校生でもちゃんと「なるほどねー」と理解できる数式なのはすばらしい。

　前置きが長くなりましたが、ざっくりと論文の内容に入っていきましょう。（どこまで書いて良いのかわかんないから、本当にざっくりと。興味を持ったら論文を購入してみてください。）

　Mathematical Marblingではマーブリング画像を生成を次の5つの手順の組み合わせで生成します

インクの滴下
直線的に絵具を引っ張る
円形に絵具を引っ張る
渦巻き状にする
波打たさせる

インクの滴下

　インクを滴下すると、インクは必ず半径rの円形に広がる物とします。
　周囲にある液体は、インクが滴下されたことによって外側に押し出されます。元々点Pにあった絵具が、インクが点Cに滴下されたことによってP'に移動したとします。この時、点P'は

$P'=C+(P-C)\sqrt{1+\frac{r^2}{|P-C|^2}}$ (1)

と表されます。（なんでこういう式になるかは論文を読もう！）
　で、これをP'の表現式からPの表現式に変換すると
$P=C+(P'-C)\sqrt{1-\frac{r^2}{|P'-C|^2}}$ (1_2)
という式になると思われます。（論文に何故か載せてくれていないんだな。私が計算したので間違いがあったらすみません。）

　(1_2)から移動後の位置から移動前の位置が分かるので、画像の変形が可能になります。

　実際に実装してみるとこんな感じになります。

直線的に絵具を引っ張る

　次に絵具を割り箸なんかでまっすぐ引っ張ったときの効果を作ります。
　この段階で

　こんな画像が生成できるようになります。楽しいね。ちなみに、この画像は3回引っ張ったよ。

　今回も点Pの位置にあった絵具が引っ張られて点P'に移動したとする。

　引っ張る直線lを、ある点AとベクトルLを用いて

$l(k) = A + kL$

という形で表現するとします。（kはただのスカラー変数)　意味としては点Aを通る方向Lに進む直線って事だね。
　で、このLの単位ベクトルをM、単位直交ベクトルをNとすると、P'は
$P'=P+\frac{\alpha\lambda}{d+\lambda}M$ 　(2)
と表される。ここでdは
$d=|(P-A)\cdot N|$
と表される。意味としては、直線Lと点Pの距離って事です。αとλはパラメータで、αが大きければ引っ張って移動する距離が大きくなり、λが大きくなると引っ張られる影響範囲が大きくなる。

　これをPの式にすると
$P=P'-\frac{\alpha\lambda}{|(P'-A)\cdot N|+\lambda}M$ 　(2_1)
となる（はず）。

円形に絵具を引っ張る

　直線に絵具を引っ張るだけでも十分楽しいんだけど、今度はぐるぐる回してみる。そしたらこんなのが作れるようになるよ〜。

　グールグール。

　さて、基本的な考え方はさっきの直線の時と同じ。（その辺の細かいことは論文読んでね）
　引っ張る軌跡の円の半径をr、中心をCとすると、P'は
$P'=C+(P-C)\left(\begin{array}{cc}\cos\theta & \sin\theta \\ -\sin\theta & \cos\theta \end{array}\right)$ 　(3)
となる。後ろの回転行列の式が普段見慣れてる式とマイナスの位置が逆だけど、これは論文の書き方に合わせました。ここで、θは
$\theta =\frac{\alpha\lambda}{(d+\lambda)|P-C|}$
で、dは
$d=||P-C|-r|$
です。これをPの式にすると
$P=C+(P'-C)\left(\begin{array}{cc}\cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end{array}\right)$ 　(3_1)
となる。θとdについてはPをP'に入れ替えれば良い。(はず)

渦巻きを作る

　なんでわざわざ分離してあるのか良くわかんないけど、さっきの円形に絵具を引っ張るにおいてr=0とする。以上

波打たせる

　こんな感じに絵具を変形させるよ。

　波の振幅Aと周波数ωと波の位相φと波打たせる方向tを用いて

$P'=P+A\sin(\omega P\cdot \left( \begin{array}{c}\sin t\\ -\cos t\end{array}\right) + \phi) \left( \begin{array}{c}\cos t\\ \sin t\end{array}\right)$ 　(4)
となる。そして、例のごとく
$P=P'-A\sin(\omega P'\cdot \left( \begin{array}{c}\sin t\\ -\cos t\end{array}\right) + \phi) \left( \begin{array}{c}\cos t\\ \sin t\end{array}\right)$ 　(4_1)
となる(と思う)。

　というわけで、駆け足でしたが、基本となるアルゴリズムを解説させていただきました。

　参考論文ではまだアンチエイリアスの方法や、ベクター画像、３Dテクスチャの話もありますが、ここでは省略させていただきます。

　私が作ったプログラムをgithubにあげておきましたので、気が向いたら弄ってみてください。
　Applicationクラスにmain関数があります。
　ちなみに、このプログラムで画像を保存しようとするとレンダリングの際に平気でCPU99%使用率を達成できます。しかも、スーパーサンプリングするのに真っ正直に全部計算してますので、メモリもすっごいたべます。貧弱なPCでは実行しないことをおすすめします。
（一応特許関連は調べてみましたが、たぶん申請していなさそうなので、載せても問題ないかと。問題があったら削除するので報告お願いします。あと、バグはもう修正しない(-_-)　）

作成したプログラムの著作権が私から研究室に移動したので公開できなくなりました。根性があれば、いつか0から作り直して公開………しないだろうなぁ。

　生成した画像