function [m] = emap(f,n)
% entropy map using x-n:x+n stencil, n=1 is probably fine
% f is either normalized iteration map or distance estimator map
% see emapOld.m for more readable snail-code

if(nargin<2)
    n=1;
end
[n1,n2] = size(f);
% subgrids
tmp = f(n+1:end-n,n+1:end-n);
g = zeros([size(tmp) (2*n+1)^2]);
k=1;
for i=-n:n
    for j=-n:n
        g(:,:,k) = f(n+1+i:end-n+i,n+1+j:end-n+j);
        k = k+1;
    end
end
[k1 k2 k3] = size(g);

g = reshape(g,k1*k2,k3);
gmin = min(g,[],2);
g = g - repmat(gmin,1,k3);
gsum = sum(g,2);
ind2 = find(gsum==0);
gsum(ind2) = 1;
g = g./repmat(gsum,1,k3);
ind = find(g(:)>0);
logg = g(:);
logg(ind) = log(g(ind));
glogg = g(:) .* logg;
glogg = reshape(glogg,k1*k2,k3);
m2 = sum(glogg,2)/log(k3);
m2 = reshape(m2,k1,k2);
% expand to original size
m = f;
m(2:end-1,2:end-1) = m2;
m(1,2:end-1) = m(2,2:end-1);
m(end,2:end-1) = m(end-1,2:end-1);
m(2:end-1,1) = m(2:end-1,2);
m(2:end-1,end) = m(2:end-1,end-1);
m(1,1) = m(2,2);
m(1,end) = m(2,end-1);
m(end,1) = m(end-1,2);
m(end,end) = m(end-1,end-1);