Date: Sat, 13 Dec 1997 15:51:52 +0800 (PST)
From: Orlando Andico <orly AT dilnet DOT upd DOT edu DOT ph>
To: "Gautam N. Lad" <gautam AT interlog DOT com>
cc: djgpp AT delorie DOT com
Subject: Re: Filtering Image Algorithm
In-Reply-To: <66supo$p8j$2@news.interlog.com>
Message-ID: <Pine.LNX.3.95.971213153935.8289B-100000@jewel.eee.upd.edu.ph>
MIME-Version: 1.0
Content-Type: TEXT/PLAIN; charset=US-ASCII
Precedence: bulk

On Sat, 13 Dec 1997, Gautam N. Lad wrote:

> Hi,
> Can someone tell me how to apply a filter (something like a 5x5 filter) to an 
> image?  I saw a site that shows it, but I don't quite understand it.  I just want to 
> modify a pixel's RGB values, and obtain a resulting RGB value, and draw that 
> coloured pixel.  I have access to 24-BIT video mode, so I don't have to care about 
> a palette, just do something like this:
> 
> putpixel(screen, 2, 4, makecol(r, g, b));	// this is put a pixel on screen at 2, 
> 			  	// 4 with colour returned by 
> 				// makecol(r,g,b);
> 
> Oh, I'm using DJGPP and Allegro 3.0 (Beta).
> Please help!

Basically what you want to do is apply a filter to your image with a 5x5
convolution kernel. There are lots of ways to do this, but the most
obvious (and slow..) is a nested for loop. 5x5 is not a good kernel size
IMHO (it's too large to implement with "hard wiring" and too small to
implement with a Fourier transform).

Anyway, here's a sample function which uses a 3x3 kernel (I'd send you my
Compleat Digital Image Processing Library (TM) source code, all written
with DJGPP for an image-recognition project a few years back, but I don't
have it on the computer I'm sending this mail from):

void convolve3x3 (int *image[], int *kernel[])
{
  int x, y, total;
  int *tmpimage[];  /* must be same size as *image[] */

  for (x = 1; x < WIDTH_OF_IMAGE; x++)
  {
    for (y = 1; y < HEIGHT_OF_IMAGE; y++)
    {
      tmpimage[x][y] = (image[x-1][y-1] * kernel[0][0] +
                        image[x][y-1] * kernel[1][0]   +
                        image[x+1][y-1] * kernel[2][0] +
                        image[x-1][y] * kernel[0][1] +
                        image[x][y] * kernel[1][1] +
                        image[x+1][y] * kernel[2][1] +
                        image[x-1][y+1] * kernel[0][2] +
                        image[x][y+1] * kernel[1][2] +
                        image[x+1][y+1] * kernel[2][2]) / 9;
    }
  }
  memcpy (image, tmpimage, sizeof(int) * WIDTH_OF_IMAGE *
          HEIGHT_OF_IMAGE);
}

The obvious problem with this code is that it leaves the borders of the
image untouched (a 3x3 will leave a 1-pixel border, a 5x5 will leave a
2-pixel border, etc).

The only way to avoid the border effect is to implement your filtering in
the frequency domain and use a Hamming window on your Fourier transform.
But, with the speed of current PC hardware, doing stuff with a Fourier is
a fool's errand: a 1024-point FFT takes about 0.1ms on a midrange Pentium,
operating on a 1024x1024 image would take you 2048 1024-point FFT's if you
use the vertical-followed-by-horizontal FFT instead of doing a brute force
(1024x1024) FFT. Even if you use the V-followed-by-F method, that's still
200 milliseconds for the entire image, which isn't fast enough for
realtime. 

According to my DSP professor, using the FFT to implement filtering starts
to get economical (compared to doing the convolution in the space domain)
when your kernel gets to the 10x10 size or so. If you're doing a 5x5,
that's too cumbersome for the manually-unrolled code I've shown above, so
you'd have to use nested for loops, which aren't terribly efficient (but
you can use -funroll-all-loops to unroll 'em anyway..)


----------------------------------------------------------------------
 L       L  Orlando Andico                             orly AT mozcom DOT com
 L L   L L  Mosaic Communications, Inc.  http://www2.mozcom.com/~orly/
 L L L L L  +63 2 848-2893 (office)
 L L L L L  +63 2 932-2385 (home)
  MosCom
----------------------------------------------------------------------