Merge remote-tracking branch 'upstream/master' into parallel_vec

Conflicts:
	astra_vc11.vcxproj.filters
	python/astra/functions.py
	src/ParallelBeamBlobKernelProjector2D.cpp
	src/ProjectionGeometry2D.cpp

15 files changed, 1032 insertions, 169 deletions

diff --git a/include/astra/AstraObjectFactory.h b/include/astra/AstraObjectFactory.h
index ba4ec11..1ed4955 100644
--- a/include/astra/AstraObjectFactory.h
+++ b/include/astra/AstraObjectFactory.h
@@ -110,7 +110,7 @@ template <typename T, typename TypeList>
 T* CAstraObjectFactory<T, TypeList>::create(const Config& _cfg)
 {
 	functor_find<T> finder = functor_find<T>();
-	finder.tofind = _cfg.self->getAttribute("type");
+	finder.tofind = _cfg.self.getAttribute("type");
 	CreateObject<TypeList>::find(finder);
 	if (finder.res == NULL) return NULL;
 	if (finder.res->initialize(_cfg))
diff --git a/include/astra/Config.h b/include/astra/Config.h
index 9893c90..c10a16e 100644
--- a/include/astra/Config.h
+++ b/include/astra/Config.h
@@ -44,13 +44,14 @@ namespace astra {
 struct _AstraExport Config {
 
 	Config();
-	Config(XMLNode* _self);
+	Config(XMLNode _self);
 	~Config();
 
 	void initialize(std::string rootname);
 
-	XMLNode* self;
-	XMLNode* global;
+	XMLNode self;
+
+	XMLDocument *_doc;
 };
 
 struct ConfigCheckData {
diff --git a/include/astra/Float32ProjectionData2D.h b/include/astra/Float32ProjectionData2D.h
index 7461491..bb99f4b 100644
--- a/include/astra/Float32ProjectionData2D.h
+++ b/include/astra/Float32ProjectionData2D.h
@@ -101,6 +101,19 @@ public:
 	 * Copy constructor
 	 */
 	CFloat32ProjectionData2D(const CFloat32ProjectionData2D& _other);
+	
+	/** Constructor. Create an instance of the CFloat32ProjectionData2D class with pre-allocated memory.
+	 *
+	 * Creates an instance of the CFloat32ProjectionData2D class. Memory 
+	 * is pre-allocated and passed via the abstract CFloat32CustomMemory handle
+	 * class. The handle will be deleted when the memory can be freed.
+	 * You should override the destructor to provide custom behaviour on free.
+	 *
+	 * @param _pGeometry Projection Geometry object.  This object will be HARDCOPIED into this class.
+	 * @param _pCustomMemory custom memory handle
+	 *
+	 */
+	CFloat32ProjectionData2D(CProjectionGeometry2D* _pGeometry, CFloat32CustomMemory* _pCustomMemory);
 
 	/**
 	 * Assignment operator
@@ -148,6 +161,18 @@ public:
 	 * @param _fScalar scalar value to be put at each index.
 	 */
 	bool initialize(CProjectionGeometry2D* _pGeometry, float32 _fScalar);
+	
+	/** Initialization. Initializes an instance of the CFloat32ProjectionData2D class with pre-allocated memory.
+	 *
+	 * Memory is pre-allocated and passed via the abstract CFloat32CustomMemory handle
+	 * class. The handle will be deleted when the memory can be freed.
+	 * You should override the destructor to provide custom behaviour on free.
+	 *
+	 * @param _pGeometry Projection Geometry object.  This object will be HARDCOPIED into this class.
+	 * @param _pCustomMemory custom memory handle
+	 *
+	 */
+	bool initialize(CProjectionGeometry2D* _pGeometry, CFloat32CustomMemory* _pCustomMemory);
 
 	/** Get the number of detectors.
 	 *
diff --git a/include/astra/Float32ProjectionData3D.h b/include/astra/Float32ProjectionData3D.h
index 79b762e..329c9a4 100644
--- a/include/astra/Float32ProjectionData3D.h
+++ b/include/astra/Float32ProjectionData3D.h
@@ -196,6 +196,11 @@ public:
 	 * @return pointer to projection geometry.
 	 */
 	virtual CProjectionGeometry3D* getGeometry() const;
+
+	/** Change the projection geometry.
+	 *  Note that this can't change the dimensions of the data.
+	 */
+	virtual void changeGeometry(CProjectionGeometry3D* pGeometry);
 };
 
 
diff --git a/include/astra/Float32VolumeData2D.h b/include/astra/Float32VolumeData2D.h
index 4f44a8c..abecebf 100644
--- a/include/astra/Float32VolumeData2D.h
+++ b/include/astra/Float32VolumeData2D.h
@@ -92,6 +92,19 @@ public:
 	 * Copy constructor
 	 */
 	CFloat32VolumeData2D(const CFloat32VolumeData2D& _other);
+	
+	/** Constructor. Create an instance of the CFloat32VolumeData2D class with pre-allocated memory.
+	 *
+	 * Creates an instance of the CFloat32VolumeData2D class. Memory 
+	 * is pre-allocated and passed via the abstract CFloat32CustomMemory handle
+	 * class. The handle will be deleted when the memory can be freed.
+	 * You should override the destructor to provide custom behaviour on free.
+	 *
+	 * @param _pGeometry Volume Geometry object.  This object will be HARDCOPIED into this class.
+	 * @param _pCustomMemory custom memory handle
+	 *
+	 */
+	CFloat32VolumeData2D(CVolumeGeometry2D* _pGeometry, CFloat32CustomMemory* _pCustomMemory);
 
 	/**
 	 * Assignment operator
@@ -132,6 +145,18 @@ public:
 	 * @param _fScalar scalar value to be put at each index.
 	 */
 	bool initialize(CVolumeGeometry2D* _pGeometry, float32 _fScalar);
+	
+	/** Initialization. Initializes an instance of the CFloat32VolumeData2D class with pre-allocated memory.
+	 *
+	 * Memory is pre-allocated and passed via the abstract CFloat32CustomMemory handle
+	 * class. The handle will be deleted when the memory can be freed.
+	 * You should override the destructor to provide custom behaviour on free.
+	 *
+	 * @param _pGeometry Volume Geometry object.  This object will be HARDCOPIED into this class.
+	 * @param _pCustomMemory custom memory handle
+	 *
+	 */
+	bool initialize(CVolumeGeometry2D* _pGeometry, CFloat32CustomMemory* _pCustomMemory);
 
 	/** Destructor.
 	 */
diff --git a/include/astra/Float32VolumeData3D.h b/include/astra/Float32VolumeData3D.h
index d8f0ae9..07df78f 100644
--- a/include/astra/Float32VolumeData3D.h
+++ b/include/astra/Float32VolumeData3D.h
@@ -214,6 +214,11 @@ public:
 	 * @return The geometry describing the data stored in this volume
 	 */
 	virtual CVolumeGeometry3D* getGeometry() const;
+
+	/** Change the projection geometry.
+	 *  Note that this can't change the dimensions of the data.
+	 */
+	virtual void changeGeometry(CVolumeGeometry3D* pGeometry);
 };
 
 //----------------------------------------------------------------------------------------
diff --git a/include/astra/Logger.h b/include/astra/Logger.h
deleted file mode 100644
index 34fd364..0000000
--- a/include/astra/Logger.h
+++ /dev/null
@@ -1,72 +0,0 @@
-/*
------------------------------------------------------------------------
-Copyright: 2010-2015, iMinds-Vision Lab, University of Antwerp
-           2014-2015, CWI, Amsterdam
-
-Contact: astra@uantwerpen.be
-Website: http://sf.net/projects/astra-toolbox
-
-This file is part of the ASTRA Toolbox.
-
-
-The ASTRA Toolbox is free software: you can redistribute it and/or modify
-it under the terms of the GNU General Public License as published by
-the Free Software Foundation, either version 3 of the License, or
-(at your option) any later version.
-
-The ASTRA Toolbox is distributed in the hope that it will be useful,
-but WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-GNU General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
-
------------------------------------------------------------------------
-$Id$
-*/
-
-#ifndef _INC_ASTRA_LOGGER
-#define _INC_ASTRA_LOGGER
-
-#include <cstdio>
-
-namespace astra
-{
-
-/**
- *  This is the first stab at a decent logger. If the file "astra_logger.txt", it will be replaced
- *  with the text sent to this logger. If the file doesn't exist when the app starts, nothing is written.
- */
-class CLogger
-{
-	static std::FILE * m_pOutFile;
-	static bool m_bInitialized;
-
-	static void _assureIsInitialized();
-
-	CLogger();
-
-public:
-
-	/**
-	 * Writes a line to the log file (newline is added). Ignored if logging is turned off.
-	 *
-	 * @param _text char pointer to text in line
-	 */
-	static void writeLine(const char * _text);
-
-	/**
-	 * Formats and writes a CUDA error to the log file. Ignored if logging is turned off.
-	 *
-	 * @param _fileName filename where error occurred (typically __FILE__)
-	 * @param _line line in file (typically __LINE__)
-	 * @param _errString string describing the error, can be output of cudaGetErrorString
-	 */
-	static void writeTerminalCUDAError(const char * _fileName, int _iLine, const char * _errString);
-};
-
-}
-
-#endif /* _INC_ASTRA_LOGGER */
-
diff --git a/include/astra/Logging.h b/include/astra/Logging.h
new file mode 100644
index 0000000..8e19ea4
--- /dev/null
+++ b/include/astra/Logging.h
@@ -0,0 +1,159 @@
+/*
+-----------------------------------------------------------------------
+Copyright: 2010-2015, iMinds-Vision Lab, University of Antwerp
+           2014-2015, CWI, Amsterdam
+
+Contact: astra@uantwerpen.be
+Website: http://sf.net/projects/astra-toolbox
+
+This file is part of the ASTRA Toolbox.
+
+
+The ASTRA Toolbox is free software: you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+The ASTRA Toolbox is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
+
+-----------------------------------------------------------------------
+$Id$
+*/
+
+#ifndef _INC_ASTRA_LOGGING
+#define _INC_ASTRA_LOGGING
+
+#include "astra/Globals.h"
+
+#define ASTRA_DEBUG(...) astra::CLogger::debug(__FILE__,__LINE__, __VA_ARGS__)
+#define ASTRA_INFO(...) astra::CLogger::info(__FILE__,__LINE__, __VA_ARGS__)
+#define ASTRA_WARN(...) astra::CLogger::warn(__FILE__,__LINE__, __VA_ARGS__)
+#define ASTRA_ERROR(...) astra::CLogger::error(__FILE__,__LINE__, __VA_ARGS__)
+
+namespace astra
+{
+
+enum log_level {
+    LOG_DEBUG,
+    LOG_INFO,
+    LOG_WARN,
+    LOG_ERROR
+};
+
+class _AstraExport CLogger
+{
+	CLogger();
+  ~CLogger();
+  static bool m_bEnabledFile;
+  static bool m_bEnabledScreen;
+  static bool m_bFileProvided;
+  static bool m_bInitialized;
+  static void _assureIsInitialized();
+  static void _setLevel(int id, log_level m_eLevel);
+
+public:
+
+	/**
+	 * Writes a line to the log file (newline is added). Ignored if logging is turned off.
+	 *
+	 * @param sfile
+   * The name of the source file making this log call (e.g. __FILE__).
+   *
+   * @param sline
+   * The line number of the call in the source code (e.g. __LINE__).
+   *
+   * @param id
+   * The id of the logger to write to.
+   *
+   * @param fmt
+   * The format string for the message (printf formatting).
+   *
+   * @param ...
+   * Any additional format arguments.
+	 */
+  static void debug(const char *sfile, int sline, const char *fmt, ...);
+  static void info(const char *sfile, int sline, const char *fmt, ...);
+  static void warn(const char *sfile, int sline, const char *fmt, ...);
+  static void error(const char *sfile, int sline, const char *fmt, ...);
+
+  /**
+	 * Sets the file to log to, with logging level.
+   *
+   * @param filename
+   * File to log to.
+	 *
+	 * @param m_eLevel
+   * Logging level (LOG_DEBUG, LOG_WARN, LOG_INFO, LOG_ERROR).
+   *
+	 */
+  static void setOutputFile(const char *filename, log_level m_eLevel);
+
+  /**
+	 * Sets the screen to log to, with logging level.
+   *
+   * @param screen_fd
+   * Screen file descriptor (1 for stdout, 2 for stderr)
+	 *
+	 * @param m_eLevel
+   * Logging level (LOG_DEBUG, LOG_WARN, LOG_INFO, LOG_ERROR).
+   *
+	 */
+  static void setOutputScreen(int fd, log_level m_eLevel);
+  
+  /**
+   * Set the format string for log messages.  Here are the substitutions you may
+   * use:
+   *
+   *     %f: Source file name generating the log call.
+   *     %n: Source line number where the log call was made.
+   *     %m: The message text sent to the logger (after printf formatting).
+   *     %d: The current date, formatted using the logger's date format.
+   *     %t: The current time, formatted using the logger's time format.
+   *     %l: The log level (one of "DEBUG", "INFO", "WARN", or "ERROR").
+   *     %%: A literal percent sign.
+   *
+   * The default format string is "%d %t %f(%n): %l: %m\n".
+   *
+   * @param fmt
+   * The new format string, which must be less than 256 bytes.
+   * You probably will want to end this with a newline (\n).
+   *
+   */
+  static void setFormatFile(const char *fmt);
+  static void setFormatScreen(const char *fmt);
+
+
+  /**
+   * Enable logging.
+   *
+   */
+  static void enable();
+  static void enableScreen();
+  static void enableFile();
+
+  /**
+   * Disable logging.
+   *
+   */
+  static void disable();
+  static void disableScreen();
+  static void disableFile();
+
+  /**
+   * Set callback function for logging to screen.
+   * @return whether callback was set succesfully.
+   *
+   */
+  static bool setCallbackScreen(void (*cb)(const char *msg, size_t len));
+
+};
+
+}
+
+#endif /* _INC_ASTRA_LOGGING */
diff --git a/include/astra/ParallelBeamBlobKernelProjector2D.inl b/include/astra/ParallelBeamBlobKernelProjector2D.inl
index 203eb51..1945cdd 100644
--- a/include/astra/ParallelBeamBlobKernelProjector2D.inl
+++ b/include/astra/ParallelBeamBlobKernelProjector2D.inl
@@ -31,21 +31,21 @@ template <typename Policy>
 void CParallelBeamBlobKernelProjector2D::project(Policy& p)
 {
 	projectBlock_internal(0, m_pProjectionGeometry->getProjectionAngleCount(),
-		                  0, m_pProjectionGeometry->getDetectorCount(), p);
+						  0, m_pProjectionGeometry->getDetectorCount(), p);
 }
 
 template <typename Policy>
 void CParallelBeamBlobKernelProjector2D::projectSingleProjection(int _iProjection, Policy& p)
 {
 	projectBlock_internal(_iProjection, _iProjection + 1,
-	                      0, m_pProjectionGeometry->getDetectorCount(), p);
+						  0, m_pProjectionGeometry->getDetectorCount(), p);
 }
 
 template <typename Policy>
 void CParallelBeamBlobKernelProjector2D::projectSingleRay(int _iProjection, int _iDetector, Policy& p)
 {
 	projectBlock_internal(_iProjection, _iProjection + 1,
-	                      _iDetector, _iDetector + 1, p);
+						  _iDetector, _iDetector + 1, p);
 }
 
 //----------------------------------------------------------------------------------------
@@ -53,12 +53,6 @@ void CParallelBeamBlobKernelProjector2D::projectSingleRay(int _iProjection, int
 template <typename Policy>
 void CParallelBeamBlobKernelProjector2D::projectBlock_internal(int _iProjFrom, int _iProjTo, int _iDetFrom, int _iDetTo, Policy& p)
 {
-	// variables
-	float32 detX, detY, x, y, c, r, update_c, update_r, offset, invBlobExtent, inv_pixelLengthX, inv_pixelLengthY, weight, d;
-	int iVolumeIndex, iRayIndex, row, col, iAngle, iDetector, colCount, rowCount, detCount;
-	int col_left, col_right, row_top, row_bottom, index;
-	const SParProjection * proj = 0;
-
 	// get vector geometry
 	const CParallelVecProjectionGeometry2D* pVecProjectionGeometry;
 	if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) {
@@ -68,16 +62,22 @@ void CParallelBeamBlobKernelProjector2D::projectBlock_internal(int _iProjFrom, i
 	}
 
 	// precomputations
-	inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
-	inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
-	colCount = m_pVolumeGeometry->getGridColCount();
-	rowCount = m_pVolumeGeometry->getGridRowCount();
-	detCount = pVecProjectionGeometry->getDetectorCount();
+	const float32 inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
+	const float32 inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
+	const int colCount = m_pVolumeGeometry->getGridColCount();
+	const int rowCount = m_pVolumeGeometry->getGridRowCount();
+	const int detCount = pVecProjectionGeometry->getDetectorCount();
 
 	// loop angles
-	for (iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+	#pragma omp parallel for
+	for (int iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+
+		// variables
+		float32 detX, detY, x, y, c, r, update_c, update_r, offset, invBlobExtent, weight, d;
+		int iVolumeIndex, iRayIndex, row, col, iDetector;
+		int col_left, col_right, row_top, row_bottom, index;
 
-		proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
+		const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
 
 		bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY);
 		if (vertical) {
diff --git a/include/astra/ParallelBeamLineKernelProjector2D.inl b/include/astra/ParallelBeamLineKernelProjector2D.inl
index dc3153d..9c58d60 100644
--- a/include/astra/ParallelBeamLineKernelProjector2D.inl
+++ b/include/astra/ParallelBeamLineKernelProjector2D.inl
@@ -55,12 +55,6 @@ void CParallelBeamLineKernelProjector2D::projectSingleRay(int _iProjection, int
 template <typename Policy>
 void CParallelBeamLineKernelProjector2D::projectBlock_internal(int _iProjFrom, int _iProjTo, int _iDetFrom, int _iDetTo, Policy& p)
 {
-	// variables
-	float32 detX, detY, S, T, I, x, y, c, r, update_c, update_r, offset;
-	float32 lengthPerRow, lengthPerCol, inv_pixelLengthX, inv_pixelLengthY, invTminSTimesLengthPerRow, invTminSTimesLengthPerCol;
-	int iVolumeIndex, iRayIndex, row, col, iAngle, iDetector, colCount, rowCount, detCount;
-	const SParProjection * proj = 0;
-
 	// get vector geometry
 	const CParallelVecProjectionGeometry2D* pVecProjectionGeometry;
 	if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) {
@@ -70,16 +64,22 @@ void CParallelBeamLineKernelProjector2D::projectBlock_internal(int _iProjFrom, i
 	}
 
 	// precomputations
-	inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
-	inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
-	colCount = m_pVolumeGeometry->getGridColCount();
-	rowCount = m_pVolumeGeometry->getGridRowCount();
-	detCount = pVecProjectionGeometry->getDetectorCount();
+	const float32 inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
+	const float32 inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
+	const int colCount = m_pVolumeGeometry->getGridColCount();
+	const int rowCount = m_pVolumeGeometry->getGridRowCount();
+	const int detCount = pVecProjectionGeometry->getDetectorCount();
 
 	// loop angles
-	for (iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+	#pragma omp parallel for
+	for (int iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+
+		// variables
+		float32 detX, detY, S, T, I, x, y, c, r, update_c, update_r, offset;
+		float32 lengthPerRow, lengthPerCol, invTminSTimesLengthPerRow, invTminSTimesLengthPerCol;
+		int iVolumeIndex, iRayIndex, row, col, iDetector;
 
-		proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
+		const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
 
 		bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY);
 		if (vertical) {
diff --git a/include/astra/ParallelBeamLinearKernelProjector2D.inl b/include/astra/ParallelBeamLinearKernelProjector2D.inl
index 9f898ce..79b82d4 100644
--- a/include/astra/ParallelBeamLinearKernelProjector2D.inl
+++ b/include/astra/ParallelBeamLinearKernelProjector2D.inl
@@ -54,12 +54,6 @@ void CParallelBeamLinearKernelProjector2D::projectSingleRay(int _iProjection, in
 template <typename Policy>
 void CParallelBeamLinearKernelProjector2D::projectBlock_internal(int _iProjFrom, int _iProjTo, int _iDetFrom, int _iDetTo, Policy& p)
 {
-	// variables
-	float32 detX, detY, x, y, c, r, update_c, update_r, offset;
-	float32 lengthPerRow, lengthPerCol, inv_pixelLengthX, inv_pixelLengthY;
-	int iVolumeIndex, iRayIndex, row, col, iAngle, iDetector, colCount, rowCount, detCount;
-	const SParProjection * proj = 0;
-
 	// get vector geometry
 	const CParallelVecProjectionGeometry2D* pVecProjectionGeometry;
 	if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) {
@@ -69,16 +63,22 @@ void CParallelBeamLinearKernelProjector2D::projectBlock_internal(int _iProjFrom,
 	}
 
 	// precomputations
-	inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
-	inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
-	colCount = m_pVolumeGeometry->getGridColCount();
-	rowCount = m_pVolumeGeometry->getGridRowCount();
-	detCount = pVecProjectionGeometry->getDetectorCount();
+	const float32 inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
+	const float32 inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
+	const int colCount = m_pVolumeGeometry->getGridColCount();
+	const int rowCount = m_pVolumeGeometry->getGridRowCount();
+	const int detCount = pVecProjectionGeometry->getDetectorCount();
 
 	// loop angles
-	for (iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+	#pragma omp parallel for
+	for (int iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+
+		// variables
+		float32 detX, detY, x, y, c, r, update_c, update_r, offset;
+		float32 lengthPerRow, lengthPerCol;
+		int iVolumeIndex, iRayIndex, row, col, iDetector;
 
-		proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
+		const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
 
 		bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY);
 		if (vertical) {
diff --git a/include/astra/ParallelBeamStripKernelProjector2D.inl b/include/astra/ParallelBeamStripKernelProjector2D.inl
index 31d222a..85faaa3 100644
--- a/include/astra/ParallelBeamStripKernelProjector2D.inl
+++ b/include/astra/ParallelBeamStripKernelProjector2D.inl
@@ -52,13 +52,6 @@ void CParallelBeamStripKernelProjector2D::projectSingleRay(int _iProjection, int
 template <typename Policy>
 void CParallelBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom, int _iProjTo, int _iDetFrom, int _iDetTo, Policy& p)
 {
-	// variables
-	float32 detLX, detLY, detRX, detRY, S, T, update_c, update_r, offsetL, offsetR, invTminS;
-	float32 inv_pixelLengthX, inv_pixelLengthY, pixelArea, res, fRxOverRy, fRyOverRx;
-	int iVolumeIndex, iRayIndex, iAngle, iDetector;
-	int row, row_top, row_bottom, col, col_left, col_right, colCount, rowCount, detCount;
-	const SParProjection * proj = 0;
-	
 	// get vector geometry
 	const CParallelVecProjectionGeometry2D* pVecProjectionGeometry;
 	if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) {
@@ -68,17 +61,24 @@ void CParallelBeamStripKernelProjector2D::projectBlock_internal(int _iProjFrom,
 	}
 
 	// precomputations
-	inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
-	inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
-	pixelArea = m_pVolumeGeometry->getPixelLengthX() * m_pVolumeGeometry->getPixelLengthY();
-	colCount = m_pVolumeGeometry->getGridColCount();
-	rowCount = m_pVolumeGeometry->getGridRowCount();
-	detCount = pVecProjectionGeometry->getDetectorCount();
+	const float32 inv_pixelLengthX = 1.0f / m_pVolumeGeometry->getPixelLengthX();
+	const float32 inv_pixelLengthY = 1.0f / m_pVolumeGeometry->getPixelLengthY();
+	const float32 pixelArea = m_pVolumeGeometry->getPixelLengthX() * m_pVolumeGeometry->getPixelLengthY();
+	const int colCount = m_pVolumeGeometry->getGridColCount();
+	const int rowCount = m_pVolumeGeometry->getGridRowCount();
+	const int detCount = pVecProjectionGeometry->getDetectorCount();
 
 	// loop angles
-	for (iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+	#pragma omp parallel for
+	for (int iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) {
+
+		// variables
+		float32 detLX, detLY, detRX, detRY, S, T, update_c, update_r, offsetL, offsetR, invTminS;
+		float32 res, fRxOverRy, fRyOverRx;
+		int iVolumeIndex, iRayIndex, iDetector;
+		int row, row_top, row_bottom, col, col_left, col_right;
 
-		proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
+		const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle];
 
 		bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY);
 		if (vertical) {
diff --git a/include/astra/XMLDocument.h b/include/astra/XMLDocument.h
index 869e1a3..eddd908 100644
--- a/include/astra/XMLDocument.h
+++ b/include/astra/XMLDocument.h
@@ -78,7 +78,7 @@ public:
 	 *
 	 * @return first XML node of the document
 	 */
-	XMLNode* getRootNode();
+	XMLNode getRootNode();
 
 	/** Save an XML DOM tree to an XML file
 	 *
diff --git a/include/astra/XMLNode.h b/include/astra/XMLNode.h
index eceffe1..4d29d5c 100644
--- a/include/astra/XMLNode.h
+++ b/include/astra/XMLNode.h
@@ -64,71 +64,68 @@ public:
 	/** Deconstructor
 	 */
 	~XMLNode();
-	
+
+	/** Check validity
+	 */
+	operator bool() const { return fDOMElement != 0; }
 
 	/** Get a single child XML node. If there are more, the first one is returned
 	 *
 	 * @param _sName tagname of the requested child node
 	 * @return first child node with the correct tagname, null pointer if it doesn't exist
 	 */
-	XMLNode* getSingleNode(string _sName);
+	XMLNode getSingleNode(string _sName) const;
 
 	/** Get all child XML nodes that have the tagname name
 	 * 
 	 * @param _sName tagname of the requested child nodes
 	 * @return list with all child nodes with the correct tagname
 	 */
-	std::list<XMLNode*> getNodes(string _sName);
+	std::list<XMLNode> getNodes(string _sName) const;
 	
 	/** Get all child XML nodes
 	 * 
 	 * @return list with all child nodes 
 	 */
-	std::list<XMLNode*> getNodes();
+	std::list<XMLNode> getNodes() const;
 
 	/** Get the name of this node
 	 * 
 	 * @return name of node
 	 */
-	std::string getName();
+	std::string getName() const;
 
 	/** Get the content of the XML node as a single string.
 	 *
 	 * @return node content
 	 */ 
-	string getContent();
+	string getContent() const;
 
 	/** Get the content of the XML node as a numerical.
 	 *
 	 * @return node content
 	 */ 
-	float32 getContentNumerical();
+	float32 getContentNumerical() const;
 
 	/** Get the content of the XML node as a boolean.
 	 *
 	 * @return node content
 	 */ 
-	bool getContentBool();
+	bool getContentBool() const;
 
 	/** Get the content of the XML node as a vector of strings.
 	 *
 	 * @return node content
 	 */ 
-	vector<string> getContentArray();
-
-	/** Get the content of the XML node as a c-array of float32 data.
-	 *
-	 * @param _pfData data array, shouldn't be initialized already.
-	 * @param _iSize number of elements stored in _pfData
-	 */ 
-	void getContentNumericalArray(float32*& _pfData, int& _iSize);
+	vector<string> getContentArray() const;
 
 	/** Get the content of the XML node as a stl container of float32 data.
+	 *  NB: A 2D matrix is returned as a linear list
 	 *
 	 * @return node content
 	 */ 
-	vector<float32> getContentNumericalArray();
-	vector<double> getContentNumericalArrayDouble();
+	vector<float32> getContentNumericalArray() const;
+	vector<double> getContentNumericalArrayDouble() const;
 
 
 
@@ -137,7 +134,7 @@ public:
 	 * @param _sName of the attribute.
 	 * @return attribute value, empty string if it doesn't exist.
 	 */ 
-	bool hasAttribute(string _sName);
+	bool hasAttribute(string _sName) const;
 
 	/** Get the value of an attribute.
 	 *
@@ -145,7 +142,7 @@ public:
 	 * @param _sDefaultValue value to return if the attribute isn't found
 	 * @return attribute value, _sDefaultValue if it doesn't exist.
 	 */ 
-	string getAttribute(string _sName, string _sDefaultValue = "");
+	string getAttribute(string _sName, string _sDefaultValue = "") const;
 
 	/** Get the value of a numerical attribute.
 	 *
@@ -153,8 +150,8 @@ public:
 	 * @param _fDefaultValue value to return if the attribute isn't found
 	 * @return attribute value, _fDefaultValue if it doesn't exist.
 	 */ 
-	float32 getAttributeNumerical(string _sName, float32 _fDefaultValue = 0);
-	double getAttributeNumericalDouble(string _sName, double _fDefaultValue = 0);
+	float32 getAttributeNumerical(string _sName, float32 _fDefaultValue = 0) const;
+	double getAttributeNumericalDouble(string _sName, double _fDefaultValue = 0) const;
 
 	/** Get the value of a boolean attribute.
 	 *
@@ -162,7 +159,7 @@ public:
 	 * @param _bDefaultValue value to return if the attribute isn't found
 	 * @return attribute value, _bDefaultValue if it doesn't exist.
 	 */ 
-	bool getAttributeBool(string _sName, bool _bDefaultValue = false);
+	bool getAttributeBool(string _sName, bool _bDefaultValue = false) const;
 
 
 
@@ -172,7 +169,7 @@ public:
 	 * @param _sKey option key
 	 * @return true if option does exist
 	 */ 
-	bool hasOption(string _sKey);
+	bool hasOption(string _sKey) const;
 
 	/** Get the value of an option within this XML Node
 	 *
@@ -180,7 +177,7 @@ public:
 	 * @param _sDefaultValue value to return if key isn't found
 	 * @return option value, _sDefaultValue if the option doesn't exist
 	 */ 
-	string getOption(string _sKey, string _sDefaultValue = "");
+	string getOption(string _sKey, string _sDefaultValue = "") const;
 
 	/** Get the value of an option within this XML Node
 	 *
@@ -188,7 +185,7 @@ public:
 	 * @param _fDefaultValue value to return if key isn't found
 	 * @return option value, _fDefaultValue if the option doesn't exist
 	 */ 
-	float32 getOptionNumerical(string _sKey, float32 _fDefaultValue = 0);
+	float32 getOptionNumerical(string _sKey, float32 _fDefaultValue = 0) const;
 
 	/** Get the value of an option within this XML Node
 	 *
@@ -196,14 +193,14 @@ public:
 	 * @param _bDefaultValue value to return if key isn't found
 	 * @return option value, _bDefaultValue if the option doesn't exist
 	 */ 
-	bool getOptionBool(string _sKey, bool _bDefaultValue = false);
+	bool getOptionBool(string _sKey, bool _bDefaultValue = false) const;
 
 	/** Get the value of an option within this XML Node
 	 *
 	 * @param _sKey option key
 	 * @return numerical array
 	 */ 
-	vector<float32> getOptionNumericalArray(string _sKey);
+	vector<float32> getOptionNumericalArray(string _sKey) const;
 
 
 
@@ -214,7 +211,7 @@ public:
 	 * @param _sNodeName the name of the new childnode
 	 * @return new child node
 	 */ 
-	XMLNode* addChildNode(string _sNodeName);
+	XMLNode addChildNode(string _sNodeName);
 
 	/** Create a new XML node as a child to this one, also add some content: 
 	 * &lt;...&gt;&lt;_sNodeName&gt;_sValue&lt;/_sNodeName>&lt;/...&gt;
@@ -223,7 +220,7 @@ public:
 	 * @param _sValue some node content
 	 * @return new child node
 	 */ 
-	XMLNode* addChildNode(string _sNodeName, string _sValue);
+	XMLNode addChildNode(string _sNodeName, string _sValue);
 
 	/** Create a new XML node as a child to this one, also add some numerical content: 
 	 * &lt;...&gt;&lt;_sNodeName&gt;_sValue&lt;/_sNodeName>&lt;/...&gt;
@@ -232,7 +229,7 @@ public:
 	 * @param _fValue some node content
 	 * @return new child node
 	 */ 
-	XMLNode* addChildNode(string _sNodeName, float32 _fValue);
+	XMLNode addChildNode(string _sNodeName, float32 _fValue);
 
 	/** Create a new XML node as a child to this one, also add a list of numerical content: 
 	 * &lt;...&gt;&lt;_sNodeName&gt;_sValue&lt;/_sNodeName>&lt;/...&gt;
@@ -242,7 +239,7 @@ public:
 	 * @param _iSize number of elements in _pfList
 	 * @return new child node
 	 */ 
-	XMLNode* addChildNode(string _sNodeName, float32* _pfList, int _iSize);
+	XMLNode addChildNode(string _sNodeName, float32* _pfList, int _iSize);
 
 	/** Add some text to the node: &lt;...&gt;_sText&lt;/...&gt;
 	 *
@@ -256,13 +253,38 @@ public:
 	 */ 
 	void setContent(float32 _fValue);
 
-	/** Add a list of numerical data to the node: &lt;...&gt;_sText&lt;/...&gt;
+	/** Add a list of numerical data to the node
 	 *
 	 * @param _pfList data
 	 * @param _iSize number of elements in the list
 	 */ 
 	void setContent(float32* _pfList, int _iSize);
 
+	/** Add a list of numerical data to the node
+	 *
+	 * @param _pfList data
+	 * @param _iSize number of elements in the list
+	 */
+	void setContent(double* _pfList, int _iSize);
+
+	/** Add a (2D) matrix of numerical data to the node
+	 *
+	 * @param _pfMatrix data
+	 * @param _iWidth width of the matrix
+	 * @param _iHeight height of the matrix
+	 * @param transposed true is C order, false is Fortran order
+	 */
+	void setContent(float32* _pfMatrix, int _iWidth, int _iHeight, bool transposed);
+
+	/** Add a (2D) matrix of numerical data to the node
+	 *
+	 * @param _pfMatrix data
+	 * @param _iWidth width of the matrix
+	 * @param _iHeight height of the matrix
+	 * @param transposed true is C order, false is Fortran order
+	 */
+	void setContent(double* _pfMatrix, int _iWidth, int _iHeight, bool transposed);
+
 	/** Add an attribute to this node: &lt;... _sName="_sValue"&gt;
 	 *
 	 * @param _sName name of the attribute
@@ -294,11 +316,11 @@ public:
 	
 	/** Print to String
 	 */
-	std::string toString();
+	std::string toString() const;
 
 	/** Print the node
 	 */ 
-	void print();
+	void print() const;
 
 protected:
 
diff --git a/include/astra/clog.h b/include/astra/clog.h
new file mode 100644
index 0000000..6ffb660
--- /dev/null
+++ b/include/astra/clog.h
@@ -0,0 +1,693 @@
+/* clog: Extremely simple logger for C.
+ *
+ * Features:
+ * - Implemented purely as a single header file.
+ * - Create multiple loggers.
+ * - Four log levels (debug, info, warn, error).
+ * - Custom formats.
+ * - Fast.
+ *
+ * Dependencies:
+ * - Should conform to C89, C++98 (but requires vsnprintf, unfortunately).
+ * - POSIX environment.
+ *
+ * USAGE:
+ *
+ * Include this header in any file that wishes to write to logger(s).  In
+ * exactly one file (per executable), define CLOG_MAIN first (e.g. in your
+ * main .c file).
+ *
+ *     #define CLOG_MAIN
+ *     #include "clog.h"
+ *
+ * This will define the actual objects that all the other units will use.
+ *
+ * Loggers are identified by integers (0 - 15).  It's expected that you'll
+ * create meaningful constants and then refer to the loggers as such.
+ *
+ * Example:
+ *
+ *  const int MY_LOGGER = 0;
+ *
+ *  int main() {
+ *      int r;
+ *      r = clog_init_path(MY_LOGGER, "my_log.txt");
+ *      if (r != 0) {
+ *          fprintf(stderr, "Logger initialization failed.\n");
+ *          return 1;
+ *      }
+ *      clog_info(CLOG(MY_LOGGER), "Hello, world!");
+ *      clog_free(MY_LOGGER);
+ *      return 0;
+ *  }
+ *
+ * The CLOG macro used in the call to clog_info is a helper that passes the
+ * __FILE__ and __LINE__ parameters for you, so you don't have to type them
+ * every time. (It could be prettier with variadic macros, but that requires
+ * C99 or C++11 to be standards compliant.)
+ *
+ * Errors encountered by clog will be printed to stderr.  You can suppress
+ * these by defining a macro called CLOG_SILENT before including clog.h.
+ *
+ * License: Do whatever you want. It would be nice if you contribute
+ * improvements as pull requests here:
+ *
+ *   https://github.com/mmueller/clog
+ *
+ * Copyright 2013 Mike Mueller <mike@subfocal.net>.
+ *
+ * As is; no warranty is provided; use at your own risk.
+ */
+
+#ifndef __CLOG_H__
+#define __CLOG_H__
+
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <errno.h>
+#include <fcntl.h>
+#include <stdarg.h>
+#include <stdlib.h>
+#include <stdio.h>
+#include <string.h>
+#include <time.h>
+#ifndef _MSC_VER
+#include <unistd.h>
+#else
+#define WIN32_LEAN_AND_MEAN
+#include <windows.h>
+#include <io.h>
+#define open _open
+#define close _close
+#define write _write
+#define snprintf _snprintf
+#endif
+
+/* Number of loggers that can be defined. */
+#define CLOG_MAX_LOGGERS 16
+
+/* Format strings cannot be longer than this. */
+#define CLOG_FORMAT_LENGTH 256
+
+/* Formatted times and dates should be less than this length. If they are not,
+ * they will not appear in the log. */
+#define CLOG_DATETIME_LENGTH 256
+
+/* Default format strings. */
+#define CLOG_DEFAULT_FORMAT "%d %t %f(%n): %l: %m\n"
+#define CLOG_DEFAULT_DATE_FORMAT "%Y-%m-%d"
+#define CLOG_DEFAULT_TIME_FORMAT "%H:%M:%S"
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+enum clog_level {
+    CLOG_DEBUG,
+    CLOG_INFO,
+    CLOG_WARN,
+    CLOG_ERROR
+};
+
+struct clog;
+
+/**
+ * Create a new logger writing to the given file path.  The file will always
+ * be opened in append mode.
+ *
+ * @param id
+ * A constant integer between 0 and 15 that uniquely identifies this logger.
+ *
+ * @param path
+ * Path to the file where log messages will be written.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_init_path(int id, const char *const path);
+
+/**
+ * Create a new logger writing to a file descriptor.
+ *
+ * @param id
+ * A constant integer between 0 and 15 that uniquely identifies this logger.
+ *
+ * @param fd
+ * The file descriptor where log messages will be written.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_init_fd(int id, int fd);
+
+/**
+ * Destroy (clean up) a logger.  You should do this at the end of execution,
+ * or when you are done using the logger.
+ *
+ * @param id
+ * The id of the logger to destroy.
+ */
+void clog_free(int id);
+
+#define CLOG(id) __FILE__, __LINE__, id
+
+/**
+ * Log functions (one per level).  Call these to write messages to the log
+ * file.  The first three arguments can be replaced with a call to the CLOG
+ * macro defined above, e.g.:
+ *
+ *     clog_debug(CLOG(MY_LOGGER_ID), "This is a log message.");
+ *
+ * @param sfile
+ * The name of the source file making this log call (e.g. __FILE__).
+ *
+ * @param sline
+ * The line number of the call in the source code (e.g. __LINE__).
+ *
+ * @param id
+ * The id of the logger to write to.
+ *
+ * @param fmt
+ * The format string for the message (printf formatting).
+ *
+ * @param ...
+ * Any additional format arguments.
+ */
+void clog_debug(const char *sfile, int sline, int id, const char *fmt, va_list ap);
+void clog_info(const char *sfile, int sline, int id, const char *fmt, va_list ap);
+void clog_warn(const char *sfile, int sline, int id, const char *fmt, va_list ap);
+void clog_error(const char *sfile, int sline, int id, const char *fmt, va_list ap);
+
+/**
+ * Set the minimum level of messages that should be written to the log.
+ * Messages below this level will not be written.  By default, loggers are
+ * created with level == CLOG_DEBUG.
+ *
+ * @param id
+ * The identifier of the logger.
+ *
+ * @param level
+ * The new minimum log level.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_level(int id, enum clog_level level);
+
+/**
+ * Set the format string used for times.  See strftime(3) for how this string
+ * should be defined.  The default format string is CLOG_DEFAULT_TIME_FORMAT.
+ *
+ * @param fmt
+ * The new format string, which must be less than CLOG_FORMAT_LENGTH bytes.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_time_fmt(int id, const char *fmt);
+
+/**
+ * Set the format string used for dates.  See strftime(3) for how this string
+ * should be defined.  The default format string is CLOG_DEFAULT_DATE_FORMAT.
+ *
+ * @param fmt
+ * The new format string, which must be less than CLOG_FORMAT_LENGTH bytes.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_date_fmt(int id, const char *fmt);
+
+/**
+ * Set the format string for log messages.  Here are the substitutions you may
+ * use:
+ *
+ *     %f: Source file name generating the log call.
+ *     %n: Source line number where the log call was made.
+ *     %m: The message text sent to the logger (after printf formatting).
+ *     %d: The current date, formatted using the logger's date format.
+ *     %t: The current time, formatted using the logger's time format.
+ *     %l: The log level (one of "DEBUG", "INFO", "WARN", or "ERROR").
+ *     %%: A literal percent sign.
+ *
+ * The default format string is CLOG_DEFAULT_FORMAT.
+ *
+ * @param fmt
+ * The new format string, which must be less than CLOG_FORMAT_LENGTH bytes.
+ * You probably will want to end this with a newline (\n).
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_fmt(int id, const char *fmt);
+
+/**
+ * Set the callback function.
+ *
+ * @param cb
+ * The new callback function.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_cb(int id, void (*cb)(const char *msg, size_t len));
+
+/**
+ * Set the file descriptor.
+ *
+ * @param id
+ * The identifier of the logger.
+ *
+ * @param fd
+ * The new file descriptor.
+ *
+ * @return
+ * Zero on success, non-zero on failure.
+ */
+int clog_set_fd(int id, int fd);
+
+
+/*
+ * No need to read below this point.
+ */
+
+/**
+ * The C logger structure.
+ */
+struct clog {
+
+    /* The current level of this logger. Messages below it will be dropped. */
+    enum clog_level level;
+
+    /* The file being written. */
+    int fd;
+
+    /* The format specifier. */
+    char fmt[CLOG_FORMAT_LENGTH];
+
+    /* Date format */
+    char date_fmt[CLOG_FORMAT_LENGTH];
+
+    /* Time format */
+    char time_fmt[CLOG_FORMAT_LENGTH];
+
+    /* Tracks whether the fd needs to be closed eventually. */
+    int opened;
+
+    /* Callback function for each log message. */
+    void (*cb)(const char *msg, size_t len);
+};
+
+void _clog_err(const char *fmt, ...);
+
+#ifdef CLOG_MAIN
+struct clog *_clog_loggers[CLOG_MAX_LOGGERS] = { 0 };
+#else
+extern struct clog *_clog_loggers[CLOG_MAX_LOGGERS];
+#endif
+
+#ifdef CLOG_MAIN
+
+const char *const CLOG_LEVEL_NAMES[] = {
+    "Debug",
+    "Info",
+    "Warning",
+    "Error",
+};
+
+int
+clog_init_path(int id, const char *const path)
+{
+    int fd = open(path, O_CREAT | O_WRONLY | O_APPEND, 0666);
+    if (fd == -1) {
+        _clog_err("Unable to open %s: %s\n", path, strerror(errno));
+        return 1;
+    }
+    if (clog_init_fd(id, fd)) {
+        close(fd);
+        return 1;
+    }
+    _clog_loggers[id]->opened = 1;
+    return 0;
+}
+
+int
+clog_init_fd(int id, int fd)
+{
+    struct clog *logger;
+
+    if (_clog_loggers[id] != NULL) {
+        _clog_err("Logger %d already initialized.\n", id);
+        return 1;
+    }
+
+    logger = (struct clog *) malloc(sizeof(struct clog));
+    if (logger == NULL) {
+        _clog_err("Failed to allocate logger: %s\n", strerror(errno));
+        return 1;
+    }
+
+    logger->level = CLOG_DEBUG;
+    logger->fd = fd;
+    logger->opened = 0;
+    strcpy(logger->fmt, CLOG_DEFAULT_FORMAT);
+    strcpy(logger->date_fmt, CLOG_DEFAULT_DATE_FORMAT);
+    strcpy(logger->time_fmt, CLOG_DEFAULT_TIME_FORMAT);
+    logger->cb = NULL;
+
+    _clog_loggers[id] = logger;
+    return 0;
+}
+
+void
+clog_free(int id)
+{
+    if (_clog_loggers[id]) {
+        if (_clog_loggers[id]->opened) {
+            close(_clog_loggers[id]->fd);
+        }
+        free(_clog_loggers[id]);
+        _clog_loggers[id]=NULL;
+    }
+}
+
+int
+clog_set_level(int id, enum clog_level level)
+{
+    if (_clog_loggers[id] == NULL) {
+        return 1;
+    }
+    if ((unsigned) level > CLOG_ERROR) {
+        return 1;
+    }
+    _clog_loggers[id]->level = level;
+    return 0;
+}
+
+int
+clog_set_fd(int id, int fd)
+{
+    if (_clog_loggers[id] == NULL) {
+        return 1;
+    }
+    _clog_loggers[id]->fd = fd;
+    return 0;
+}
+
+int
+clog_set_time_fmt(int id, const char *fmt)
+{
+    struct clog *logger = _clog_loggers[id];
+    if (logger == NULL) {
+        _clog_err("clog_set_time_fmt: No such logger: %d\n", id);
+        return 1;
+    }
+    if (strlen(fmt) >= CLOG_FORMAT_LENGTH) {
+        _clog_err("clog_set_time_fmt: Format specifier too long.\n");
+        return 1;
+    }
+    strcpy(logger->time_fmt, fmt);
+    return 0;
+}
+
+int
+clog_set_date_fmt(int id, const char *fmt)
+{
+    struct clog *logger = _clog_loggers[id];
+    if (logger == NULL) {
+        _clog_err("clog_set_date_fmt: No such logger: %d\n", id);
+        return 1;
+    }
+    if (strlen(fmt) >= CLOG_FORMAT_LENGTH) {
+        _clog_err("clog_set_date_fmt: Format specifier too long.\n");
+        return 1;
+    }
+    strcpy(logger->date_fmt, fmt);
+    return 0;
+}
+
+int
+clog_set_fmt(int id, const char *fmt)
+{
+    struct clog *logger = _clog_loggers[id];
+    if (logger == NULL) {
+        _clog_err("clog_set_fmt: No such logger: %d\n", id);
+        return 1;
+    }
+    if (strlen(fmt) >= CLOG_FORMAT_LENGTH) {
+        _clog_err("clog_set_fmt: Format specifier too long.\n");
+        return 1;
+    }
+    strcpy(logger->fmt, fmt);
+    return 0;
+}
+
+int
+clog_set_cb(int id, void (*cb)(const char *msg, size_t len))
+{
+    struct clog *logger = _clog_loggers[id];
+    if (logger == NULL) {
+        _clog_err("clog_set_cb: No such logger: %d\n", id);
+        return 1;
+    }
+    logger->cb = cb;
+    return 0;
+}
+
+/* Internal functions */
+
+size_t
+_clog_append_str(char **dst, char *orig_buf, const char *src, size_t cur_size)
+{
+    size_t new_size = cur_size;
+
+    while (strlen(*dst) + strlen(src) >= new_size) {
+        new_size *= 2;
+    }
+    if (new_size != cur_size) {
+        if (*dst == orig_buf) {
+            *dst = (char *) malloc(new_size);
+            strcpy(*dst, orig_buf);
+        } else {
+            *dst = (char *) realloc(*dst, new_size);
+        }
+    }
+
+    strcat(*dst, src);
+    return new_size;
+}
+
+size_t
+_clog_append_int(char **dst, char *orig_buf, long int d, size_t cur_size)
+{
+    char buf[40]; /* Enough for 128-bit decimal */
+    if (snprintf(buf, 40, "%ld", d) >= 40) {
+        return cur_size;
+    }
+    return _clog_append_str(dst, orig_buf, buf, cur_size);
+}
+
+size_t
+_clog_append_time(char **dst, char *orig_buf, struct tm *lt,
+                  const char *fmt, size_t cur_size)
+{
+    char buf[CLOG_DATETIME_LENGTH];
+    size_t result = strftime(buf, CLOG_DATETIME_LENGTH, fmt, lt);
+
+    if (result > 0) {
+        return _clog_append_str(dst, orig_buf, buf, cur_size);
+    }
+
+    return cur_size;
+}
+
+const char *
+_clog_basename(const char *path)
+{
+    const char *slash = strrchr(path, '/');
+    if (slash) {
+        path = slash + 1;
+    }
+#ifdef _WIN32
+    slash = strrchr(path, '\\');
+    if (slash) {
+        path = slash + 1;
+    }
+#endif
+    return path;
+}
+
+char *
+_clog_format(const struct clog *logger, char buf[], size_t buf_size,
+             const char *sfile, int sline, const char *level,
+             const char *message)
+{
+    size_t cur_size = buf_size;
+    char *result = buf;
+    enum { NORMAL, SUBST } state = NORMAL;
+    size_t fmtlen = strlen(logger->fmt);
+    size_t i;
+    time_t t = time(NULL);
+    struct tm *lt = localtime(&t);
+
+    sfile = _clog_basename(sfile);
+    result[0] = 0;
+    for (i = 0; i < fmtlen; ++i) {
+        if (state == NORMAL) {
+            if (logger->fmt[i] == '%') {
+                state = SUBST;
+            } else {
+                char str[2] = { 0 };
+                str[0] = logger->fmt[i];
+                cur_size = _clog_append_str(&result, buf, str, cur_size);
+            }
+        } else {
+            switch (logger->fmt[i]) {
+                case '%':
+                    cur_size = _clog_append_str(&result, buf, "%", cur_size);
+                    break;
+                case 't':
+                    cur_size = _clog_append_time(&result, buf, lt,
+                                                 logger->time_fmt, cur_size);
+                    break;
+                case 'd':
+                    cur_size = _clog_append_time(&result, buf, lt,
+                                                 logger->date_fmt, cur_size);
+                    break;
+                case 'l':
+                    cur_size = _clog_append_str(&result, buf, level, cur_size);
+                    break;
+                case 'n':
+                    cur_size = _clog_append_int(&result, buf, sline, cur_size);
+                    break;
+                case 'f':
+                    cur_size = _clog_append_str(&result, buf, sfile, cur_size);
+                    break;
+                case 'm':
+                    cur_size = _clog_append_str(&result, buf, message,
+                                                cur_size);
+                    break;
+            }
+            state = NORMAL;
+        }
+    }
+
+    return result;
+}
+
+void
+_clog_log(const char *sfile, int sline, enum clog_level level,
+          int id, const char *fmt, va_list ap)
+{
+    /* For speed: Use a stack buffer until message exceeds 4096, then switch
+     * to dynamically allocated.  This should greatly reduce the number of
+     * memory allocations (and subsequent fragmentation). */
+    char buf[4096];
+    size_t buf_size = 4096;
+    char *dynbuf = buf;
+    char *message;
+    int result;
+    struct clog *logger = _clog_loggers[id];
+
+    if (!logger) {
+        _clog_err("No such logger: %d\n", id);
+        return;
+    }
+
+    if (level < logger->level) {
+        return;
+    }
+
+    /* Format the message text with the argument list. */
+    result = vsnprintf(dynbuf, buf_size, fmt, ap);
+    if ((size_t) result >= buf_size) {
+        buf_size = result + 1;
+        dynbuf = (char *) malloc(buf_size);
+        result = vsnprintf(dynbuf, buf_size, fmt, ap);
+        if ((size_t) result >= buf_size) {
+            /* Formatting failed -- too large */
+            _clog_err("Formatting failed (1).\n");
+            free(dynbuf);
+            return;
+        }
+    }
+
+    /* Format according to log format and write to log */
+    {
+        char message_buf[4096];
+        message = _clog_format(logger, message_buf, 4096, sfile, sline,
+                               CLOG_LEVEL_NAMES[level], dynbuf);
+        if (!message) {
+            _clog_err("Formatting failed (2).\n");
+            if (dynbuf != buf) {
+                free(dynbuf);
+            }
+            return;
+        }
+        result = write(logger->fd, message, strlen(message));
+        if (logger->cb) logger->cb(message,strlen(message));
+        if (result == -1) {
+            _clog_err("Unable to write to log file: %s\n", strerror(errno));
+        }
+        if (message != message_buf) {
+            free(message);
+        }
+        if (dynbuf != buf) {
+            free(dynbuf);
+        }
+#ifndef _MSC_VER
+        fsync(logger->fd);
+#else
+        HANDLE h = (HANDLE) _get_osfhandle(logger->fd);
+        if (h != INVALID_HANDLE_VALUE) {
+            // This call will fail on a console fd, but that's ok.
+            FlushFileBuffers(h);
+        }
+#endif
+    }
+}
+
+void
+clog_debug(const char *sfile, int sline, int id, const char *fmt, va_list ap)
+{
+    _clog_log(sfile, sline, CLOG_DEBUG, id, fmt, ap);
+}
+
+void
+clog_info(const char *sfile, int sline, int id, const char *fmt, va_list ap)
+{
+    _clog_log(sfile, sline, CLOG_INFO, id, fmt, ap);
+}
+
+void
+clog_warn(const char *sfile, int sline, int id, const char *fmt, va_list ap)
+{
+    _clog_log(sfile, sline, CLOG_WARN, id, fmt, ap);
+}
+
+void
+clog_error(const char *sfile, int sline, int id, const char *fmt, va_list ap)
+{
+    _clog_log(sfile, sline, CLOG_ERROR, id, fmt, ap);
+}
+
+void
+_clog_err(const char *fmt, ...)
+{
+#ifdef CLOG_SILENT
+    (void) fmt;
+#else
+    va_list ap;
+
+    va_start(ap, fmt);
+    vfprintf(stderr, fmt, ap);
+#endif
+}
+
+#endif /* CLOG_MAIN */
+
+#ifdef __cplusplus
+} /* extern "C" */
+#endif
+
+#endif /* __CLOG_H__ */