凹みTips

C++、JavaScript、Unity、ガジェット等の Tips について雑多に書いています。

Unity でボリュームレンダリングをしてみる - vol.4 プロシージャル

はじめに

本記事はボリュームレンダリングシリーズで以下の記事の続きです(2年ぶり)。

tips.hecomi.com

前回までは魚の 3D テクスチャを使いましたが、今回はノイズ関数を使って雲のようなものを作ってみたいと思います。前回同様キューブポリゴン内に生成する形を取り、ノイズ関数によるプロシージャルな密度関数をトレースするところを解説します。

参考

雲の生成は以下の Shadertoy の方法を参考に実装します。

www.shadertoy.com

手法ついては edo_m18さんが解説記事を書かれており、補足コメント付きの Shadertoy でのサンプルも公開されています。

qiita.com

www.shadertoy.com

これらは以下の Shader Bits さんの記事を参考に作られています。

shaderbits.com

この記事に関しても edo_m18さんが日本語へ翻訳されています。

qiita.com

ありがたいですね。。これらを色々改変しつつ Unity でのサンプルを作ってみようと思います。ただ今回は Shadertoy のものを再現するまでで、他のオブジェクトとの整合性やライティングの一貫性までは手を付けません。それっぽく見えるのを目標にします。

コード

7. Procedural が今回のサンプルになります。

github.com

密度関数

3D テクスチャ部分を式にしてみる感覚を掴むのがまず早いと思いますので、簡単な球体の中心に行くほど密度が高くなるようなケースを考えてみましょう。中心に行くほど 1.0 を返すような直径 1 の球の密度関数は次のように書くことが出来ます。

inline float densityFunction(float3 p)
{   
    return 0.5 - length(p);
}

要はレイマーチングでの数式(SDF)の符号を反転させた形にすればよいわけですね。

これを前回までと同じようにキューブポリゴンのオブジェクト空間でレイを通過させるコードは以下のようになります。

Shader "VolumeSphereUnlit"
{

Properties
{
    _Color("Color", Color) = (1, 1, 1, 1)
    _Intensity("Intensity", Range(0, 1)) = 0.1
    [IntRange] _Loop("Loop", Range(0, 128)) = 32
}

CGINCLUDE

#include "UnityCG.cginc"

struct appdata
{
    float4 vertex : POSITION;
};

struct v2f
{
    float4 vertex   : SV_POSITION;
    float3 worldPos : TEXCOORD1;
};

#define MAX_LOOP 100
float4 _Color;
float _Intensity;
int _Loop;

inline float densityFunction(float3 p)
{   
    return 0.5 - length(p);
}

v2f vert(appdata v)
{
    v2f o;
    o.vertex = UnityObjectToClipPos(v.vertex);
    // ポリゴン表面の座標がフラグメントシェーダで使えるようにする
    o.worldPos = mul(unity_ObjectToWorld, v.vertex);
    return o;
}

float4 frag(v2f i) : SV_Target
{
    // ワールド空間でのポリゴン表面座標とそこへのカメラからの向き
    float3 worldPos = i.worldPos;
    float3 worldDir = normalize(worldPos - _WorldSpaceCameraPos);

    // オブジェクト空間に変換
    float3 localPos = mul(unity_WorldToObject, float4(worldPos, 1.0));
    float3 localDir = UnityWorldToObjectDir(worldDir);

    // オブジェクト空間でのレイのステップ長
    float step = 1.0 / _Loop;
    float3 localStep = localDir * step;

    // レイを通過させて得られる透過率
    float alpha = 0.0;

    for (int i = 0; i < _Loop; ++i)
    {
        // ポリゴン中心ほど大きな値が返ってくる
        float density = densityFunction(localPos);

        // 球の外側ではマイナスの値が返ってくるのでそれを弾く
        if (density > 0.001)
        {
            // 透過率の足し合わせ
            alpha += (1.0 - alpha) * density * _Intensity;
        }

        // ステップを進める
        localPos += localStep;

        // ポリゴンの外に出たら終わり
        if (!all(max(0.5 - abs(localPos), 0.0))) break;
    }

    float4 color = _Color;
    color.a *= alpha;
    return color;
}

ENDCG

SubShader
{

Tags 
{ 
    "Queue" = "Transparent"
    "RenderType" = "Transparent" 
}

Pass
{
    Cull Back
    ZWrite Off
    Blend SrcAlpha OneMinusSrcAlpha 
    Lighting Off

    CGPROGRAM
    #pragma vertex vert
    #pragma fragment frag
    ENDCG
}

}

}

これを表示すると次のような形になります。

f:id:hecomi:20200425203517p:plain

ノイズ関数

さて、シンプルな球の密度関数を今度は fBm で置き換えてみます。

thebookofshaders.com

簡単に言えば、今の球の密度関数とノイズ関数を組み合わせて雲っぽい見た目を作ります。fBm は以下の Shadertoy のコードを利用させてもらいました。

www.shadertoy.com

Properties
{
    ...
    _NoiseScale("NoiseScale", Range(0, 100)) = 5
    _Radius("Radius", Range(0, 2)) = 1.0 
}

...

float _NoiseScale;
float _Radius;

// ref. https://www.shadertoy.com/view/lss3zr
inline float hash(float n)
{
    return frac(sin(n) * 43758.5453);
}

inline float noise(float3 x)
{
    float3 p = floor(x);
    float3 f = frac(x);
    f = f * f * (3.0 - 2.0 * f);
    float n = p.x + p.y * 57.0 + 113.0 * p.z;
    float res = 
        lerp(lerp(lerp(hash(n +   0.0), hash(n +   1.0), f.x),
                  lerp(hash(n +  57.0), hash(n +  58.0), f.x), f.y),
             lerp(lerp(hash(n + 113.0), hash(n + 114.0), f.x),
                  lerp(hash(n + 170.0), hash(n + 171.0), f.x), f.y), f.z);
    return res;
}

inline float fbm(float3 p)
{
    float3x3 m = float3x3(
        +0.00, +0.80, +0.60,
        -0.80, +0.36, -0.48,
        -0.60, -0.48, +0.64);
    float f = 0.0;
    f += 0.5 * noise(p); p = mul(m, p) * 2.02;
    f += 0.3 * noise(p); p = mul(m, p) * 2.03;
    f += 0.2 * noise(p);
    return f;
}

inline float densityFunction(float3 p)
{   
    return fbm(p * _NoiseScale) - length(p / _Radius);
}

...

fBm のコードは複雑ですが要は色々な周波数のノイズを組み合わせたノイズを生成しています。密度関数では fBM を返すだけでは一様に濃い場になってしまうので、先ほどと同じように遠くなるほど密度が小さくなっていくように length(p) を返しています。ただし、大きさが制御できるように _Radius を導入しています(値が小さいほどすぐにこの第 2 項が小さくなるのをイメージしてください)。

f:id:hecomi:20200425214822g:plain

アーティファクトの除去

よくよく画を見てみると繰り返しのアーティファクトが見えてしまいます。

f:id:hecomi:20200425220835p:plain

これはレイの開始地点がボックスの表面なので各面の繰り返しパターンがステップ長ごとに見えてしまうからです。これを防ぐには冒頭にも貼りました Shader Bits さんで解説があるように、レイの開始地点を少しボックス内部に移動させ、カメラから見たプレーンが開始地点となるようにします。

float4 frag(v2f i) : SV_Target
{
    float step = 1.0 / _Loop;

    float3 worldPos = i.worldPos;
    float3 worldDir = normalize(worldPos - _WorldSpaceCameraPos);
    float3 camToWorldPos = worldPos - _WorldSpaceCameraPos;
    worldPos += (step - fmod(length(camToWorldPos), step)) * worldDir;
    ...

f:id:hecomi:20200425220815p:plain

または開始点をジッターを掛けて位置・時間でずらしてしまっても良いです。

float4 frag(v2f i) : SV_Target
{
    ...
    float jitter = step * hash(
        worldPos.x + 
        worldPos.y * 10 + 
        worldPos.z * 100 + 
        _Time.x);
    worldPos += jitter * worldDir;
    ...

プレーンを揃える方式ではポリゴンの境界面でまだ破綻があるため、今回はジッターを採用します。

透過吸収

さて、このままでは見た目に陰影が無いので、雲の透過吸収の物理挙動に沿うように密度の反映の仕方を変えてみます。alpha の代わりに透過率 transmittance および光がどれだけ吸収されるのかの係数 _Absorption、そして雲の透明度を表す _Opacity を導入して次のように計算を行います。

Properties
{
    ...
    _Absorption("Absorption", Range(0, 1000)) = 50
    _Opacity("Opacity", Range(0, 100)) = 100
    ...
}

...
float _Absorption;
float _Opacity;
...

float4 frag(v2f i) : SV_Target
{
    ...
    float4 color = float4(_Color.rgb, 0.0);
    float transmittance = 1.0;
    ...

    for (int i = 0; i < _Loop; ++i)
    {
        float density = densityFunction(localPos);

        if (density > 0.01)
        {
            float d = density * step;
            transmittance *= 1.0 - d * _Absorption;
            if (transmittance < 0.01) break;

            color.a +=  _Opacity * d * transmittance;
        }
        ...
    }

    return color;
}

ENDCG

f:id:hecomi:20200425232540g:plain

見た目がそれっぽくなりました。

ライティング

ライティングは Shader Bits さんに書かれているように重いですが柔軟性のある方法を取ります。現在のレイの点から光源方向にレイを伸ばし、それぞれの点で密度を計算、それを見て光方向の透過具合を調べる、というものです。ループの中で更にループでレイを伸ばすので重いのは想像がつくと思いますが、コード自体はとても簡単になります。

Properties
{
    ...
    _AbsorptionLight("AbsorptionLight", Range(0, 100)) = 50
    _OpacityLight("OpacityLight", Range(0, 100)) = 50
    _LightStepScale("LightStepScale", Range(0, 1)) = 0.5
    [IntRange] _LoopLight("LoopLight", Range(0, 128)) = 6
}

...
float _AbsorptionLight;
float _OpacityLight;
int _LoopLight;
float _LightStepScale;
float4 _LightColor0;
...

float4 frag(v2f i) : SV_Target
{
    ...
    float lightStep = 1.0 / _LoopLight;
    float3 localLightDir = UnityWorldToObjectDir(_WorldSpaceLightPos0.xyz);
    float3 localLightStep = localLightDir * lightStep * _LightStepScale;

    float4 color = float4(_Color.rgb, 0.0);
    float transmittance = 1.0;

    for (int i = 0; i < _Loop; ++i)
    {
        float density = densityFunction(localPos);
        if (density > 0.0)
        {
            ...
            float transmittanceLight = 1.0;
            float3 lightPos = localPos;

            for (int j = 0; j < _LoopLight; ++j)
            {
                float densityLight = densityFunction(lightPos);
                if (densityLight > 0.0)
                {
                    float dl = densityLight * lightStep;
                    transmittanceLight *= 1.0 - dl * _AbsorptionLight;
                    if (transmittanceLight < 0.01) 
                    {
                        transmittanceLight = 0.0;
                        break;
                    }
                }
                lightPos += localLightStep;
            }
            ...
            color.rgb += 
                _LightColor0 * 
                (_OpacityLight * d * transmittance * transmittanceLight);
        }
        ...
    }

    return color;
}

f:id:hecomi:20200426183713g:plain

アニメーション

もう一度、密度関数 densityFunction() に注目してみます。最初の方にこれはレイマーチングの距離関数の符号反転したものだと述べました。ということはレイマーチングと同じように距離関数の合成などの操作も可能、ということになります。おさらいは以下の記事をご参照ください:

tips.hecomi.com

www.iquilezles.org

ではこの関数を次のようにしてみます。

inline float sphere(float3 pos, float radius)
{
    return length(pos) - radius;
}

inline float torus(float3 pos, float2 radius)
{
    float2 r = float2(length(pos.xz) - radius.x, pos.y);
    return length(r) - radius.y;
}

inline float densityFunction(float3 p)
{   
    float f = fbm(p * _NoiseScale);

    // 半径 0.2 の急に 0.3 くらいのノイズを混ぜたもの
    float d1 = -sphere(p, 0.2) + f * 0.3;

    // 内径 0.36、外径 0.46 のトーラスに 0.2 くらいのノイズを混ぜたもの
    float d2 = -torus(p, float2(0.36, 0.1)) + f * 0.2;

    // 2 つのブレンド比率を 0 ~ 1 でなめらかに
    float blend = 0.5 + 0.5 * sin(_Time.z);

    // 2 つの密度関数を合成
    return lerp(d1, d2, blend);
}

f:id:hecomi:20200426184806g:plain

これでポリゴンの位置移動・回転・スケールも効くアニメーション付き雲が出来上がりました。なおこのトーラス雲は以下の かねたさんのデモからヒントを得たものになります。

www.shadertoy.com

おわりに

オブジェクトとのインターセクションやライティングの変更、最適化周りについて続けて書きたいと思います。