request: simple spatial blur code

[E-Male]

i got some problems getting started on a filter
so i ask for another simple code sample:
spatial 3x3 blur
every pixel gets componentwise replaced by the average of itself and the 8 spatial neighboars
if someone could code this for yv12 and without any mmx and similar stuf
it would be a great help
thx
e-male

[MfA]

Im not gonna do a complete filter, just look at the simpleresize example code for a complete filter you can build on.


#define CLIP(x, low, high) ((x) < (low) ? (low) : ((x) > (high) ? (high) : (x)))

void blur(unsigned char *src, unsigned char *dst, int height, int width, int pitchSrc, int pitchDst)
{
for (int y = 0; y < height; y++)
for (int x = 0; x < width; x++)
{
int sum = 0;

for (int i = -1; i < 2; i++)
for (int j = -1; j < 2; j++)
sum += src[CLIP(y + i, 0, height) * pitchSrc + CLIP(x + j, 0, width)];

sum /= 9;

dst[y * pitchDst + x] = sum;
}
}


Just pass it the correct pointers and height/width/pitch data you get from GetReadPtr, call it for Y/U/V seperately.

This is a inefficient code, but very simple. It is inefficient for a couple of reasons :

Handling border effects inside the main loop is inefficient, the CLIPping, it is better to either do edge extension or handle them with seperate code.

A 3x3 uniform blur is seperable, that means you can first do it vertically and then horizontally.

A moving average can be computed in O(1), and this code is O(N) ... not really relevant for a 3x3 filter, but any larger and it starts to be.

I think you should learn to walk before you try to run BTW :)

[sh0dan]

First: a 3x3 blur in MMX is a b*tch, so if you don't know it, I'd start with something simpler.

Simple-not-tested-at-all code:

PVideoFrame MyFilter::GetFrame(int n, IScriptEnvironment* env) {
// Blur cannot be done inplace.
PVideoFrame src = child->GetFrame(n, env);
PVideoFrame dst = env->NewVideoFrame(vi);

BlurMe_Planar(dst->GetWritePtr(PLANAR_Y),dst->GetPitch(PLANAR_Y),
src->GetReadPtr(PLANAR_Y),src->GetPitch(PLANAR_Y), vi.width, vi.height, env);

BlurMe_Planar(dst->GetWritePtr(PLANAR_U),dst->GetPitch(PLANAR_U),
src->GetReadPtr(PLANAR_U),src->GetPitch(PLANAR_U), src->GetWidth(PLANAR_U), src->GetHeight(PLANAR_U), env);

BlurMe_Planar(dst->GetWritePtr(PLANAR_V),dst->GetPitch(PLANAR_V),
src->GetReadPtr(PLANAR_V),src->GetPitch(PLANAR_V), , src->GetWidth(PLANAR_V), src->GetHeight(PLANAR_V), env);

return dst;
}

void MyFilter::BlurMe_Planar(BYTE* dst_plane, int dst_pitch,
BYTE* src_plane, int src_pitch, int width, int height, IScriptEnvironment* env) {

// Do special case for top-line. For now we just copy
env->BitBlt(dst_plane, dst_pitch, src_plane, src_pitch, width,1);

src_plane+=src_pitch;
dst_plane+=dst_pitch;
for (int y = 1; y < (height-1); y++) {
dstplane[0] = src_plane[0]; // Special case for x = 0 - just copy.
for (int x = 1; x < (width-1); x++) {
int pixel = 0;
for (int i=-1; i<2; i++) {
for (int o=-1; o<2; o++) {
pixel+=src_plane[x+(i*src_pitch)+o];
}
}
dst_plane[x] = pixel/9;
}
dst_plane[width-1] = src_plane[width-1];
src_plane+=src_pitch;
dst_plane+=dst_pitch;
}
// Do special case for bottom-line. For now we just copy
env->BitBlt(dst_plane, dst_pitch, src_plane, src_pitch, width,1);
}

This is the basic algorithm. Problems:
1) VERY slow.
2) No coefficients.
3) Border pixels not treated properly.
4) Nearly impossible to MMX-optimize.
5) Rounding problems at the division.

You should do vertical and horizontal filtering separately. It is much easier to optimize, since you can avoid loops, or excessive unrolling. Also the division should be replaced by a multiplication + bitshift.

For examples of this have a look at focus.cpp (http://cvs.sourceforge.net/viewcvs.py/avisynth2/avisynth/src/filters/focus.cpp?view=auto&rev=1.4&sortby=date&only_with_tag=MAIN) (the internal Blur() filter) - it implements a 3x3 blur with optional weights, and contains MMX versions.

There are however still some border issues with Blur(), and strange luma flickering has been reported (https://sourceforge.net/tracker/?func=detail&aid=748265&group_id=57023&atid=482673), when changing the blur amount between frames.

Edit: Found another problem. :)

[Didée]

By the way:

Shouldn't SpatialSoften() support YV12 also? It's the only core filter available for thresholded smoothing, and does only YUY2 ...

- Didée

[sh0dan]

@Didée: SpatialSoften() bl*ws, IMO. It's incredible slow, and doesn't really serve much purpose. There are a lot of better filters for that kind of stuff.

I never really spent any time on it, as filters like this should be external plugins, IMO.

[Mug Funky]

hmm. blur filters are quite generic. i'm not sure if something like a gaussian blur should be a plugin or not, as they are extremely useful.

of course, it can be scripted quite adequately with masktools (or chained blur(1) if absolutely necessary).

i'd love to see a plugin of the "box filter" that comes in virtualdub (but for yv12 of course). it's extremely fast even in RGB.

*sigh* i wish i was a coder.

[sh0dan]

w = width()
h = height()
BicubicResize(w/4,h/4)
BicubicResize(w,h)

;)

PS. Yes - I know it isn't quite the same. ;)
PPS. A generic blur would be fine. A threshold-based selective spacial blur, I'm not sure.

[Manao]

Instead of simply divide by 9 ( pix = pix / 9; ), do the following :pix = (2 * pix + 9) / 18;This will allow the average luma of the picture to be unchanged.I'd love to see a plugin of the "box filter" that comes in virtualdub (but for yv12 of course).You can do it with YV12Concolution ( a gaussian blur is separable, hence YV12Convolution can handle it ), it *should* be as fast.

[Mug Funky]

manao - that's what i've done so far. it's not as fast as virtualdub's filters (in spite of not being RGB).

of course, it might be my implementation. basically i generate a string to feed to YV12convo, which is a simple sin^2 function of variable period.

it's not the same as a genuine low pass filter as it does not ring (i don't have the maths for something that fancy) but rather is a single cycle of a sine curve.

however, it is quite fast. i just wish there was an easier way :)

[Leak]

Originally posted by sh0dan
First: a 3x3 blur in MMX is a b*tch, so if you don't know it, I'd start with something simpler.

Oh great... I'm just working on a 3x3-15x15 averaging blur for my BlendBob filter, and now you tell me... :D

Actually, it's already done in C++ and I think it's quite fast (though I'm only using it for YV12 chroma, so at double the resolution it's likely that you'll take a larger performance hit). For 15x15, it uses about 5 times the memory (but only 12 additions in total) that you need for a single plane, as it keeps several buffers for several stages of addition - but as it is, the code should be quite straightforward to port to MMX or SSE assembler code, which is what is next on my list.

(snip code)

This is the basic algorithm. Problems:
1) VERY slow.
2) No coefficients.


Of course, if I'd add coefficients the optimization of using several buffers and trading memory usage for speed goes straight out the window... :(

Still, a straight area averaging is sufficient for my purposes.

3) Border pixels not treated properly.


It took me some tinkering, but I finally got that working properly.

You should do vertical and horizontal filtering separately. It is much easier to optimize, since you can avoid loops, or excessive unrolling. Also the division should be replaced by a multiplication + bitshift.


That's my experience also. I probably won't replace the division for the C++ version, but surely for the MMX/SSE version.

Edit: Found another problem. :)

Ain't life grand? ;)

By the way, I've noticed a few strange things in the AviSynth code when looking for the YUV2RGB conversion used by AviSynth:

In color.cpp's Color::YUV2RGB(...), there's this code:


*b = ((scaled_y + (u-128) * cbu) + 2^15) >> 16; // blue
*g = ((scaled_y - (u-128) * cgu - (v-128) * cgv) + 2^15) >> 16; // green
*r = ((scaled_y + (v-128) * crv) + 2^15) >> 16; // red


I don't think that those instances of "2^15" do what whoever wrote that code meant... but does it ever get used anyway?

Also, in layer.cpp's ColorKeyMask::GetFrame(...), the non-MMX code looks like this:


for (int y=0; y<vi.height; y++) {
for (int x=0; x<rowsize; x+=4) {
if (IsClose(pf[x],B,tol) && IsClose(pf[x+1],B,tol) && IsClose(pf[x+2],R,tol))
pf[x+3]=0;
}
pf += pitch;
}


I guess one of the B's in the IsClose calls really wants to be a G...

np: Radiohead - Subterranean Homesick Alien (OK Computer)

[sh0dan]

ColorYUV only uses "YUV2RGB" for printing out some debugging information, when it is compiled in debug mode.

Regarding ColorKeyMask, it is wrong. However it is "thank god" only Non-MMX fallback code, which has probably never been used, since most of AviSynth doesn't work on non-MMX machines. The calculation if "R" was also wrong, if alpha was already present. Both issues have been fixed.

Unless it is completely unavoidable, always replace divisions - one division per pixel is quite a large overhead.


Thanks!