Android Tutorial OpenGL 1.7 RayPlane Intersection

So ive been messing round with OpenGL for a while now as discovered that OpenGL ES doesnt include a 'GluUnproject' function, which is rather pivotal for discovering world coordinates with Ray Plane intersection. So I did some research and ported some Mesa code for the function and its related matrix manipulation functions. Heres what I discovered (hope it helps someone because i spent ages trying to get it to work)

To get the OpenGL world coordinates from a point we touch on screen we must make a "ray" going from the place we touch on the near clipping plane to the far clipping plane and see if it intersects anything along the way. Since this is just for the purpose of world coordinates the plane we intersect the ray with will just have a Y value of zero(so essentially we are just trying to find the point on a flat plane(or gameboard if youd prefer to think about it like that) lying down on the Z axis.

To make a ray we need the GluUnproject function, heres my implementation of it.

Sub UnProject(winx As Float, winy As Float, winz As Float, modelview() As Float, projection() As Float, viewport() As Int) As Int
Dim m(16) As Float
Dim A(16) As Float
Dim in(4) As Float
Dim out(4) As Float
Dim testmatrix As Matrix
Dim testmult As Matrix
testmult.MultiplyMM(A,0,projection,0,modelview,0)
'A = MultiplyMatrices4by4OpenGL_FLOAT(projection, modelview)

If testmatrix.InvertM(m,0,A,0) = False Then '  glhInvertMatrixf2(A,m) = 0 Then
   Return 0
End If

in(0) = (winx - viewport(0)) / viewport(2) * 2.0 - 1.0
in(1) = (winy - viewport(1)) / viewport(3) * 2.0 - 1.0
in(2) = 2.0 * winz - 1.0
in(3) = 1.0

out = MultiplyMatrixByVector4by4OpenGL_FLOAT(m, in)

If out(3) = 0.0 Then
   Return 0
End If

out(3) = 1.0 / out(3)

WorldXYZNear.x = out(0) * out(3)
WorldXYZNear.y = out(1) * out(3)
WorldXYZNear.z = out(2) * out(3)

Return 1

Unproject uses 'MultiplyMatrixByVector4by4OpenGL_FLOAT', so heres my port of it from Mesa.

Sub MultiplyMatrices4by4OpenGL_FLOAT(matrix1() As Float, Matrix2() As Float) As Float()
     Dim Result(16) As Float
    Result(0) = matrix1(0) * Matrix2(0) + matrix1(4) * Matrix2(1) + matrix1(8) * Matrix2(2) + matrix1(12) * Matrix2(3)
    Result(4) = matrix1(0) * Matrix2(4) + matrix1(4) * Matrix2(5) + matrix1(8) * Matrix2(6) + matrix1(12) * Matrix2(7)
    Result(8) = matrix1(0) * Matrix2(8) + matrix1(4) * Matrix2(9) + matrix1(8) * Matrix2(10)+ matrix1(12) * Matrix2(11)
    Result(12)= matrix1(0) * Matrix2(12)+ matrix1(4) * Matrix2(13)+ matrix1(8) * Matrix2(14)+ matrix1(12) * Matrix2(15)
    Result(1) = matrix1(1) * Matrix2(0) + matrix1(5) * Matrix2(1) + matrix1(9) * Matrix2(2) + matrix1(13) * Matrix2(3)
    Result(5) = matrix1(1) * Matrix2(4) + matrix1(5) * Matrix2(5) + matrix1(9) * Matrix2(6) + matrix1(13) * Matrix2(7)
    Result(9) = matrix1(1) * Matrix2(8) + matrix1(5) * Matrix2(9) + matrix1(9) * Matrix2(10)+ matrix1(13) * Matrix2(11)
    Result(13)= matrix1(1) * Matrix2(12)+ matrix1(5) * Matrix2(13)+ matrix1(9) * Matrix2(14)+ matrix1(13) * Matrix2(15)
    Result(2) = matrix1(2) * Matrix2(0) + matrix1(6) * Matrix2(1) + matrix1(10)* Matrix2(2) + matrix1(14) * Matrix2(3)
    Result(6) = matrix1(2) * Matrix2(4) + matrix1(6) * Matrix2(5) + matrix1(10)* Matrix2(6) + matrix1(14) * Matrix2(7)
    Result(10)= matrix1(2) * Matrix2(8) + matrix1(6) * Matrix2(9) + matrix1(10)* Matrix2(10)+ matrix1(14) * Matrix2(11)
    Result(14)= matrix1(2) * Matrix2(12)+ matrix1(6) * Matrix2(13)+ matrix1(10)* Matrix2(14)+ matrix1(14) * Matrix2(15)
    Result(3) = matrix1(3) * Matrix2(0) + matrix1(7) * Matrix2(1) + matrix1(11)* Matrix2(2) + matrix1(15) * Matrix2(3)
    Result(7) = matrix1(3) * Matrix2(4) + matrix1(7) * Matrix2(5) + matrix1(11)* Matrix2(6) + matrix1(15) * Matrix2(7)
    Result(11)= matrix1(3) * Matrix2(8) + matrix1(7) * Matrix2(9) + matrix1(11)* Matrix2(10)+ matrix1(15) * Matrix2(11)
    Result(15)= matrix1(3) * Matrix2(12)+ matrix1(7) * Matrix2(13)+ matrix1(11)* Matrix2(14)+ matrix1(15) * Matrix2(15)
   Return Result

End Sub

Im aware I could have used the inbuilt Matrix manipulation in the OpenGL library but I wanted to keep things relatively native to Mesa.

To create the ray we need to Unproject twice, so I have another Unproject function renamed slightly differently so it outputs to :
WorldXYZFar.x = out(0) * out(3)
WorldXYZFar.y = out(1) * out(3)
WorldXYZFar.z = out(2) * out(3)

I created a few custom types to help with implementing code on the net - they are:

Type Vector3(x As Float, y As Float, z As Float)
Type Vector2(x As Float, y As Float)
Type Ray(Pos As Vector3, Direction As Vector3)
Type Plane(Normal As Vector3, d As Float)

We have our 2 points we create a ray out of through unprojecting twice. We can now intersect that with the Y = 0 plane to discover the world coordinate of where we touched.

Sub CalculateWorldCoordinates(Position As Vector2) As Vector3
Dim result As Vector3

isworking = UnProject(Position.x,m_viewport(3) - Position.y,-1,m_modelview,m_projection,m_viewport)
isworking = Unproject1(Position.x,m_viewport(3) - Position.y,1,m_modelview,m_projection,m_viewport)

Dim direction As Vector3
direction.x = WorldXYZFar.x - WorldXYZNear.x
direction.y = WorldXYZFar.y - WorldXYZNear.y
direction.z = WorldXYZFar.z - WorldXYZNear.z

Dim raylength As Float
raylength = Sqrt((WorldXYZNear.x * WorldXYZNear.x) + (WorldXYZNear.y * WorldXYZNear.y) + (WorldXYZNear.z * WorldXYZNear.z))


'normalize
direction.x = direction.x / raylength
direction.y = direction.y / raylength
direction.z = direction.z / raylength

Dim r As Ray : r.Pos = WorldXYZNear : r.direction = direction

Dim n As Vector3 : n.x = 0 : n.y = 1: n.z = 0

Dim p As Plane : p.Normal = n : p.d = 0

Dim denominator, numerator As Float
denominator = Dot(p.normal,r.direction) + p.d
numerator = Dot(p.normal,r.Pos)

Dim t As Float
t = -(numerator / denominator)

result.x = WorldXYZNear.x + direction.x * t
result.y = WorldXYZNear.y + direction.y * t
result.z = WorldXYZNear.z + direction.z * t

Return result

I am not a mathematician so dont ask me the nitty gritty details behind the math, in fact if anyone could enlighten me on what the hell a Dot Product computes and its meaning then id be forever in your debt!

Sub Dot(A As Vector3, B As Vector3) As Float
Return A.x * B.x + A.y * B.y + A.z * B.z
End Sub

The algorithm to intersect with the plane is just a port of some XNA code, it can be found on the internet as well as other techniques.


Some pointers:
You have to pass the viewport/projection/modelview matrices to the Unproject function, you MUST grab the viewport/projection matrices at the point of creation, and grabbing the modelview matrix should be done in the main drawing routine as I 'think' it modelview matrix is cumulative so if you apply transformations to it then the values contained within the matrix will change.
I grab my viewport/projection matrices at their point of creation.


Sub glsv_SurfaceChanged(gl As GL1, width As Int, height As Int)
gl.glViewport(0, 0, width, height)
gl.glGetIntegerv(gl.GL_VIEWPORT,m_viewport,0)
gl.glMatrixMode(gl.GL_PROJECTION)
gl.glLoadIdentity
Dim ratio As Float
ratio = width/height
gl.gluPerspective(45,screenratio,0.5,100)
gl.gluLookAt(0,4,5,0,0,0,0,1,0)
gl.glgetfloatv(gl.GL_PROJECTION_MATRIX,m_projection,0)
gl.glMatrixMode(gl.GL_MODELVIEW)
gl.glLoadIdentity
End Sub


The text in bold is how you grab the respective matrices, the viewport matrix needs an array of integers with 4 elements, the projection needs an array of floats with 16 elements.

The text in italics - For a long time my code seemed to not be working, and then when I put in gl.gluLookAt(0,4,5,0,0,0,0,1,0) it all seemed to work for some reason. I put it down to beginners luck.

Anyway! Just thought id share my discoveries as it took me a while to work this out as everything on the net is in java.

Thanks!

 

[Erel]

Thank you for sharing.