大津の手法による閾値自動決定アルゴリズムを用いた二値化をGPUを用いて実装しました
Boltzmanです. 前回の記事から1か月以上が経過しました. 極めて難航しましたが, 進捗が上がりましたので記事にします. (遅すぎ)
前回, 大津の手法による自動閾値決定2値化アルゴリズムをCPUで実装してあまりにも重くなったため, GPUを用いて, 具体的にはComputeShaderを用いて実装しました.
実行環境は私のノートPC(Intel Core i5-5200U/RAM 4.00GB/Intel HD Graphics 5500)です. 決して良いとは言えないスペック.
実行結果は
CPUのみで実装 (v0.1.0):MonoBehaviourが平均330msくらい
GPUを取り入れて実装(v0.2.3):Gfx.WaitForPresentが平均16~20msくらい?(ガタつきが多い)
改めてこうしてみるとCPUのみで実装したほうが理不尽に重いので, どこかでコーディングミスをしている可能性が大いにあります.
また大津の手法による自動閾値決定により2値化画像がどのように変化するかが以下です.
上を見るとわかる通り, マーカーは濃淡の変化が顕著になるように作られているためか, あまり見栄えは変わりません.
しかしその他の背景画像は大津の手法を実装したときのほうが見やすくなっています.
以下にComputeShaderのコードを記載します
| #pragma kernel CreateHistogram | |
| #pragma kernel MergeHistogram | |
| #pragma kernel OtsuMethod | |
| #pragma kernel Threshold | |
| Texture2D<float4> Input; | |
| float thresh; | |
| RWTexture2D<float4> Result; | |
| RWStructuredBuffer<uint> histogram : register(u0); | |
| RWStructuredBuffer<uint> hbuffer; | |
| groupshared uint sbuffer[256]; | |
| RWStructuredBuffer<float2> otsubuffer; //x:thresh, y:variable | |
| RWStructuredBuffer<float4> data;//x:num1, y:num2, z:sum1, w:sum2 | |
| [numthreads(32,32,1)] | |
| void CreateHistogram (uint3 thid : SV_DispatchThreadID, uint gridx : SV_GroupIndex, uint3 grid : SV_GroupID) | |
| { | |
| float w, h; | |
| int i=0; | |
| uint grayscale = 0; | |
| uint value=0; | |
| uint id; | |
| w = 0; h = 0; | |
| Input.GetDimensions(w, h); | |
| if(gridx == 256) for(i=0; i<256; i++) sbuffer[i] = 0; | |
| GroupMemoryBarrierWithGroupSync(); | |
| grayscale = (uint)floor((Input[thid.xy].r+Input[thid.xy].g+Input[thid.xy].b)/3*255); | |
| value = step( distance(thid.xy, float2(w/2, h/2)), w*26/63.5 );//魚眼レンズ使用時に使用 | |
| if((thid.x < (uint)w)&&(thid.y < (uint)h)){InterlockedAdd(sbuffer[grayscale], /*value*/1); } | |
| GroupMemoryBarrierWithGroupSync(); | |
| if(gridx < 256) InterlockedAdd(histogram[gridx], sbuffer[gridx]); | |
| } | |
| [numthreads(256,1,1)] | |
| void OtsuMethod(uint3 id : SV_DispatchThreadID){ //id = dim of histogram | |
| float w=0; | |
| float2 w2 = float2(0, 0); | |
| float4 w4 = float4(0, 0, 0, 0); | |
| uint i=0; //ループ用 | |
| uint thid; //threadID格納用 | |
| int offset; //Scan用変数 | |
| int ai, bi; //Scan用変数 | |
| float m1, m2; //平均 | |
| float num1, num2; | |
| m1=0; | |
| m2=0; | |
| thid = id.x; | |
| data[thid].x = (float)histogram[thid]; | |
| data[255-thid].y = (float)histogram[thid]; | |
| data[thid].z = (float)histogram[thid] * thid; | |
| data[255-thid].w = (float)histogram[thid] * thid; | |
| AllMemoryBarrierWithGroupSync(); | |
| //Scanでクラス内要素数とクラス内和を算出 | |
| //上にsweep(Reduce) | |
| offset = 1; | |
| for(i=256>>1; i>0; i>>=1){ | |
| AllMemoryBarrierWithGroupSync(); | |
| if(thid<i){ | |
| ai = offset*(2*thid+1)-1; | |
| bi = offset*(2*thid+2)-1; | |
| data[bi]+=data[ai]; | |
| } | |
| offset *= 2; | |
| } | |
| if(thid == 0){data[255]=0;} | |
| //下にsweep | |
| for(i=1; i<256; i*=2){ | |
| offset >>= 1; | |
| AllMemoryBarrierWithGroupSync(); | |
| if(thid<i){ | |
| ai = offset*(2*thid+1)-1; | |
| bi = offset*(2*thid+2)-1; | |
| w4 = data[ai]; | |
| data[ai] = data[bi]; | |
| data[bi] += w4; | |
| } | |
| } | |
| AllMemoryBarrierWithGroupSync(); | |
| //クラス内平均を算出 | |
| ai = thid; | |
| m1 = data[ai].z / data[ai].x; | |
| num1 = data[ai].x; | |
| AllMemoryBarrierWithGroupSync(); | |
| bi = 256-thid; | |
| m2 = data[bi].w / data[bi].y; | |
| num2 = data[bi].y; | |
| AllMemoryBarrierWithGroupSync(); | |
| //分離度を算出 | |
| w = num1 * num2 * (m1 - m2) * (m1 - m2); | |
| //0割をはじく | |
| if(isfinite(w)){ | |
| otsubuffer[thid] = float2((float)thid, w); | |
| } | |
| else{ | |
| otsubuffer[thid] = float2((float)thid, 0); | |
| } | |
| AllMemoryBarrierWithGroupSync(); | |
| //Scanで最大値算出 | |
| //上にsweep(Reduce) | |
| offset = 1; | |
| for(i = 256>>1; i>0; i >>= 1){ | |
| AllMemoryBarrierWithGroupSync(); | |
| if(thid<i){ | |
| ai = offset*(2*thid+1)-1; | |
| bi = offset*(2*thid+2)-1; | |
| otsubuffer[bi] = lerp(otsubuffer[ai], otsubuffer[bi], step(otsubuffer[ai].y, otsubuffer[bi].y)); | |
| } | |
| offset *= 2; | |
| } | |
| } | |
| [numthreads(8,8,1)] | |
| void Threshold (uint3 id : SV_DispatchThreadID) | |
| { | |
| //uint w, h; | |
| thresh = (float)otsubuffer[255].x/255; | |
| Result[id.xy] = lerp(float4(1, 1, 1, 0), float4(0, 0, 0, 0), step( (Input[id.xy].r+Input[id.xy].g+Input[id.xy].b)/3, thresh) ); | |
| //ケラレ範囲確認用 | |
| //Input.GetDimensions(w, h); | |
| //if(distance(id.xy, float2(w/2, h/2)) < w*26/63.5) Result[id.xy] = float4(1, 0, 0, 0); | |
| //元映像確認用 | |
| //Result[id.xy] = Input[id.xy]; | |
| } |
by 前回からの変更点としては,
・大津の手法の計算と最大値算出にScanを使用した
(参考URL:http://http.developer.nvidia.com/GPUGems3/gpugems3_ch39.html)
・ヒストグラムの生成時にSharedMemoryを利用してアクセスの集中を軽減した
(参考URL:https://devblogs.nvidia.com/parallelforall/gpu-pro-tip-fast-histograms-using-shared-atomics-maxwell/)
といった感じです.
参考元はNvidiaなのでCUDAで書いてありますが, 実装はHLSLで行いました.
いずれもまだまだ最適には程遠いかと思います.
もう少し実行速度が安定するまでコードを練るか, それともGPUを用いたラベリングの実装に移るか少し迷っています.
あと魚眼レンズはまだまだ遠いので一度置いておこうと思います.
そもそも魚眼レンズから取得したゆがんだ画像からマーカーを検出できるのかという問題を一切考えていないのでできる保証はないです.
次回はどちらになるでしょうか…
それではまた.