clear
IT = 1000; % Number of iterations
R = 3:7; % Sparsity of the signals
SN = 5; % Number of trials
M = 20; % Signal space dimension
N = 40; % Ambient space dimension
L = 32*N; % Number of training signals

ParamGen(M,N,L,R,SN) % Generating generic problems

%%%% Dictionary recovery algorithms
for sn = 1:SN    
    for r = R
        mmdlcn_exactRec_demo(IT,r,sn,'un')
        mapcn_exactRec_demo(IT,r,sn,'bn')
        ksvd_exactRec_demo(IT,r,sn)
        modcn_exactRec_demo(IT,r,sn)
    end
end
%%%%%%%
Tre = .99;
for k = R
    for i=1:SN
        load(['RDLl1',num2str(M),'t',num2str(IT),'ikiun',num2str(k),'v1d',num2str(i),'.mat'])
        nrRDLu(i,k-R(1)+1) = sum(max(abs((Phid'*Phi)))>=Tre);
        %%%%%%
        load(['MAPl1',num2str(M),'t',num2str(IT),'ikiun',num2str(k),'v1d',num2str(i),'.mat'])
        nrMAP(i,k-R(1)+1) = sum(max(abs((Phid'*Phi)))>=Tre);        
        %%%%%%
        load(['KSVDl1',num2str(M),'t',num2str(IT),'ikiun',num2str(k),'v1d',num2str(i),'.mat'])
        nrKSVD(i,k-R(1)+1) = sum(max(abs((Phid'*Phi)))>=Tre);
        %%%%%%
        load(['MODl1',num2str(M),'t',num2str(IT),'ikiun',num2str(k),'v1d',num2str(i),'.mat'])
        nrMOD(i,k-R(1)+1) = sum(max(abs((Phid'*Phi)))>=Tre);
    end
end
clf
errorbar(R+.01,10*mean(nrRDLu,1)/4,std(nrRDLu,0,1),'k-.')
hold on
errorbar(R-.01,10*mean(nrKSVD,1)/4,std(nrKSVD,0,1),'r*-')
errorbar(R+.01,10*mean(nrMOD,1)/4,std(nrMOD,0,1),'b-v')
errorbar(R-.01,10*mean(nrMAP,1)/4,std(nrMAP,0,1),'b-^')
axis([2.5 7.5 15 105]);
title('Constrained column-norm dictionaries')
xlabel('Sparsity (# of non-zero elements in each coefficient vector)')
ylabel(['Average percents of exact recovery',sprintf('\n'),'after ',num2str(IT),' iterations'])
grid on
legend('MM Unit norm','K-SVD','MOD','MAPbased-DL');
axis([2.8 7.2 -5 105])