2017-11-11

Compute Shader でFFTと畳み込み演算でブラー

はじめに

前回,ComputeShaderを利用してSummed Area Tableを作りましたが, 今回はFFT(高速フーリエ変換)をしてみたいと思います。フーリエ変換ができれば,大きなカーネルの畳み込み計算を現実的な時間で出来て, Unreal Engineで言うところの,Convolution Bloomのようなエフェクトを実装できます。

https://docs.unrealengine.com/latest/INT/Engine/Rendering/PostProcessEffects/Bloom/

高速フーリエ変換

フーリエ変換の表式はいろいろありますが,今回は連続の空間だと以下のもので表されるとします。

${ \displaystyle \begin{eqnarray} F(k) = \int f(x) e^{-ikx} dx \end{eqnarray} }$

離散空間でも,いろいろ表式はありますが, ${N}$ 個のデータがあるとして

${ \displaystyle \begin{eqnarray} F(t) = \sum_{x=0}^{N-1} f(x) e^{-i\frac{2\pi x}{N}t} \end{eqnarray} }$

${N\times N}$ サイズの画像にフーリエ変換をかける場合は

${ \displaystyle \begin{eqnarray} F(s,t) = \sum_{x=0}^{N-1}\sum_{y=0}^{N-1} f(x,y) e^{-i\frac{2\pi x}{N}t} e^{-i\frac{2\pi y}{N}s} \end{eqnarray} }$

これはX方向に1次元のフーリエ変換を行ってから,Y方向にフーリエ変換を行えばよいことになるのですが,
何も考えずに計算をすると,1ピクセルあたり ${O(N^2)}$ の計算量となります。

これを,フーリエ変換の係数の重みの周期性に着目して加算の仕方を工夫したものがFFT(高速フーリエ変換)です。 1ピクセルあたりの計算量が ${O(N\log_2{(N)})}$ まで減ります。

https://ja.wikipedia.org/wiki/%E9%AB%98%E9%80%9F%E3%83%95%E3%83%BC%E3%83%AA%E3%82%A8%E5%A4%89%E6%8F%9B

FFTの実装

ComputeShaderでFFTを行う実装は,以下のサイトを参考にして実装しました。

https://software.intel.com/en-us/articles/fast-fourier-transform-for-image-processing-in-directx-11

コードはこんな感じです。

クリックで展開

// #define INVERSE      :  逆フーリエ変換
// #define PASS_1ST        :  1pass目, X方向変換
// #define PASS_2ND        :  2pass目, Y方向変換

#if defined(PASS_1ST) && !defined(INVERSE)
#define INPUT_REALPART     // 入力は実数空間のみ
#endif // #if defined(PASS_1ST) && !defined(INVERSE)
#if defined(PASS_2ND) && defined(INVERSE)
#define OUTPUT_REALPART    // 出力は実数空間のみ
#endif // #if defined(PASS_1ST) && !defined(INVERSE)

// 共有メモリが許すならば,ダブルバッファを行う
#if LENGTH<=512
#define DOUBLE_BUFFER
#endif // #ifLENGTH<=512

Texture2D<float3> g_SrcTextureR   : register(t0);
#if !defined(INPUT_REALPART)
Texture2D<float3> g_SrcTextureI   : register(t1);
#endif // #if !defined(INPUT_REALPART)
RWTexture2D<float3> g_DstTextureR : register(u0);
#if !defined(OUTPUT_REALPART)
RWTexture2D<float3> g_DstTextureI : register(u1);
#endif // #if !defined(OUTPUT_REALPART)

static const uint REAL_PART = 0;
static const uint IMAGINARY_PART = 1;

///< @brief バタフライ演算用の計算バッファ
#ifdef DOUBLE_BUFFER
groupshared float3 s_aaTempBuffer[2][2][LENGTH];
#define TempBuffer(tempId, x, ri) (s_aaTempBuffer[(tempId)][(ri)][(x)])
#else
groupshared float3 s_aaTempBuffer[2][LENGTH];
#define TempBuffer(tempId, x, ri) (s_aaTempBuffer[(ri)][(x)])
#endif

/// @brief バタフライ演算に必要な重みと加算インデックスを取得
void GetButterflyValues(uint passIndex, uint x, out uint2 indices, out float2 weights)
{
    int sectionWidth = 2 << passIndex;
    int halfSectionWidth = sectionWidth / 2;

    int sectionStartOffset = x & ~(sectionWidth - 1);
    int halfSectionOffset = x & (halfSectionWidth - 1);
    int sectionOffset = x & (sectionWidth - 1);

    sincos(2.0f*PI*sectionOffset / (float)sectionWidth, weights.y, weights.x);
    weights.y = -weights.y;

    indices.x = sectionStartOffset + halfSectionOffset;
    indices.y = sectionStartOffset + halfSectionOffset + halfSectionWidth;

    if (passIndex == 0)
    {
        indices = reversebits(indices) >> (32 - BUTTERFLY_COUNT) & (LENGTH - 1);
    }
}

/// @brief バタフライ演算
void ButterflyPass(int passIndex, uint x, uint tempId, out float3 resultR, out float3 resultI)
{
    // 重みとインデックスの取得
    uint2 indices;
    float2 weights;
    GetButterflyValues(passIndex, x, indices, weights);
    // 入力値の取得
    const float3 inputR1   = TempBuffer(tempId, indices.x, REAL_PART);
    const float3 inputI1   = TempBuffer(tempId, indices.x, IMAGINARY_PART);
    const float3 inputR2   = TempBuffer(tempId, indices.y, REAL_PART);
    const float3 inputI2   = TempBuffer(tempId, indices.y, IMAGINARY_PART);

    // ダブルバッファを持っていない場合は待機
#ifndef DOUBLE_BUFFER
    GroupMemoryBarrierWithGroupSync();
#endif // #ifndef DOUBLE_BUFFER

    // 加算
#ifdef INVERSE
    resultR = (inputR1 + weights.x * inputR2 + weights.y * inputI2) * 0.5f;
    resultI = (inputI1 - weights.y * inputR2 + weights.x * inputI2) * 0.5f;
#else
    resultR = inputR1 + weights.x * inputR2 - weights.y * inputI2;
    resultI = inputI1 + weights.y * inputR2 + weights.x * inputI2;
#endif
}

#ifdef OUTPUT_REALPART
/// @brief 逆フーリエ変換出力用のバタフライ演算パス
void InverseButterflyPassRealPart(int passIndex, uint x, uint tempId, out float3 resultR)
{
    // 重みとインデックスの取得
    uint2 indices;
    float2 weights;
    GetButterflyValues(passIndex, x, indices, weights);
    // 入力値の取得
    const float3 inputR1 = TempBuffer(tempId, indices.x, REAL_PART);
    const float3 inputR2 = TempBuffer(tempId, indices.y, REAL_PART);
    const float3 inputI2 = TempBuffer(tempId, indices.y, IMAGINARY_PART);
    // 加算
    resultR = (inputR1 + weights.x * inputR2 + weights.y * inputI2) * 0.5f; // // 逆フーリエ変換でしか利用しない
}
#endif // OUTPUT_REALPART


/// @brief CSエントリポイント
[numthreads(LENGTH, 1, 1)]
void CS_Transform(uint3 position : SV_DispatchThreadID)
{
    const uint bufferId = position.x;
    // フェッチするテクスチャ座標の取得
#ifdef PASS_1ST
    // X方向
    uint2 texturePos = position.xy;
#endif // #ifdef PASS_1ST
#ifdef PASS_2ND
    // Y方向
    uint2 texturePos = position.yx;
#endif // #ifdef PASS_2ND
    texturePos = (texturePos + LENGTH / 2) % LENGTH; // 中心を周波数空間原点とする処理

    // 入力テクスチャを格納
    TempBuffer(0, bufferId, REAL_PART) = g_SrcTextureR[texturePos];
#if defined(INPUT_REALPART)
    TempBuffer(0, bufferId, IMAGINARY_PART) = 0.0f; // FFT入力は実部のみ
#else 
    TempBuffer(0, bufferId, IMAGINARY_PART) = g_SrcTextureI[texturePos];
#endif // #if defined(INPUT_REALPART)

    // バタフライ演算
    // 最後の1回のパスの計算は行われない,(直接RWTextureに書き出すため)
    for (uint i = 0; i < BUTTERFLY_COUNT - 1; i++)
    {
        GroupMemoryBarrierWithGroupSync();
        const uint srcTempId = i % 2;
        const uint dstTempId = (i+1) % 2;
        ButterflyPass(i, bufferId, srcTempId, TempBuffer(dstTempId, bufferId, REAL_PART), TempBuffer(dstTempId, bufferId, IMAGINARY_PART));
    }

    // バタフライ演算最後の1回
    // 出力へ書き出し
    {
        GroupMemoryBarrierWithGroupSync();
        const uint srcTempId = (BUTTERFLY_COUNT - 1)% 2;
#ifdef OUTPUT_REALPART
        InverseButterflyPassRealPart(BUTTERFLY_COUNT - 1, bufferId, srcTempId, g_DstTextureR[texturePos]); // IFFT出力は実部のみ
#else
        ButterflyPass(BUTTERFLY_COUNT - 1, bufferId, srcTempId, g_DstTextureR[texturePos], g_DstTextureI[texturePos]);
#endif // #ifdef OUTPUT_REALPART
    }
}

前回のSATと同様に,共有メモリをうまく使って高速化をしています。インテルの実装だと,バタフライ演算でダブルバッファを利用しているので,共有メモリを2倍消費しますが,無駄なバリアが発生せずに高速です。
(何故バタフライ演算のインデックスの計算でreversebits命令を使うとうまくいくのか,あんまりよくわかりません。)

FFTの適用結果

512x512サイズの画像に対して,FFTを適用しました。 FFTを適用すると,実数成分の他に虚数成分も出てくるので,パワースペクトルとして ${|F(s,t)|^2}$ を出力しています。

元画像	FFT適用後パワースペクトル

512x512サイズのFFTのGPU実行時間は,GTX970の環境で0.32ms程度でした。

畳み込み演算

次に畳み込み演算について考えてみます。畳み込み演算の表式は

${ \displaystyle \begin{eqnarray} (f \ast g)(x) = \int f(x) g(x-t) dt \end{eqnarray} }$

ブルーム表現などで行っていたガウスブラーなどの処理は, この畳み込み演算を離散的に行っていたという事になります。

ところで,畳み込み演算に対してフーリエ変換を行うと,以下のように畳み込む関数をフーリエ変換した関数の積になります。

${ \displaystyle \begin{eqnarray} \int (f \ast g)(x) e^{-ikx} dx &=& \int dx \int dt f(x) g(x-t) e^{-ikx} \\\ &=& \int f(x) e^{-ikx} dx \int g(x) e^{-ikx} dx \\\ &=& F(k)G(k) \end{eqnarray} }$

つまり,入力テクスチャのフーリエ変換結果と,フィルタ画像のフーリエ変換結果を用意して,ピクセル同士の積を行った結果を逆フーリエ変換することで,畳み込み演算の結果を取得できます。
まともにカーネルを作って畳み込み計算をおこなうと,カーネルサイズに比例して計算速度が下がりますが, フーリエ変換を行えばどのようなサイズのカーネルでも一定の時間で計算が可能で,縮小バッファも用いないのでキレイです。

畳み込みブラーの適用結果

まずはガウスブラーをかけてみます。計算するとわかりますが,ガウス関数はフーリエ変換してもガウス関数です。

${ \displaystyle \begin{eqnarray} \frac{1}{\sqrt{2\pi}\sigma} e^{-\frac{x^2}{2\sigma^2}} \rightarrow e^{-\frac{\sigma^2x^2}{2}} \end{eqnarray} }$

標準偏差が ${\sigma}$ のガウスブラーをかけるときは,入力テクスチャをフーリエ変換した結果の中央からの距離 ${r}$ の位置に, ${e^{-\frac{\sigma^2r^2}{2}}}$ の重みを掛けてから逆フーリエ変換することで,どのようなサイズのガウスブラーもかけることができます。

${\sigma=32}$	${\sigma=128}$

512x512サイズのテクスチャに対して,FFTと畳み込みブラーを実行した結果です。見やすくなるように,トーンマップとガンマ補正を適用しています。 512x512サイズの畳み込み演算のGPU実行時間は,GTX970の環境で0.13ms程度でした。

大雑把に計算して,FFTに0.3ms,畳み込みに0.13ms,逆FFTで0.3msで0.73ms程度のGPU負荷で,どんなサイズのブラーも掛けられるということになります。実際は入力のFFT結果はブルーム処理やDoF処理で使いまわすので,実用上では1枚の縮小バッファブラーあたり0.5ms程度となるのでしょうか。決して速くはありあせんが,綺麗な結果を得られることを考えると,良い方法な気がします。

ついでに,ディスクブラーもFFTとの畳込みでやってみました。円形開口のマスクのフーリエ変換はベッセル関数とかになった記憶があるんですが,面倒なのでこれはそのまま円形マスクテクスチャをFFTした結果を利用します。

${r=16}$	${r=32}$

大きなディスクブラーもキレイに出来ています。普通にカーネルを作ってディスクブラーを行おうとすると, 13x13のディスクブラーでも49サンプルぐらいしないとダメだったのでそれを考えるとFFTの方が良いのかもしれません。

また,ディスクブラーやガウスブラーみたいな定まった形ではなく,好きなカーネルを使えるのもメリットです。

f:id:hikita12312:20171111174923p:plain:w128

例えば,星型のカーネルを用意して,このカーネルをフーリエ変換したもので畳み込めば...

f:id:hikita12312:20171111175642p:plain:w600

ギザギザとした感じのフィルタがかけられます。どんなに複雑な形状でも,計算速度は変わりません。

なんだか,学校の課題みたいなネタでしたが,畳み込みブラーをうまく使ったら高品質なブルームやDoFができそうです。

2017-11-05

Compute ShaderでSummed Area Table

Summed Area Table

コンピュートシェーダーに興味が出てきたので,練習がてらSummed Area Table作成してブラーをかけてみます。

https://docs.nvidia.com/gameworks/content/gameworkslibrary/graphicssamples/d3d_samples/d3dcomputefiltersample.htm

https://www.slideshare.net/egorodet/cpu-is-in-focus-again-implementing-dof-on-cpu

Summed Area Tableとは,テクスチャのテクセル ${t_{xy}}$ と表したときに,以下のようにX方向,Y方向に順番に和を取っていったテーブルのことです

${ \displaystyle \begin{eqnarray} T_{xy} = \sum_{i\leq x} \sum_{j\leq y} t_{ij} \end{eqnarray} }$

f:id:hikita12312:20171105202149p:plain:w400

Summed Area Tableを作成しておくと,テーブルから4点フェッチして重みを工夫してサンプルするだけで, 元のテクスチャの任意の矩形範囲内部のピクセルの積分が計算できます。

f:id:hikita12312:20171105202836p:plain:w300

${ \displaystyle \begin{eqnarray} I(x_1,y_1,x_2,y_2) = \frac{T_{x_2y_2} - T_{x_1y_2} - T_{x_2y_1} + T_{x_1y_1}}{(x_2-x_1)(y_2-y_1)} \end{eqnarray} }$

テーブルさえ作ってしまえば,どんな積分範囲であっても高速に計算ができるので, 被写界深度処理のようなピクセルによってブラー範囲(=積分範囲)を変えるというケースで特に有効かもしれません。

ただし,Summed Area Table自体の作成するときに, 全ピクセルに対して依存するピクセルの総和と計算しておく必要があるので,この計算コストが問題となります。

Summed Area Tableの作成

Summed Area Tableの作成には,Compute Shaderを利用しました。
実際のシェーダーコードとCPU側のソースは次のような感じです。

//// HLSL側のCompute Shader
Texture2D<float4> g_SrcTexture : register(t0);
SamplerState g_LinearSampler : register(s0);
RWTexture2D<float4> g_DstTexture : register(u0);
#define    PIXEL_SIZE g_aTemp[0].xy
groupshared float3 s_aTempBuffer[LENGTH];

[numthreads(LENGTH, 1, 1)]
void CS_CreateTable(uint3 position : SV_DispatchThreadID)
{
    const uint bufferId = position.x;
#ifdef PASS_1ST
    // 水平パス
    const uint2 texturePos = position.xy;
#endif
#ifdef PASS_2ND
    // 垂直パス
    const uint2 texturePos = position.yx;
#endif
        // テクスチャの読み込み
    s_aTempBuffer[bufferId] = SampleLod0(g_SrcTexture, g_LinearSampler, float2(texturePos)*PIXEL_SIZE).rgb;

    // イテレーション
    float3 addBuf;
    for (uint i = 0; i < ITERATION; ++i)
    {
        GroupMemoryBarrierWithGroupSync();
        const uint delta = (0x1 << i);
        const float3 buffer = (bufferId >= delta) ? (s_aTempBuffer[bufferId - delta]) : (0.0f);
        GroupMemoryBarrierWithGroupSync();
        s_aTempBuffer[bufferId] += buffer;

    }
    float4 ret = float4(0.0f, 0.0f, 0.0f, 1.0f);
    ret.rgb = s_aTempBuffer[bufferId];
    g_DstTexture[texturePos] = ret;
}

//// CPP側の擬似コード
// 水平方向加算
auto pPass = rPfxEffect.GetPass("SummedAreaTableCreate", "L512_1ST").lock();
auto pTempBuffer = rPfxEffect.GetTempTexturePool().Get(dstSize, dstSize, tableFormat, bindFlags).GetTextureView();
pPass->CSBindConstantBuffer(0, tempCB);
pPass->CSBindTexture(0, pSrcTextureView);
pPass->CSBindSampler(0, rPfxEffect.GetLinearSampler());
pPass->CSBindRWTexture(0, pTempBuffer);
result = rPfxEffect.Dispatch(pDeviceContext, pPass, 1, dstSize, 1);
// 垂直方向加算
pPass = rPfxEffect.GetPass("SummedAreaTableCreate", "L512_2ND").lock();
auto pResultTextureView = rPfxEffect.GetTempTexturePool().Get(dstSize, dstSize, tableFormat, bindFlags).GetTextureView();
pPass->CSBindConstantBuffer(0, tempCB);
pPass->CSBindTexture(0, pTempBuffer);
pPass->CSBindSampler(0, rPfxEffect.GetLinearSampler());
pPass->CSBindRWTexture(0, pResultTextureView);
result = rPfxEffect.Dispatch(pDeviceContext, pPass, 1, dstSize, 1);

このコードの要点をまとめると

X方向の加算とY方向の加算を別パスに分離
再帰的に加算をして計算量をテクスチャ幅のlogに落とす
シェーダーのグループ共有メモリを利用して,テクスチャフェッチ回数を最小限に減らす

という感じです。

何も考えずに加算をすると,テクスチャ幅 ${w}$ ]として,一つのピクセルあたり ${w}$ 回のループを回すことになります。
これを,以下の図のように2べき飛ばしで再帰的に加算をしていくことで, ${\log{(w)}}$ 回まで減らすことが可能です。

f:id:hikita12312:20171105204230p:plain:w300

この再帰計算を行うときに便利なのがgroupsharedで宣言されたグループ共有メモリです。 ${w}$ 個のスレッドを走らせるグループを ${h(=w)}$ 個作成し,それぞれのグループで ${w}$ 個のfloat3を計算バッファとして確保しておきます。
最初の一回だけこのバッファに対してテクスチャ読み出しをおこなって,あとの再帰的な加算を行うときはバッファを利用することで,テクスチャ読み出しを必要最低限に抑えます。

f:id:hikita12312:20171105205815p:plain:w600

Summed Area Tableの作成自体はVertex ShaderとPixel Shaderでも工夫をすれば可能ですが,再帰的に加算をする度にテクスチャ読み出しをするので, Compute Shaderほど高速には行なえません。 DirectX11ではグループ共有メモリは32KBまでしか利用できないので,float4で考えるとテクスチャサイズは2048が限界ということに注意が必要です。また,加算時の誤差が激しいので,32bit浮動小数テクスチャを利用しています。

ブラー結果

512x512ピクセルのSummed Area Tableを作成すると,次のような感じになりました。

入力テクスチャ	Summed Area Table

Summed Area Tableの右端の値は9000ぐらいで,普通に表示すると飽和して見えないので,適当にスケールしてあります。
このテーブルを利用して,ただ矩形状に平均値を計算するボックスフィルタを適用してみます。

f:id:hikita12312:20171105210816p:plain:w400

だいぶカクカクとした結果になりました。今度は矩形を3つ重ねて少し滑らかになるようなフィルタを適用してみます。

f:id:hikita12312:20171105211251p:plain:w150

f:id:hikita12312:20171105211313p:plain:w400

まずまずではないでしょうか。
このパス全体のGPU実行時間は,GTX970の環境で

512x512テクスチャのTable作成 : 0.14 ms
1280x720のテクスチャへのブラー書き出し : 0.34 ms

程度のGPU時間でした。なかなか高速だと思います。

Log Summed Area Table

おまけです。
Summed Area Tableを作成するときに,対数平均を取るようにして,ブラー時にexpで元に戻すようにしたら,誤差が減るのではないかと思ってやってみました。

${ \displaystyle \begin{eqnarray} T_{xy} = \sum_{i\leq x} \sum_{j\leq y} \log{(t_{ij})} \end{eqnarray} }$

Summed Area Table	Log Summed Area Table

Summed Area Table	Log Summed Area Table

少しブラーの角が取れたような気もします。

2017-10-25

プロファイラの作成

ImGuiを使って,フレームのCPU時間とGPU時間を可視化するプロファイラを作りました。

f:id:hikita12312:20171025222821p:plain

それだけです。

プロファイラを作ってわかったのですが,被写界深度処理が重いので,やはり最適化の余地ありという気がします。

2017-10-15

ImGuiの導入

ImGuiとは

ImGuiとは,楽に組み込みができて,簡易的なGUIを提供するライブラリです。

https://github.com/ocornut/imgui

例えば,このようなコードを書くと,

ImGui::Begin("Stats", nullptr, window_flags);
ImGui::LabelText("", "FPS:%4.2f (%4.2f ms)", s_Fps);
ImGui::LabelText("", "Camera Pos:(%.2f, %.2f, %.2f )", cameraPosition.x, cameraPosition.y, cameraPosition.z);
ImGui::Spacing();
ImGui::SetNextTreeNodeOpen(true, ImGuiSetCond_Once);
if (ImGui::CollapsingHeader("PostEffect"))
{
    if (ImGui::TreeNode("Tonemap"))
    {
        ImGui::SliderFloat("Gamma", &pfxDesc.Tonemap.Gamma, 0.0f, 4.0f);
        // ...
        ImGui::TreePop();
    }
    // ...
    ImGui::TreePop();
}
//...
ImGui::End();

ImGui::Render();

こんな感じでGUIが作れます。

f:id:hikita12312:20171015113931p:plain:w600

手続き的にImGuiのnamespaceのAPIをコールするだけで,それっぽいGUIが出来上がるのでとても便利です。

ImGuiの組み込み方法

ソースコードを持ってくる

ImGuiのリポジトリからソースを持ってきて,自分のプロジェクトに必要なソースとヘッダを組み込みます。必要なファイルは

imgui.h
imgui.cpp
imgui_draw.cpp
imgui_internal.h
stb_rect_pack.h
stb_textedit.h
stb_truetype.h

です。自分のライブラリに組み込むような場合は,imgui.hのみ公開ヘッダとすれば大丈夫かと思います。

プラットフォームに合わせた実装をする

(簡単) examplesを利用する方法

ImGuiを利用するためには,描画部分とマウスやキーボードの入力部分を実装しなければいけません。
その実装例がImGuiのリポジトリのexamplesディレクトリ以下に格納されているので,面倒ならばこのサンプルをそのまま使っても良いと思います。

例としてDX11版の実装を見てみます。 imgui_impl_dx11.hを見ると,

IMGUI_API bool        ImGui_ImplDX11_Init(void* hwnd, ID3D11Device* device, ID3D11DeviceContext* device_context);
IMGUI_API void        ImGui_ImplDX11_Shutdown();
IMGUI_API void        ImGui_ImplDX11_NewFrame();
IMGUI_API void        ImGui_ImplDX11_InvalidateDeviceObjects();
IMGUI_API bool        ImGui_ImplDX11_CreateDeviceObjects();

というAPIが定義されていて,Init()で初期化を行い,デバイスの作成・破棄が行われたらInvalidateDeviceObjects(),CreateDeviceObjects()をコールします。 ImGui名前空間のAPIをコールする前には,NewFrame()を呼ぶ必要があります。また,WindowsAPIのイベントを取得するために,

IMGUI_API LRESULT   ImGui_ImplDX11_WndProcHandler(HWND hWnd, UINT msg, WPARAM wParam, LPARAM lParam);

ImGui_ImplDX11_WndProcHandler()をグローバルな名前空間にexternで定義しておいて,組み込み先のアプリケーションのWNDPROCのコールバック内部でコールするようにします。

これだけで,とりあえずは自分のプロジェクトでImGuiを動作させるところまでは持っていくことができます。 DX11以外のプラットフォームでも同じような感じで組み込めるかと思います。

真面目に組み込む方法

examplesに自分の使いたいプラットフォームが無い
独自に抽象化されたフレームワークを利用している,もしくは作成している

という場合には,自分で真面目にImGui::GetIO()を利用してImGuiの描画部分と入力部分を実装する必要があります。これもexampleのサンプルコードの中身を参考にすれば良いです。 ImGuiIO::KeyMapやImGuiIO::MousePosなどに入力値を格納しておけばImGui内部で利用する入力として利用してくれます。

描画部分は,ImGuiIO::RenderDrawListsFnにImGui内部で利用するImDrawData構造体のを引数に,描画処理を行うコールバックを実装する必要があります。 ImDrawDataには描画に利用する頂点データや行列変換の定数バッファが格納されていて,これを利用して描画を行うことになります。そのために,ImGui利用前にプラットフォーム固有のシェーダーや,Font用のテクスチャを作成する必要があります。

ひっかかり
シェーダーリフレクションを利用してシェーダーコードからInputLayoutを作っていたので, ImGuiの組み込みをしたときに,渡される頂点バッファのColor部分のデータが8bitUNORMのフォーマットだったことに少し躓きました。 DX11の場合はリフレクションで8bitフォーマットかどうかはわからないので,手動でInputLayoutを作成する必要があります。

2017-10-12

スキャッターベースのボケ

ギャザーベースとスキャッターベース

前回までの被写界深度処理は,ギャザーベースと呼ばれているカテゴリの方法です。周囲のサンプリングして集めて一つのピクセルの色を決定するという意味でギャザーベースと言われているのだと思います。
ただし,ギャザーベースでは"注目するピクセル"に対する処理なので, 前ボケのブラーを作るときにも問題になったように, 前ボケが焦点の合っている場所の上に重なるように生じる現象はできません。

https://gpuopen.com/gdc2017-cinematic-depth-of-field/

また,大きなサイズのボケを作ると,大きなカーネルが必要になったり,ボケの細かい形状の制御が出来なかったりする問題があります。

そこで,今回はスキャッターベースと呼ばれる,あるピクセルの情報を周囲のピクセルに書き込む方法でボケの処理を実装してみようと思います。

スキャッターベースの概略

以下のサイトを参考にしました。

https://www.slideshare.net/siliconstudio/cedec-2008imagire-day

https://www.slideshare.net/siliconstudio/cedec-2012-34760454

http://hexadrive.jp/lab/demo/726/

今回やってみた方法の基本的な流れはこんな感じです。

f:id:hikita12312:20171010102539p:plain:w600

画面全体のピクセルに対応する位置に縮退した矩形を作れるように頂点として入力します。その4頂点をVertexShaderでボケサイズに応じて膨らませて,PixelShaderでボケの形状を決定します。周囲のピクセルからサンプリングするわけではなく,ボケ自体が膨らむのがポイントです。
膨らませる処理にGeometryShaderを利用する方法もあると思いますが,今回は試していません。

マスク処理

画面の解像度が1280x720である場合,そのままピクセルとボケの矩形を対応させると, 一つの矩形あたり2個の三角形を書くとして約164万ポリゴンのレンダリングとなります。こうなってしまうとリアルタイムで行うことが難しいので,例えば縮小バッファを利用して,1/2サイズのバッファを利用します。バッファのサイズを半分にすれば4倍高速化できますが,それでも40万ポリゴン近くです。

そこで,今回は目立つボケ以外はスキャッター処理を行わないというマスクを作成しました。

元画像	マスク

マスクの条件は,

輝度のエッジとなっている
輝度が高い
ボケの直径(CoC)がある程度の大きさ以上

となっていて,輝度が高く周囲よりも明るいボケが目立つ領域であると仮定したマスクとなっています。エッジの計算はSobelフィルターを利用しています。マスクされる領域のスキャッターボケは,VertexShaderで頂点を膨らませないで縮退させておくと,ラスタライズの時点で無視されるので高速化出来ます。

ボケの形状の作成

VertexShaderで膨らんだスキャッターボケの矩形領域に対して,実際に五角形や六角形のボケを書いていきます。今回書いたボケ形状用のPixelShaderのコード片は以下です。

float4 PS_Create(SCATTER_BOKSH_VS_OUT In) : SV_TARGET0
{
    const float coc = In.Color.a;
    // 表示範囲が前ボケ,後ボケ境界からはみ出すならば描画終了(後述)
    if ((coc > 0.0f) && (In.UV.z >= 0.0f)) discard;
    if ((coc < 0.0f) && (In.UV.z <= 0.0f)) discard;

    float3 color = In.Color.rgb;
    float2 rPos = In.UV.xy;
    if (coc < 0.0f) rPos *= 1.0f; // 前ボケのときの形状反転

    // ボケ形状の計算
    const float deltaTheta = 2 * PI / APERTURE_BLADE;
    // 描画点の光軸中心からの角度 (0~2pi)
    float theta = atan2(rPos.y, rPos.x) + PI + APERTURE_ROTATE;
    // APERTURE_BLADE角形のi番目頂点の角度を基準にした相対theta (0~2pi/APERTURE_BLADE)
    const float modTheta = frac(theta / deltaTheta)*deltaTheta;
    // bokehRadius=1.0としたときの,APERTURE_BLADE角形の境界となる境界
    const float effBound = 1.0f - (1.0f - APERTURE_CIRCULARITY)*(1.0f - cos(deltaTheta / 2.0f))*sin(APERTURE_BLADE / 2.0f*modTheta); 
    const float bokehRadius2 = rPos.x*rPos.x + rPos.y*rPos.y; // (0~1)
    const float effBokehRadius = sqrt(bokehRadius2) / effBound;
    // ボケ形状範囲外ならば描画しない
    if (effBokehRadius > 1.0f) discard;

    // ボケの形状による減衰(明度)を計算(0.0~1.0)
    float attenuation = saturate(1.0f - pow(effBokehRadius, 10.0f)); // 後ボケは中心ほど明るい
    // 前ボケならば明度を反転,
    if (coc < 0.0f) attenuation = 1.0f - attenuation;
    // CoCの大きさによる減衰
    // 本来は面積に反比例するが,マスクによって,積分されるボケ減っているので,面積反比例だと,あまりに薄い
    // 厳密な計算をするべきだが,coc^-1の重みでなんとなく良さそうだったのでこれで
    const float cocAttenuationWeight = min(1.0f, 1.0f / coc);
    attenuation *= cocAttenuationWeight;
    // 明度の計算
    color.rgb *= attenuation*APERTURE_BRIGHTNESS_SCALE;

    return float4(color, coc);
}

ボケの内部の色味だとか収差だとかは今回は適当に作っています。ボケ形状を作成する部分は,APERTURE_BLADEで ${n}$ 角形を指定し,APERTURE_CIRCULARITYでボケの円形度合いを決めています。

f:id:hikita12312:20171011104928p:plain:w300

円弧を多角形の辺を軸と見たときのsinカーブに近似できるとして計算しています。実際にPixelShaderから出力されたボケ形状は以下の画像のようになりました。適当な近似でもそれっぽくなっていると思います。

f:id:hikita12312:20171011105944p:plain

本来ボケの面積に従って明度が反比例していくべきですが,今回の場合ですとマスクによってボケの総量が減っているので,単純に面積に反比例だと色が薄くなりすぎます。アドホックに直径に反比例という形で明度を減衰しています。

スキャッターボケテクスチャの作成

PixelShaderで作成したボケテクスチャを加算合成してテクスチャに書き出します。このとき,入力テクスチャと同じサイズのテクスチャではなく,入力サイズの高さを2倍にしたアスペクト比で,異なる倍率の縮小バッファを3枚用意しました。

f:id:hikita12312:20171012102919p:plain:w400

VertexShaderで頂点を膨らませるときに細工をして,前ボケならばテクスチャの上半分,後ボケならばテクスチャの下半分に書き出すようにします。カメラのボケは前ボケと後ボケで反転したり収差が変わったりして異なる形状になりますが,それを分離して表現することが可能となります。

f:id:hikita12312:20171012104625p:plain:w600

さらに,縮小バッファとして用意した異なるサイズ毎に,VerterxShaderでの膨らましとPixelShaderでの書き込みの一連のパスをそれぞれ行います。大きなサイズのボケを書き出すときに,膨らんだスキャッター矩形の内部のピクセル数が多ければ多いほど負荷が高くなりますが, 大きいサイズのボケほど小さい縮小バッファを使うようにすれば負荷を抑えることが出来ます。大きいボケには,それほどディティールが必要ないので縮小バッファで十分です。ドローコールは3倍に増えますが,大きいボケの負荷を下げることが出来るので全体としての負荷は下がります。さらに,ボケのサイズ毎に3つのテクスチャに分類して書き出すことで,ボケの前後関係を分けて書き出すことができることもメリットです。

結果

出来上がったスキャッターボケテクスチャをギャザーボケの結果に加算合成することで,最終出力が得られます。スキャッターボケテクスチャのアルファ値にそのボケの生じた距離を格納することで,合成を行うときに後ボケをマスクをすることもできます。
(2017/10/12現在,合成周りがまだ甘いです...)

以下が結果です。 f:id:hikita12312:20171012105802p:plain:w600 こちら五角形の後ボケ,ある程度大きいボケでもキレイに表示されています。

ギャザーボケ六角形	スキャッターボケ六角形

フラットなブラーの結果と比較しました。ギャザーベースのボケよりハッキリボケの輪郭が出ていて,大きいボケが表現できています。

f:id:hikita12312:20171012110415p:plain:w600

今度は円形ボケにして接写的に焦点を手前に合わせてみました。だいたい良いんですが,一部の後ボケがマスクされていなかったり, ギャザーボケとスキャッターボケの合成が甘くて,ギャザーボケの結果が薄っすら見えてしまったりと改善の余地はありそうです。

2017-10-06

倍率色収差の実装

倍率色収差

今回は簡単に出来そうな倍率色収差を実装してみます。収差とは,レンズを通る光の波長の違いによるズレによる現象のことですが,倍率色収差は画面の周辺部にかけて色がズレていく現象のことだそうです。

http://cweb.canon.jp/eos/special/dlo/factor/index02.html
https://docs.unity3d.com/jp/540/Manual/script-VignettingAndChromaticAberration.html

実装

基本的な実装は,元となるテクスチャに対して幾つかのチャンネルでマスクをかけたものをそれぞれ拡大して加算合成すれば実現できます。

f:id:hikita12312:20171005235956p:plain:w600

#define CHROMATIC_ABERRATION_ARRAY_NUM 5
static const float3 s_aChoromaticAberrationFilter[CHROMATIC_ABERRATION_ARRAY_NUM] =
{
    float3(0.0f, 0.0f, 0.5f),
    float3(0.0f, 0.25f, 0.5f),
    float3(0.0f, 0.5f, 0.0f),
    float3(0.5f, 0.25f, 0.0f),
    float3(0.5f, 0.0f, 0.0f)
};

float4 PS_ChromaticAberration(PPFX_OUT In) : SV_TARGET0
{
    const float2 uv = In.UV.xy;
    float4 outputColor = float4(0.0f, 0.0f, 0.0f, 1.0f);

    float2 screenPos = TexCoordToScreenSpace(uv);
    for (uint i = 0; i < CHROMATIC_ABERRATION_ARRAY_NUM; ++i)
    {
        const float2 tempUV = ScreenSpaceToTexCoord(screenPos * (1.0f - CRHOMATIC_ABERRATION_SCALE * (i + 1) / CHROMATIC_ABERRATION_ARRAY_NUM));
        const float3 mask = s_aChoromaticAberrationFilter[i];
        const float3 srcColor = SampleLod0(g_SrcTexture, g_LinearSampler, tempUV).rgb;
        outputColor.rgb += srcColor * mask;
    }

    return outputColor;
}

以下が実装結果です。

倍率色収差無し	倍率色収差有り

f:id:hikita12312:20171006000534p:plain:w600

倍率色収差をかけると,画面の端の方ほど,色がズレていることが分かるかと思います。

2017-09-26

フラットなブラー

ボケの形

被写界深度表現を実装しましたが, なんとなくガウスブラーで画面をボカシていました。実際はカメラのボケというようのは,ガウシアンのボケではなく,カメラの絞りの形状による円形だったり六角形のフラットなボケとなるはずです。今回はガウスブラーではない,カメラの絞りの形状のような形となるブラーの作り方を考えます。

ガウスブラー

まずは,ブルームや被写界深度を実装したときに作った,ガウスブラーの画像を比較までに置いておきます。

f:id:hikita12312:20170926223258p:plain:w600

ブラー無し	ガウスブラー

なんだかボケの輪郭がボヤッとしてますね。

ディスクブラー

円形のボケを実装します。以下のスライドの134ページ目あたりを参考にしました。

https://www.slideshare.net/siliconstudio/cedec-2009imagire-day-2009

f:id:hikita12312:20170926223603p:plain:w300

上の図が9x9のディスクブラーを作成するときのサンプリングの概略図です。
シェーダーのBilinear補間を考慮してサンプル位置とサンプル重みを決定します。円に完全にフェッチする4x4(2x1,1x2)のテクセルの矩形が含まれている場合には,目的のテクセルの中央を等しい重みでサンプルしますが, 円周がフェッチするテクセルの矩形に含まれる場合には,含まれる面積を考慮してサンプリング位置を円の内側へずらし,重みも面積に比例するように調整します。また,ブラーをかけて元の画像よりも明度が変わるのを避けるために,サンプリングする重みは正規化して総和が1となるように調整します。

それで出来上がった結果がこちらです。

f:id:hikita12312:20170926224131p:plain:w600

ガウスブラー	ディスクブラー

輝度が高い部分のボケが円形のサンプリングに従ってキレイに丸くなっています。

六角形ブラー

今度は六角形のボケを実装します。以下の記事を参考にしました。

https://colinbarrebrisebois.com/2017/04/18/hexagonal-bokeh-blur-revisited/

実現方法として,まずは1回目のパスで六角形の垂直方向の対角線のブラーと対角120度方向のブラーのMRTを行います。 2回目のパスで1回目の結果を六角形になるように延ばすことで,六角形を作ります。

f:id:hikita12312:20170926224800p:plain:w600

詳しいシェーダーコードなどは,参考にした記事を参照してください。　　 1パス目のMRTを使って作った2本の方向のブラーが頭が良くて,2回めのパスで六角形の余った領域を240度方向への1回のブラーでキレイに埋めています。実は上の図で表した方法は,参考にした記事から少しだけ実装が変えていて, 垂直方向に関連するブラーをするときに+1のオフセットを与えて重なる領域を作らないようにすることで,六角形の対角が重なって見える問題を解決しています。

半径 ${r}$ の六角形を作るときに, ${2r+2r}$ 回程度のテクセルフェッチで実現できます。 MRTを使っている分,コスト面では劣りますが,実質ガウスブラーと同程度のサンプル数で実現できます。

f:id:hikita12312:20170926225548p:plain:w600

ガウスブラー	六角形ブラー

実際に適用した結果です。六角形が出来ていることがわかります。