This repository has been archived on 2023-08-20. You can view files and clone it, but cannot push or open issues or pull requests.
yap-6.3/packages/udi/rtree.c

525 lines
11 KiB
C
Raw Normal View History

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <float.h>
#include "rtree.h"
static node_t RTreeNewNode (void);
static void RTreeDestroyNode (node_t);
static void RTreeNodeInit (node_t);
static int RTreeSearchNode (node_t, rect_t, SearchHitCallback, void *);
static int RTreeInsertNode (node_t, int, rect_t,void *,node_t *);
static int RTreePickBranch (rect_t, node_t);
static int RTreeAddBranch(node_t, branch_t, node_t *);
static void RTreeSplitNode (node_t, branch_t, node_t *);
static void RTreePickSeeds(partition_t *, node_t, node_t);
static void RTreeNodeAddBranch(rect_t *, node_t, branch_t);
static void RTreePickNext(partition_t *, node_t, node_t);
static rect_t RTreeNodeCover(node_t);
static double RectArea (rect_t);
static rect_t RectCombine (rect_t, rect_t);
static int RectOverlap (rect_t, rect_t);
static void RectPrint (rect_t);
static partition_t PartitionNew (void);
static void PartitionPush (partition_t *, branch_t);
static branch_t PartitionPop (partition_t *);
static branch_t PartitionGet (partition_t *, int);
rtree_t RTreeNew (void)
{
rtree_t t;
t = RTreeNewNode();
t->level = 0; /*leaf*/
return t;
}
void RTreeDestroy (rtree_t t)
{
if (t)
RTreeDestroyNode (t);
}
static node_t RTreeNewNode (void)
{
node_t n;
n = (node_t) malloc (sizeof(*n));
assert(n);
RTreeNodeInit(n);
return n;
}
static void RTreeDestroyNode (node_t node)
{
int i;
if (node->level == 0) /* leaf level*/
{
for (i = 0; i < MAXCARD; i++)
if (node->branch[i].child)
;/* allow user free data*/
else
break;
}
else
{
for (i = 0; i < MAXCARD; i++)
if (node->branch[i].child)
RTreeDestroyNode (node->branch[i].child);
else
break;
}
free (node);
}
static void RTreeNodeInit (node_t n)
{
memset((void *) n,0, sizeof(*n));
n->level = -1;
}
int RTreeSearch (rtree_t t, rect_t s, SearchHitCallback f, void *arg)
{
assert(t);
return RTreeSearchNode(t,s,f,arg);
}
static int RTreeSearchNode (node_t n, rect_t s, SearchHitCallback f, void *arg)
{
int i;
int c = 0;
if (n->level > 0)
{
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child &&
RectOverlap (s,n->branch[i].mbr))
c += RTreeSearchNode ((node_t) n->branch[i].child, s, f, arg);
}
else
{
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child &&
RectOverlap (s,n->branch[i].mbr))
{
c ++;
if (f)
if ( !f(n->branch[i].mbr,n->branch[i].child,arg))
return c;
}
}
return c;
}
void RTreeInsert (rtree_t *t, rect_t r, void *data)
{
node_t n2;
node_t new_root;
branch_t b;
assert(t && *t);
if (RTreeInsertNode(*t, 0, r, data, &n2))
/* deal with root split */
{
new_root = RTreeNewNode();
new_root->level = (*t)->level + 1;
b.mbr = RTreeNodeCover(*t);
b.child = (void *) *t;
RTreeAddBranch(new_root, b, NULL);
b.mbr = RTreeNodeCover(n2);
b.child = (void *) n2;
RTreeAddBranch(new_root, b, NULL);
*t = new_root;
}
}
static int RTreeInsertNode (node_t n, int level,
rect_t r, void *data,
node_t *new_node)
{
int i;
node_t n2;
branch_t b;
assert(n && new_node);
assert(level >= 0 && level <= n->level);
if (n->level > level)
{
i = RTreePickBranch(r,n);
if (!RTreeInsertNode((node_t) n->branch[i].child, level,
r, data,&n2)) /* not split */
{
n->branch[i].mbr = RectCombine(r,n->branch[i].mbr);
return FALSE;
}
else /* node split */
{
n->branch[i].mbr = RTreeNodeCover(n->branch[i].child);
b.child = n2;
b.mbr = RTreeNodeCover(n2);
return RTreeAddBranch(n, b, new_node);
}
}
else /*insert level*/
{
b.mbr = r;
b.child = data;
return RTreeAddBranch(n, b, new_node);
}
}
static int RTreeAddBranch(node_t n, branch_t b, node_t *new_node)
{
int i;
assert(n);
if (n->count < MAXCARD) /*split not necessary*/
{
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child == NULL)
{
n->branch[i] = b;
n->count ++;
break;
}
return FALSE;
}
else /*needs to split*/
{
assert(new_node);
RTreeSplitNode (n, b, new_node);
return TRUE;
}
}
static int RTreePickBranch (rect_t r, node_t n)
{
int i;
double area;
double inc_area;
rect_t tmp;
int best_i;
double best_inc;
double best_i_area;
best_i = 0;
best_inc = DBL_MAX; /* double Max value */
best_i_area = DBL_MAX;
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child)
{
area = RectArea (n->branch[i].mbr);
tmp = RectCombine (r, n->branch[i].mbr);
inc_area = RectArea (tmp) - area;
if (inc_area < best_inc)
{
best_inc = inc_area;
best_i = i;
best_i_area = area;
}
else if (inc_area == best_inc && best_i_area > area)
{
best_inc = inc_area;
best_i = i;
best_i_area = area;
}
}
else
break;
return best_i;
}
static void RTreeSplitNode (node_t n, branch_t b, node_t *new_node)
{
partition_t p;
int level;
int i;
assert(n);
assert(new_node);
p = PartitionNew();
for (i = 0; i < MAXCARD; i ++)
PartitionPush(&p,n->branch[i]);
PartitionPush(&p,b);
level = n->level;
RTreeNodeInit(n);
n->level = level;
*new_node = RTreeNewNode();
(*new_node)->level = level;
RTreePickSeeds(&p, n, *new_node);
while (p.n)
if (n->count + p.n <= MINCARD)
/* first group (n) needs all entries */
RTreeNodeAddBranch(&(p.cover[0]), n, PartitionPop(&p));
else if ((*new_node)->count + p.n <= MINCARD)
/* second group (new_node) needs all entries */
RTreeNodeAddBranch(&(p.cover[1]), *new_node, PartitionPop(&p));
else
RTreePickNext(&p, n, *new_node);
}
static void RTreePickNext(partition_t *p, node_t n1, node_t n2)
/* linear version */
{
branch_t b;
double area[2], inc_area[2];
rect_t tmp;
b = PartitionPop(p);
area[0] = RectArea (p->cover[0]);
tmp = RectCombine (p->cover[0], b.mbr);
inc_area[0] = RectArea (tmp) - area[0];
area[1] = RectArea (p->cover[1]);
tmp = RectCombine (p->cover[1], b.mbr);
inc_area[1] = RectArea (tmp) - area[1];
if (inc_area[0] < inc_area[1] ||
(inc_area[0] == inc_area[1] && area[0] < area[1]))
RTreeNodeAddBranch(&(p->cover[0]),n1,b);
else
RTreeNodeAddBranch(&(p->cover[1]),n2,b);
}
static void RTreePickSeeds(partition_t *p, node_t n1, node_t n2)
/* puts in index 0 of each node the resulting entry, forming the two
groups
This is the linear version
*/
{
int dim,high, i;
int highestLow[NUMDIMS], lowestHigh[NUMDIMS];
double width[NUMDIMS];
int seed0, seed1;
double sep, best_sep;
assert(p->n == MAXCARD + 1);
for (dim = 0; dim < NUMDIMS; dim++)
{
high = dim + NUMDIMS;
highestLow[dim] = lowestHigh[dim] = 0;
for (i = 1; i < MAXCARD +1; i++)
{
if (p->buffer[i].mbr.coords[dim] >
p->buffer[highestLow[dim]].mbr.coords[dim])
highestLow[dim] = i;
if (p->buffer[i].mbr.coords[high] <
p->buffer[lowestHigh[dim]].mbr.coords[high])
lowestHigh[dim] = i;
}
width[dim] = p->cover_all.coords[high] - p->cover_all.coords[dim];
assert(width[dim] >= 0);
}
seed0 = lowestHigh[0];
seed1 = highestLow[0];
best_sep = 0;
for (dim = 0; dim < NUMDIMS; dim ++)
{
high = dim + NUMDIMS;
sep = (p->buffer[highestLow[dim]].mbr.coords[dim] -
p->buffer[lowestHigh[dim]].mbr.coords[high]) / width[dim];
if (sep > best_sep)
{
seed0 = lowestHigh[dim];
seed1 = highestLow[dim];
best_sep = sep;
}
}
/* assert (seed0 != seed1); */
if (seed0 > seed1)
{
RTreeNodeAddBranch(&(p->cover[0]),n1,PartitionGet(p,seed0));
RTreeNodeAddBranch(&(p->cover[1]),n2,PartitionGet(p,seed1));
}
else if (seed0 < seed1)
{
RTreeNodeAddBranch(&(p->cover[0]),n1,PartitionGet(p,seed1));
RTreeNodeAddBranch(&(p->cover[1]),n2,PartitionGet(p,seed0));
}
}
static void RTreeNodeAddBranch(rect_t *r, node_t n, branch_t b)
{
int i;
assert(n);
assert(n->count < MAXCARD);
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child == NULL)
{
n->branch[i] = b;
n->count ++;
break;
}
*r = RectCombine(*r,b.mbr);
}
void RTreePrint(node_t t)
{
int i;
/* printf("rtree([_,_,_,_,_]).\n"); */
printf("rtree(%p,%d,[",t,t->level);
for (i = 0; i < MAXCARD; i++)
{
if (t->branch[i].child != NULL)
{
printf("(%p,",t->branch[i].child);
RectPrint(t->branch[i].mbr);
printf(")");
}
else
{
printf("nil");
}
if (i < MAXCARD-1)
printf(",");
}
printf("]).\n");
if (t->level != 0)
for (i = 0; i < MAXCARD; i++)
if (t->branch[i].child != NULL)
RTreePrint((node_t) t->branch[i].child);
else
break;
}
/*
* Partition related
*/
static partition_t PartitionNew (void)
{
partition_t p;
memset((void *) &p,0, sizeof(p));
p.cover[0] = p.cover[1] = p.cover_all = RectInit();
return p;
}
static void PartitionPush (partition_t *p, branch_t b)
{
assert(p->n < MAXCARD + 1);
p->buffer[p->n] = b;
p->n ++;
p->cover_all = RectCombine(p->cover_all,b.mbr);
}
static branch_t PartitionPop (partition_t *p)
{
assert(p->n > 0);
p->n --;
return p->buffer[p->n];
}
static branch_t PartitionGet (partition_t *p, int n)
{
branch_t b;
assert (p->n > n);
b = p->buffer[n];
p->buffer[n] = PartitionPop(p);
return b;
}
/*
* Rect related
*/
rect_t RectInit (void)
{
rect_t r = {{DBL_MAX, DBL_MAX, DBL_MIN, DBL_MIN}};
return (r);
}
static double RectArea (rect_t r)
{
int i;
double area;
for (i = 0,area = 1; i < NUMDIMS; i++)
area *= r.coords[i+NUMDIMS] - r.coords[i];
/* area = (r.coords[1] - r.coords[0]) * */
/* (r.coords[3] - r.coords[2]); */
return area;
}
static rect_t RectCombine (rect_t r, rect_t s)
{
int i;
rect_t new_rect;
for (i = 0; i < NUMDIMS; i++)
{
new_rect.coords[i] = MIN(r.coords[i],s.coords[i]);
new_rect.coords[i+NUMDIMS] = MAX(r.coords[i+NUMDIMS],s.coords[i+NUMDIMS]);
}
return new_rect;
}
static int RectOverlap (rect_t r, rect_t s)
{
int i;
for (i = 0; i < NUMDIMS; i++)
if (r.coords[i] > s.coords[i + NUMDIMS] ||
s.coords[i] > r.coords[i + NUMDIMS])
return FALSE;
return TRUE;
}
static rect_t RTreeNodeCover(node_t n)
{
int i;
rect_t r = RectInit();
for (i = 0; i < MAXCARD; i++)
if (n->branch[i].child)
{
r = RectCombine (r, n->branch[i].mbr);
}
else
break;
return r;
}
static void RectPrint (rect_t r)
{
int i;
printf("[");
for (i = 0; i < 2*NUMDIMS; i++)
{
printf("%f",r.coords[i]);
if ( i < 2*NUMDIMS - 1)
printf(",");
}
printf("]");
}