ml_detect.c

feret人脸图象数据库处理代码

/*----------------------------------------------------------------------
PROGRAM: ml_detect.c
DATE:    7/11/94
AUTHOR:  Baback Moghaddam, [e-mail]
------------------------------------------------------------------------


  Feature Geometry Validation using Maximum Likelihood Minimization

  This routine reads in BF format files: dmf%dt%.bf from the indir
  and a feature model matrix (in BF) format as well as a geometry_musigma.bf
  file

  And dumps the coods of best feature matches in output file

	       
---------------------------------------------------------------------- */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <limits.h>
#include <float.h>
#include "util.h"
#include "io.h"


struct Points {
  int n;
  int *x;
  int *y;
  float *f;
};

struct Point {
  int x;
  int y;
  float f;
};

#include "functions.h"

/* ----------- CONFIGURES ---------------- */


#define MAX_NUM_TEMPLATES         6
#define MAX_NUM_MINIMA          100    /* max num of minima for array storage */
#define MAX_CHARS               512


/* ----------- Command-Line Parsing Stuff ------- */

extern int optind;
extern char *optarg;
char *progname;          /* used to store the name of the program  */
char comline[170];       /* used to store the entire command line  */

#define OPTIONS "i:m:g:l:o:s:t:f"

char *usage = "\t-i indir -l list -m model.bf -g geometry_musigma.bf \n\
\t\t\t[-o outfile] [-s mean_scale] [-t top_N_max] [-f]\n";


char *help="\
Eigenfeature Detection with Maximum Likelihood Matching\n\n\
-i indir    \t input directory where distancemaps (or minima) reside\n\
-l list     \t ASCII list of filenames to process (one per line)\n\
-m model.bf \t BF-format file containing canonical feature locations\n\
-g priors.bf\t geometry_musigma.bf file containing geometry priors\n\
-o outfile  \t output file\n\
-s mean_scale\t float, expected scale of feature constellation (default = 1.0)\n\
-t top_N_max\t integer, # of local minima for each feature (default = 10)\n\
-f          \t read minima files instead of distancemaps\n";



/*----------------------------------------------------------------------*/
/* ---------------------------- MAIN ---------------------------------- */

main(int argc, char *argv[])
{

  register int i,j,k,l,ii,jj;
  int f, c, nframe, nfeatures, sets, bytes_pixel;
  int nrow, ncol, image_nrow, image_ncol, max_N, xc, yc;
  char command[MAX_CHARS],indir[MAX_CHARS],listfile[MAX_CHARS],
       modelfile[MAX_CHARS],outfile[MAX_CHARS],line[MAX_CHARS],
       infile[MAX_CHARS],filename[MAX_CHARS],geometryfile[MAX_CHARS];

  float **image;
  float fval1, fval2, fval;
  float maxerror;
  FILE *fp, *fp1, *fp2;         /* for output values dump */


  int i_min, i_max, j_min, j_max, N_minima;

  struct Points myPoints_array[MAX_NUM_TEMPLATES];

  int top_N[MAX_NUM_TEMPLATES];
  int **Combination_matrix;
  int num_combinations_max;
  int num_combinations;
  int best_c;
  float **Model;
  float **Match;
  float **Stiffness;
  float match_scale, match_energy;

  int Geom_Dim;
  float *Geom_Mean;
  float **Geom_Cov;
  float **Geom_ICov;

  /* required input flags */
  
  int errflag   = 0;
  int inflag    = 0;
  int outflag   = 0;
  int modelflag = 0;
  int listflag  = 0;
  int geomflag  = 0;

  /* command line defaults */

  float mean_scale = 1.0;
  int top_N_max = 10;
  int do_read_minima = 0;      /* default is to read detection maps */


  /* setup program name and command line strings */
  
  progname = argv[0];
  for (i=0; i<argc; i++)
    strcat(comline, argv[i]),strcat(comline, " ");
  sprintf(outfile,"%s.out", progname);            /* default output file */  
  

  /* ----------------------  Command Line Parse ------------------------ */
  
  while ((c = getopt(argc, argv, OPTIONS)) != EOF)
    switch (c) {
      
    case 'i':
      strcpy(indir, optarg);
      inflag = 1;
      break;

    case 'l':
      strcpy(listfile, optarg);
      listflag = 1;
      break;

    case 'm':
      strcpy(modelfile, optarg);
      modelflag = 1;
      break;

    case 'g':
      strcpy(geometryfile, optarg);
      geomflag = 1;
      break;

    case 'o':
      strcpy(outfile, optarg);
      outflag = 1;
      break;

    case 's':
      mean_scale = atof(optarg);
      break;

    case 't':
      top_N_max = atoi(optarg);
      break;

    case 'f':
      do_read_minima = 1;
      break;

    case '?':
      errflag = 1;
      break;
      
    }
  
  
  /* command line error check */
  
  if (errflag || !inflag || !modelflag || !listflag || !geomflag) {
    fprintf(stderr,"\nUSAGE: %s %s\n%s\n", progname, usage, help);
    exit(1);
  }
  


  /* ----  read indir descriptor file -------- */
  
  if (do_read_minima == 0) {
    read_descriptor(indir, &nframe, &sets, &bytes_pixel, 
		    &image_ncol, &image_nrow);
    if (sets>1) 
      myerror("Input files must be single-set DAT files!");
    if (bytes_pixel != 4) {
      fprintf(stderr,"ERROR: detection maps in %s must be float!\n", indir);
      exit(1);
    }
    image = matrix(1, image_nrow, 1, image_ncol);
  }


  /* ---- read model file ----- */

  Model = read_BIN(modelfile, &nrow, &ncol);
  if (nrow<2) {
    fprintf(stderr,"ERROR: Model file %s must have 2 rows!\n\n",
	    modelfile);
    exit(1);
  }
  nfeatures = ncol;


  /* ---- read geometry file ----- */
  { 
    float **M;

    M = read_BIN(geometryfile, &nrow, &ncol);

    Geom_Dim = 0;
    for (i=1; i<= nfeatures; i++)
      for (j=i+1; j<=nfeatures; j++)
	Geom_Dim++;
    Geom_Dim = 2 * Geom_Dim;             /* inter-dist + inter-angles */

    if (ncol != Geom_Dim) {
      fprintf(stderr,"ERROR: Geometry file %s must have %d cols not %d!\n\n",
	      geometryfile, Geom_Dim, ncol);
      exit(2);
    }
    if (nrow != Geom_Dim+1) {
     fprintf(stderr,"ERROR: Geometry file %s must have %d rows not %d!\n\n",
	      geometryfile, Geom_Dim+1, nrow);
      exit(2);
    }
    
    Geom_Mean = vector(1, Geom_Dim);
    Geom_Cov  = matrix(1, Geom_Dim, 1, Geom_Dim);
    Geom_ICov = matrix(1, Geom_Dim, 1, Geom_Dim);
    
    for (i=1; i<=Geom_Dim; i++)
      Geom_Mean[i] = M[1][i];

    for (i=1; i<=Geom_Dim; i++)
      for (j=1; j<=Geom_Dim; j++)
	Geom_Cov[i][j] = M[i+1][j];

    matrix_inverse(Geom_Cov, Geom_Dim, Geom_ICov);

    free_matrix(M, 1, nrow, 1, ncol);
  }



  /* ---- now determine max # of minima and allocate space ------ */

  N_minima = MAX_NUM_MINIMA;
  if (do_read_minima == 0) {
    i_max = image_nrow;
    j_max = image_ncol;
    i_min = j_min = 1;
    N_minima = (i_max - i_min) * (j_max - j_min) / 9;
  }


  for (k=1; k<=nfeatures; k++) {
    myPoints_array[k].n = N_minima;
    myPoints_array[k].x = ivector(1, N_minima);
    myPoints_array[k].y = ivector(1, N_minima);
    myPoints_array[k].f =  vector(1, N_minima);
  }



  /* ----  set up minimization stuff ---- */

  Stiffness = matrix(1, nfeatures, 1, nfeatures);
  for (i=1; i<=nfeatures; i++)
    for (j=1; j<=nfeatures; j++)
      Stiffness[i][j] = 1.0;

  num_combinations_max = 1;
  for (i=1; i<=nfeatures; i++) {
    top_N[i] = top_N_max;
    num_combinations_max *= top_N_max;
  }
  Combination_matrix = imatrix(1, num_combinations_max, 1, nfeatures);
  Match = matrix(1, 3, 1, nfeatures);


  

  /* ----  open output data file in outdir ---- */

  sprintf(filename,"%s",outfile);
  if ((fp = fopen(filename, "w")) == NULL) {
    fprintf(stderr,"ERROR Could not open output file %s \n\n", filename);
    exit(1);
  }
    

  fprintf(stdout,"%s\n\n",comline);     /* echo command line */
  fprintf(fp,"# %s\n\n",comline); 


  /* ------------------------------------------------------ */
  /* -------------- MAIN Sequence Loop -------------------- */
  /* ------------------------------------------------------ */
 
  /* ---- loop over input list file and warp ------- */

  if ((fp2 = fopen(listfile, "r")) == NULL) {
    fprintf(stderr,"ERROR: Could not open list file %s \n\n", listfile);
    exit(1);
  }

  nframe = 0;
  while (fgets(line, MAX_CHARS, fp2)) {

    if (strncmp(line, "#", 1) != 0 && strlen(line)>1) {

      nframe++;


      /* --- read filename  --- */

      sscanf(line, "%s", infile);
      
      
      /* ---- find/read the local minima in each map -------- */
      
      if (do_read_minima) {
	for (k=1; k<=nfeatures; k++) {
	  sprintf(filename,"%s/%s_t%d_minima", indir, infile, k);
	  if ((fp1 = fopen(filename, "r")) == NULL) {
	    fprintf(stderr,
		    "ERROR: Could not open minima file %s \n\n", filename);
	    exit(1);
	  }
	  i=1;
	  while (fscanf(fp1, "%d %d %f", 
			myPoints_array[k].x+i, myPoints_array[k].y+i,
			myPoints_array[k].f+i) != EOF  && i<N_minima)
	    i++;
	  myPoints_array[k].n = i-1;
	  fclose(fp1);
	}
      }
      else {
	for (k=1; k<=nfeatures; k++) {
	  sprintf(filename,"%s/%s_t%d", indir, infile, k);
	  read_RAW_float(filename, image, image_nrow, image_ncol);
	  minima(image, 1, image_nrow, 1, image_ncol, &myPoints_array[k]);
	}
      }
      
      
      
      /* ---- minimzation search ------- */
      
      
      for (k=1; k<=nfeatures; k++)
	if (myPoints_array[k].n<top_N_max)
	  top_N[k] = myPoints_array[k].n;
      i = 1;
      for (k=1; k<=nfeatures; k++)   /* this is superfluous, cuz inside */
	i *= top_N[k];               /* enumerate, the same is computed */
      
      num_combinations = enumerate(top_N,
				   nfeatures,
				   Combination_matrix);
      
      best_c = best_ml_match(myPoints_array,
			     Combination_matrix, num_combinations, nfeatures,
			     Model,
			     Match,
			     Geom_Mean,
			     Geom_ICov,
			     Geom_Dim,
			     &match_scale,
			     &match_energy);
      
    
      
      /* ---------- Dump information to stdout --------- */
      
      fprintf(stdout,"FILE : %s\n", infile);
      fprintf(stdout,"Minima: ");
      for (k=1; k<=nfeatures; k++) 
	fprintf(stdout,"%3d ", myPoints_array[k].n);
      fprintf(stdout,"\n");
      fprintf(stdout,"Best  : ");
      for (k=1; k<=nfeatures; k++)
	fprintf(stdout,"%3d ", top_N[k]);
      fprintf(stdout," -->  %d  combinations\n", num_combinations);
      fprintf(stdout,"Rank  : ");
      for (k=1; k<=nfeatures; k++)
	fprintf(stdout,"%3d ", Combination_matrix[best_c][k]);
      fprintf(stdout,"\nScale : %1.6f \nEnergy: %1.6e \n",
	      match_scale, match_energy);
      for (k=1; k<=nfeatures; k++) {
	fprintf(stdout,"%3d %3d\n", 
		(int) Match[1][k],  (int) Match[2][k]);
      }
      fprintf(stdout,"\n");
      
      
      
      
      /* --------- dump matches to output file ------------- */
      
      fprintf(fp,"%s \t", infile);
      for (j=1; j<=2; j++) 
	for (k=1; k<=nfeatures; k++) 
	  fprintf(fp,"%3d ",  (int) Match[j][k]);
      fprintf(fp,"  %1.6f  %1.6e", match_scale, match_energy);
      fprintf(fp,"\n");
      fflush(fp);

    } /* end of valid list filename */

  } /* end of frame sequence loop */
  
  fclose(fp); 
  fclose(fp2);
  
    
  /* --- free up allocated matrices ----- */


  for (k=1; k<=nfeatures; k++) {
    free_ivector(myPoints_array[k].x, 1, N_minima);
    free_ivector(myPoints_array[k].y, 1, N_minima);
    free_vector( myPoints_array[k].f, 1, N_minima);
  }

  free_matrix(Stiffness, 1, nfeatures, 1, nfeatures); 
  free_imatrix(Combination_matrix, 1, num_combinations_max, 1, nfeatures);
  free_matrix(Match, 1, 3, 1, nfeatures);  

  if (do_read_minima == 0)
    free_matrix(image, 1, image_nrow, 1, image_ncol);

  
  free_vector(Geom_Mean, 1, Geom_Dim);
  free_matrix(Geom_Cov,  1, Geom_Dim, 1, Geom_Dim);
  free_matrix(Geom_ICov, 1, Geom_Dim, 1, Geom_Dim);
  
  return 0;

}



/*------------------------------------------------------------------- */
/* -------------------- end of main() ------------------------------- */
/*------------------------------------------------------------------- */