NxM DCT matrix generation
type dctmtx


you will get:

function c = dctmtx(n) 
[cc,rr] = meshgrid(0:n-1); 
c = sqrt(2 / n) * cos(pi * (2*cc + 1) .* rr / (2 * n)); 
c(1,:) = c(1,:) / sqrt(2); 


So it's not that hard to generate an non-square DCT matrix, just with minor revision:

function c = dctmtx1(n,m) 
[cc,rr] = meshgrid(0:n-1, 0:m-1); 
c = sqrt(2 / n) * cos(pi * (2*cc + 1) .* rr / (2 * n)); 
c(1,:) = c(1,:) / sqrt(2);