stdvqe.c

Vector Quantization压缩算法

  while (fread(tempvector, sizeof(DATA), dimension, inputfile) == dimension
	 && !feof(inputfile) && !ferror(inputfile) ) {
    bestdistortion = HUGE; /* keep convention that ties go to lower index */
    bestcodeword = 0;
    /* find the best codeword */
    for (i = 0; i < codebooksize; i++) {

      tempdistortion = 0.0;
      for (j = 0; j < dimension; j++) {
	temp = ( (double) tempvector[j]) - codebook[i][j];
	tempdistortion += temp*temp;
	if (tempdistortion > bestdistortion) j = dimension;
      }

      if (tempdistortion < bestdistortion) {
	bestdistortion = tempdistortion;
	bestcodeword = i;
      }

      /* if bestdistortion is 0 then the best codeword has been found */
      if (bestdistortion == 0.0) i = codebooksize;
    } /* the best codeword has been found */

    count[bestcodeword]++;
    celldistortion[bestcodeword] += bestdistortion;
    distortion += bestdistortion;
    /* round the codeword */
    for (i = 0; i < dimension; i++) { /* for double data remove floor */
      tempvector[i] = (DATA) (floor(codebook[bestcodeword][i] + 0.5));
      /* tempvector[i] = (DATA) codebook[bestcodeword][i]; */
    }
    /* write the data to the output file */
    if (fwrite(tempvector,sizeof(DATA),dimension,outputfile) != dimension) {
      fprintf(stderr,"%s: %s: %s\n",programname,outputname,NOWRITE);
      return(-1.0);
    }
  } /* all vectors have been encoded */

  return(distortion);
}

double annulus(codebook,count,celldistortion)
     double **codebook;
     long   *count;
     double *celldistortion;
{
  DATA    *tempvector;    /* vector to be encoded */
  double  *norm;          /* array of norms for each codeword */
  double  tempdouble;     /* temporary variable */
  double  tempnorm;       /* temporary variable used to sort codebook */
  double  distortion;     /* total image distortion */
  double  bestdistortion; /* distortion between vector and best codeword */
  double  temp;           /* temporary variable */
  long    *index;         /* array to store sorted codeword order */
  long    bestcodeword;   /* index of closest codeword */
  long    i,j;            /* counters and indices */
  long    imin,imax;      /* indices */
  long    bestnorm;       /* used as an index for sorting norm */

  /* allocate memory for tempvector, norm, and index */
  if(!(tempvector = (DATA *) calloc(dimension,sizeof(DATA))) ||
     !(norm = (double *) calloc(codebooksize, sizeof(double))) ||
     !(index = (long *) calloc(codebooksize, sizeof(long)))) {
    fprintf(stderr,"%s: %s\n",programname,NOMEMORY);
    return(-1.0);
  }

  /* initialize index order */
  for (i = 0; i < codebooksize; i++) {
    index[i] = i;
  }

  /* find the norm of each codeword and store in an array */
  for (i = 0; i < codebooksize; i++) {
    norm[i] = 0.0;
    for (j = 0; j < dimension; j++) {
      norm[i] += codebook[i][j]*codebook[i][j];
    }
    norm[i] = sqrt(norm[i]);
  }

  /* reorder the codebook by ascending norm */
  for (i = 0; i < codebooksize; i++) {
    bestnorm = i;
    /* find the lowest norm */
    for ( j = i; j < codebooksize; j++) {
      if (norm[j] < norm[bestnorm]) bestnorm = j;
    }
    /* record index change */
    j = index[i];
/*    index[i] = bestnorm; this in in error and is changed 1/3/96 */
    index[i] = index[bestnorm];
    index[bestnorm] = j;

    /* reorder the norm arrray */
    tempdouble = norm[i];
    norm[i] = norm[bestnorm];
    norm[bestnorm] = tempdouble;
  }

  /* encode the data and write them to the output file, also find statistics */
  rewind(inputfile);
  distortion = 0.0;
  while (fread(tempvector, sizeof(DATA), dimension, inputfile) == dimension
	 && !feof(inputfile) && !ferror(inputfile) ) {

    /* compute the norm of the input vector */
    tempnorm = 0.0;
    for (j = 0; j < dimension; j++) {
      tempnorm += (double) tempvector[j]*tempvector[j];
    }
    tempnorm = sqrt(tempnorm);

    /* find the codeword with a norm closest to the input vector */
    bestcodeword = min_sup((long) 0,codebooksize-1,tempnorm,norm);

    /* compute the distance between the codeword and the input vector that
       has the closest norm */
    tempdouble = 0.0;
    for (j = 0; j < dimension; j++) {
      temp = ( (double) tempvector[j]) - codebook[index[bestcodeword]][j];
      tempdouble += temp*temp;
    }
    tempdouble = sqrt(tempdouble);

    /* identify the subset of the codewords to search such that
       norm(input)-distance <= norm(codeword) <= norm(input)+distance */
    imin = min_sup((long) 0, codebooksize-1, tempnorm-tempdouble, norm);
    imax = max_inf(imin, codebooksize-1, tempnorm+tempdouble, norm);

    /* find the bestcodeword using partial distortion method */
    bestdistortion = HUGE; /* keep convention that ties go to lower index */
    /* find the best codeword */
    for (i = imin; i <= imax; i++) {
      tempdouble = 0.0;
      for (j = 0; j < dimension; j++) {
	temp = ( (double) tempvector[j]) - codebook[index[i]][j];
	tempdouble += temp*temp;
	if (tempdouble > bestdistortion) j = dimension; /* abort loop */
      }

      if (tempdouble < bestdistortion) {
	bestdistortion = tempdouble;
	bestcodeword = i;
      }

      /* if bestdistortion is 0 then the best codeword has been found */
      if (bestdistortion == 0.0) i = codebooksize;
    } /* the best codeword has been found */

    count[index[bestcodeword]]++;
    celldistortion[index[bestcodeword]] += bestdistortion;
    distortion += bestdistortion;
    /* round the codeword */
    for (i = 0; i < dimension; i++) { /* for double data remove floor */
/* tempvector[i] = (DATA) (floor(codebook[bestcodeword][i] + 0.5));
   this in in error and is changed 9/4/97 */
      tempvector[i] = (DATA) (floor(codebook[index[bestcodeword]][i] + 0.5));
      /* tempvector[i] = (DATA) codebook[index[bestcodeword]][i]; */
    }
    /* write the data to the output file */
    if (fwrite(tempvector,sizeof(DATA),dimension,outputfile) != dimension) {
      fprintf(stderr,"%s: %s: %s\n",programname,outputname,NOWRITE);
      return(-1.0);
    }
  } /* all vectors have been encoded */

  return(distortion);
}


double sphere_annulus(codebook,count,celldistortion)
     double **codebook;
     long   *count;
     double *celldistortion;
{
  DATA    *tempvector;    /* vector to be encoded */
  double  **dist_matrix;  /* array listing distances between all codewords */
  double  *norm;          /* array of norms for each codeword */
  double  tempdouble;     /* temporary variable */
  double  tempnorm;       /* temporary variable used to sort codebook */
  double  distortion;     /* total image distortion */
  double  bestdistortion; /* distortion between vector and best codeword */
  double  temp;           /* temporary variable */
  long    *index;         /* array to store sorted codeword order */
  long    bestcodeword;   /* index of closest codeword */
  long    i,j,k;          /* counters and indices */
  long    imin,imax;      /* indices */
  long    bestnorm;       /* used as an index for sorting norm */

  /* allocate memory for tempvector, norm, index and dist_matrix */
  if(!(tempvector = (DATA *) calloc(dimension,sizeof(DATA))) ||
     !(norm = (double *) calloc(codebooksize, sizeof(double))) ||
     !(index = (long *) calloc(codebooksize, sizeof(long))) ||
    !(dist_matrix  = (double **) calloc(codebooksize,sizeof(double *)))) {
    fprintf(stderr,"%s: %s\n",programname,NOMEMORY);
    return(-1.0);
  }
  for (i = 0; i < codebooksize; i++) {
    if (!(dist_matrix[i] = (double *) calloc(codebooksize,sizeof(double)))) {
      fprintf(stderr,"%s: %s\n",programname,NOMEMORY);
      return(-1);
    }
  }

  /* initialize index order */
  for (i = 0; i < codebooksize; i++) {
    index[i] = i;
  }

  /* find the norm of each codeword and store in an array */
  for (i = 0; i < codebooksize; i++) {
    norm[i] = 0.0;
    for (j = 0; j < dimension; j++) {
      norm[i] += codebook[i][j]*codebook[i][j];
    }
    norm[i] = sqrt(norm[i]);
  }

  /* reorder the codebook by ascending norm */
  for (i = 0; i < codebooksize; i++) {
    bestnorm = i;
    /* find the lowest norm */
    for ( j = i; j < codebooksize; j++) {
      if (norm[j] < norm[bestnorm]) bestnorm = j;
    }
    /* record index change */
    j = index[i];
/*    index[i] = bestnorm; this in in error and is changed 1/3/96 */
    index[i] = index[bestnorm];
    index[bestnorm] = j;

    /* reorder the norm arrray */
    tempdouble = norm[i];
    norm[i] = norm[bestnorm];
    norm[bestnorm] = tempdouble;
  }

  /* find the table of distortions between all of the codewords */
  for (i = 0; i < codebooksize; i++) {
    for (j = i; j < codebooksize; j++) {
      tempdouble = 0.0;
      for (k = 0; k < dimension; k++) {
	temp = ( (double) codebook[index[i]][k]) - codebook[index[j]][k];
	tempdouble += temp*temp;
      }
      tempdouble = sqrt(tempdouble);
      dist_matrix[i][j] = tempdouble;
      dist_matrix[j][i] = tempdouble;
    }
  }

  /* encode the data and write them to the output file, also find statistics */
  distortion = 0.0;
  rewind(inputfile);
  while (fread(tempvector, sizeof(DATA), dimension, inputfile) == dimension
	 && !feof(inputfile) && !ferror(inputfile) ) {
    /* compute the norm of the input vector */
    tempnorm = 0.0;
    for (j = 0; j < dimension; j++) {
      tempnorm += (double) tempvector[j]*tempvector[j];
    }
    tempnorm = sqrt(tempnorm);

    /* find the codework with a norm closest to the input vector */
    bestcodeword = min_sup((long) 0,codebooksize-1,tempnorm,norm);

    /* compute the distance between the codeword and the input vector that
       has the closest norm */
    tempdouble = 0.0;
    for (j = 0; j < dimension; j++) {
      temp = ( (double) tempvector[j]) - codebook[index[bestcodeword]][j];
      tempdouble += temp*temp;
    }
    tempdouble = sqrt(tempdouble);

    /* identify the subset of the codewords to search such that
       norm(input)-distance <= norm(codeword) <= norm(input)+distance */
    imin = min_sup((long) 0, codebooksize-1, tempnorm-tempdouble, norm);
    imax = max_inf(imin, codebooksize-1, tempnorm+tempdouble, norm);

    /* find the bestcodeword using partial distortion method */
    bestdistortion = HUGE; /* keep convention that ties go to lower index */
    /* find the best codeword */
    for (i = imin; i <= imax; i++) {
      if (dist_matrix[i][bestcodeword] <= 2*sqrt(bestdistortion)) {
	tempdouble = 0.0;
	for (j = 0; j < dimension; j++) {
	  temp = ( (double) tempvector[j]) - codebook[index[i]][j];
	  tempdouble += temp*temp;
	  if (tempdouble > bestdistortion) j = dimension; /* abort loop */
	}

	if (tempdouble < bestdistortion) {
	  bestdistortion = tempdouble;
	  bestcodeword = i;
	}

	/* if bestdistortion is 0 then the best codeword has been found */
	if (bestdistortion == 0.0) i = codebooksize;
      }
    } /* the best codeword has been found */

    count[index[bestcodeword]]++;
    celldistortion[index[bestcodeword]] += bestdistortion;
    distortion += bestdistortion;
    /* round the codeword */
    for (i = 0; i < dimension; i++) { /* for double data remove floor */
/* tempvector[i] = (DATA) (floor(codebook[bestcodeword][i] + 0.5));
   this in in error and is changed 9/4/97 */
      tempvector[i] = (DATA) (floor(codebook[index[bestcodeword]][i] + 0.5));
      /* tempvector[i] = (DATA) codebook[index[bestcodeword]][i]; */
    }
    /* write the data to the output file */
    if (fwrite(tempvector,sizeof(DATA),dimension,outputfile) != dimension) {
      fprintf(stderr,"%s: %s: %s\n",programname,outputname,NOWRITE);
      return(-1.0);
    }
  } /* all vectors have been encoded */

  return(distortion);
}