View Full Version : QuickHull and OrientedBBox

Offtopic, but some of you might be interested.

I just finished implementing the algorithm presented in Real-time rendering to calculate a tight fitting bounding box.

It uses quickhull to create the convex hull and then compute the covariance matrix and it's eigenvectors to then find the OBB.

I tested it with all sorts of models and it looks pretty robust. If one of you ever use it and make it crash or throw, just send me the points you used.

You can pick it up here.

http://members.rogers.com/deseric/hullbbox.zip

I believe it is documented well enough(at least the convex hull part!). I know it compiles under gcc3.1, gcc3.0.4, and Borland 5.

MSVC++ users should only need to fix for-loop index variable repetitions.

Well, I hope it will be of use to someone.

I forgot to mention, the point list can contain duplicated points. So that way you just need to dump your vertices in a vector.

[This message has been edited by Gorg (edited 06-02-2002).]

painterb

06-02-2002, 10:31 AM

Thanks Gorg ... that's just the thing I was looking for. http://www.opengl.org/discussion_boards/ubb/smile.gif

jwatte

06-02-2002, 03:12 PM

> MSVC++ users should only need to fix for-

> loop index variable repetitions.

Put this in StdAfx.h (or some other always-included file, or maybe even on the command line):

#define for if(0); else for

This moves the variables declared in the scope of the for clause into the scope of the else clause, so they won't leak out outside the end of the for.

As you know, you CAN turn off language extensions to make the for variable scoping be correct in MSVC. It's a pity the Microsoft headers don't compile correctly when you turn these off, though.

Originally posted by jwatte:

> MSVC++ users should only need to fix for-

> loop index variable repetitions.

Put this in StdAfx.h (or some other always-included file, or maybe even on the command line):

#define for if(0); else for

This moves the variables declared in the scope of the for clause into the scope of the else clause, so they won't leak out outside the end of the for.

As you know, you CAN turn off language extensions to make the for variable scoping be correct in MSVC. It's a pity the Microsoft headers don't compile correctly when you turn these off, though.

Yes it is a pity.

But I am not going to write a fix. I don`t use MSVC so, I prefer to let people take the decision on the solution they want to use. Besides, It`s only in one function anyway.

[This message has been edited by Gorg (edited 06-02-2002).]

SirKnight

06-02-2002, 04:35 PM

Cool, I have been wanting some code to do this. gg! http://www.opengl.org/discussion_boards/ubb/smile.gif

-SirKnight

I just uploaded a new version. I found some problems with highly tessaled models with very uniform distributio n of points(like a sphere of 60000 vertices ).

I added a scale to reduce floating point issue. If you ever get an exception that says "Found only trivial solution. Did you try to pass a plane". Go in ConvexHull.cpp and raise the SCALE variable. For now it seems good.

As side effect, this changes make the algorithm find tighter bounding box because the coefficients of the covariance matrix are bigger.

And I also optimised the initialisation where the redundant points are removed. I got about a 25% speed up.

http://members.rogers.com/hullbbox.h

On a side notes, does anybody knows how to find a good hash value for floats? If I get that, I might be able to speed it up a little more.

Kilam Malik

06-03-2002, 01:01 AM

I use hashing for building up a half-edge datastructure out of triangles. For this I have to find duplicate points fast. My hashing function is:

hVal = (unsigned int)fmod(fabs(100 * p.x() + 10000 * p.y() + 1000000 * p.z()), MESH_BUCKETS);

I use MESH_BUCKTES = 16381.

The problem with floats is only, that for one point you get 1.0000 and for the other you get 0.9999 for a coordinate. As I have a tolerance for which two coordinates are considered the same, I wrote a function to find a point in the hash table which searches the bucket which is returned by the hash function and also the bucket before and after this. But it's still faster... I search three buckets instead all points.

After this I could read a 35 MByte STL Model into a half edge structure in 15 seconds instead of 5 Minutes!

Kilam.

[This message has been edited by Kilam Malik (edited 06-03-2002).]

Thank you Kilam. That's much better than just hashing each float separatly.!!!

I have tried a few formula like that, but I got more chaining than I should and I could not find one the web.

*update*

Great I just implemented it with the hashing of point got another minor speed up.

I guess one day I do some key repartition tests, but for now it is good for me.

The new version is uploaded for anybody interested. This will be last version until a while(unless someone finds a way to break it),

[This message has been edited by Gorg (edited 06-03-2002).]

[This message has been edited by Gorg (edited 06-03-2002).]

WhatEver

06-03-2002, 03:57 PM

Kilam Malik, could you explain the hashing in more depth? I just did some searches and they resulted in only descriptions of what it does and not how to do it. I'm not quite undertanding how the method works.

I need it to remove redundent vertices. A faster method than a modified balloon sort would be great.

WhatEver

06-03-2002, 04:03 PM

Here's my first question: what do you do with hVal?

I'm quite excited to know that there's an algorithm out there that's really fast compaired to my modified balloon sort.

WhatEver

06-03-2002, 04:28 PM

I found this, but I don't understand it: http://www.concentric.net/~Ttwang/tech/inthash.htm

I don't know what to do with the Key, or even how to generate it.

Mr.Freeze

06-03-2002, 05:43 PM

Hi Grog,

I didn't compile your code yet, but I have a question:

I just implemented the same thing a few days ago and it gives me the results I was expecting (tight fitted bounding boxes) except if the shape is a cube (xSize=ySize=zSize)... then the bounding box is randomly oriented and doesn't tightly fit anymore. How did (if you did) you solve that problem?

Thanks,

Marc

(In case I try to get the convex hull of a planar shape, I just add an artificial point outside of the plane and compute the convex hull. Then I remove all triangles containing the added point. This gives me the convex planar shape)

[This message has been edited by Mr.Freeze (edited 06-03-2002).]

WhatEver, just download my code and have look at ConvexHull.cpp.

Look for "PointHashTable".

Mr.Freeze :

Not too be annoying but it's Gorg http://www.opengl.org/discussion_boards/ubb/wink.gif grog sounds grouchy.

Anyway, for the cube. In that case the covariance matrix should be diagonal ==> all the values are 0 aside from those on the diagonal(which are the eigenvalues). If the matrix is diagonal, then the vector are simply ( i, j, k ) the standard basis.

This is true for all object with a symetrical distribution. So spheres and boxes fall in that category.

So before trying to compute the eigenvectors(which you cannot actually because if your matrix is diagonal then you get only the trivial solution(0,0,0)), I simply check if the matrix is diagonal or not.

I stupidly update a version where the hash table doesn't have a prime number as a size. I updated the new version if anybody is interested.

Hopefully it will be the last time I update it for a while., I passed 5 consecutive days doing doing the #%$#!% implementation.

[This message has been edited by Gorg (edited 06-03-2002).]

WhatEver

06-03-2002, 07:05 PM

Thx Gorg! I will check that out asap!

Mr. Freeze, I forgot to mention : If you have a rotated cube, then your covariance matrix will not be diagonal and you will have to solve for the eigenvalues and eigenvectors.

And everythings comes down to how you do that. If you actually test my code with a rotated box, you will see that OBB it finds is slightly off. This is because of the errors caused by finding the eigenvectors. I find the error very small so I don't really care. You can try it and see for yourself I guess.

If you want to know this how I find the eigenvectors.

1. CreateCovariance Matrix

2a. If matrix is diagonal, assign i,j,k to eigenvectors.

and we are done

2b. Find the characteristic polynomial. Solve it. This gives you the eigenvalues. Now solve the 3 following system. If the eigenvalues are e1,e2,e3

(e1 * I - A)v = 0

The nice thing about our matrices, is that they always have 3 eigenvectors.

Now this is a 3x3 system. 2x2 are much more easier to solve. So we do a trick where we set one the component of v to 1. This might change the results of the true eigenvectors by a scale only, so you will still get the correct vectors.

Which component you set to 1. Well any that is not really 0!! a component will really be 0 if all coefficient of this component in the matrix are 0.

So in my code I simply check if x can be 1. So I look for mat[0][0] != 0 | | mat[0][1] != 0 | | mat[0][2] != 0.

If all the coefficients of X are zero, then I check those of Y. If they are all 0, then I check for those of Z. If they are all 0 I quit, because it means you can only find the trivial solution.

But if you found a correct convexhull, this should never happen and you should always find one of the coordinate to set to 1.

[This message has been edited by Gorg (edited 06-03-2002).]

Marc alias Mr. Freeze http://www.opengl.org/discussion_boards/ubb/smile.gif

About you planar shapes. this will not give you 3d hull will it?

I guess you could take the remaining triangles and change their winding and offset them by a bit... hmm.

Anyway, I don't support or plan to support planar meshes anyway, but I will keep it in my mind if I ever need to.

WhatEver : Comment of the Insert function is wrong. It should say "return true if the point is inserted and false if the point was already inserted". But I guess you can figure this out by yourself http://www.opengl.org/discussion_boards/ubb/wink.gif

[This message has been edited by Gorg (edited 06-03-2002).]

Mr.Freeze

06-03-2002, 07:52 PM

Gorg,

Sorry for my wrong spelling of your name... http://www.opengl.org/discussion_boards/ubb/wink.gif

About the rotated cube, well, I guess there must be a subtle mistake in my code and I don't get the correct bounding box. But any other rectangle (where sizeX!=sizeY!=sizeZ!=sizeX) works perfectly fine as for other shapes... ?!? I'll have a look at your code http://www.opengl.org/discussion_boards/ubb/smile.gif

But basically I compute it the exact same way as you.

My convex hull routine doesn't really return the convex hull of a planar shape, just it's corresponding planar convex shape, but as you said, from there the hull is quickly found.

Cheers,

Marc

the best advice I can give you is to output the matrices and eigen values and eigen vectors for each case. (normal cube and rotated cube) and check the difference.

And If you spot anything that can be improved in my code, please let me know.

[This message has been edited by Gorg (edited 06-03-2002).]

Mr.Freeze

06-03-2002, 08:58 PM

Gorg,

it seems very strange to me that for a rotated cube you don't get a diagonal matrix for the covariance matrix.

I checked and check, and I always get a diagonal matrix with a cube... ??

Marc

I get a symetric matrix(like in all case pretty much) but it is not diagonal. Which is normal, Otherwise the only possible eigenvectors would be i,j,k and this is obviously not correct for a rotated cube.

I have a question about your quickhull implementation. What do you do when a point is directly on the plane of the triangle?

What I currenttly do is say that it's outside. And if it ever gives me degenerate triangles, I just don't create then when filling the hole.

Mr.Freeze

06-04-2002, 04:00 PM

Hi Gorg,

I had a look how you compute the covariance matrix and I have a few questions:

1. The covariance matrix shouldn't depend on anything else than vertices if you want to get a tight bounding box. Why do you include triangle information in the calculation? Take a collection of points and find a tight bounding box. Now take the same collection of points and add some random triangles between these points (you don't add more points). The bounding box should be the same as in the first calculation -> it doesn't depend on triangle information.

2. The way I compute the covariance matrix is:

I compute the mean position of all vertices of my convex hull: mu[3]

Then, for each vertice (p[3]) of my convex hull (c[3][3] which is the covariance matrix is initialized to zero):

for (i=0;i<n;i++)

{

float p[3];

p[0]=convHull[3*i+0]; // x component of n-th vertice

p[1]=convHull[3*i+1]; // y component of n-th vertice

p[2]=convHull[3*i+2]; // z component of n-th vertice

for (int j=0;j<3;j++)

for (int k=0;k<3;k++)

c[j][k]=c[j][k]+(p[j]-mu[j])*(p[k]-mu[k]);

}

for (j=0;j<3;j++)

for (int k=0;k<3;k++)

c[j][k]=c[j][k]/n;

3. By computing it my way you can clearly see that the covariance matrix of a rotated cube is always diagonal. Where did I do a mistake?

4. The bounding box one can compute is never the tightest because to compute the covariance in a exact way, one should need to know in advance the center ((max+min)/2) of the bounding box (that means with max and min computed along the axis of the bounding box). Since me don't know it in advance, we compute the mean (mu) which is an approximation of the real center.

What if we compute the bounding box in more than one passes?

First pass: mu=mean position of all vertices of the convex hull

Second pass: mu=center of the obtained bounding box (mu=(max+min)/2)

Third pass: mu=center of the obtained bounding box

...

I think that in that way we could find the tightest bounding box. What do you think?

Cheers,

Marc

(About your question: I implemented a incremental convex hull algorithm and when a point lies on a plane containing a triangle of my conv. hull, I just ignore it (I decide that it is just inside))

I got the covariance matrix from this paper.

http://www.cs.unc.edu/~geom/OBB/paper.html

You can look at the correct equation. But I do not use this equation, because for some reason, I cannot find the eigenvalues of the matrix created by it. So I use the formula from the Book Real-time rendering.

for the equations you use, if the convex hull you find happens to have a huge set of points in a corner, the bbox might get aligned with it and it means it might not be the correct orientation.

the equations I use try to resolve this by "weighting" the covariance with area of the triangles so that the box aligns with the longest triangles.

The explanation for using the triangle is not in the paper and I cannot find where I read it. hmm I am starting to doubt I ever did!! In that case, then I don't why my equation or your equation is better.

[This message has been edited by Gorg (edited 06-04-2002).]

Mr.Freeze

06-04-2002, 08:33 PM

I was inspired by the same paper, but I don't agree with their triangle "trick" (weighting with the triangles). I think it's a "trick" (approximation) because it is impossible to get the main axis in one unique pass and that's why I proposed a second pass (or more) to avoid the problem you mentioned (with a huge set of points in a corner...). In the second pass the "huge set of points" don't weight anymore in the balance since we just take the (min+max)/2 along the axis found in the first pass.

When you say that you want your box to get aligned with the longest triangles... again, I have to say that triangles don't matter at all in the calculation of the bounding box, only vertices matter (if all vertices are inside of the box, then all triangles are inside too).

Hum, I hope you can understand my poor explanation http://www.opengl.org/discussion_boards/ubb/wink.gif

I misread your last post where you proposed multiple passes. It seems it could work well, though I have to sit down and do some calculations and tests.

How well does it work for you for something else than a box?

I have tried to find an explanation for taking the area of the triangle. In the paper, it mention that a better way to calculate the covariance matrix is integrating over the surface of the triangle and then they give the integrated formula.

My guess is that the formula I used is a tweaked formula or simply modified formula of the integral.

But I could not find more details on the web.

Statistics is not my forte, so those things are starting to be out of my limited knowledge.

Mr.Freeze

06-05-2002, 03:41 PM

For something else than a box (with sizeX=sizeY=sizeZ), the results are what one expects: tight bounding boxes (at least they look tight).

I compared bounding boxes generated in one pass and in two passes (as I proposed) and the difference is really negligible for normal shapes (quite equilibrated). For heavily non-equilibrated shapes, the difference is also very small (I got an improvement of 3% with 2 passes). This surprised me a bit.

About your question if a point lies on the plane of a triangle of your convex hull... I think I missunderstood you. What I do is:

If the point lies inside or on that plane, fine.

If the point lies outside of that plane, I remove the triangle associated with that plane and all other triangles of the convex hull which are in the same plane (triangles in the same planes are tagged with a same identifier). By doing so I don't get overlapping triangles in my convex hull (because of precision errors). Then I close the hull with the new point and the edges resulting from the triangles removal.

In your quickhull implementation (I implemented it in an incremental algogithm) it might be a bit more difficult to implement I guess...

[This message has been edited by Mr.Freeze (edited 06-05-2002).]

I believe you will find this interesting :

I tried your method, just the first pass, and I got the exact same boxes. the difference in eigenvectors were usually in the second or third decimal place or they were simply scaled.

So I believe it is actually the same technique, but my version uses an interface where the convex hull comes described in triangles instead of only points(or the triangle contains indexes to the points).

As for the small difference you see when you use 2 or 3 passes, I not sure if it suprising or not, because for most of the models they are actually quite close. But I don't really know how to interpret what you see.

[This message has been edited by Gorg (edited 06-05-2002).]

[This message has been edited by Gorg (edited 06-05-2002).]

Mr.Freeze

06-05-2002, 05:21 PM

Now I understand... :

First I used to compute the bounding box out of a collection of vertices which were not filtered by the convex hull function (I was too lazy to implement it http://www.opengl.org/discussion_boards/ubb/wink.gif ). That's why I didn't use the proposed method in the paper. But up until now (even after implementing the convex hull function) I somehow kept in mind that my bounding box function should work with any kind of vertices... My fault m(__)m. Sure, if you weight your triangles of your convex hull, it's perfectly correct!

Well, this one's over! http://www.opengl.org/discussion_boards/ubb/smile.gif

GoWithTheFlow

06-17-2010, 02:24 AM

Hey guys,

that is exactly what I'm looking for at the moment. But the link on the first site is down... :(

Gorg, maybe you can upload the code somewhere for me? Or can someone else provide the software? That would be great!

Thanks!

Flo

Powered by vBulletin® Version 4.2.3 Copyright © 2017 vBulletin Solutions, Inc. All rights reserved.