/********************************************************************/ 
/* File name: MINING 10.0                                           */
/* Synopsis : Performs context-free image retrieval.                */
/*								    */
/* Processes:		    					    */
/*          1) RGB to rgb color conversion.                         */
/*          2) Lowpass filter images.                               */
/*          3) Generate ROI matched filter.                         */
/*          4) Correlate with data set.                             */
/*          5) Find nearest match.                                  */    
/********************************************************************/
#include <stdio.h> 
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <sys/types.h>
#include <sys/wait.h>
#include <signal.h>

#define VERSION "10.0"
#define OUTPUT_FILE_NAME "output-10.0.txt"
#define IMAGES 5000
#define BUFFERSIZE 1024

#define unchar unsigned char

static double COSINE[BUFFERSIZE], SINE[BUFFERSIZE];

static double rDataR[BUFFERSIZE][BUFFERSIZE]; // Red buffer.
static double cDataR[BUFFERSIZE][BUFFERSIZE]; // Red buffer.

static unchar iDataR[BUFFERSIZE][BUFFERSIZE]; // Red buffer.
static unchar iDataG[BUFFERSIZE][BUFFERSIZE]; // Green or BW buffer.
static unchar iDataB[BUFFERSIZE][BUFFERSIZE]; // Blue buffer.

static unchar hDataR[BUFFERSIZE][BUFFERSIZE]; // Red buffer.
static unchar pDataR[BUFFERSIZE][BUFFERSIZE]; // Red buffer.
  
/********************************************************************/
/* Unix PGM access routines                                         */
/********************************************************************/
/****************************************/
/* PGM filename prompt.                 */  
/****************************************/
void pgmFileName(char FileName[80], char *Prompt)
{
   printf("* Enter a %s filename: ", Prompt);
   scanf("%s", FileName);
   printf("\n");
}

/****************************************/
/* Clear comments from PGM file.	*/
/****************************************/
void pgmCommentClear(FILE *disk)
{
int    data;
int    err;
long   l;
unsigned char ch;
fpos_t *pos;
   pos = &l;
   fread(&ch, 1, 1, disk);
   err = fgetpos(disk, pos);
   fread(&ch, 1, 1, disk);
   if (ch == '#') {
      while (ch == '#') {
         while (ch != '\n') fread(&ch, 1, 1, disk);
         err = fgetpos(disk, pos);
         fread(&ch, 1, 1, disk);
      }
      err = fsetpos(disk, pos);      
   } else {
      err = fsetpos(disk, pos);
   }
}

/****************************************/
/* PGM file load.                       */  
/****************************************/
void pgmLoad(char FileName[80], int *xSize, int *ySize, int *color)
{ 
int    x, y;
int    xmax, ymax;
int    pmax;
int    cFlag;
unchar *ptrR, *ptrG, *ptrB;
char   ch;
char   type[3];
FILE   *disk;

   // Open file for reading.
   printf("* Loading PGM file: %s \n", FileName);
   if ((disk = fopen(FileName, "rb")) == NULL) {
      printf(">> ERROR - Cannot open PGM file.\n");
      return;
   }

   // Read in file header.
   fscanf(disk, "%s", &type);
   pgmCommentClear(disk);
   fscanf(disk, "%d", &xmax);
   fscanf(disk, "%d", &ymax);
   fscanf(disk, "%d", &pmax);
   fread(&ch, 1, 1, disk);

   // Check file.
   if (strcmp(type, "P5") == 0) {
      cFlag = 0;
   } else {
      if (strcmp(type, "P6") == 0) {
         cFlag = 1;
      } else {
         printf(">> ERROR - Incorrect PGM format.\n");
         fclose(disk);
         return;
      }
   }  
   if (xmax > BUFFERSIZE) {
      printf(">> ERROR - Image is too big.\n");
      fclose(disk);
      return;
   }  
   if (ymax > BUFFERSIZE) {
      printf(">> ERROR - Image is too big.\n");
      fclose(disk);
      return;
   }  
   if (pmax > 255) {
      printf(">> ERROR - Image is too deep.\n");
      fclose(disk);
      return;
   }  

   // Clear buffers.
   for (y = 0; y < BUFFERSIZE; y++) {
       ptrR = iDataR[y];
       ptrG = iDataG[y];
       ptrB = iDataB[y];
       for (x = 0; x < BUFFERSIZE; x++) {
           *ptrR = 0;
           *ptrG = 0;
           *ptrB = 0;
           ptrR++;
           ptrG++;
           ptrB++;
       }
   }

   if (cFlag == 0) {
      // Load b/w image.
      for (y = 0; y < ymax; y++) {
          ptrG = iDataG[y];
          for (x = 0; x < xmax; x++) {
              *ptrG = getc(disk);
              ptrG++;
          }
      }
   } else {
      // Load color image.
      for (y = 0; y < ymax; y++) {
          ptrR = iDataR[y];
          ptrG = iDataG[y];
          ptrB = iDataB[y];
          for (x = 0; x < xmax; x++) {
              *ptrR = getc(disk);
              *ptrG = getc(disk);
              *ptrB = getc(disk);
              ptrR++;
              ptrG++;
              ptrB++;
          }
      }
   }
   fclose(disk);

   // Return image size and color.;
   *xSize = xmax;
   *ySize = ymax;
   *color = cFlag;
}

/****************************************/
/* PGM file save.                       */  
/****************************************/
void pgmSave(char FileName[80], int xSize, int ySize, int cFlag)
{
int    x, y;
unchar *ptrR, *ptrG, *ptrB; 
FILE   *disk;

   // Open file for writing.
   printf("* Saving PGM file: %s \n", FileName);
   disk = fopen(FileName, "wb");

   // Write file header.
   if (cFlag == 0) fprintf(disk, "P5\n");
              else fprintf(disk, "P6\n");
   fprintf(disk, "%d %d\n", xSize, ySize);
   fprintf(disk, "255\n");

   if (cFlag == 0) {
      // Write BW data.
      for (y = 0; y < ySize; y++) {
          ptrG = iDataG[y];
          for (x = 0; x < xSize; x++) {
              putc(*ptrG, disk);
              ptrG++;
          }
      }  
   } else {
      // Write color data.
      for (y = 0; y < ySize; y++) {
          ptrR = iDataR[y];
          ptrG = iDataG[y];
          ptrB = iDataB[y];
          for (x = 0; x < xSize; x++) {
              putc(*ptrR, disk);
              putc(*ptrG, disk);
              putc(*ptrB, disk);
              ptrR++;
              ptrG++;
              ptrB++;
          }
      }
   }
   fclose(disk);
}

/********************************************************************/
/* Display generation routines                              	    */
/********************************************************************/
/****************************************/
/* Display image.        	        */  
/****************************************/
void pgmDisplay(char FileName[24], int *Pid)
{
int  exitCode;
char UnixCmd[255];
char CmdArg1[255];
char CmdArg2[255];

   strcpy(UnixCmd, "xv");
   strcpy(CmdArg1, "-name");
   strcpy(CmdArg2, FileName);
   if (*Pid != 0) kill(*Pid, 9);
   *Pid = fork();
   if (*Pid == 0) {
      execlp(UnixCmd, CmdArg1, CmdArg2, NULL);
   }
}

/********************************************************************/
/* Filter routines.                                                 */
/********************************************************************/
/****************************************/
/* Spatial filter image.       	     	*/
/****************************************/
void pgmFilterConvolve(char buf, int xSize, int ySize, int filter)
{
int x, y;
int i, j;
int sum1, sum2;
unsigned char Buf1[xSize], Buf2[xSize], Buf3[xSize], BufN[xSize];
int m;
int kern1[3][3];
int kern2[3][3];
int lopass[3][3] = { 1,  1,  1,  1,  1,  1,  1,  1,  1};
int hipass[3][3] = {-1, -1, -1, -1,  9, -1, -1, -1, -1};
int laplac[3][3] = {-1, -1, -1, -1,  8, -1, -1, -1, -1};
int sobelV[3][3] = {-1, -2, -1,  0,  0,  0,  1,  2,  1};
int sobelH[3][3] = {-1,  0,  1, -2,  0,  2, -1,  0,  1};

    // Set scale term.
    m = 1;
    if (filter == 0) m = 9;

    // Initialize kernels.
    for (i = 0; i < 3; i++) {
        for (j = 0; j < 3; j++) {
            kern1[i][j] = 0;
            kern2[i][j] = 0;
         }
    }

    // Set filter kernel.
    for (i = 0; i < 3; i++) {
        for (j = 0; j < 3; j++) {
            if (filter == 0) kern1[i][j] = lopass[i][j];
            if (filter == 1) kern1[i][j] = hipass[i][j];
            if (filter == 2) kern1[i][j] = laplac[i][j];
            if (filter == 3) kern1[i][j] = sobelH[i][j];
            if (filter == 3) kern2[i][j] = sobelV[i][j];
        }
    }

    // Convolve image with kernel.
    for (x = 0; x < xSize; x++) {
        if (buf == 'R') {
           Buf1[x] = iDataR[0][x];
           Buf2[x] = iDataR[1][x];
        }     
        if (buf == 'G') {
           Buf1[x] = iDataG[0][x];
           Buf2[x] = iDataG[1][x];
        }     
        if (buf == 'B') {
           Buf1[x] = iDataB[0][x];
           Buf2[x] = iDataB[1][x];
        }   
    }
    for (y = 2; y < ySize-1; y++) {
        for (x = 0; x < xSize; x++) {
            if (buf == 'R') Buf3[x] = iDataR[y][x];
            if (buf == 'G') Buf3[x] = iDataG[y][x];
            if (buf == 'B') Buf3[x] = iDataB[y][x];
        }
        for (x = 1; x < xSize-1; x++) {
            sum1 = 0;
            sum2 = 0;
            for (i = 0; i < 3; i++) {
                sum1 += (Buf1[x+i-1] * kern1[i][0]);
                sum1 += (Buf2[x+i-1] * kern1[i][1]);
                sum1 += (Buf3[x+i-1] * kern1[i][2]); 
                sum2 += (Buf1[x+i-1] * kern2[i][0]);
                sum2 += (Buf2[x+i-1] * kern2[i][1]);
                sum2 += (Buf3[x+i-1] * kern2[i][2]);
            }
	    //sum1 = abs(sum1 / m) | abs(sum2 / m);
            sum1 = abs(sum1 / m);
            sum2 = abs(sum2 / m);
            sum1 = sqrt(sum1*sum1 + sum2*sum2);
            if (sum1 > 255) sum1 = 255;
            BufN[x] = sum1;
	} 
        for (x = 1; x < xSize-1; x++) {
            if (buf == 'R') iDataR[y-1][x] = BufN[x];
            if (buf == 'G') iDataG[y-1][x] = BufN[x];
            if (buf == 'B') iDataB[y-1][x] = BufN[x];
        }
        for (x = 0; x < xSize; x++) {
	    Buf1[x] = Buf2[x];
	    Buf2[x] = Buf3[x];
	}
    }
}

/********************************************************************/
/* Color routines.                                                  */
/********************************************************************/
/****************************************/
/* RGB to rgb.        	        	*/  
/****************************************/
void RGB2rgb(int xSize, int ySize)
{
int x, y;
int I;
    for (y = 0; y < ySize; y++) {
        for (x = 0; x < xSize; x++) {
            I = (iDataR[y][x] + iDataG[y][x] + iDataB[y][x]);
            if (I == 0) I = 1; 
            iDataR[y][x] = (255 * iDataR[y][x]) / I;
            iDataG[y][x] = (255 * iDataG[y][x]) / I;
            iDataB[y][x] = (255 * iDataB[y][x]) / I;
        }
    }
}

/****************************************************************/
/* FOURIER TRANSFORM subroutines.			        */
/****************************************************************/
/********************************************************/
/* Generates Fourier transform twiddle table.		*/
/********************************************************/
void MAKE_TWIDDLE_TABLE(void)
{
int    i;
int    Size2 = (BUFFERSIZE / 2);
double arg1, arg2;
   arg1  = (6.283185307 / (double)BUFFERSIZE);
   for (i = 0; i < BUFFERSIZE; i++) {
       arg2      = (double)(i - Size2) * arg1;
       COSINE[i] = cos(arg2);
       SINE[i]   = sin(arg2);
   }
}

/********************************************************/
/* Fast Fourier transform. 				*/
/********************************************************/
void FFT(double RData[BUFFERSIZE], double CData[BUFFERSIZE], double Flag)
{
int    i, z, zz, j, jj, j2, k, r, d;
int    pwr, nd2, nml, arg, twf;
double tempr, tempc, c, s;
   for (i = 0; i < BUFFERSIZE; i+=2) {
       RData[i] *= -1.0 * Flag;
       CData[i] *= -1.0 * Flag;
   }
   j  = 1;
   for (z = 1; z < BUFFERSIZE; z++) {
       if (z < j) {
          jj        = j - 1;
          zz        = z - 1;
          tempr     = RData[jj];
          RData[jj] = RData[zz];
          RData[zz] = tempr; 
          tempc     = CData[jj];
          CData[jj] = CData[zz];
          CData[zz] = tempc; 
       }
       k = (BUFFERSIZE / 2);
       while (k < j) {
          j -= k;
          k  = (k / 2);
       }
       j += k;
   }
   if (BUFFERSIZE ==   64) pwr = 6;
   if (BUFFERSIZE ==  128) pwr = 7;
   if (BUFFERSIZE ==  256) pwr = 8;
   if (BUFFERSIZE ==  512) pwr = 9;
   if (BUFFERSIZE == 1024) pwr = 10;
   for (z = 1; z <= pwr; z++) {
       r   = pow(2, z);
       d   = (r / 2);
       arg = (BUFFERSIZE / r);
       for (zz = 0; zz < d; zz++) {
           twf = (BUFFERSIZE / 2) + (zz * arg);
           c   = COSINE[twf];
           s   = (SINE[twf] * Flag);
           k   = zz;
           while (k < BUFFERSIZE) {
              j2        = (k + d);
              tempr     = (RData[j2] * c) + (CData[j2] * s);
              tempc     = (CData[j2] * c) - (RData[j2] * s);
              RData[j2] = (RData[k] - tempr);
              CData[j2] = (CData[k] - tempc);
              RData[k]  = (RData[k] + tempr);
              CData[k]  = (CData[k] + tempc);
              k        += r;
           }
       }
   }
}

/********************************************************/
/* 2D Fast Fourier transform. 	        		*/
/********************************************************/
void TWO_D_FFT(double Flag)
{
int    x, y;
double *ptrR;
double *ptrC;
double RData[BUFFERSIZE], CData[BUFFERSIZE];
    /* X axis FFT. */
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            RData[x] = *ptrR;
            CData[x] = *ptrC;
            ptrR++;
            ptrC++;
        } 
        FFT(RData, CData, Flag);
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrR = RData[x] / (double)BUFFERSIZE;
            *ptrC = CData[x] / (double)BUFFERSIZE;
            ptrR++;
            ptrC++;
        }
    }

    /* Y axis FFT. */
    for (x = 0; x < BUFFERSIZE; x++) {
        for (y = 0; y < BUFFERSIZE; y++) {
            RData[y] = rDataR[y][x];
            CData[y] = cDataR[y][x];
        }            
        FFT(RData, CData, Flag);
        for (y = 0; y < BUFFERSIZE; y++) {
            rDataR[y][x] = RData[y];
            cDataR[y][x] = CData[y];
        }
    }
}

/********************************************************/
/* Copies TARGET to FILTER. 				*/
/********************************************************/
void COPY_TO_FILTER(void)
{
int    x, y;
unchar *ptriR, *ptriG, *ptriB;
unchar *ptrhR;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptriR = iDataR[y];
        ptriG = iDataG[y];
        ptriB = iDataB[y];
        ptrhR = hDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrhR = (*ptriR + *ptriG + *ptriB) / 3;
            ptriR++;
            ptriG++;
            ptriB++;
            ptrhR++;
        } 
    }           
}

/********************************************************/
/* Copies FILTER TO TARGET. 				*/
/********************************************************/
void COPY_TO_TARGET(void)
{
int    x, y;
unchar *ptriR, *ptriG, *ptriB;
unchar *ptrhR;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptriR = iDataR[y];
        ptriG = iDataG[y];
        ptriB = iDataB[y];
        ptrhR = hDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptriR = *ptrhR;
            *ptriG = *ptrhR;
            *ptriB = *ptrhR;
            ptriR++;
            ptriG++;
            ptriB++;
            ptrhR++;
        } 
    }           
}

/********************************************************/
/* Quarters object in image window (for FFT). 		*/
/********************************************************/
void SECTION_FILTER(void)
{
int    x, y, c;
unchar junkR;
    c = (BUFFERSIZE / 2) - 1;
    for (y = 0; y < BUFFERSIZE; y++) {
        for (x = 0; x <= c; x++) {
            junkR            = hDataR[y][x];
            hDataR[y][x]     = hDataR[y][x+c+1];
            hDataR[y][x+c+1] = junkR;     
        } 
    }           
    for (x = 0; x < BUFFERSIZE; x++) {
        for (y = 0; y <= c; y++) {
            junkR            = hDataR[y][x];
            hDataR[y][x]     = hDataR[y+c+1][x];
            hDataR[y+c+1][x] = junkR;     
        } 
    }         
}

/********************************************************/
/* Copies integer to complex arrays. 			*/
/********************************************************/
void COMPLEX_FILTER(void)
{
int    x, y;
unchar *ptrF;
double *ptrR, *ptrC;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrF = hDataR[y];
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrR = (double)(*ptrF);
            *ptrC = 0.0L;
            ptrF++;
            ptrR++;
            ptrC++;
        } 
    }              
}

/********************************************************/
/* Generates 2D BPOHF from TARGET image.	 	*/
/********************************************************/
void FILTER_MAKE_BPOHF(void)
{
int    x, y;
unchar *ptrH;
double *ptrR, *ptrC;
    COPY_TO_FILTER();
    SECTION_FILTER();
    COMPLEX_FILTER(); 
    TWO_D_FFT(1.0); 
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrH = hDataR[y];
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            if ((*ptrR - *ptrC) > 0.0) *ptrH = 255;
                                  else *ptrH = 0;
            ptrH++;
            ptrR++;
            ptrC++;
        }
    }
}

/********************************************************/
/* Copies TARGET to CORRELATE. 				*/
/********************************************************/
void COPY_TO_CORRELATE(void)
{
int    x, y;
unchar *ptriR, *ptriG, *ptriB;
unchar *ptrpR;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptriR = iDataR[y];
        ptriG = iDataG[y];
        ptriB = iDataB[y];
        ptrpR = pDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrpR = (*ptriR + *ptriG + *ptriB) / 3;
            ptriR++;
            ptriG++;
            ptriB++;
            ptrpR++;
        } 
    }           
}

/********************************************************/
/* Copies FILTER TO TARGET. 				*/
/********************************************************/
void CORRELATION_TO_TARGET(void)
{
int    x, y;
unchar *ptriR, *ptriG, *ptriB;
unchar *ptrpR;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptriR = iDataR[y];
        ptriG = iDataG[y];
        ptriB = iDataB[y];
        ptrpR = pDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptriR = *ptrpR;
            *ptriG = *ptrpR;
            *ptriB = *ptrpR;
            ptriR++;
            ptriG++;
            ptriB++;
            ptrpR++;
        } 
    }           
}

/********************************************************/
/* Quarters object in image window (for FFT). 		*/
/********************************************************/
void SECTION_CORRELATE(void)
{
int    x, y, c;
unchar junkR;
    c = (BUFFERSIZE / 2) - 1;
    for (y = 0; y < BUFFERSIZE; y++) {
        for (x = 0; x <= c; x++) {
            junkR            = pDataR[y][x];
            pDataR[y][x]     = pDataR[y][x+c+1];
            pDataR[y][x+c+1] = junkR;     
        } 
    }           
    for (x = 0; x < BUFFERSIZE; x++) {
        for (y = 0; y <= c; y++) {
            junkR            = pDataR[y][x];
            pDataR[y][x]     = pDataR[y+c+1][x];
            pDataR[y+c+1][x] = junkR;     
        } 
    }         
}

/********************************************************/
/* Copies integer to complex arrays. 			*/
/********************************************************/
void COMPLEX_CORRELATE(void)
{
int    x, y;
unchar *ptrP;
double *ptrR, *ptrC;
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrP = pDataR[y];
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrR = (double)(*ptrP);
            *ptrC = 0.0L;
            ptrP++;
            ptrR++;
            ptrC++;
        } 
    }              
}

/********************************************************/
/* Generates 2D FFT correlation.	                */
/********************************************************/
double CORRELATE_MAKE_FFT(void)
{
int    x, y;
int    xLoc, yLoc;
unchar *ptrH;
unchar *ptrP;
double *ptrR, *ptrC;
double pixel;
double temp;
double peak = 0.0;

    COPY_TO_CORRELATE();
    SECTION_CORRELATE();
    COMPLEX_CORRELATE(); 
    TWO_D_FFT(1.0);
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrH = hDataR[y];
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            *ptrR *= (double)(*ptrH - 128) / 127.0;
            *ptrC *= (double)(*ptrH - 128) / 127.0;
            ptrH++;
            ptrR++;
            ptrC++;
        }
    }
    TWO_D_FFT(1.0);
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            pixel = sqrt((*ptrR * *ptrR)+(*ptrC * *ptrC));
            if (peak < pixel) {
               peak = pixel; 
               xLoc = x;
               yLoc = y;
            }
            ptrR++;
            ptrC++;
        }
    }
    printf("Peak: %lf (X: %d Y: %d)\n", peak, xLoc, yLoc);

    temp = (255.0 / (double)peak);
    for (y = 0; y < BUFFERSIZE; y++) {
        ptrP = pDataR[y];
        ptrR = rDataR[y];
        ptrC = cDataR[y];
        for (x = 0; x < BUFFERSIZE; x++) {
            pixel = sqrt((*ptrR * *ptrR)+(*ptrC * *ptrC));
            *ptrP = (unchar)(pixel * temp);
            ptrP++;
            ptrR++;
            ptrC++;
        }
    }
    return peak;
}

/********************************************************************/
/* Search routines                                	    	    */
/********************************************************************/
/****************************************/
/* Bubble sort.       	     		*/
/****************************************/
void bubbleSort(long lData[IMAGES], char cData[IMAGES][80], int count)
{
int  i, j;
long temp;
char text[80];
    for (i = 0; i < count-1; i++) {
        for (j = 0; j < count-1-i; j++) {
            if (lData[j] > lData[j+1]) { 
               temp       = lData[j];
               lData[j]   = lData[j+1];
               lData[j+1] = temp;      
               strcpy(text,       cData[j]);
               strcpy(cData[j],   cData[j+1]);
               strcpy(cData[j+1], text);
            }
        }
    }
}

/****************************************/
/* Search for like image(s).            */  
/****************************************/
pgmFind(char RefName[80], int xSize, int ySize, int cFlag, char MatchName[80])
{
int    i;
int    cntr = 0;
int    TxSize, TySize, TcFlag;
char   FileName[80];
char   BestName[IMAGES][80];
long   Diffs[IMAGES];
FILE   *data;
int    Pid[3];

    RGB2rgb(xSize, ySize);
    pgmFilterConvolve('R', xSize, ySize, 0);
    pgmFilterConvolve('G', xSize, ySize, 0);
    pgmFilterConvolve('B', xSize, ySize, 0);
    FILTER_MAKE_BPOHF();

    data = fopen("images.dat", "r");
    while (!feof(data)) {
        fscanf(data, "%s", FileName); 
        pgmLoad(FileName, &TxSize, &TySize, &TcFlag);
        RGB2rgb(TxSize, TySize);
        if (cFlag == TcFlag) {
           pgmFilterConvolve('R', TxSize, TySize, 0);
           pgmFilterConvolve('G', TxSize, TySize, 0);
           pgmFilterConvolve('B', TxSize, TySize, 0);
           Diffs[cntr] = (long)CORRELATE_MAKE_FFT();
           strcpy(BestName[cntr], FileName);
           cntr++;
        }
    }
    fclose(data);
    cntr--;

    bubbleSort(Diffs, BestName, cntr);
    strcpy(MatchName, BestName[1]);

    // Open output file for writing.
    printf("* Saving output file: %s \n", OUTPUT_FILE_NAME);
    data = fopen(OUTPUT_FILE_NAME, "ab");
    cntr--;

    // Write results.
    fprintf(data, "\n\n");
    fprintf(data, "+--------------------------------------------------+\n");
    fprintf(data, "  Image Mining (v %s)\n", VERSION);
    fprintf(data, "+--------------------------------------------------+\n");
    fprintf(data, "  Target : (%s)\n", RefName);
    fprintf(data, "  Count  : (%d)\n", cntr);
    fprintf(data, "  Best   : (%ld)\n", Diffs[1]);
    fprintf(data, "  Worst  : (%ld)\n", Diffs[cntr]);
    fprintf(data, "  Match  : (%s)\n", MatchName);
    fprintf(data, "\n");
    for (i = 1; i <= 5; i++) {
        fprintf(data, "  %s : (%ld)\n", BestName[i], Diffs[i]);
    }
    fprintf(data, "+--------------------------------------------------+\n");
    fclose(data);

    // Display results.
    printf("Target : (%s)\n", RefName);
    printf("Count  : (%d)\n", cntr);
    printf("Best   : (%ld)\n", Diffs[1]);
    printf("Worst  : (%ld)\n", Diffs[cntr]);
    printf("Match  : (%s)\n", MatchName);
}

/********************************************************************/
/* Main routine                                                     */
/********************************************************************/
void main(int argc, char *argv[])
{
int    xSize, ySize;
int    cFlag;
int    Pid1, Pid2;
char   choice[2];
char   band[2];
char   FileName[80];
char   BestName[80];
FILE   *data;

    MAKE_TWIDDLE_TABLE();

   /******************************/
   /* Command mode.              */
   /******************************/
   if (argc == 2) {
         data = fopen("testset.dat", "r");
         while (1) {
             fscanf(data, "%s", FileName);
             if (feof(data)) break;
             pgmLoad(FileName, &xSize, &ySize, &cFlag);
             pgmFind(FileName, xSize, ySize, cFlag, BestName);
         }
         fclose(data);
         exit(0);
   }

   /******************************/
   /* Menu mode.                 */
   /******************************/
   Pid1 = 0;
   do {
      printf("  Image Mining (v %s)              \n", VERSION);
      printf("+----------------------------------+\n");
      printf("| Press:                           |\n");
      printf("| L) Load file.                    |\n"); 		
      printf("| S) Save file.                    |\n");
      printf("| D) Display file.                 |\n");
      printf("| H) Hartley filter display.       |\n");
      printf("| C) Correlation display (test).   |\n");
      printf("| F) Find like file(s).            |\n");
      printf("|                                  |\n");
      printf("| 1) auto-find 1                   |\n");
      printf("|                                  |\n");
      printf("| X) eXit.                         |\n");
      printf("+----------------------------------+\n");
      printf("Enter option: ");
      scanf("%s", &choice);
      printf("\n\n");
      if (toupper(choice[0]) == 'L') {
         pgmFileName(FileName, "Load");
         pgmLoad(FileName, &xSize, &ySize, &cFlag);
         pgmDisplay(FileName, &Pid1);
      }
      if (toupper(choice[0]) == 'S') {
         pgmFileName(FileName, "Save");
         pgmSave(FileName, xSize, ySize, cFlag);
      }
      if (toupper(choice[0]) == 'D') {
         pgmSave("temp", xSize, ySize, cFlag);
         pgmDisplay("temp", &Pid1);
      }
      if (toupper(choice[0]) == 'H') {
         FILTER_MAKE_BPOHF();
         COPY_TO_TARGET();
         pgmSave("Hartley", BUFFERSIZE, BUFFERSIZE, 0);
         pgmDisplay("Hartley", &Pid2);
      }
      if (toupper(choice[0]) == 'C') {
         FILTER_MAKE_BPOHF();
         pgmFileName(FileName, "Load");
         pgmLoad(FileName, &xSize, &ySize, &cFlag);
         pgmDisplay(FileName, &Pid1);
         CORRELATE_MAKE_FFT();
         CORRELATION_TO_TARGET();
         pgmSave("Correlation", BUFFERSIZE, BUFFERSIZE, 0);
         pgmDisplay("Correlation", &Pid2);
      }
      if (toupper(choice[0]) == 'F') {
         pgmFind(FileName, xSize, ySize, cFlag, BestName);
         pgmDisplay(BestName, &Pid2);
      }
      if (choice[0] == '1') {
         data = fopen("testset.dat", "r");
         while (!feof(data)) {
             fscanf(data, "%s", FileName);
             pgmLoad(FileName, &xSize, &ySize, &cFlag);
             pgmFind(FileName, xSize, ySize, cFlag, BestName);
         }
         fclose(data);
      }
   } while (choice[0] != 'X');
   if (Pid1 != 0) kill(Pid1, 9);
   if (Pid2 != 0) kill(Pid2, 9);
}