diff --git a/apps/specular_estimation/CMakeFiles/specular_estimation.dir/src/specular_estimation.cc.o b/apps/specular_estimation/CMakeFiles/specular_estimation.dir/src/specular_estimation.cc.o
index 59809c72e6c2d7826035549e2235d96de11d102c..b6aa8bee2599dffe58502ea20650860108c1b484 100644
Binary files a/apps/specular_estimation/CMakeFiles/specular_estimation.dir/src/specular_estimation.cc.o and b/apps/specular_estimation/CMakeFiles/specular_estimation.dir/src/specular_estimation.cc.o differ
diff --git a/apps/specular_estimation/src/CookTorrance.fragmentshader b/apps/specular_estimation/src/CookTorrance.fragmentshader
new file mode 100644
index 0000000000000000000000000000000000000000..ecfae888279923a0c70f65dddd5ab80dfdd2cb2e
--- /dev/null
+++ b/apps/specular_estimation/src/CookTorrance.fragmentshader
@@ -0,0 +1,188 @@
+#version 330 core
+// https://stackoverflow.com/a/39684336
+
+// Interpolated values from the vertex shaders
+in vec2 UV;
+in vec3 Normal_cameraspace;
+in vec3 EyeDirection_cameraspace;
+in vec3 LightDirection_cameraspace;
+in vec4 ShadowCoord;
+//in vec3 Position_worldspace;
+
+// Ouput data
+layout(location = 0) out vec3 color;
+
+// Values that stay constant for the whole mesh.
+uniform sampler2D myTextureSampler;
+//uniform sampler2D mySpecularSampler;
+uniform mat4 MV;
+uniform sampler2DShadow shadowMap;
+//uniform vec3 LightPosition_worldspace;
+uniform float specularIntensity;
+uniform float specularPower;
+
+vec2 poissonDisk[16] = vec2[]( 
+   vec2( -0.94201624,  -0.39906216 ), 
+   vec2(  0.94558609,  -0.76890725 ), 
+   vec2( -0.094184101, -0.92938870 ), 
+   vec2(  0.34495938,   0.29387760 ), 
+   vec2( -0.91588581,   0.45771432 ), 
+   vec2( -0.81544232,  -0.87912464 ), 
+   vec2( -0.38277543,   0.27676845 ), 
+   vec2(  0.97484398,   0.75648379 ), 
+   vec2(  0.44323325,  -0.97511554 ), 
+   vec2(  0.53742981,  -0.47373420 ), 
+   vec2( -0.26496911,  -0.41893023 ), 
+   vec2(  0.79197514,   0.19090188 ), 
+   vec2( -0.24188840,   0.99706507 ), 
+   vec2( -0.81409955,   0.91437590 ), 
+   vec2(  0.19984126,   0.78641367 ), 
+   vec2(  0.14383161,  -0.14100790 ) 
+);
+
+// Returns a random number based on a vec3 and an int.
+float random(vec3 seed, int i) {
+	vec4 seed4 = vec4(seed, i);
+	float dot_product = dot(seed4, vec4(12.9898,78.233,45.164,94.673));
+	return fract(sin(dot_product) * 43758.5453);
+}
+
+
+
+float D(float m, float c) {
+    float c2 = c*c, m2 = m*m, c2m2 = c2*m2;
+    return exp((c2 - 1)/c2m2)/(3.14159*c2m2*c2);
+}
+float F(float R0, float NV) { return R0 + (1 - R0)*pow(1 - NV, 5); }
+float G(float NL, float NV, float NH, float HV) { return min(1, 2*NH*min(NV, NL)/HV); }
+
+
+void accumulate_light(special3 tangent, vec3 viewDir, float roughness, PerLight light, inout vec3 diffuse, inout vec3 specular) {
+    vec3 lightDir = quat_apply(tangent.q, light.position);
+
+    if(lightDir.z > 0)
+    {
+        lightDir = normalize(lightDir);
+        float NL = lightDir.z;
+        diffuse += NL * light.color;
+
+        float NV = viewDir.z;
+        if(NV > 0)
+        {
+            vec3 halfDir = normalize(lightDir + viewDir);
+            float NH = halfDir.z;
+            float HV = dot(halfDir, viewDir);
+            specular += D(roughness, NH)*F(0.034, NV)*G(NL, NV, NH, HV)/(4*NV*NL) * light.color;
+        }
+    }
+}
+
+void main() {
+    special3 tangent = {
+        vec3(0),
+        texture(u_bump, IN.UV.xy)
+    };
+
+    tangent = special_mul(IN.pos, tangent);
+    tangent = special_inverse(tangent);
+
+    const vec3 viewDir = normalize(tangent.v);
+
+    vec3 diffuse = texture( myTextureSampler, UV ).rgb;
+    vec3 specular = vec3(0);
+    float roughness = max(0.215 + texture(u_roughness, IN.IV.xy).r - .5, 0.001);
+    for(int i = 0; i < LIGHTS.nlights; ++i)
+        accumulate_light(tangent, viewDir, roughness, LIGHTS.lights[i], diffuse, specular);
+
+    OUT = vec4(diffuse*IN.diffuse.rgb + specular, IN.diffuse.a);
+}
+
+
+/*
+void main() {
+
+	// Light emission properties
+	vec3 LightColor = vec3(1, 1, 1);
+
+	// The power and colour of each light could be obtained from either the chrome sphere reflection or the Macbeth colour chart
+	float LightPower = 1.0f;
+	
+	// Material properties
+	float ambientPower = 0.0f;
+	vec3 MaterialDiffuseColor  = texture( myTextureSampler, UV ).rgb;
+	vec3 MaterialAmbientColor  = vec3(ambientPower, ambientPower, ambientPower) * MaterialDiffuseColor;
+	
+	vec3 MaterialSpecularColor = vec3(specularIntensity, specularIntensity, specularIntensity);
+	//vec3 MaterialSpecularColor = texture( myTextureSampler, UV ).rgb;
+
+	// Distance to the light
+	//float distance = length( LightPosition_worldspace - Position_worldspace );
+
+	// Normal of the computed fragment, in camera space
+	vec3 n = normalize( Normal_cameraspace );
+	// Direction of the light (from the fragment to the light)
+	vec3 l = normalize( LightDirection_cameraspace );
+	// Cosine of the angle between the normal and the light direction, 
+	// clamped above 0
+	//  - light is at the vertical of the triangle -> 1
+	//  - light is perpendiular to the triangle -> 0
+	//  - light is behind the triangle -> 0
+	float cosTheta = clamp( dot( n, l ), 0, 1 );
+	
+	// Eye vector (towards the camera)
+	vec3 E = normalize(EyeDirection_cameraspace);
+	// Direction in which the triangle reflects the light
+	vec3 R = reflect(-l, n);
+	// Cosine of the angle between the Eye vector and the Reflect vector, clamped to 0
+	//  - Looking into the reflection -> 1
+	//  - Looking elsewhere -> <1
+	float cosAlpha = clamp( dot( E, R ), 0, 1 );
+	
+	float visibility = 1.0f;
+
+	// Fixed bias
+	//float bias = 0.005;
+
+	// Variable bias, which reduces shadow acne by modifying the bias according to the slope.
+	float bias = 0.005*tan(acos(cosTheta));
+	bias = clamp(bias, 0, 0.01);
+
+	int samples = 16;
+	float darkness = 1.0 - ambientPower;
+
+	// Sample the shadow map 16 times
+	for (int i = 0; i < samples; i++) {
+		// use either :
+		//  - Always the same samples. Gives a fixed pattern in the shadow, but no noise
+		// int index = i;
+		//  - A random sample, based on the pixel's screen location. No banding, but the shadow moves with the camera, causing shimmering.
+		// int index = int(16.0*random(gl_FragCoord.xyy, i))%16;
+		//  - A random sample, based on the pixel's position in world space. The position is rounded to the millimetre to avoid too much aliasing
+		// int index = int(16.0 * random(floor(Position_worldspace.xyz * 1000.0), i)) % 16;
+		
+		// being fully in the shadow will eat up 4*0.2 = 0.8
+		// 0.2 potentially remain, which is quite dark.
+		//visibility -= 0.2*(1.0-texture( shadowMap, vec3(ShadowCoord.xy + poissonDisk[index]/700.0,  (ShadowCoord.z-bias)/ShadowCoord.w) ));
+
+		// Orthographic Shadow Projection
+		//visibility -= (darkness/samples)*(1.0-texture( shadowMap, vec3(ShadowCoord.xy + poissonDisk[i]/700.0,  (ShadowCoord.z-bias)/ShadowCoord.w) ));
+
+		// Perspective Shadow Projection
+		// texture( shadowMap, ShadowCoord.xy ).z is the distance between the light and the nearest occluder
+		// ShadowCoord.z is the distance between the light and the current fragment
+		visibility -= (darkness/samples)*(1.0-texture( shadowMap, vec3(ShadowCoord.xy/ShadowCoord.w + poissonDisk[i]/700.0,  (ShadowCoord.z-bias)/ShadowCoord.w) ));
+	}
+
+	// For spot lights, use either one of these lines instead.
+	// if ( texture( shadowMap, (ShadowCoord.xy/ShadowCoord.w) ).z  <  (ShadowCoord.z-bias)/ShadowCoord.w )
+	// if ( textureProj( shadowMap, ShadowCoord.xyw ).z  <  (ShadowCoord.z-bias)/ShadowCoord.w )
+	
+	color = 
+		// Ambient: simulates indirect lighting
+		MaterialAmbientColor +
+		// Diffuse: "colour" of the object
+		visibility * MaterialDiffuseColor  * LightColor * LightPower * cosTheta +
+		// Specular: reflective highlight, like a mirror
+		visibility * MaterialSpecularColor * LightColor * LightPower * pow(cosAlpha, specularPower);
+}
+*/
\ No newline at end of file
diff --git a/apps/specular_estimation/src/OpenGL.h b/apps/specular_estimation/src/OpenGL.h
index 03ad8add82de58016fe0efe427e663d00dcfe740..799d5b60120ff818c959e083a47347fce88fba50 100644
--- a/apps/specular_estimation/src/OpenGL.h
+++ b/apps/specular_estimation/src/OpenGL.h
@@ -52,7 +52,7 @@ void renderPolygons(int windowWidth, int windowHeight, GLuint programID, glm::ve
 void renderShadows(GLuint FramebufferName, int shadowResolution, GLuint depthProgramID, glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, GLuint depthMatrixID, GLuint vertexbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, glm::mat4 depthModelMatrix, glm::mat4 depthMVP);
 void renderShadowMap(GLuint quad_programID, GLuint depthTexture, GLuint texID, GLuint quad_vertexbuffer);
 //void meanSquaredError(int imageNumber, int windowHeight, int windowWidth, std::string modelPath, std::vector < cv::Mat > textureImages, cv::Mat& residual, cv::Vec3d& residuals, double& sum, double SpecularIntensity, double SpecularPower);
-void meanSquaredError(int imageNumber, int windowHeight, int windowWidth, std::vector <cv::Mat> textureImages, double& sum);
+void meanSquaredError(int imageNumber, int windowHeight, int windowWidth, std::vector<cv::Mat> textureImages, double& sum);
 double returnResidual(int imageNumber, int windowHeight, int windowWidth, std::vector<cv::Mat> textureImages, cv::Mat& residual, cv::Vec3d& residuals, double SpecularIntensity, double SpecularPower);
 //void computeResiduals(cv::Vec3d& residual, std::vector< cv::Vec3d >& residuals, double& totalResidual, cv::Mat& residualImage, bool calculateResidual, glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat> textureImages, int lightNumber, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath);
 void viewModel(cv::Vec3d& residual, double& totalResidual, cv::Mat& residualImage, bool calculateResidual, glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat> textureImages, int lightNumber, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath);
@@ -60,6 +60,12 @@ void RenderSynthetic(glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix,
 void checkForErrors(std::string stage, bool printIfOK);
 void calibrateLights(int& windowWidth, int& windowHeight, cv::Mat cv_texture);
 
+void initialiseOpenGLCT(cv::Mat cv_depth, cv::Mat cv_normals, cv::Mat cv_texture, cv::Mat lightDirections, int& windowWidth, int& windowHeight, glm::mat4& depthProjectionMatrix, glm::mat4& depthViewMatrix, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, GLuint& programID, GLuint& MatrixID, GLuint& ModelMatrixID, GLuint& ViewMatrixID, GLuint& DepthBiasID, GLuint& lightInvDirID, GLuint& Texture, GLuint& TextureID, GLuint& depthTexture, GLuint& ShadowMapID, GLuint& vertexbuffer, GLuint& uvbuffer, GLuint& normalbuffer, GLuint& elementbuffer, std::vector<unsigned int> & indices, GLuint& depthProgramID, GLuint& quad_programID, GLuint& FramebufferName, GLuint& quad_vertexbuffer, GLuint& VertexArrayID, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, int numberOfLights, bool calculateResidual, glm::mat4& depthMVP, glm::mat4& depthModelMatrix, glm::mat4& MVP, glm::mat4& ProjectionMatrix, glm::mat4& ViewMatrix, glm::mat4& ModelMatrix, glm::mat4& depthBiasMVP, glm::mat4& biasMatrix);
+void renderPolygonsCT(int windowWidth, int windowHeight, GLuint programID, glm::vec3 lightInvDir, GLuint MatrixID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, glm::mat4 MVP, glm::mat4 ModelMatrix, glm::mat4 ViewMatrix, glm::mat4 depthBiasMVP, GLuint SpecularIntensityID, double SpecularIntensity, GLuint SpecularPowerID, double SpecularPower);
+void terminateOpenGLCT(GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, GLuint programID, GLuint depthProgramID, GLuint quad_programID, GLuint Texture, GLuint FramebufferName, GLuint depthTexture, GLuint quad_vertexbuffer, GLuint VertexArrayID);
+void viewModelCT(cv::Vec3d& residual, double& totalResidual, cv::Mat& residualImage, bool calculateResidual, glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat> textureImages, int lightNumber, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath);
+void RenderSyntheticCT(glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat>& textureImages, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath);
+
 void initialiseOpenGL(cv::Mat cv_depth, cv::Mat cv_normals, cv::Mat cv_texture, cv::Mat lightDirections, int& windowWidth, int& windowHeight, glm::mat4& depthProjectionMatrix, glm::mat4& depthViewMatrix, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, GLuint& programID, GLuint& MatrixID, GLuint& ModelMatrixID, GLuint& ViewMatrixID, GLuint& DepthBiasID, GLuint& lightInvDirID, GLuint& Texture, GLuint& TextureID, GLuint& depthTexture, GLuint& ShadowMapID, GLuint& vertexbuffer, GLuint& uvbuffer, GLuint& normalbuffer, GLuint& elementbuffer, std::vector<unsigned int> & indices, GLuint& depthProgramID, GLuint& quad_programID, GLuint& FramebufferName, GLuint& quad_vertexbuffer, GLuint& VertexArrayID, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, int numberOfLights, bool calculateResidual, glm::mat4& depthMVP, glm::mat4& depthModelMatrix, glm::mat4& MVP, glm::mat4& ProjectionMatrix, glm::mat4& ViewMatrix, glm::mat4& ModelMatrix, glm::mat4& depthBiasMVP, glm::mat4& biasMatrix) {
 	// Initialise the GLFW library first
 	if (!glfwInit()) {
@@ -1508,6 +1514,24 @@ void calibrateLights(int& windowWidth, int& windowHeight, cv::Mat cv_texture) {
 	//depthProjectionMatrix = glm::perspective<float>(45.0f, 1.0f, 1.0f, 100.0f);
 	//depthViewMatrix = glm::lookAt(glm::vec3(0, 0, 1), glm::vec3(0, 0, 0), glm::vec3(0, 1, 0));
 
+	// Initial position: on +Z
+	glm::vec3 position = glm::vec3(0, 0, 1);
+	// Initial horizontal angle: toward -Z
+	float horizontalAngle = M_PI;
+	// Initial vertical angle: none
+	float verticalAngle = 0.0f;
+	// Initial Field of View
+	float FoV = 60.0f;
+	// Mouse Speed
+	float mouseSpeed = 0.005f;
+		
+	bool perspectiveProjection = false, shadowControl = false, calculateResidual = false;
+	int lightNumber = 0;
+
+	glm::mat4 ViewMatrix, ProjectionMatrix;
+
+	double SpecularIntensity = 0.5, SpecularPower = 2;
+
 	do {
 		// Render to our framebuffer
 		glBindFramebuffer(GL_FRAMEBUFFER, FramebufferName);
@@ -1582,23 +1606,7 @@ void calibrateLights(int& windowWidth, int& windowHeight, cv::Mat cv_texture) {
 		// Use our shader
 		glUseProgram(programID);
 
-		// Initial position: on +Z
-		glm::vec3 position = glm::vec3(0, 0, 1);
-		// Initial horizontal angle: toward -Z
-		float horizontalAngle = M_PI;
-		// Initial vertical angle: none
-		float verticalAngle = 0.0f;
-		// Initial Field of View
-		float FoV = 60.0f;
-		// Mouse Speed
-		float mouseSpeed = 0.005f;
 		
-		bool perspectiveProjection = false, shadowControl = false, calculateResidual = false;
-		int lightNumber = 0;
-
-		glm::mat4 ViewMatrix, ProjectionMatrix;
-
-		double SpecularIntensity = 0.5, SpecularPower = 2;
 
 		// Compute the MVP matrix from keyboard and mouse input
 		computeMatricesFromInputs(windowWidth, windowHeight, position, horizontalAngle, verticalAngle, FoV, mouseSpeed, ProjectionMatrix, ViewMatrix, lightInvDir, depthProjectionMatrix, depthViewMatrix, perspectiveProjection, shadowControl, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower, calculateResidual);
@@ -1729,8 +1737,7 @@ void calibrateLights(int& windowWidth, int& windowHeight, cv::Mat cv_texture) {
 		glfwPollEvents();
 
 	} // Check if the ESC key was pressed or the window was closed
-	while( glfwGetKey(window, GLFW_KEY_ESCAPE ) != GLFW_PRESS &&
-		   glfwWindowShouldClose(window) == 0 );
+	while( glfwGetKey(window, GLFW_KEY_ESCAPE ) != GLFW_PRESS && glfwWindowShouldClose(window) == 0 );
 
 	// Cleanup VBO and shader
 	glDeleteBuffers(1, &vertexbuffer);
@@ -2513,4 +2520,454 @@ class OpenGL {
 	return fin_pt;
 }*/
 
+
+void initialiseOpenGLCT(cv::Mat cv_depth, cv::Mat cv_normals, cv::Mat cv_texture, cv::Mat lightDirections, int& windowWidth, int& windowHeight, glm::mat4& depthProjectionMatrix, glm::mat4& depthViewMatrix, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, GLuint& programID, GLuint& MatrixID, GLuint& ModelMatrixID, GLuint& ViewMatrixID, GLuint& DepthBiasID, GLuint& lightInvDirID, GLuint& Texture, GLuint& TextureID, GLuint& depthTexture, GLuint& ShadowMapID, GLuint& vertexbuffer, GLuint& uvbuffer, GLuint& normalbuffer, GLuint& elementbuffer, std::vector<unsigned int> & indices, GLuint& depthProgramID, GLuint& quad_programID, GLuint& FramebufferName, GLuint& quad_vertexbuffer, GLuint& VertexArrayID, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, int numberOfLights, bool calculateResidual, glm::mat4& depthMVP, glm::mat4& depthModelMatrix, glm::mat4& MVP, glm::mat4& ProjectionMatrix, glm::mat4& ViewMatrix, glm::mat4& ModelMatrix, glm::mat4& depthBiasMVP, glm::mat4& biasMatrix) {
+	// Initialise the GLFW library first
+	if (!glfwInit()) {
+		std::cerr << "Failed to initialise GLFW.\n";
+		getchar();
+		return;
+	}
+	std::cout << "Initialised GLFW.\n";
+
+	glfwWindowHint(GLFW_SAMPLES, 4); // 4x antialiasing
+	glfwWindowHint(GLFW_RESIZABLE, GL_FALSE);
+	glfwWindowHint(GLFW_CONTEXT_VERSION_MAJOR, 3); // OpenGL 3.3
+	glfwWindowHint(GLFW_CONTEXT_VERSION_MINOR, 3);
+	glfwWindowHint(GLFW_OPENGL_FORWARD_COMPAT, GL_TRUE); // Only needed for macOS
+	glfwWindowHint(GLFW_OPENGL_PROFILE, GLFW_OPENGL_CORE_PROFILE); // Disable old OpenGL
+
+	// Set the shadow resolution to 4096x4096
+	int shadowResolution = 4096;
+
+	// Open a window and create its OpenGL context
+	window = glfwCreateWindow(windowWidth, windowHeight, "Window", NULL, NULL);
+	if (window == NULL) {
+		std::cerr << "Failed to open GLFW window. Older Intel GPUs are not OpenGL 3.3 compatible.\n";
+		getchar();
+		glfwTerminate();
+		return;
+	}
+	glfwMakeContextCurrent(window);
+	std::cout << "Opened GLFW window.\n";
+
+	// On macOS with a retina screen it will be windowWidth*2 and windowHeight*2, so the actual framebuffer size is:
+	glfwGetFramebufferSize(window, &windowWidth, &windowHeight);
+
+	// Initialize GLEW
+	glewExperimental = true; // Needed for core profile
+	if (glewInit() != GLEW_OK) {
+		std::cerr << "Failed to initialise GLEW.\n";
+		getchar();
+		glfwTerminate();
+		return;
+	}
+	std::cout << "Initialised GLEW.\n";
+
+	// Ensure that the escape key being pressed can be captured
+	glfwSetInputMode(window, GLFW_STICKY_KEYS, GL_TRUE);
+	// Hide the mouse and enable unlimited mouvement
+	//glfwSetInputMode(window, GLFW_CURSOR, GLFW_CURSOR_DISABLED);
+
+	// Set the mouse at the center of the screen
+	glfwPollEvents();
+	//glfwSetCursorPos(window, windowWidth / 2, windowHeight / 2);
+
+	// Dark blue background
+	glClearColor(0.0f, 0.0f, 0.4f, 0.0f);
+
+	// Enable depth test, which stores fragments in the Z-buffer
+	glEnable(GL_DEPTH_TEST);
+
+	// Accept fragment if it is closer to the camera than the former one
+	glDepthFunc(GL_LESS);
+
+	// Cull triangles which do not have a normal facing towards the camera
+	glEnable(GL_CULL_FACE);
+
+	// Create a Vertex Array Object (VAO) and set it as the current one
+	//GLuint VertexArrayID;
+	glGenVertexArrays(1, &VertexArrayID);
+	glBindVertexArray(VertexArrayID);
+
+	// Create and compile the GLSL program from the shaders
+	depthProgramID = LoadShaders("/home/thomas/Documents/Minimisation/apps/specular_estimation/src/DepthRTT.vertexshader", "/home/thomas/Documents/Minimisation/apps/specular_estimation/src/DepthRTT.fragmentshader");
+
+	// Get a handle for the "MVP" uniform
+	GLuint depthMatrixID = glGetUniformLocation(depthProgramID, "depthMVP");
+
+	// Load the texture
+	//GLuint Texture = loadDDS("C:/Users/tdwal/Documents/repos/OpenGLTutorials/tutorial16_shadowmaps/uvmap.DDS");
+	Texture = loadMat(cv_texture);
+	std::cout << "Texture loaded.\n\n";
+
+	// Read the .obj file
+	std::vector<glm::vec3> vertices;
+	std::vector<glm::vec2> uvs;
+	std::vector<glm::vec3> normals;
+	//std::vector<unsigned int> indices;
+	//bool res = loadOBJ("/home/thomas/Documents/Minimisation/apps/specular_estimation/src/sphere.obj", vertices, uvs, normals);
+	createMesh(cv_depth, cv_normals, vertices, uvs, normals, indices);
+	//createSphere(100, 100);
+	std::cout << "\nModel created.\n\n";
+
+	// Load it into a VBO
+	//GLuint vertexbuffer; // This will identify the vertex buffer
+	glGenBuffers(1, &vertexbuffer); // Generate 1 buffer, and put the resulting identifier in vertexbuffer
+	glBindBuffer(GL_ARRAY_BUFFER, vertexbuffer); // The following commands will talk about the 'vertexbuffer' buffer
+	glBufferData(GL_ARRAY_BUFFER, vertices.size() * sizeof(glm::vec3), &vertices[0], GL_STATIC_DRAW); // Give the vertices to OpenGL
+
+	//GLuint uvbuffer;
+	glGenBuffers(1, &uvbuffer);
+	glBindBuffer(GL_ARRAY_BUFFER, uvbuffer);
+	glBufferData(GL_ARRAY_BUFFER, uvs.size() * sizeof(glm::vec2), &uvs[0], GL_STATIC_DRAW);
+
+	//GLuint normalbuffer;
+	glGenBuffers(1, &normalbuffer);
+	glBindBuffer(GL_ARRAY_BUFFER, normalbuffer);
+	glBufferData(GL_ARRAY_BUFFER, normals.size() * sizeof(glm::vec3), &normals[0], GL_STATIC_DRAW);
+
+	// Generate a buffer for the indices as well
+	//GLuint elementbuffer;
+	glGenBuffers(1, &elementbuffer);
+	glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, elementbuffer);
+	glBufferData(GL_ELEMENT_ARRAY_BUFFER, indices.size() * sizeof(unsigned int), &indices[0], GL_STATIC_DRAW);
+
+
+	// ---------------------------------------------
+	// Render to Texture - specific code begins here
+	// ---------------------------------------------
+
+	// The framebuffer, which regroups 0, 1 or more textures, and 0 or 1 depth buffers
+	FramebufferName = 0;
+	glGenFramebuffers(1, &FramebufferName);
+	glBindFramebuffer(GL_FRAMEBUFFER, FramebufferName);
+
+	// Depth texture. Slower than a depth buffer, but can be sampled later in the shader
+	//GLuint depthTexture;
+	glGenTextures(1, &depthTexture);
+	glBindTexture(GL_TEXTURE_2D, depthTexture);
+	glTexImage2D(GL_TEXTURE_2D, 0, GL_DEPTH_COMPONENT16, shadowResolution, shadowResolution, 0, GL_DEPTH_COMPONENT, GL_FLOAT, 0);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_COMPARE_FUNC, GL_LEQUAL);
+	glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_COMPARE_MODE, GL_COMPARE_R_TO_TEXTURE);
+
+	glFramebufferTexture(GL_FRAMEBUFFER, GL_DEPTH_ATTACHMENT, depthTexture, 0);
+
+	// No colour output in the bound framebuffer, only depth.
+	glDrawBuffer(GL_NONE);
+
+	// Always check if there are any errors in the framebuffer
+	if (glCheckFramebufferStatus(GL_FRAMEBUFFER) != GL_FRAMEBUFFER_COMPLETE)
+		return;
+
+	// The quad's FBO. Used only for visualising the shadowmap.
+	static const GLfloat g_quad_vertex_buffer_data[] = {
+		-1.0f, -1.0f, 0.0f,
+		 1.0f, -1.0f, 0.0f,
+		-1.0f,  1.0f, 0.0f,
+		-1.0f,  1.0f, 0.0f,
+		 1.0f, -1.0f, 0.0f,
+		 1.0f,  1.0f, 0.0f,
+	};
+
+	//GLuint quad_vertexbuffer;
+	glGenBuffers(1, &quad_vertexbuffer);
+	glBindBuffer(GL_ARRAY_BUFFER, quad_vertexbuffer);
+	glBufferData(GL_ARRAY_BUFFER, sizeof(g_quad_vertex_buffer_data), g_quad_vertex_buffer_data, GL_STATIC_DRAW);
+
+	// Create and compile the GLSL program from the shaders
+	quad_programID = LoadShaders("/home/thomas/Documents/Minimisation/apps/specular_estimation/src/Passthrough.vertexshader", "/home/thomas/Documents/Minimisation/apps/specular_estimation/src/SimpleTexture.fragmentshader");
+	std::cout << "Shaders loaded.\n";
+	GLuint texID = glGetUniformLocation(quad_programID, "texture");
+	std::cout << "Uniform Location of texture found.\n";
+
+	// Create and compile the GLSL program from the shaders
+	programID = LoadShaders("/home/thomas/Documents/Minimisation/apps/specular_estimation/src/ShadowMapping.vertexshader", "/home/thomas/Documents/Minimisation/apps/specular_estimation/src/CookTorrance.fragmentshader");
+
+	// Get a handle for the "myTextureSampler" uniform
+	TextureID = glGetUniformLocation(programID, "myTextureSampler");
+	//SpecularID = glGetUniformLocation(programID, "mySpecularSampler");
+	
+	SpecularIntensityID = glGetUniformLocation(programID, "specularIntensity");
+	SpecularPowerID = glGetUniformLocation(programID, "specularPower");
+
+	// Get a handle for the "MVP" uniform
+	MatrixID = glGetUniformLocation(programID, "MVP");
+	ViewMatrixID = glGetUniformLocation(programID, "V");
+	ModelMatrixID = glGetUniformLocation(programID, "M");
+	DepthBiasID = glGetUniformLocation(programID, "DepthBiasMVP");
+	ShadowMapID = glGetUniformLocation(programID, "shadowMap");
+
+	// Get a handle for the "LightPosition" uniform
+	lightInvDirID = glGetUniformLocation(programID, "LightInvDirection_worldspace");
+
+	// Compute the MVP matrix from the light's point of view
+	//lightInvDir = glm::vec3(0, 0, 1);
+	depthProjectionMatrix = glm::perspective<float>(45.0f, 1.0f, 1.0f, 100.0f);
+	depthViewMatrix = glm::lookAt(glm::vec3(0, 0, 1), glm::vec3(0, 0, 0), glm::vec3(0, 1, 0));
+
+	// Initial position: on +Z
+	position = glm::vec3(0, 0, 1);
+	// Initial horizontal angle: toward -Z
+	horizontalAngle = M_PI;
+	// Initial vertical angle: none
+	verticalAngle = 0.0f;
+	// Initial Field of View
+	FoV = 60.0f;
+	// Mouse Speed
+	float mouseSpeed = 0.005f;
+	
+	depthMVP = depthProjectionMatrix * depthViewMatrix * depthModelMatrix;
+	bool perspectiveProjection = false;
+	int lightNumber = 0;
+
+	// Compute the MVP matrix
+	computeMatricesFromLights(windowWidth, windowHeight, position, horizontalAngle, verticalAngle, FoV, ProjectionMatrix, ViewMatrix, lightInvDir, depthProjectionMatrix, depthViewMatrix, perspectiveProjection, lightDirections, lightNumber, numberOfLights, calculateResidual, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower);
+	
+	MVP = ProjectionMatrix * ViewMatrix * ModelMatrix;
+	depthBiasMVP = biasMatrix * depthMVP;
+}
+
+
+void renderPolygonsCT(int windowWidth, int windowHeight, GLuint programID, glm::vec3 lightInvDir, GLuint MatrixID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, glm::mat4 MVP, glm::mat4 ModelMatrix, glm::mat4 ViewMatrix, glm::mat4 depthBiasMVP, GLuint SpecularIntensityID, double SpecularIntensity, GLuint SpecularPowerID, double SpecularPower) {
+
+	// Render to the screen
+	glBindFramebuffer(GL_FRAMEBUFFER, 0);
+	glViewport(0, 0, windowWidth, windowHeight); // Render on the whole framebuffer, complete from the lower left corner to the upper right
+
+	glEnable(GL_CULL_FACE);
+	glCullFace(GL_BACK); // Cull back-facing triangles -> draw only front-facing triangles
+
+	// Clear the screen
+	glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
+
+	// Use the shader
+	glUseProgram(programID);
+	
+	// Send the transformation to the currently bound shader, in the "MVP" uniform
+	glUniformMatrix4fv(MatrixID,      1, GL_FALSE, &MVP[0][0]);
+	glUniformMatrix4fv(ModelMatrixID, 1, GL_FALSE, &ModelMatrix[0][0]);
+	glUniformMatrix4fv(ViewMatrixID,  1, GL_FALSE, &ViewMatrix[0][0]);
+	glUniformMatrix4fv(DepthBiasID,   1, GL_FALSE, &depthBiasMVP[0][0]);
+	
+	glUniform1f(SpecularIntensityID, SpecularIntensity);
+	glUniform1f(SpecularPowerID,     SpecularPower);
+
+	glUniform3f(lightInvDirID, lightInvDir.x, lightInvDir.y, lightInvDir.z);
+
+	// Bind the texture in Texture Unit 0
+	glActiveTexture(GL_TEXTURE0);
+	glBindTexture(GL_TEXTURE_2D, Texture);
+	// Set the "myTextureSampler" sampler to use Texture Unit 0
+	glUniform1i(TextureID, 0);
+
+	glActiveTexture(GL_TEXTURE1);
+	glBindTexture(GL_TEXTURE_2D, depthTexture);
+	glUniform1i(ShadowMapID, 1);
+
+	// 1st attribute buffer: vertices
+	glEnableVertexAttribArray(0);
+	glBindBuffer(GL_ARRAY_BUFFER, vertexbuffer);
+	glVertexAttribPointer(
+		0,			// attribute
+		3,			// size
+		GL_FLOAT,	// type
+		GL_FALSE,	// normalized?
+		0,			// stride
+		(void*)0	// array buffer offset
+	);
+
+	// 2nd attribute buffer: UVs
+	glEnableVertexAttribArray(1);
+	glBindBuffer(GL_ARRAY_BUFFER, uvbuffer);
+	glVertexAttribPointer(
+		1,			// attribute
+		2,			// size
+		GL_FLOAT,	// type
+		GL_FALSE,	// normalized?
+		0,			// stride
+		(void*)0	// array buffer offset
+	);
+
+	// 3rd attribute buffer: normals
+	glEnableVertexAttribArray(2);
+	glBindBuffer(GL_ARRAY_BUFFER, normalbuffer);
+	glVertexAttribPointer(
+		2,			// attribute
+		3,			// size
+		GL_FLOAT,	// type
+		GL_FALSE,	// normalized?
+		0,			// stride
+		(void*)0	// array buffer offset
+	);
+
+	// Index buffer
+	glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, elementbuffer);
+
+	// Draw the triangles
+	glDrawElements(
+		GL_TRIANGLES,				// mode
+		GLsizei(indices.size()),	// count
+		GL_UNSIGNED_INT,			// type
+		(void*)0					// element array buffer offset
+	);
+
+	glDisableVertexAttribArray(0);
+	glDisableVertexAttribArray(1);
+	glDisableVertexAttribArray(2);
+
+}
+
+void terminateOpenGLCT(GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, GLuint programID, GLuint depthProgramID, GLuint quad_programID, GLuint Texture, GLuint FramebufferName, GLuint depthTexture, GLuint quad_vertexbuffer, GLuint VertexArrayID) {
+	
+	// Cleanup VBO and shader
+	glDeleteBuffers(1, &vertexbuffer);
+	glDeleteBuffers(1, &uvbuffer);
+	glDeleteBuffers(1, &normalbuffer);
+	glDeleteBuffers(1, &elementbuffer);
+	glDeleteProgram(programID);
+	glDeleteProgram(depthProgramID);
+	glDeleteProgram(quad_programID);
+	glDeleteTextures(1, &Texture);
+
+	glDeleteFramebuffers(1, &FramebufferName);
+	glDeleteTextures(1, &depthTexture);
+	glDeleteBuffers(1, &quad_vertexbuffer);
+	glDeleteVertexArrays(1, &VertexArrayID);
+
+	// Close OpenGL window and terminate GLFW
+	glfwTerminate();
+}
+
+void viewModelCT(cv::Vec3d& residual, double& totalResidual, cv::Mat& residualImage, bool calculateResidual, glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat> textureImages, int lightNumber, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath) {
+	do {
+		// Render the shadows
+		glm::mat4 depthModelMatrix = glm::mat4(1.0);
+		glm::mat4 depthMVP = depthProjectionMatrix * depthViewMatrix * depthModelMatrix;
+		//renderShadows(FramebufferName, shadowResolution, depthProgramID, depthProjectionMatrix, depthViewMatrix, depthMatrixID, vertexbuffer, elementbuffer, indices, depthModelMatrix, depthMVP);
+
+		// Compute the MVP matrix from keyboard and mouse input
+		glm::mat4 ModelMatrix = glm::mat4(1.0);
+		glm::mat4 MVP, ViewMatrix, depthBiasMVP, ProjectionMatrix;
+		bool perspectiveProjection = false, shadowControl;
+		//float mouseSpeed = 0.005f;
+		//computeMatricesFromInputs(width, height, position, horizontalAngle, verticalAngle, FoV, mouseSpeed, ProjectionMatrix, ViewMatrix, lightInvDir, depthProjectionMatrix, depthViewMatrix, perspectiveProjection, shadowControl, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower, calculateResidual);
+		computeMatricesFromLights(width, height, position, horizontalAngle, verticalAngle, FoV, ProjectionMatrix, ViewMatrix, lightInvDir, depthProjectionMatrix, depthViewMatrix, perspectiveProjection, lightDirections, lightNumber, numberOfLights, calculateResidual, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower);
+		MVP = ProjectionMatrix * ViewMatrix * ModelMatrix;
+
+		glm::mat4 biasMatrix(
+			0.5, 0.0, 0.0, 0.0,
+			0.0, 0.5, 0.0, 0.0,
+			0.0, 0.0, 0.5, 0.0,
+			0.5, 0.5, 0.5, 1.0
+		);
+
+		depthBiasMVP = biasMatrix * depthMVP;
+		
+			
+		// Render the polygons
+		renderPolygons(width, height, programID, lightInvDir, MatrixID, ModelMatrixID, ViewMatrixID, DepthBiasID, lightInvDirID, Texture, TextureID, depthTexture, ShadowMapID, vertexbuffer, uvbuffer, normalbuffer, elementbuffer, indices, MVP, ModelMatrix, ViewMatrix, depthBiasMVP, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower);
+		
+		// Optionally render the shadowmap (for debug only)
+		//void renderShadowMap(GLuint quad_programID, GLuint depthTexture, GLuint texID, GLuint quad_vertexbuffer);
+				
+		if (calculateResidual) {
+			double sum = 0;
+			meanSquaredError(lightNumber, height, width, textureImages, sum);
+			//computeResiduals(residual, residuals, totalResidual, residualImage, calculateResidual, depthProjectionMatrix, depthViewMatrix, width, height, position, horizontalAngle, verticalAngle, FoV, lightInvDir, lightDirections, textureImages, lightNumber, numberOfLights, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower, programID, ModelMatrixID, ViewMatrixID, DepthBiasID, lightInvDirID, Texture, TextureID, depthTexture, ShadowMapID, vertexbuffer, uvbuffer, normalbuffer, elementbuffer, indices, MatrixID, modelPath);
+			std::cout << "Light number: " << lightNumber << ", per-pixel residual = " << totalResidual << ", specular intensity = " << SpecularIntensity << ", specular power = " << SpecularPower << std::endl;
+		}
+		
+		
+		/*
+		// TODO: Check if the viewpoint has been changed
+		if (calculateResidual) {
+			// Compute the residual image and sum of differences
+			cv::Mat residualImage;
+			meanSquaredError(lightNumber, height, width, textureImages, residualImage, residuals, residual, SpecularIntensity, SpecularPower);
+			
+			specularMinimisation(SpecularIntensity, SpecularPower, residual);
+			calculateResidual = false;
+		}
+		*/
+		
+		//totalResidual = returnResidual(lightNumber, height, width, textureImages, residualImage, residuals, SpecularIntensity, SpecularPower);
+		//meanSquaredError(lightNumber, height, width, textureImages, residualImage, residual, totalResidual, SpecularIntensity, SpecularPower);
+		//residuals.push_back(residual);
+		
+		// Set up the only cost function (also known as residual). This uses auto-differentiation to obtain the derivative (jacobian).
+		// AutoDiffCostFunction<CostFunctor, (Dimensions of Residual), (Dimensions of Variables)>
+		//CostFunction* cost_function = new AutoDiffCostFunction<CostFunctor, 1, 2>(new CostFunctor);
+		//problem.AddResidualBlock(cost_function, NULL, &SpecularIntensity, &SpecularPower);
+		
+		// Swap buffers
+		glfwSwapBuffers(window);
+		glfwPollEvents();
+	} while (glfwGetKey(window, GLFW_KEY_ESCAPE) != GLFW_PRESS && glfwWindowShouldClose(window) == 0); // Check if the ESC key was pressed or the window was closed
+}
+
+void RenderSyntheticCT(glm::mat4 depthProjectionMatrix, glm::mat4 depthViewMatrix, int width, int height, glm::vec3& position, float& horizontalAngle, float& verticalAngle, float& FoV, glm::vec3& lightInvDir, cv::Mat lightDirections, std::vector<cv::Mat>& textureImages, int numberOfLights, GLuint& SpecularIntensityID, double& SpecularIntensity, GLuint& SpecularPowerID, double& SpecularPower, GLuint programID, GLuint ModelMatrixID, GLuint ViewMatrixID, GLuint DepthBiasID, GLuint lightInvDirID, GLuint Texture, GLuint TextureID, GLuint depthTexture, GLuint ShadowMapID, GLuint vertexbuffer, GLuint uvbuffer, GLuint normalbuffer, GLuint elementbuffer, std::vector<unsigned int> indices, GLuint MatrixID, std::string modelPath) {
+	for(int lightNumber = 0; lightNumber < numberOfLights; lightNumber++) {
+		// Render the shadows
+		glm::mat4 depthModelMatrix = glm::mat4(1.0);
+		glm::mat4 depthMVP = depthProjectionMatrix * depthViewMatrix * depthModelMatrix;
+		
+		// Compute the MVP matrix from keyboard and mouse input
+		glm::mat4 ModelMatrix = glm::mat4(1.0);
+		glm::mat4 MVP, ViewMatrix, depthBiasMVP, ProjectionMatrix;
+		bool perspectiveProjection = false, shadowControl, calculateResidual = false;
+		computeMatricesFromLights(width, height, position, horizontalAngle, verticalAngle, FoV, ProjectionMatrix, ViewMatrix, lightInvDir, depthProjectionMatrix, depthViewMatrix, perspectiveProjection, lightDirections, lightNumber, numberOfLights, calculateResidual, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower);
+		MVP = ProjectionMatrix * ViewMatrix * ModelMatrix;
+
+		glm::mat4 biasMatrix(
+			0.5, 0.0, 0.0, 0.0,
+			0.0, 0.5, 0.0, 0.0,
+			0.0, 0.0, 0.5, 0.0,
+			0.5, 0.5, 0.5, 1.0
+		);
+
+		depthBiasMVP = biasMatrix * depthMVP;
+		
+		// Render the polygons
+		renderPolygons(width, height, programID, lightInvDir, MatrixID, ModelMatrixID, ViewMatrixID, DepthBiasID, lightInvDirID, Texture, TextureID, depthTexture, ShadowMapID, vertexbuffer, uvbuffer, normalbuffer, elementbuffer, indices, MVP, ModelMatrix, ViewMatrix, depthBiasMVP, SpecularIntensityID, SpecularIntensity, SpecularPowerID, SpecularPower);
+		
+		// Create an empty Mat the same size as the image
+		cv::Mat temp = cv::Mat(height, width, CV_8UC3);
+		
+		// Read the image into the Mat as an 8-bit colour image
+		glReadPixels(0, 0, width, height, GL_BGR, GL_UNSIGNED_BYTE, temp.data);
+		
+		// 3D models have the origin in the bottom-left corner, while images have the origin in the top-left corner. The image is flipped to convert between the two.
+		cv::flip(temp, temp, 0);
+
+		//temp = addNoise(temp, 0.1);
+		temp = AddGaussianNoise_Opencv(temp, 0, 0.001);
+		textureImages.push_back(temp);
+
+		std::ostringstream a, b, c;
+		a << lightNumber;
+		b << SpecularIntensity;
+		c << SpecularPower;
+		std::string imageNumber = a.str();
+		std::string specularIntensity = b.str();
+		std::string specularPower = c.str();
+
+		//cv::imwrite(modelPath + "Intensity" + specularIntensity + "Power" + specularPower + imageNumber + ".png", temp);
+		cv::imwrite(modelPath + imageNumber + ".png", temp);
+
+		// Swap buffers
+		glfwSwapBuffers(window);
+		glfwPollEvents();
+	} //while (glfwGetKey(window, GLFW_KEY_ESCAPE) != GLFW_PRESS && glfwWindowShouldClose(window) == 0); // Check if the ESC key was pressed or the window was closed
+}
+
+
+
+
+
+
 #endif
diff --git a/apps/specular_estimation/src/ShadowMapping.fragmentshader b/apps/specular_estimation/src/ShadowMapping.fragmentshader
index d63f26fe6251040a4343027955b2721bdc7f6b71..17a335c3258e6fa7bf8ce7fcb8e84b460a696494 100644
--- a/apps/specular_estimation/src/ShadowMapping.fragmentshader
+++ b/apps/specular_estimation/src/ShadowMapping.fragmentshader
@@ -51,14 +51,11 @@ void main() {
 	// Light emission properties
 	vec3 LightColor = vec3(1, 1, 1);
 
-	// The power and colour of each light could be obtained from either the chrome sphere reflection
-	// or the Macbeth colour chart
+	// The power and colour of each light could be obtained from either the chrome sphere reflection or the Macbeth colour chart
 	float LightPower = 1.0f;
 	
 	// Material properties
 	float ambientPower = 0.0f;
-	//specularPower = 5.0f;
-	//specularIntensity = 0.5f;
 	vec3 MaterialDiffuseColor  = texture( myTextureSampler, UV ).rgb;
 	vec3 MaterialAmbientColor  = vec3(ambientPower, ambientPower, ambientPower) * MaterialDiffuseColor;
 	
diff --git a/apps/specular_estimation/src/specular_estimation.cc b/apps/specular_estimation/src/specular_estimation.cc
index 9f05e043e8677dad9642f4e73d40916ce1edc339..5a248f4ca5676114b292f6f524bcac3959036250 100644
--- a/apps/specular_estimation/src/specular_estimation.cc
+++ b/apps/specular_estimation/src/specular_estimation.cc
@@ -7,9 +7,9 @@ void denseSyntheticSample(std::string imageName, std::string calibration, double
 
 int main(int argc, char** argv) {
 	
-	int windowWidth = 1092, windowHeight = 728;
+	/*int windowWidth = 1092, windowHeight = 728;
 	cv::Mat black = cv::Mat(windowHeight, windowWidth, CV_8UC3, cv::Scalar::all(0));
-	calibrateLights(windowWidth, windowHeight, black);
+	calibrateLights(windowWidth, windowHeight, black);*/
 
 	// Required for the Ceres solver
 	google::InitGoogleLogging(argv[0]);
diff --git a/bin/specular_estimation b/bin/specular_estimation
index dd728f6f6c8af3dcea6eec57052514a615149049..1970896bb6a1733e85e36fa87ad7c8c540bc0e4f 100755
Binary files a/bin/specular_estimation and b/bin/specular_estimation differ