% Usage: Y = blkdct2(y, m, n).
% Block DCT routine on array y, using m by n blocks. m must divide
% number of rows in y, n divide number of colums. Transform returned in
% y. For small blocks it's faster to use matrix multiplies than
% to repeatedly call the "dct2" command.
% If you don't have the signal processing toolbox substitute kdct() and
% kidct() for dct() and idct().
function Y = blkdct2(y, m, n)
% for small blocks the dct2 is computed faster by matrix multiplication
% that by using the dct2 command
[r,s] = size(y);
Y = zeros(r,s);
C1 = dct(eye(m));
C2 = dct(eye(n))';
for i=1:(r/m)
for j=1:(s/n)
e = (i-1)*m+1;
f = (j-1)*n+1;
Z = double(y(e:(e+m-1),(f:f+n-1)));
Y(e:(e+m-1),(f:f+n-1)) = C1*Z*C2;
end;
end;