Quadtrees

[std10093]

I think it's time to get a look at your code.

I had implemented GRAHAM-SCAN algorithm for finding convex hulls, that runs in O(nlogn).

[std10093]

Your test point seems to be outside the polytope bounding box, so it should not even descend into the tree.

Recall that the index of the child for descending into the tree is:

child = 1 * (testpoint.x >= split.x)
      + 2 * (testpoint.y >= split.y)
      + 4 * (testpoint.z >= split.z);

so if we compute the index, it is correct to be five, since the Test point is 1.1, 0.5, 0.7

But, we should probably first check if the test point is inside the bounding box.

This seems to do the trick:

/**
 * Checks if a point lies inside the bounding box of the polytope.
 * @param box - the bounding box of the polytope (of the root)
 * @param p - test point
 * @return - result
 */
bool Quadtree::insidePolytopeBox(const BoundingBox* box, const Point3d& p)
{
    for(int i = 0 ; i < D ; ++i) // Check all coordinates
        if(box->getMinBox()[i] > p[i] || box->getMaxBox()[i] < p[i])
            return false;;
    return true;        
}

[Nominal Animal]

Recall that the index of the child for descending into the tree is:

child = 1 * (testpoint.x >= split.x)
      + 2 * (testpoint.y >= split.y)
      + 4 * (testpoint.z >= split.z);

so if we compute the index, it is correct to be five, since the Test point is 1.1, 0.5, 0.7

No; you are forgetting that the halfspaces retained for inclusion tests in each node (or its child nodes) are only relevant to points within the bounding box of that node.

Remember, we specifically store only the information relevant to that bounding box. So, descending into the tree with a point outside the original bounding box yields undefined results.

In practice, it should not matter, as every direction from inside the halfspace should be bounded by a halfspace.

Our implementations don't do that, though, because of the way we define "interesting halfspace": we ignore halfspaces that are perpendicular to the X, Y, or Z axes, on the lower-coordinate side of the polytope (the inside==0, outside==0 case; i.e. one or more corners of the box are on the halfspace plane, but no corner is clearly inside or outside the halfspace). This does not mean there should be any possibility of infinite recursion or stuff like that; it just means the tree constructed is only valid for points inside the original bounding box.

Take the (0,0,0),(1,0,0),(0,1,0),(0,0,1) tetrahedron as an example. If you test point (-1,-1,-1), ignoring the original bounding box, you get result inside. That is because the three axis-perpendicular faces are coplanar with the bounding box, and found "uninteresting" -- as they don't do any further limiting to point inclusion compared to the bounding box alone. The tree is really just one leaf with the diagonal halfspace, as the bounding box coincides with the other sides. So, negative coordinates will always yield "inside", because they are "inside" if you only consider the diagonal halfspace.

If you were to consider the halfspaces that are coplanar with a bounding box corner, edge, or side interesting, then the initial bounding box does not matter, and it would be correct to descend into child five in the case you showed.

[std10093]

Ok, but having like the one I posted above ( insidePolytopeBox ), has a really low cost, so I think I am gonna use it.

I got rid of the BoundingBox and I am calculating the boxes locally.. The problem is that my tree gets very very big and that is because I find many interesting halfspaces, until I have created many levels...

[Nominal Animal]

I think it's time to get a look at your code.

The any-dimensional code is way too spaghetti for me to post it; I cringe every time I look at it.

Besides, the algorithm does not truly scale to "any" number of dimensions D, because each inner node will need to contain 2^D child pointers. On 32-bit systems, a single 28-dimensional inner node takes a gigabyte!

The relationship between the number of dimensions D, and leaf nodes having up to t halfspaces, is such that on 64-bit systems where pointers and doubles are the same size, inner and leaf nodes are about the same size if

t ≃ 2^D / (D + 1)

For 20-dimensional data, an inner node is the same size as a leaf node with about 50,000 halfspaces.

In other words, I strongly believe this algorithm is suitable for small number of dimensions only.

Also, I haven't had time to implement a nice convex hull algorithm, either; partly because the methods get pretty hairy (to fully grok) when there are more than three spatial dimensions. (I personally would need to take the time to develop the two and three dimensional cases in parallel first, then compare the two to develop the generic case.)

I think I'll just clean up my 3D convex hull case, post it here, and call it done.