function regions = image2regions(imgFile, tmpDir, overlap, maxEig)
%%%%%%%%%%%%%%
%%%
%%% Computes the UCM regions of an image in a format compatible with the
%%% Java annotation tool. The output is a matrix of the same dimensions as
%%% the original image and type uint16. Pixels internal to a region have a
%%% unique ID of that region (starting with 1), and pixels on the contour
%%% edges have a value of 16384 + edge strength where edge strength is
%%% 1..255).
%%%
%%%%%%%%%%%%%
addpath debug/BSR/grouping/lib;

%maxImg = 10000*10000;

if nargin<4 maxEig=480*320; end
if nargin<3 overlap=50; end
if nargin<2 tmpDir='tmp'; end

if ~exist(tmpDir,'file')
   mkdir(tmpDir); 
end


%% prepare data
[foo,stem] = fileparts(imgFile);

img=imread(imgFile);

[tx,ty,nch]=size(img);
img=imresize(img,2.0,'nearest');

num_pix = size(img,1)*size(img,2);
rsz = 1/sqrt(num_pix / maxEig);

fprintf('rsz=%4.2f\n',rsz);

imwrite(img(1:tx+overlap, 1:ty+overlap,:),          [tmpDir '/' stem '_img1.bmp']);
imwrite(img(1:tx+overlap, ty-overlap+1:end,:),      [tmpDir '/' stem '_img2.bmp']);
imwrite(img(tx-overlap+1:end, 1:ty+overlap,:),      [tmpDir '/' stem '_img3.bmp']);
imwrite(img(tx-overlap+1:end, ty-overlap+1:end,:), 	[tmpDir '/' stem '_img4.bmp']);


%% gPb in quadrants

for i = 1 : 4
    gPb2_orient_all{i} = globalPb(sprintf('%s/%s_img%d.bmp',tmpDir,stem,i), sprintf('%s/%s_img%d.mat',tmpDir,stem,i), rsz); %#ok<*AGROW>
end
%save([outFile(1:end-4) '_all.mat'],'gPb2_orient_all');

%% assemble output with blending

mask = zeros(tx+overlap,ty+overlap);

gPb2_orient = zeros(2*tx,2*ty,8);
mask_all = gPb2_orient;
for o = 1: 8
    mask(end,:)=1;mask(:,end)=1; dst=bwdist(mask); dst(dst>2*overlap)=2*overlap; dst=dst+1; dst=dst/max(dst(:));
    gPb2_orient(1:tx+overlap, 1:ty+overlap,o)= dst.*double(gPb2_orient_all{1}(:,:,o))+gPb2_orient(1:tx+overlap, 1:ty+overlap,o);
    mask_all(1:tx+overlap, 1:ty+overlap,o)= dst+mask_all(1:tx+overlap, 1:ty+overlap,o);
    
    mask(:)=0;  mask(:,1)=1;mask(end,:)=1; dst=bwdist(mask); dst(dst>2*overlap)=2*overlap; dst=dst/max(dst(:));
    gPb2_orient(1:tx+overlap, ty-overlap+1:end,o)=dst.*double(gPb2_orient_all{2}(:,:,o))+gPb2_orient(1:tx+overlap, ty-overlap+1:end,o);
    mask_all(1:tx+overlap, ty-overlap+1:end,o)=dst+mask_all(1:tx+overlap, ty-overlap+1:end,o);
    
    mask(:)=0;  mask(1,:)=1;mask(:,end)=1; dst=bwdist(mask); dst(dst>2*overlap)=2*overlap; dst=dst/max(dst(:));
    gPb2_orient(tx-overlap+1:end, 1:ty+overlap,o)=dst.*double(gPb2_orient_all{3}(:,:,o))+gPb2_orient(tx-overlap+1:end, 1:ty+overlap,o);
    mask_all(tx-overlap+1:end, 1:ty+overlap,o)=dst+mask_all(tx-overlap+1:end, 1:ty+overlap,o);
    
    mask(:)=0;  mask(:,1)=1;mask(1,:)=1; dst=bwdist(mask); dst(dst>2*overlap)=2*overlap; dst=dst/max(dst(:));
    gPb2_orient(tx-overlap+1:end, ty-overlap+1:end,o)=dst.*double(gPb2_orient_all{4}(:,:,o))+gPb2_orient(tx-overlap+1:end, ty-overlap+1:end,o);
    mask_all(tx-overlap+1:end, ty-overlap+1:end,o)=dst+mask_all(tx-overlap+1:end, ty-overlap+1:end,o);
    
end

gPb2_orient = gPb2_orient./mask_all;

for o = 1: 8
    gPb_orient(:,:,o)=imresize(gPb2_orient(:,:,o),0.5); %#ok<AGROW>
end

ucm = uint8(contours2ucm(gPb_orient)*255);

regions = uint16(bwlabel(~ucm)+1);
regions(ucm>0) = uint16(ucm(ucm>0))+16384;

for i =1:4
    delete(sprintf('%s/%s_img%d.bmp',tmpDir,stem,i));
    delete(sprintf('%s/%s_img%d.mat',tmpDir,stem,i));
end