Dark Bit Factory & Gravity
PROGRAMMING => General coding questions => Topic started by: zawran on January 25, 2008
-
Do anyone have a working example, basic language variant if possible, of a gaussian blur function?
-
if you know how to make a simple box blur (simple average of pixels in a region), then apply the algorthim 3 times, each pass reuse the previous result => fast approximation to gaussian blur. The more passes you use, the closer you get to a perfect gaussian blur...
-
@zawran:
possible one of this nice sites will help you:
http://www.catenary.com/howto/gaussian.html (http://www.catenary.com/howto/gaussian.html)
http://www.jhlabs.com/ip/blurring.html (http://www.jhlabs.com/ip/blurring.html)
http://canopus.physik.uni-potsdam.de/~axm/bgsmooth.html (http://canopus.physik.uni-potsdam.de/~axm/bgsmooth.html)
http://en.wikipedia.org/wiki/Gaussian_blur (http://en.wikipedia.org/wiki/Gaussian_blur)
http://www.gamasutra.com/features/20010209/evans_01.htm (http://www.gamasutra.com/features/20010209/evans_01.htm)
http://www.pages.drexel.edu/~nk752/cannyTut2.html (http://www.pages.drexel.edu/~nk752/cannyTut2.html)
http://www.gamedev.net/community/forums/topic.asp?topic_id=427696 (http://www.gamedev.net/community/forums/topic.asp?topic_id=427696)
http://ciri.be/blog/?p=14 (http://ciri.be/blog/?p=14)
-
I found some blitzbasic code, which I have converted to BlitzMax. It seems to be working like it is supposed to, so I will use this for now.
' pixmap gaussian blur function
' based on blitzbasic code by Elias_T
' maskWidth & maskHeight should always be uneven numbers [3,4,7,9]
Function gaussianBlur:TPixmap(pMap:TPixmap,maskWidth:Int,maskHeight:Int)
Local tmp:TPixmap = CopyPixmap(pMap)
Local width:Int = pMap.width
Local height:Int = pMap.height
Local texel:Float[width,height,3]
Local result:Float[width,height,3]
Local maskData:Float[width*height]
Local x:Int,y:Int,ym:Int,xm:Int
Local cy:Float,cx:Float,rt:Float
Local r1:Float,g1:Float,b1:Float
Local rr:Float,gg:Float,bb:Float
Local mult:Float = 0.0
For x = 0 To width-1
For y = 0 To height-1
rgb = ReadPixel(pMap,x,y)
texel[x,y,0] = (rgb Shr 16) & 255
texel[x,y,1] = (rgb Shr 8) & 255
texel[x,y,2] = rgb & 255
Next
Next
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = Float(xm - (maskWidth - 1) / 2.0)
cy = Float(ym - (maskHeight - 1) / 2.0)
rt = cx*cx + cy*cy
mult :+ Exp(-0.35 * rt)
Next
Next
mult = 1.0 / mult
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - (maskWidth-1) / 2.0
cy = ym - (maskHeight-1) / 2.0
rt = cx*cx + cy*cy
maskData[ym * maskWidth + xm] = mult * Exp(-0.35 * rt)
Next
Next
For ym = 0 To height-1
For xm = 0 To width-1
rr = 0.0
gg = 0.0
bb = 0.0
For yy = 0 To maskHeight-1
For xx = 0 To maskWidth-1
If (xm+xx-Floor(maskWidth/2.0)<0) Or (ym+yy-Floor(maskHeight/2.0)<0) Or (xm+xx-Floor(maskWidth/2.0)>width-1) Or (ym+yy-Floor(maskHeight/2.0)>height-1) Then
r1 = 0.0
g1 = 0.0
b1 = 0.0
Else
r1 = texel[xm + xx - Floor(maskWidth/2.0), ym + yy - Floor(maskHeight/2.0),0]
g1 = texel[xm + xx - Floor(maskWidth/2.0), ym + yy - Floor(maskHeight/2.0),1]
b1 = texel[xm + xx - Floor(maskWidth/2.0), ym + yy - Floor(maskHeight/2.0),2]
End If
rr :+ r1 * maskData[xx + yy * maskWidth]
gg :+ g1 * maskData[xx + yy * maskWidth]
bb :+ b1 * maskData[xx + yy * maskWidth]
Next
Next
result[xm,ym,0] = rr
result[xm,ym,1] = gg
result[xm,ym,2] = bb
Next
Next
For x = 0 To width-1
For y = 0 To height-1
WritePixel(tmp,x,y,Int(result[x,y,0]) Shl 16 + Int(result[x,y,1]) Shl 8 + Int(result[x,y,2]))
Next
Next
Return tmp
End Function
Its not the fastest thing around, and I should probably try and convert it into using pointers instead of read-/write-pixels. And it just might be faster to do a regular box blur a couple of times with almost the same end result. But for now, this will do fine for what I needed it for, where it isn't too critical if its fast or not, at least not at this point.
-
first thing to optimize in your code is the reading of pixels. Take advantage of the fact that a gaussian blur is separable, which means you can process lines and columns separately. I won't go in the math details here why it is separable but it is, so do the following:
1) blur the lines first (blur on x axis)
2) take the result and feed it into the next blur for the columns (blur on y axis)
this will reduce the pixel access by a great factor and definitely improve the speed of the algorithm.
-
Which reduces the bluring to one add, one sub and one div/lookup per pixel and per pass.
-
http://blackpawn.com/texts/blur/default.html
fastest one I've ever used. :*)
-
Hell yes!
Nice one Rel, that's very useful. I came home from work to find this forum full of lovely useful replies today :D
-
relsoft: ;) this is exactly what I meant by one add + one sub + one div/lookUp per pixel per pass.
-
@relsoft, thanks for the link, I will take a look.
-
I did one small change to the blur code, which knocked a lot of the processing time of it. Don't know, why I didn't do that in the first place:
Graphics 1024,768
' test image = 256x256 pixels
' first version = 4956ms
' second version = 108ms
Cls
Local pixm:TPixmap = LoadPixmap("testimage.png")
Local time:Int = MilliSecs()
Local pixm2:TPixmap = gaussianBlur(pixm,9,9)
Print MilliSecs()-time
DrawPixmap(pixm,0,0)
DrawPixmap(pixm2,256,0)
Flip
WaitKey()
End
' pixmap gaussian blur function
' based on blitzbasic code by Elias_T
' maskWidth & maskHeight should always be uneven numbers [3,5,7,9]
Function gaussianBlur:TPixmap(pMap:TPixmap,maskWidth:Int,maskHeight:Int)
Local tmp:TPixmap = CopyPixmap(pMap)
Local width:Int = pMap.width
Local height:Int = pMap.height
Local texel:Float[width,height,3]
Local result:Float[width,height,3]
Local maskData:Float[width*height]
Local x:Int,y:Int,ym:Int,xm:Int
Local cy:Float,cx:Float,rt:Float
Local r1:Float,g1:Float,b1:Float
Local rr:Float,gg:Float,bb:Float
Local mult:Float = 0.0
For x = 0 To width-1
For y = 0 To height-1
rgb = ReadPixel(pMap,x,y)
texel[x,y,0] = (rgb Shr 16) & 255
texel[x,y,1] = (rgb Shr 8) & 255
texel[x,y,2] = rgb & 255
Next
Next
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = Float(xm - (maskWidth - 1) / 2.0)
cy = Float(ym - (maskHeight - 1) / 2.0)
rt = cx*cx + cy*cy
mult :+ Exp(-0.35 * rt)
Next
Next
mult = 1.0 / mult
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - (maskWidth-1) / 2.0
cy = ym - (maskHeight-1) / 2.0
rt = cx*cx + cy*cy
maskData[ym * maskWidth + xm] = mult * Exp(-0.35 * rt)
Next
Next
Local mhhalf:Int = Floor(maskHeight/2.0)
Local mwhalf:Int = Floor(maskWidth/2.0)
For ym = 0 To height-1
For xm = 0 To width-1
rr = 0.0
gg = 0.0
bb = 0.0
For yy = 0 To maskHeight-1
For xx = 0 To maskWidth-1
If (xm+xx-mwhalf<0) Or (ym+yy-mhhalf<0) Or (xm+xx-mwhalf>width-1) Or (ym+yy-mhhalf>height-1) Then
r1 = 0.0
g1 = 0.0
b1 = 0.0
Else
r1 = texel[xm + xx - mwhalf, ym + yy - mhhalf,0]
g1 = texel[xm + xx - mwhalf, ym + yy - mhhalf,1]
b1 = texel[xm + xx - mwhalf, ym + yy - mhhalf,2]
End If
rr :+ r1 * maskData[xx + yy * maskWidth]
gg :+ g1 * maskData[xx + yy * maskWidth]
bb :+ b1 * maskData[xx + yy * maskWidth]
Next
Next
result[xm,ym,0] = rr
result[xm,ym,1] = gg
result[xm,ym,2] = bb
Next
Next
For x = 0 To width-1
For y = 0 To height-1
WritePixel(tmp,x,y,Int(result[x,y,0]) Shl 16 + Int(result[x,y,1]) Shl 8 + Int(result[x,y,2]))
Next
Next
Return tmp
End Function
Difference is that I calculate that / 2.0 value before the loop instead of each time within the loop. And it went from almost 5000ms to about 110ms, so thats a big saving just for that.
-
since maskwidth and maskheight are ints you could use shift right instead of divide, also since xm and ym are ints you could make cx and cy ints and convert to float in the line:
rt = float(cx*cx + cy*cy)
so you only have one int to float in each of the loops.
Another thing I would do is something like this:
y_index=ym-mhhalf
For yy = 0 To maskHeight-1
If (y_index>=0) and (y_index<height) Then
mask_index=yy * maskwidth
x_index=xm-mwhalf
For xx = 0 To maskWidth-1
If (x_index>=0) and (x_index<width) Then
rr :+ texel[x_index, y_index,0] * maskData[mask_index]
gg :+ texel[x_index, y_index,1] * maskData[mask_index]
bb :+ texel[x_index, y_index,2] * maskData[mask_index]
End If
mask_index :+ 1
x_index :+ 1
Next
End If
y_index :+ 1
Next
to take as much of the calculations as possible out of the loops (although no idea if i have that all right as i just edited it in notepad), it should even be possible to remove the If from the inner loop completely.
-
By making the changes you suggested, it came down another 20ms.
Graphics 1024,768
' test image = 256x256 pixels
' first version = 4956ms
' second version = 108ms
' third version = 82ms
Cls
Local pixm:TPixmap = LoadPixmap("testimage.png")
Local time:Int = MilliSecs()
Local pixm2:TPixmap = gaussianBlur(pixm,9,9)
Print "processtime was: "+(MilliSecs()-time)+"ms"
DrawPixmap(pixm,0,0)
DrawPixmap(pixm2,256,0)
Flip
WaitKey()
End
' pixmap gaussian blur function
' based on blitzbasic code by Elias_T
' maskWidth & maskHeight should always be uneven numbers [3,5,7,9]
Function gaussianBlur:TPixmap(pMap:TPixmap,maskWidth:Int,maskHeight:Int)
Local tmp:TPixmap = CopyPixmap(pMap)
Local width:Int = pMap.width
Local height:Int = pMap.height
Local texel:Float[width,height,3]
Local result:Float[width,height,3]
Local maskData:Float[width*height]
Local x:Int,y:Int,ym:Int,xm:Int
Local cy:Int,cx:Int,rt:Float
Local r1:Float,g1:Float,b1:Float
Local rr:Float,gg:Float,bb:Float
Local mult:Float = 0.0
Local mwhalf:Int,mhhalf:Int
Local y_index:Int,x_index:Int,mask_index:Int
mwhalf = Floor((maskWidth-1) / 2)
mhhalf = Floor((maskHeight-1) / 2)
For x = 0 To width-1
For y = 0 To height-1
rgb = ReadPixel(pMap,x,y)
texel[x,y,0] = (rgb Shr 16) & 255
texel[x,y,1] = (rgb Shr 8) & 255
texel[x,y,2] = rgb & 255
Next
Next
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - mwhalf
cy = ym - mhhalf
rt = Float(cx*cx + cy*cy)
mult :+ Exp(-0.35 * rt)
Next
Next
mult = 1.0 / mult
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - mwhalf
cy = ym - mhhalf
rt = Float(cx*cx + cy*cy)
maskData[ym * maskWidth + xm] = mult * Exp(-0.35 * rt)
Next
Next
mhhalf = Floor(maskHeight/2.0)
mwhalf = Floor(maskWidth/2.0)
For ym = 0 To height-1
For xm = 0 To width-1
rr = 0.0
gg = 0.0
bb = 0.0
y_index = ym-mhhalf
For yy = 0 To maskHeight-1
If y_index => 0 And y_index < height Then
mask_index = yy * maskwidth
x_index = xm-mwhalf
For xx = 0 To maskWidth-1
If x_index => 0 And x_index < width Then
rr :+ texel[x_index, y_index,0] * maskData[mask_index]
gg :+ texel[x_index, y_index,1] * maskData[mask_index]
bb :+ texel[x_index, y_index,2] * maskData[mask_index]
End If
mask_index :+ 1
x_index :+ 1
Next
End If
y_index :+ 1
Next
result[xm,ym,0] = rr
result[xm,ym,1] = gg
result[xm,ym,2] = bb
Next
Next
For x = 0 To width-1
For y = 0 To height-1
WritePixel(tmp,x,y,Int(result[x,y,0]) Shl 16 + Int(result[x,y,1]) Shl 8 + Int(result[x,y,2]))
Next
Next
Return tmp
End Function
I tried removing the if/then check in the inner loop, but it crashes, so its still needed. But at least it took off a good chunk of the process time, about 24%, which isn't half bad. Thanks for the tips. :)
-
Sorry, yeah it'd need a few other changes before you can remove the If but good to hear the other changes made some difference.
-
I went ahead and changed it from read/write pixels to pointers, and shave another 8ms of the process time, not a huge difference, but still :)
Graphics 1024,768
' test image = 256x256 pixels
' first version = 4956ms
' second version = 108ms
' third version = 82ms
' forth version = 74ms
Cls
Local pixm:TPixmap = LoadPixmap("testimage.png")
Local time:Int = MilliSecs()
Local pixm2:TPixmap = gaussianBlur(pixm,9,9)
Print "processtime was: "+(MilliSecs()-time)+"ms"
DrawPixmap(pixm,0,0)
DrawPixmap(pixm2,256,0)
Flip
WaitKey()
End
' pixmap gaussian blur function
' based on blitzbasic code by Elias_T
' maskWidth & maskHeight should always be uneven numbers [3,5,7,9]
Function gaussianBlur:TPixmap(pMap:TPixmap,maskWidth:Int,maskHeight:Int)
Local tmp:TPixmap = CopyPixmap(pMap)
Local pMapPtr:Byte Ptr = PixmapPixelPtr(pMap,0,0)
Local tmpPtr:Byte Ptr = PixmapPixelPtr(tmp,0,0)
Local width:Int = pMap.width
Local height:Int = pMap.height
Local texel:Float[width,height,3]
Local result:Float[width,height,3]
Local maskData:Float[width*height]
Local x:Int,y:Int,ym:Int,xm:Int
Local cy:Int,cx:Int,rt:Float
Local r1:Float,g1:Float,b1:Float
Local rr:Float,gg:Float,bb:Float
Local mult:Float = 0.0
Local mwhalf:Int,mhhalf:Int
Local y_index:Int,x_index:Int,mask_index:Int
mwhalf = Floor((maskWidth-1) / 2)
mhhalf = Floor((maskHeight-1) / 2)
For x = 0 To width-1
For y = 0 To height-1
texel[x,y,0] = pMapPtr[x*4+y*4*width]
texel[x,y,1] = pMapPtr[x*4+y*4*width+1]
texel[x,y,2] = pMapPtr[x*4+y*4*width+2]
Next
Next
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - mwhalf
cy = ym - mhhalf
rt = Float(cx*cx + cy*cy)
mult :+ Exp(-0.35 * rt)
Next
Next
mult = 1.0 / mult
For ym = 0 To maskHeight-1
For xm = 0 To maskWidth-1
cx = xm - mwhalf
cy = ym - mhhalf
rt = Float(cx*cx + cy*cy)
maskData[ym * maskWidth + xm] = mult * Exp(-0.35 * rt)
Next
Next
mhhalf = Floor(maskHeight/2.0)
mwhalf = Floor(maskWidth/2.0)
For ym = 0 To height-1
For xm = 0 To width-1
rr = 0.0
gg = 0.0
bb = 0.0
y_index = ym-mhhalf
For yy = 0 To maskHeight-1
If y_index => 0 And y_index < height Then
mask_index = yy * maskwidth
x_index = xm-mwhalf
For xx = 0 To maskWidth-1
If x_index => 0 And x_index < width Then
rr :+ texel[x_index, y_index,0] * maskData[mask_index]
gg :+ texel[x_index, y_index,1] * maskData[mask_index]
bb :+ texel[x_index, y_index,2] * maskData[mask_index]
End If
mask_index :+ 1
x_index :+ 1
Next
End If
y_index :+ 1
Next
result[xm,ym,0] = rr
result[xm,ym,1] = gg
result[xm,ym,2] = bb
Next
Next
For x = 0 To width-1
For y = 0 To height-1
tmpPtr[x*4+y*4*width] = Int(result[x,y,0])
tmpPtr[x*4+y*4*width+1] = Int(result[x,y,1])
tmpPtr[x*4+y*4*width+2] = Int(result[x,y,2])
Next
Next
Return tmp
End Function
I am going to move on to other filters now. I am slowly working my way through a lot of the different filters used by a lot of art packages. Eventually they are to be used in my own paint program, but its a long way ahead and its taking its time. But I will get there eventually. Thanks for the help so far. If I run into a wall on some of the other filters, I will be back posing what I got to see if it can be fixed/improved.