tsvq_util.c

Vector Quantization压缩算法

      node->right_child->data[j] = rightcentroid[j];
      rightcentroid[j] = 0.0;
    }
    leftcount = 0;
    rightcount = 0;
    leftdist = 0.0;
    rightdist = 0.0;

    /* do a lloyd iteration */
    for (i = 0; i < node->count; i++) {
      current_leftdist = dist(node->trainseq[i],node->left_child->data);
      current_rightdist = dist(node->trainseq[i],node->right_child->data);
      if (current_leftdist < current_rightdist) {
	for (j = 0; j < dim; j++) {
	  leftcentroid[j] += ((DISTTYPE) node->trainseq[i][j]);
	}
	leftcount += 1;
	leftdist += current_leftdist;
      }
      else {
	for (j = 0; j < dim; j++) {
	  rightcentroid[j] += ((DISTTYPE) node->trainseq[i][j]);
	}
	rightcount += 1;
	rightdist += current_rightdist;
      }
    }

    /* find average mse distortion */
    if (leftcount != 0)
      leftdist /= ((DISTTYPE) leftcount);
    if (rightcount != 0)
      rightdist /= ((DISTTYPE) rightcount);

    /* save old distortion, and find new distortion */
    past_dist = current_dist;
    current_dist = (leftdist * ((DISTTYPE) leftcount)/((DISTTYPE) node->count))
      + (rightdist * ((DISTTYPE) rightcount) / ((DISTTYPE) node->count));

    /* find the centroids for the next iteration */
    for (j = 0; j < dim; j++) {
      if (leftcount != 0)
	leftcentroid[j] /= ((DISTTYPE) leftcount);
      if (rightcount != 0)
	rightcentroid[j] /= ((DISTTYPE) rightcount);
    }

    /* if the distortion is 0, stop trying to split the node */
    if (current_dist == 0.0) break;

    /* if alternate_flag is true, a forced split has already been attempted,
       do not try it again since an infinite loop would occur */
    if ((leftcount == 0 || rightcount == 0) && alternate_flag) {
      /* do nothing, will exit when current_dist = past_dist */
    }

    /* else if the distortion is not 0 and the split was not successful,
       force a better split */
    else if (leftcount == 0 || rightcount == 0) {
      split_alternate(node,leftcentroid,rightcentroid);
      alternate_flag = TRUE;
      past_dist = HUGE;
    }

  }  while ( !(((past_dist - current_dist)/past_dist) <= thresh) );
  /* finished with iterations */

  /* save the avmse and count of the node's children */
  node->left_child->avmse = leftdist;
  node->right_child->avmse = rightdist;
  node->left_child->count = leftcount;
  node->right_child->count = rightcount;

  /* prepare to store the training sequence among the children */
  if (!(node->left_child->trainseq =
	(DATATYPE **)calloc(leftcount,sizeof(DATATYPE *))) ||
      !(node->right_child->trainseq =
	(DATATYPE **)calloc(rightcount,sizeof(DATATYPE *)))) {
    fprintf(stderr,"%s: %s\n",programname,NOMEMORY);
    return(FALSE);
  }

  /* divide the training set among the node's children */
  leftcount = 0;
  rightcount = 0;
  for (i = 0; i < node->count; i++) {
    if (dist(node->trainseq[i],node->left_child->data) <
	dist(node->trainseq[i],node->right_child->data)) {
      node->left_child->trainseq[leftcount] = node->trainseq[i];
      leftcount += 1;
    }
    else {
      node->right_child->trainseq[rightcount] = node->trainseq[i];
      rightcount += 1;
    }
  }

  /* When the convergence threshold thresh is greater than zero,the codeword
     of a given voronoi region is not the average of the training vectors in
     the voronoi region.  Hence it is possible that the Lloyd iteration may
     stop with a voronoi region that has one or more identical training
     vectors, but that the codeword may not be the same as the training
     vector.  If a tree that completely resolves the training data is to be
     grown, i.e. one leaf for each distinct codeword, then it may not be
     possible to create a tree with zero distortion.  This is solved by
     allowing any node with a non-zero distortion to split, even if
     the node's count is only one, since there is still a positive slope
     that represents the goodness of the split.  Allowing the node
     to split may cause the growing algorithm to add an empty cell to
     the tree.

     There may be another problem when the convergence threshold is greater
     than zero.  It is possible that a node with a non-zero avmse has failed
     to split and has a child node with the parent node's training vectors
     and a lower avmse. This would result in a positive slope which may cause
     the growing algorithm to add an empty cell to the tree. This
     would happen rarely and only if the threshold is too large.

     Since the option of a changing threshold is still allowed, a test
     is inserted here to check for an empty cell problem.  If an
     empty cell is found, the user is notified of the problem and the
     use of a smaller threshold is recommended. An empty cell problem
     does not cause the program to terminate. */

  if ((rightcount == 0 && node->left_child->avmse < node->avmse) ||
      (leftcount == 0 && node->right_child->avmse < node->avmse)) {
    fprintf(stderr,"%s: %s: %s\n",programname,EMPTY_CELL,SMALL_THRESH);
  }

  free((char *) node->trainseq);
  free((char *) leftcentroid);
  free((char *) rightcentroid);
  return(TRUE);
}

void split_node(node,leftcentroid,rightcentroid)
     TreeNode *node;
     DISTTYPE *leftcentroid;
     DISTTYPE *rightcentroid;
{
  int i;

  for (i = 0; i < dim; i++) {
    leftcentroid[i] = node->data[i];
    rightcentroid[i] = (node->data[i]) * (1 + mult_offset);
  }
}

void split_alternate(node,leftcentroid,rightcentroid)
     TreeNode *node;
     DISTTYPE *leftcentroid;
     DISTTYPE *rightcentroid;
{
  long   i, j;
  double maxdist = 0.0;
  double tempdist;

  /* find the training vector that is the furthest from the centroid */
  j = 0;
  for (i = 0; i < node->count; i++) {
    if ((tempdist = dist(node->trainseq[i],node->data)) > maxdist) {
      j = i;
      maxdist = tempdist;
    }
  }

  /* copy new codewords */
  for (i = 0; i < dim; i++) {
    leftcentroid[i] = node->data[i];
    rightcentroid[i] = (DISTTYPE) node->trainseq[j][i];
  }
}

long tree_clean(root)
     TreeNode *root;
{
  TreeNode *node;
  long     n;

  if(!(node = newnode())) {
    fprintf(stderr,"%s: %s\n",programname,NOMEMORY);
    return(0);
  }
  node = root;

  /* count the number of designed nodes in the tree, and get rid of those
     nodes that are not designed */
  n = 0;
  for ( ; ; ) {
    /* determine node position and from this find the next node to use */
    if (node->left_child == NULL) { /* node is terminal */

      if(node->right_child != NULL) {  /* test tree fidelity */
        fprintf(stderr,"%s: %s\n",programname,NOTREE);
        return(0);
      }

      if(node->designed) n++; /* increment count of designed nodes */
      else { /* node is not designed, get rid of it and its sibling */
	node = node->parent;

	node->left_child = NULL;
	node->right_child = NULL;
      }

      /* find the next node to examine */
      while ((node != root) && (node == node->parent->right_child)) {
        node = node->parent; /* continue with the leaf node's parent */
      }

      /* tree search complete */
      if (node == root) break;

      else {
        node = node->parent->right_child;
      }
    }

    else { /* node has a child, go to left child */

      /* test tree fidelity, a non-terminal node should be designed and
         should have two childen */
      if(node->right_child == NULL || !(node->designed)) {
        fprintf(stderr,"%s: %s\n",programname,NOTREE);
        return(0);
      }
      n++; /* increment count of designed nodes */
      node = node->left_child; /* search the next node of the tree */
    }
  }

  /* return a count of the number of designed nodes */
  return(n);
}