diff --git a/stim/math/matrix.h b/stim/math/matrix.h
index 9ef0c48..0aca681 100644
--- a/stim/math/matrix.h
+++ b/stim/math/matrix.h
@@ -4,6 +4,7 @@
 //#include "rts/vector.h"
 #include <string.h>
 #include <iostream>
+#include <fstream>
 #include <stim/math/vector.h>
 #include <stim/math/vec3.h>
 //#include <stim/cuda/cudatools/callable.h>
@@ -31,6 +32,36 @@ namespace stim{
 			default:           return 0;
 		}
 	}
+
+	void save_mat4(char* data, std::string filename, std::string varname, size_t sx, size_t sy, mat4Format format){
+		//save the matrix file here (use the mat4 function above)
+		//data format: https://maxwell.ict.griffith.edu.au/spl/matlab-page/matfile_format.pdf (page 32)
+		
+		int MOPT = 0;									//initialize the MOPT type value to zero
+		int m = 0;										//little endian
+		int o = 0;										//reserved, always 0
+		int p = format;
+		int t = 0;
+		MOPT = m * 1000 + o * 100 + p * 10 + t;			//calculate the type value
+		int mrows = sx;
+		int ncols = sy;
+		int imagf = 0;									//assume real (for now)
+		varname.push_back('\0');									//add a null to the string
+		int namlen = varname.size();						//calculate the name size
+
+		size_t bytes = sx * sy * mat4Format_size(format);
+		std::ofstream outfile(filename, std::ios::binary);
+		outfile.write((char*)&MOPT, 4);
+		outfile.write((char*)&mrows, 4);
+		outfile.write((char*)&ncols, 4);
+		outfile.write((char*)&imagf, 4);
+		outfile.write((char*)&namlen, 4);
+		outfile.write((char*)&varname[0], namlen);
+		outfile.write((char*)data, bytes);				//write the matrix data
+		outfile.close();
+
+	}
+
 template <class T>
 class matrix {
 	//the matrix will be stored in column-major order (compatible with OpenGL)
@@ -380,8 +411,17 @@ public:
 	}
 
 	// saves the matrix as a Level-4 MATLAB file
-	void mat4(std::string filename, std::string name = std::string("unknown"), mat4Format = mat4_float) {
-		int MOPT = 0;									//initialize the MOPT type value to zero
+	void mat4(std::string filename, std::string name = std::string("unknown"), mat4Format format = mat4_float) {
+		if (format == mat4_float) {
+			if (sizeof(T) == 4) format = mat4_float32;
+			else if (sizeof(T) == 8) format = mat4_float64;
+			else {
+				std::cout << "stim::matrix ERROR - incorrect format specified" << std::endl;
+				exit(1);
+			}
+		}
+		stim::save_mat4((char*)M, filename, name, rows(), cols(), format);
+		/*int MOPT = 0;									//initialize the MOPT type value to zero
 		int m = 0;										//little endian
 		int o = 0;										//reserved, always 0
 		int p = mat4_float;
@@ -409,14 +449,10 @@ public:
 		outfile.write((char*)&namlen, 4);
 		outfile.write((char*)&name[0], namlen);
 		outfile.write((char*)M, bytes());				//write the matrix data
-		outfile.close();
+		outfile.close();*/
 	}
 };
 
-void save_mat4(char* data, std::string filename, std::string varname, size_t sx, size_t sy, mat4Format format){
-	//save the matrix file here (use the mat4 function above)
-	//data format: https://maxwell.ict.griffith.edu.au/spl/matlab-page/matfile_format.pdf (page 32)
-}
 }	//end namespace rts
 
 
--
libgit2 0.21.4