function auto_regions(reg_x, reg_y, Xm)
% function auto_regions(reg_x, reg_y, Xm)
%
% A function that takes the number of regions in the x-direction
Xm=round(Xm);
R = uint8(pts2binim(Xm));
xr = range(Xm(:,2));
yr = range(Xm(:,1));
dx = round((xr(2)-xr(1))/reg_x);
dy = round((yr(2)-yr(1))/reg_y);
layer=1;
sx = xr(1);
corner_pts = [];
regions = {};

idx = linspace(xr(1), xr(2), reg_x+1);
idy = linspace(yr(1), yr(2), reg_y+1);

[pts_x, pts_y] = meshgrid(idx, idy);
pts_x = round(pts_x); pts_y = round(pts_y);
%G = latticegrid(size(pts_x,1), size(pts_x,2), 8) - latticegrid(size(pts_x,1), size(pts_x,2), 4);
%G = latticegrid(size(pts_x,1), size(pts_x,2), 8);
G = latticegrid(size(pts_x,1), size(pts_x,2), 4);

figure(1); clf;
imagesc(R)
hold on
plot(pts_x(:), pts_y(:), 'wo');
%drawgrid([pts_x(:), pts_y(:)], G, gca);
axis image

for i=1:reg_x
    sy = yr(1);
    for j=1:reg_y
        
        
        s = R(rect(2):rect(2)+rect(4), rect(1):rect(1)+rect(3));

        R(sy:sy+dy, sx:sx+dx) = R(sy:sy+dy, sx:sx+dx)*layer;
        %        subim = R(sy:sy+dy, sx:sx+dx);
        %        if sum(subim(:))>0
        %            edge = bwboundaries(subim);
        %            edge = fliplr(edge{1});
        %            edge = edge + ones(size(edge,1),1)*[sx, sy];
        %            region.edge = edge;
        %            regions = cat(1, regions, region);
        %        end
        %        cpts = getRegionCornerPts(R(sy:sy+dy, sx:sx+dx));
        %        cpts = cpts + ones(size(cpts,1),1)*[sx, sy];
        %        corner_pts = cat(1, corner_pts, cpts);
        sy = sy + dy + 1;
        layer = layer + 1;
    end
    sx = sx + dx + 1;
end
% remove the duplicates
d = matdistance(corner_pts, corner_pts);
idx = [];
for i=1:size(d,1)
    for j=i:size(d,1)
        if d(i, j)<5 && d(i,j)>0
            idx = cat(2, idx, j);
        end
    end
end
corner_pts(idx,:) = [];
imagesc(R); hold on
plot(corner_pts(:,1), corner_pts(:,2), 'ro')

return;




% function [R, E2, E, edge_petal, dx, dy]=auto_regions(reg_x, reg_y, Z)
% Z=round(Z);
% p = uint8(pts2binim(Z));
% xr=range(find(any(p)));
% yr=range(find(any(p')));
% dx=ceil(diff(xr)/reg_x);
% dy=ceil(diff(yr)/reg_y);
% edge_petal = [];
% waitbar_h = waitbar(0, 'Calculating grid...');
% for i=1:reg_x,
%     for j=1:reg_y,
%         r=logical(zeros([size(p,1) size(p,2)]));
%         r(yr(1)+dy*(j-1)+(1:dy),xr(1)+dx*(i-1)+(1:dx))=logical(ones(dy,dx));
%         r=r&p;
%         if any(r(:))
%             layer=layer+1;
%             R(:,:,layer) = r;
%             p=xor(r,p);
%             edge=[];
%             vertex=[];
%             mask = imdilate(r, ones(10));
%             [x, y] = find(bwperim(mask));
%             in = find(inpolygon(Eptl(:,1), Eptl(:,2), y, x)==1);
%             edge = [edge; Eptl(in,1:2)];
%             edge_petal = [edge_petal; [Eptl(in,1:2) ones(length(in),1)*layer]];
%             for k=0:dx,
%                 if ~ r(yr(1)+dy*(j-1)+(1),xr(1)+dx*(i-1)+(k))== r(yr(1)+dy*(j-1)+(1),xr(1)+dx*(i-1)+(k+1)),
%                     if r(yr(1)+dy*(j-1)+(1),xr(1)+dx*(i-1)+(k))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(k),yr(1)+dy*(j-1)+(1)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(k),yr(1)+dy*(j-1)+(1)]);
%                     elseif r(yr(1)+dy*(j-1)+(1),xr(1)+dx*(i-1)+(k+1))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(k+1),yr(1)+dy*(j-1)+(1)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(k+1),yr(1)+dy*(j-1)+(1)]);
%                     end;
%                 end;
%                 if ~ r(yr(1)+dy*(j-1)+(dy),xr(1)+dx*(i-1)+(k))== r(yr(1)+dy*(j-1)+(dy),xr(1)+dx*(i-1)+(k+1)),
%                     if r(yr(1)+dy*(j-1)+(dy),xr(1)+dx*(i-1)+(k))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(k),yr(1)+dy*(j-1)+(dy)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(k),yr(1)+dy*(j-1)+(dy)]);
%                     elseif r(yr(1)+dy*(j-1)+(dy),xr(1)+dx*(i-1)+(k+1))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(k+1),yr(1)+dy*(j-1)+(dy)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(k+1),yr(1)+dy*(j-1)+(dy)]);
%                     end
%                 end;
%             end;
%
%             for k=0:dy,
%                 if ~ r(yr(1)+dy*(j-1)+(k),xr(1)+dx*(i-1)+(1))== ...
%                         r(yr(1)+dy*(j-1)+(k+1),xr(1)+dx*(i-1)+(1)),
%                     if r(yr(1)+dy*(j-1)+(k),xr(1)+dx*(i-1)+(1))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(1),yr(1)+dy*(j-1)+(k)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(1),yr(1)+dy*(j-1)+(k)]);
%                     elseif r(yr(1)+dy*(j-1)+(k+1),xr(1)+dx*(i-1)+(1))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(1),yr(1)+dy*(j-1)+(k+1)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(1),yr(1)+dy*(j-1)+(k+1)]);
%                     end;
%                 end;
%
%                 if ~ r(yr(1)+dy*(j-1)+(k),xr(1)+dx*(i-1)+(dx))== ...
%                         r(yr(1)+dy*(j-1)+(k+1),xr(1)+dx*(i-1)+(dx)),
%                     if r(yr(1)+dy*(j-1)+(k),xr(1)+dx*(i-1)+(dx))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(dx),yr(1)+dy*(j-1)+(k)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(dx),yr(1)+dy*(j-1)+(k)]);
%                     elseif r(yr(1)+dy*(j-1)+(k+1),xr(1)+dx*(i-1)+(dx))==1,
%                         edge=cat(1,edge,[xr(1)+dx*(i-1)+(dx),yr(1)+dy*(j-1)+(k+1)]);
%                         vertex=cat(1,vertex,[xr(1)+dx*(i-1)+(dx),yr(1)+dy*(j-1)+(k+1)]);
%                     end;
%                 end;
%             end;
%             E{layer} = unique(edge, 'rows');
%             E2{layer} = unique(vertex, 'rows');
%             waitbar((i*j)/(reg_x*reg_y), waitbar_h);
%         end
%     end
%
% end
% close(waitbar_h);