commit

MATLAB/calib_meth/EMNeNMF/emnenmf.m

@@ -1,4 +1,13 @@
-function [ T , RMSE ] = emnenmf( W , X , G , F , Omega_G, Omega_F, Phi_G, Phi_F , InnerMinIter , InnerMaxIter , Tmax , v, F_theo, delta_measure)
+function [ T , RMSE , MERs ] = emnenmf( X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, G, F, config)
+
+%% loading the config parameters 
+Tmax = config.Tmax;
+delta_measure = config.delta_measure;
+InnerMinIter = config.InnerMinIter;
+InnerMaxIter = config.InnerMaxIter;
+M_loop = config.M_loop;
+
+%% 
 
 X0 =X;
 Omega_G = (Omega_G == 1); % Logical mask is faster than indexing in matlab.
@@ -6,11 +15,11 @@ Omega_F = (Omega_F == 1); % Logical mask is faster than indexing in matlab.
 nOmega_G = ~Omega_G; % Logical mask is faster than indexing in matlab.
 nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
 
-[~, num_sensor] = size(F);
-num_sensor = num_sensor-1;
+num_sensor = config.numSensor;
 em_iter_max = round(Tmax / delta_measure) ;
 T = nan(1,em_iter_max);
-RMSE = nan(2,em_iter_max);
+RMSE = nan(1+config.numSubSensor,em_iter_max);
+MERs = nan(1+config.numSubSensor,em_iter_max);
 
 X = G*F+W.*(X0-G*F);
 GG = G'*G;
@@ -28,7 +37,7 @@ StoppingCritG = StoppingCritF;
 tic
 i = 1;
 T(i) = toc;
-RMSE(:,i) = vecnorm(F(:,1:end-1)- F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
+RMSE(:,i) = vecnorm(F(:,1:num_sensor)- F_theo(:,1:num_sensor),2,2)/sqrt(num_sensor);
 niter = 0;
 T_E = [];
 T_M = [];
@@ -38,7 +47,7 @@ while toc<Tmax
     X = G*F+W.*(X0-G*F);
     T_E = cat(1,T_E,toc - t_e);
     
-    for j =1:v
+    for j =1:M_loop
         
         t_m = toc;
         FF = F*F';
@@ -68,6 +77,8 @@ while toc<Tmax
             if i > em_iter_max
                 break
             end
+            [MER,~]=bss_eval_mix(F_theo',F');
+            MERs(:,i) = MER;
             T(i) = toc;
             RMSE(:,i) = vecnorm(F(:,1:end-1) - F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
         end
@@ -77,8 +88,7 @@ while toc<Tmax
   
 end
 niter  
-disp(['em E step : ',num2str(median(T_E))])
-disp(['em M step : ',num2str(median(T_M))])
-end
-
-
+disp(['em E step : ',num2str(mean(T_E))])
+disp(['em M step : ',num2str(mean(T_M))])
+disp(['RMSE 2ème colonne de G : ',num2str(norm(G(:,2) - X_theo(:,end),2)/sqrt(config.sceneWidth*config.sceneLength))])
+end

MATLAB/calib_meth/INCAL/IN_Cal.m

@@ -1,4 +1,10 @@
-function [ T , RMSE ] = IN_Cal( W , X , G , F , Omega_G , Omega_F , Phi_G , Phi_F , F_theo , Tmax, delta_measure )
+function [ T , RMSE , MERs ] = IN_Cal( X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, G, F, config )
+
+%% loading the config parameters 
+Tmax = config.Tmax;
+delta_measure = config.delta_measure;
+
+%% 
 
 [~, num_sensor] = size(F);
 num_sensor = num_sensor-1;
@@ -11,7 +17,8 @@ Delta_F = F.*Omega_Bar_F;
 
 iter_max = round(Tmax / delta_measure) ;
 T = nan(1,iter_max);
-RMSE = nan(2,iter_max);
+RMSE = nan(1+config.numSubSensor,iter_max);
+MERs = nan(1+config.numSubSensor,iter_max);
 
 tic
 i = 1;
@@ -33,6 +40,8 @@ while toc<Tmax
         if i > iter_max
             break
         end
+        [MER,~]=bss_eval_mix(F_theo',F');
+        MERs(:,i) = MER;
         T(i) = toc;
         RMSE(:,i) = vecnorm(F(:,1:end-1) - F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
     end

MATLAB/calib_meth/REMNeNMF/remnenmf.m

@@ -1,4 +1,14 @@
-function [ T , RMSE ] = remnenmf( W , X , G , F , Omega_G, Omega_F, Phi_G, Phi_F , InnerMinIter , InnerMaxIter , Tmax , v, F_theo, delta_measure, nu)
+function [ T , RMSE ] = remnenmf( X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, G, F, config)
+
+%% loading the config parameters 
+Tmax = config.Tmax;
+delta_measure = config.delta_measure;
+InnerMinIter = config.InnerMinIter;
+InnerMaxIter = config.InnerMaxIter;
+M_loop = config.M_loop;
+nu = config.nu;
+
+%% 
 
 r = nu;
 X0 = X;
@@ -11,7 +21,7 @@ nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
 num_sensor = num_sensor-1;
 em_iter_max = round(Tmax / delta_measure) ;
 T = nan(1,em_iter_max);
-RMSE = nan(2,em_iter_max);
+RMSE = nan(1+config.numSubSensor,em_iter_max);
 
 nW = (1-W);
 % X = G*F+W.*(X0-G*F);
@@ -54,7 +64,7 @@ while toc<Tmax
     T_E = cat(1,T_E,toc - t_e);
     
     % Maximization step
-    for j =1:v
+    for j =1:M_loop
         
         t_m = toc;
 %         F_R = F * R;
@@ -100,8 +110,8 @@ while toc<Tmax
     
 end
 niter
-disp(['rem E step : ',num2str(median(T_E))])
-disp(['rem M step : ',num2str(median(T_M))])
+disp(['rem E step : ',num2str(mean(T_E))])
+disp(['rem M step : ',num2str(mean(T_M))])
 end
 
 function [ L,R ] = RSI_compression(X,r,varargin)

MATLAB/calib_meth/test/remnenmf_not_full.m

@@ -1,4 +1,14 @@
-function [ T , RMSE ] = remnenmf_not_full( W , X , G , F , Omega_G, Omega_F, Phi_G, Phi_F , InnerMinIter , InnerMaxIter , Tmax , v, F_theo, delta_measure, nu)
+function [ T , RMSE , MERs ] = remnenmf_not_full( X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, G, F, config)
+
+%% loading the config parameters 
+Tmax = config.Tmax;
+delta_measure = config.delta_measure;
+InnerMinIter = config.InnerMinIter;
+InnerMaxIter = config.InnerMaxIter;
+M_loop = config.M_loop;
+nu = config.nu;
+
+%% 
 
 r = nu;
 X0 = X;
@@ -11,7 +21,8 @@ nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
 num_sensor = num_sensor-1;
 em_iter_max = round(Tmax / delta_measure) ;
 T = nan(1,em_iter_max);
-RMSE = nan(2,em_iter_max);
+RMSE = nan(1+config.numSubSensor,em_iter_max);
+MERs = nan(1+config.numSubSensor,em_iter_max);
 
 nW = (1-W);
 % X = G*F+W.*(X0-G*F);
@@ -54,7 +65,7 @@ while toc<Tmax
     T_E = cat(1,T_E,toc - t_e);
     
     % Maximization step
-    for j =1:v
+    for j =1:M_loop
         
         t_m = toc;
 %         F_R = F * R;
@@ -91,6 +102,8 @@ while toc<Tmax
             if i > em_iter_max
                 break
             end
+            [MER,~]=bss_eval_mix(F_theo',F');
+            MERs(:,i) = MER;
             T(i) = toc;
             RMSE(:,i) = vecnorm(F(:,1:end-1) - F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
         end
@@ -100,8 +113,8 @@ while toc<Tmax
     
 end
 niter
-disp(['rem E step : ',num2str(median(T_E))])
-disp(['rem M step : ',num2str(median(T_M))])
+disp(['rem E step : ',num2str(mean(T_E))])
+disp(['rem M step : ',num2str(mean(T_M))])
 end
 
 function [ L,R ] = RSI_compression(X,r,varargin)

MATLAB/common/bss_eval_mix.m

@@ -0,0 +1,63 @@
+function [MER,perm]=bss_eval_mix(Ae,A)
+
+% BSS_EVAL_MIX Ordering and measurement of the quality of an estimated
+% (possibly frequency-dependent) mixing matrix
+%
+% [MER,perm]=bss_eval_mix(Ae,A)
+%
+% Inputs:
+% Ae: either a nchan x nsrc estimated mixing matrix (for instantaneous
+% mixtures) or a nchan x nsrc x nbin estimated frequency-dependent mixing
+% matrix (for convolutive mixtures)
+% A: the true nchan x nsrc or nchan x nsrc x nbin mixing matrix
+%
+% Outputs:
+% MER: nsrc x 1 vector of Mixing Error Ratios (SNR-like criterion averaged
+% over frequency and expressed in decibels, allowing arbitrary scaling for
+% each source in each frequency bin)
+% perm: nsrc x 1 vector containing the best ordering of estimated sources
+% in the maximum MER sense (estimated source number perm(j) corresponds to
+% true source number j)
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Copyright 2008 Emmanuel Vincent
+% This software is distributed under the terms of the GNU Public License
+% version 3 (http://www.gnu.org/licenses/gpl.txt)
+% If you find it useful, please cite the following reference:
+% Emmanuel Vincent, Shoko Araki and Pau Bofill, "The 2008 Signal Separation
+% Evaluation Campaign: A community-based approach to large-scale
+% evaluation," In Proc. Int. Conf. on Independent Component Analysis and
+% Signal Separation (ICA), pp. 734-741, 2009.
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+%%% Errors %%%
+if nargin<2, error('Not enough input arguments.'); end
+[nchan,nsrc,nbin]=size(Ae);
+[nchan2,nsrc2,nbin2]=size(A);
+if ~((nchan2==nchan)&&(nsrc2==nsrc)&&(nbin2==nbin)), error('The estimated and true mixing matrix must have the same size.'); end
+
+%%% Performance criterion %%%
+% Computation of the criterion for all possible pair matches
+MER=zeros(nsrc,nsrc,nbin);
+for f=1:nbin,
+    for jest=1:nsrc,
+        for jtrue=1:nsrc,
+            Aproj=A(:,jtrue,f)'*Ae(:,jest,f)/sum(abs(A(:,jtrue,f)).^2)*A(:,jtrue,f);
+            MER(jest,jtrue,f)=10*log10(sum(abs(Aproj).^2)/sum(abs(Ae(:,jest,f)-Aproj).^2));
+        end
+    end
+end
+MER=mean(MER,3);
+% Selection of the best ordering
+perm=perms(1:nsrc);
+nperm=size(perm,1);
+meanMER=zeros(nperm,1);
+for p=1:nperm,
+    meanMER(p)=mean(MER((0:nsrc-1)*nsrc+perm(p,:)));
+end
+[meanMER,popt]=max(meanMER);
+perm=perm(popt,:).';
+MER=MER((0:nsrc-1).'*nsrc+perm);
+
+return;

MATLAB/common/data_gen.m

@@ -1,34 +1,88 @@
-function [X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit] = data_gen(s_width, s_length, run, N_Ref, N_Cpt, Mu_beta, Mu_alpha, Bound_beta, Bound_alpha, MV, RV, var_n) 
+function [X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit] = data_gen(config, run) 
 
 rng(run+1306) % Random seed
 
+s_width         = config.sceneWidth;
+s_length        = config.sceneLength;
+N_Ref           = config.numRef;
+N_Cpt           = config.numSensor;
+N_sousCpt       = config.numSubSensor;
+Bound_phen      = config.Bound_phen;
+Mu_beta         = config.Mu_beta;
+Mu_alpha        = config.Mu_alpha;
+Bound_beta      = config.Bound_beta;
+Bound_alpha     = config.Bound_alpha;
+gamma           = config.gamma;
+MV              = config.mvR;
+RV              = config.rdvR;
+var_n           = config.var_n;
+
 %% Scene simulation 
 
-phenLowerBound = 0.05;
-phenUpperBound = 0.15;
 n_pic = 15;
 s_n = s_width*s_length; % Total number of areas in the scene
 [xx,yy] = meshgrid((-1:2/(s_width-1):1),(-1:2/(s_length-1):1));
 xxyy = cat(2,xx(:),yy(:));
-
-g = zeros(1,s_n);
-for pic = 1:n_pic
-    mu = 2*(rand(1,2)-0.5);
-    sig = diag([ phenLowerBound + (phenUpperBound - phenLowerBound)*abs(randn()+0.5) , phenLowerBound + (phenUpperBound - phenLowerBound)*abs(randn()+0.5) ]);
-    g = g + mvnpdf(xxyy,mu,sig);
+G_theo = ones(s_n,1);
+for sensor = 1:N_sousCpt
+    g = zeros(s_n,1);
+    for pic = 1:n_pic
+        mu = 2*(rand(1,2)-0.5);
+        sig = diag([ Bound_phen((sensor-1)*2+1) + (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))*abs(randn()+0.5) , ...
+            Bound_phen((sensor-1)*2+1) + (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))*abs(randn()+0.5) ]);
+        g = g + mvnpdf(xxyy,mu,sig);
+    end
+    g = (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))/(max(g)-min(g))*(g-min(g)) + Bound_phen((sensor-1)*2+1);
+    % .5*(g/max(g))+1e-5
+    G_theo = cat(2, G_theo, g);
 end
 
-g = g-min(g);
-g = .5*(g/max(g))+1e-5;
-G_theo = [ones(s_n,1),g]; % Theoretical matrix G (see eq.(3) of [1])
-
 
 %% Sensors simulation
 
-F_theo = [max(Bound_beta(1),min(Bound_beta(2),Mu_beta+.5*randn(1,N_Cpt)));...
-     max(Bound_alpha(1),min(Bound_alpha(2),Mu_alpha+.5*randn(1,N_Cpt)))];
-F_theo = [F_theo,[0;1]]; % Theoretical matrix F (see eq.(3) of [1])
+F_theo = zeros(N_sousCpt+1, N_sousCpt*(N_Cpt+1));
+
+for sensor = 1:N_sousCpt
+    F_theo(1,(sensor-1)*(N_Cpt+1)+1:sensor*(N_Cpt+1)) = ...
+        cat(2, max(Bound_beta((sensor-1)*2+1), min(Bound_beta((sensor-1)*2+2), ...
+        Mu_beta(sensor)+(Bound_beta((sensor-1)*2+2)-Bound_beta((sensor-1)*2+1))*0.55*randn(1,N_Cpt))), 0);
+end
+
+for sen = 1:N_sousCpt
+    f_theo = cat(2, max(Bound_alpha((sen-1)*2+1),...
+        min( Bound_alpha((sen-1)*2+2), ...
+        Mu_alpha(sen)+(Bound_alpha((sen-1)*2+2)-Bound_alpha((sen-1)*2+1))*0.25*randn(1,N_Cpt))), 1);
+    F_theo(sen+1, (sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)) = f_theo;
+    
+    C = (1-gamma)/gamma*(F_theo(1,(sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)-1)+Bound_phen((sen-1)*2+1)*f_theo(1:end-1));
+    list_nosen = 1:N_sousCpt;
+    list_nosen(sen) = [];
+    maxPhen_nosen = norm(Bound_phen(2*list_nosen));
+    for sor = list_nosen
+        f_theo_nosen = rand(1,N_Cpt).*C/(sqrt(N_sousCpt)*maxPhen_nosen);
+        other_f_theo = cat(2, f_theo_nosen, 0);
+        F_theo(sor+1, (sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)) = other_f_theo;
+    end
+end
 
+% F_theo = [];
+% 
+% for sensor = 1:N_sousCpt
+%     F_theo = cat(2, F_theo, cat(2, max(Bound_beta(1), min(Bound_beta(2), Mu_beta+.5*randn(1,N_Cpt))), 0));
+% end
+% 
+% 
+% for sen = 1:N_sousCpt
+%     f_theo = [];
+%     for sor = 1:N_sousCpt
+%         if sen==sor
+%             f_theo = cat(2, f_theo, cat(2, max(Bound_alpha(1), min(Bound_alpha(2), Mu_alpha+.5*randn(1,N_Cpt))), sen==sor));
+%         else
+%             f_theo = cat(2, f_theo, cat(2, 0*max(Bound_alpha(1), min(Bound_alpha(2), Mu_alpha+.5*randn(1,N_Cpt))), sen==sor));
+%         end
+%     end
+%     F_theo = cat(1, F_theo, f_theo);
+% end
 
 %% Data simulation
 
@@ -60,8 +114,10 @@ for i = 1 : N_data
     W(xx(idx_data(i)),yy(idx_data(i))) = 1; % Sensor measurement placement
 end
 
-N = var_n*randn(s_n,N_Cpt+1); % Noise simulation
-N(:,end) = 0;
+W = repmat(W, 1, N_sousCpt);
+
+N = var_n*randn(s_n,N_sousCpt*(N_Cpt+1)); % Noise simulation
+N(:,(N_Cpt+1)*(1:N_sousCpt)) = 0;
 N = max(N,-X_theo);
 
 X = W.*(X_theo+N); % Data matrix X
@@ -69,15 +125,54 @@ X = W.*(X_theo+N); % Data matrix X
 %% Calibration parameters
 
 % % Common parameters
-Omega_G = [ones(s_n,1),W(:,end)]; % Mask on known values in G (see eq.(14) of [1])
-Omega_F = [zeros(2,N_Cpt),[1;1]]; % Mask on known values in F (see eq.(15) of [1])
-Phi_G = [ones(s_n,1),X(:,end)]; % Known values in G (see eq.(14) of [1])
-Phi_F = [zeros(2,N_Cpt),[0;1]]; % Known values in F (see eq.(15) of [1])
-Ginit = abs(randn(s_n,2)+mean(Phi_G(idx_Ref,end))); % Initial matrix G : randn + mean of known ref values
+Omega_G = [ones(s_n,1),W(:,(N_Cpt+1)*(1:N_sousCpt))]; % Mask on known values in G (see eq.(14) of [1])
+Omega_F = zeros(N_sousCpt+1,N_sousCpt*(N_Cpt+1)); 
+Omega_F(:,(N_Cpt+1)*(1:N_sousCpt)) = 1; % Mask on known values in F (see eq.(15) of [1])
+Phi_G = [ones(s_n,1),X(:,(N_Cpt+1)*(1:N_sousCpt))]; % Known values in G (see eq.(14) of [1])
+Phi_F = F_theo .* Omega_F; % Known values in F (see eq.(15) of [1])
+Ginit = ones(s_n,1);
+for sensor = 1:N_sousCpt
+    g = zeros(s_n,1);
+    for pic = 1:n_pic
+        mu = 2*(rand(1,2)-0.5);
+        sig = diag([ Bound_phen((sensor-1)*2+1) + (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))*abs(randn()+0.5) , ...
+            Bound_phen((sensor-1)*2+1) + (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))*abs(randn()+0.5) ]);
+        g = g + mvnpdf(xxyy,mu,sig);
+    end
+    g = (Bound_phen((sensor-1)*2+2) - Bound_phen((sensor-1)*2+1))/(max(g)-min(g))*(g-min(g)) + Bound_phen((sensor-1)*2+1);
+    % .5*(g/max(g))+1e-5
+    Ginit = cat(2, Ginit, g);
+end
 Ginit = (1-Omega_G).*Ginit+Phi_G;
-Finit = [max(Bound_beta(1),min(Bound_beta(2),Mu_beta+.5*randn(1,N_Cpt)));...
-     max(Bound_alpha(1),min(Bound_alpha(2),Mu_alpha+.5*randn(1,N_Cpt)))];
-Finit = [Finit,[0;1]]; % Initial matrix F
-Finit = (1-Omega_F).*Finit+Phi_F;
+% Ginit = G_theo;
+
+Finit = zeros(N_sousCpt+1, N_sousCpt*(N_Cpt+1));
+
+for sensor = 1:N_sousCpt
+    Finit(1,(sensor-1)*(N_Cpt+1)+1:sensor*(N_Cpt+1)) = ...
+        cat(2, max(Bound_beta((sensor-1)*2+1), min(Bound_beta((sensor-1)*2+2), ...
+        Mu_beta(sensor)+(Bound_beta((sensor-1)*2+2)-Bound_beta((sensor-1)*2+1))*0.25*randn(1,N_Cpt))), 0);
+end
+
+for sen = 1:N_sousCpt
+    finit = cat(2, max(Bound_alpha((sen-1)*2+1),...
+        min( Bound_alpha((sen-1)*2+2), ...
+        Mu_alpha(sen)+(Bound_alpha((sen-1)*2+2)-Bound_alpha((sen-1)*2+1))*0.25*randn(1,N_Cpt))), 1);
+    Finit(sen+1, (sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)) = finit;
+    
+%     C = (1-gamma)/gamma*(Finit(1,(sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)-1)+Bound_phen((sen-1)*2+1)*finit(1:end-1));
+%     list_nosen = 1:N_sousCpt;
+%     list_nosen(sen) = [];
+%     maxPhen_nosen = norm(Bound_phen(2*list_nosen));
+    for sor = 1:N_sousCpt
+        if sen~=sor
+%             finit_nosen = rand(1,N_Cpt).*C/(sqrt(N_sousCpt)*maxPhen_nosen);
+%             other_finit = cat(2, finit_nosen, 0);
+%             Finit(sor+1, (sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)) = other_finit;
+            Finit(sor+1, (sen-1)*(N_Cpt+1)+1:sen*(N_Cpt+1)) = zeros(1,N_Cpt+1);
+        end
+    end
+end
+% Finit = F_theo;
 
 end

MATLAB/config.json

@@ -0,0 +1,27 @@
+{
+    "sceneWidth" : 10,
+    "sceneLength" : 10,
+    "numSensor" : 100,
+    "numSubSensor" : 1,
+    "CalibratedSensor" : [1],
+    "numRef" : 4,
+    "rdvR" : 0.3,
+    "mvR" : 0.5,
+    "calibrationMethods" : ["IN_Cal"],
+    "numRuns" : 1,
+    "Bound_phen" : [0.05,0.35 , 10,20],
+    "Mu_beta" : [2 , 5],
+    "Mu_alpha" : [30 , 0.5],
+    "Bound_beta" : [1,3 , 4.5,5.5],
+    "Bound_alpha" : [20,40 , 0.4,0.6],
+    "gamma" : 0.8,
+    "var_n" : 0,
+    "Tmax" : 30,
+    "nu" : 12,
+    "save2dat" : 1,
+    "display" : 1,
+    "delta_measure" : 1,
+    "M_loop" : 50,
+    "InnerMinIter" : 5,
+    "InnerMaxIter" : 20
+}

MATLAB/config2.json

@@ -0,0 +1,27 @@
+{
+    "sceneWidth" : 10,
+    "sceneLength" : 10,
+    "numSensor" : 100,
+    "numSubSensor" : 2,
+    "CalibratedSensor" : [1],
+    "numRef" : 4,
+    "rdvR" : 0.3,
+    "mvR" : 0.5,
+    "calibrationMethods" : ["emnenmf"],
+    "numRuns" : 1,
+    "Bound_phen" : [0.05,0.35 , 10,20],
+    "Mu_beta" : [2 , 5],
+    "Mu_alpha" : [30 , 0.5],
+    "Bound_beta" : [1,3 , 4.5,5.5],
+    "Bound_alpha" : [20,40 , 0.4,0.6],
+    "gamma" : 0.8,
+    "var_n" : 0,
+    "Tmax" : 30,
+    "nu" : 12,
+    "save2dat" : 1,
+    "display" : 1,
+    "delta_measure" : 1,
+    "M_loop" : 50,
+    "InnerMinIter" : 5,
+    "InnerMaxIter" : 20
+}

MATLAB/main.m

@@ -0,0 +1,8 @@
+clear all
+
+%% Adding every files in the path
+
+addpath(genpath(pwd))
+
+run_config("config2.json")
+% run_config_forget_other_sen("config2.json")

MATLAB/run_config.m

@@ -0,0 +1,73 @@
+function [res] = run_config(varargin)
+    if nargin == 0
+        file_name = "default_config.json";
+    elseif nargin == 1
+        file_name = varargin{1};
+    else
+        file_name = varargin{1};
+        warning("run_config only accepts one argument. Every other argument is ignored")
+    end
+    config = jsondecode(fileread(file_name));
+    methods = config.calibrationMethods;
+    iter_max = round(config.Tmax / config.delta_measure);
+    RMSEs = cell(length(methods),1+config.numSubSensor);
+    times = cell(length(methods),1);
+    
+    for run = 1:config.numRuns
+        % data generation
+        [X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit] = data_gen(config, run);
+        seed = floor(rand*10000);
+        
+        % We want to only calibrate the sensor of interest
+        X = X(:,1:config.numSensor+1);
+        X_theo = X_theo(:,1:config.numSensor+1);
+        W = W(:,1:config.numSensor+1);
+        F_theo = F_theo(:,1:config.numSensor+1);
+        Omega_F = Omega_F(:,1:config.numSensor+1);
+        Phi_F = Phi_F(:,1:config.numSensor+1);
+        Finit = Finit(:,1:config.numSensor+1);
+
+%         X = X(:,1:config.numSensor+1);
+%         X_theo = X_theo(:,1:config.numSensor+1);
+%         W = W(:,1:config.numSensor+1);
+%         F_theo = F_theo(1:2,1:config.numSensor+1);
+%         Omega_F = Omega_F(1:2,1:config.numSensor+1);
+%         Phi_F = Phi_F(1:2,1:config.numSensor+1);
+%         Finit = Finit(1:2,1:config.numSensor+1);
+%         Omega_G = Omega_G(:,1:2);
+%         Phi_G = Phi_G(:,1:2);
+%         Ginit = Ginit(:,1:2);
+%         config.numSubSensor = 1;
+        
+        % Testing each method for the same data
+        for it_meth = 1:length(methods)
+            rng(seed) % resetting the seed before starting the method
+            fct = str2func(methods{it_meth});
+            [T, RMSE , MERs] = fct(X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit, config);
+            
+            times{it_meth} = T; % overwritting previous record times for simplicity's sake
+%             figure
+            for measure = 1:(1+config.numSubSensor)
+                RMSEs{it_meth,measure} = cat(1, RMSEs{it_meth,measure}, RMSE(measure,:)); % saving the RMSE 
+%                 % RMSEs{it_meth,1} is for the offsets, RMSEs{it_meth,>1} is
+%                 % for the gains
+%                 plot(MERs(measure,:))
+%                 hold on
+            end
+        end
+    end
+    if config.display
+        figure
+        for it_meth = 1:length(methods)
+            for measure = 1:(config.numSubSensor+1)
+                subplot(length(methods),config.numSubSensor+1,(it_meth-1)*(config.numSubSensor+1)+measure)
+                for k = 1:size(RMSEs{it_meth,measure},1)
+                    semilogy(RMSEs{it_meth,measure}(k,:))
+                    hold on
+                end
+                hold off
+            end
+        end
+    end
+end
+

MATLAB/run_config_forget_other_sen.m

@@ -0,0 +1,73 @@
+function [res] = run_config(varargin)
+    if nargin == 0
+        file_name = "default_config.json";
+    elseif nargin == 1
+        file_name = varargin{1};
+    else
+        file_name = varargin{1};
+        warning("run_config only accepts one argument. Every other argument is ignored")
+    end
+    config = jsondecode(fileread(file_name));
+    methods = config.calibrationMethods;
+    iter_max = round(config.Tmax / config.delta_measure);
+    RMSEs = cell(length(methods),1+config.numSubSensor);
+    times = cell(length(methods),1);
+    
+    for run = 1:config.numRuns
+        % data generation
+        [X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit] = data_gen(config, run);
+        seed = floor(rand*10000);
+        
+        % We want to only calibrate the sensor of interest
+%         X = X(:,1:config.numSensor+1);
+%         X_theo = X_theo(:,1:config.numSensor+1);
+%         W = W(:,1:config.numSensor+1);
+%         F_theo = F_theo(:,1:config.numSensor+1);
+%         Omega_F = Omega_F(:,1:config.numSensor+1);
+%         Phi_F = Phi_F(:,1:config.numSensor+1);
+%         Finit = Finit(:,1:config.numSensor+1);
+
+        X = X(:,1:config.numSensor+1);
+        X_theo = X_theo(:,1:config.numSensor+1);
+        W = W(:,1:config.numSensor+1);
+        F_theo = F_theo(1:2,1:config.numSensor+1);
+        Omega_F = Omega_F(1:2,1:config.numSensor+1);
+        Phi_F = Phi_F(1:2,1:config.numSensor+1);
+        Finit = Finit(1:2,1:config.numSensor+1);
+        Omega_G = Omega_G(:,1:2);
+        Phi_G = Phi_G(:,1:2);
+        Ginit = Ginit(:,1:2);
+        config.numSubSensor = 1;
+        
+        % Testing each method for the same data
+        for it_meth = 1:length(methods)
+            rng(seed) % resetting the seed before starting the method
+            fct = str2func(methods{it_meth});
+            [T, RMSE , MERs] = fct(X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit, config);
+            
+            times{it_meth} = T; % overwritting previous record times for simplicity's sake
+%             figure
+            for measure = 1:(1+config.numSubSensor)
+                RMSEs{it_meth,measure} = cat(1, RMSEs{it_meth,measure}, RMSE(measure,:)); % saving the RMSE 
+%                 % RMSEs{it_meth,1} is for the offsets, RMSEs{it_meth,>1} is
+%                 % for the gains
+%                 plot(MERs(measure,:))
+%                 hold on
+            end
+        end
+    end
+    if config.display
+        figure
+        for it_meth = 1:length(methods)
+            for measure = 1:(config.numSubSensor+1)
+                subplot(length(methods),config.numSubSensor+1,(it_meth-1)*(config.numSubSensor+1)+measure)
+                for k = 1:size(RMSEs{it_meth,measure},1)
+                    semilogy(RMSEs{it_meth,measure}(k,:))
+                    hold on
+                end
+                hold off
+            end
+        end
+    end
+end
+

README.md

@@ -0,0 +1,22 @@
+Une version Python et une version MATLAB sont disponibles. La version MATLAB est la plus à jour. La version Python a été mise de côté compte tenu de la lenteur non expliquée de la méthode EM par rapport à la méthode MU.
+
+###### Version MATLAB ######
+
+Le code fonctionne sur la version R2018a.
+
+L'organisation est la suivante :
+
+MATLAB
+	├─ calib_meth				# Contient les différentes méthodes d'étalonnage
+	│	├─ common_nesterov		# Contient des fonctions communes à EMNeNMF et REMNeNMF
+	│	├─ EMNeNMF				# Contient la méthode EMNeNMF
+	│	├─ INCAL				# Contient la méthode INCAL
+	│	└─ REMNeNMF				# Contient la méthode REMNeNMF
+	│
+	├─ common
+	│	└─ data_gen.m			# Fonction permettant de générer les données
+	│
+	├─ config.json				# Fichier permettant de configurer les données à générer et les méthodes
+	└─ main.m					# Fichier exécutant le test en fonction de config.json
+	
+Pour configurer un test, les variables à paramétrer sont les suivantes :