srtree_build.c

srtree算法实现

void cal_centroid(double centroid[], node_type *over_node, node_type *extra_node)
{
 int i, j;
 int total_size;

 if (extra_node->attribute == LEAF) {
 	for (j=0; j < dim; j++) {
		centroid[j] = extra_node->a[j];
		total_size = 1;
	
		for (i=0; i < M; i++) {
			centroid[j] += over_node->ptr[i]->a[j]; 
			total_size++;
		}
	}

	centroid[j] = centroid[j] * 1.0 / total_size;

 } else {

 	for (j=0; j < dim; j++) {
		centroid[j] = extra_node->centroid[j];

		total_size = 1;
		for (i=0; i < M; i++) {
			centroid[j] += over_node->ptr[i]->centroid[j] * over_node->ptr[i]->total_size;
			total_size += over_node->ptr[i]->total_size;
		}

		centroid[j] = centroid[j] * 1.0 / total_size;
	}
 }

 return;

} /* cal_centroid */


void reinsert(node_type *over_node, int over_level, node_type *extra_node, 
node_type *root) 
{
  node_type *node_found;
  node_type **overnode;
  double *value; 
  int *sorted_index;
  int i, start, stop;
  double centroid[dim];

  //printf("h5\n");

  value = (double *)malloc(sizeof(double) * M+1);
  sorted_index = (int *)malloc(sizeof(int) * M+1);

  overnode = (node_type **)malloc((M+1) * sizeof(node_type *));
  for (i=0; i<M; i++) {
        overnode[i] = over_node->ptr[i];
  }

  overnode[M] = extra_node;
  overnode[M]->parent = over_node;

  cal_centroid(centroid, over_node, extra_node);

  if (extra_node->attribute == LEAF) {
  	for (i=0; i<M+1; i++) {
		value[i] = dist_centroid(overnode[i]->a, centroid);
		sorted_index[i] = i;
	}
  } else {
  	for (i=0; i<M+1; i++) {
		//value[i] = Dist2(overnode[i], over_node);
		value[i] = dist_centroid(overnode[i]->centroid, centroid);
		sorted_index[i] = i;
	}
  }

  quicksort(sorted_index, value, 0, M);

  for (i=0; i<dim; i++) {
	over_node->a[i] = overnode[sorted_index[0]]->a[i];
	over_node->b[i] = overnode[sorted_index[0]]->b[i];
  }
  over_node->ptr[0] = overnode[sorted_index[0]];


  stop = M+1-reinsert_p;

  for (i=1; i<stop; i++) {
	over_node->ptr[i] = overnode[sorted_index[i]];
	cal_MBR_node_node(over_node->a, over_node->b, over_node, overnode[sorted_index[i]]);
  }
  for (i=stop; i<M; i++)
	over_node->ptr[i] = NULL;

  over_node->vacancy = reinsert_p - 1;

  adjust_MBR_delete(over_node);

  update_centroid_node(over_node);
  update_radius_node(over_node);

  update_centroid(over_node);

  start = M + 1 - reinsert_p;

  for (i=start; i<M+1; i++) {

	node_found = root;
	
	choose_leaf_level(&node_found, root, 0, overnode[sorted_index[i]], over_level+extra_level);
 
	/* Test whether the node has room or not */
	if(node_found->vacancy!=0) {
 
		/* have room to insert the entry */

		overnode[sorted_index[i]]->parent = node_found;

	        node_found->ptr[M - node_found->vacancy] = overnode[sorted_index[i]];

        	node_found->vacancy--;

	        adjust_MBR(overnode[sorted_index[i]]);

  		update_centroid(overnode[sorted_index[i]]);

	}
	else {
		overflow(node_found, over_level, over_level, overnode[sorted_index[i]], root);
	}

  }        

  //printf("h6\n");

  return;
  
} /* reinsert */




void overflow(node_type *over_node, int over_level, int old_level, node_type
*extra_node, node_type *root)
{
  node_type *node1, *node2;

  //printf("h3\n");

  if (over_level < old_level && over_level != 0) {
	reinsert(over_node, over_level, extra_node, root);

  }
  else {

	tree_node_allocate(&node1);
	tree_node_allocate(&node2);

        split(over_node, extra_node, node1, node2);

        adjust_tree(over_node, over_level, old_level, node1, node2, root);

	//printf("overflow, go out adjust, extra_node->id %d\n", extra_node->id);


  }
        
  //printf("h4\n");
  
  return;
  
}



void insert_node(node_type *root, double *data, int seq_id)
{
  node_type *node_found, *new_node, *data_node;
  int level_found;


  int i;

  /******/
  /* I1 */
  /******/

  //printf("//////////////////////////////////////\n");

  node_found = root;
  extra_level=0;

  tree_node_allocate(&data_node);

  for (i=0; i<dim; i++) {
	data_node->a[i] = data[i];
	data_node->b[i] = data[i];
  }
  data_node->attribute = LEAF;

  level_found = choose_leaf(&node_found, root, 0, data_node);
  /* Now, node_found is the pointer to the leaf chosen */


  /******/
  /* I2 */
  /******/
  
  /********************************************************/
  /* Add record to leaaf node:                            */
  /* If L has room for another entry, install the entry   */
  /* Otherwise invoke split() to split the node */
  /********************************************************/  
  
  /* Make a leaf node */

  tree_leaf_node_allocate(&new_node);


  for(i=0; i < dim; i++) { 
       	new_node->a[i]=data[i];
       	new_node->b[i]=data[i];
  }

  tree_node_deallocate(data_node);

  new_node->id = seq_id;		//data[dim] is the seq_id   
  new_node->attribute=LEAF;
  new_node->vacancy=M;
  new_node->parent=node_found;
  new_node->radius = 0.0;		// Ling
  new_node->total_size = 0;		// Ling
  new_node->ptr = NULL;
  new_node->centroid = NULL;
  /*
  for(i=0; i < M; i++) {
     new_node->ptr[i]=NULL;
     //new_node->child_radius[i] = 0.0;
     //new_node->child_size[i] = 0;
  }
  */

  /* Test whether the node has room or not */
  if(node_found->vacancy!=0) { 


	//printf("insert new->id %d\n", seq_id);

	//printf("h11\n");

	/* have room to insert the entry */

	node_found->ptr[M - node_found->vacancy] = new_node;

      	node_found->vacancy--;
     
	adjust_MBR(new_node);

	update_centroid(new_node);
	//update_centroid(new_node, M-node_found->vacancy);

	//printf("h12\n");
  }
  else {

	//printf("h13\n");

	overflow(node_found, level_found, level_found+1, new_node, root);

	//printf("h14\n");
  }

  //printf("root->radius = %f, vacancy = %d, size = %d\n", root->radius, root->vacancy, root->total_size);

  return;

} /* insert_node */
  




/********************/
/* build_tree():    */
/* build the R-tree */
/********************/
      
void build_tree(node_type **root, double **data, int num_data)
{ 
  int i, j;


  /* make tree root node first */

  tree_node_allocate(root);

  for(i=0; i < dim; i++) { 
	(*root)->b[i] = (double)(-1 * INT_MAX);
       	(*root)->a[i] = (double)(INT_MAX);
  }
  (*root)->id = NO_ID;		//NO_ID
  (*root)->attribute=ROOT;
  (*root)->vacancy=M;   
  (*root)->parent=NULL; 
  (*root)->radius=0.0;
  (*root)->total_size=0;

  for(j=0; j < M; j++) {
     	(*root)->ptr[j]=NULL;
	//(*root)->child_radius[j]=0.0;
	//(*root)->child_size[j]=0;
 }


  /* add data to the tree */
  for(i=0; i<num_data; i++) { 

       #ifdef DEBUG
  	if ((i+1)%10==0)
		printf("insert data %d\n", i+1);
       #endif

      	insert_node(*root, data[i], i);  // i is the seq id.

       #ifndef DEBUG
	free(data[i]);
       #endif
  }

 #ifndef DEBUG
  free(data);
 #endif

} /*build_tree */


int make_data(char *datafile, double ***data)
{ 
  int num_data;
  int i,j;
  FILE *fp_position;

  fp_position=fopen(datafile,"r");  

  num_data = no_histogram;

  (*data) = (double **)malloc(sizeof(double*) * num_data);
  for(i=0; i<num_data; i++)
        (*data)[i] = (double *)malloc(sizeof(double) * dim);

  for(i=0; i<num_data; i++)
	for (j=0; j<dim; j++)
		fscanf(fp_position,"%lf",&((*data)[i][j]));

  fclose(fp_position);

  return(num_data);

} /* make_data */



void write_leaf_node(node_type *node, FILE *fp)
{
  int i;
                
  for (i = 0; i<dim; i++)
        fprintf(fp, "%f\n", (node->a)[i]);
 
  for (i = 0; i<dim; i++)
        fprintf(fp, "%f\n", (node->b)[i]);
 
  fprintf(fp, "%d\n", node->attribute);
  fprintf(fp, "%d\n", node->id);
  fprintf(fp, "%d\n", node->vacancy);

  return;

}


void write_inter_node(node_type *node, FILE *fp)
{
  int i, count;

  for (i = 0; i<dim; i++)
	fprintf(fp, "%f\n", (node->a)[i]);

  for (i = 0; i<dim; i++)  
        fprintf(fp, "%f\n", (node->b)[i]);

  fprintf(fp, "%d\n", node->attribute);

  for (i = 0; i<dim; i++)  
        fprintf(fp, "%f\n", (node->centroid)[i]);

  fprintf(fp, "%d\n", node->vacancy);

  fprintf(fp, "%f\n", node->radius);
  fprintf(fp, "%d\n", node->total_size);

  count = M - node->vacancy;
  for (i=0; i<count; i++) {
	if (node->ptr[i]->attribute != LEAF)
		write_inter_node(node->ptr[i], fp);
	else
		write_leaf_node(node->ptr[i], fp);
  }

  return;

}
  
void check_leaf_node(node_type *node, FILE *fp, int level, char point[])
{
  int i, j;

  for (j=0; j < level; j++)
                fprintf(fp, "  ");

  for (i = 0; i<dim; i++) {
        fprintf(fp, "[%.3f] ", (node->a)[i]);
  }
  fprintf(fp,"\n");

  for (j=0; j < level; j++)
                fprintf(fp, "  ");

  fprintf(fp, "attribute: %d  ", node->attribute);
  fprintf(fp, "ID: %d  ", node->id);
  fprintf(fp, "vacancy: %d\n", node->vacancy);

  point[node->id] = '1';

  return;

}

int within_rectangle(double parent_a[], double parent_b[], double child_a[], double child_b[])
{
 int i;

  for (i=0; i < dim; i++) {
	if (child_a[i] < parent_a[i] || child_a[i] > parent_b[i] ||
		child_b[i] > parent_b[i]) 
		return FALSE;
  }
  return TRUE;
}

/*
int within_sphere(double parent_centroid[], double parent_radius, double child_centroid[], double child_radius)
{
 int i;
 double centroid_centroid = 0.0;

 for (i=0; i < dim; i++) {
	centroid_centroid += pow(parent_centroid[i] - child_centroid[i], 2.0);
 }

 if (centroid_centroid + child_radius > parent_radius)
	return FALSE;
 else
	return TRUE;
}
*/

int within_sphere(double parent_centroid[], double parent_radius, double point[])
{
 int i=dim;
 double centroid_centroid = 0.0;

 while (i--)
        centroid_centroid += pow(parent_centroid[i] - point[i], 2.0);

  centroid_centroid = sqrt(centroid_centroid);

 if (centroid_centroid - _EPSILON > parent_radius) {
        printf("Error: centroid_centroid = %f, parent_radius = %f\n", centroid_centroid, parent_radius);
        fprintf(stderr, "Error: centroid_centroid = %f, parent_radius = %f\n", centroid_centroid, parent_radius);
        return FALSE;
 } else
        return TRUE;
}

void check_inter_node(node_type *node, FILE *fp, int level, char point[], double **data, int total_level)
{
  char local_point[no_histogram];
  int i, count;
  int j;

  memset(local_point, '0', sizeof(local_point));

  for (j=0; j < level; j++)
                fprintf(fp, "  ");
  fprintf(fp, " [a] ");

  for (i = 0; i<dim; i++) {
        fprintf(fp, "[%.3f,", (node->a)[i]);
        fprintf(fp, "%.3f]  ", (node->b)[i]);
  }
  fprintf(fp, "\n");

  for (j=0; j < level; j++)
                fprintf(fp, "  ");
  fprintf(fp, " [centroid] ");

  for (i = 0; i<dim; i++)
        fprintf(fp, "%.3f ", (node->centroid)[i]);

  fprintf(fp, "\n");
  for (j=0; j < level; j++)
                fprintf(fp, "  ");

  fprintf(fp, "attribute: %d  ", node->attribute);
  fprintf(fp, "vacancy: %d  ", node->vacancy);

  fprintf(fp, "radius: %.3f  ", node->radius);
  fprintf(fp, "size: %d\n", node->total_size);

  count = M - node->vacancy;
  for (i=0; i<count; i++) {
        fprintf(fp, "\n");

        if (!within_rectangle(node->a, node->b, node->ptr[i]->a, node->ptr[i]->b)) {
                fprintf(fp, "\nError!  Child %d [%p] of [%p] is out of rectangle bound\n", i, node->ptr[i], node);
                fprintf(stderr, "\nError!  Child %d [%p] of [%p] is out of rectangle bound\n", i, node->ptr[i], node);
        }

        for (j=0; j < level+1; j++)
                fprintf(fp, "--");

        if (node->ptr[i]->attribute != LEAF) {
                fprintf(fp, "[%p] NODE level %d\n", node->ptr[i], level+1);
                check_inter_node(node->ptr[i], fp, level+1, local_point, data, total_level);
        } else {
                fprintf(fp, "[%p] LEAF level %d\n", node->ptr[i], level+1);
                check_leaf_node(node->ptr[i], fp, level+1, local_point);
                if (level+1 != total_level-1) {
                        fprintf(fp, "Error!  leaf level = %d, total level = %d, node = %p\n", level+2, total_level, node);
                        fprintf(stderr, "Error!  leaf level = %d, total level = %d, node = %p\n", level+2, total_level, node);
                }
        }
  }

  for (i=0; i < no_histogram; i++) {
        if (local_point[i] == '1') {
                within_sphere(node->centroid, node->radius, data[i]);
                point[i] = '1';
        }
  }

  return;

}


void print_node_address(node_type *node, FILE *fp)
{
  int count = M - node->vacancy;
  int i;

  fprintf(fp, "[%p] : ", node);
  for (i=0; i<count; i++)
        fprintf(fp, "[%p][%d]%d ", node->ptr[i], node->ptr[i]->id, node->ptr[i]->attribute);
  fprintf(fp, "\n");

  if (node->ptr[0]->attribute != LEAF)
        for (i=0; i < count; i++)
                print_node_address(node->ptr[i], fp);

  return;
}

void check_srtree(node_type *root, double **data)
{
  int i;
  char point[no_histogram];
  int total_level=1;
  FILE *fp;
  node_type *node;

  memset(point, '0', sizeof(point));

  node = root;
  while (node->attribute != LEAF) {
        total_level++;
        node = node->ptr[0];
  }

  if (!(fp=fopen(CHECK_TREE_LOG, "w"))) {
        printf("Can't open file, %s\n", CHECK_TREE_LOG);
        exit(EXIT_FAILURE);
  }

  fprintf(fp, "\nTree height: %d\n", total_level);

  fprintf(fp, "\nROOT [%p]\n", root);
  check_inter_node(root, fp, 0, point, data, total_level);

  for (i=0; i < no_histogram; i++) {
        if (point[i] == '1') {
                within_sphere(root->centroid, root->radius, data[i]);
        }
  }

  fprintf(fp, "\n");

  print_node_address(root, fp);

  fclose(fp);

  return;

}


void save_srtree(node_type *root, char save_tree_file[])
{
  FILE *fp;

  //fp = fopen(SAVE_SRTREE_FILE, "w");
  fp = fopen(save_tree_file, "w");
  if (!fp) {
	printf("Can't write to %s\n", save_tree_file);
	exit(EXIT_FAILURE);
  }
  
  write_inter_node(root, fp);

  fclose(fp);
 
  return;

}


void free_tree(node_type *node)
{
  int i;

  if (!node) return;

  free(node->a);
  free(node->b);

  if (node->attribute != LEAF)
  	free(node->centroid);

  for (i=0; i < M-node->vacancy; i++)
	free_tree(node->ptr[i]);

  free(node);

  return;
}


void destroy(double **data)
{
 int i;

 for (i=0; i < no_histogram; i++)
	free(data[i]);
 free(data);

 return;
}

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

  int num_data;
  double **data;
  config_type config;
  node_type *root;
  float  userTime, sysTime;
  struct rusage myTime1, myTime2;

  ////////////////////////////
  if (argc != 2) {
	printf("Build SR-tree\n\n");
        printf("Usage: %s <config>\n", argv[0]);
        exit(EXIT_FAILURE);
  }
  ////////////////////////////

  getrusage(RUSAGE_SELF,&myTime1);

  printf("initialize\n");
  initialize(&config, argv[1]);   

  printf("make_data\n");
  num_data=make_data(config.datafile, &data);

  printf("build_tree\n");
  build_tree(&root, data, num_data);

  getrusage(RUSAGE_SELF,&myTime2);

  ////////////////////////////

  #ifdef DEBUG
  printf("check_srtree\n");
  check_srtree(root, data);
  destroy(data);
  #endif

  printf("save_srtree\n");
  save_srtree(root, config.save_tree_file);

  printf("free_tree\n");
  free_tree(root);

  ////////////////////////////

  userTime =
                ((float) (myTime2.ru_utime.tv_sec  - myTime1.ru_utime.tv_sec)) +
                ((float) (myTime2.ru_utime.tv_usec - myTime1.ru_utime.tv_usec)) * 1e-6;

   sysTime =
                ((float) (myTime2.ru_stime.tv_sec  - myTime1.ru_stime.tv_sec)) +
                ((float) (myTime2.ru_stime.tv_usec - myTime1.ru_stime.tv_usec)) * 1e-6;

  fprintf(stdout, "User time : %f seconds\n",userTime);
  fprintf(stdout, "System time : %f seconds\n",sysTime);
  fprintf(stdout, "Total time : %f seconds\n",userTime+sysTime);

  return(0);  
} /* main */