Can someone check this triangle intersection code?

Blitz3D Forums/Blitz3D Programming/Can someone check this triangle intersection code?

sswift(Posted 2003) [#1]
I'm not sure if the math at the end of the ray intersect triangle function is correct for determining if the point is beyond the endpoint of the ray defined by Dxyz.

I'm attempting to make a custom linepick for that lightmapper I've beem improving on, because if I can use a linepick instead of entityvisible then not only will the code be a lot faster, but I can create realistic soft shadows with light sources that have area.

Then again maybe I shouild use a different method... Hm....

But I would like this function to work anyhow cause it seems really fast, though it also is broken somehow.


; -------------------------------------------------------------------------------------------------------------------
; This function returns TRUE if a ray intersects a triangle.
; 	It also calculates the UV coordinates of said colision as part of the intersection test,
; 	but does not return them.
;
; V0xyz, V1xyz, and V2xyz are the locations of the three vertices of the triangle.
;
; V0 is location of UV(0,0).  V1 is the location of UV(1, 0).  V2 is the location of UV(0,1)
; 	These are important to know if you want to return the exact location in texture space of the collision, but
; 	you don't have to worry about them if you only want to find out if a collision occured.
;
; Pxyz is a the start of the line.
;	Triangles which are "behind" this point are ignored. 
; 	"Behind" is defined as the direction opposite that which Dxyz points in.
;
; Dxyz is a vector providing the slope of the line.
; 	(Does this have to be a normal?  I'm not sure.)
;
; If Extend_To_Infinity is set to false, then the length of Dxyz is how far the ray extends. 
; 	So if you want an endpoint on your ray beyond which no triangles will be detected, subtract the position 
; 	of Pxyz from your endpoint's position, and pass that ato the function as Dxyz.  Ie: (Dx = P2x-P1x)
;
; If Cull_Backface is set to true, then if the specified ray passes through the triangle from it's back side, then
; it will not register that it hit that triangle.
; -------------------------------------------------------------------------------------------------------------------
Function Ray_Intersect_Triangle(Px#, Py#, Pz#, Dx#, Dy#, Dz#, V0x#, V0y#, V0z#, V1x#, V1y#, V1z#, V2x#, V2y#, V2z#, Extend_To_Infinity=True, Cull_Backfaces=False)


	; To do:
	; Cull triangles behind Pxyz.
	; Cull triangles beyond Dxyz If Extend_To_Infinity = False
	
	; crossproduct(b,c) =
	; ax = (by * cz) - (cy * bz) 
	; ay = (bz * cx) - (cz * bx) 	
	; az = (bx * cy) - (cx * by)

	; dotproduct(v,q) =
	; (vx * qx) + (vy * qy) + (vz * qz)	
	; DP =  1 = Vectors point in same direction.          (  0 degrees of seperation)
	; DP =  0 = Vectors are perpendicular to one another. ( 90 degrees of seperation)
	; DP = -1 = Vectors point in opposite directions.     (180 degrees of seperation) 
	;
	; The dot product is also reffered to as "the determinant" or "the inner product"


	; Calculate the vector that represents the first side of the triangle.
	E1x# = V1x# - V0x#
	E1y# = V1y# - V0y#
	E1z# = V1z# - V0z#


	; Calculate the vector that represents the second side of the triangle.
	E2x# = V2x# - V0x#
	E2y# = V2y# - V0y#
	E2z# = V2z# - V0z#


	; Calculate the triangle's normal by finding a vector which is perpendicular to two of it's sides.
	; (This "normal" is not actually normalized.)
	; Hxyz = Crossproduct(Dxyz, E2xyz)
	Hx# = (Dy# * E2z#) - (E2y# * Dz#)
	Hy# = (Dz# * E2x#) - (E2z# * Dx#)
	Hz# = (Dx# * E2y#) - (E2x# * Dy#)


	; Calculate the dotproduct of the ray and the triangle's normal.
	A# = (E1x# * Hx#) + (E1y# * Hy#) + (E1z# * Hz#)


	; If we should ignore triangles the ray passes through the back side of,
	; and the ray points in the same direction as the normal of the plane,
	; then the ray passed through the back side of the plane,  
	; and the ray does not intersect the plane the triangle lies in.
	If (Cull_Backface = True) And (A# >= 0) Then Return False
	
		
	; If the ray is almost parralel to the plane,
	; then the ray does not intersect the plane the triangle lies in.
	If (A# > -0.00001) And (A# < 0.00001) Then Return False
	
	
	; Inverse Determinant(Dot Product). 
	; (Scaling factor for UV's?)
	F# = 1.0 / A#


	; Calculate a vector between the starting point of our ray, and the first point of the triangle,
	; which is at UV(0,0)
	; (Note:  I'm not sure why this is done!)
	Sx# = Px# - V0x#
	Sy# = Py# - V0y#
	Sz# = Pz# - V0z#
	
	
	; Calculate the U coordinate of the intersection point.
	;
	;	Sxyz is the vector between the start of our ray and the first point of the triangle.
	;	Hxyz is the normal of our triangle.
	;	
	; U# = F# * (DotProduct(Sxyz, Hxyz))
	U# = F# * ((Sx# * Hx#) + (Sy# * Hy#) + (Sz# * Hz#))
	
	
	; Is the U coordinate outside the range of values inside the triangle?
	If (U# < 0.0) Or (U# > 1.0)

		; The ray has intersected the plane outside the triangle.
		Return False
	
	EndIf


	; I think this is the intersection point?
	;
	;	Sxyz is the vector from the starting point of the ray to the first corner of the triangle.
	;	E1xyz is the vector which represents the first side of the triangle.
	;	The crossproduct of these two would be a vector which is perpendicular to both.
	;
	; Qxyz = CrossProduct(Sxyz, E1xyz)
	Qx# = (Sy# * E1z#) - (E1y# * Sz#)
	Qy# = (Sz# * E1x#) - (E1z# * Sx#)
	Qz# = (Sx# * E1y#) - (E1x# * Sy#)

	
	; Calculate the V coordinate of the intersection point.
	;	
	;	Dxyz is the vector which represents the direction the ray is pointing in.
	;	Qxyz is the intersection point I think?
	;
	; V# = F# * DotProduct(Dxyz, Qxyz)
	V# = F# * ((Dx# * Qx#) + (Dy# * Qy#) + (Dz# * Qz#))

	
	; Is the V coordinate outside the range of values inside the triangle?	
	; Does U+V exceed 1.0?  
	If (V# < 0.0) Or ((U# + V#) > 1.0)

		; The ray has intersected the plane outside the triangle.		
		Return(False)

		; The reason we check U+V is because if you imagine the triangle as half a square, U=1 V=1 would
		; be in the lower left hand corner which would be in the lower left triangle making up the square.
		; We are looking for the upper right triangle, and if you think about it, U+V will always be less
		; than or equal to 1.0 if the point is in the upper right of the triangle.

	EndIf


	; Calculate the distance of the intersection point from the starting point of the ray, Pxyz.
	;
	; I don't know if this is the true distance or some scaled distance because none of the vectors are normalized,
	; but I think it is safe to assume that if it is negative, then the plane is behind Pxyz on the line.
	T# = F# * ((E2x# * Qx#) + (E2y# * Qy#) + (E2z# * Qz#))


	; If the distance to the plane is negative, then the triangle is behind Pxyz. 
	; Because we want a directional ray, if the triangle is behind Pxyz, then we ignore this intersection.
	If (T# < 0) Then Return False
	

	; Ray-triangle intersection point equation:	
	; p + t * d = (1-u-v) * V0 + u * V1 + v * V2


	If Extend_To_Infinity = False

		; Calculate T for the intersection point UV on the triangle defined by V0xyz, V1xyz, and V2xyz.
		; There's gotta be a better way to do this...
		If Dx# <> 0 
			T2# = (((1.0 - (U# + V#)) * V0x#) + (U# * V1x#) + (V# * V2x#) - Px#) / Dx#
		Else
			If Dy# <> 0 	
				T2# = (((1.0 - (U# + V#)) * V0y#) + (U# * V1y#) + (V# * V2y#) - Py#) / Dy#
			Else
				T2# = (((1.0 - (U# + V#)) * V0z#) + (U# * V1z#) + (V# * V2z#) - Pz#) / Dz#
			EndIf
		EndIf
			
		; If the intersection point is farther along the line than the line's endpoint, then no intersection occured.	
		If T2# > T# Then Return False
	
	EndIf
		

	; The intersection point lies within the triangle.
	; The ray intersects the triangle!		
	Return True

End Function


; -------------------------------------------------------------------------------------------------------------------
; This function returns true if a ray intersects a mesh.
;
; This function differs from LinePick in that the specified mesh does not need to have a pickmode set,
; and this function can optionally ignore backfacing polygons in the mesh.  So if you do a pick from 
; inside a mesh, it will not register a hit.
; -------------------------------------------------------------------------------------------------------------------
Function Ray_Intersect_Mesh(Mesh, Px#, Py#, Pz#, Dx#, Dy#, Dz#, Extend_To_Infinity=True, Cull_Backfaces=False)

	; Make sure there's a surface, because the mesh might be empty.
	If CountSurfaces(Mesh) > 0

		Surface = GetSurface(Mesh, 1)
	
		; Examine all triangles in this surface.	
		Tris  = CountTriangles(Surface)
		For TriLoop = 0 To Tris-1
	
			V0 = TriangleVertex(Surface, TriLoop, 0)
			V1 = TriangleVertex(Surface, TriLoop, 1)
			V2 = TriangleVertex(Surface, TriLoop, 2)
		
			V0x# = VertexX#(Surface, V0)
			V0y# = VertexY#(Surface, V0)
			V0z# = VertexZ#(Surface, V0)

			V1x# = VertexX#(Surface, V1)
			V1y# = VertexY#(Surface, V1)
			V1z# = VertexZ#(Surface, V1)

			V2x# = VertexX#(Surface, V2)
			V2y# = VertexY#(Surface, V2)
			V2z# = VertexZ#(Surface, V2)
			
			Intersected = Ray_Intersect_Triangle(Px#, Py#, Pz#, Dx#, Dy#, Dz#, V0x#, V0y#, V0z#, V1x#, V1y#, V1z#, V2x#, V2y#, V2z#, Extend_To_Infinity, Cull_Backfaces)
			If Intersected Then Return True
		
		Next
		
	EndIf
	
	Return False

End Function



JoshK(Posted 2003) [#2]
Good function to have. I have a RayIntersectsPlane function in the maths3d.dll lib, but no triangle yet.

If all you need to do now is check the range of the ray, just find the distance between the start of the ray and the intersected point(?)


Bot Builder(Posted 2003) [#3]
ohh. I was going to say that there's a line-tri dll, but I didn't realize it was giving UVs. looks good though. perhaps, instead of checking the distance of the ray with a threshold value, you could have an equation giving the intensity of the shadow(farther away, less intense shadow). while I realize that this is not a method that truely represents reality, it's a good way to give the effect of rays bouncing around.

[blabber]
I always thought it'd be cool to make a lightmapper. Maybe realtime radiosity (sounds like an oxymoron eh?). If I do try, I'll take it on the same way I'm taking on my physics system. simply as an educational and possibly useful challenge. I'm just blabbering. Maybe I'll do a sophesticated particle system... been done, but It'd be another fun challenge. hehehe. maybe sometime I could take all these "fun challenges" and make a game. Shadows-physics-gui system-particlew system-ohhhh maybe even a cutscene scripter! hehehe. That'd make a pretty slow game unless I optimize my systems sevierly which I intend to...........
[/blabber]


sswift(Posted 2003) [#4]
"If all you need to do now is check the range of the ray, just find the distance between the start of the ray and the intersected point(?)"

And T# is supposed to be that distance. When you do this function you get TUV, not XYZ.

That's the source of my difficulty. I'm not sure if T is in world coordinates or some other coordinate system that's scaled by some amount. If it's in world coordinates, things are easy. I suppose I'll have to write a test program. I just thought I'd ask last night before I went to bed.

Also I want to find a way to convert TUV to XYZ so that I can return that too.

Hm. Maybe I should have just re-written this intersection code from scratch and used plane math instead of this UV stuff, but it is useful to know the UV.

Bot Builder:
It's not that simple. I need to constrain the ray because I am casting rays from the lumels TO the light source. So they have to start at the lumel and stop at the light source, otherwise they'll detect the walls beyond the light source and falsely assume that that lumel is in shadow.

And radiosity is something I've been considering myself. I might take a whack at it if I can come up with a method which is easy and fast enough for my tastes.


sswift(Posted 2003) [#5]
Okay I figured out what sort of distance T represents with some test code.

T is scaled so that at Pxyz, the start of the ray, T=0, and at Dxyz, the end of the ray, T=1. If the intersection point is behind Pxyz, then T will be negative, and if the intersection point is beyond Dxyz then T will be greater than 1.

So this makes it really simple to know if the triangle is beyond either end of the ray.


patmaba(Posted 2006) [#6]
I think you have a small error in your Function Ray_Intersect_Triangle with the last parameters Cullback_faces.

The function used a condition with the parameter i think Cullback_faces but you have write "Cullback_face" without S

Conclusion, you cannnot test this parameter, i think.

Is it Just ?


Moreover I have rename your code into Inc_Ray.bb

Include "Inc_Ray.bb"

; CreateCube Example
; ------------------

Graphics3D 640,480
SetBuffer BackBuffer()

camera=CreateCamera()

light=CreateLight()
RotateEntity light,90,0,0

; Create cube
cube=CreateCube()
PositionEntity cube,0,0,5

PositionEntity camera,0,2,0
PointEntity camera,cube


While Not KeyDown( 1 )
  If Ray_Intersect_Mesh(cube, EntityX(camera), EntityY(camera), EntityZ(camera), 0, 0, 1, True, False)
        EntityColor cube, Rnd(255),0,0
  Else 
        EntityColor cube, 255,255,255
  EndIf 
  TurnEntity cube, 1, 1,-1
RenderWorld
Flip
Wend

End


Normally my cube must be red, but it's always white.

My condition is it correct with your ray function ?


Beaker(Posted 2006) [#7]
This thread is 3 years old. Look at the latest version in the code archive.


Defoc8(Posted 2006) [#8]
lol..i didnt even notice :] - well spotted :p ;)


patmaba(Posted 2006) [#9]
Yes, your latest version use TFormPoint and Tformedxyz function.

I need a solution without this function! why ? To implement a similar function in BlitzMax for MiniB3d. :)

Minib3d for Bmax does provide TFormPoint, TFormedxyz and LinePick !:(

So, i need a mathematical function Ray_Intersect_Triangle to simulate a linepick-like function.