Wesley Faler is a senior Manufacturing Systems Engineering student at GMI Engineering & Management Institute in Flint, MI. He has worked with C for 2 years now, and interests include user interface coding, AI, and scientific computing. He is presently on the programming team that is writing the GMI CIM Lab system code. All comments should be sent to 1409 New York Ave., Flint, MI 48506.
Several months ago, a project led me to a particularly interesting image manipulation algorithm called convolution, which is capable of finding edges, horizontal lines, vertical lines, and even diagonal lines in an image.
Convolution is a reasonably simple process whereby we interpret an image's pixels in terms of how they fit with their neighbors. In other words, the image's pixels after convolution represent some function of their original neighbors. The following pseudocode makes the process a little clearer:
For every pixel in the original image, redraw the pixel with a color that is a function of the colors of the original pixel's neighbors.
Continue the first step until the original image has been completely "convolved" into the new image. This algorithm doesn't choose the function involving the original neighbors. I chose a simple and effective function: a weighted sum, including the original pixel.
In a simple, one-dimentional example, a row of pixels is represented by (1 2 3 4), with each integer defining a pixel's color. The pixel colored (2) is the first that can be redrawn since it is the first (left to right) with left and right neighbors. Assume we apply a "convolution matrix" of (-5 0 6) to the pixel and its neighbors. The pixel's new color will be
(1)(-5)+(2)(0)+(3)(6)or (13). The convolved value for the original pixel colored (3) will be
(2)(-5)+(3)(0)+(4)(6)or (14). The two end pixels (1) and (4) cannot be convolved since there are no left or right points to which to apply the convolution matrix. They are assigned a value of (0). The new row of pixels is now (0 13 14 0).The two-dimentional case is similar, involving 2-D regions of pixels and a 2-D convolution matrix.
Choosing A Convolution Matrix
How you choose the convolution matrix determines the final picture's appearance. The program in Listing 1 accepts a file named matrix.dat that has the following format:
<matrix width> <matrix height> <row 1,col 1 factor> ... <row 1,col w factor> ... ... ... <row h,col 1 factor> ... <row h,col w factor>Note that both the width and height must be odd, since the matrix must have a center point.Some sample matrices will give you a feel for what convolution can do. For example, the matrix file
3 3 -1 0 1 -1 0 1 -1 0 1will light pixels that were on only the left edge of a dark-to-light transition on the screen (see Figure 1) . Reversing the 1s and -1s will detect a light-to-dark transition (i.e. lighting the "right" side of the object). Likewise, the matrix file:
3 3 1 1 1 0 0 0 -1 -1 -1will detect a dark-to-light transition from top to bottom (Figure 2) . Switching the 1s and -1s here allows you to detect top-side lines. To generate a picture of only an object's edges (like the image on this issue's cover), use the matrix file:
3 3 -1 -1 -1 -1 8 -1 -1 -1 -1Diagonal lines can be detected by using a matrix such as
3 3 0 1 0 -1 0 1 0 -1 0In choosing your own matrix, observe the following simple rules for best results:
Make sure the sum of the factors in the matrix is positive or zero. A zero sum will eliminate extraneous pixels, while a small net positive sum can seem to fill in missing pixels. A negative sum usually filters out too many pixels (especially on a monochrome screen).
Use a zero in matrix locations where you are ambivalent about the color.
Use negative numbers for locations where you wish not to have color.
Use positive numbers for locations where you want a non-zero value.
The Program
Listing 1 contains a program that accepts images in the CUT file format used by the Dr. Halo products. This format was easy to use with a hand scanner, but you may substitute any other means of getting a picture onto the screen. If you should replace the present routine for loading images (load_cut (...) ), be sure to set the global variables minx, miny, maxx, maxy to the image's location on the screen, but do not change the display page. To keep the convolution code reasonably clean, the program also requires that you leave at least <matrix height+1>/2 and <matrix width+1 >/2 margin on the top and left side of the image, respectively.
CONVOLVE.C was written in Turbo C 2.0 and the graphics calls should be fairly portable to other Cs. I frequently used pointers to minimize the execution time that several hundred thousand array accesses would normally take. The program could be further optimized by forcing it to ignore matrix entries of zero.
When you run the program, it asks for the name of a CUT file to load and some information about image cleanup. If the file is not in the same directory, enter the path and the filename (including the .CUT extension).
The image cleanup option will delete stray pixels before convolution. A single pixel with zero neighbors usually does not belong to a surface of interest. You can delete such pixels before convolution by answering 0 at the cleanup prompt. To disable the cleanup process, respond with -1.
Once the image has been convolved, press the space bar to switch between the original image (or cleaned image) and the convolved image. The effects of convolution become very apparent by simply holding down the space bar! Should you be fortunate enough to use the program on a color system, you can toggle colors by pressing the letters A-P (A=black, B=blue, etc...). Finally, hitting the ESC key will get you out (as will pressing any key during convolution).
Conclusion
The program works fine now, though you may wish to add the capability for single page display modes, developing a method to save convolved images, further optimization (assembly?), and the use of other image file formats.
Convolution is an amazingly straightforward algorithm for image manipulation and can serve as the starting point for AI programs, contour analysis, or image recognition.
Figure 1
Figure 2
Listing 1
/**********************************************************/ /* CONVOLVE.C - Turbo C 2.0 implementation of image */ /* convolution */ /* ---------- by Wesley G. Faler. All code is "as is". */ /* There is NO copyright. Use this code as you will, and */ /* if you make money at it, good for you. */ /**********************************************************/ #include<stdlib.h> #include<stdio.h> #include<graphics.h> #include<alloc.h> #include<ctype.h> int load_cut(char *fname); int load_convolution_matrix(char *fname); int convolve_image(void); int swap_pictures(void); int minx,maxx,miny,maxy; int LOADPAGE=0; int ENHANCEPAGE=1; int *cmat, *pmat, *vmat; int cmx,cmy,cmnum; struct palettetype palette,newpal; int driver,mode; int cleancut=-1; int init_graphics(void) { driver=DETECT; mode=0; detectgraph(&driver,&mode); if(driver==VGA) mode=VGAMED; initgraph(&driver,&mode,""); getpalette(&palette); getpalette(&newpal); } int cleanup_image(void) { int i,j,num,x,y,k; if(cleancut<0) return; setactivepage(LOADPAGE); setvisualpage(ENHANCEPAGE); for(x=minx;x<maxx;x++) { for(y=miny;y<maxy;y++) { if(getpixel(x,y)!=0) num=-1; else num=0; for(j=-1;j<2;j++) { for(i=-1;i<2;i++) { if(getpixel(x+i,y+j)!=0) num++; } } if(num>cleancut) { k=getpixel(x,y); setactivepage(ENHANCEPAGE); putpixel(x,y,k); setactivepage(LOADPAGE); } } } k=ENHANCEPAGE; ENHANCEPAGE=LOADPAGE; LOADPAGE=k; } void show_test_image(void) { int i; minx-cmx; miny=cmy; maxx=100+minx; maxy=100+miny; setcolor(1); moveto(minx,miny); randomize(); for(i=0;i<20;i++) lineto(random(100)+minx,random(100)+miny); for(i=0;i<10;i++) fillellipse(random(50)+25+minx,random(50)+25+miny, random(25),random(25)); } main( ) { char fname[50]; int flag=0; load_convolution_matrix("matrix.dat"); printf(".CUT file (1) or test image (0)?"); scanf("%d",&flag); flag= flag? 1:0; if(flag) { fflush(stdin); printf("filename to process:"); gets(fname); } printf("Delete pixels with x or fewer neighbors. x="); scanf("%d",&cleancut); if(cleancut>8) cleancut=8; init_graphics(); setactivepage(1); cleardevice(); setactivepage(0); cleardevice(); setactivepage(LOADPAGE); setvisualpage(LOADPAGE); if(flag) load_cut(fname); else show_test_image(); cleanup_image(); setvisualpage(ENHANCEPAGE); convolve_image(); swap_pictures(); restorecrtmode(); } int toggle_colors(char c) { c=tolower(c); c=c-'a'; if(c<0 || c>=palette.size) return 0; newpal.colors[c]= palette.colors[c]-newpal.colors[c]; setpalette(c,newpal.colors[c]); return 1; } int swap_pictures(void) { int mode=0; char a; setvisualpage(LOADPAGE); for (;;) { a=getch(); if(a==27) return; if(toggle colors(a)) continue; if(mode==0) setvisualpage(ENHANCEPAGE); if(mode==1) setvisualpage(LOADPAGE); mode=1-mode; } } int convolve_image(void) { int i,j,k,nval; int *vx, *vy, *c; int colmax,offset,end,midy; char **lines=NULL; char *temp=NULL; offset =-minx+(mcx/2); end=cmy-1; midy=cmy/2; lines=(char **)malloc(cmy*sizeof(char *)); for(i=0;i<cmy;i++) lines[i]-(char *)malloc(sizeof(char)*(maxx-minx+cmx+1)); setactivepage(LOADPAGE); for(j=-cmy/2;j<cmy/2;j++) { for(i=minx-cmx/2;i<(maxx+cmx/2+1);i++) { lines[j+midy][i+offset]-getpixel(i,j+miny); } } colmax=getmaxcolor(); for(j = miny;j<maxy;j++) { setactivepage(LOADPAGE); for(i=j+cmy/2,k=minx-cmx/2,nval=maxx+cmx/2;k<nval;k++) lines[end][k+offset] = getpixel(k,i); for(i=minx;i<maxx;i++) { /* Load & multiply neighbors into matrix */ setactivepage(LOADPAGE); vx=vmat; vy=vmat+1; c=cmat; nval=0; for(k=0;k<cmnum;k++) { if(*c) nval+= lines[(*vy)+midy][i+(*vx)+offset]*(*c); /* if(*c) nval+= getpixel(i+(*vx),j+(*vy)) * (*c); */ c++; vx+=2; vy+=2; } /* Cut off values too high or too low */ if(nval<0) nval=0; if(nval>colmax) nval=colmax; /* Place new pixel value */ setactivepage(ENHANCEPAGE); putpixel(i,j,nval); } if(kbhit()) { getch(); break; } /* rotate line pointers */ temp=lines[0]; for(i=1;i<cmy;i++) lines[i-1]=lines[i]; lines[end]=temp; } for(i=0;i<cmy;i++) { if(lines[i]!=NULL) free(lines[i]); } if(lines!=NULL) { free(lines); } return; } int build_offset_vectors(void) { int *t; int il,im,jl,jm,i,j; il=-cmx/2; im=cmx+il; jl=-cmy/2; jm=cmy+jl; t=vmat; for(j = jl;j<jm;j++) { for(i=il;i<im;i++) { *t++=i; *t++=j; } } } int load_convolution_matrix(char *fname) { /* Layout of matrix file: #x #y x0y0 x1y0 ... xny1 .... .... ... .... x0ym x1ym ... xnym */ FILE *mf; int *t; int i,j,im,jm; if( (mf=fopen(fname,"rt"))==NULL) { printf("Cannot load matrix file.\n"); abort(); } fscanf(mf,"%d%d",&im,&jm); if( (im&1)==0 || (jm&1)==0) { printf("Convolution matrix MUST have a center point.\n"); abort(); } if( (cmat=(int *)calloc(im*jm,sizeof(int)))==NULL ) { printf("Unable to calloc convolution matrix.\n"); abort(); } cmx=im; cmy=jm; cmnum=im*jm; t=cmat; for(j=0;j<jm;j++) { for(i=0;i<im;i++) { if( fscanf(mf,"%d",t++)!=1 ) { printf("Unable to read matrix.\n"); abort(); } } } fclose(mf); build_offset_vectors(); } int load_cut(char *fname) { static unsigned char st[3000]; char *sp=st,*spend; int stp=0; int width,height; FILE *fp; int x,y,xl,yl; int i,n,len,d,j; fp=fopen(fname,"rb"); width=getw(fp); height=getw(fp); xl =cmx; yl=cmy; minx=xl; miny=yl; maxx=xl+width; maxy=yl+height; if(maxy>(getmaxy()-cmy)) { maxy=getmaxy()-cmy; height=maxy-yl; } getw(fp); y=yl-1; for(sp=st,n=0;n<height;n++) { stp=getw(fp); for(sp=st,spend=st+stp;sp<spend;) *sp++=getc(fp); sp=st; spend=sp+stp; x=xl; y++; while(sp<spend) { if(*((unsigned char *)sp)>0x80) { len=(*sp++) & 0x7f; if(!(*sp)) { x+=len; continue; } setcolor(*sp++); moveto(x,y); linerel(len,0); x+=len; continue; } else { len=*sp++; for(j-0;j<len;j++) putpixel(x++,y,*sp++); continue; } } } fclose(fp); }