Applying a "keep top N values" mask to image, in blocks

Question

I have troubles understanding how to apply a function on a block of data / submatrix of a matrix.

My task is to take an image, divide it into 8x8 blocks, and from each block pick 8 biggest values and set the rest to zero. I know the way could be through the for loops, but I would like to learn how to apply function on block of data.

Currently, I am applying a function to whole 256x256 matrix, but I would need to find a way how to apply it on just each block 8x8. Here is my commented code:

%% DCT transformation
I = imread('cameraman.tif');
I = im2double(I);
T = dctmtx(8); %returns the 8-by-8 DCT transform matrix
dct = @(block_struct) T * block_struct.data * T';
B = blockproc(I,[8 8],dct);
% Here I want to apply my function applyMask to blocks of 8x8 individualy
%this function will take a block 8x8, sort it, pick 8 biggest ones, save
%them and set rest to zero
f = @applyMask;
b = f(B)
function y = applyMask(x)
vector = x(:); %retransform matrix to be sorted
y=zeros(8,8)   %prepare matrix where 8 biggest values will be stored (rest is zero)
sorted = transpose(sort(vector,'descend')) %sort vecotr in descending order
pick = sorted(1:1, 1:8) %takes just first 8 biggest values 
for k=1 : 8
for i=1 : 8
    for j=1 : 8
        if ((x(i,j)==pick(1,k)) && nnz(y)<8 ) %if there is one of the 8 biggest - keep
           y(i,j)= pick(1,k) %put 8 biggest values to matrix
        end
    end
end

end
end

Answer 1

Your code can benefit from vectorization (i.e. removal of for loops).

function C = q52688681
%% DCT transformation
I = imread('cameraman.tif');
I = im2double(I);
T = dctmtx(8); %returns the 8-by-8 DCT transform matrix
B = blockproc(I,[8 8], @(block_struct) T * block_struct.data * T.');
C = blockproc(I,[8 8], @applyMask);

function out = applyMask(img)
NMAX = 8;
out = zeros(size(img.data));
[~,idx] = maxk(img.data(:), NMAX);
out(idx) = img.data(idx);

If your MATLAB version is >= R2017b, you can use maxk, otherwise:

function out = applyMask(img)
NMAX = 8;
out = zeros(size(img.data));
[~,idx] = sort(img.data(:), 'descend');
out( idx(1:NMAX) ) = img.data( idx(1:NMAX) );

And you can further reduce the amount of code and computations by doing this:

function B = q52688681
NMAX = 8;
I = im2double(imread('cameraman.tif'));
B = blockproc(I, [NMAX NMAX], @(x)applyMask(x, NMAX, dctmtx(NMAX)) );

function out = applyMask(blk, nmax, T)
img = T * blk.data * T.';
out = zeros(size(img));
[~,idx] = sort(img(:), 'descend');
out( idx(1:nmax) ) = img( idx(1:nmax) );

Let me know in the comments if there's anything unclear about this code, and I'll try to explain it.