Houdiniでアルゴリズミック・デザイン～Isosurface編

# Houdiniでアルゴリズミック・デザイン～Isosurface編

この記事はHoudiniアドベントカレンダー2017　10日目の記事です。


*本文最下部に記事で使っているHoudiniファイルのダウンロード先のリンクがあります*

## Isosurfaceとは

**[Isosurface](https://en.wikipedia.org/wiki/Isosurface)** は建築界隈ではおなじみのワードで、空間上に数値が与えられた点群があると考えた時（**[ボクセル空間](https://ja.wikipedia.org/wiki/%E3%83%9C%E3%82%AF%E3%82%BB%E3%83%AB)**）、その点群上の近い数値同士を結んだ時にできる面のことです。

![Isosurface.png](https://qiita-image-store.s3.amazonaws.com/0/42909/9be775c4-4f68-c655-effe-a658cef6d2e4.png)
* Pinterest上のisosurfaceで検索すると出て来る画像

このボクセル空間から面をを作るためのアルゴリズムとして**[Marching Cube](https://en.wikipedia.org/wiki/Marching_cubes)** というものがありますが、ありがたいことに **[Houdini](https://www.sidefx.com/)** を使うことで簡単にこのボクセル空間が作れ、またそのボクセル空間からMarching Cubeと思われるアルゴリズムを使ったノードでメッシュを作成することが可能です。

この記事ではその方法を説明したいと思います。

## 数理的にボクセル空間をコントロールするための関数を集める

なにはともあれボクセル空間上の点に数値を与えていくことが必要になりますが、今回はお題として **[アルゴリズミック・デザイン](http://code.algorithmicdesign.net/)** とよんでいることもあり、ボクセル空間上の点の位置情報（XYZ座標）を利用してsineやcosineといった関数を使って数理的に数値を決定してみたいと思います。

ただ闇雲に関数を作ってもきれいな面ができるとは限らないので、とりあえず先例はないかなと探します。ざっと検索した結果、すぐにでも使えそうな関数が載っているページが２つ見つかったので、今回はそちらに載っている関数を拝借して、一部 **[パラメトリック](https://en.wikipedia.org/wiki/Parametric_design)** にすることを前提に利用してみたいと思います。

- [Metamagic "3D Printing with Isosurfaces"](https://metamagics.wordpress.com/2012/10/16/3d-printing-part-1/)
- [stackoverflow "How do I reproduce this heart-shaped mesh in MATLAB?"](https://stackoverflow.com/questions/1526898/how-do-i-reproduce-this-heart-shaped-mesh-in-matlab)

## Houdiniのノードのセットアップ

ボクセル空間をコントロールするための関数が集まったところで、次にHoudiniのノードのセットアップです。今回は関数を入力することがメインなので、ノードのネットワーク自体は簡素なものにしたいと思います。

### Step 1: Boxを作成

![step-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/9ba1ecec-2907-65ec-9a7d-508dfd96dcf8.png)
ボクセル空間の領域となるBoxを作成します。
（Boxノード：Uniform Scale=5)

### Step 2: BoxからFog Volumeを作成

![step-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/d9175d42-ca88-015a-1e7e-cc3b8ddbe906.png)
Houdini上のボクセル形式の一つであるFog Volumeをボックスから作成します。
(Isooffsetノード：Uniform Sampling Divs=100)

### Step 3: 空のVolume Wrangleを６個つなげる
![step-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/8c6fb6b6-0a7a-f757-8955-1e63e8c6156e.png)
関数の参照先の名称にのっとり、Volume Wrangleノードの名前をそれぞれ"sphere", "wonder_tree", "space_egg", "orbis_gravis", "invader", "heart"にします。最後にその６つのVolume WrangleをSwitchノードで切り替えられるようにします。

### Step 4: FogをSDFに変換する
![step-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/92434413-d60c-c244-9f45-b0b473c90ea5.png)
メッシュに変換するために、Fogを一度Houdini上のもう一つのボクセル形式であるSDFに変換します。
(Convert VDBノード：Convert to=VDB, VDB Class=Convert Fog to SDF)

### Step 5: SDFをPolygonに変換する
![step-5.png](https://qiita-image-store.s3.amazonaws.com/0/42909/62ecd7ff-36d5-5202-ee55-0a5d5b1c6671.png)
これでボクセル空間がメッシュに変換されました。
(Convert VDBノード：Convert to: Polygons)

### Step 6: 色を付ける
![step-6.png](https://qiita-image-store.s3.amazonaws.com/0/42909/f20fd998-f409-9772-d3c4-4cb16ed549c1.png)
ポイントの位置情報を利用して色をつけます。
(Colorノード：Color Type=Ramp from Attribute, Range= -2.5, 2.5)

## 関数をVolume Wrangleに入力

これで準備はできたので、早速６つあるそれぞれのVolume Wrangleに参照先URLから拝借してきた関数を入力してみようと思います。

### Sphere

まずはMetamagicで一個目で説明しているシンプルな球体から。サイトでは下記のような関数で、半径４の球が表現できるとしています。

>The application lets you visualize an implicit, 3D surface, i.e., the solution to an equation in three variables. A sphere centered on origin with a radius of four, for example, would be the solution to the equation
> `x^2+y^2+z^2-4^2`

ではこれをVolume Wrangleをつかって、この関数を導入してみたいと思います。

![sphere-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/0269f2bc-992e-d301-3594-4e7466814514.png)


``` c
float radius = chf("radius");

float x = @P.x;
float y = @P.y;
float z = @P.z;

float base = pow(x,2) + pow(y,2) + pow(z,2) - pow(radius,2);
@density = fit(base,-1,1,1,0);

```

ここでは、radiusというfloatのパラメータを作成して、半径を決定していると考えられる4.0の部分をradiusで置き換えています。結果としてパラメトリックに半径を変えられる球体のフォグができました。
これをメッシュ化すると・・・

![sphere-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/4243b4ee-8ac8-d468-c111-8d5f95b7e83c.png)
(radius=2.5)

![sphere-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/45b26f36-4b6f-591d-d97d-5a8d0bfbd37f.png)
(radius=1)

![sphere-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/67fc25b3-7d53-9a37-68ac-d7a3050cd0b2.png)
(radius=3.5)


球じゃあまりおもしろくないですね。


### Wonder Tree

同じMetamagicで説明しているWonder Treeと呼んでいる関数は下記のとおり。

> Next, try the formula

> ``` c
  cos(4.0*x/(x*x+y*y+z*z+0.0001))*sin(4.0*y/(x*x+y*y+z*z+0.0001))
+ cos(4.0*y/(x*x+y*y+z*z+0.0001))*sin(4.0*z/(x*x+y*y+z*z+0.0001))
+ cos(4.0*z/(x*x+y*y+z*z+0.0001))*sin(4.0*x/(x*x+y*y+z*z+0.0001))
+ exp(0.1*(x*x+y*y+z*z-0.2)) 
- exp(-10.0*(x*x+y*y+z*z-0.15))
```
> inside the bounds -2 to 2.

これをHoudini用に変換すると・・・

![wonder_tree-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/7478f9b0-79ac-db23-876d-4842df5810ea.png)

``` c
float x = @P.x;
float y = @P.y;
float z = @P.z;

float parm = chf("parm");

float base = cos(parm*x/(x*x+y*y+z*z+0.0001))*sin(parm*y/(x*x+y*y+z*z+0.0001))
+ cos(parm*y/(x*x+y*y+z*z+0.0001))*sin(parm*z/(x*x+y*y+z*z+0.0001))
+ cos(parm*z/(x*x+y*y+z*z+0.0001))*sin(parm*x/(x*x+y*y+z*z+0.0001))
+ exp(0.1*(x*x+y*y+z*z-0.2)) 
- exp(-10.0 *(x*x+y*y+z*z-0.15));

@density = fit(base,-1,1,1,0);
```

ここではparmというfloatパラメータを作り、模様を変形するのに最適と感じた部分（4.0の部分）をparmで置き換えています。メッシュ化した結果は・・・

![wonder_tree-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/8684c6d6-931d-ff78-4562-3f21f8e4ba08.png)
(parm=4.0)

![wonder_tree-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/417c1fed-8b6a-6667-d02e-4d67ca6228ec.png)
(parm=3.16)

![wonder_tree-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/62674ead-736a-6547-e711-c3e8b9d8657d.png)
(parm=9.0)

気持ち悪い！いい感じ！

### Space Egg

同じMetamagicで説明しているSpace Eggと呼んでいる関数は下記のとおり。

> XYZ range: [-8.5, 8.5] x [-8.5, 8.5] x [-17, 17].

> ``` c
-(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
*(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
+0.05-exp(100.0*(x*x/64+y*y/64 + z*z/(1.6*64)*exp(-0.4*z/8) - 1))
```

Houdiniで書くと・・・

![space_egg-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/7d444200-85dd-6e16-a451-38a24fbfde35.png)

``` c
float size = chf("size");
float thickness = chf("thickness");

float x = @P.x*size;
float y = @P.y*size;
float z = @P.z*size;

float base = -(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
*(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
+thickness-exp(100.0*(x*x/64+y*y/64 + z*z/(1.6*64)*exp(-0.4*z/8) - 1));

@density = fit(base,-1,1,0,1);
```

ここではsizeとthicknessという２つのパラメータを設定しています。これをメッシュ化すると・・・

![space_egg-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/966a1dd7-52be-ec3b-05d7-d556688d6ab9.png)
(size=5, thickness=0.2)

![space_egg-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/dbbb7366-9d96-7859-f9e0-7f0369bbd816.png)
(size=1.5, thickness=0.1)

![space_egg-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/b02c96d3-924f-6ee8-82ef-c9166f872c95.png)
(size=3, thickness=0.8)

エイリアン感あっていい！

### Orbis Gravis

同じMetamagicで説明しているOrbis Gravisと呼んでいる関数は下記のとおり。

> XYZ range: [-6, 6] x [-6, 6] x [-6, 6].

> ``` c
-(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
*(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
+0.02
+exp((x*x+y*y+z*z)-32.0)
```

Houdiniに翻訳・・・

![gravis_orbis-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/5032a891-8568-588c-2760-425c62ef591e.png)

``` c
float size = chf("size");
float parm = chf("parm");

float x = @P.x*size;
float y = @P.y*size;
float z = @P.z*size;

float base = -(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
*(cos(x) * sin(y) + cos(y) * sin(z) + cos(z) * sin(x))
+0.02
+exp((x*x+y*y+z*z)-parm);

@density = fit(base,-1,1,1,0);
```

パラメータとしてサイズをコントロールするsizeと模様を変えるparm（他に言い方思いつかなかった）の二種類を作っています。

これをメッシュ化すると・・・

![gravis_orbis-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/851eb5d0-ac65-6a69-c7c4-a0a37664f22b.png)
(size=2.5, parm=32)

![gravis_orbis-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/ee61cb9a-cfd7-9645-c3a1-82279b8e7f05.png)
(size=1.5, parm=12)

![gravis_orbis-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/85fa0af4-639d-14ac-a021-b63de08481e3.png)
(size=3, parm=100)

未来球技のボールデザイン！


### Invader

同じMetamagicで説明しているInvaderと呼んでいる関数は下記のとおり。

>XYZ range: [-1, 1] x [-1, 1] x [-1, 1].

> ``` c
  cos(4.0*x/(x*x+y*y+z*z+0.2)) * sin(4.0*y/(x*x+y*y+z*z+0.2)) 
+ cos(4.0*y/(x*x+y*y+z*z+0.2)) * sin(4.0*z/(x*x+y*y+z*z+0.2)) 
+ cos(4.0*z/(x*x+y*y+z*z+0.2)) * sin(4.0*x/(x*x+y*y+z*z+0.2))
+ exp(0.5*(x*x+y*y+z*z - 0.1))
```

Houdiniで書き直すと・・・

![invader-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/3cd41589-863e-56c3-126d-b2b6c743e19a.png)

``` c
float size = chf("size");
float parm = chf("parm");

float x = @P.x*size;
float y = @P.y*size;
float z = @P.z*size;



float base = cos(parm*x/(x*x+y*y+z*z+0.2)) * sin(parm*y/(x*x+y*y+z*z+0.2)) 
+ cos(parm*y/(x*x+y*y+z*z+0.2)) * sin(parm*z/(x*x+y*y+z*z+0.2)) 
+ cos(parm*z/(x*x+y*y+z*z+0.2)) * sin(parm*x/(x*x+y*y+z*z+0.2))
+ exp(0.5*(x*x+y*y+z*z - 0.1));

@density = fit(base,-1,1,1,0);
```

ここでもsize（サイズコントロール）とparm（模様コントロール）という二種類のパラメータを作成しています。パラメータを変えてメッシュ化すると・・・

![invader-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/81d5e74f-cca6-17c5-a71d-a2a80209cbb6.png)
(size=0.35, parm=4)

![invader-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/e48c988c-89f8-4c5e-1ffa-a08272fb6eb5.png)
(size=0.3, parm=5)

![invader-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/f9ec6641-8e40-5f7a-3952-d2a8eca39b7a.png)
(size=0.15, parm=10)

宇宙に咲く花といった感じ。かっこいい！



### Heart

最後にハートマークを作ります。stackoverflowで解説されていたハートマークを作れる関数は下記のとおり。MATLABを使っているみたいです。

> ``` matlab
% volume data
step = 0.05;
[X,Y,Z] = meshgrid(-3:step:3, -3:step:3, -3:step:3);
F = (-(X.^2).*(Z.^3)-(9/80).*(Y.^2).*(Z.^3))+((X.^2)+(9/4).*(Y.^2)+(Z.^2)-1).^3;
```

これもHoudiniで書いてみると・・・

![heart-1.png](https://qiita-image-store.s3.amazonaws.com/0/42909/1a613d40-8e0f-a6ec-3387-8da98c79b624.png)

``` c
float size = chf("size");
float parm = chf("parm");

float x = @P.x*size;
float y = @P.y*size;
float z = @P.z*size;



float base = (-pow(x,2)*pow(z,3)-(parm/80.0)*pow(y,2)*pow(z,3))+pow((pow(x,2)+(9.0/4.0)*pow(y,2)+pow(z,2)-1),3);

@density = fit(base,-1,1,1,0);
```


ここでのパラメータはsize（サイズコントロール）とparm（模様のコントロール）。parm=9.0の時にハートの形になります。それ以外だと・・・

![heart-2.png](https://qiita-image-store.s3.amazonaws.com/0/42909/3982200b-a1df-dc13-1f98-49e3d2b04474.png)
(size=0.5, parm=9.0)

![heart-3.png](https://qiita-image-store.s3.amazonaws.com/0/42909/a36d87a0-1a3a-221d-7146-e9618d3f3948.png)
(size=0.5, parm=300.0)

![heart-4.png](https://qiita-image-store.s3.amazonaws.com/0/42909/a460c863-e1fc-a9b2-aa25-9528c0024886.png)
(size=0.7, parm=10000.0)

ハートを成長させるとこうなるという見方でみると面白い。


## おわりに

以上がHoudiniをツールとして用いて、数理的に**Isosurface**を作る方法となります。

数理的にボクセル空間をいじると簡単に非常に面白い造形ができるということは感じていただけたかと思います。この手法にプラスでマニュアルコントロール（ramp parameterとか、参照点・面との関係性によるコントロール、vdb combine...etc）を加えることで**もっと面白い造形**が望めると思います。

こんな操作を簡単にできるHoudini、**素敵**です。

## ファイルダウンロード先

**データの動作環境:** 
Houdini Apprentice 16.5.268

**ダウンロード先リンク(Github):**
https://github.com/jhorikawa/ProgrammingArticleFiles/tree/master/20171210_Algorithmic%20Design%20with%20Houdini%20-%20Isosurface