function [D,X] = ksvd_upravene(param)
%KSVDM simple implementation od K-SVD algorithm
%   [D] = ksvdm(param)
%     reguired parame ters is:
%           param.data - matrix with training data (Matice Y)
%           param.initdic - matrix with training dictionary (maticeD
%           startovaci)
%           param.Tdata - sparsity threshold (ridkost X)
%           param.N - size of dictionary (pocet sloupcu v D)
%           param.iternum - number of iterations

%t = clock;
%% initialization - nacte, trenovaci signaly, normalizuje slovnik, 
Y = param.data;
D = normcols(param.initdic);
T = param.Tdata; 
K = param.N;
[a,K] = size(D); %ulozi rozmery D do promenych a - radky,K sloupce
iter = param.iternum;
%  par.StopMaxerror=1e-2;
par.StopMaxiter=T; %max pocet iteraci je pocet nulovych hodnot x (ridkost)


%% main iteration loop
for j=1:iter %j od 1 po pocet iteraci
    %% 1. step - sparse coding stage
    % G is Gramm matrix, for OMP
    %G = D'*D;
    %X = omp1(D,Y,G,T);
    % Gamma = omp1(D,X,[],T); %faster OMP, more memory
%%  
    if j~=1 %j ruzno od 1 
        Xp = X; %previous X used to compute the error 
    end;
    X=zeros(K,size(Y,2)); %prazdna matice X o velikosti nxp (kde n je K a p je pocet sloupcu Y) 
    for i=1:size(Y,2) %odhad X pomoci OMP - pro i od 1 do poctu sloupcu v Y
        X(:,i)=omp(D,Y(:,i),par); % i-ty sloupec z X pouzij OMP
    end
    
        %error computing - pocitani chyby
    abs_error(j)=norm(Y-D*X,'fro'); %celkov? chyba pro jty krok frobeinova norma rozdilu
    if j==1 
        rel_error(j)=0; % pro prvni j-ty krok je 0 
    else
        rel_error(j)=(norm(Y-Dp*Xp,'fro')-norm(Y-D*X,'fro'))/norm(Y-Dp*Xp,'fro'); % eror v procentech proti posledni iteraci j
    end;


    %% 2. step - dictionary update stage
    for k=1:K %column loop - iterace po sloupcich slovniku
%         vekt_dkxt = D(:,k)*Gamma(k,:);

        %fixing different sparsity of each column in X, that could make problems
        b=length(find(X(:,k))); %%difference between wanted and real sparsity b = pocet nenulovych prvku v sloupci
        if b<T %b je mensi jako pozadovana ridkost
            for c=1:T-b %c od 1 do ridkost - pocet nenulovych prvku ve sloupci
                for i=c:K %pristup k prvku ve sloupci      
                    if X(i,k)==0    %%replacement of necessary number of zeros with very small numbers
                        X(i,k)=1e-15; %pokud neni nulovy zmeni X na velice male cislo
                        break;
                    end;
                end;
            end;
        end;
        
        gamma_k = find(X(:,k)); %gamma_k index na kolikatem miste v k-tym sloupci jsou hodnoty
        
        %mask of size N x |omega_k|
        Mask = zeros(size(Y,2),length(gamma_k)); %mask = nulov? matice o velikosti sloupcu Y x delky gamma_k
        
        %
        for iii=1:length(gamma_k) %od 1 do velikosti gamma_k = pos sloupcich 
            Mask(gamma_k(iii),iii) = 1; 
        end
        
        %remove k-th row from gamma - odstrani kty radek z X_j
        X_j = X; 
        X_j(k,:) = []; 
        
        %remove k-th column from dictionary
        D_j = D;
        D_j(:,k) = []; %odstrani kty sloupec z D_j
        d_jxt = D_j * X_j; %d_jxt je soucin ochuzeneho slovniku a ochuzeneho X

        Ek = Y - d_jxt; %
        %masking by omega_k
        Ek = Ek*Mask; % pokd je mask v kazdem kroku eye tak proc ?

        %SVD
        [U,S,V]= svd(Ek, 0); %svd rozklad Ek 
     %   U(U(:,1)<0)=0;   U(:,1)=U(:,1)/norm(U(:,1));    %%zku?ebn? upraven? atom? slovn?ku
        %update
        D_next(:,k) = U(:,1);
%         V = transpose(V); %%???
        Gamma_next(k,:) = V(:,1)*S(1,1);
    end
    %end of column loop

    %update dictionary *??
    Dp=D; %previous dictionary used to compute the error
    D = D_next;
    
    D_next = [];
    Gamma_next = [];
%     clearvars D_next Gamma_next;
    
    %disp(sprintf('%d-th iteration',j));
    
    %error computing
    abs_error(j)=norm(Y-D*X);
    if j==1 
        rel_error(j)=0;
    else
        rel_error(j)=(norm(Y-Dp*Xp)-norm(Y-D*X))/norm(Y-Dp*Xp);
    end;
  
    
    
end
%plot(abs_error)
%title('Celkova zavislost chyby v zavislosti na iteracich od pocatku')
%xlabel('Pocet iteraci')
%ylabel('Velikost chyby [-]')
%figure, plot(rel_error)
%title('Celkova zavislost chyby v zavislosti na iteracich k predchozi iteraci')
%xlabel('Pocet iteraci')
%ylabel('Velikost chyby [%]')
%axis tight
%time = etime(clock,t)/60;
end

