UT/UT_RTree.C

Go to the documentation of this file.

/*
 * PROPRIETARY INFORMATION.  This software is proprietary to
 * Side Effects Software Inc., and is not to be reproduced,
 * transmitted, or disclosed in any way without written permission.
 *
 * Produced by:
 *      Michiel Hagedoorn
 *      Side Effects Software Inc
 *      123 Front Street West, Suite 1401
 *      Toronto, Ontario
 *      Canada   M5J 2M2
 *      416-504-9876
 *
 * NAME:        UT_RTreeGeneric.C (CLO library, C++)
 *
 * COMMENTS:
 */

#define NOMINMAX

#include "UT_RTree.h"
#include "UT_Assert.h"
#include <ostream>
#include <algorithm>

#undef min
#undef max

template< typename T >
inline
UT_BoxT<T>::UT_BoxT()
{
    myMin[0] = std::numeric_limits<T>::max();
    myMin[1] = std::numeric_limits<T>::max();
    myMin[2] = std::numeric_limits<T>::max();

    myMax[0] = -std::numeric_limits<T>::max();
    myMax[1] = -std::numeric_limits<T>::max();
    myMax[2] = -std::numeric_limits<T>::max();
}

template< typename T >
inline
UT_BoxT<T>::UT_BoxT(const T p[3])
{
    myMin[0] = myMax[0] = p[0];
    myMin[1] = myMax[1] = p[1];
    myMin[2] = myMax[2] = p[2];
}

template< typename T >
inline bool
UT_BoxT<T>::isEmpty() const
{
    return(
        (myMin[0] > myMax[0]) ||
        (myMin[1] > myMax[1]) ||
        (myMin[2] > myMax[2])
    );
}

template< typename T >
inline T
UT_BoxT<T>::getMin(const int c) const
{
    UT_ASSERT( (0 <= c) && (c < 3) );
    UT_ASSERT( !isEmpty() );

    return myMin[c];
}

template< typename T >
inline T
UT_BoxT<T>::getMax(const int c) const
{
    UT_ASSERT( (0 <= c) && (c < 3) );
    UT_ASSERT( !isEmpty() );

    return myMax[c];
}

template< typename T >
inline void
UT_BoxT<T>::assignPoint(const T*const p)
{
    myMin[0] = p[0];
    myMin[1] = p[1];
    myMin[2] = p[2];

    myMax[0] = p[0];
    myMax[1] = p[1];
    myMax[2] = p[2];
}

template< typename T >
inline void
UT_BoxT<T>::absorbPoint(const T*const p)
{
    myMin[0] = (myMin[0] <= p[0]) ? myMin[0] : p[0];
    myMin[1] = (myMin[1] <= p[1]) ? myMin[1] : p[1];
    myMin[2] = (myMin[2] <= p[2]) ? myMin[2] : p[2];

    myMax[0] = (myMax[0] >= p[0]) ? myMax[0] : p[0];
    myMax[1] = (myMax[1] >= p[1]) ? myMax[1] : p[1];
    myMax[2] = (myMax[2] >= p[2]) ? myMax[2] : p[2];
}

template< typename T >
inline void
UT_BoxT<T>::absorbBox(const UT_BoxT<T>& b)
{
    myMin[0] = (myMin[0] <= b.myMin[0]) ? myMin[0] : b.myMin[0];
    myMin[1] = (myMin[1] <= b.myMin[1]) ? myMin[1] : b.myMin[1];
    myMin[2] = (myMin[2] <= b.myMin[2]) ? myMin[2] : b.myMin[2];

    myMax[0] = (myMax[0] >= b.myMax[0]) ? myMax[0] : b.myMax[0];
    myMax[1] = (myMax[1] >= b.myMax[1]) ? myMax[1] : b.myMax[1];
    myMax[2] = (myMax[2] >= b.myMax[2]) ? myMax[2] : b.myMax[2];
}

template< typename T >
inline void
UT_BoxT<T>::expandDistance(const T l)
{
    if( isEmpty() )
        return;

    myMin[0] -= l;
    myMin[1] -= l;
    myMin[2] -= l;

    myMax[0] += l;
    myMax[1] += l;
    myMax[2] += l;
}

template< typename T >
inline bool
UT_BoxT<T>::contains(const T*const p) const
{
    return ( (myMin[0] <= p[0]) && (p[0] <= myMax[0]) &&
             (myMin[1] <= p[1]) && (p[1] <= myMax[1]) &&
             (myMin[2] <= p[2]) && (p[2] <= myMax[2]) );
}

template< typename T >
inline bool
UT_BoxT<T>::intersects(const UT_BoxT<T>& b) const
{
    return (
        (UTmax(myMin[0], b.myMin[0]) <= UTmin(myMax[0], b.myMax[0])) &&
        (UTmax(myMin[1], b.myMin[1]) <= UTmin(myMax[1], b.myMax[1])) &&
        (UTmax(myMin[2], b.myMin[2]) <= UTmin(myMax[2], b.myMax[2]))
    );
}

template< typename T >
inline int
UT_BoxT<T>::getLargestAxis() const
{
    if( isEmpty() )
        return 0;

    int max_axis(0);
    T max_length(myMax[0] - myMin[0]);

    T length(0);

    length = myMax[1] - myMin[1];
    if( length > max_length )
    {
        max_length = length;
        max_axis = 1;
    }

    length = myMax[2] - myMin[2];
    if( length > max_length )
    {
        max_length = length;
        max_axis = 2;
    }
        
    return max_axis;
}

template< typename T >
inline ostream& operator<<(ostream& os, const UT_BoxT<T>& b)
{
    cout << "Box{ min = ";
    cout << b.myMin[0] << ", " << b.myMin[1] << ", " << b.myMin[2];
    cout << ", max = ";
    cout << b.myMax[0] << ", " << b.myMax[1] << ", " << b.myMax[2];
    cout << " }";
    return os;
}

template< typename T>
struct UT_BoxedItemT
{
    typedef UT_BoxT<T> Box;

    Box myBox;
    int myItem;

    UT_BoxedItemT() :
        myBox(),
        myItem(-1)
    {
    }
};

template< typename T >
inline ostream& operator<<(ostream& os, const UT_BoxedItemT<T>& bf)
{
    cout << "BoxedItem{ ";
    cout << "box = ";
    cout << bf.myBox << ", ";
    cout << "item = ";
    cout << bf.myItem;
    cout << " }";
    return os;
}

// Compute the bounding box for a range of boxes [begin, end)
template< typename T>
inline void
computeBoundingBox(
    UT_BoxT<T>& bounding_box,
    const UT_BoxT<T>*const begin, const UT_BoxT<T>*const end
)
{
    UT_ASSERT( end != begin );

    bounding_box = *begin;

    for( const UT_BoxT<T>* i = begin + 1; i != end; ++i )
    {
        bounding_box.absorbBox(*i);
    }
}

// Compute the bounding box for all the boxes in 
// a range of boxed items [begin, end)
template< typename T >
inline void
computeBoundingBox(
    UT_BoxT<T>& bounding_box,
    const UT_BoxedItemT<T>*const begin, const UT_BoxedItemT<T>*const end
)
{
    UT_ASSERT( end != begin );

    bounding_box = begin->myBox;

    for( const UT_BoxedItemT<T>* i = begin + 1; i != end; ++i )
    {
        bounding_box.absorbBox(i->myBox);
    }
}

// Sort boxed items in the direction of an axis, by their centers
template< typename T >
struct UT_BoxedItemAxisOrderT
{
    int myAxis;

    explicit UT_BoxedItemAxisOrderT(const int axis) :
        myAxis(axis)
    {
    }

    ~UT_BoxedItemAxisOrderT()
    {
    }

    bool operator()(const UT_BoxedItemT<T>& a, const UT_BoxedItemT<T>& b)
    {
        return ( a.myBox.getMin(myAxis) + a.myBox.getMax(myAxis) <
                 b.myBox.getMin(myAxis) + b.myBox.getMax(myAxis) );
    }
};

template< int ORDER >
class UT_RNode
{
public:
    typedef uint32 Index;

    UT_RNode()
    {
        clear();
    }

    ~UT_RNode()
    {
    }

    void clear()
    {
        setInternal(false);
        for( int s = 0; s < ORDER; ++s )
        {
            mySlots[s] = SLOT_VALUE_EMPTY;          
        }
    }

    void setInternal(const bool isInternal)
    {
        mySlots[0] = ( (mySlots[0] & SLOT_MASK_INDEX) | 
                       (isInternal ? SLOT_FLAG_INTERNAL : 0) );
    }

    // Assuming this is an internal node, assign a node index to slot s
    void assignSubtree(const int s, const int index_root_subtree)
    {
        UT_ASSERT( isInternal() );
        UT_ASSERT( (0 <= s) && (s < ORDER) );
        UT_ASSERT( 0 <= index_root_subtree );
        UT_ASSERT(
            (static_cast<uint32>(index_root_subtree) & SLOT_MASK_INDEX) ==
            static_cast<uint32>(index_root_subtree) 
        );

        mySlots[s] = index_root_subtree | SLOT_FLAG_INTERNAL;  
    }

    // Assign an item to slot s
    void assignItem(const int s, const int item)
    {
        UT_ASSERT( !isInternal() );
        UT_ASSERT( (0 <= s) && (s < ORDER) );
        UT_ASSERT( 0 <= item );
        UT_ASSERT(
            (static_cast<uint32>(item) & SLOT_MASK_INDEX) ==
            static_cast<uint32>(item) 
        );

        mySlots[s] = item;
    }

    bool isInternal() const
    {
        return ((mySlots[0] & SLOT_FLAG_INTERNAL) != 0);
    }

    // Get node stored in slot.
    // This assumes that the last assignment to this slot was an node!
    // (not an item)
    int getSubtree(const int s) const
    {
        UT_ASSERT( isInternal() );
        UT_ASSERT( (0 <= s) && (s < ORDER) );
        UT_ASSERT( (mySlots[s] & SLOT_MASK_INDEX) != SLOT_VALUE_EMPTY );

        return (mySlots[s] & SLOT_MASK_INDEX);
    }

    int getNumItems() const
    {
        UT_ASSERT( !isInternal() );

        int num_items(0);
        for( int s = 0; s < ORDER; ++s )
        {
            if( (mySlots[s] & SLOT_MASK_INDEX) == SLOT_VALUE_EMPTY )
                break;

            ++num_items;
        }

        return num_items;
    }

    // Get item stored in slot.
    // This assumes that the last assignment to this slot was an item!
    // (not a node)
    int getItem(const int s) const
    {
        UT_ASSERT( !isInternal() );
        UT_ASSERT( (0 <= s) && (s < ORDER) );
        UT_ASSERT( (mySlots[s] & SLOT_MASK_INDEX) != SLOT_VALUE_EMPTY );

        return (mySlots[s] & SLOT_MASK_INDEX);
    }

private:
    typedef Index Slot;
    static const Slot SLOT_FLAG_INTERNAL = 0x80000000;
    static const Slot SLOT_MASK_INDEX    = 0x7fffffff;
    static const Slot SLOT_VALUE_EMPTY   = SLOT_MASK_INDEX;

    // mySlots[0] has a bit that defines whether the node is an internal
    // node. This is the SLOT_FLAG_INTERNAL bit.
    // If the node is internal, then all the slots represent indices
    // in the shared node array (each slot points to a subtree root).
    // An internal node is always full: each slot contains a valid node index.
    // If the node is a leaf (not internal), then the slots represent
    // items. A leaf may not always be full. When it's not full, then
    // the empty slots form a contiguous range after all the nonempty slots.
    // The full range consists of the slots s in the half-open range
    // [ 0, getNumItems() )
    Slot mySlots[ORDER];
};

// Create a subtree for the range of boxed items [begin, end).
// Return the index of the substree root in shared_nodes
// In addition, bounding_box will be a bounding box for the subtree.
template< int ORDER, typename T >
inline UT_RNode<ORDER>*
subtreeCreate(
    UT_BoxT<T>& bounding_box,
    UT_BoxT<T> shared_boxes[],
    UT_BoxedItemT<T>*const begin, UT_BoxedItemT<T>*const end,
    UT_RNode<ORDER>*const shared_nodes,
    UT_RNode<ORDER>*& shared_nodes_end
)
{
    UT_ASSERT_COMPILETIME( ORDER >= 2 );

    UT_RNode<ORDER>* node(shared_nodes_end++);
    node->clear();

    // node_boxes is the sub-array that contains the ORDER boxes
    // for the node that we're building
    UT_BoxT<T>*const node_boxes(
        shared_boxes +
        ((node - shared_nodes) * ORDER)
    );

    const int size(end - begin);

    UT_ASSERT( size > 0 );

    computeBoundingBox(bounding_box, begin, end);
    
    if( size > ORDER )
    {
        // Create an internal node with ORDER nonempty slots.
        // slot s will have the index of the root of a subtree,
        // which is contained in node_boxes[s]
        node->setInternal(true);
        
        // TODO: replace sort by O(n) way of partitioning
        // TODO: perhaps splitting along a single axis is not the best
        UT_BoxedItemAxisOrderT<T> order(bounding_box.getLargestAxis());
        std::sort(begin, end, order);

        const int n(end - begin);
        const int m(n / ORDER);

        UT_ASSERT( m >= 1 );

        UT_BoxedItemT<T>* begin_sub(begin);

        for( int s = 0; s < ORDER - 1; ++s )
        {
            UT_BoxedItemT<T>*const end_sub(begin_sub + m);

            node->assignSubtree(
                s,
                subtreeCreate(
                    node_boxes[s], 
                    shared_boxes,
                    begin_sub, end_sub,
                    shared_nodes,
                    shared_nodes_end
                ) - shared_nodes //FIXME: ugly pointer arithmetic
            );

            begin_sub = end_sub;
        }

        UT_ASSERT( begin_sub < end );

        node->assignSubtree(
            ORDER - 1,
            subtreeCreate(
                node_boxes[ORDER - 1],
                shared_boxes,
                begin_sub, end,
                shared_nodes,
                shared_nodes_end
            ) - shared_nodes // FIXME: ugly pointer arithmetic
        );
    }
    else
    {
        // Create a leaf with size nonempty slots.
        // Each slot has the index of a shared node box.
        node->setInternal(false);

        int s(0);
        for( const UT_BoxedItemT<T>* i = begin; i != end; ++i )
        {
            node_boxes[s] = i->myBox;
            node->assignItem(s, i->myItem);
            ++s;
        }
    }
 
    return node;
}

template< int ORDER >
inline int
UT_RTreeGeneric<ORDER>::getNumItems() const
{
    return myNumItems;
}

template< int ORDER, typename T >
inline void
subtreeCreateBoxAssignment(
    UT_BoxT<T>& bounding_box,
    UT_BoxT<T>*const shared_boxes,
    const UT_RNode<ORDER>*const shared_nodes,
    const UT_RNode<ORDER>*const node, 
    const UT_BoxT<T> item_boxes[]
)
{
    UT_ASSERT( node != 0 );

    // node_boxes is the sub-array that contains the ORDER boxes
    // for the node that we're visiting
    UT_BoxT<T>*const node_boxes(
        shared_boxes +
        ((node - shared_nodes) * ORDER)
    );

    if( node->isInternal() )
    {
        for( int s = 0; s < ORDER; ++s )
        {       
            subtreeCreateBoxAssignment(
                node_boxes[s],
                shared_boxes,
                shared_nodes,
                shared_nodes + node->getSubtree(s),
                item_boxes
            );
        }
    }
    else
    {
        const int num_items(node->getNumItems());
        for( int s = 0; s < num_items; ++s )
        {
            node_boxes[s] = item_boxes[node->getItem(s)];
        }
        for( int s = num_items; s < ORDER; ++s )
        {
            node_boxes[s] = UT_BoxT<T>(); // Empty set
        }
    }

    computeBoundingBox(bounding_box, node_boxes, node_boxes + ORDER);
}

template< int ORDER, typename T >
inline int*
subtreeGetIntersectingItems(
    const UT_RNode<ORDER>*const shared_nodes,
    const UT_BoxT<T>*const shared_boxes,
    const UT_RNode<ORDER>*const node,
    const UT_BoxT<T>& query_box,
    int*const items_begin
)
{
    UT_ASSERT( node != 0 );

    // node_boxes is the sub-array that contains the ORDER boxes
    // for the node that we're visiting
    const UT_BoxT<T>*const node_boxes(
        shared_boxes +
        ((node - shared_nodes) * ORDER)
    );

    int* i(items_begin);

    if( node->isInternal() )
    {
        for( int s = 0; s < ORDER; ++s )
        {
            if( node_boxes[s].intersects(query_box) )
            {
                i = subtreeGetIntersectingItems(
                    shared_nodes, shared_boxes,
                    shared_nodes + node->getSubtree(s),
                    query_box, i
                );
            }
        }
    }
    else
    {
        for( int s = 0; s < ORDER; ++s )
        {
            if( node_boxes[s].intersects(query_box) )
            {
                *i++ = node->getItem(s);
            }
        }
    }

    return i;
}

template< int ORDER >
template< typename T >
inline UT_RTreeGeneric<ORDER>::UT_RTreeGeneric(
    const UT_BoxT<T> item_boxes[], const int size
) :
    mySharedNodesCapacity(0),
    mySharedNodes(0),
    myRoot(0),
    myNumItems(size)
{
    UT_ASSERT_COMPILETIME( ORDER >= 2 );

    UT_ASSERT( size >= 0 );

    if( size <= 0 )
        return;

    mySharedNodesCapacity = 2 * size; // enough for any tree!
    mySharedNodes = new UT_RNode<ORDER>[mySharedNodesCapacity];

    UT_RNode<ORDER>* shared_nodes_end(mySharedNodes);

    // Put item #i in box[i], resulting in boxes_items
    std::vector< UT_BoxedItemT<T> > boxed_items;
    boxed_items.resize(size);
    for( int i = 0; i < size; ++i )
    {
        boxed_items[i].myBox = item_boxes[i];
        boxed_items[i].myItem = i;
    }

    // Create the R-tree nodes, based on boxed_items
    UT_BoxT<T> bounding_box_root;
    std::vector< UT_BoxT<T> > shared_boxes(
        ORDER * mySharedNodesCapacity
    );
    myRoot = subtreeCreate(
        bounding_box_root,
        &(shared_boxes[0]),
        &(boxed_items[0]), &(boxed_items[0]) + boxed_items.size(),
        mySharedNodes,
        shared_nodes_end
    );

    UT_ASSERT( shared_nodes_end - mySharedNodes < mySharedNodesCapacity );
}

template< int ORDER >
inline UT_RTreeGeneric<ORDER>::~UT_RTreeGeneric()
{
    delete[] mySharedNodes;
    mySharedNodes = 0;
}

template< int ORDER >
template< typename T>
inline void
UT_RTreeGeneric<ORDER>::createBoxAssignment(
    UT_RTreeBoxAssignmentT<T>& assignment, 
    const UT_BoxT<T> item_boxes[]
) const
{
    if( myRoot == 0 )
        return;

    // Allocate room for ORDER boxes for each shared node
    // Box c of node n is stored in slot (ORDER * n) + c
    assignment.myBoxesForSharedNodes.resize(
        ORDER * mySharedNodesCapacity
    );

    UT_BoxT<T> bounding_box_root;
    subtreeCreateBoxAssignment(
        bounding_box_root,
        &(assignment.myBoxesForSharedNodes[0]),
        mySharedNodes,
        myRoot, 
        item_boxes
    );
}

template< int ORDER >
template< typename T>
inline int* 
UT_RTreeGeneric<ORDER>::getIntersectingItems(
    const UT_BoxT<T>& query_box, 
    const UT_RTreeBoxAssignmentT<T>& assignment, 
    int*const items_begin
) const
{
    if( myRoot == 0 )
        return items_begin;

#if 0
    // Return all items, for debugging
    int* items_end(items_begin);
    for( int i = 0; i < myNumItems; ++i )
    {
        *items_end++ = i;
    }
    return items_end;
#else
    return subtreeGetIntersectingItems(
        mySharedNodes, 
        &(assignment.myBoxesForSharedNodes[0]),
        myRoot, query_box, items_begin
    );
#endif // 0
}

---