% This function creates a decimation operator
function decimate = matrix_decimation_func(SuS)
% SuS - size of the image boundariess
    smp_size =(SuS - 1)/2;
    [ out_row, out_col ] = ...
        meshgrid(1:smp_size(1), 1 : smp_size(2));
    %create grid of vector 1:smp_size(1), to make the base for image
    in_col = out_col + out_col;
    %for faster pipeline performances
    in_row = out_row+out_row;   
     %for faster pipeline performances
     
    outputInd = sub2ind(smp_size, out_row, out_col);
    
    inputInd = sub2ind(SuS, in_row, in_col);
    calc_assist=ones(numel(outputInd), 1);
    decimate = sparse(outputInd, inputInd, ...
     calc_assist , prod(smp_size), prod(SuS));

end