Light probe volume generation

Source/BansheeUtility/Include/BsVector3.h

@@ -14,22 +14,22 @@ namespace bs
 	 */
 
 	/** A three dimensional vector. */
-    class BS_UTILITY_EXPORT Vector3
-    {
-    public:
+	class BS_UTILITY_EXPORT Vector3
+	{
+	public:
 		float x, y, z;
 
-    public:
+	public:
 		Vector3()
 		{ }
 
-        Vector3(BS_ZERO zero)
+		Vector3(BS_ZERO zero)
 			:x(0.0f), y(0.0f), z(0.0f)
-        { }
+		{ }
 
-        Vector3(float x, float y, float z)
-            :x(x), y(y), z(z)
-        { }
+		Vector3(float x, float y, float z)
+			:x(x), y(y), z(z)
+		{ }
 
 		explicit Vector3(const Vector4& vec);
 
@@ -42,18 +42,18 @@ namespace bs
 		}
 
 		float operator[] (UINT32 i) const
-        {
-            assert(i < 3);
+		{
+			assert(i < 3);
 
-            return *(&x+i);
-        }
+			return *(&x + i);
+		}
 
 		float& operator[] (UINT32 i)
-        {
-            assert(i < 3);
+		{
+			assert(i < 3);
 
-            return *(&x+i);
-        }
+			return *(&x + i);
+		}
 
 		/** Pointer accessor for direct copying. */
 		float* ptr()
@@ -67,214 +67,214 @@ namespace bs
 			return &x;
 		}
 
-        Vector3& operator= (const Vector3& rhs)
-        {
-            x = rhs.x;
-            y = rhs.y;
-            z = rhs.z;
-
-            return *this;
-        }
-
-        Vector3& operator= (float rhs)
-        {
-            x = rhs;
-            y = rhs;
-            z = rhs;
-
-            return *this;
-        }
-
-        bool operator== (const Vector3& rhs) const
-        {
-            return (x == rhs.x && y == rhs.y && z == rhs.z);
-        }
-
-        bool operator!= (const Vector3& rhs) const
-        {
-            return (x != rhs.x || y != rhs.y || z != rhs.z);
-        }
-
-        Vector3 operator+ (const Vector3& rhs) const
-        {
-            return Vector3(x + rhs.x, y + rhs.y, z + rhs.z);
-        }
-
-        Vector3 operator- (const Vector3& rhs) const
-        {
-            return Vector3(x - rhs.x, y - rhs.y, z - rhs.z);
-        }
-
-        Vector3 operator* (float rhs) const
-        {
-            return Vector3(x * rhs, y * rhs, z * rhs);
-        }
-
-        Vector3 operator* (const Vector3& rhs) const
-        {
-            return Vector3(x * rhs.x, y * rhs.y, z * rhs.z);
-        }
-
-        Vector3 operator/ (float val) const
-        {
-            assert(val != 0.0);
-
-            float fInv = 1.0f / val;
-            return Vector3(x * fInv, y * fInv, z * fInv);
-        }
-
-        Vector3 operator/ (const Vector3& rhs) const
-        {
-            return Vector3(x / rhs.x, y / rhs.y, z / rhs.z);
-        }
-
-        const Vector3& operator+ () const
-        {
-            return *this;
-        }
-
-        Vector3 operator- () const
-        {
-            return Vector3(-x, -y, -z);
-        }
-
-        friend Vector3 operator* (float lhs, const Vector3& rhs)
-        {
-            return Vector3(lhs * rhs.x, lhs * rhs.y, lhs * rhs.z);
-        }
-
-        friend Vector3 operator/ (float lhs, const Vector3& rhs)
-        {
-            return Vector3(lhs / rhs.x, lhs / rhs.y, lhs / rhs.z);
-        }
-
-        friend Vector3 operator+ (const Vector3& lhs, float rhs)
-        {
-            return Vector3(lhs.x + rhs, lhs.y + rhs, lhs.z + rhs);
-        }
-
-        friend Vector3 operator+ (float lhs, const Vector3& rhs)
-        {
-            return Vector3(lhs + rhs.x, lhs + rhs.y, lhs + rhs.z);
-        }
-
-        friend Vector3 operator- (const Vector3& lhs, float rhs)
-        {
-            return Vector3(lhs.x - rhs, lhs.y - rhs, lhs.z - rhs);
-        }
-
-        friend Vector3 operator- (float lhs, const Vector3& rhs)
-        {
-            return Vector3(lhs - rhs.x, lhs - rhs.y, lhs - rhs.z);
-        }
-
-        Vector3& operator+= (const Vector3& rhs)
-        {
-            x += rhs.x;
-            y += rhs.y;
-            z += rhs.z;
-
-            return *this;
-        }
-
-        Vector3& operator+= (float rhs)
-        {
-            x += rhs;
-            y += rhs;
-            z += rhs;
-
-            return *this;
-        }
-
-        Vector3& operator-= (const Vector3& rhs)
-        {
-            x -= rhs.x;
-            y -= rhs.y;
-            z -= rhs.z;
-
-            return *this;
-        }
-
-        Vector3& operator-= (float rhs)
-        {
-            x -= rhs;
-            y -= rhs;
-            z -= rhs;
-
-            return *this;
-        }
-
-        Vector3& operator*= (float rhs)
-        {
-            x *= rhs;
-            y *= rhs;
-            z *= rhs;
-
-            return *this;
-        }
-
-        Vector3& operator*= (const Vector3& rhs)
-        {
-            x *= rhs.x;
-            y *= rhs.y;
-            z *= rhs.z;
-
-            return *this;
-        }
-
-        Vector3& operator/= (float rhs)
-        {
-            assert(rhs != 0.0f);
-
-            float inv = 1.0f / rhs;
-
-            x *= inv;
-            y *= inv;
-            z *= inv;
-
-            return *this;
-        }
-
-        Vector3& operator/= (const Vector3& rhs)
-        {
-            x /= rhs.x;
-            y /= rhs.y;
-            z /= rhs.z;
-
-            return *this;
-        }
-
-        /** Returns the length (magnitude) of the vector. */
+		Vector3& operator= (const Vector3& rhs)
+		{
+			x = rhs.x;
+			y = rhs.y;
+			z = rhs.z;
+
+			return *this;
+		}
+
+		Vector3& operator= (float rhs)
+		{
+			x = rhs;
+			y = rhs;
+			z = rhs;
+
+			return *this;
+		}
+
+		bool operator== (const Vector3& rhs) const
+		{
+			return (x == rhs.x && y == rhs.y && z == rhs.z);
+		}
+
+		bool operator!= (const Vector3& rhs) const
+		{
+			return (x != rhs.x || y != rhs.y || z != rhs.z);
+		}
+
+		Vector3 operator+ (const Vector3& rhs) const
+		{
+			return Vector3(x + rhs.x, y + rhs.y, z + rhs.z);
+		}
+
+		Vector3 operator- (const Vector3& rhs) const
+		{
+			return Vector3(x - rhs.x, y - rhs.y, z - rhs.z);
+		}
+
+		Vector3 operator* (float rhs) const
+		{
+			return Vector3(x * rhs, y * rhs, z * rhs);
+		}
+
+		Vector3 operator* (const Vector3& rhs) const
+		{
+			return Vector3(x * rhs.x, y * rhs.y, z * rhs.z);
+		}
+
+		Vector3 operator/ (float val) const
+		{
+			assert(val != 0.0);
+
+			float fInv = 1.0f / val;
+			return Vector3(x * fInv, y * fInv, z * fInv);
+		}
+
+		Vector3 operator/ (const Vector3& rhs) const
+		{
+			return Vector3(x / rhs.x, y / rhs.y, z / rhs.z);
+		}
+
+		const Vector3& operator+ () const
+		{
+			return *this;
+		}
+
+		Vector3 operator- () const
+		{
+			return Vector3(-x, -y, -z);
+		}
+
+		friend Vector3 operator* (float lhs, const Vector3& rhs)
+		{
+			return Vector3(lhs * rhs.x, lhs * rhs.y, lhs * rhs.z);
+		}
+
+		friend Vector3 operator/ (float lhs, const Vector3& rhs)
+		{
+			return Vector3(lhs / rhs.x, lhs / rhs.y, lhs / rhs.z);
+		}
+
+		friend Vector3 operator+ (const Vector3& lhs, float rhs)
+		{
+			return Vector3(lhs.x + rhs, lhs.y + rhs, lhs.z + rhs);
+		}
+
+		friend Vector3 operator+ (float lhs, const Vector3& rhs)
+		{
+			return Vector3(lhs + rhs.x, lhs + rhs.y, lhs + rhs.z);
+		}
+
+		friend Vector3 operator- (const Vector3& lhs, float rhs)
+		{
+			return Vector3(lhs.x - rhs, lhs.y - rhs, lhs.z - rhs);
+		}
+
+		friend Vector3 operator- (float lhs, const Vector3& rhs)
+		{
+			return Vector3(lhs - rhs.x, lhs - rhs.y, lhs - rhs.z);
+		}
+
+		Vector3& operator+= (const Vector3& rhs)
+		{
+			x += rhs.x;
+			y += rhs.y;
+			z += rhs.z;
+
+			return *this;
+		}
+
+		Vector3& operator+= (float rhs)
+		{
+			x += rhs;
+			y += rhs;
+			z += rhs;
+
+			return *this;
+		}
+
+		Vector3& operator-= (const Vector3& rhs)
+		{
+			x -= rhs.x;
+			y -= rhs.y;
+			z -= rhs.z;
+
+			return *this;
+		}
+
+		Vector3& operator-= (float rhs)
+		{
+			x -= rhs;
+			y -= rhs;
+			z -= rhs;
+
+			return *this;
+		}
+
+		Vector3& operator*= (float rhs)
+		{
+			x *= rhs;
+			y *= rhs;
+			z *= rhs;
+
+			return *this;
+		}
+
+		Vector3& operator*= (const Vector3& rhs)
+		{
+			x *= rhs.x;
+			y *= rhs.y;
+			z *= rhs.z;
+
+			return *this;
+		}
+
+		Vector3& operator/= (float rhs)
+		{
+			assert(rhs != 0.0f);
+
+			float inv = 1.0f / rhs;
+
+			x *= inv;
+			y *= inv;
+			z *= inv;
+
+			return *this;
+		}
+
+		Vector3& operator/= (const Vector3& rhs)
+		{
+			x /= rhs.x;
+			y /= rhs.y;
+			z /= rhs.z;
+
+			return *this;
+		}
+
+		/** Returns the length (magnitude) of the vector. */
 		float length() const
 		{
 			return std::sqrt(x * x + y * y + z * z);
 		}
 
-        /** Returns the square of the length(magnitude) of the vector. */
-        float squaredLength() const
-        {
-            return x * x + y * y + z * z;
-        }
+		/** Returns the square of the length(magnitude) of the vector. */
+		float squaredLength() const
+		{
+			return x * x + y * y + z * z;
+		}
 
-        /**	Returns the distance to another vector. */
-        float distance(const Vector3& rhs) const
-        {
-            return (*this - rhs).length();
-        }
+		/**	Returns the distance to another vector. */
+		float distance(const Vector3& rhs) const
+		{
+			return (*this - rhs).length();
+		}
 
-        /** Returns the square of the distance to another vector. */
-        float squaredDistance(const Vector3& rhs) const
-        {
-            return (*this - rhs).squaredLength();
-        }
+		/** Returns the square of the distance to another vector. */
+		float squaredDistance(const Vector3& rhs) const
+		{
+			return (*this - rhs).squaredLength();
+		}
 
-        /** Calculates the dot (scalar) product of this vector with another. */
-        float dot(const Vector3& vec) const
-        {
-            return x * vec.x + y * vec.y + z * vec.z;
-        }
+		/** Calculates the dot (scalar) product of this vector with another. */
+		float dot(const Vector3& vec) const
+		{
+			return x * vec.x + y * vec.y + z * vec.z;
+		}
 
-        /** Normalizes the vector. */
+		/** Normalizes the vector. */
 		float normalize()
 		{
 			float len = length();
@@ -291,60 +291,60 @@ namespace bs
 		}
 
 
-        /** Calculates the cross-product of 2 vectors, that is, the vector that lies perpendicular to them both. */
-        Vector3 cross(const Vector3& other) const
-        {
-            return Vector3(
-                y * other.z - z * other.y,
-                z * other.x - x * other.z,
-                x * other.y - y * other.x);
-        }
+		/** Calculates the cross-product of 2 vectors, that is, the vector that lies perpendicular to them both. */
+		Vector3 cross(const Vector3& other) const
+		{
+			return Vector3(
+				y * other.z - z * other.y,
+				z * other.x - x * other.z,
+				x * other.y - y * other.x);
+		}
 
-        /** Sets this vector's components to the minimum of its own and the ones of the passed in vector. */
-        void floor(const Vector3& cmp)
-        {
-            if(cmp.x < x) x = cmp.x;
-            if(cmp.y < y) y = cmp.y;
-            if(cmp.z < z) z = cmp.z;
-        }
+		/** Sets this vector's components to the minimum of its own and the ones of the passed in vector. */
+		void floor(const Vector3& cmp)
+		{
+			if (cmp.x < x) x = cmp.x;
+			if (cmp.y < y) y = cmp.y;
+			if (cmp.z < z) z = cmp.z;
+		}
 
-        /** Sets this vector's components to the maximum of its own and the ones of the passed in vector. */
-        void ceil(const Vector3& cmp)
-        {
-            if(cmp.x > x) x = cmp.x;
-            if(cmp.y > y) y = cmp.y;
-            if(cmp.z > z) z = cmp.z;
-        }
+		/** Sets this vector's components to the maximum of its own and the ones of the passed in vector. */
+		void ceil(const Vector3& cmp)
+		{
+			if (cmp.x > x) x = cmp.x;
+			if (cmp.y > y) y = cmp.y;
+			if (cmp.z > z) z = cmp.z;
+		}
 
-        /** Generates a vector perpendicular to this vector. */
-        Vector3 perpendicular() const
-        {
-            static const float squareZero = (float)(1e-06 * 1e-06);
+		/** Generates a vector perpendicular to this vector. */
+		Vector3 perpendicular() const
+		{
+			static const float squareZero = (float)(1e-06 * 1e-06);
 
-            Vector3 perp = this->cross(Vector3::UNIT_X);
+			Vector3 perp = this->cross(Vector3::UNIT_X);
 
-            if(perp.squaredLength() < squareZero)
-                perp = this->cross(Vector3::UNIT_Y);
+			if (perp.squaredLength() < squareZero)
+				perp = this->cross(Vector3::UNIT_Y);
 
 			perp.normalize();
-            return perp;
-        }
+			return perp;
+		}
 
 		/** Gets the angle between 2 vectors. */
 		inline Radian angleBetween(const Vector3& dest) const;
 
-        /** Returns true if this vector is zero length. */
-        bool isZeroLength() const
-        {
-            float sqlen = (x * x) + (y * y) + (z * z);
-            return (sqlen < (1e-06 * 1e-06));
-        }
+		/** Returns true if this vector is zero length. */
+		bool isZeroLength() const
+		{
+			float sqlen = (x * x) + (y * y) + (z * z);
+			return (sqlen < (1e-06 * 1e-06));
+		}
 
-        /** Calculates a reflection vector to the plane with the given normal. */
-        Vector3 reflect(const Vector3& normal) const
-        {
-            return Vector3(*this - (2 * this->dot(normal) * normal));
-        }
+		/** Calculates a reflection vector to the plane with the given normal. */
+		Vector3 reflect(const Vector3& normal) const
+		{
+			return Vector3(*this - (2 * this->dot(normal) * normal));
+		}
 
 		/** Calculates two vectors orthonormal to the current vector, and normalizes the current vector if not already. */
 		void orthogonalComplement(Vector3& a, Vector3& b)
@@ -364,7 +364,7 @@ namespace bs
 		{
 			vec0.normalize();
 
-			float dot0 = vec0.dot(vec1); 
+			float dot0 = vec0.dot(vec1);
 			vec1 -= dot0*vec0;
 			vec1.normalize();
 
@@ -374,7 +374,7 @@ namespace bs
 			vec2.normalize();
 		}
 
-        /** Calculates the dot (scalar) product of two vectors. */
+		/** Calculates the dot (scalar) product of two vectors. */
 		static float dot(const Vector3& a, const Vector3& b)
 		{
 			return a.x * b.x + a.y * b.y + a.z * b.z;
@@ -385,16 +385,16 @@ namespace bs
 
 		/** Calculates the cross-product of 2 vectors, that is, the vector that lies perpendicular to them both. */
 		static Vector3 cross(const Vector3& a, const Vector3& b)
-        {
-            return Vector3(
+		{
+			return Vector3(
 				a.y * b.z - a.z * b.y,
 				a.z * b.x - a.x * b.z,
 				a.x * b.y - a.y * b.x);
-        }
+		}
 
-		/** 
+		/**
 		 * Linearly interpolates between the two vectors using @p t. t should be in [0, 1] range, where t = 0 corresponds
-		 * to the left vector, while t = 1 corresponds to the right vector. 
+		 * to the left vector, while t = 1 corresponds to the right vector.
 		 */
 		static Vector3 lerp(float t, const Vector3& a, const Vector3& b)
 		{
@@ -416,13 +416,13 @@ namespace bs
 			return Vector3(std::max(a.x, b.x), std::max(a.y, b.y), std::max(a.z, b.z));
 		}
 
-        static const Vector3 ZERO;
+		static const Vector3 ZERO;
 		static const Vector3 ONE;
 		static const Vector3 INF;
-        static const Vector3 UNIT_X;
-        static const Vector3 UNIT_Y;
-        static const Vector3 UNIT_Z;
-    };
+		static const Vector3 UNIT_X;
+		static const Vector3 UNIT_Y;
+		static const Vector3 UNIT_Z;
+	};
 
 	/** @} */
 

Source/BansheeUtility/Source/BsTriangulation.cpp

@@ -41,6 +41,7 @@ namespace bs
 			memcpy(volume.tetrahedra[i].neighbors, &output.neighborlist[i * 4], sizeof(INT32) * 4);
 		}
 
+		// Generate boundary faces
 		UINT32 numFaces = (UINT32)output.numberoftrifaces;
 		for (UINT32 i = 0; i < numFaces; ++i)
 		{
@@ -59,6 +60,33 @@ namespace bs
 			face.tetrahedron = tetIdx;
 		}
 
+		// Ensure that vertex at the specified location points a neighbor opposite to it
+		for(UINT32 i = 0; i < numTetrahedra; ++i)
+		{
+			INT32 neighbors[4];
+			memcpy(neighbors, volume.tetrahedra[i].neighbors, sizeof(INT32) * 4);
+
+			for(UINT32 j = 0; j < 4; ++j)
+			{
+				INT32 vert = volume.tetrahedra[i].vertices[j];
+
+				for (UINT32 k = 0; k < 4; ++k)
+				{
+					INT32 neighborIdx = neighbors[k];
+					if (neighborIdx == -1)
+						continue;
+
+					Tetrahedron& neighbor = volume.tetrahedra[neighborIdx];
+					if (vert != neighbor.vertices[0] && vert != neighbor.vertices[1] &&
+						vert != neighbor.vertices[2] && vert != neighbor.vertices[3])
+					{
+						volume.tetrahedra[i].neighbors[j] = neighborIdx;
+						break;
+					}
+				}
+			}
+		}
+
 		return volume;
 	}
 }

Source/RenderBeast/CMakeSources.cmake

@@ -16,6 +16,7 @@ set(BS_RENDERBEAST_INC_NOFILTER
 	"Include/BsShadowRendering.h"
 	"Include/BsRendererScene.h"
 	"Include/BsStandardDeferredLighting.h"
+	"Include/BsLightProbes.h"
 )
 
 set(BS_RENDERBEAST_SRC_NOFILTER
@@ -35,6 +36,7 @@ set(BS_RENDERBEAST_SRC_NOFILTER
 	"Source/BsShadowRendering.cpp"
 	"Source/BsRendererScene.cpp"
 	"Source/BsStandardDeferredLighting.cpp"
+	"Source/BsLightProbes.cpp"
 )
 
 source_group("Header Files" FILES ${BS_RENDERBEAST_INC_NOFILTER})

Source/RenderBeast/Include/BsLightProbes.h

@@ -0,0 +1,29 @@
+//********************************** Banshee Engine (www.banshee3d.com) **************************************************//
+//**************** Copyright (c) 2016 Marko Pintera ([email protected]). All rights reserved. **********************//
+#pragma once
+
+#include "BsRenderBeastPrerequisites.h"
+#include "BsTriangulation.h"
+#include "BsMatrix4.h"
+
+namespace bs { namespace ct
+{
+	/** @addtogroup RenderBeast
+	 *  @{
+	 */
+
+	/** Handles any pre-processing for light (irradiance) probe lighting. */
+	class LightProbes
+	{
+		struct LightTetrahedron
+		{
+			Tetrahedron volume;
+			Matrix4 transform;
+		};
+
+		// TODO - This should accept LightProbe types as input probably
+		void updateProbes(const Vector<Vector3>& positions);
+	};
+
+	/** @} */
+}}

Source/RenderBeast/Source/BsLightProbes.cpp

@@ -0,0 +1,365 @@
+//********************************** Banshee Engine (www.banshee3d.com) **************************************************//
+//**************** Copyright (c) 2016 Marko Pintera ([email protected]). All rights reserved. **********************//
+#include "BsLightProbes.h"
+#include "../../External/XShaderCompiler/inc/Xsc/Reflection.h"
+
+namespace bs { namespace ct 
+{
+	/** Hash value generator for std::pair<INT32, INT32>. */
+	struct pair_hash
+	{
+		size_t operator()(const std::pair<INT32, INT32>& key) const
+		{
+			size_t hash = 0;
+			bs::hash_combine(hash, key.first);
+			bs::hash_combine(hash, key.second);
+
+			return hash;
+		}
+	};
+
+	void LightProbes::updateProbes(const Vector<Vector3>& positions)
+	{
+		bs_frame_mark();
+		{
+			TetrahedronVolume volume = Triangulation::tetrahedralize(positions);
+
+			// Put outer faces into the Tetrahedron structure, for convenience
+			UINT32 outerFaceOffset = (UINT32)volume.tetrahedra.size();
+			FrameVector<Tetrahedron> outerTetrahedrons(volume.outerFaces.size());
+
+			for (UINT32 i = 0; i < (UINT32)volume.outerFaces.size(); ++i)
+			{
+				Tetrahedron outerTetrahedron;
+				memcpy(outerTetrahedron.vertices, volume.outerFaces[i].vertices, sizeof(INT32) * 3);
+				memset(outerTetrahedron.neighbors, -1, sizeof(INT32) * 3);
+
+				outerTetrahedron.vertices[4] = -1; // Marks the tetrahedron as an outer face
+				outerTetrahedron.neighbors[4] = volume.outerFaces[i].tetrahedron;
+
+				outerTetrahedrons[i] = outerTetrahedron;
+			}
+
+			// Connect boundary tetrahedrons with these new outer tetrahedrons
+			for (UINT32 i = 0; i < (UINT32)volume.outerFaces.size(); ++i)
+			{
+				Tetrahedron& tet = volume.tetrahedra[volume.outerFaces[i].tetrahedron];
+				for (UINT32 j = 0; j < 4; j++)
+				{
+					if (tet.neighbors[j] == -1)
+						tet.neighbors[j] = outerFaceOffset + i;
+				}
+			}
+
+			// Make a map between outer edges and faces, used in the following algorithms
+			struct Edge
+			{
+				INT32 faces[2];
+				INT32 oppositeVerts[2];
+			};
+
+			FrameUnorderedMap<std::pair<INT32, INT32>, Edge, pair_hash> edgeMap;
+			for (UINT32 i = 0; i < (UINT32)volume.outerFaces.size(); ++i)
+			{
+				for (UINT32 j = 0; j < 3; ++j)
+				{
+					INT32 v0 = volume.outerFaces[i].vertices[j];
+					INT32 v1 = volume.outerFaces[i].vertices[(j + 1) % 3];
+
+					// Keep the same ordering so other faces can find the same edge
+					if (v0 > v1)
+						std::swap(v0, v1);
+
+					auto iterFind = edgeMap.find(std::make_pair(v0, v1));
+					if (iterFind != edgeMap.end())
+					{
+						iterFind->second.faces[1] = i;
+						iterFind->second.oppositeVerts[1] = (j + 2) % 3;
+					}
+					else
+					{
+						Edge edge;
+						edge.faces[0] = i;
+						edge.oppositeVerts[0] = (j + 2) % 3;
+
+						edgeMap.insert(std::make_pair(std::make_pair(v0, v1), edge));
+					}
+				}
+			}
+
+			// Form connections between outer tetrahedrons
+			for (auto& entry : edgeMap)
+			{
+				const Edge& edge = entry.second;
+
+				Tetrahedron& tet0 = outerTetrahedrons[outerFaceOffset + edge.faces[0]];
+				tet0.neighbors[edge.oppositeVerts[0]] = outerFaceOffset + edge.faces[1];
+
+				Tetrahedron& tet1 = outerTetrahedrons[outerFaceOffset + edge.faces[1]];
+				tet1.neighbors[edge.oppositeVerts[1]] = outerFaceOffset + edge.faces[0];
+			}
+
+			// Generate face normals
+			FrameVector<Vector3> faceNormals(volume.outerFaces.size());
+			for (UINT32 i = 0; i < (UINT32)volume.outerFaces.size(); ++i)
+			{
+				const Vector3& v0 = positions[volume.outerFaces[i].vertices[0]];
+				const Vector3& v1 = positions[volume.outerFaces[i].vertices[1]];
+				const Vector3& v2 = positions[volume.outerFaces[i].vertices[2]];
+				
+				Vector3 e0 = v1 - v0;
+				Vector3 e1 = v2 - v0;
+
+				faceNormals[i] = Vector3::normalize(e1.cross(e0));
+			}
+
+			// Generate vertex normals
+			struct VertexAccum
+			{
+				Vector3 sum;
+				float weight;
+			};
+
+			FrameUnorderedMap<INT32, Vector3> vertexNormals;
+			for (auto& entry : edgeMap)
+			{
+				const Edge& edge = entry.second;
+
+				auto accumulateNormalForEdgeVertex = [&](UINT32 v0Idx, UINT32 v1Idx)
+				{
+					auto iter = vertexNormals.insert(std::make_pair(v0Idx, Vector3(BsZero)));
+
+					Vector3& accum = iter.first->second;
+					const Vector3& v0 = positions[v0Idx];
+
+					auto accumulateNormalForFace = [&](INT32 faceIdx, INT32 v2LocIdx)
+					{
+						const TetrahedronFace& face = volume.outerFaces[faceIdx];
+
+						// Vertices on the face, that aren't the vertex we're calculating the normal for
+						const Vector3& v1 = positions[v1Idx];
+						const Vector3& v2 = positions[face.vertices[v2LocIdx]];
+
+						// Weight the contribution to the normal based on the angle spanned by the triangle
+						Vector3 e0 = Vector3::normalize(v1 - v0);
+						Vector3 e1 = Vector3::normalize(v2 - v0);
+
+						float weight = acos(e0.dot(e1));
+						accum += weight * faceNormals[faceIdx];
+					};
+
+					accumulateNormalForFace(edge.faces[0], entry.second.oppositeVerts[0]);
+					accumulateNormalForFace(edge.faces[1], entry.second.oppositeVerts[1]);
+				};
+
+				accumulateNormalForEdgeVertex(entry.first.first, entry.first.second);
+				accumulateNormalForEdgeVertex(entry.first.second, entry.first.first);
+			}
+
+			for (auto& entry : vertexNormals)
+				entry.second.normalize();
+
+			// Generate matrices
+			Vector<LightTetrahedron> output;
+			output.reserve(volume.tetrahedra.size() + volume.outerFaces.size());
+
+			// Insert innert tetrahedrons, generate matrices
+			for(UINT32 i = 0; i < (UINT32)volume.tetrahedra.size(); ++i)
+			{
+				LightTetrahedron entry;
+				entry.volume = volume.tetrahedra[i];
+
+				// Generate a matrix that can be used for calculating barycentric coordinates
+				// To determine a point within a tetrahedron, using barycentric coordinates, we use:
+				// P = (P1 - P4) * a + (P2 - P4) * b + (P3 - P4) * c + P4
+				//
+				// Where P1, P2, P3, P4 are the corners of the tetrahedron.
+				//
+				// Expanded for each coordinate this is:
+				// x = (x1 - x4) * a + (x2 - x4) * b + (x3 - x4) * c + x4
+				// y = (y1 - y4) * a + (y2 - y4) * b + (y3 - y4) * c + y4
+				// z = (z1 - z4) * a + (z2 - z4) * b + (z3 - z4) * c + z4
+				//
+				// In matrix form this is:
+				//                                      a
+				// P = [P1 - P4, P2 - P4, P3 - P4, P4] [b]
+				//                                      c
+				//                                      1
+				//
+				// Solved for barycentric coordinates:
+				//  a
+				// [b] = Minv * P 
+				//  c
+				//  1
+				//
+				// Where Minv is the inverse of the matrix above.
+
+				const Vector3& P1 = positions[volume.tetrahedra[i].vertices[0]];
+				const Vector3& P2 = positions[volume.tetrahedra[i].vertices[1]];
+				const Vector3& P3 = positions[volume.tetrahedra[i].vertices[2]];
+				const Vector3& P4 = positions[volume.tetrahedra[i].vertices[3]];
+
+				Matrix4 mat;
+				mat.setColumn(0, Vector4(P1 - P4, 0.0f));
+				mat.setColumn(1, Vector4(P2 - P4, 0.0f));
+				mat.setColumn(2, Vector4(P3 - P4, 0.0f));
+				mat.setColumn(3, Vector4(P4, 1.0f));
+
+				entry.transform = mat.inverse();
+
+				output.push_back(entry);
+			}
+
+			// Insert outer tetrahedrons, generate matrices
+			for(UINT32 i = 0; i < (UINT32)outerTetrahedrons.size(); ++i)
+			{
+				LightTetrahedron entry;
+				entry.volume = outerTetrahedrons[i];
+
+				// We need a way to project a point outside the tetrahedron volume onto an outer face, then calculate
+				// triangle's barycentric coordinates. Use use the per-vertex normals to extrude the triangle face into
+				// infinity.
+
+				// Our point can be represented as:
+				// p == a (p0 + t*v0) + b (p1 + t*v1) + c (p2 + t*v2)
+				//
+				// where a, b and c are barycentric coordinates,
+				// p0, p1, p2 are the corners of the face
+				// v0, v1, v2 are the vertex normals, per corner
+				// t is the distance from the triangle to the point
+				//
+				// Essentially we're calculating the corners of a bigger triangle that's "t" units away from the
+				// face, and its corners lie along the per-vertex normals. Point "p" will lie on that triangle, for which
+				// we can then calculate barycentric coordinates normally.
+				//
+				// First we substitute: c = 1 - a - b
+				// p == a (p0 + t v0) + b (p1 + t v1) + (1 - a - b) (p2 + t v2)
+				// p == a (p0 + t v0) + b (p1 + t v1) + (p2 + t v2) - a (p2 + t v2) - b (p2 + t v2)
+				// p == a (p0 - p2 + t v0 - t v2) + b (p1 - p2 + t v1 - t v2) + (p2 + t v2)
+				//
+				// And move everything to one side:
+				// p - p2 - t v2 == a (p0 - p2 + t ( v0 - v2)) + b (p1 - p2 + t ( v1 - v2))
+				// a (p0 - p2 + t ( v0 - v2)) + b (p1 - p2 + t ( v1 - v2)) - (p - p2 - t v2) == 0
+				//
+				// We rewrite it using:
+				// Ap = p0 - p2
+				// Av = v0 - v2
+				// Bp = p1 - p2
+				// Bv = v1 - v2
+				// Cp = p - p2
+				// Cv = -v2
+				//
+				// Which yields:
+				// a (Ap + t Av) + b (Bp + t Bv) - (Cp + t Cv) == 0
+				//
+				// Which can be written in matrix form:
+				//
+				// M = {Ap + t Av, Bp + t Bv, Cp + t Cv}
+				//       a      0
+				// M * [ b ] = [0]
+				//      -1      0
+				//
+				// From that we can tell that matrix M cannot be inverted, because if we multiply the zero vector with the
+				// inverted matrix the result would be zero, and not [a, b, -1]. Since the matrix cannot be inverted
+				// det(M) == 0.
+				//
+				// We can use that fact to calculate "t". After we have "t" we can calculate barycentric coordinates
+				// normally.
+				//
+				// Solving equation det(M) == 0 yields a cubic in form:
+				// p t^3 + q t^2 + r t + s = 0
+				//
+				// We'll convert this to monic form, by dividing by p:
+				// t^3 + q/p t^2 + r/p t + s/p = 0
+				//
+				// Or if p ends up being zero, we end up with a quadratic instead:
+				// q t^2 + r t + s = 0
+				// 
+				// We want to create a matrix that when multiplied with the position, yields us the three coefficients,
+				// which we can then use to solve for "t". For this we create a 4x3 matrix, where each row represents
+				// a solution for one of the coefficients. We factor contributons to each coefficient whether they depend on
+				// position x, y, z, or don't depend on position (row columns, in that order respectively).
+
+				const Vector3& p0 = positions[entry.volume.vertices[0]];
+				const Vector3& p1 = positions[entry.volume.vertices[1]];
+				const Vector3& p2 = positions[entry.volume.vertices[2]];
+
+				const Vector3& v0 = vertexNormals[entry.volume.vertices[0]];
+				const Vector3& v1 = vertexNormals[entry.volume.vertices[1]];
+				const Vector3& v2 = vertexNormals[entry.volume.vertices[2]];
+
+				float p =
+					v2.x * v1.y * v0.z -
+					v1.x * v2.y * v0.z -
+					v2.x * v0.y * v1.z +
+					v0.x * v2.y * v1.z +
+					v1.x * v0.y * v2.z -
+					v0.x * v1.y * v2.z;
+				
+				float qx = -v1.y * v0.z + v2.y * v0.z + v0.y * v1.z - v2.y * v1.z - v0.y * v2.z + v1.y * v2.z;
+				float qy = v1.x * v0.z - v2.x * v0.z - v0.x * v1.z + v2.x * v1.z + v0.x * v2.z - v1.x * v2.z;
+				float qz = -v1.x * v0.y + v2.x * v0.y + v0.x * v1.y - v2.x * v1.y - v0.x * v2.y + v1.x * v2.y;
+				float qw = v2.y * v1.z * p0.x - v1.y * v2.z * p0.x - v2.y * v0.z * p1.x + v0.y * v2.z * p1.x + 
+					v1.y * v0.z * p2.x - v0.y * v1.z * p2.x - v2.x * v1.z * p0.y + v1.x * v2.z * p0.y + 
+					v2.x * v0.z * p1.y - v0.x * v2.z * p1.y - v1.x * v0.z * p2.y + v0.x * v1.z * p2.y + 
+					v2.x * v1.y * p0.z - v1.x * v2.y * p0.z - v2.x * v0.y * p1.z + v0.x * v2.y * p1.z + 
+					v1.x * v0.y * p2.z - v0.x * v1.y * p2.z;
+
+				float rx = v1.z * p0.y - v2.z * p0.y - v0.z * p1.y + v2.z * p1.y + v0.z * p2.y - v1.z * p2.y -
+					v1.y * p0.z + v2.y * p0.z + v0.y * p1.z - v2.y * p1.z - v0.y * p2.z + v1.y * p2.z;
+				float ry = -v1.z * p0.x + v2.z * p0.x + v0.z * p1.x - v2.z * p1.x - v0.z * p2.x + v1.z * p2.x +
+					v1.x * p0.z - v2.x * p0.z - v0.x * p1.z + v2.x * p1.z + v0.x * p2.z - v1.x * p2.z;
+				float rz = v1.y * p0.x - v2.y * p0.x - v0.y * p1.x + v2.y * p1.x + v0.y * p2.x - v1.y * p2.x -
+					v1.x * p0.y + v2.x * p0.y + v0.x * p1.y - v2.x * p1.y - v0.x * p2.y + v1.x * p2.y;
+				float rw = v2.z * p1.x * p0.y - v1.z * p2.x * p0.y - v2.z * p0.x * p1.y + v0.z * p2.x * p1.y +
+					v1.z * p0.x * p2.y - v0.z * p1.x * p2.y - v2.y * p1.x * p0.z + v1.y * p2.x * p0.z +
+					v2.x * p1.y * p0.z - v1.x * p2.y * p0.z + v2.y * p0.x * p1.z - v0.y * p2.x * p1.z -
+					v2.x * p0.y * p1.z + v0.x * p2.y * p1.z - v1.y * p0.x * p2.z + v0.y * p1.x * p2.z +
+					v1.x * p0.y * p2.z - v0.x * p1.y * p2.z;
+
+				float sx = -p1.y * p0.z + p2.y * p0.z + p0.y * p1.z - p2.y * p1.z - p0.y * p2.z + p1.y * p2.z;
+				float sy = p1.x * p0.z - p2.x * p0.z - p0.x * p1.z + p2.x * p1.z + p0.x * p2.z - p1.x * p2.z;
+				float sz = -p1.x * p0.y + p2.x * p0.y + p0.x * p1.y - p2.x * p1.y - p0.x * p2.y + p1.x * p2.y;
+				float sw = p2.x * p1.y * p0.z - p1.x * p2.y * p0.z - p2.x * p0.y * p1.z + 
+					p0.x * p2.y * p1.z + p1.x * p0.y * p2.z - p0.x * p1.y * p2.z;
+
+				entry.transform[0][0] = qx;
+				entry.transform[0][1] = qy;
+				entry.transform[0][2] = qz;
+				entry.transform[0][3] = qw;
+
+				entry.transform[1][0] = rx;
+				entry.transform[1][1] = ry;
+				entry.transform[1][2] = rz;
+				entry.transform[1][3] = rw;
+
+				entry.transform[2][0] = sx;
+				entry.transform[2][1] = sy;
+				entry.transform[2][2] = sz;
+				entry.transform[2][3] = sw;
+
+				// Unused
+				entry.transform[3][0] = 0.0f;
+				entry.transform[3][1] = 0.0f;
+				entry.transform[3][2] = 0.0f;
+				entry.transform[3][3] = 0.0f;
+
+				if (fabs(p) > 0.00001f)
+					entry.transform = entry.transform * (1.0f / p);
+				else // Quadratic
+					entry.volume.neighbors[3] = -2;
+
+				output.push_back(entry);
+			}
+		}
+		bs_frame_clear();
+	}
+}}
+
+/** @cond STDLIB */
+
+namespace std
+{
+}
+
+/** @endcond */