From: Benjamin Auder <benjamin.auder@somewhere>
Date: Wed, 11 Jan 2017 16:49:31 +0000 (+0100)
Subject: add test folder
X-Git-Url: https://git.auder.net/images/assets/js/img/pieces/cn.svg?a=commitdiff_plain;h=7b2720733e9aebe177c211119a9ec160c7e7117b;p=valse.git

add test folder
---

diff --git a/src/test/.test.EMGLLF.c.swp b/src/test/.test.EMGLLF.c.swp
new file mode 100644
index 0000000..176ac00
Binary files /dev/null and b/src/test/.test.EMGLLF.c.swp differ
diff --git a/src/test/.utils.c.swp b/src/test/.utils.c.swp
new file mode 100644
index 0000000..61f8f0a
Binary files /dev/null and b/src/test/.utils.c.swp differ
diff --git a/src/test/Makefile b/src/test/Makefile
new file mode 100644
index 0000000..c4931ad
--- /dev/null
+++ b/src/test/Makefile
@@ -0,0 +1,37 @@
+CC = gcc
+CFLAGS = -g -std=gnu99 -Wno-implicit-function-declaration
+LDFLAGS = -Lobj/ -lvalse_core
+LIB = valse_core.so
+LIB_SRC = $(wildcard ../sources/*.c)
+LIB_OBJ = $(LIB_SRC:.c=.o)
+
+all: $(LIB) test.EMGLLF test.EMGrank test.constructionModelesLassoMLE test.EMGrank test.constructionModelesLassoRank test.selectionTotale
+
+$(LIB): $(LIB_OBJ)
+	$(CC) -o $@ $^
+
+test.EMGLLF: test.EMGLLF.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.constructionModelesLassoMLE: test.constructionModelesLassoMLE.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.EMGrank: test.EMGrank.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.constructionModelesLassoRank: test.constructionModelesLassoRank.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.selectionTotale: test.selectionTotale.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+%.o: %.c
+	$(CC) -o $@ -c $< $(CFLAGS)
+
+clean:
+	rm -f *.o ../sources/*.o
+
+cclean: clean
+	rm -f $(LIB) $(TEST)
+
+.PHONY: all clean cclean
diff --git a/src/test/OLD_TEST_MATLAB/TODO b/src/test/OLD_TEST_MATLAB/TODO
new file mode 100644
index 0000000..d4b05cf
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/TODO
@@ -0,0 +1,4 @@
+1) transformer les .m en .R
+2) sauvegarder les résultats + entrée/sorties à partir du code R uniquement
+3) dans test.Truc.c, prendre en entrée les paramètres de data/ juste sauvegardés,
+   puis comparer les sorties aux sorties enregistrées
diff --git a/src/test/OLD_TEST_MATLAB/basicInitParameters.m b/src/test/OLD_TEST_MATLAB/basicInitParameters.m
new file mode 100644
index 0000000..50410f5
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/basicInitParameters.m
@@ -0,0 +1,19 @@
+function[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n,p,m,k)
+
+	phiInit = zeros(p,m,k);
+	
+	piInit = (1.0/k) * ones(1,k);
+	
+	rhoInit = zeros(m,m,k);
+	for r=1:k
+		rhoInit(:,:,r) = eye(m,m);
+	end
+	
+	gamInit = 0.1 * ones(n,k);
+	R = random('unid',k,n,1);
+	for i=1:n
+		gamInit(i,R(i)) = 0.9;
+	end
+	gamInit = gamInit / (sum(gamInit(1,:)));
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/checkOutput.m b/src/test/OLD_TEST_MATLAB/checkOutput.m
new file mode 100644
index 0000000..c56ed0b
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/checkOutput.m
@@ -0,0 +1,11 @@
+function[]=checkOutput(varName, matrix, refMatrix, tol)
+
+	fprintf('Checking %s\n',varName);
+	maxError = max(max(max(max(abs(matrix - refMatrix)))));
+	if maxError >= tol
+		fprintf('    Inaccuracy: max(abs(error)) = %g >= %g\n',maxError,tol);
+	else
+		fprintf('    OK\n');
+	end
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/covariance.m b/src/test/OLD_TEST_MATLAB/covariance.m
new file mode 100644
index 0000000..839c33c
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/covariance.m
@@ -0,0 +1,9 @@
+%covariance matrix for tests on synthetic data: A(i,j) = a^|i-j|
+function[A] = covariance(p,a)
+
+	A = a*ones(p,p);
+	for i=1:p
+		A(i,:) = A(i,:) .^ abs(i-(1:p));
+	end
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateIO.m b/src/test/OLD_TEST_MATLAB/generateIO.m
new file mode 100644
index 0000000..c677fd2
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateIO.m
@@ -0,0 +1,37 @@
+%X is generated following a gaussian mixture \sum pi_r N(meanX_k, covX_r)
+%Y is generated then, with Y_i ~ \sum pi_r N(Beta_r.X_i, covY_r)
+function[X,Y,Z] = generateIO(meanX, covX, covY, pi, beta, n)
+
+	[p, ~, k] = size(covX);
+	[m, ~, ~] = size(covY);
+	if exist('octave_config_info')
+		%Octave statistics package	doesn't have gmdistribution()
+		X = zeros(n, p);
+		Z = zeros(n);
+		cs = cumsum(pi);
+		for i=1:n
+			%TODO: vectorize ? http://stackoverflow.com/questions/2977497/weighted-random-numbers-in-matlab
+			tmpRand01 = rand();
+			[~,Z(i)] = min(cs - tmpRand01 >= 0);
+			X(i,:) = mvnrnd(meanX(Z(i),:), covX(:,:,Z(i)), 1);
+		end
+	else
+		gmDistX = gmdistribution(meanX, covX, pi);
+		[X, Z] = random(gmDistX, n);
+	end
+	
+	Y = zeros(n, m);
+	BX = zeros(n,m,k);
+	for i=1:n
+		for r=1:k
+			%compute beta_r . X_i
+			BXir = zeros(1, m);
+			for mm=1:m
+				BXir(mm) = dot(X(i,:), beta(:,mm,r));
+			end
+			%add pi(r) * N(beta_r . X_i, covY) to Y_i
+			Y(i,:) = Y(i,:) + pi(r) * mvnrnd(BXir, covY(:,:,r), 1);
+		end
+	end
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateIOdefault.m b/src/test/OLD_TEST_MATLAB/generateIOdefault.m
new file mode 100644
index 0000000..f4c3c1f
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateIOdefault.m
@@ -0,0 +1,23 @@
+%call generateIO with default parameters (random means, covariances = identity, equirepartition)
+function[X,Y,Z] = generateIOdefault(n, p, m, k)
+
+	rangeX = 100;
+	meanX = rangeX * (1 - 2*rand(k, p));
+	covX = zeros(p,p,k);
+	covY = zeros(m,m,k);
+	for r=1:k
+		covX(:,:,r) = eye(p);
+		covY(:,:,r) = eye(m);
+	end
+	pi = (1/k) * ones(1,k);
+	
+	%initialize beta to a random number of non-zero random value
+	beta = zeros(p,m,k);
+	for j=1:p
+		nonZeroCount = ceil(m*rand(1));
+		beta(j,1:nonZeroCount,:) = rand(nonZeroCount, k);
+	end
+	
+	[X,Y,Z] = generateIO(meanX, covX, covY, pi, beta, n);
+	
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m
new file mode 100644
index 0000000..bf0badf
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m
@@ -0,0 +1,47 @@
+function[] = generateRunSaveTest_EMGLLF(n, p, m, k, mini, maxi, gamma, lambda, varargin)
+
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 0.5};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, lambda] = optargs{:};
+	tau = 1e-6;
+	mini = int64(mini);
+	maxi = int64(maxi);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'lambda'), lambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k], delimiter);
+
+	[phi,rho,pi,LLF,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau);
+
+	%save output
+	dlmwrite(strcat(testFolder,'phi'), reshape(phi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rho'), reshape(rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'pi'), pi, delimiter);
+	dlmwrite(strcat(testFolder,'LLF'), LLF, delimiter);
+	dlmwrite(strcat(testFolder,'S'), reshape(S,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m
new file mode 100644
index 0000000..6e48d45
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m
@@ -0,0 +1,62 @@
+function[] = generateRunSaveTest_constructionModelesLassoMLE(n, p, m, k, mini, maxi, gamma, glambda, varargin)
+	
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 [0.0,0.01,0.02,0.03,0.05,0.1,0.2,0.3,0.5,0.7,0.85,0.99]};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, glambda] = optargs{:};
+	tau = 1e-6;
+	seuil = 1e-15;
+	mini = int64(mini);
+	maxi = int64(maxi);
+	L = length(glambda);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	A2 = zeros(p,m+1,L,'int64');
+    for i=1:L
+        A2(:,1,i) = 1:p;
+        A2(1:5,2,i) = 1:5;
+    end
+	A1 = zeros(p,m+1,L,'int64');
+    for i=1:L
+        A1(:,1,i) = 1:p;
+        A1(1:5,2,i) = 1:5;
+    end
+    
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'glambda'), glambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'seuil'), seuil, delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'A1'), reshape(A1,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'A2'), reshape(A2,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k,L], delimiter);
+
+	[phi,rho,pi,lvraisemblance] = constructionModelesLassoMLE(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2);
+
+	%save output
+	dlmwrite(strcat(testFolder,'phi'), reshape(phi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rho'), reshape(rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'pi'), reshape(pi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'lvraisemblance'), reshape(lvraisemblance,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m
new file mode 100644
index 0000000..6caee15
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m
@@ -0,0 +1,49 @@
+function[] = generateRunSaveTest_selectiontotale(n, p, m, k, mini, maxi, gamma, glambda, varargin)
+
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 [0.0,0.01,0.02,0.03,0.05,0.1,0.2,0.3,0.5,0.7,0.85,0.99]};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, glambda] = optargs{:};
+	tau = 1e-6;
+	seuil = 1e-15;
+	mini = int64(mini);
+	maxi = int64(maxi);
+	L = length(glambda);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+	
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'glambda'), glambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'seuil'), seuil, delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k,L], delimiter);
+
+	[A1,A2,Rho,Pi] = selectiontotale(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau);
+
+	%save output
+	dlmwrite(strcat(testFolder,'A1'), reshape(A1,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'A2'), reshape(A2,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Rho'), reshape(Rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Pi'), reshape(Pi,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m b/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m
new file mode 100644
index 0000000..27c1208
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m
@@ -0,0 +1,46 @@
+function[] = testConstructionModelesLassoMLE()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+	L = dimensions(5);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	glambda = dlmread(strcat(testFolder,'glambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	seuil = dlmread(strcat(testFolder,'seuil'), delimiter);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+	A1 = int64(reshape(dlmread(strcat(testFolder,'A1'), delimiter), p, m+1, L));
+	A2 = int64(reshape(dlmread(strcat(testFolder,'A2'), delimiter), p, m+1, L));
+
+	%run constructionModelesLassoMLE.m
+	[phi,rho,pi,lvraisemblance] = constructionModelesLassoMLE(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2);
+
+	%get all stored outputs
+	ref_phi = reshape(dlmread(strcat(testFolder,'phi'), delimiter), p, m, k, L);
+	ref_rho = reshape(dlmread(strcat(testFolder,'rho'), delimiter), m, m, k, L);
+	ref_pi = reshape(dlmread(strcat(testFolder,'pi'), delimiter), k, L);
+	ref_lvraisemblance = reshape(dlmread(strcat(testFolder,'lvraisemblance'), delimiter), L, 2);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('phi',phi,ref_phi,tol);
+	checkOutput('rho',rho,ref_rho,tol);
+	checkOutput('pi',pi,ref_pi,tol);
+	checkOutput('lvraisemblance',lvraisemblance,ref_lvraisemblance,tol);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testEMGLLF.m b/src/test/OLD_TEST_MATLAB/testEMGLLF.m
new file mode 100644
index 0000000..0dd9db8
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testEMGLLF.m
@@ -0,0 +1,44 @@
+function[] = testEMGLLF()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	lambda = dlmread(strcat(testFolder,'lambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+
+	%run EMGLLF.m
+	[phi,rho,pi,LLF,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau);
+
+	%get all stored outputs
+	ref_phi = reshape(dlmread(strcat(testFolder,'phi'), delimiter), p, m, k);
+	ref_rho = reshape(dlmread(strcat(testFolder,'rho'), delimiter), m, m, k);
+	ref_pi = dlmread(strcat(testFolder,'pi'), delimiter);
+	ref_LLF = dlmread(strcat(testFolder,'LLF'), delimiter);
+	ref_S = reshape(dlmread(strcat(testFolder,'S'), delimiter), p, m, k);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('phi',phi,ref_phi,tol);
+	checkOutput('rho',rho,ref_rho,tol);
+	checkOutput('pi',pi,ref_pi,tol);
+	checkOutput('LLF',LLF,ref_LLF,tol);
+	checkOutput('S',S,ref_S,tol);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testSelectiontotale.m b/src/test/OLD_TEST_MATLAB/testSelectiontotale.m
new file mode 100644
index 0000000..996017f
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testSelectiontotale.m
@@ -0,0 +1,44 @@
+function[] = testSelectiontotale()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+	L = dimensions(5);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	glambda = dlmread(strcat(testFolder,'glambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	seuil = dlmread(strcat(testFolder,'seuil'), delimiter);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+
+	%run constructionModelesLassoMLE.m
+	[A1,A2,Rho,Pi] = selectiontotale(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau);
+
+	%get all stored outputs
+	ref_A1 = int64(reshape(dlmread(strcat(testFolder,'A1'), delimiter), p, m+1, L));
+	ref_A2 = int64(reshape(dlmread(strcat(testFolder,'A2'), delimiter), p, m+1, L));
+	ref_Rho = reshape(dlmread(strcat(testFolder,'Rho'), delimiter), m, m, k, L);
+	ref_Pi = reshape(dlmread(strcat(testFolder,'Pi'), delimiter), k, L);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('A1',A1,ref_A1,tol);
+	checkOutput('A2',A2,ref_A2,tol);
+	checkOutput('Rho',Rho,ref_Rho,tol);
+	checkOutput('Pi',Pi,ref_Pi,tol);
+
+end
diff --git a/src/test/test.ConstructionModelesLassoMLE.c b/src/test/test.ConstructionModelesLassoMLE.c
new file mode 100644
index 0000000..a6c0cce
--- /dev/null
+++ b/src/test/test.ConstructionModelesLassoMLE.c
@@ -0,0 +1,143 @@
+#include "constructionModelesLassoMLE.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// phiInit
+	const mwSize dimPhiInit[] = {p, m, k};
+	Real* phiInit = readArray_real("phiInit",dimPhiInit,3);
+
+	// rhoInit
+	const mwSize dimRhoInit[] = {m, m, k};
+	Real* rhoInit = readArray_real("rhoInit",dimRhoInit,3);
+
+	// piInit
+	Real* piInit = readArray_real("piInit",&k,1);
+
+	// gamInit
+	const mwSize dimGamInit[] = {n, k};
+	Real* gamInit = readArray_real("gamInit",dimGamInit,2);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+	
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// gamma
+	Real* pgamma = readArray_real("gamma",&lengthOne,1);
+	Real gamma = *pgamma;
+	free(pgamma);
+
+	// lambda
+	Real* glambda = readArray_real("glambda",&L,1);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// seuil
+	Real* pseuil = readArray_real("seuil",&lengthOne,1);
+	Real seuil = *pseuil;
+	free(pseuil);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+	
+	// A1
+	const mwSize dimA[] = {p, m+1, L};
+	Int* A1 = readArray_int("A1",dimA,3);
+	
+	// A2
+	Int* A2 = readArray_int("A2",dimA,3);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	const mwSize dimPhi[] = {dimPhiInit[0], dimPhiInit[1], dimPhiInit[2], L};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]*sizeof(Real));
+
+	// rho
+	const mwSize dimRho[] = {dimRhoInit[0], dimRhoInit[1], dimRhoInit[2], L};
+	Real* rho = (Real*)malloc(dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]*sizeof(Real));
+
+	// pi
+	const mwSize dimPi[] = {k, L};
+	Real* pi = (Real*)malloc(dimPi[0]*dimPi[1]*sizeof(Real));
+
+	// lvraisemblance
+	const mwSize dimLvraisemblance[] = {L, 2};
+	Real* lvraisemblance = (Real*)malloc(dimLvraisemblance[0]*dimLvraisemblance[1]*sizeof(Real));
+
+	/////////////////////////////////////////
+	// Call to constructionModelesLassoMLE //
+	/////////////////////////////////////////
+
+	constructionModelesLassoMLE(
+		phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2,
+		phi,rho,pi,lvraisemblance,
+		n,p,m,k,L);
+	
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(X);
+	free(Y);
+	free(A1);
+	free(A2);
+	free(glambda);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi,4);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]);
+	free(phi);
+	free(ref_phi);
+	
+	// rho
+	Real* ref_rho = readArray_real("rho",dimRho,4);
+	compareArray_real("rho", rho, ref_rho, dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]);
+	free(rho);
+	free(ref_rho);
+	
+	// pi
+	Real* ref_pi = readArray_real("pi",dimPi,2);
+	compareArray_real("pi", pi, ref_pi, dimPi[0]*dimPi[1]);
+	free(pi);
+	free(ref_pi);
+	
+	// lvraisemblance
+	Real* ref_lvraisemblance = readArray_real("lvraisemblance",dimLvraisemblance,2);
+	compareArray_real("lvraisemblance", lvraisemblance, ref_lvraisemblance, dimLvraisemblance[0]*dimLvraisemblance[1]);
+	free(lvraisemblance);
+	free(ref_lvraisemblance);
+
+	return 0;
+}
diff --git a/src/test/test.EMGLLF.c b/src/test/test.EMGLLF.c
new file mode 100644
index 0000000..60516bc
--- /dev/null
+++ b/src/test/test.EMGLLF.c
@@ -0,0 +1,81 @@
+#include "EMGLLF.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	////////////
+	// INPUTS //
+	////////////
+
+	Int* dimensions = readArray_int("dimensions");
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	free(dimensions);
+
+	Real* phiInit = readArray_real("phiInit");
+	Real* rhoInit = readArray_real("rhoInit");
+	Real* piInit = readArray_real("piInit");
+	Real* gamInit = readArray_real("gamInit");
+	Int mini = read_int("mini");
+	Int maxi = read_int("maxi");
+	Real gamma = read_real("gamma");
+	Real lambda = read_real("lambda");
+	Real* X = readArray_real("X");
+	Real* Y = readArray_real("Y");
+	Real tau = read_real("tau");
+
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	Real* phi = (Real*)malloc(p*m*k*sizeof(Real));
+	Real* rho = (Real*)malloc(m*m*k*sizeof(Real));
+	Real* pi = (Real*)malloc(k*sizeof(Real));
+	Real* LLF = (Real*)malloc(maxi*sizeof(Real));
+	Real* S = (Real*)malloc(p*m*k*sizeof(Real));
+
+	////////////////////
+	// Call to EMGLLF //
+	////////////////////
+
+	EMGLLF_core(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau,
+		phi,rho,pi,LLF,S,
+		n,p,m,k);
+
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(X);
+	free(Y);
+
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi");
+	compareArray_real("phi", phi, ref_phi, p*m*k);
+	free(phi);
+	free(ref_phi);
+
+	Real* ref_rho = readArray_real("rho");
+	compareArray_real("rho", rho, ref_rho, m*m*k);
+	free(rho);
+	free(ref_rho);
+
+	Real* ref_pi = readArray_real("pi");
+	compareArray_real("pi", pi, ref_pi, k);
+	free(pi);
+	free(ref_pi);
+
+	Real* ref_LLF = readArray_real("LLF", maxi);
+	compareArray_real("LLF", LLF, ref_LLF);
+	free(LLF);
+	free(ref_LLF);
+
+	Real* ref_S = readArray_real("S");
+	compareArray_real("S", S, ref_S, p*m*k);
+	free(S);
+	free(ref_S);
+
+	return 0;
+}
diff --git a/src/test/test.EMGrank.c b/src/test/test.EMGrank.c
new file mode 100644
index 0000000..1227738
--- /dev/null
+++ b/src/test/test.EMGrank.c
@@ -0,0 +1,92 @@
+#include "EMGrank.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 4;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// Rho
+	const mwSize dimRho[] = {m, m, k};
+	Real* Rho = readArray_real("Rho",dimRho,3);
+
+	// Pi
+	Real* Pi = readArray_real("Pi",&k,1);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX, 2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY, 2);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+
+	// tau
+	Int* rank = readArray_int("rank",&k,1);
+
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	const mwSize dimPhi[] = {p, m, k};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*sizeof(Real));
+
+	// LLF
+	Real* LLF = (Real*)malloc(1*sizeof(Real));
+
+	//////////////////////////
+	// Main call to EMGrank //
+	//////////////////////////
+
+	EMGrank(Pi,Rho,mini,maxi,X,Y,tau,rank,
+		phi,LLF,
+		n,p,m,k);
+	
+	// free input pointers
+	free(Rho);
+	free(Pi);
+	free(X);
+	free(Y);
+	free(rank);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi, 3);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]);
+	free(phi);
+	free(ref_phi);
+	
+	// LLF
+	Real* ref_LLF = readArray_real("LLF",&lengthOne,1);
+	compareArray_real("LLF", LLF, ref_LLF, 1);
+	free(LLF);
+	free(ref_LLF);
+	
+	return 0;
+}
diff --git a/src/test/test.constructionModelesLassoRank.c b/src/test/test.constructionModelesLassoRank.c
new file mode 100644
index 0000000..2232e9e
--- /dev/null
+++ b/src/test/test.constructionModelesLassoRank.c
@@ -0,0 +1,107 @@
+#include "constructionModelesLassoRank.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// piInit
+	const mwSize dimPi[] = {k, L};
+	Real* Pi = readArray_real("Pi",dimPi,2);
+
+	// rhoInit
+	const mwSize dimRho[] = {m, m, k, L};
+	Real* Rho = readArray_real("Rho",dimRho,4);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+
+	// A1
+	const mwSize dimA[] = {p, L};
+	Int* A1 = readArray_int("A1",dimA,2);
+
+	// rangmin
+	Int* prangmin = readArray_int("rangmin",&lengthOne,1);
+	Int rangmin = *prangmin;
+	free(prangmin);
+
+	// rangmax
+	Int* prangmax = readArray_int("rangmax",&lengthOne,1);
+	Int rangmax = *prangmax;
+	free(prangmax);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	mwSize Size = (mwSize)pow(rangmax-rangmin+1, k);
+	const mwSize dimPhi[] = {p, m, k, L*Size};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]*sizeof(Real));
+
+	// lvraisemblance
+	const mwSize dimLvraisemblance[] = {L*Size, 2};
+	Real* lvraisemblance = (Real*)malloc(dimLvraisemblance[0]*dimLvraisemblance[1]*sizeof(Real));
+
+	//////////////////////////////////////////////
+	// Main call to constructionModelesLassoMLE //
+	//////////////////////////////////////////////
+
+	constructionModelesLassoRank(
+		Pi,Rho,mini,maxi,X,Y,tau,A1,rangmin,rangmax,
+		phi,lvraisemblance,
+		n,p,m,k,L);
+	
+	free(Rho);
+	free(Pi);
+	free(X);
+	free(Y);
+	free(A1);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi, 4);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]);
+	free(phi);
+	free(ref_phi);
+	
+	// lvraisemblance
+	Real* ref_lvraisemblance = readArray_real("lvraisemblance",dimLvraisemblance,2);
+	compareArray_real("lvraisemblance", lvraisemblance, ref_lvraisemblance, dimLvraisemblance[0]*dimLvraisemblance[1]);
+	free(lvraisemblance);
+	free(ref_lvraisemblance);
+	
+	return 0;
+}
diff --git a/src/test/test.selectiontotale.c b/src/test/test.selectiontotale.c
new file mode 100644
index 0000000..45cb461
--- /dev/null
+++ b/src/test/test.selectiontotale.c
@@ -0,0 +1,134 @@
+#include "selectiontotale.h"
+#include "ioutils.h"
+
+mwSize main(mwSize argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+	
+	////////////
+	// INPUTS //
+	////////////
+
+	// phiInit
+	const mwSize dimPhiInit[] = {p, m, k};
+	Real* phiInit = readArray_real("phiInit",dimPhiInit,3);
+
+	// rhoInit
+	const mwSize dimRhoInit[] = {m, m, k};
+	Real* rhoInit = readArray_real("rhoInit",dimRhoInit,3);
+
+	// piInit
+	Real* piInit = readArray_real("piInit",&k,1);
+
+	// gamInit
+	const mwSize dimGamInit[] = {n, k};
+	Real* gamInit = readArray_real("gamInit",dimGamInit, 2);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+	
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+	
+	// gamma
+	Real* pgamma = readArray_real("gamma",&lengthOne,1);
+	Real gamma = *pgamma;
+	free(pgamma);
+	
+	// lambda
+	Real* glambda = readArray_real("glambda",&L,1);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// seuil
+	Real* pseuil = readArray_real("seuil",&lengthOne,1);
+	Real seuil = *pseuil;
+	free(pseuil);
+	
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// A1
+	const mwSize dimA[] = {p, m+1, L};
+	Int* A1 = (Int*)malloc(dimA[0]*dimA[1]*dimA[2]*sizeof(Int));
+
+	// A2
+	Int* A2 = (Int*)malloc(dimA[0]*dimA[1]*dimA[2]*sizeof(Int));
+
+	// rho
+	const mwSize dimRho[] = {m, m, k, L};
+	Real* Rho = (Real*)malloc(dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]*sizeof(Real));
+
+	// pi
+	const mwSize dimPi[] = {k, L};
+	Real* Pi = (Real*)malloc(dimPi[0]*dimPi[1]*sizeof(Real));
+
+	//////////////////////////////////////////////
+	// Main call to constructionModelesLassoMLE //
+	//////////////////////////////////////////////
+
+	selectiontotale(
+		phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,
+		A1,A2,Rho,Pi,
+		n,p,m,k,L);
+	
+	// free input pointers
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(glambda);
+	free(X);
+	free(Y);
+	
+	// Compare to reference outputs
+	Int* ref_A1 = readArray_int("A1",dimA, 3);
+	compareArray_int("A1", A1, ref_A1, dimA[0]*dimA[1]*dimA[2]);
+	free(A1);
+	free(ref_A1);
+	
+	// A2
+	Int* ref_A2 = readArray_int("A2",dimA, 3);
+	compareArray_int("A2", A2, ref_A2, dimA[0]*dimA[1]*dimA[2]);
+	free(A2);
+	free(ref_A2);
+	
+	// Rho
+	Real* ref_Rho = readArray_real("Rho",dimRho, 4);
+	compareArray_real("Rho", Rho, ref_Rho, dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]);
+	free(Rho);
+	free(ref_Rho);
+	
+	// Pi
+	Real* ref_Pi = readArray_real("Pi",dimPi, 2);
+	compareArray_real("Pi", Pi, ref_Pi, dimPi[0]*dimPi[1]);
+	free(Pi);
+	free(ref_Pi);
+	
+	return 0;
+}
diff --git a/src/test/utils.c b/src/test/utils.c
new file mode 100644
index 0000000..7cdf740
--- /dev/null
+++ b/src/test/utils.c
@@ -0,0 +1,81 @@
+// Check if array == refArray
+void compareArray(const char* ID, const void* array, const void* refArray, int size,
+	int isInteger)
+{
+	Real EPS = 1e-5; //precision
+	printf("Checking %s\n",ID);
+	Real maxError = 0.0;
+	for (mwSize i=0; i<size; i++)
+	{
+		Real error = isInteger
+			? fabs(((Int*)array)[i] - ((Int*)refArray)[i])
+			: fabs(((Real*)array)[i] - ((Real*)refArray)[i]);
+		if (error >= maxError)
+			maxError = error;
+	}
+	if (maxError >= EPS)
+		printf("    Inaccuracy: max(abs(error)) = %g >= %g\n",maxError,EPS);
+	else
+		printf("    OK\n");
+}
+
+void compareArray_real(const char* ID, const void* array, const void* refArray, int size)
+{
+	return compareArray(ID, array, refArray, size, 0);
+}
+
+void compareArray_int(const char* ID, const void* array, const void* refArray, int size)
+{
+	return compareArray(ID, array, refArray, size, 1);
+}
+
+// Read array by columns (as in MATLAB) and return by-rows encoding
+void* readArray(const char* fileName, int isInteger)
+{
+	// need to prepend '../data/' (not really nice code...)
+	char* fullFileName = (char*)calloc(5+strlen(fileName)+1,sizeof(char));
+	strcat(fullFileName, "data/");
+	strcat(fullFileName, fileName);
+	FILE* file = fopen(fullFileName, "r");
+	free(fullFileName);
+
+	// first pass to know how many elements to allocate
+	// /////................... TODO
+
+	int d = 1;
+	for (int i=0; i<nbDims; i++)
+		totalDim *= dimensions[i];
+	size_t elementSize = isInteger
+		? sizeof(Int)
+		: sizeof(Real);
+	
+	// read all values, and convert them to by-rows matrices format
+	void* matlabArray = malloc(totalDim*elementSize);
+	char curChar = ' ';
+	char bufferNum[64];
+	for (mwSize u=0; u<totalDim; u++)
+	{
+		// position to next non-separator character
+		while (!feof(file) && (curChar==' ' || curChar=='\n' || curChar=='\t' || curChar==','))
+			curChar = fgetc(file);
+		// read number (as a string)
+		int bufferIndex = 0;
+		while (!feof(file) && curChar!=' ' && curChar!='\n' && curChar!='\t' && curChar!=',')
+		{
+			bufferNum[bufferIndex++] = curChar;
+			curChar = fgetc(file);
+		}
+		bufferNum[bufferIndex] = 0;
+		// transform string into Real, and store it at appropriate location
+		if (isInteger)
+			((Int*)matlabArray)[u] = atoi(bufferNum);
+		else
+			((Real*)matlabArray)[u] = atof(bufferNum);
+	}
+	fclose(file);
+	
+	void* brArray = matlabToBrArray(matlabArray, dimensions, nbDims, isInteger);
+	free(matlabArray);
+	return brArray;
+}
+