c++ - Oriented Bounding Box is Misshapen and the Wrong Size in OpenGL

Question

I'm writing a program in OpenGL to load a mesh and draw an oriented bounding box around said mesh. The mesh loads correctly but when I draw the bounding box the box is the wrong shape and far too small.

The process I used to calculate this box was to use principle component analysis to find a covariance matrix. I then got the eigenvectors of that matrix and used those as a local co-ordinate system to find the 8 vertices of the cube. Then I calculated a transformation to move the cube from the local co-ordinate system to the global co-ordinate system.

The code for calculating the covariance is here:

std::array<std::array<double, 3>, 3> covarianceCalc2()
{
 std::array<std::array<double, 3>, 3> sum = {{{0, 0, 0}, {0, 0, 0}, {0, 0, 0,}}};
 std::array<double, 3> tempVec;
 double mean = 0;
 for(int i = 0; i < meshVertices.size(); i++)
 {
     mean += meshVertices[i].x;
     mean += meshVertices[i].y;
     mean += meshVertices[i].z;
 }
 mean = mean/(meshVertices.size() * 3);

 for(int i = 0; i < meshVertices.size(); i++)
 {
     //mean = (meshVertices[i].x + meshVertices[i].y + meshVertices[i].z)/3;
     tempVec[0] = meshVertices[i].x - mean;
     tempVec[1] = meshVertices[i].y - mean;
     tempVec[2] = meshVertices[i].z - mean;
     sum = matrixAdd(sum, vectorTranposeMult(tempVec));
 }
 sum = matrixMultNum(sum,(double) 1/(meshVertices.size()));
 return sum;
 }

The code for calculating the eigenvectors is here:

void Compute_EigenV(std::array<std::array<double, 3>, 3> covariance, double eigenValues[3], double eigenVectors_1[3], double eigenVectors_2[3], double eigenVectors_3[3])
{

    printf("Matrix Stuff
");
    MatrixXd m(3, 3);
    m << covariance[0][0], covariance[0][1], covariance[0][2],
         covariance[1][0], covariance[1][1], covariance[1][2],
         covariance[2][0], covariance[2][1], covariance[2][2];

    // volving SVD
    printf("EigenSolver
");
    EigenSolver<MatrixXd> solver(m);
    MatrixXd all_eigenVectors = solver.eigenvectors().real();
    MatrixXd all_eigenValues = solver.eigenvalues().real();

    // find the max index
    printf("Find Max Index
");
    int INDEX[3];
    double max;
    max=all_eigenValues(0,0);
    int index=0;
    for (int i=1;i<3;i++){
        if (max<all_eigenValues(i,0)){
            max=all_eigenValues(i,0);
            index=i;
        }
    }
    INDEX[0]=index;

    // find the min index
    printf("Find Min Index
");
    double min;
    min=all_eigenValues(0,0);

    index=0;
    for (int i=1;i<3;i++){
        if (min>all_eigenValues(i,0)){
            min=all_eigenValues(i,0);
            index=i;
        }
    }
    INDEX[1]=3-index-INDEX[0];
    INDEX[2]=index;

    // giave eigenvalues and eien vectors to matrix
    printf("Give values and vector to matrix
");
    eigenValues[0]=all_eigenValues(INDEX[0],0);
    printf("1");
    eigenValues[1]=all_eigenValues(INDEX[1],0);
    printf("1
");
    eigenValues[2]=all_eigenValues(INDEX[2],0);

    printf("Vector 1
");
    VectorXd featureVector_1 = all_eigenVectors.col(INDEX[0]);
    eigenVectors_1[0]=featureVector_1(0);
    eigenVectors_1[1]=featureVector_1(1);
    eigenVectors_1[2]=featureVector_1(2);

    printf("Vector 2
");
    VectorXd featureVector_2 = all_eigenVectors.col(INDEX[1]);
    eigenVectors_2[0]=featureVector_2(0);
    eigenVectors_2[1]=featureVector_2(1);
    eigenVectors_2[2]=featureVector_2(2);

    printf("Vector 3
");
    VectorXd featureVector_3 = all_eigenVectors.col(INDEX[2]);
    eigenVectors_3[0]=featureVector_3(0);
    eigenVectors_3[1]=featureVector_3(1);
    eigenVectors_3[2]=featureVector_3(2);

}

The code that finds the global co-ordinates is this:

std::array<double, 3> localToGlobal(std::array<double, 3> vec, double eigenVector1[3], double eigenVector2[3], double eigenVector3[3], double mean)
{
std::array<double, 3> tempVec;
std::array<std::array<double, 3>, 3> eigenArray;
eigenArray[0][0] = eigenVector1[0]; eigenArray[0][1] = eigenVector2[0]; eigenArray[0][2] = eigenVector3[0];
eigenArray[1][0] = eigenVector1[1]; eigenArray[1][1] = eigenVector2[1]; eigenArray[1][2] = eigenVector3[1];
eigenArray[2][0] = eigenVector1[2]; eigenArray[2][1] = eigenVector2[2]; eigenArray[2][2] = eigenVector3[2];

tempVec = matrixVectorMult(eigenArray, vec);
tempVec[0] += mean;
tempVec[1] += mean;
tempVec[2] += mean;

return tempVec;
}

The code that calls all of these and draws the box is:

void obbBoundingBox()
{

double eigenValues[3] = {0, 0, 0};
double eigenVectors_1[3] = {0, 0, 0}, eigenVectors_2[3] = {0, 0, 0}, eigenVectors_3[3] = {0, 0, 0};

Compute_EigenV(covarianceCalc2(), eigenValues, eigenVectors_1, eigenVectors_2, eigenVectors_3);


std::array<double, 3> point1 = {findVectorMax(eigenVectors_1), findVectorMax(eigenVectors_2), findVectorMax(eigenVectors_3)};
std::array<double, 3> point2 = {findVectorMax(eigenVectors_1), findVectorMax(eigenVectors_2), findVectorMin(eigenVectors_3)};
std::array<double, 3> point3 = {findVectorMax(eigenVectors_1), findVectorMin(eigenVectors_2), findVectorMin(eigenVectors_3)};
std::array<double, 3> point4 = {findVectorMax(eigenVectors_1), findVectorMin(eigenVectors_2), findVectorMin(eigenVectors_3)};
std::array<double, 3> point5 = {findVectorMin(eigenVectors_1), findVectorMax(eigenVectors_2), findVectorMax(eigenVectors_3)};
std::array<double, 3> point6 = {findVectorMin(eigenVectors_1), findVectorMax(eigenVectors_2), findVectorMin(eigenVectors_3)};
std::array<double, 3> point7 = {findVectorMin(eigenVectors_1), findVectorMin(eigenVectors_2), findVectorMax(eigenVectors_3)};
std::array<double, 3> point8 = {findVectorMin(eigenVectors_1), findVectorMin(eigenVectors_2), findVectorMin(eigenVectors_3)};


 double mean = 0;
 for(int i = 0; i < meshVertices.size(); i++)
 {
     mean += meshVertices[i].x;
     mean += meshVertices[i].y;
     mean += meshVertices[i].z;
 }
 mean = mean/(meshVertices.size() * 3);

point1 = localToGlobal(point1, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point2 = localToGlobal(point2, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point3 = localToGlobal(point3, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point4 = localToGlobal(point4, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point5 = localToGlobal(point5, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point6 = localToGlobal(point6, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point7 = localToGlobal(point7, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);
point8 = localToGlobal(point8, eigenVectors_1, eigenVectors_2, eigenVectors_3, mean);


glEnable(GL_BLEND);
 glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
 glColor4f(1.0f, 1.0f, 1.0f, 0.5f);

glBegin(GL_QUADS);
//Front Face
    glVertex3f(point1[0], point1[1], point1[2]);
    glVertex3f(point3[0], point3[1], point3[2]);
    glVertex3f(point7[0], point7[1], point7[2]);
    glVertex3f(point5[0], point5[1], point5[2]);
glEnd();

glBegin(GL_QUADS);
//Left Face
    glVertex3f(point5[0], point5[1], point5[2]);
    glVertex3f(point7[0], point7[1], point7[2]);
    glVertex3f(point8[0], point8[1], point8[2]);
    glVertex3f(point6[0], point6[1], point6[2]);
glEnd();

glBegin(GL_QUADS);
//Back Face
    glVertex3f(point6[0], point6[1], point6[2]);
    glVertex3f(point8[0], point8[1], point8[2]);
    glVertex3f(point4[0], point4[1], point4[2]);
    glVertex3f(point2[0], point2[1], point2[2]);
glEnd();

glBegin(GL_QUADS);
//Right Face
    glVertex3f(point2[0], point2[1], point2[2]);
    glVertex3f(point4[0], point4[1], point4[2]);
    glVertex3f(point3[0], point3[1], point3[2]);
    glVertex3f(point1[0], point1[1], point1[2]);
glEnd();

glBegin(GL_QUADS);
//Top Face
    glVertex3f(point2[0], point2[1], point2[2]);
    glVertex3f(point1[0], point1[1], point1[2]);
    glVertex3f(point5[0], point5[1], point5[2]);
    glVertex3f(point6[0], point6[1], point6[2]);
glEnd();

glBegin(GL_QUADS);
//Bottom Face
    glVertex3f(point4[0], point4[1], point4[2]);
    glVertex3f(point3[0], point3[1], point3[2]);
    glVertex3f(point7[0], point7[1], point7[2]);
    glVertex3f(point8[0], point8[1], point8[2]);
glEnd();




 }

Answer 1

The code looks right (I haven't checked all the functions like matrixmult and like that but I assume you have tested and checked them).

Problem: you have a small misunderstanding of what are you doing there.

So, to help you a bit, but not code your coursework by myself, as that would get us both into trouble, I decided to make a small tutorial in Matlab (as I know you have access to) of what you want to do and why. Some functions are missing, but you should be able to understand what's going on:

clear;clc;
%% your "bunny"
vtx=    [ 1 0
         1 1
         2 2
         3 3 
         1 3];

% Lets draw it  

hold on
plot(vtx(:,1),vtx(:,2),'.')

% lets draw XY axis also
plot([0 4],[0 0],'k')
plot([0 0],[0 4],'k')
axis([-1 5 -1 5])
axis square

%% Draw abb
maxX=max(vtx(:,1));
maxY=max(vtx(:,2));
minX=min(vtx(:,1));
minY=min(vtx(:,2));

% Mising: Create a square and draw it 
cub=createcube(maxX,maxY,minX,minY)
drawabb(cub);

%% Create obb

C=cov(vtx);
vtxmean=mean(vtx);

[eVect,~]=eig(C);

% Draw new local coord system
plot([vtxmean(1) vtxmean(1)+eVect(1,1)],[vtxmean(2) vtxmean(2)+eVect(1,2)],'k')
plot([vtxmean(1) vtxmean(1)+eVect(2,1)],[vtxmean(2) vtxmean(2)+eVect(2,2)],'k')

Now you can see that if we get the max and min of the eigenvector, it doesnt make too much sense. Thats NOT the obb. So what are we suposed to do? Well, we can create an abb, but alligned to the NEW axis, not the XY axis!

What would we need for that? well we need to know the values of our points in the new coordinate axis, dont we?

Localvtx=fromGlobalToLocal(vtx,eVect,vtxmean);

% get the max and min of the points IN THE NEW COORD SYSTEM!
maxX=max(Localvtx(:,1));
maxY=max(Localvtx(:,2));
minX=min(Localvtx(:,1));
minY=min(Localvtx(:,2));

Fantastic!!! Now , we can create a square in this coord system, and using fromLocalToGlobal, draw it in XY!!

obbcub=createcube(maxX,maxY,minX,minY);
obbcubXY=fromLocalToGlobal(obbcub,eVect,vtxmean);
drawcube(obbcubXY);

Logical question: WHY ARE WE DOING ALL THIS!?!?

Well its an interesting question indeed. Do you play videogames? Have you ever took "an arrow in the knee?". How does a computer know if you have shooted the guy in the head or in the leg with you sniper rifle if the guy was jumping and crouching and lying all the time!!

What about an oriented bounding box! If you know the bounding box of the leg, or the head, independently of the geometrical position of the model, you can compute if the shot went inside that box or not! (dont take everything literally, this is a huge world and there are infinite ways of doing things like this).

See example:

Sidenote: Dont use my words or image or code in you report, as that will be considered cheating! (just in case)