function h = get_height_map(n, mask)
mask1 = mask(:,:,1);

%nz = non_zero_seq_num(mask);


%[rr cc] = find(mask);
%pixels = size(rr,1)

rows = size(n, 1)
cols = size(n, 2)
pixels = rows*cols

indices = zeros(rows, cols);
sum = 1;
for r=1:rows
    for c=1:cols
        indices(r,c) = sum;
        sum = sum+1;
    end
end

non_zero = prod(size(mask1));
eq_zero = pixels - non_zero;

A = sparse([],[],[], pixels + eq_zero+2, pixels, 3*pixels);
B = sparse([],[],[], pixels + eq_zero+2, 1, 2*pixels);

pb = waitbar(0,'starting');

rrrr = 0;
r = 1;
for x=1:cols
    waitbar(x/cols, pb, sprintf('%d % complete', 100*x/cols));
    for y=1:rows
        if (mask1(y,x) < 255)
            A(r, indices(y,x)) = 1;
            B(r, 1) = 0;
            r = r+1;
        else
            % has right neighbor
            if (x+1 <= cols && mask1(y, x+1) == 255)
                A(r, indices(y, x+1)) =  n(y,x,3);
                A(r, indices(y, x))   = -n(y,x,3);
                B(r, 1)               = -n(y,x,1);
                r = r+1;
                rrrr = rrrr + 1;
                if rrrr == 1000
                    A(r, indices(y, x)) = 1;
                    B(r, 1) =             40;
                    r = r + 1
                end
            end
            
            % has upper neighbor
            if (y-1 > 0 && mask(y-1, x) == 255)
                A(r, indices(y-1, x)) =  n(y,x,3);
                A(r, indices(y,   x)) = -n(y,x,3);
                B(r, 1)               = -n(y,x,2);
                r = r+1;
            end
        end
    end
end

2*pixels - r
close(pb);

hh=A\B;

h = zeros(340,512);
for r=1:340
    h(r,1:512)=hh((r-1)*512+1:r*512);
end


