int reinsertCount = numClusters * info.getReinsertFactor() / 100;
    HnPoint centroid;
    int i, axis;
    HnSRTreeReinsertionEntry *entries;

    /* initialize flags */
    if ( flags != NULL ) {
	delete[] flags;
    }
    flags = new ReinsertFlag[numClusters];

    for ( i=0; i<numClusters; i++ ) {
	flags[i] = REINSERT_NONE;
    }

    /* calculate centroid */
    centroid = new_HnPoint(dimension);

    for ( axis=0; axis<dimension; axis++ ) {
	double sum = 0;
	int total = 0;

	for ( i=0; i<numClusters; i++ ) {
	    double coord = clusters[i].getCentroid().getCoordAt(axis);
	    int weight = clusters[i].getWeight();

	    sum += coord * weight;
	    total += weight;
	}

	centroid.setCoordAt(sum / total, axis);
    }

    /* initialize entries */
    entries = (HnSRTreeReinsertionEntry *)
	HnMalloc(sizeof(HnSRTreeReinsertionEntry) * numClusters);

    for ( i=0; i<numClusters; i++ ) {
	entries[i].index = i;
	entries[i].distance = clusters[i].getCentroid().getDistance(centroid);
    }

    /* sort clusters by distance in descent order */
    qsort(entries, numClusters, sizeof(HnSRTreeReinsertionEntry),
	  HnSRTreeCompareReinsertionEntries);

    /* make reinsert flags */
    i=0;
    while ( i<reinsertCount ) {
	flags[entries[i].index] = REINSERT;
	i++;
    }
    while ( i<numClusters ) {
	flags[entries[i].index] = STAY;
	i++;
    }

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::selectClusters:\n");
	fprintf(stderr, "    REINSERT\n");
	for ( i=0; i<numClusters; i++ ) {
	    if ( flags[i] == REINSERT ) {
		fprintf(stderr,
			"    clusters[%d].getCentroid()."
			"getDistance(centroid) = %g\n",
			i, clusters[i].getCentroid().getDistance(centroid));
	    }
	}
	fprintf(stderr, "    STAY\n");
	for ( i=0; i<numClusters; i++ ) {
	    if ( flags[i] == STAY ) {
		fprintf(stderr,
			"    clusters[%d].getCentroid()."
			"getDistance(centroid) = %g\n",
			i, clusters[i].getCentroid().getDistance(centroid));
	    }
	}
    }

    *flags_return = flags;

    HnFree(entries, sizeof(HnSRTreeReinsertionEntry) * numClusters);
}

/*
 * Split
 */

void
HnSRTreeFileObj::splitLeaf(HnSRTreeStack &stack,
			   const HnPoint &newPoint,
			   const HnDataItem &newDataItem)
{
    HnSRTreeLeaf leaf, left, right;
    int i;
    HnPointArray points;
    SplitFlag *flags;

    leaf = stack.getTopLeaf();

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::splitLeaf: leaf = %s\n",
		(char *)leaf.toString());
    }

    /* put points into an array */
    points = new_HnPointArray(leaf.getCount() + 1);
    for ( i=0; i<leaf.getCount(); i++ ) {
	points[i] = leaf.getPointAt(i);
    }
    points[i] = newPoint;

    /* split points */
    splitPoints(points, &flags);

    /* put points into the left and right leaves */
    left = new_HnSRTreeLeaf(info, leaf.getOffset());
    right = allocateLeaf();

    for ( i=0; i<leaf.getCount(); i++ ) {
	switch ( flags[i] ) {
	case LEFT:
	    left.addElement(leaf.getPointAt(i), leaf.getDataItemAt(i));
	    break;
	case RIGHT:
	    right.addElement(leaf.getPointAt(i), leaf.getDataItemAt(i));
	    break;
	default:
	    HnAbort("split flag is incorrectly assigned.");
	}
    }
    switch ( flags[i] ) {
    case LEFT:
	left.addElement(newPoint, newDataItem);
	break;
    case RIGHT:
	right.addElement(newPoint, newDataItem);
	break;
    default:
	HnAbort("split flag is incorrectly assigned.");
    }

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::splitLeaf: left  = %s\n",
		(char *)left.toString());
	fprintf(stderr, "HnSRTreeFileObj::splitLeaf: right = %s\n",
		(char *)right.toString());
    }

    writeLeaf(left);
    writeLeaf(right);

    stack.pop();

    updateNode(stack,
	       left.getCluster(), left.getOffset(),
	       right.getCluster(), right.getOffset());

    profile->numLeafSplits ++;
}

void
HnSRTreeFileObj::splitNode(HnSRTreeStack &stack,
			   const HnSRTreeCluster &newCluster, long newOffset)
{
    HnSRTreeNode node, left, right;
    int i;
    HnSRTreeClusterArray clusters;
    SplitFlag *flags;

    node = stack.getTopNode();

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::splitNode: node = %s\n",
		(char *)node.toString());
    }

    /* put clusters into an array */
    clusters = new_HnSRTreeClusterArray(node.getCount() + 1);
    for ( i=0; i<node.getCount(); i++ ) {
	clusters[i] = node.getClusterAt(i);
    }
    clusters[i] = newCluster;

    /* split clusters */
    splitClusters(clusters, &flags);

    /* put clusters into the left and right nodes */
    left = new_HnSRTreeNode(info, node.getOffset());
    right = allocateNode();

    for ( i=0; i<node.getCount(); i++ ) {
	switch ( flags[i] ) {
	case LEFT:
	    left.addElement(node.getClusterAt(i), node.getOffsetAt(i));
	    break;
	case RIGHT:
	    right.addElement(node.getClusterAt(i), node.getOffsetAt(i));
	    break;
	default:
	    HnAbort("split flag is incorrectly assigned.");
	}
    }
    switch ( flags[i] ) {
    case LEFT:
	left.addElement(newCluster, newOffset);
	break;
    case RIGHT:
	right.addElement(newCluster, newOffset);
	break;
    default:
	HnAbort("split flag is incorrectly assigned.");
    }

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::splitNode: left  = %s\n",
		(char *)left.toString());
	fprintf(stderr, "HnSRTreeFileObj::splitNode: right = %s\n",
		(char *)right.toString());
    }

    writeNode(left);
    writeNode(right);

    stack.pop();

    updateNode(stack,
	       left.getCluster(), left.getOffset(),
	       right.getCluster(), right.getOffset());

    profile->numNodeSplits ++;
}

void
HnSRTreeFileObj::updateNode(HnSRTreeStack &stack,
			    const HnSRTreeCluster &leftCluster,
			    long leftOffset,
			    const HnSRTreeCluster &rightCluster,
			    long rightOffset)
{
    HnSRTreeNode node;
    int cursor;

    if ( stack.getDepth() == 0 ) {
	extendTree(leftCluster, leftOffset, rightCluster, rightOffset);
    }
    else {
	node = stack.getTopNode();
	cursor = stack.getCursor();

	node.setElementAt(leftCluster, leftOffset, cursor);

	if ( node.isFull() ) {
	    int reinsertCount =
		(node.getCount() + 1) * info.getReinsertFactor() / 100;

	    if ( stack.getDepth() == 1 ) {
		/* the parent is the root node */
		splitNode(stack, rightCluster, rightOffset);
	    }
	    else if ( reinsertCount <= 0 ) {
		/* no entry needs to be reinserted */
		splitNode(stack, rightCluster, rightOffset);
	    }
	    else {
		/* reinsert some of entries */
		int index = -1;

		if ( (index =
		      reinsertedBlocks.indexOf(node.getOffset())) < 0 ) {
		    reinsertedBlocks.addElement(node.getOffset());
		    reinsertNode(stack, rightCluster, rightOffset);
		}
		else {
		    reinsertedBlocks.removeElementAt(index);
		    splitNode(stack, rightCluster, rightOffset);
		}
	    }
	}
	else {
	    node.addElement(rightCluster, rightOffset);

	    writeNode(node);
	    updateCluster(stack);
	}
    }
}

void
HnSRTreeFileObj::extendTree(const HnSRTreeCluster &leftCluster,
			    long leftOffset,
			    const HnSRTreeCluster &rightCluster,
			    long rightOffset)
{
    HnSRTreeNode node;

    node = allocateNode();
    node.addElement(leftCluster, leftOffset);
    node.addElement(rightCluster, rightOffset);
    writeNode(node);

    info.setRootOffset(node.getOffset());
    info.setHeight(info.getHeight() + 1);
    writeSuperBlock();
}

/*
 * Split
 */

struct HnSRTreeAxisSortEntry {
    int index;
    double coord;
};

static int
HnSRTreeCompareAxisSortEntries(const void *ptr1, const void *ptr2)
{
    HnSRTreeAxisSortEntry *entry1 = (HnSRTreeAxisSortEntry *)ptr1;
    HnSRTreeAxisSortEntry *entry2 = (HnSRTreeAxisSortEntry *)ptr2;

    if ( entry1->coord == entry2->coord ) {
	return 0;
    }
    else {
	if ( entry1->coord < entry2->coord ) {
	    return -1;
	}
	else {
	    return 1;
	}
    }
}

/*
 * Split (leaf)
 */

void
HnSRTreeFileObj::splitPoints(const HnPointArray &points,
			     SplitFlag **flags_return)
{
    static SplitFlag *flags = NULL;
    int numPoints = points.length;
    int dimension = info.getDimension();
    int minCount, maxCount;
    struct {
	int axis;
	double var;
    } max;
    struct {
	int count;
	double leftVar, rightVar;
    } min;
    int axis, count;
    int i;
    HnSRTreeAxisSortEntry *entries;

    /* initialize minCount and maxCount */
    if ( (minCount = numPoints * info.getSplitFactor() / 100) == 0 ) {
	minCount = 1;
    }
    maxCount = numPoints - minCount;

    /* initialize flags */
    if ( flags != NULL ) {
	delete[] flags;
    }
    flags = new SplitFlag[numPoints];

    for ( i=0; i<numPoints; i++ ) {
	flags[i] = SPLIT_NONE;
    }

    max.axis = -1;
    max.var = -1;

    /* calculate variance */
    for ( axis=0; axis<dimension; axis++ ) {
	double avg, var;
	double sum = 0;
	double sum2 = 0;

	for ( i=0; i<numPoints; i++ ) {
	    double coord = points[i].getCoordAt(axis);
	    sum += coord;
	    sum2 += coord * coord;
	}

	avg = sum / numPoints;
	var = sum2 / numPoints - avg * avg;

	if ( max.axis == -1 ) {
	    max.axis = axis;
	    max.var = var;
	}
	else {
	    if ( var > max.var ) {
		max.axis = axis;
		max.var = var;
	    }
	}
    }

    /* sort points along the max variance axis */
    entries = (HnSRTreeAxisSortEntry *)
	HnMalloc(sizeof(HnSRTreeAxisSortEntry) * numPoints);

    for ( i=0; i<numPoints; i++ ) {
	entries[i].index = i;
	entries[i].coord = points[i].getCoordAt(max.axis);
    }

    qsort(entries, numPoints, sizeof(HnSRTreeAxisSortEntry),
	  HnSRTreeCompareAxisSortEntries);

    /* choose split point */
    min.count = -1;
    min.leftVar = 0;
    min.rightVar = 0;

    for ( count=minCount; count<=maxCount; count++ ) {
	double leftVar, rightVar;

	/* calculate the sum of variances on the left side */
	leftVar = 0;

	for ( axis=0; axis<dimension; axis++ ) {
	    double sum, sum2, avg, var;

	    sum = 0;
	    sum2 = 0;

	    for ( i=0; i<count; i++ ) {
		double coord = points[entries[i].index].getCoordAt(axis);
		sum += coord;
		sum2 += coord * coord;
	    }

	    avg = sum / count;
	    var = sum2 / count - avg * avg;

	    leftVar += var;
	}

	/* calculate the sum of variances on the right side */
	rightVar = 0;

	for ( axis=0; axis<dimension; axis++ ) {
	    double sum, sum2, avg, var;

	    sum = 0;
	    sum2 = 0;

	    for ( i=count; i<numPoints; i++ ) {
		double coord = points[entries[i].index].getCoordAt(axis);
		sum += coord;
		sum2 += coord * coord;
	    }
	    avg = sum / (numPoints - count);
	    var = sum2 / (numPoints - count) - avg * avg;

	    rightVar += var;
	}

	if ( min.count < 0 ) {
	    min.count = count;
	    min.leftVar = leftVar;
	    min.rightVar = rightVar;
	}
	else {
	    if ( leftVar + rightVar < min.leftVar + min.rightVar ) {
		min.count = count;
		min.leftVar = leftVar;
		min.rightVar = rightVar;
	    }
	}
    }

    /* make split flags */
    i=0;
    while ( i<min.count ) {
	flags[entries[i].index] = LEFT;
	i++;
    }
    while ( i<numPoints ) {
	flags[entries[i].index] = RIGHT;
	i++;
    }

    if ( debug ) {
	fprintf(stderr, "HnSRTreeFileObj::splitPoints\n");
	fprintf(stderr, "    LEFT\n");
	for ( i=0; i<numPoints; i++ ) {
	    if ( flags[i] == LEFT )
		fprintf(stderr,
			"    points[%d].getCoordAt(%d) = %g\n",
			i, max.axis, points[i].getCoordAt(max.axis));
	}
	fprintf(stderr, "    RIGHT\n");
	for ( i=0; i<numPoints; i++ ) {
	    if ( flags[i] == RIGHT )
		fprintf(stderr,
			"    points[%d].getCoordAt(%d) = %g\n",
			i, max.axis, points[i].getCoordAt(max.axis));
	}
    }

    *flags_return = flags;

    HnFree(entries, sizeof(HnSRTreeAxisSortEntry) * numPoints);
}

/*
 * Split (node)
 */

void
HnSRTreeFileObj::splitClusters(const HnSRTreeClusterArray &clusters,
			       SplitFlag **flags_return)
{
    static SplitFlag *flags = NULL;
    int numClusters = clusters.length;
    int dimension = info.getDimension();
    int minCount, maxCount;
    struct {
	int axis;
	double var;
    } max;
    struct {
	int count;
	double leftVar, rightVar;
    } min;
    int axis, count;
    int i;
    HnSRTreeAxisSortEntry *entries;

    /* initialize minCount and maxCount */
    if ( (minCount = numClusters * info.getSplitFactor() / 100) == 0 ) {
	minCount = 1;
    }