この記事はアドカレに参加しています。

luajitでガウシアンぼかし

前回の記事の延長線上にあることをやります。

ガウシアンぼかしとは

ガウシアンぼかしとは、ガウス関数で重みを付けたぼかし処理のことです。

ガウス関数は中心近いほど値が大きくなり、中心から離れるほど値が0に近くなっていきます。これをぼかし処理に適用すると、近くにあるピクセルほど重みが強くなり、遠くにあるピクセルの影響は微小になります。
　

アルゴリズム

先人の方達が素晴らしいアルゴリズムを共有してくださっています(take2
言語不問で実装する高速な移動平均処理
フィルタ処理の高速化アルゴリズム（縦横に処理を分ける）
gaussian フィルタ
画像フィルター：ガウスぼかし←これすき

コード

.anmとして使用できます。

--[[

blur.anm
	ガウシアンぼかしです。

]]

--track0:size,1,50,10,1

local ffi=require"ffi"

pcall(ffi.cdef,[[
    typedef struct Pixel_ {
      uint8_t b,g,r,a;
    } Pixel;
  ]])

local work,data,w,h=obj.getpixeldata"work",obj.getpixeldata()
local c=ffi.cast("Pixel*",data)
local r=ffi.cast("Pixel*",work)

local function f(a)--0~255の範囲で四捨五入する関数。
 if(a<0)then return 0
 elseif(a>255)then return 255
 else
  return math.floor(a+0.5)
 end
end

local t0=obj.track0/3
local ha=math.floor(obj.track0)
local ha2=ha*2+1
local ta={}
local total_ta=0
for i=1,ha2 do--taにガウス関数で重みを付ける
 local s=(i-1)-ha
 ta[i]=math.exp(-(s*s)/(2*t0*t0))/(math.sqrt(2*math.pi)*t0)
 total_ta=total_ta+ta[i]
end

for i=1,ha2 do--taの合計値が1.0になるようにする
 ta[i]=ta[i]/total_ta
end

local qt={}--ぼかす際の変数を入れる
local af=ffi.new("Pixel",{0})

--横方向
for y=0,h-1 do
 local ix=y*w
 local s=-ha

 for i=1,ha2 do
  if(s>=0 and s<w)then
   qt[i]=c[ix+s]
  else
   qt[i]=af
  end
  s=s+1
 end

 for x=0,w-1 do
  local to={0,0,0,0}

  for i=1,ha2 do
   to[1]=to[1]+qt[i].r*ta[i]
   to[2]=to[2]+qt[i].g*ta[i]
   to[3]=to[3]+qt[i].b*ta[i]
   to[4]=to[4]+qt[i].a*ta[i]
  end

  r[ix].r=f(to[1])
  r[ix].g=f(to[2])
  r[ix].b=f(to[3])
  r[ix].a=f(to[4])

  for i=1,ha2-1 do
   qt[i]=qt[i+1]
  end

  if(x+ha+1<w)then
   qt[ha2]=c[ix+ha+1]
  else qt[ha2]=af end
  ix=ix+1
 end
end

--縦方向
for x=0,w-1 do
 local ix=x
 local s=-ha

 for i=1,ha2 do
  if(s>=0 and s<h)then
   qt[i]=r[ix+s*w]
  else
   qt[i]=af
  end
  s=s+1
 end

 for y=0,h-1 do
  local s=-ha
  local to={0,0,0,0}

  for i=1,ha2 do
   to[1]=to[1]+qt[i].r*ta[i]
   to[2]=to[2]+qt[i].g*ta[i]
   to[3]=to[3]+qt[i].b*ta[i]
   to[4]=to[4]+qt[i].a*ta[i]
  end

  c[ix].r=f(to[1])
  c[ix].g=f(to[2])
  c[ix].b=f(to[3])
  c[ix].a=f(to[4])

  for i=1,ha2-1 do
   qt[i]=qt[i+1]
  end

  if(y+ha+1<h)then
   qt[ha2]=r[ix+(ha+1)*w]
  else qt[ha2]=af end
  ix=ix+w
 end
end

obj.putpixeldata(data)

おわりに

実はこれのgpu ver.をc++ ampで書いてみたのですが、上手く高速化できませんでした。gpu、難しいですね。