function [signal_in, signal_tex, signal_car, signal_out, psnr,  psnr_t, psnr_c] = decomposition_2D_impl(image, ridkost_k1, ridkost_k2, rozmer,original_picture, baze_1, baze_2, wavelet_1, wavelet_2, hloubka_decom_1, hloubka_decom_2)
%
global M;
global N;
disp(strcat('image: ', image));
disp(['Ridkost k1: ' num2str(ridkost_k1)]);
disp(['Ridkost k2: ' num2str(ridkost_k2)]);
disp(['Rozmer: ' num2str(rozmer)]);

disp(strcat('Real: ', original_picture));
disp(strcat('Baze: ', baze_1));
disp(strcat('Baze: ', baze_2));
disp(strcat('Wavelet: ', wavelet_1));
disp(strcat('Wavelet: ', wavelet_2));
disp(['Hloubka Dekompozice: ' num2str(hloubka_decom_1)]);
disp(['Hloubka Dekompozice: ' num2str(hloubka_decom_2)]);




J_1 = hloubka_decom_1;
J_2 = hloubka_decom_2;

L = rozmer;
k1 = ridkost_k1;
k2 = ridkost_k2;

%% devise of picture
velikost_deleni = 16;
switch  L < velikost_deleni 
    case true
        velikost_deleni = L;
        m_ = 1;
        n_ = 1;
    case false
        m_ = M/velikost_deleni;
        n_ = N/velikost_deleni;
      
end

%% Base
 
switch baze_1
    case 'DCT'
    U_1 = dctmtx(velikost_deleni)';                 % baze pro zpracovani 
    U_1_or = dctmtx(L)';                            % baze pro vygenerovani ridke komponenty
     
    case 'DWT'
    [Lo_D,Hi_D,Lo_R,Hi_R] = w_filters(wavelet_1); 
     U_1 = dwtmatrix(Lo_D,Hi_D,J_1,velikost_deleni)';               % baze pro zpracovani 
     U_1_or = dwtmatrix(Lo_D,Hi_D,J_1,L)';                          % baze pro vygenerovani ridke komponenty
   
    
end

switch baze_2
    case 'DCT'
    U_2 = dctmtx(velikost_deleni)';         % baze pro zpracovani 
    U_2_or = dctmtx(L)';                    % baze pro vygenerovani ridke komponenty
     
    case 'DWT'
    [Lo_D,Hi_D,Lo_R,Hi_R] = w_filters(wavelet_2); 
     U_2 = dwtmatrix(Lo_D,Hi_D,J_2,velikost_deleni)';               % baze pro zpracovani 
     U_2_or = dwtmatrix(Lo_D,Hi_D,J_2,L)';                          % baze pro vygenerovani ridke komponenty
    
end

%% Original_signal
switch image
    case 'nahodny'
        N = L;
        M = L;
        S = N*N;
        % vstupni obraz je souctem dvou komponent
        %1. komponenta
        x_signal1 = zeros(S, 1);
        perm = randperm(S); %random positions of nonzeros
        %vector of coordinates contains random values from [-1,1] at nonzero positions
        indexes = perm(1:k1);
        x_signal1(indexes) = randi(255,k1,1);    %2*rand(k1,1)-1;
        % max_freq = max(indexes) %maximum 'frequency'
        koef1 = reshape(x_signal1,N,N);

        signal1_in = U_1_or' * koef1 * U_1_or; 
                
       %2. komponenta
        x_signal2 = zeros(S, 1);
        perm = randperm(S); %random positions of nonzeros
                                                     %vector of coordinates contains random values from [-1,1] at nonzero positions
        indexes = perm(1:k2);
        x_signal2(indexes) = randi(255,k2,1);         %2*rand(k2,1)-1;
        % max_freq = max(indexes) %maximum 'frequency'
        koef2 = reshape(x_signal2,N,N);

        signal2_in = U_2_or' * koef2 * U_2_or;

        original_picture_n = signal1_in + signal2_in;
         
        signal_in = original_picture_n(:);
    
    case 'real'
    
            switch size(size(imread(original_picture)),2)
    
                case 3
     
                original_picture_n = im2double(rgb2gray(imread(original_picture)));

                case 2
                 original_picture_n = im2double(imread(original_picture));
            end

            original_picture_n = imresize(original_picture_n,[32,32]);

            M = size(original_picture_n,1);   % size of picture M x N 
            N = size(original_picture_n,2);
 
            signal_in = original_picture_n(:);
end



%% Dictionary
 A_1 = kron(U_1',U_1');
 A_2 = kron(U_2',U_2');
 
%% Decomposition

%inicializace    
signal_tex = zeros(M,N);
signal_car = zeros(M,N);
for i = 1:m_
    for j = 1:n_
        
    
        m_dilci = ((i-1)*velikost_deleni + 1):i*velikost_deleni;
        n_dilci = ((j-1)*velikost_deleni + 1):j*velikost_deleni;
        
        K = velikost_deleni*velikost_deleni; 
        
        part_of_picture = original_picture_n(m_dilci,n_dilci);
        y = part_of_picture(:);
        
           
        %% Initialisation  - UnLockBox
         % Author: Nathanael Perraudin
            % Date: 2 October 2013

%             % Loading toolbox
%             init_unlocbox();

            verbose = 1;    % verbosity level

        %% General parameter for the problem
            mu1=1;
            mu2=1;

        %% Creation of the Problem
           

            % Creation of the operator Psi
            Psi = @(x) A_1*x;            %1/sqrt(N)*fft(x);
            Psit = @(x) A_1'*x;          %sqrt(N)*ifft(x);
            nu_Psi = 1; % Frame bound

            % Creation of th operator Phi
            Phi = @(x) A_2*x;              %M.*Psi(x);
            Phit = @(x) A_2'*x;              %Psit(M.*x);
            nu_Phi = 1; % Frame bound


           % Creation of the measurements
            s = y;

        %% Define the UNLocBoX proximity operators and gradient

            %Define mask functions
            Mc = @(x) x(1:K);         % select the variable c            Mc = @(x) x(1:N);         % select the variable c
            Mb = @(x) x(K+1:end);     % select the variable b             Mb = @(x) x(N+1:end);     % select the variable b
            Cc = @(x,c) [c;Mb(x)];    % recreate all opt variable from c                Cc = @(x,c) [c;Mb(x)];    % recreate all opt variable from c
            Cb = @(x,b) [Mc(x);b];    % recreate all opt variable from b                Cb = @(x,b) [Mc(x);b];    % recreate all opt variable from b
            
        
            % Functional mu1*||c||_1
            param_l1_1.verbose = verbose - 1;
            f1.prox = @(x,T) Cc(x,prox_l1(Mc(x),T*mu1,param_l1_1));
            f1.eval = @(x)   mu1*norm(Mc(x),1);

            % Functional mu2*||c||_1
            param_l1_2.verbose = verbose - 1;
            f2.prox = @(x,T) Cb(x,prox_l1(Mb(x),T*mu1,param_l1_2));
            f2.eval = @(x)   mu2*norm(Mb(x),1);

            % Define the gradient of ||s - Psi c - Phi b||_2^2 
            g.grad = @(x)  [ Psit( Psi(Mc(x)) + Phi(Mb(x)) - s ) ; ...
                 Phit( Psi(Mc(x)) + Phi(Mb(x)) - s ) ];
            g.eval = @(x)    norm((s-(Psi(Mc(x))+Phi(Mb(x))))/sqrt(2));


        %% Define solver parameter

            paramfb.maxit=1000;
            paramfb.gamma=1/(nu_Phi+nu_Psi);
            paramfb.verbose=verbose;

        % Solve the problem
            x0= eps*ones(2*K,1);
            sol=generalized_forward_backward(x0,{f1,f2},g,paramfb);
            
            uv = reshape(sol,K,2);
sol_1 = uv(:, 1);
sol_2 = uv(:, 2);
            
signal_tex_v = A_1*sol_1;
signal_car_v = A_2*sol_2;      %reconstructed sparse vector of coordinates  




signal_tex(m_dilci,n_dilci) = reshape(signal_tex_v,velikost_deleni,velikost_deleni);
signal_car(m_dilci,n_dilci) = reshape(signal_car_v,velikost_deleni,velikost_deleni);
     
        

    end
    
     

end
disp('konec');

  
 %% PSNR
signal_out = signal_tex(:) + signal_car(:);
signal_in    = double(signal_in);
signal_out = double(signal_out(:));
absD = abs(signal_in-signal_out);
moc    = absD.^2;
mse  = sum(moc(:))/numel(signal_in);
psnr = 10*log10((max(signal_in)-min(signal_in))^2/mse);

switch image
    case 'nahodny'
        
 %% PSNR_T

signal1_in = double(signal1_in(:));
signal_tex = double(signal_tex(:));
absD = abs(signal1_in-signal_tex);
moc    = absD.^2;
mse  = sum(moc(:))/numel(signal1_in);
psnr_t = 10*log10((max(signal1_in)-min(signal1_in))^2/mse);

 %% PSNR_C

signal2_in    = double(signal2_in(:));
signal_car = double(signal_car(:));
absD = abs(signal2_in-signal_car);
moc    = absD.^2;
mse  = sum(moc(:))/numel(signal2_in);
psnr_c = 10*log10((max(signal2_in)-min(signal2_in))^2/mse);

    case 'real'
        psnr_t = 0;
        psnr_c = 0;
end

