3

I have a set of points (3D) taken from a range scanner. Sample data can be found here: http://pastebin.com/RBfQLm56

I also have the following parameters for the scanner:

camera matrix
[3871.88184, 0, 950.736938;
  0, 3871.88184, 976.1383059999999;
  0, 0, 1]



distortion coeffs
[0.020208003; -1.41251862; -0.00355229038; -0.00438868301; 6.55825615]



camera to reference point (transform)

[0.0225656671, 0.0194614234, 0.9995559233, 1.2656986283;

  -0.9994773883, -0.0227084301, 0.0230060289, 0.5798922567;

  0.0231460759, -0.99955269, 0.0189388219, -0.2110195758;

  0, 0, 0, 1]

I am trying to render these points properly using opengl but the rendering does not look right. What is the correct way to set openGL projection and modelview matrix? This is what I currently do -

znear = 0.00001
zfar =  100
K = array([[3871.88184, 0, 950.736938],[0, 3871.88184, 976.1383059999999],[0, 0, 1]])
Rt =array([[0.0225656671, 0.0194614234, 0.9995559233, 1.2656986283],[-0.9994773883, -0.0227084301, 0.0230060289, 0.5798922567],[0.0231460759, -0.99955269, 0.0189388219, -0.2110195758]])
ren.set_projection(K,zfar,znear)
ren.set_projection_from_camera(Rt)

The function being used are:

def set_projection(self,K,zfar,znear):
    glMatrixMode(GL_PROJECTION);
    glLoadIdentity();
    f_x = K[0,0]
    f_y = K[1,1]
    c_x = K[0,2]
    c_y = K[1,2]
    fovY = 1/(float(f_x)/height * 2);
    aspectRatio = (float(width)/height) * (float(f_y)/f_x);
    near = zfar
    far = znear
    frustum_height = near * fovY;
    frustum_width = frustum_height * aspectRatio;
    offset_x = (width/2 - c_x)/width * frustum_width * 2;
    offset_y = (height/2 - c_y)/height * frustum_height * 2;
    glFrustum(-frustum_width - offset_x, frustum_width - offset_x, -frustum_height - offset_y, frustum_height - offset_y, near, far);


def set_modelview_from_camera(self,Rt):
    glMatrixMode(GL_MODELVIEW)
    glLoadIdentity()
    Rx = array([[1,0,0],[0,0,-1],[0,1,0]])
    R = Rt[:,:3]
    U,S,V = linalg.svd(R)
    R = dot(U,V)
    R[0,:]=-R[0,:]
    t=Rt[:,3]
    M=eye(4)
    M[:3,:3]=dot(R,Rx)
    M[:3,3]=t
    M=M.T
    m=M.flatten()
    glLoadMatrixf(m)

Then I just render points (pasting snippet):

def renderLIDAR(self,filename):
    glClear(GL_COLOR_BUFFER_BIT|GL_DEPTH_BUFFER_BIT)
    glPushMatrix();

    glEnable(GL_DEPTH_TEST)
    glClear(GL_DEPTH_BUFFER_BIT)
    glPointSize(1.0)
    f = open(filename,'r')
    f.readline() #Contains number of particles
    for line in f:
        line = line.split(' ')
        glBegin(GL_POINTS)
        glColor3f (0.0,1.0,0.0); 
        x = float(line[0])
        y = float(line[1])
        z = float(line[2])
        glVertex3f(x,y,z)
        #print x,y,z
        glEnd()

    glPopMatrix();
mrkulk
  • 33
  • 1
  • 4

1 Answers1

3

The matrices you get back, most notably the last one in your question are what in OpenGL is the composition of projection and modelview, also called Modelviewprojection, i.e.

MVP = P · M

As long as you're not interested in performing illumination calculations, you can use just that in a vertex shader, i.e.

#version 330

uniform mat4 MVP;
in vec3 position;

void main()
{
    gl_Position = MVP * vec4(position, 1);
}

BTW, OpenGL and probably the library you're using as well, are using column major order, i.e. the order of the elements in memory is

0 4 8 c
1 5 9 d
2 6 a e
3 7 b f

so what's written in source code must be thought as "transposed" (of course it is not). Since the matrix you wrote follows the same scheme you can just put it into the uniform as it is. The only question that remains are the boundaries of the NDC space used by the range scanner. But that could be taken care of with an additional matrix applied. OpenGL uses the range [-1, 1]^3 so the worst thing that can happen is, that if it's in the other popular NDC range [0, 1]^3, you'll see your geometry just squeezed into the upper left hand corner of your window, and maybe turned "inside out" if the Z axis goes into the other direction. Just try it, I'd say it already matches OpenGL.

Anyway, if you want to use it with illumination, you have to decompose it into a projection and a modelview part. Easier said than done, but a good starting point is to orthonormalize the upper left 3×3 submatrix, which yields the rotational part of the modelview 'M'. You then have to find a matrix P, that, when left multiplied with M yields the original matrix. That's an overdetermined set of linear equations, so a Gauss-Jordan scheme can do it. And if I'm not entirely mistaken, what you already got in form of that camera matrix is either the decomposed M or P (I'd go for M).

Once you got that you may want to get the translational part (the 4th column) into the modelview matrix as well.

datenwolf
  • 159,371
  • 13
  • 185
  • 298
  • Thanks a lot datenwolf. As far as the matrices are concerned, the camera matrix (first one) is intrinsic camera parameter, the third (camera to reference) is extrinsic camera parameter. I assume with fixed openGL pipeline, what is the easiest way to set the matrices? I suppose one option is to multiply these in the shader but what if with illumination? – mrkulk Feb 11 '13 at 21:05
  • If you want working illumination you must decompose the matrices as outlined. Do you have some documentation about those matrices I could read? Technically for illumination to work you just need to extract the orthonormalized upper left 3×3 to transform the normals with. Unfortunately this works only with the programmable pipeline, as the fixed function pipeline make certain assumptions about the matrices. So for fixed function pipeline you must decompose it. The orthonormalized part is, what goes into modelview then, and the complementing part would be the projection. – datenwolf Feb 11 '13 at 22:07
  • Also missing are properly set near and far clip planes. – datenwolf Feb 11 '13 at 22:07
  • After decomposition, it seems to work properly. thanks for your help – mrkulk Feb 12 '13 at 09:05