Gaussian Blur

## 1. Gaussian Blur (一) ―― 基礎アルゴリズムとその実装


### 1.1 基礎知識
ガウシアンぼかし(Gaussian Blur)は、Gaussian Smoothingとも呼ばれ、画像をぼかすための典型的なアルゴリズムです。ガウシアンぼかしは、簡単に言うと、画像全体を加重平均する処理で、各画素値は自身とフィールド内の他の画素値の加重平均で得られます。

ガウシアンぼかしのアルゴリズムに入る前に、畳み込み(Convolution)演算を理解する必要があります。ここでは、行列の畳み込み演算プロセスを簡単に説明します。具体的な畳み込み演算については、関連資料を調べてください。

処理する行列があります。
![14.1611111467927.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/044e7371-632c-f6d8-189c-804f229cf4ea.png)



この行列が3*3の畳み込みカーネルを用いて処理されると仮定すると、このカーネルは次のようになります。

![15.1611111505375.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/8ac99b1a-7579-9e18-db3f-bb7b59992cf9.png)




まず、畳み込みカーネルを180°反転し、元の畳み込みカーネルと同じ結果を得ます。次に、畳み込みカーネルの中心を処理対象の行列の対応する要素、例えば左上の最初の要素に合わせ、要素のないところには0を加え、対応する要素を掛け合わせて次の結果を得ます。

![16.1611111556573.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/27119429-0dd6-b2b2-2a28-386fed95d531.png)




0*0.1+0*0.1+0*0.1+0*0.1+0*0.1+1*0.2+2*0.1+0*0.1+5*0.1+6*0.1=1.5



畳み込み演算後の最初の要素の値は1.5です。同様に、行列の他の値も最終結果まで処理します。
![17.1611113313051.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/be044292-0d70-943f-f4cc-99b0a7e42815.png)





これは、SWA（Sliding Window Algorithm）の応用で、簡単な行列の畳み込み演算プロセスです。

ガウシアンぼかしはこのような畳み込み演算を利用したもので、畳み込みカーネルの計算がガウシアン分布関数を適用していることから「ガウシアン」ぼかしと呼ばれています。




![18.1611111916122.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/70fe9f7c-5ff0-b7b4-25be-2b6eed7a8c16.png)

(一次元ガウシアン正規分布の密度関数）

![19.1611111954342.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/f32f7fcb-9bce-1350-461d-7307d42bbd9b.png)


(一次元ガウシアン正規分布の密度関数グラフ）。



さらに、2次元のガウシアン分布の式は次のように求められます。
![20.1611111614910.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/94715e4f-ce26-e5cd-7c44-0eaf5174f1a6.png)



これを正規化すると、次のようになります。
![21.1611112175805.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/87899896-56a6-7399-1821-caddd20a682b.png)



ここで、uは行、vは列を表し、u, v∈｛-ω, -ω+1, .........｝とします。ω-1, ω}

Sは正規化定数。
![22.1611111643013.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/c5335ccd-b0e6-b065-4f3a-2a0315dc8f9c.png)



概略図は次のようになります。

![23.1611111678843.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/d5a9224c-ad58-c41a-9579-d369f72af69b.png)




3×3の畳み込みカーネルの場合、その個々の位置の相対座標は次のように表すことができます。

![24.1611111698440.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/bac9a79f-44ad-dfdc-a72b-0ba5c89a6cb6.png)




式に代入して演算すると、正規化された3*3ガウシアン畳み込みカーネルが得られます。

![25.1611111723104.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/5c8a8144-402e-e259-d6b2-7780fe3f1811.png)




画像に対して、ガウシアン畳み込みカーネルを用いて全画素を畳み込み演算します。 最終的に得られる画像は、画像がぼかされたように見えるので、ガウシアンぼかしと呼ばれています。

出力画像をY、入力画像をXとし、i行j列のデータをX（i,j）、Y（i,j）とすると、サイズが（2ω+1）*（2ω+1）、標準偏差がσのガウシアンカーネルを用いた計算後の結果は、次のようになります。
![26.1611112100752.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/9931cbc8-d5f6-ff37-b8de-c3557a0870cf.png)



この式に従って、ガウシアンカーネルの中心を入力画像の位置(i,j)に置き、ガウシアンカーネルの各値と入力画像の対応する位置の値を掛け合わせて、（2ω+1）*（2ω+1）回の乗算計算を行い、（2ω+1）*（2ω+1）-1回の加算計算を行って、時間複雑度はO（ω^2）である。

また、サイズ1024*1024の画像に対して、サイズ33*33の畳み込みカーネルを使って画像全体に1回のぼかし処理を行うには、1024*1024*33*33≒11億4000万回のマッピングリードが必要となり、非常に時間がかかるのは間違いない。



### 1.2 Unityの実装
上記のアルゴリズムを元に、Unityで実装すると以下のようになります。

Shaderを作成し、 2 次元ガウシアンカーネルの計算式を実装します。

```
float gauss(float x, float y, float sigma)
{
	return  1.0f / (2.0f * PI * sigma * sigma) * exp(-(x * x + y * y) / (2.0f * sigma * sigma));
}
```

頂点シェーダは UnityCG.cginc に実装された vert_img 関数を選択し、その主な機能は、頂点をモデル空間からクロップ空間へ変換することです。
```
v2f_img vert_img( appdata_img v )
{
    v2f_img o;
    UNITY_INITIALIZE_OUTPUT(v2f_img, o);
    UNITY_SETUP_INSTANCE_ID(v);
    UNITY_INITIALIZE_VERTEX_OUTPUT_STEREO(o);
    o.pos = UnityObjectToClipPos (v.vertex);
    o.uv = v.texcoord;
    return o;
}
```
ガウシアンカーネルのサイズを設定します。
```
#ifdef MEDIUM_KERNEL
#define KERNEL_SIZE 35
#elif BIG_KERNEL
#define KERNEL_SIZE 127
#else
#define KERNEL_SIZE 7
#endif
```
各ピクセルとそのフィールドのピクセルに重みをつけるフラグメントシェーダを実装します。重みは上記のgauss関数によって得られます。
```
float4 frag(v2f_img i) : COLOR
{
	float4 o = 0;
	float sum = 0;
	float2 uvOffset;
	float weight;
	for (int x = -KERNEL_SIZE / 2; x <= KERNEL_SIZE / 2; ++x)
		for (int y = -KERNEL_SIZE / 2; y <= KERNEL_SIZE / 2; ++y)
		{
            //オフセットを計算する
			uvOffset = i.uv;
			uvOffset.x += x * _MainTex_TexelSize.x;
			uvOffset.y += y * _MainTex_TexelSize.y;
            //重みを確認する
			weight = gauss(x, y, _Sigma);
			o += tex2D(_MainTex, uvOffset) * weight;
			sum += weight;
		}
	o *= (1.0f / sum);
	return o;
}
```
以下のような結果が得られます。
![27.1611111799580.jpg](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/35d6f87c-25df-9eb0-339c-4f5db0bb78ca.jpeg)

原画像


![28.1611111771747.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/d0843fb0-f844-912f-df02-beff294a519d.png)

(7*7のガウシアンカーネルを使用したもの)



![1.1609836981604.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/1577059e-324d-9f1b-ad88-45012246a684.png)

(35*35のガウシアンカーネルを使用したもの)

DEMOリンク

https://github.com/UWA-MakeItSimple/Course-PostProcessingEffect/blob/main/Assets/Shaders/Blur/GaussianBlur/SecondDimensionGaussianBlur.shader

---

## 2. Gaussian Blur（二）―― 2段階の1次元演算アルゴリズムとその実装
### 2.1 基礎知識
2次元ガウシアン関数を用いた畳み込みカーネルによる画像のぼかし処理は，リアルタイム計算を行うゲームに適用するとオーバーヘッドが大きくなります。加速を最適化し、テクスチャ読み取り操作と算術操作の数を減らすために、いくつかの方法を検討する必要があります。最適化は、2次元ガウシアン分布関数の分離可能な性質を利用することによって実行できます。

2次元ガウシアン分布関数の分離可能性とは、2次元ガウシアン関数を同じ1次元ガウシアン関数に分解して2回処理しても、どちらも等価であり、得られる結果は同じだが処理時間が大幅に短縮されることを意味します。数学的な証明式は以下の通りです。

![30.1611112495702.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/83bdaf9e-d905-3b8d-727a-059f37918fc1.png)


とすると、

![31.1611112523903.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/fbe2b716-c829-b05d-64cd-e9dbf6556ad5.png)




ガウシアンカーネル行列を、行ベクトルと列ベクトルを掛けた正規化定数として書き直します。
![32.1611112549278.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/8afd6f30-8473-4408-9718-d86855bc27fc.png)



このときの正規化定数Sは、2つの1次元ガウシアンカーネルの正規化係数の逆数乗算に分割できます。
![33.1611112580384.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/e594f178-d656-b248-941d-720789f231d1.png)





したがって、Gを2つのベクトルの積に分割することもできます。
![34.1611112611726.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/a9f71599-504f-f8c6-a87a-63a1fc684c35.png)



この中には、
![35.1611112636310.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/c2e895e5-55aa-d55f-f1ea-21fe3c675a4d.png)



したがって、ガウシアンぼかした画像Yは次のように表すこともできます。

![36.1611112660992.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/a617d4b7-6180-2ec9-306a-28990282fe8b.png)




上記の式は、最初に1次元ガウシアンカーネルと入力画像Xを使用して水平畳み込み演算を実行し、中間結果Zを取得し、次に1次元ガウシアンカーネルと中間結果を使用して垂直畳み込み演算を実行することを示しています。結果Zを使用して、最終結果Yを取得します。

この方法でピクセルのぼかし効果を計算するには、（4ω+ 2）乗算演算と4ω加算演算が必要で、時間計算量はO（ω）であることがわかります。

1024 * 1024サイズの画像の場合、33 * 33の畳み込みカーネルを使用して画像全体にぼかし操作を実行するには、1024 * 1024 *33*2≈6.9kwのテクスチャ読み取りを実行する必要があります。これは2次元ガウシアン畳み込みカーネルを使用した直接処理と比較して大幅に低下します。

### 2.2 Unityの実装
上記のアルゴリズムに従って、Unityで実装します。

まず、1次元ガウシアンカーネル計算式を実装します。
```
float gauss(float x, float sigma)
{
	return  1.0f / (sqrt(2.0f * PI) * sigma) * exp(-(x * x) / (2.0f * sigma * sigma));
}
```
1次元ぼかしアルゴリズムを実装します。
```
float4 GaussianBlur(pixel_info pinfo, float sigma, float2 dir)
{
	float4 o = 0;
	float sum = 0;
	float2 uvOffset;
	float weight;

	for (int kernelStep = -KERNEL_SIZE / 2; kernelStep <= KERNEL_SIZE / 2; ++kernelStep)
	{
		uvOffset = pinfo.uv;
		uvOffset.x += ((kernelStep)* pinfo.texelSize.x) * dir.x;
		uvOffset.y += ((kernelStep)* pinfo.texelSize.y) * dir.y;
		weight = gauss(kernelStep, sigma);
		o += tex2D(pinfo.tex, uvOffset) * weight;
		sum += weight;
	}
	o *= (1.0f / sum);
	return o;
}
```
最初のpassでの水平方向の1次元ぼかしを実行する：
```
struct pixel_info
{
	sampler2D tex;
	float2 uv;
	float4 texelSize;
};

float4 frag_horizontal(v2f_img i) : COLOR
{
	pixel_info pinfo;
	pinfo.tex = _MainTex;
	pinfo.uv = i.uv;
	pinfo.texelSize = _MainTex_TexelSize;
	return GaussianBlur(pinfo, _Sigma, float2(1,0));
}

Pass
{
	CGPROGRAM
	#pragma target 3.0
	#pragma vertex vert_img
	#pragma fragment frag_horizontal
	ENDCG
}
```
Unityの内蔵関数GrabPass{}を使用して、最初のpassの実行後に画面イメージをキャプチャします。この関数は、デフォルトで_GrabTextureという名前の変数に画像を保存し、次に垂直方向のぼかしを実現します。
```
float4 frag_vertical(v2f_img i) : COLOR
{
	pixel_info pinfo;
	pinfo.tex = _GrabTexture;
	pinfo.uv = i.uv;
	pinfo.texelSize = _GrabTexture_TexelSize;
	return GaussianBlur(pinfo, _Sigma, float2(0,1));
}
uniform sampler2D _GrabTexture;
uniform float4 _GrabTexture_TexelSize;
```
GrabPass {}関数によって取得された画像のUV座標は通常の画像と一致しないため、ぼかしアルゴリズム後の結果が反転しないようにUV座標を処理する必要があります。したがって、頂点シェーダーは次のように実装されます。
```
v2f_img vert_img_grab(appdata_img v)
{
	v2f_img o;
	UNITY_INITIALIZE_OUTPUT(v2f_img, o);
	UNITY_SETUP_INSTANCE_ID(v);
	UNITY_INITIALIZE_VERTEX_OUTPUT_STEREO(o);
	o.pos = UnityObjectToClipPos(v.vertex);
	o.uv = half2(v.vertex.x, 1-v.vertex.y);
	return o;
}
```
垂直ぼかしを実装するpass
```
GrabPass{}

Pass
{
	CGPROGRAM
	#pragma target 3.0
	#pragma vertex vert_img_grab
	#pragma fragment frag_vertical
	ENDCG
```
このように、ぼかし処理を行うには、2つのパスを使用する必要があって、結果として2つのDrawCallが必要となり、これはパフォーマンスに影響を与える重要な要素です。

ここでの実装では、比較的に時間がかかるGrabPass{}関数を使用します。OnRenderImage()関数で2つのblitを使用して、それぞれ水平操作と垂直操作を実行することを検討してください。

効果は次のとおりです。
![23.1611122001943.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/72f3e11b-cdcd-05d9-aab9-c33931daa570.png)
（35 * 35のガウシアンカーネルを使用した）

---

## 3. ガウシアンぼかし(三)-線形サンプリングを使用したアルゴリズムと実装
### 3.1 基礎知識
ガウシアン分布関数の特性を利用することで、テクスチャの読み取りと算術演算の数を減らすことに成功しました。GPUを活用することで、さらなる高速化が可能です。これまでのアルゴリズムでは、1回のテクスチャ読み取りで1ピクセルの情報しか得られないことを前提としていましたが、GPUの場合は必ずしもそうとは限りません。画像の読み取りにバイリニア補間サンプリング（Bilinear Sampling）が有効になっている場合、GPUは複数のピクセル情報を一度に読み取ることができます。バイリニア補間を使用しても、GPUへ追加の負担をかけません。

バイリニア補間サンプリングについて、簡単に紹介します：
![38.1611112941173.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/a3c49a4d-7def-84dd-5039-09250f110606.png)

（バイリニア補間の概略図）



図に示すように、最初にx方向に2つの単線形補間を取得して、2つの一時点R1（x、y1）とR2（x、y2）を取得し、次にy方向に1つの単線形補間を計算してP（x、 y）を得ます。

つまり、テクスチャピクセルを読み取る場合、テクセルの中心にあるテクスチャを読み取らないように選択し、読み取るのに適した位置を選択することで、GPUバイリニア補間によって2つのピクセルの情報を取得できるようになります。

離散サンプリングと同じ効果を維持するには、位置の重みが2つのテクセルの重みの合計に等しくなるように、読み取り座標を調整する必要があります。つまり、簡略化されたガウシアンカーネルの重みと位置は次の要件を満たす必要があります：
![39.1611112980998.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/e8f8120f-6204-73f9-e1bb-3ab09bb56dd7.png)





たとえば、5*5のガウシアンカーネルの場合、2つの1次元ガウシアンカーネルに分割します。ここで、横方向のガウシアンカーネルは次のようになります。

位置：

![40.1611112995563.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/05cde4a2-c157-05c4-d019-453f5bf85482.png)




重み：
![41.1611113032657.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/d16f0fe0-56e4-a1de-5413-fa69de52d59b.png)





次のように置き換えることができます：

位置：

![42.1611113070103.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/c4c304bb-72c7-56c5-ec95-c2cbdd19c6f1.png)




重み：
![43.1611113109500.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/01e8d167-74f2-5a05-3424-3b603ee1e180.png)





ハードウェア実装による丸め誤差を無視すると、線形サンプリングによって得られた結果は、離散サンプリングによって得られた結果と同じになります。このようにして、（2ω+ 1）*（2ω+ 1）のガウシアンカーネルは（ω+ 1）*（ω+ 1）のガウシアンカーネルに簡略化されます。ピクセルをぼかすには、（2ω+ 2）回の乗算、2ω回の加算、時間計算量はO（ω）になります。

1024 * 1024サイズの画像の場合、33 * 33サイズの畳み込みカーネルを使用して画像全体にぼかし操作を実行する場合は、1024 * 1024 *17*2≈3.56kwのテクスチャ読み取りを実行する必要があります。以前の離散サンプリングより、テクスチャ読み取りの数が大幅に減少しました。

### 3.2 Unityの実装
上記のアルゴリズムに従って、Unityで実装します。

操作を容易にするために、選択したサンプリングポイントは直接中間点として設定されるため、得られた結果には前の結果と一定の誤差があります。基本的なプロセスは前の方法と同様であり、ぼかしアルゴリズムが実装されています。この時点では、ガウシアンカーネルを1つずつ読み取る必要はなく、ステップサイズを2に変更します。
```
float4 GaussianBlurLinearSampling(pixel_info pinfo, float sigma, float2 dir)
{
	float4 o = 0;
	float sum = 0;
	float2 uvOffset;
	float weight;

	for (int kernelStep = -KERNEL_SIZE / 2; kernelStep <= KERNEL_SIZE / 2; kernelStep += 2)
	{
		uvOffset = pinfo.uv;
		uvOffset.x += ((kernelStep + 0.5f) * pinfo.texelSize.x) * dir.x;
		uvOffset.y += ((kernelStep + 0.5f) * pinfo.texelSize.y) * dir.y;
		weight = gauss(kernelStep, sigma) + gauss(kernelStep + 1, sigma);
		o += tex2D(pinfo.tex, uvOffset) * weight;
		sum += weight;
	}
	o *= (1.0f / sum);
	return o;
}
```
効果は次のとおりです。
![26.1611122422195.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/653855/539afa31-457c-bace-f310-ff8c64696e77.png)



（35 * 35のガウシアンカーネルを使用）

### 3.3 まとめ

ガウシアンぼかしアルゴリズムによって実装されたぼかし効果の品質は良好ですが、パフォーマンスの点からみると思いどおりになりません。畳み込みアルゴリズムでは画像のサンプリング回数が多過ぎて、線形サンプリングアルゴリズムでは、2回のパスが必要のためDrawCallの数を2倍にしました。これらはすべて、ガウシアンぼかしアルゴリズムの実際のアプリケーションパフォーマンスに影響します。これにより、開発者はぼかしアルゴリズムを最適化し、ガウシアンぼかし効果の品質を可能な限り達成しながら、実際のニーズを満たすようにパフォーマンスを向上させることができます。次の章では、画像ぼかしアルゴリズムのいくつかの改善されたアルゴリズムを学習します。

DEMOリンク：https://github.com/UWA-MakeItSimple/Course-PostProcessingEffect/blob/main/Assets/Shaders/Blur/GaussianBlur/FirstDimensionGaussianBlur.shader

参考資料：

Gaussian Kernel Calculator：http://dev.theomader.com/gaussian-kernel-calculator/

デモアドレス：https://github.com/VestLee/VestImageEffect



https://zhuanlan.zhihu.com/p/110754637

https://zhuanlan.zhihu.com/p/58182228

https://blog.csdn.net/qq_36359022/article/details/80188873

https://zhuanlan.zhihu.com/p/120102048

http://www.ruanyifeng.com/blog/2017/12/image-and-wave-filters.html

https://blog.csdn.net/yangtrees/article/details/8740933

https://www.zhihu.com/question/54918332

https://zhuanlan.zhihu.com/p/106587858

https://yemi.me/2018/04/30/super-fast-blur-gpu/

https://lab.uwa4d.com/lab/5b65e90ad7f10a201ff92a5d

https://lab.uwa4d.com/lab/5b5cfe50d7f10a201fea2686

---
UWA Technologyは、モバイル/VRなど様々なゲーム開発者向け、**パフォーマンス分析**と**最適化ソリューション**及び**コンサルティングサービス**を提供している会社でございます。

今なら、UWA GOTローカルツールが**15日間に無償試用**できます！！
よければ、ぜひ！

UWA公式サイト：https://jp.uwa4d.com
UWA GOT OnlineレポートDemo：https://jp.uwa4d.com/u/got/demo.html
UWA公式ブログ：https://blog.jp.uwa4d.com