Initial commit

.gitattributes

+# Auto detect text files and perform LF normalization
+* text=auto
+
+# Explicitly declare text files we want to always be normalized and converted 
+# to native line endings on checkout.
+*.c text
+*.cpp text
+*.h text
+*.hxx text
+*.hpp text
+
+# Denote all files that are truly binary and should not be modified.
+*.png binary
+*.jpg binary
+
+# Standard to msysgit
+*.doc	 diff=astextplain
+*.DOC	 diff=astextplain
+*.docx diff=astextplain
+*.DOCX diff=astextplain
+*.dot  diff=astextplain
+*.DOT  diff=astextplain
+*.pdf  diff=astextplain
+*.PDF	 diff=astextplain
+*.rtf	 diff=astextplain
+*.RTF	 diff=astextplain

.gitignore

+# Ignore the user-specific initial_cache file
+initial_cache.cmake

CMakeLists.txt

+project(eos)
+cmake_minimum_required(VERSION 2.8.7)
+set(eos_VERSION_MAJOR 0)
+set(eos_VERSION_MINOR 1)
+set(eos_VERSION_PATCH 0)
+
+# g++ needs a compiler flag to enable C++11 support
+include(CheckCXXCompilerFlag)
+check_cxx_compiler_flag(-std=c++11 HAS_CXX11_FLAG)
+if (HAS_CXX11_FLAG)
+	set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+endif()
+
+# All the options for building the library. Can be changed on the command-line or in initial_cache.cmake.
+message(STATUS "Options:")
+option(BUILD_EXAMPLES "Build the example applications." ON)
+message(STATUS "BUILD_EXAMPLES: ${BUILD_EXAMPLES}")
+option(BUILD_DOCUMENTATION "Build the library documentation." ON)
+message(STATUS "BUILD_DOCUMENTATION: ${BUILD_DOCUMENTATION}")
+
+# Build a CPack driven installer package:
+include(InstallRequiredSystemLibraries) # This module will include any runtime libraries that are needed by the project for the current platform
+set(CPACK_RESOURCE_FILE_LICENSE "${CMAKE_CURRENT_SOURCE_DIR}/LICENSE")
+set(CPACK_PACKAGE_VERSION_MAJOR "${eos_VERSION_MAJOR}")
+set(CPACK_PACKAGE_VERSION_MINOR "${eos_VERSION_MINOR}")
+set(CPACK_PACKAGE_VERSION_PATCH "${eos_VERSION_PATCH}")
+include(CPack)
+
+# Find dependencies:
+find_package(OpenCV 2.4.3 REQUIRED core)
+message(STATUS "OpenCV include dir found at ${OpenCV_INCLUDE_DIRS}")
+message(STATUS "OpenCV lib dir found at ${OpenCV_LIB_DIR}")
+
+find_package(Boost 1.48.0 COMPONENTS system REQUIRED)
+if(Boost_FOUND)
+  message(STATUS "Boost found at ${Boost_INCLUDE_DIRS}")
+else(Boost_FOUND)
+  message(FATAL_ERROR "Boost not found")
+endif()
+
+# Source and header files:
+set(SOURCES
+	src/eos/core/LandmarkMapper.cpp
+	src/eos/morphablemodel/PcaModel.cpp
+	src/eos/morphablemodel/MorphableModel.cpp
+	src/eos/morphablemodel/io/cvssp.cpp
+	src/eos/fitting/AffineCameraEstimation.cpp
+	src/eos/fitting/LinearShapeFitting.cpp
+	src/eos/render/Mesh.cpp
+	src/eos/render/utils.cpp
+)
+
+set(HEADERS
+	include/eos/core/LandmarkMapper.hpp
+	include/eos/morphablemodel/PcaModel.hpp
+	include/eos/morphablemodel/MorphableModel.hpp
+	include/eos/morphablemodel/io/cvssp.hpp
+	include/eos/fitting/AffineCameraEstimation.hpp
+	include/eos/fitting/LinearShapeFitting.hpp
+	include/eos/render/Mesh.hpp
+	include/eos/render/utils.hpp
+)
+
+# Add header includes:
+include_directories("include")
+include_directories(${Boost_INCLUDE_DIRS})
+include_directories(${OpenCV_INCLUDE_DIRS})
+
+# Make the library:
+add_library(eos ${SOURCES} ${HEADERS})
+target_link_libraries(eos ${OpenCV_LIBS} ${Boost_LIBRARIES})
+
+# The install target:
+install(TARGETS eos DESTINATION lib)
+install(DIRECTORY ${CMAKE_SOURCE_DIR}/include/ DESTINATION include)
+install(DIRECTORY ${CMAKE_SOURCE_DIR}/share/ DESTINATION share)
+
+if(BUILD_EXAMPLES)
+	add_subdirectory(examples)
+endif()
+
+if(BUILD_DOCUMENTATION)
+	add_subdirectory(doc)
+endif()

README.md

+# Eos: A 3D Morphable Model fitting library in modern C++11/14.
+
+
+Eos is a tiny 3D Morphable Model fitting library that provides just the bare minimum to load a model and perform camera and shape fitting. It's written in modern C++11/14.
+
+At the moment, it mainly provides the following functionality:
+
+* MorphableModel class to represent a 3DMM (using OpenCVs `cv::Mat`)
+* Camera estimation, implementation of the _Gold Standard Algorithm_ from Multiple View Geometry, Hartley & Zisserman
+* Shape fitting, implementation of the linear shape-to-landmarks fitting of O. Aldrian & W. Smith, Inverse Rendering of Faces with a 3D Morphable Model, PAMI 2013.
+
+## Usage
+
+* Tested with the following compilers: gcc-4.8.2, clang-3.5, Visual Studio 2013
+* Needed dependencies: CMake, Boost system filesystem program_options (1.54.0), OpenCV core, imgproc, highgui (2.4.3). Older versions might work as well.
+
+### Build the examples and tests
+
+* copy `initial_cache.cmake.template` to `initial_cache.cmake`, edit the necessary paths
+* create a build directory next to the `eos` folder: `mkdir build; cd build`
+* `cmake -C ../eos/initial_cache.cmake -G "<your favourite generator>" ../eos -DCMAKE_INSTALL_PREFIX=install/`
+* build using your favourite tools, e.g. `make; make install` or open the solution in Visual Studio.
+
+
+## Sample code
+
+See examples/fit_model.cpp.
+
+
+## Documentation
+
+The code is fully documented with Doxygen and you can also build the `doc` target.
+
+
+## License & contributions
+
+This code is licensed under the Apache License, Version 2.0
+
+Contributions are very welcome! (best in the form of pull requests.) Please use Github issues for any bug reports, ideas, and discussions.
+
+If you use this code in your own work, please cite the following paper: _Random Cascaded-Regression Copse for Robust Facial Landmark Detection_, Z. Feng, P. Huber, J. Kittler, W. Christmas, X.J. Wu, IEEE Signal Processing Letters, Vol: 22(1), 2015 (http://ieeexplore.ieee.org/xpl/articleDetails.jsp?tp=&arnumber=6877655).
+
+(_We are working on publishing a paper that more closely resembles this library, and we will replace the previously mentioned paper here with the new one as soon as it is published._)

cmake/FindEigen3.cmake

+# - Try to find Eigen3 lib
+#
+# This module supports requiring a minimum version, e.g. you can do
+#   find_package(Eigen3 3.1.2)
+# to require version 3.1.2 or newer of Eigen3.
+#
+# Once done this will define
+#
+#  EIGEN3_FOUND - system has eigen lib with correct version
+#  EIGEN3_INCLUDE_DIR - the eigen include directory
+#  EIGEN3_VERSION - eigen version
+
+# Copyright (c) 2006, 2007 Montel Laurent, <montel@kde.org>
+# Copyright (c) 2008, 2009 Gael Guennebaud, <g.gael@free.fr>
+# Copyright (c) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
+# Redistribution and use is allowed according to the terms of the 2-clause BSD license.
+
+if(NOT Eigen3_FIND_VERSION)
+  if(NOT Eigen3_FIND_VERSION_MAJOR)
+    set(Eigen3_FIND_VERSION_MAJOR 2)
+  endif(NOT Eigen3_FIND_VERSION_MAJOR)
+  if(NOT Eigen3_FIND_VERSION_MINOR)
+    set(Eigen3_FIND_VERSION_MINOR 91)
+  endif(NOT Eigen3_FIND_VERSION_MINOR)
+  if(NOT Eigen3_FIND_VERSION_PATCH)
+    set(Eigen3_FIND_VERSION_PATCH 0)
+  endif(NOT Eigen3_FIND_VERSION_PATCH)
+
+  set(Eigen3_FIND_VERSION "${Eigen3_FIND_VERSION_MAJOR}.${Eigen3_FIND_VERSION_MINOR}.${Eigen3_FIND_VERSION_PATCH}")
+endif(NOT Eigen3_FIND_VERSION)
+
+macro(_eigen3_check_version)
+  file(READ "${EIGEN3_INCLUDE_DIR}/Eigen/src/Core/util/Macros.h" _eigen3_version_header)
+
+  string(REGEX MATCH "define[ \t]+EIGEN_WORLD_VERSION[ \t]+([0-9]+)" _eigen3_world_version_match "${_eigen3_version_header}")
+  set(EIGEN3_WORLD_VERSION "${CMAKE_MATCH_1}")
+  string(REGEX MATCH "define[ \t]+EIGEN_MAJOR_VERSION[ \t]+([0-9]+)" _eigen3_major_version_match "${_eigen3_version_header}")
+  set(EIGEN3_MAJOR_VERSION "${CMAKE_MATCH_1}")
+  string(REGEX MATCH "define[ \t]+EIGEN_MINOR_VERSION[ \t]+([0-9]+)" _eigen3_minor_version_match "${_eigen3_version_header}")
+  set(EIGEN3_MINOR_VERSION "${CMAKE_MATCH_1}")
+
+  set(EIGEN3_VERSION ${EIGEN3_WORLD_VERSION}.${EIGEN3_MAJOR_VERSION}.${EIGEN3_MINOR_VERSION})
+  if(${EIGEN3_VERSION} VERSION_LESS ${Eigen3_FIND_VERSION})
+    set(EIGEN3_VERSION_OK FALSE)
+  else(${EIGEN3_VERSION} VERSION_LESS ${Eigen3_FIND_VERSION})
+    set(EIGEN3_VERSION_OK TRUE)
+  endif(${EIGEN3_VERSION} VERSION_LESS ${Eigen3_FIND_VERSION})
+
+  if(NOT EIGEN3_VERSION_OK)
+
+    message(STATUS "Eigen3 version ${EIGEN3_VERSION} found in ${EIGEN3_INCLUDE_DIR}, "
+                   "but at least version ${Eigen3_FIND_VERSION} is required")
+  endif(NOT EIGEN3_VERSION_OK)
+endmacro(_eigen3_check_version)
+
+if (EIGEN3_INCLUDE_DIR)
+
+  # in cache already
+  _eigen3_check_version()
+  set(EIGEN3_FOUND ${EIGEN3_VERSION_OK})
+
+else (EIGEN3_INCLUDE_DIR)
+
+  find_path(EIGEN3_INCLUDE_DIR NAMES signature_of_eigen3_matrix_library
+      PATHS
+      ${CMAKE_INSTALL_PREFIX}/include
+      ${KDE4_INCLUDE_DIR}
+      PATH_SUFFIXES eigen3 eigen
+    )
+
+  if(EIGEN3_INCLUDE_DIR)
+    _eigen3_check_version()
+  endif(EIGEN3_INCLUDE_DIR)
+
+  include(FindPackageHandleStandardArgs)
+  find_package_handle_standard_args(Eigen3 DEFAULT_MSG EIGEN3_INCLUDE_DIR EIGEN3_VERSION_OK)
+
+  mark_as_advanced(EIGEN3_INCLUDE_DIR)
+
+endif(EIGEN3_INCLUDE_DIR)
+

doc/CMakeLists.txt

+find_package(Doxygen)
+if(NOT DOXYGEN_FOUND)
+	message(FATAL_ERROR "Doxygen is needed to build the documentation.")
+endif()
+
+set(doxyfile_in ${CMAKE_CURRENT_SOURCE_DIR}/Doxyfile.in)
+set(doxyfile ${CMAKE_CURRENT_BINARY_DIR}/Doxyfile)
+
+configure_file(${doxyfile_in} ${doxyfile} @ONLY)
+
+add_custom_target(doc
+	COMMAND ${DOXYGEN_EXECUTABLE} ${doxyfile}
+	WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}
+	COMMENT "Generating API documentation with Doxygen"
+	VERBATIM
+)
+
+install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/html DESTINATION doc)

doc/Doxyfile.in

+PROJECT_NAME           = "@CMAKE_PROJECT_NAME@"
+PROJECT_NUMBER         = @eos_VERSION_MAJOR@.@eos_VERSION_MINOR@.@eos_VERSION_PATCH@
+STRIP_FROM_PATH        = @PROJECT_SOURCE_DIR@
+INPUT                  = @PROJECT_SOURCE_DIR@/README.md \
+                         @PROJECT_SOURCE_DIR@/include \
+						 @PROJECT_SOURCE_DIR@/src
+RECURSIVE              = YES
+USE_MDFILE_AS_MAINPAGE = README.md
+
+USE_MATHJAX            = YES
+MATHJAX_RELPATH        = http://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML

examples/CMakeLists.txt

+project(examples)
+cmake_minimum_required(VERSION 2.8.7)
+
+# The examples need a few additional dependencies (e.g. boost filesystem and OpenCV highgui):
+find_package(OpenCV 2.4.3 REQUIRED core imgproc highgui)
+message(STATUS "OpenCV include dir found at ${OpenCV_INCLUDE_DIRS}")
+message(STATUS "OpenCV lib dir found at ${OpenCV_LIB_DIR}")
+
+find_package(Boost 1.48.0 COMPONENTS system filesystem program_options REQUIRED)
+if(Boost_FOUND)
+  message(STATUS "Boost found at ${Boost_INCLUDE_DIRS}")
+else(Boost_FOUND)
+  message(FATAL_ERROR "Boost not found")
+endif()
+
+# Model fitting (affine cam & shape to landmarks) example:
+add_executable(fit_model fit_model.cpp)
+target_link_libraries(fit_model eos ${OpenCV_LIBS} ${Boost_LIBRARIES})
+# TODO Add dependency to library-target??
+
+# install target:
+install(TARGETS fit_model DESTINATION bin)
+install(DIRECTORY ${CMAKE_SOURCE_DIR}/examples/data DESTINATION bin)

examples/data/image_0001.pts

+version: 1
+n_points:  68
+{
+115.167660 220.807529
+116.164839 245.721357
+120.208690 270.389841
+126.818237 293.156059
+135.630967 313.719094
+148.115667 336.485312
+162.069155 357.782742
+175.244288 375.646752
+196.629641 382.645477
+218.251194 381.353984
+240.107363 364.628222
+258.624908 347.083661
+276.182758 326.369112
+287.933064 300.889178
+292.640683 273.327417
+296.864703 246.434902
+300.377962 218.496828
+123.880661 214.575887
+139.490678 215.929559
+153.266075 219.788449
+167.524839 225.304335
+181.178785 232.537082
+209.926602 232.365554
+228.769747 226.989098
+245.270669 222.183443
+262.623085 220.119501
+278.396024 221.284493
+195.438309 247.230838
+193.878987 266.412131
+192.715650 283.702654
+191.413450 299.777554
+175.394672 303.212546
+183.252244 309.321263
+192.860805 313.386123
+203.029501 309.992086
+211.167139 305.865789
+143.511472 240.605877
+152.873610 234.904027
+165.885257 234.655325
+175.515413 243.565778
+165.006732 245.420441
+152.026189 244.447496
+222.758288 246.624566
+233.246081 239.604592
+245.639459 238.315050
+255.057142 244.458884
+245.790086 249.092411
+234.108232 249.773901
+162.766194 318.246815
+173.239365 319.463029
+184.243161 320.665411
+193.282615 323.761650
+203.914072 322.478920
+217.708849 323.720545
+234.039780 323.266218
+216.911898 345.507827
+203.746197 349.087438
+192.056977 349.252320
+181.960552 346.648746
+170.687582 338.371143
+169.898521 322.774535
+183.940980 324.969545
+193.291010 327.191247
+203.605936 327.453212
+227.819411 327.220701
+203.605936 336.539751
+193.356785 336.194279
+183.164758 334.765621
+}
\ No newline at end of file

examples/data/notes.txt

+The data in this folder was downloaded from http://ibug.doc.ic.ac.uk/resources/300-W/ on 15 December 2014.

examples/fit_model.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: examples/fit_model.cpp
+ *
+ * Copyright 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/core/LandmarkMapper.hpp"
+#include "eos/fitting/AffineCameraEstimation.hpp"
+#include "eos/fitting/LinearShapeFitting.hpp"
+#include "eos/morphablemodel/io/cvssp.hpp"
+#include "eos/render/utils.hpp"
+
+#include "opencv2/core/core.hpp"
+#include "opencv2/highgui/highgui.hpp"
+
+#ifdef WIN32
+	#define BOOST_ALL_DYN_LINK	// Link against the dynamic boost lib. Seems to be necessary because we use /MD, i.e. link to the dynamic CRT.
+	#define BOOST_ALL_NO_LIB	// Don't use the automatic library linking by boost with VS2010 (#pragma ...). Instead, we specify everything in cmake.
+#endif
+#include "boost/program_options.hpp"
+#include "boost/filesystem.hpp"
+
+#include <vector>
+#include <iostream>
+#include <fstream>
+
+using namespace eos;
+namespace po = boost::program_options;
+namespace fs = boost::filesystem;
+using cv::Mat;
+using cv::Vec2f;
+using cv::Vec3f;
+using cv::Vec4f;
+using std::cout;
+using std::endl;
+using std::vector;
+using std::string;
+
+/**
+ * Reads an ibug .pts landmark file and returns an ordered vector with
+ * the 68 2D landmark coordinates.
+ *
+ * @param[in] filename Path to a .pts file.
+ * @return An ordered vector with the 68 ibug landmarks.
+ */
+vector<Vec2f> readPtsLandmarks(std::string filename)
+{
+	using std::getline;
+	vector<Vec2f> landmarks;
+	landmarks.reserve(68);
+
+	std::ifstream file(filename);
+	if (!file.is_open()) {
+		throw std::runtime_error(string("Unable to open landmark file: " + filename));
+	}
+	
+	string line;
+	// Skip the first 3 lines, they're header lines:
+	getline(file, line); // 'version: 1'
+	getline(file, line); // 'n_points : 68'
+	getline(file, line); // '{'
+
+	while (getline(file, line))
+	{
+		if (line == "}") { // end of the file
+			break;
+		}
+		std::stringstream lineStream(line);
+		Vec2f landmark(0.0f, 0.0f);
+		if (!(lineStream >> landmark[0] >> landmark[1])) {
+			throw std::runtime_error(string("Landmark format error while parsing the line: " + line));
+		}
+		// From the iBug website:
+		// "Please note that the re-annotated data for this challenge are saved in the Matlab convention of 1 being 
+		// the first index, i.e. the coordinates of the top left pixel in an image are x=1, y=1."
+		// ==> So we shift every point by 1:
+		landmark[0] -= 1.0f;
+		landmark[1] -= 1.0f;
+		landmarks.emplace_back(landmark);
+	}
+	return landmarks;
+};
+
+/**
+ * This app demonstrates estimation of the camera and fitting of the shape
+ * model of a 3D Morphable Model from an ibug LFPW image with its landmarks.
+ *
+ * First, the 68 ibug landmarks are loaded from the .pts file and converted
+ * to vertex indices using the LandmarkMapper. Then, an affine camera matrix
+ * is estimated, and then, using this camera matrix, the shape is fitted
+ * to the landmarks.
+ */
+int main(int argc, char *argv[])
+{
+	fs::path modelfile, isomapfile, imagefile, landmarksfile, mappingsfile;
+	try {
+		po::options_description desc("Allowed options");
+		desc.add_options()
+			("help,h",
+				"display the help message")
+			("model,m", po::value<fs::path>(&modelfile)->required(),
+				"a CVSSP .scm Morphable Model file")
+			("isomap,t", po::value<fs::path>(&isomapfile),
+				"optional isomap containing the texture mapping coordinates")
+			("image,i", po::value<fs::path>(&imagefile)->required()->default_value("data/image_0001.png"),
+				"an input image")
+			("landmarks,l", po::value<fs::path>(&landmarksfile)->required()->default_value("data/image_0001.pts"),
+				"2D landmarks for the image, in ibug .pts format")
+			("mapping,p", po::value<fs::path>(&mappingsfile)->required()->default_value("../share/ibug2did.txt"),
+				"landmark identifier to model vertex number mapping")
+			;
+		po::variables_map vm;
+		po::store(po::command_line_parser(argc, argv).options(desc).run(), vm);
+		if (vm.count("help")) {
+			cout << "Usage: fit_model [options]" << endl;
+			cout << desc;
+			return EXIT_SUCCESS;
+		}
+		po::notify(vm);
+	}
+	catch (po::error& e) {
+		cout << "Error while parsing command-line arguments: " << e.what() << endl;
+		cout << "Use --help to display a list of options." << endl;
+		return EXIT_SUCCESS;
+	}
+
+	// Load the image, landmarks, LandmarkMapper and the Morphable Model:
+	Mat image = cv::imread(imagefile.string());
+	auto landmarks = readPtsLandmarks(landmarksfile.string());
+	morphablemodel::MorphableModel morphableModel = morphablemodel::loadScmModel(modelfile, isomapfile);
+	morphablemodel::LandmarkMapper landmarkMapper = mappingsfile.empty() ? morphablemodel::LandmarkMapper() : morphablemodel::LandmarkMapper(mappingsfile);
+
+	// Draw the loaded landmarks:
+	Mat outimg = image.clone();
+	for (auto&& lm : landmarks) {
+		cv::rectangle(outimg, cv::Point2f(lm[0] - 2.0f, lm[1] - 2.0f), cv::Point2f(lm[0] + 2.0f, lm[1] + 2.0f), { 255, 0, 0 });
+	}
+	
+	// Convert the landmarks to clip-space:
+	std::transform(begin(landmarks), end(landmarks), begin(landmarks), [&image](const Vec2f& lm) { return render::screenToClipSpace(lm, image.cols, image.rows); });
+
+	// These will be the final 2D and 3D points used for the fitting:
+	vector<Vec4f> modelPoints; ///< the points in the 3D shape model
+	vector<int> vertexIndices; ///< their vertex indices
+	vector<Vec2f> imagePoints; ///< the corresponding 2D landmark points
+
+	// Sub-select all the landmarks which we have a mapping for (i.e. that are defined in the 3DMM):
+	int ibugId = 1;
+	for (int i = 0; i < landmarks.size(); ++i) {
+		try {
+			int vertexIdx = boost::lexical_cast<int>(landmarkMapper.convert(std::to_string(ibugId)));
+			Vec4f vertex = morphableModel.getShapeModel().getMeanAtPoint(vertexIdx);
+			modelPoints.emplace_back(vertex);
+			vertexIndices.emplace_back(vertexIdx);
+			imagePoints.emplace_back(landmarks[i]);
+		}
+		catch (std::out_of_range& e) {
+			// just continue if the point isn't defined in the mapping
+		}
+		++ibugId;
+	}
+	
+	// Estimate the camera from the 2D - 3D point correspondences
+	Mat affineCam = fitting::estimateAffineCamera(imagePoints, modelPoints);
+
+	// Draw the mean-face landmarks projected using the estimated camera:
+	for (auto&& vertex : modelPoints) {
+		Vec2f screenPoint = fitting::projectAffine(vertex, affineCam, image.cols, image.rows);
+		cv::circle(outimg, cv::Point2f(screenPoint), 5.0f, { 0.0f, 255.0f, 0.0f });
+	}
+
+	// Estimate the shape coefficients by fitting the shape to the landmarks:
+	float lambda = 5.0f; ///< the regularisation parameter
+	vector<float> fittedCoeffs = fitting::fitShapeToLandmarksLinear(morphableModel, affineCam, imagePoints, vertexIndices, lambda);
+
+	// Obtain the full mesh and draw it using the estimated camera:
+	render::Mesh mesh = morphableModel.drawSample(fittedCoeffs, vector<float>());
+	render::writeObj(mesh, "out.obj"); // save the mesh as obj
+
+	// Draw the projected points again, this time using the fitted model shape:
+	for (auto&& idx : vertexIndices) {
+		Vec4f modelPoint(mesh.vertices[idx][0], mesh.vertices[idx][1], mesh.vertices[idx][2], mesh.vertices[idx][3]);
+		Vec2f screenPoint = fitting::projectAffine(modelPoint, affineCam, image.cols, image.rows);
+		cv::circle(outimg, cv::Point2f(screenPoint), 3.0f, { 0.0f, 0.0f, 255.0f });
+	}
+
+	// Save the output image:
+	cv::imwrite("out.png", outimg);
+	std::cout << "Finished fitting and wrote result to out.png." << std::endl;
+
+	return EXIT_SUCCESS;
+}

include/eos/core/LandmarkMapper.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/core/LandmarkMapper.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef LANDMARKMAPPER_HPP_
+#define LANDMARKMAPPER_HPP_
+
+#ifdef WIN32
+	#define BOOST_ALL_DYN_LINK	// Link against the dynamic boost lib. Seems to be necessary because we use /MD, i.e. link to the dynamic CRT.
+	#define BOOST_ALL_NO_LIB	// Don't use the automatic library linking by boost with VS2010 (#pragma ...). Instead, we specify everything in cmake.
+#endif
+#include "boost/filesystem/path.hpp"
+
+#include <string>
+#include <map>
+
+namespace eos {
+	namespace morphablemodel {
+
+/**
+ * Represents a mapping from one kind of landmarks
+ * to a different format. Mappings are stored in a
+ * file (see share/ for an example for ibug landmarks).
+ */
+class LandmarkMapper {
+public:
+	/**
+	 * Constructs a new landmark mapper that performs an identity mapping,
+	 * that is, its output is the same as the input.
+	 *
+	 */
+	LandmarkMapper() = default;
+
+	/**
+	 * Constructs a new landmark mapper from a mappings-file.
+	 *
+	 * @param[in] filename A file with landmark mappings.
+	 */
+	LandmarkMapper(boost::filesystem::path filename);
+
+	/**
+	 * Converts the given landmark name to the mapped name.
+	 *
+	 * @param[in] landmarkName A landmark name to convert.
+	 * @return The mapped landmark name.
+	 * @throws out_of_range exception if there is no mapping
+	 *         for the given landmarkName.
+	 */
+	std::string convert(std::string landmarkName);
+
+private:
+	std::map<std::string, std::string> landmarkMappings; ///< Mapping from one landmark name to a name in a different format.
+};
+
+	} /* namespace morphablemodel */
+} /* namespace eos */
+
+#endif /* LANDMARKMAPPER_HPP_ */

include/eos/fitting/AffineCameraEstimation.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/fitting/AffineCameraEstimation.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef AFFINECAMERAESTIMATION_HPP_
+#define AFFINECAMERAESTIMATION_HPP_
+
+#include "opencv2/core/core.hpp"
+
+#include <vector>
+
+namespace eos {
+	namespace fitting {
+
+/**
+ * The Gold Standard Algorithm for estimating an affine
+ * camera matrix from world to image correspondences.
+ * See Algorithm 7.2 in Multiple View Geometry, Hartley &
+ * Zisserman, 2nd Edition, 2003.
+ *
+ * Requires >= 4 corresponding points.
+ *
+ * @param[in] imagePoints A list of 2D image points.
+ * @param[in] modelPoints Corresponding points of a 3D model.
+ * @return A 3x4 affine camera matrix (the third row is [0, 0, 0, 1]).
+ */
+cv::Mat estimateAffineCamera(std::vector<cv::Vec2f> imagePoints, std::vector<cv::Vec4f> modelPoints);
+
+/**
+ * Projects a point from world coordinates to screen coordinates.
+ * First, an estimated affine camera matrix is used to transform
+ * the point to clip space. Second, the point is transformed to
+ * screen coordinates using the window transform. The window transform
+ * also flips the y-axis (the image origin is top-left, while in
+ * clip space top is +1 and bottom is -1).
+ *
+ * @param[in] vertex A vertex in 3D space. vertex[3] = 1.0f.
+ * @param[in] affineCameraMatrix A 3x4 affine camera matrix.
+ * @param[in] screenWidth Width of the screen or window used for projection.
+ * @param[in] screenHeight Height of the screen or window used for projection.
+ * @return A vector with x and y coordinates transformed to screen coordinates.
+ */
+cv::Vec2f projectAffine(cv::Vec4f vertex, cv::Mat affineCameraMatrix, int screenWidth, int screenHeight);
+
+	} /* namespace fitting */
+} /* namespace eos */
+
+#endif /* AFFINECAMERAESTIMATION_HPP_ */

include/eos/fitting/LinearShapeFitting.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/fitting/LinearShapeFitting.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef LINEARSHAPEFITTING_HPP_
+#define LINEARSHAPEFITTING_HPP_
+
+#include "eos/morphablemodel/MorphableModel.hpp"
+
+#include "opencv2/core/core.hpp"
+
+#include "boost/optional.hpp"
+
+#include <vector>
+
+namespace eos {
+	namespace fitting {
+
+/**
+ * Fits the shape of a Morphable Model to given 2D landmarks (i.e. estimates the maximum likelihood solution of the shape coefficients) as proposed in [1].
+ * It's a linear, closed-form solution fitting of the shape, with regularisation (prior towards the mean).
+ *
+ * [1] O. Aldrian & W. Smith, Inverse Rendering of Faces with a 3D Morphable Model, PAMI 2013.
+ *
+ * Note: Using less than the maximum number of coefficients to fit is not thoroughly tested yet and may contain an error.
+ * Note: Returns coefficients following standard normal distribution (i.e. all have similar magnitude). Why? Because we fit using the normalised basis?
+ * Note: The standard deviations given should be a vector, i.e. different for each landmark. This is not implemented yet.
+ *
+ * @param[in] morphableModel The Morphable Model whose shape (coefficients) are estimated.
+ * @param[in] affineCameraMatrix A 3x4 affine camera matrix from world to clip-space (should probably be of type CV_32FC1 as all our calculations are done with float).
+ * @param[in] landmarks 2D landmarks from an image, given in clip-coordinates.
+ * @param[in] vertexIds The vertex ids in the model that correspond to the 2D points.
+ * @param[in] lambda The regularisation parameter (weight of the prior towards the mean).
+ * @param[in] numCoefficientsToFit How many shape-coefficients to fit (all others will stay 0). Not tested thoroughly.
+ * @param[in] detectorStandardDeviation The standard deviation of the 2D landmarks given (e.g. of the detector used).
+ * @param[in] modelStandardDeviation The standard deviation of the 3D vertex points in the 3D model.
+ * @return The estimated shape-coefficients (alphas).
+ */
+std::vector<float> fitShapeToLandmarksLinear(morphablemodel::MorphableModel morphableModel, cv::Mat affineCameraMatrix, std::vector<cv::Vec2f> landmarks, std::vector<int> vertexIds, float lambda=20.0f, boost::optional<int> numCoefficientsToFit=boost::optional<int>(), boost::optional<float> detectorStandardDeviation=boost::optional<float>(), boost::optional<float> modelStandardDeviation=boost::optional<float>());
+
+	} /* namespace fitting */
+} /* namespace eos */
+
+#endif /* LINEARSHAPEFITTING_HPP_ */

include/eos/morphablemodel/MorphableModel.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/morphablemodel/MorphableModel.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef MORPHABLEMODEL_HPP_
+#define MORPHABLEMODEL_HPP_
+
+#include "eos/morphablemodel/PcaModel.hpp"
+
+#include "eos/render/Mesh.hpp"
+
+#include <vector>
+
+namespace eos {
+	namespace morphablemodel {
+
+/**
+ * A class representing a 3D Morphable Model.
+ * It consists of a shape- and albedo (texture) PCA model.
+ * 
+ * For the general idea of 3DMMs see T. Vetter, V. Blanz,
+ * 'A Morphable Model for the Synthesis of 3D Faces', SIGGRAPH 1999
+ */
+class MorphableModel
+{
+public:
+
+	/**
+	 * Create a Morphable Model from a shape and a color PCA model, and optional
+	 * texture coordinates.
+	 *
+	 * @param[in] shapeModel A PCA model over the shape.
+	 * @param[in] colorModel A PCA model over the color (albedo).
+	 * @param[in] textureCoordinates Optional texture coordinates for every vertex.
+	 */
+	MorphableModel(PcaModel shapeModel, PcaModel colorModel, std::vector<cv::Vec2f> textureCoordinates = std::vector<cv::Vec2f>());
+
+	/**
+	 * Returns the PCA shape model of this Morphable Model.
+	 * as a Mesh.
+	 *
+	 * @return The shape model.
+	 */
+	PcaModel getShapeModel() const;
+	
+	/**
+	 * Returns the PCA color or albedo model of this Morphable Model.
+	 *
+	 * @return The color model.
+	 */	
+	PcaModel getColorModel() const;
+
+	/**
+	 * Returns the mean of the shape- and color model as a Mesh.
+	 *
+	 * @return An mesh instance of the mean of the Morphable Model.
+	 */
+	render::Mesh getMean() const;
+
+	/**
+	 * Draws a random sample from the model, where the coefficients
+	 * for the shape- and color models are both drawn from a standard
+	 * normal (or with the given standard deviation).
+	 *
+	 * @param[in] shapeSigma The shape model standard deviation.
+	 * @param[in] colorSigma The color model standard deviation.
+	 * @return A random sample from the model.
+	 */
+	render::Mesh drawSample(float shapeSigma = 1.0f, float colorSigma = 1.0f);
+
+	/**
+	 * Returns a sample from the model with the given shape- and
+	 * color PCA coefficients. 
+	 *
+	 * If one of the given vectors is empty, the mean is used.
+	 * The coefficient vectors should contain normalised, i.e. standard normal distributed coefficients.
+	 *
+	 * @param[in] shapeCoefficients The PCA coefficients used to generate the shape sample.
+	 * @param[in] colorCoefficients The PCA coefficients used to generate the shape sample.
+	 * @return A model instance with given coefficients.
+	 */
+	render::Mesh drawSample(std::vector<float> shapeCoefficients, std::vector<float> colorCoefficients);
+
+private:
+	PcaModel shapeModel; ///< A PCA model of the shape
+	PcaModel colorModel; ///< A PCA model of vertex color information
+	std::vector<cv::Vec2f> textureCoordinates; ///< uv-coordinates for every vertex
+
+	/**
+	 * Returns whether the model has texture mapping coordinates, i.e.
+	 * coordinates for every vertex.
+	 *
+	 * @return True if the model contains texture mapping coordinates.
+	 */
+	bool hasTextureCoordinates() const {
+		return textureCoordinates.size() > 0 ? true : false;
+	};
+
+};
+
+	} /* namespace morphablemodel */
+} /* namespace eos */
+
+#endif /* MORPHABLEMODEL_HPP_ */

include/eos/morphablemodel/PcaModel.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/morphablemodel/PcaModel.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef PCAMODEL_HPP_
+#define PCAMODEL_HPP_
+
+#include "opencv2/core/core.hpp"
+
+#ifdef WIN32
+	#define BOOST_ALL_DYN_LINK	// Link against the dynamic boost lib. Seems to be necessary because we use /MD, i.e. link to the dynamic CRT.
+	#define BOOST_ALL_NO_LIB	// Don't use the automatic library linking by boost with VS2010 (#pragma ...). Instead, we specify everything in cmake.
+#endif
+#include "boost/filesystem/path.hpp"
+
+#include <string>
+#include <vector>
+#include <array>
+#include <random>
+
+namespace eos {
+	namespace morphablemodel {
+
+/**
+ * This class represents a PCA-model that consists of:
+ *   - a mean vector (y x z)
+ *   - a PCA basis matrix (unnormalised and normalised)
+ *   - a PCA variance vector.
+ *
+ * It also contains a list of triangles to built a mesh as well as a mapping
+ * from landmark points to the corresponding vertex-id in the mesh.
+ * It is able to return instances of the model as meshes.
+ */
+class PcaModel
+{
+public:
+	
+	/**
+	 * Construct a PCA model from given mean, normalised PCA basis, eigenvalues
+	 * and triangle list.
+	 *
+	 * See the documentation of the member variables for how the data should
+	 * be arranged.
+	 *
+	 * @param[in] mean The mean used to build the PCA model.
+	 * @param[in] pcaBasis The PCA basis (eigenvectors), normalised (multiplied by the eigenvalues).
+	 * @param[in] eigenvalues The eigenvalues used to build the PCA model.
+	 * @param[in] triangleList An index list of how to assemble the mesh.
+	 */
+	PcaModel(cv::Mat mean, cv::Mat pcaBasis, cv::Mat eigenvalues, std::vector<std::array<int, 3>> triangleList);
+
+	/**
+	 * Returns the number of principal components in the model.
+	 *
+	 * @return The number of principal components in the model.
+	 */
+	int getNumberOfPrincipalComponents() const;
+
+	/**
+	 * Returns the dimension of the data, i.e. the number of shape dimensions.
+	 *
+	 * As the data is arranged in a [x y z x y z ...] fashion, dividing this by
+	 * three yields the number of vertices in the model.
+	 *
+	 * @return The dimension of the data.
+	 */
+	int getDataDimension() const;
+
+	/**
+	 * Returns a list of triangles on how to assemble the vertices into a mesh.
+	 *
+	 * @return The list of triangles to build a mesh.
+	 */
+	std::vector<std::array<int, 3>> getTriangleList() const;
+
+	/**
+	 * Returns the mean of the model.
+	 *
+	 * @return The mean of the model.
+	 */
+	cv::Mat getMean() const;
+
+	/**
+	 * Return the value of the mean at a given vertex id.
+	 *
+	 * @param[in] vertexIndex A vertex id.
+	 * @return A homogeneous vector containing the values at the given vertex id.
+	 */
+	cv::Vec4f getMeanAtPoint(int vertexIndex) const;
+
+	/**
+	 * Draws a random sample from the model, where the coefficients are drawn
+	 * from a standard normal (or with the given standard deviation).
+	 *
+	 * @param[in] sigma The standard deviation.
+	 * @return A random sample from the model.
+	 */
+	cv::Mat drawSample(float sigma = 1.0f);
+
+	/**
+	 * Returns a sample from the model with the given PCA coefficients.
+	 * The given coefficients should follow a standard normal distribution, i.e.
+	 * not be "normalised" with their eigenvalues/variances.
+	 *
+	 * @param[in] coefficients The PCA coefficients used to generate the sample.
+	 * @return A model instance with given coefficients.
+	 */
+	cv::Mat drawSample(std::vector<float> coefficients);
+
+	/**
+	 * Returns the PCA basis matrix, i.e. the eigenvectors.
+	 * Each column of the matrix is an eigenvector.
+	 * The returned basis is normalised, i.e. every eigenvector
+	 * is normalised by multiplying it with the square root of its eigenvalue.
+	 * 
+	 * Returns a clone of the matrix so that the original cannot
+	 * be modified. TODO: No, don't return a clone.
+	 *
+	 * @return Returns the normalised PCA basis matrix.
+	 */
+	cv::Mat getNormalisedPcaBasis() const;
+	
+	/**
+	 * Returns the PCA basis for a particular vertex.
+	 * The returned basis is normalised, i.e. every eigenvector
+	 * is normalised by multiplying it with the square root of its eigenvalue.
+	 *
+	 * @param[in] vertexId A vertex index. Make sure it is valid.
+	 * @return A Mat that points to the rows in the original basis.
+	 */
+	cv::Mat getNormalisedPcaBasis(int vertexId) const;
+
+	/**
+	 * Returns the PCA basis matrix, i.e. the eigenvectors.
+	 * Each column of the matrix is an eigenvector.
+	 * The returned basis is unnormalised, i.e. not scaled by their eigenvalues.
+	 *
+	 * Returns a clone of the matrix so that the original cannot
+	 * be modified. TODO: No, don't return a clone.
+	 *
+	 * @return Returns the unnormalised PCA basis matrix.
+	 */
+	cv::Mat getUnnormalisedPcaBasis() const;
+
+	/**
+	 * Returns the PCA basis for a particular vertex.
+	 * The returned basis is unnormalised, i.e. not scaled by their eigenvalues.
+	 *
+	 * @param[in] vertexId A vertex index. Make sure it is valid.
+	 * @return A Mat that points to the rows in the original basis.
+	 */
+	cv::Mat getUnnormalisedPcaBasis(int vertexId) const;
+
+	/**
+	 * Returns an eigenvalue.
+	 *
+	 * @param[in] index The index of the eigenvalue to return.
+	 * @return The eigenvalue.
+	 */
+	float getEigenvalue(int index) const;
+
+private:
+	std::mt19937 engine; ///< Random number engine used to draw random coefficients.
+	
+	cv::Mat mean; ///< A 3m x 1 col-vector (xyzxyz...)', where m is the number of model-vertices.
+	cv::Mat normalisedPcaBasis; ///< The normalised PCA basis matrix. m x n (rows x cols) = numShapeDims x numShapePcaCoeffs, (=eigenvector matrix V). Each column is an eigenvector.
+	cv::Mat unnormalisedPcaBasis; ///< The unnormalised PCA basis matrix. m x n (rows x cols) = numShapeDims x numShapePcaCoeffs, (=eigenvector matrix V). Each column is an eigenvector.
+	cv::Mat eigenvalues; ///< A col-vector of the eigenvalues (variances in the PCA space).
+
+	std::vector<std::array<int, 3>> triangleList; ///< List of triangles that make up the mesh of the model.
+};
+
+/**
+ * Takes an unnormalised PCA basis matrix (a matrix consisting
+ * of the eigenvectors and normalises it, i.e. multiplies each
+ * eigenvector by the square root of its corresponding
+ * eigenvalue.
+ *
+ * @param[in] unnormalisedBasis An unnormalised PCA basis matrix.
+ * @param[in] eigenvalues A row or column vector of eigenvalues.
+ * @return The normalised PCA basis matrix.
+ */
+cv::Mat normalisePcaBasis(cv::Mat unnormalisedBasis, cv::Mat eigenvalues);
+
+/**
+ * Takes a normalised PCA basis matrix (a matrix consisting
+ * of the eigenvectors and denormalizes it, i.e. multiplies each
+ * eigenvector by 1 over the square root of its corresponding
+ * eigenvalue.
+ *
+ * @param[in] normalisedBasis A normalised PCA basis matrix.
+ * @param[in] eigenvalues A row or column vector of eigenvalues.
+ * @return The unnormalised PCA basis matrix.
+ */
+cv::Mat unnormalisePcaBasis(cv::Mat normalisedBasis, cv::Mat eigenvalues);
+
+	} /* namespace morphablemodel */
+} /* namespace eos */
+
+#endif /* PCAMODEL_HPP_ */

include/eos/morphablemodel/io/cvssp.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/morphablemodel/io/cvssp.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef IO_CVSSP_HPP_
+#define IO_CVSSP_HPP_
+
+#include "eos/morphablemodel/MorphableModel.hpp"
+
+#include "opencv2/core/core.hpp"
+
+#ifdef WIN32
+	#define BOOST_ALL_DYN_LINK	// Link against the dynamic boost lib. Seems to be necessary because we use /MD, i.e. link to the dynamic CRT.
+	#define BOOST_ALL_NO_LIB	// Don't use the automatic library linking by boost with VS2010 (#pragma ...). Instead, we specify everything in cmake.
+#endif
+#include "boost/filesystem/path.hpp"
+
+#include <vector>
+
+namespace eos {
+	namespace morphablemodel {
+
+/**
+ * Load a shape or color model from a .scm file containing
+ * a Morphable Model in the Surrey format.
+ *
+ * Note on multi-resolution models: The landmarks to vertex-id mapping is
+ * always the same. The lowest resolution model has all the landmarks defined
+ * and for the higher resolutions, the mesh is divided from that on.
+ * Note: For new landmarks we add, this might not be the case if we add them
+ * in the highest resolution model, so take care!
+ *
+ * - The pcaBasis matrix stored in the file and loaded is the orthogonal PCA basis, i.e. it is not normalised by the eigenvalues.
+ *
+ * @param[in] modelFilename A binary .scm-file containing the model.
+ * @param[in] isomapFile An optional path to an isomap containing texture coordinates.
+ * @return The Morphable Model loaded from the file.
+ */
+MorphableModel loadScmModel(boost::filesystem::path modelFilename, boost::filesystem::path isomapFile = boost::filesystem::path());
+
+/**
+ * Load a set of 2D texture coordinates pre-generated by the isomap algorithm.
+ * After loading, we rescale the coordinates to [0, 1] x [0, 1].
+ *
+ * @param[in] isomapFile Path to an isomap file containing texture coordinates.
+ * @return The 2D texture coordinates for every vertex.
+ */
+std::vector<cv::Vec2f> loadIsomap(boost::filesystem::path isomapFile);
+
+	} /* namespace morphablemodel */
+} /* namespace eos */
+
+#endif /* IO_CVSSP_HPP_ */

include/eos/render/Mesh.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/render/Mesh.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef MESH_HPP_
+#define MESH_HPP_
+
+#include "opencv2/core/core.hpp"
+
+#include <vector>
+#include <array>
+#include <string>
+
+namespace eos {
+	namespace render {
+
+/**
+ * This class represents a 3D mesh consisting of vertices and vertex color
+ * information. Additionally it stores the indices that specify which vertices
+ * to use to generate the triangle mesh out of the vertices.
+ */
+class Mesh
+{
+public:
+	std::vector<cv::Vec4f> vertices; ///< 3D vertex positions.
+	std::vector<cv::Vec3f> colors; ///< Color information for each vertex. Expected to be in RGB order.
+	std::vector<cv::Vec2f> texcoords; ///< Texture coordinates for each vertex.
+
+	std::vector<std::array<int, 3>> tvi; ///< Triangle vertex indices
+	std::vector<std::array<int, 3>> tci; ///< Triangle color indices
+};
+
+/**
+ * Writes an obj file of the given Mesh that can be read by e.g. Meshlab.
+ *
+ * @param[in] mesh The mesh to save as obj.
+ * @param[in] filename Output filename.
+ */
+void writeObj(Mesh mesh, std::string filename);
+
+	} /* namespace render */
+} /* namespace eos */
+
+#endif /* MESH_HPP_ */

include/eos/render/utils.hpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/render/utils.hpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#pragma once
+
+#ifndef RENDER_UTILS_HPP_
+#define RENDER_UTILS_HPP_
+
+#include "opencv2/core/core.hpp"
+
+namespace eos {
+	namespace render {
+
+/**
+ * Transforms a point from clip space ([-1, 1] x [-1, 1]) to
+ * image (screen) coordinates, i.e. the window transform.
+ * Note that the y-coordinate is flipped because the image origin
+ * is top-left while in clip space top is +1 and bottom is -1.
+ * No z-division is performed.
+ * Note: It should rather be called from NDC to screen space?
+ *
+ * Exactly conforming to the OpenGL viewport transform, except that
+ * we flip y at the end.
+ * Qt: Origin top-left. OpenGL: bottom-left. OCV: top-left.
+ *
+ * @param[in] clipCoordinates A point in clip coordinates.
+ * @param[in] screenWidth Width of the screen or window.
+ * @param[in] screenHeight Height of the screen or window.
+ * @return A vector with x and y coordinates transformed to screen space.
+ */
+cv::Vec2f clipToScreenSpace(cv::Vec2f clipCoordinates, int screenWidth, int screenHeight);
+
+/**
+ * Transforms a point from image (screen) coordinates to
+ * clip space ([-1, 1] x [-1, 1]).
+ * Note that the y-coordinate is flipped because the image origin
+ * is top-left while in clip space top is +1 and bottom is -1.
+ *
+ * @param[in] screenCoordinates A point in screen coordinates.
+ * @param[in] screenWidth Width of the screen or window.
+ * @param[in] screenHeight Height of the screen or window.
+ * @return A vector with x and y coordinates transformed to clip space.
+ */
+cv::Vec2f screenToClipSpace(cv::Vec2f screenCoordinates, int screenWidth, int screenHeight);
+
+	} /* namespace render */
+} /* namespace eos */
+
+#endif /* RENDER_UTILS_HPP_ */

initial_cache.cmake.template

+# Mechanism via FindLIB.cmake
+# ==============================
+# Boost:
+# -------
+# Windows: Download the pre-built binaries from http://sourceforge.net/projects/boost/files/boost-binaries/ for VS2013 (msvc-12 64bit).
+# Set the windows PATH variable to "<YOUR_BOOST_DIRECTORY>\lib64-msvc-12.0" and CMake will find it.
+# Linux: Boost can usually be installed via a package manager (e.g. apt-get install boost-all-dev) and this variable can be left uncommented.
+#set(BOOST_ROOT     "/home/user/boost/install"    CACHE STRING "Boost search location" FORCE)
+
+# Mechanism via ConfigLIB.cmake in 3rd party library directory
+# ==============================
+# OpenCV:
+# -------
+# Windows: Download the package from opencv.org, use 2.4.7.2 or never. It includes binaries for VS2013. Set this path accordingly.
+#set(OpenCV_DIR   "C:\\opencv\\install"   CACHE STRING "OpenCV config dir, where OpenCVConfig.cmake can be found" FORCE)
+# Linux: Usually can be left blank but it can be used if OpenCV is not found.
+#set(OpenCV_DIR   "/home/user/opencv/install/share/OpenCV"   CACHE STRING "OpenCV config dir, where OpenCVConfig.cmake can be found" FORCE)
+
+# Configuration options
+# ==============================
+set(BUILD_EXAMPLES ON CACHE BOOL "Build the example applications." FORCE)
+set(BUILD_DOCUMENTATION ON CACHE BOOL "Build the library documentation." FORCE)

share/ibug2did.txt

+; Mapping from the 68-point iBug annotations to the Surrey DID format (3DMM vertex indices).
+; Partial mapping
+; Created by: Patrik Huber, 19.04.2014
+
+landmarkMappings ; A mapping from input landmarks (iBug, lhs) to output landmarks (DID, rhs)
+{
+;	 9    35 ; 9 is not the chin tip (4) - the DID point is much higher up
+	18	 225 ; right eyebrow outer-corner (18)
+	20   233 ; right eyebrow middle (top-of-eyebrow?) (20)
+	22	 157 ; right eyebrow inner-corner (top-of-eyebrow?) (19)
+	23	 590 ; left eyebrow inner-corner (top-of-eyebrow?) (23)
+	25   666 ; left eyebrow middle (top-of-eyebrow?) (24)
+	27	 658 ; left eyebrow outer-corner (22)
+	31   114 ; nose-tip (3)
+	34   270 ; nose-lip junction (28)
+	37   177 ; right eye outer-corner (1)
+	40   181 ; right eye inner-corner (5)
+	43   614 ; left eye inner-corner (8)
+	46   610 ; left eye outer-corner (2)
+	49   398 ; right mouth corner (12)
+	52   329 ; upper lip middle top (14)
+	55   812 ; left mouth corner (13)
+	58   411 ; lower lip middle bottom (17)
+}

src/eos/core/LandmarkMapper.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/core/LandmarkMapper.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/core/LandmarkMapper.hpp"
+
+#include "boost/property_tree/ptree.hpp"
+#include "boost/property_tree/info_parser.hpp"
+
+#include <iostream>
+
+using boost::property_tree::ptree;
+using std::string;
+
+namespace eos {
+	namespace morphablemodel {
+
+LandmarkMapper::LandmarkMapper(boost::filesystem::path filename)
+{
+	ptree configTree;
+	try {
+		boost::property_tree::info_parser::read_info(filename.string(), configTree);
+	}
+	catch (const boost::property_tree::ptree_error& error) {
+		throw std::runtime_error(string("LandmarkMapper: Error reading landmark-mappings file: ") + error.what());
+	}
+
+	try {
+		ptree ptLandmarkMappings = configTree.get_child("landmarkMappings");
+		for (auto&& mapping : ptLandmarkMappings) {
+			landmarkMappings.insert(make_pair(mapping.first, mapping.second.get_value<string>()));
+		}
+		std::cout << "Loaded a list of " << landmarkMappings.size() << " landmark mappings." << std::endl;
+	}
+	catch (const boost::property_tree::ptree_error& error) {
+		throw std::runtime_error(string("LandmarkMapper: Error while parsing the mappings file: ") + error.what());
+	}
+	catch (const std::runtime_error& error) {
+		throw std::runtime_error(string("LandmarkMapper: Error while parsing the mappings file: ") + error.what());
+	}
+}
+
+string LandmarkMapper::convert(string landmarkName)
+{
+	if (landmarkMappings.empty()) {
+		// perform identity mapping, i.e. return the input
+		return landmarkName;
+	}
+	else {
+		return landmarkMappings.at(landmarkName); // throws an out_of_range exception if landmarkName does not match the key of any element in the map
+	}
+}
+
+	} /* namespace morphablemodel */
+} /* namespace eos */

src/eos/fitting/AffineCameraEstimation.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/fitting/AffineCameraEstimation.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/fitting/AffineCameraEstimation.hpp"
+
+#include "eos/render/utils.hpp"
+
+#include "opencv2/core/core_c.h" // for CV_REDUCE_AVG
+
+#include <iostream>
+#include <exception>
+#include <cassert>
+
+using cv::Mat;
+using cv::Vec2f;
+using cv::Vec3f;
+using cv::Vec4f;
+using std::vector;
+
+namespace eos {
+	namespace fitting {
+
+Mat estimateAffineCamera(vector<Vec2f> imagePoints, vector<Vec4f> modelPoints)
+{
+	assert(imagePoints.size() == modelPoints.size());
+
+	const auto numCorrespondences = imagePoints.size();
+	if (numCorrespondences < 4) {
+		std::string errorMsg("AffineCameraEstimation: Number of points given needs to be equal to or larger than 4.");
+		std::cout << errorMsg << std::endl;
+		throw std::runtime_error(errorMsg);
+	}
+
+	Mat matImagePoints; // will be numCorrespondences x 2, CV_32FC1
+	Mat matModelPoints; // will be numCorrespondences x 3, CV_32FC1
+	
+	for (int i = 0; i < imagePoints.size(); ++i) {
+		Mat imgPoint(1, 2, CV_32FC1);
+		imgPoint.at<float>(0, 0) = imagePoints[i][0];
+		imgPoint.at<float>(0, 1) = imagePoints[i][1];
+		matImagePoints.push_back(imgPoint);
+
+		Mat mdlPoint(1, 3, CV_32FC1);
+		mdlPoint.at<float>(0, 0) = modelPoints[i][0];
+		mdlPoint.at<float>(0, 1) = modelPoints[i][1];
+		mdlPoint.at<float>(0, 2) = modelPoints[i][2];
+		matModelPoints.push_back(mdlPoint);
+	}
+
+	// translate the centroid of the image points to the origin:
+	Mat imagePointsMean; // use non-homogeneous coords for the next few steps? (less submatrices etc overhead)
+	cv::reduce(matImagePoints, imagePointsMean, 0, CV_REDUCE_AVG);
+	imagePointsMean = cv::repeat(imagePointsMean, matImagePoints.rows, 1); // get T_13 and T_23 from imagePointsMean
+	matImagePoints = matImagePoints - imagePointsMean;
+	// scale the image points such that the RMS distance from the origin is sqrt(2):
+	// 1) calculate the average norm (root mean squared distance) of all vectors
+	float averageNorm = 0.0f; // TODO change to double!
+	for (int row = 0; row < matImagePoints.rows; ++row) {
+		averageNorm += cv::norm(matImagePoints.row(row), cv::NORM_L2);
+	}
+	averageNorm /= matImagePoints.rows;
+	// 2) multiply every vectors coordinate by sqrt(2)/avgnorm
+	float scaleFactor = std::sqrt(2)/averageNorm;
+	matImagePoints *= scaleFactor; // add unit homogeneous component here
+	// The points in matImagePoints now have a RMS distance from the origin of sqrt(2).
+	// The normalisation matrix so that the 2D points are mean-free and their norm is as described above.
+	Mat T = Mat::zeros(3, 3, CV_32FC1);
+	T.at<float>(0, 0) = scaleFactor; // s_x
+	T.at<float>(1, 1) = scaleFactor; // s_y
+	T.at<float>(0, 2) = -imagePointsMean.at<float>(0, 0) * scaleFactor; // t_x
+	T.at<float>(1, 2) = -imagePointsMean.at<float>(0, 1) * scaleFactor; // t_y
+	T.at<float>(2, 2) = 1;
+	
+	// center the model points to the origin:
+	Mat tmpOrigMdlPoints = matModelPoints.clone(); // Temp for testing: Save the original coordinates.
+	// translate the centroid of the model points to the origin:
+	Mat modelPointsMean; // use non-homogeneous coords for the next few steps? (less submatrices etc overhead)
+	cv::reduce(matModelPoints, modelPointsMean, 0, CV_REDUCE_AVG);
+	modelPointsMean = cv::repeat(modelPointsMean, matModelPoints.rows, 1);
+	matModelPoints = matModelPoints - modelPointsMean;
+	// scale the model points such that the RMS distance from the origin is sqrt(3):
+	// 1) calculate the average norm (root mean squared distance) of all vectors
+	averageNorm = 0.0f;
+	for (int row = 0; row < matModelPoints.rows; ++row) {
+		averageNorm += cv::norm(matModelPoints.row(row), cv::NORM_L2);
+	}
+	averageNorm /= matModelPoints.rows;
+	// 2) multiply every vectors coordinate by sqrt(3)/avgnorm
+	scaleFactor = std::sqrt(3) / averageNorm;
+	matModelPoints *= scaleFactor; // add unit homogeneous component here
+	// The points in matModelPoints now have a RMS distance from the origin of sqrt(3).
+	// The normalisation matrix so that the 3D points are mean-free and their norm is as described above.
+	Mat U = Mat::zeros(4, 4, CV_32FC1);
+	U.at<float>(0, 0) = scaleFactor; // s_x
+	U.at<float>(1, 1) = scaleFactor; // s_y
+	U.at<float>(2, 2) = scaleFactor; // s_z
+	U.at<float>(0, 3) = -modelPointsMean.at<float>(0, 0) * scaleFactor; // t_x
+	U.at<float>(1, 3) = -modelPointsMean.at<float>(0, 1) * scaleFactor; // t_y
+	U.at<float>(2, 3) = -modelPointsMean.at<float>(0, 2) * scaleFactor; // t_z
+	U.at<float>(3, 3) = 1;
+
+	// Estimate the normalised camera matrix (C tilde).
+	// We are solving the system $A_8 * p_8 = b$
+	// The solution is obtained by the pseudo-inverse of A_8:
+	// $p_8 = A_8^+ * b$
+	Mat A_8 = Mat::zeros(numCorrespondences * 2, 8, CV_32FC1);
+	//Mat p_8(); // p_8 is 8 x 1. We are solving for it.
+	Mat b(numCorrespondences * 2, 1, CV_32FC1);
+	for (int i = 0; i < numCorrespondences; ++i) {
+		A_8.at<float>(2*i, 0) = matModelPoints.at<float>(i, 0); // could maybe made faster by assigning the whole row/col-range if possible?
+		A_8.at<float>(2*i, 1) = matModelPoints.at<float>(i, 1);
+		A_8.at<float>(2*i, 2) = matModelPoints.at<float>(i, 2);
+		A_8.at<float>(2*i, 3) = 1;
+		A_8.at<float>((2*i)+1, 4) = matModelPoints.at<float>(i, 0);
+		A_8.at<float>((2*i)+1, 5) = matModelPoints.at<float>(i, 1);
+		A_8.at<float>((2*i)+1, 6) = matModelPoints.at<float>(i, 2);
+		A_8.at<float>((2*i)+1, 7) = 1;
+		b.at<float>(2*i, 0) = matImagePoints.at<float>(i, 0);
+		b.at<float>((2*i)+1, 0) = matImagePoints.at<float>(i, 1);
+	}
+	Mat p_8 = A_8.inv(cv::DECOMP_SVD) * b;
+	Mat C_tilde = Mat::zeros(3, 4, CV_32FC1);
+	C_tilde.row(0) = p_8.rowRange(0, 4).t(); // The first row of C_tilde consists of the first 4 entries of the col-vector p_8
+	C_tilde.row(1) = p_8.rowRange(4, 8).t(); // Second row = last 4 entries
+	C_tilde.at<float>(2, 3) = 1; // the last row is [0 0 0 1]
+
+	Mat P_Affine = T.inv() * C_tilde * U;
+	return P_Affine;
+}
+
+cv::Vec2f projectAffine(cv::Vec4f vertex, cv::Mat affineCameraMatrix, int screenWidth, int screenHeight)
+{
+	// Transform to clip space:
+	Mat clipCoords = affineCameraMatrix * Mat(vertex);
+	// Take the x and y coordinates in clip space and apply the window transform:
+	cv::Vec2f screenCoords = render::clipToScreenSpace(cv::Vec2f(clipCoords.rowRange(0, 2)), screenWidth, screenHeight);
+	return screenCoords;
+}
+
+	} /* namespace fitting */
+} /* namespace eos */

src/eos/fitting/LinearShapeFitting.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/fitting/LinearShapeFitting.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/fitting/LinearShapeFitting.hpp"
+
+#include "eos/render/utils.hpp"
+
+#include <cassert>
+
+using eos::morphablemodel::MorphableModel;
+using cv::Mat;
+using std::vector;
+
+namespace eos {
+	namespace fitting {
+
+vector<float> fitShapeToLandmarksLinear(MorphableModel morphableModel, Mat affineCameraMatrix, vector<cv::Vec2f> landmarks, std::vector<int> vertexIds, float lambda/*=20.0f*/, boost::optional<int> numCoefficientsToFit/*=boost::optional<int>()*/, boost::optional<float> detectorStandardDeviation/*=boost::optional<float>()*/, boost::optional<float> modelStandardDeviation/*=boost::optional<float>()*/)
+{
+	assert(landmarks.size() == vertexIds.size());
+	
+	int numCoeffsToFit = numCoefficientsToFit.get_value_or(morphableModel.getShapeModel().getNumberOfPrincipalComponents());
+
+	// $\hat{V} \in R^{3N\times m-1}$, subselect the rows of the eigenvector matrix $V$ associated with the $N$ feature points
+	// And we insert a row of zeros after every third row, resulting in matrix $\hat{V}_h \in R^{4N\times m-1}$:
+	Mat V_hat_h = Mat::zeros(4 * landmarks.size(), numCoeffsToFit, CV_32FC1);
+	int rowIndex = 0;
+	for (int i = 0; i < landmarks.size(); ++i) {
+		Mat basisRows = morphableModel.getShapeModel().getNormalisedPcaBasis(vertexIds[i]); // In the paper, the not-normalised basis might be used? I'm not sure, check it. It's even a mess in the paper. PH 26.5.2014: I think the normalised basis is fine/better.
+		//basisRows.copyTo(V_hat_h.rowRange(rowIndex, rowIndex + 3));
+		basisRows.colRange(0, numCoeffsToFit).copyTo(V_hat_h.rowRange(rowIndex, rowIndex + 3));
+		rowIndex += 4; // replace 3 rows and skip the 4th one, it has all zeros
+	}
+	// Form a block diagonal matrix $P \in R^{3N\times 4N}$ in which the camera matrix C (P_Affine, affineCam) is placed on the diagonal:
+	Mat P = Mat::zeros(3 * landmarks.size(), 4 * landmarks.size(), CV_32FC1);
+	for (int i = 0; i < landmarks.size(); ++i) {
+		Mat submatrixToReplace = P.colRange(4 * i, (4 * i) + 4).rowRange(3 * i, (3 * i) + 3);
+		affineCameraMatrix.copyTo(submatrixToReplace);
+	}
+	// The variances: Add the 2D and 3D standard deviations.
+	// If the user doesn't provide them, we choose the following:
+	// 2D (detector) variance: Assuming the detector has a standard deviation of 3 pixels, and the face size (IED) is around 80px. That's 3.75% of the IED. Assume that an image is on average 512x512px so 80/512 = 0.16 is the size the IED occupies inside an image.
+	//                         Now we're in clip-coords ([-1, 1]) and take 0.16 of the range [-1, 1], 0.16/2 = 0.08, and then the standard deviation of the detector is 3.75% of 0.08, i.e. 0.0375*0.08 = 0.003.
+	// 3D (model) variance: 0.0f. It only makes sense to set it to something when we have a different variance for different vertices.
+	float sigma_2D_3D = detectorStandardDeviation.get_value_or(0.003f) + modelStandardDeviation.get_value_or(0.0f);
+	// Note: Isn't it a bit strange to add these as they have different units/normalisations? Check the paper.
+	Mat Sigma = Mat::zeros(3 * landmarks.size(), 3 * landmarks.size(), CV_32FC1);
+	for (int i = 0; i < 3 * landmarks.size(); ++i) {
+		Sigma.at<float>(i, i) = 1.0f / sigma_2D_3D; // the higher the sigma_2D_3D, the smaller the diagonal entries of Sigma will be
+	}
+	Mat Omega = Sigma.t() * Sigma; // just squares the diagonal
+	// The landmarks in matrix notation (in homogeneous coordinates), $3N\times 1$
+	Mat y = Mat::ones(3 * landmarks.size(), 1, CV_32FC1);
+	for (int i = 0; i < landmarks.size(); ++i) {
+		y.at<float>(3 * i, 0) = landmarks[i][0];
+		y.at<float>((3 * i) + 1, 0) = landmarks[i][1];
+		//y.at<float>((3 * i) + 2, 0) = 1; // already 1, stays (homogeneous coordinate)
+	}
+	// The mean, with an added homogeneous coordinate (x_1, y_1, z_1, 1, x_2, ...)^t
+	Mat v_bar = Mat::ones(4 * landmarks.size(), 1, CV_32FC1);
+	for (int i = 0; i < landmarks.size(); ++i) {
+		cv::Vec4f modelMean = morphableModel.getShapeModel().getMeanAtPoint(vertexIds[i]);
+		v_bar.at<float>(4 * i, 0) = modelMean[0];
+		v_bar.at<float>((4 * i) + 1, 0) = modelMean[1];
+		v_bar.at<float>((4 * i) + 2, 0) = modelMean[2];
+		//v_bar.at<float>((4 * i) + 3, 0) = 1; // already 1, stays (homogeneous coordinate)
+		// note: now that a Vec4f is returned, we could use copyTo?
+	}
+
+	// Bring into standard regularised quadratic form with diagonal distance matrix Omega
+	Mat A = P * V_hat_h; // camera matrix times the basis
+	Mat b = P * v_bar - y; // camera matrix times the mean, minus the landmarks.
+	//Mat c_s; // The x, we solve for this! (the variance-normalised shape parameter vector, $c_s = [a_1/sigma_{s,1} , ..., a_m-1/sigma_{s,m-1}]^t$
+	//int numShapePc = morphableModel.getShapeModel().getNumberOfPrincipalComponents();
+	int numShapePc = numCoeffsToFit;
+	Mat AtOmegaA = A.t() * Omega * A;
+	Mat AtOmegaAReg = AtOmegaA + lambda * Mat::eye(numShapePc, numShapePc, CV_32FC1);
+	// Invert using OpenCV:
+	Mat AtOmegaARegInv = AtOmegaAReg.inv(cv::DECOMP_SVD); // DECOMP_SVD calculates the pseudo-inverse if the matrix is not invertible.
+	// Invert (and perform some sanity checks) using Eigen:
+	/*
+	Eigen::Map<Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>> AtOmegaAReg_Eigen(AtOmegaAReg.ptr<float>(), AtOmegaAReg.rows, AtOmegaAReg.cols);
+	Eigen::FullPivLU<Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>> luOfAtOmegaAReg(AtOmegaAReg_Eigen); // Calculate the full-pivoting LU decomposition of the regularized AtA. Note: We could also try FullPivHouseholderQR if our system is non-minimal (i.e. there are more constraints than unknowns).
+	float rankOfAtOmegaAReg = luOfAtOmegaAReg.rank();
+	bool isAtOmegaARegInvertible = luOfAtOmegaAReg.isInvertible();
+	Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> AtARegInv_EigenFullLU = luOfAtOmegaAReg.inverse();
+	Mat AtOmegaARegInvFullLU(AtARegInv_EigenFullLU.rows(), AtARegInv_EigenFullLU.cols(), CV_32FC1, AtARegInv_EigenFullLU.data()); // create an OpenCV Mat header for the Eigen data
+	*/
+	Mat AtOmegatb = A.t() * Omega.t() * b;
+	Mat c_s = -AtOmegaARegInv * AtOmegatb; // Note/Todo: We get coefficients ~ N(0, sigma) I think. They are not multiplied with the eigenvalues.
+	
+	return vector<float>(c_s);
+}
+
+	} /* namespace fitting */
+} /* namespace eos */

src/eos/morphablemodel/MorphableModel.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/morphablemodel/MorphableModel.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/morphablemodel/MorphableModel.hpp"
+
+#include "opencv2/core/core.hpp"
+
+#include <iostream>
+
+using cv::Mat;
+using cv::Vec2f;
+using cv::Vec3f;
+using cv::Vec4f;
+using std::vector;
+using std::array;
+
+namespace {
+/**
+ * Internal helper function that creates a Mesh from given shape and color
+ * PCA instances. Needs the vertex index lists as well to assemble the mesh -
+ * and optional texture coordinates.
+ *
+ * @param[in] shape PCA shape model instance.
+ * @param[in] color PCA color model instance.
+ * @param[in] tvi Triangle vertex indices.
+ * @param[in] tci Triangle color indices (usually identical to the vertex indices).
+ * @param[in] textureCoordinates Optional texture coordinates for each vertex.
+ * @return A mesh created from given parameters.
+ */
+eos::render::Mesh sampleToMesh(Mat shape, Mat color, vector<array<int, 3>> tvi, vector<array<int, 3>> tci, vector<Vec2f> textureCoordinates = vector<Vec2f>())
+{
+	assert(shape.rows == color.rows); // The number of vertices (= model.getDataDimension() / 3) has to be equal for both models.
+	
+	auto numVertices = shape.rows / 3;
+	
+	eos::render::Mesh mesh;
+
+	// Construct the mesh vertices and vertex color information:
+	mesh.vertices.resize(numVertices);
+	mesh.colors.resize(numVertices);
+	for (auto i = 0; i < numVertices; ++i) {
+		mesh.vertices[i] = Vec4f(shape.at<float>(i * 3 + 0), shape.at<float>(i * 3 + 1), shape.at<float>(i * 3 + 2), 1.0f);
+		mesh.colors[i] = Vec3f(color.at<float>(i * 3 + 0), color.at<float>(i * 3 + 1), color.at<float>(i * 3 + 2));        // order in hdf5: RGB. Order in OCV: BGR. But order in vertex.color: RGB
+	}
+
+	// Assign the triangle lists:
+	mesh.tvi = tvi;
+	mesh.tci = tci;
+
+	// Texture coordinates, if the model has them:
+	if (!textureCoordinates.empty()) {
+		mesh.texcoords.resize(numVertices);
+		for (auto i = 0; i < numVertices; ++i) {
+			mesh.texcoords[i] = textureCoordinates[i];
+		}
+	}
+
+	return mesh;
+}
+} /* unnamed namespace */
+
+namespace eos {
+	namespace morphablemodel {
+
+MorphableModel::MorphableModel(PcaModel shapeModel, PcaModel colorModel, std::vector<cv::Vec2f> textureCoordinates /*= std::vector<cv::Vec2f>()*/) : shapeModel(shapeModel), colorModel(colorModel), textureCoordinates(textureCoordinates)
+{
+}
+
+
+PcaModel MorphableModel::getShapeModel() const
+{
+	return shapeModel;
+}
+
+PcaModel MorphableModel::getColorModel() const
+{
+	return colorModel;
+}
+
+render::Mesh MorphableModel::getMean() const
+{
+	assert(shapeModel.getDataDimension() == colorModel.getDataDimension()); // The number of vertices (= model.getDataDimension() / 3) has to be equal for both models.
+
+	Mat shape = shapeModel.getMean();
+	Mat color = colorModel.getMean();
+
+	render::Mesh mesh;
+	if (hasTextureCoordinates()) {
+		mesh = sampleToMesh(shape, color, shapeModel.getTriangleList(), colorModel.getTriangleList(), textureCoordinates);
+	}
+	else {
+		mesh = sampleToMesh(shape, color, shapeModel.getTriangleList(), colorModel.getTriangleList());
+	}
+	return mesh;
+}
+
+render::Mesh MorphableModel::drawSample(float shapeSigma /*= 1.0f*/, float colorSigma /*= 1.0f*/)
+{
+	assert(shapeModel.getDataDimension() == colorModel.getDataDimension()); // The number of vertices (= model.getDataDimension() / 3) has to be equal for both models.
+
+	Mat shapeSample = shapeModel.drawSample(shapeSigma);
+	Mat colorSample = colorModel.drawSample(colorSigma);
+
+	render::Mesh mesh;
+	if (hasTextureCoordinates()) {
+		mesh = sampleToMesh(shapeSample, colorSample, shapeModel.getTriangleList(), colorModel.getTriangleList(), textureCoordinates);
+	}
+	else {
+		mesh = sampleToMesh(shapeSample, colorSample, shapeModel.getTriangleList(), colorModel.getTriangleList());
+	}
+	return mesh;
+}
+
+render::Mesh MorphableModel::drawSample(vector<float> shapeCoefficients, vector<float> colorCoefficients)
+{
+	assert(shapeModel.getDataDimension() == colorModel.getDataDimension()); // The number of vertices (= model.getDataDimension() / 3) has to be equal for both models.
+
+	Mat shapeSample;
+	Mat colorSample;
+
+	if (shapeCoefficients.empty()) {
+		shapeSample = shapeModel.getMean();
+	} else {
+		shapeSample = shapeModel.drawSample(shapeCoefficients);
+	}
+	if (colorCoefficients.empty()) {
+		colorSample = colorModel.getMean();
+	} else {
+		colorSample = colorModel.drawSample(colorCoefficients);
+	}
+
+	render::Mesh mesh;
+	if (hasTextureCoordinates()) {
+		mesh = sampleToMesh(shapeSample, colorSample, shapeModel.getTriangleList(), colorModel.getTriangleList(), textureCoordinates);
+	}
+	else {
+		mesh = sampleToMesh(shapeSample, colorSample, shapeModel.getTriangleList(), colorModel.getTriangleList());
+	}
+	return mesh;
+}
+
+	} /* namespace morphablemodel */
+} /* namespace eos */

src/eos/morphablemodel/PcaModel.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: include/eos/morphablemodel/PcaModel.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/morphablemodel/PcaModel.hpp"
+
+#include "boost/lexical_cast.hpp"
+#include "boost/algorithm/string.hpp"
+
+#include <fstream>
+
+using cv::Mat;
+using cv::Vec4f;
+using boost::lexical_cast;
+using boost::filesystem::path;
+using std::string;
+using std::vector;
+using std::array;
+
+namespace eos {
+	namespace morphablemodel {
+
+PcaModel::PcaModel(Mat mean, Mat pcaBasis, Mat eigenvalues, vector<array<int, 3>> triangleList) : mean(mean), normalisedPcaBasis(pcaBasis), eigenvalues(eigenvalues), triangleList(triangleList)
+{
+	const auto seed = std::random_device()();
+	engine.seed(seed);
+	unnormalisedPcaBasis = unnormalisePcaBasis(normalisedPcaBasis, eigenvalues);
+}
+
+int PcaModel::getNumberOfPrincipalComponents() const
+{
+	// Note: we could assert(normalisedPcaBasis.cols==unnormalisedPcaBasis.cols)
+	return normalisedPcaBasis.cols;
+}
+
+int PcaModel::getDataDimension() const
+{
+	// Note: we could assert(normalisedPcaBasis.rows==unnormalisedPcaBasis.rows)
+	return normalisedPcaBasis.rows;
+}
+
+std::vector<std::array<int, 3>> PcaModel::getTriangleList() const
+{
+	return triangleList;
+}
+
+Mat PcaModel::getMean() const
+{
+	return mean;
+}
+
+Vec4f PcaModel::getMeanAtPoint(int vertexIndex) const
+{
+	vertexIndex *= 3;
+	if (vertexIndex >= mean.rows) {
+		throw std::out_of_range("The given vertex id is larger than the dimension of the mean.");
+	}
+	return Vec4f(mean.at<float>(vertexIndex), mean.at<float>(vertexIndex+1), mean.at<float>(vertexIndex+2), 1.0f);
+}
+
+Mat PcaModel::drawSample(float sigma/*=1.0f*/)
+{
+	std::normal_distribution<float> distribution(0.0f, sigma); // this constructor takes the stddev
+
+	vector<float> alphas(getNumberOfPrincipalComponents());
+
+	for (auto&& a : alphas) {
+		a = distribution(engine);
+	}
+
+	return drawSample(alphas);
+}
+
+Mat PcaModel::drawSample(vector<float> coefficients)
+{
+	// Fill the rest with zeros if not all coefficients are given:
+	if (coefficients.size() < getNumberOfPrincipalComponents()) {
+		coefficients.resize(getNumberOfPrincipalComponents());
+	}
+	Mat alphas(coefficients);
+
+	Mat modelSample = mean + normalisedPcaBasis * alphas;
+
+	return modelSample;
+}
+
+cv::Mat PcaModel::getNormalisedPcaBasis() const
+{
+	return normalisedPcaBasis.clone();
+}
+
+cv::Mat PcaModel::getNormalisedPcaBasis(int vertexId) const
+{
+	vertexId *= 3; // the basis is stored in the format [x y z x y z ...]
+	return normalisedPcaBasis.rowRange(vertexId, vertexId + 3);
+}
+
+cv::Mat PcaModel::getUnnormalisedPcaBasis() const
+{
+	return unnormalisedPcaBasis.clone();
+}
+
+cv::Mat PcaModel::getUnnormalisedPcaBasis(int vertexId) const
+{
+	vertexId *= 3; // the basis is stored in the format [x y z x y z ...]
+	return unnormalisedPcaBasis.rowRange(vertexId, vertexId + 3);
+}
+
+float PcaModel::getEigenvalue(int index) const
+{
+	return eigenvalues.at<float>(index);
+}
+
+cv::Mat normalisePcaBasis(cv::Mat unnormalisedBasis, cv::Mat eigenvalues)
+{
+	Mat normalisedPcaBasis(unnormalisedBasis.size(), unnormalisedBasis.type()); // empty matrix with the same dimensions
+	Mat sqrtOfEigenvalues = eigenvalues.clone();
+	for (int i = 0; i < eigenvalues.rows; ++i)	{
+		sqrtOfEigenvalues.at<float>(i) = std::sqrt(eigenvalues.at<float>(i));
+	}
+	// Normalise the basis: We multiply each eigenvector (i.e. each column) with the square root of its corresponding eigenvalue
+	for (int basis = 0; basis < unnormalisedBasis.cols; ++basis) {
+		Mat normalisedEigenvector = unnormalisedBasis.col(basis).mul(sqrtOfEigenvalues.at<float>(basis));
+		normalisedEigenvector.copyTo(normalisedPcaBasis.col(basis));
+	}
+	
+	return normalisedPcaBasis;
+}
+
+cv::Mat unnormalisePcaBasis(cv::Mat normalisedBasis, cv::Mat eigenvalues)
+{
+	Mat unnormalisedBasis(normalisedBasis.size(), normalisedBasis.type()); // empty matrix with the same dimensions
+	Mat oneOverSqrtOfEigenvalues = eigenvalues.clone();
+	for (int i = 0; i < eigenvalues.rows; ++i)	{
+		oneOverSqrtOfEigenvalues.at<float>(i) = 1.0f / std::sqrt(eigenvalues.at<float>(i));
+	}
+	// De-normalise the basis: We multiply each eigenvector (i.e. each column) with 1 over the square root of its corresponding eigenvalue
+	for (int basis = 0; basis < normalisedBasis.cols; ++basis) {
+		Mat unnormalisedEigenvector = normalisedBasis.col(basis).mul(oneOverSqrtOfEigenvalues.at<float>(basis));
+		unnormalisedEigenvector.copyTo(unnormalisedBasis.col(basis));
+	}
+
+	return unnormalisedBasis;
+}
+
+	} /* namespace morphablemodel */
+} /* namespace eos */

src/eos/morphablemodel/io/cvssp.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/morphablemodel/io/cvssp.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/morphablemodel/io/cvssp.hpp"
+
+#include "boost/lexical_cast.hpp"
+#include "boost/algorithm/string.hpp"
+
+#include <fstream>
+#include <iostream>
+
+namespace fs = boost::filesystem;
+using cv::Mat;
+using cv::Vec2f;
+using cv::Vec3f;
+using boost::lexical_cast;
+using std::string;
+using std::vector;
+using std::array;
+
+namespace eos {
+	namespace morphablemodel {
+
+MorphableModel loadScmModel(fs::path modelFilename, fs::path isomapFile)
+{
+	if (sizeof(unsigned int) != 4) { // note/todo: maybe use uint32 or similar instead? Yep, but still we could encounter endianness-trouble.
+		std::cout << "Warning: We're reading 4 Bytes from the file but sizeof(unsigned int) != 4. Check the code/behaviour." << std::endl;
+	}
+	if (sizeof(double) != 8) {
+		std::cout << "Warning: We're reading 8 Bytes from the file but sizeof(double) != 8. Check the code/behaviour." << std::endl;
+	}
+
+	std::ifstream modelFile(modelFilename.string(), std::ios::binary);
+	if (!modelFile.is_open()) {
+		std::string msg("Unable to open model file: " + modelFilename.string());
+		std::cout << msg << std::endl;
+		throw std::runtime_error(msg);
+	}
+
+	// Reading the shape model
+	// Read (reference?) num triangles and vertices
+	unsigned int numVertices = 0;
+	unsigned int numTriangles = 0;
+	modelFile.read(reinterpret_cast<char*>(&numVertices), 4); // 1 char = 1 byte. uint32=4bytes. float64=8bytes.
+	modelFile.read(reinterpret_cast<char*>(&numTriangles), 4);
+
+	// Read triangles
+	vector<array<int, 3>> triangleList;
+
+	triangleList.resize(numTriangles);
+	unsigned int v0, v1, v2;
+	for (unsigned int i=0; i < numTriangles; ++i) {
+		v0 = v1 = v2 = 0;
+		modelFile.read(reinterpret_cast<char*>(&v0), 4);	// would be nice to pass a &vector and do it in one
+		modelFile.read(reinterpret_cast<char*>(&v1), 4);	// go, but didn't work. Maybe a cv::Mat would work?
+		modelFile.read(reinterpret_cast<char*>(&v2), 4);
+		triangleList[i][0] = v0;
+		triangleList[i][1] = v1;
+		triangleList[i][2] = v2;
+	}
+
+	// Read number of rows and columns of the shape projection matrix (pcaBasis)
+	unsigned int numShapePcaCoeffs = 0;
+	unsigned int numShapeDims = 0;	// dimension of the shape vector (3*numVertices)
+	modelFile.read(reinterpret_cast<char*>(&numShapePcaCoeffs), 4);
+	modelFile.read(reinterpret_cast<char*>(&numShapeDims), 4);
+
+	if (3 * numVertices != numShapeDims) {
+		std::cout << "Warning: Number of shape dimensions is not equal to three times the number of vertices. Something will probably go wrong during the loading." << std::endl;
+	}
+
+	// Read shape projection matrix
+	Mat unnormalisedPcaBasisShape(numShapeDims, numShapePcaCoeffs, CV_32FC1); // m x n (rows x cols) = numShapeDims x numShapePcaCoeffs
+	std::cout << "Loading shape PCA basis matrix with " << unnormalisedPcaBasisShape.rows << " rows and " << unnormalisedPcaBasisShape.cols << " cols." << std::endl;
+	for (unsigned int col = 0; col < numShapePcaCoeffs; ++col) {
+		for (unsigned int row = 0; row < numShapeDims; ++row) {
+			double var = 0.0;
+			modelFile.read(reinterpret_cast<char*>(&var), 8);
+			unnormalisedPcaBasisShape.at<float>(row, col) = static_cast<float>(var);
+		}
+	}
+
+	// Read mean shape vector
+	unsigned int numMean = 0; // dimension of the mean (3*numVertices)
+	modelFile.read(reinterpret_cast<char*>(&numMean), 4);
+	if (numMean != numShapeDims) {
+		std::cout << "Warning: Number of shape dimensions is not equal to the number of dimensions of the mean. Something will probably go wrong during the loading." << std::endl;
+	}
+	Mat meanShape(numMean, 1, CV_32FC1);
+	unsigned int counter = 0;
+	double vd0, vd1, vd2;
+	for (unsigned int i=0; i < numMean/3; ++i) {
+		vd0 = vd1 = vd2 = 0.0;
+		modelFile.read(reinterpret_cast<char*>(&vd0), 8);
+		modelFile.read(reinterpret_cast<char*>(&vd1), 8);
+		modelFile.read(reinterpret_cast<char*>(&vd2), 8);
+		meanShape.at<float>(counter, 0) = static_cast<float>(vd0);
+		++counter;
+		meanShape.at<float>(counter, 0) = static_cast<float>(vd1);
+		++counter;
+		meanShape.at<float>(counter, 0) = static_cast<float>(vd2);
+		++counter;
+	}
+
+	// Read shape eigenvalues
+	unsigned int numEigenValsShape = 0;
+	modelFile.read(reinterpret_cast<char*>(&numEigenValsShape), 4);
+	if (numEigenValsShape != numShapePcaCoeffs) {
+		std::cout << "Warning: Number of coefficients in the PCA basis matrix is not equal to the number of eigenvalues. Something will probably go wrong during the loading." << std::endl;
+	}
+	Mat eigenvaluesShape(numEigenValsShape, 1, CV_32FC1);
+	for (unsigned int i=0; i < numEigenValsShape; ++i) {
+		double var = 0.0;
+		modelFile.read(reinterpret_cast<char*>(&var), 8);
+		eigenvaluesShape.at<float>(i, 0) = static_cast<float>(var);
+	}
+
+	// We read the unnormalised basis from the file. Now let's normalise it and store the normalised basis separately.
+	Mat normalisedPcaBasisShape = normalisePcaBasis(unnormalisedPcaBasisShape, eigenvaluesShape);
+	PcaModel shapeModel(meanShape, normalisedPcaBasisShape, eigenvaluesShape, triangleList);
+
+	// Reading the color model
+	// Read number of rows and columns of projection matrix
+	unsigned int numTexturePcaCoeffs = 0;
+	unsigned int numTextureDims = 0;
+	modelFile.read(reinterpret_cast<char*>(&numTexturePcaCoeffs), 4);
+	modelFile.read(reinterpret_cast<char*>(&numTextureDims), 4);
+	// Read color projection matrix
+	Mat unnormalisedPcaBasisColor(numTextureDims, numTexturePcaCoeffs, CV_32FC1);
+	std::cout << "Loading color PCA basis matrix with " << unnormalisedPcaBasisColor.rows << " rows and " << unnormalisedPcaBasisColor.cols << " cols." << std::endl;
+	for (unsigned int col = 0; col < numTexturePcaCoeffs; ++col) {
+		for (unsigned int row = 0; row < numTextureDims; ++row) {
+			double var = 0.0;
+			modelFile.read(reinterpret_cast<char*>(&var), 8);
+			unnormalisedPcaBasisColor.at<float>(row, col) = static_cast<float>(var);
+		}
+	}
+
+	// Read mean color vector
+	unsigned int numMeanColor = 0; // dimension of the mean (3*numVertices)
+	modelFile.read(reinterpret_cast<char*>(&numMeanColor), 4);
+	Mat meanColor(numMeanColor, 1, CV_32FC1);
+	counter = 0;
+	for (unsigned int i=0; i < numMeanColor/3; ++i) {
+		vd0 = vd1 = vd2 = 0.0;
+		modelFile.read(reinterpret_cast<char*>(&vd0), 8); // order in hdf5: RGB. Order in OCV: BGR. But order in vertex.color: RGB
+		modelFile.read(reinterpret_cast<char*>(&vd1), 8);
+		modelFile.read(reinterpret_cast<char*>(&vd2), 8);
+		meanColor.at<float>(counter, 0) = static_cast<float>(vd0);
+		++counter;
+		meanColor.at<float>(counter, 0) = static_cast<float>(vd1);
+		++counter;
+		meanColor.at<float>(counter, 0) = static_cast<float>(vd2);
+		++counter;
+	}
+
+	// Read color eigenvalues
+	unsigned int numEigenValsColor = 0;
+	modelFile.read(reinterpret_cast<char*>(&numEigenValsColor), 4);
+	Mat eigenvaluesColor(numEigenValsColor, 1, CV_32FC1);
+	for (unsigned int i=0; i < numEigenValsColor; ++i) {
+		double var = 0.0;
+		modelFile.read(reinterpret_cast<char*>(&var), 8);
+		eigenvaluesColor.at<float>(i, 0) = static_cast<float>(var);
+	}
+
+	// We read the unnormalised basis from the file. Now let's normalise it and store the normalised basis separately.
+	Mat normalisedPcaBasisColor = normalisePcaBasis(unnormalisedPcaBasisColor, eigenvaluesColor);
+	PcaModel colorModel(meanColor, normalisedPcaBasisColor, eigenvaluesColor, triangleList);
+
+	modelFile.close();
+
+	// Load the isomap with texture coordinates if a filename has been given:
+	vector<Vec2f> texCoords;
+	if (!isomapFile.empty()) {
+		texCoords = loadIsomap(isomapFile);
+		if (shapeModel.getDataDimension() / 3.0f != texCoords.size()) {
+			std::string errorMessage("Error, wrong number of texture coordinates. Don't have the same number of texcoords than the shape model has vertices.");
+			std::cout << errorMessage << std::endl;
+			throw std::runtime_error(errorMessage);
+		}
+	}
+	
+	return MorphableModel(shapeModel, colorModel, texCoords);
+}
+
+vector<Vec2f> loadIsomap(fs::path isomapFile)
+{
+	vector<float> xCoords, yCoords;
+	string line;
+	std::ifstream myfile(isomapFile.string());
+	if (!myfile.is_open()) {
+		string logMessage("The Isomap file could not be opened. Did you specify a correct filename? " + isomapFile.string());
+		std::cout << logMessage << std::endl;
+		throw std::runtime_error(logMessage);
+	}
+	else {
+		while (getline(myfile, line))
+		{
+			std::istringstream iss(line);
+			string x, y;
+			iss >> x >> y;
+			xCoords.push_back(lexical_cast<float>(x));
+			yCoords.push_back(lexical_cast<float>(y));
+		}
+		myfile.close();
+	}
+	// Process the coordinates: Find the min/max and rescale to [0, 1] x [0, 1]
+	auto minMaxX = std::minmax_element(begin(xCoords), end(xCoords)); // minMaxX is a pair, first=min, second=max
+	auto minMaxY = std::minmax_element(begin(yCoords), end(yCoords));
+
+	vector<Vec2f> texCoords;
+	float divisorX = *minMaxX.second - *minMaxX.first;
+	float divisorY = *minMaxY.second - *minMaxY.first;
+	for (int i = 0; i < xCoords.size(); ++i) {
+		texCoords.push_back(Vec2f((xCoords[i] - *minMaxX.first) / divisorX, 1.0f - (yCoords[i] - *minMaxY.first) / divisorY)); // We rescale to [0, 1] and at the same time flip the y-coords (because in the isomap, the coordinates are stored upside-down).
+	}
+
+	return texCoords;
+}
+
+	} /* namespace morphablemodel */
+} /* namespace eos */

src/eos/render/Mesh.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/render/Mesh.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/render/Mesh.hpp"
+
+#include <cassert>
+#include <fstream>
+
+namespace eos {
+	namespace render {
+
+void writeObj(Mesh mesh, std::string filename)
+{
+	std::ofstream objFile(filename);
+
+	assert(mesh.vertices.size() == mesh.colors.size());
+	
+	for (std::size_t i = 0; i < mesh.vertices.size(); ++i) {
+		objFile << "v " << mesh.vertices[i][0] << " " << mesh.vertices[i][1] << " " << mesh.vertices[i][2] << " " << mesh.colors[i][0] << " " << mesh.colors[i][1] << " " << mesh.colors[i][2] << " " << std::endl;
+	}
+
+	for (auto&& v : mesh.tvi) {
+		// Add one because obj starts counting triangle indices at 1
+		objFile << "f " << v[0] + 1 << " " << v[1] + 1 << " " << v[2] + 1 << std::endl;
+	}
+
+	return;
+}
+
+	} /* namespace render */
+} /* namespace eos */

src/eos/render/utils.cpp

+/*
+ * Eos - A 3D Morphable Model fitting library written in modern C++11/14.
+ *
+ * File: src/eos/render/utils.cpp
+ *
+ * Copyright 2014, 2015 Patrik Huber
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+#include "eos/render/utils.hpp"
+
+using cv::Vec2f;
+
+namespace eos {
+	namespace render {
+
+Vec2f clipToScreenSpace(Vec2f clipCoordinates, int screenWidth, int screenHeight)
+{
+	// Window transform:
+	float x_ss = (clipCoordinates[0] + 1.0f) * (screenWidth / 2.0f);
+	float y_ss = screenHeight - (clipCoordinates[1] + 1.0f) * (screenHeight / 2.0f); // also flip y; Qt: Origin top-left. OpenGL: bottom-left.
+	return Vec2f(x_ss, y_ss);
+	/* Note: What we do here is equivalent to
+	   x_w = (x *  vW/2) + vW/2;
+	   However, Shirley says we should do:
+	   x_w = (x *  vW/2) + (vW-1)/2;
+	   (analogous for y)
+	   Todo: Check the consequences.
+	*/
+}
+
+Vec2f screenToClipSpace(Vec2f screenCoordinates, int screenWidth, int screenHeight)
+{
+	float x_cs = screenCoordinates[0] / (screenWidth / 2.0f) - 1.0f;
+	float y_cs = screenCoordinates[1] / (screenHeight / 2.0f) - 1.0f;
+	y_cs *= -1.0f;
+	return Vec2f(x_cs, y_cs);
+}
+
+	} /* namespace render */
+} /* namespace eos */