根据位置将两个列向量排序为 3D 矩阵

Question

使用 MATLAB 中的 imfindcircles 函数跟踪两个图像中的圆。我从大约变形的圆网格开始。我正在尝试将 imfindcircles 中的两列向量排序为矩阵，以便相邻的圆是矩阵中的相邻元素。第一张图片，圆圈符合网格，以下代码有效：

[centXsort,IX] = sortrows(centres1,1); %sort by x 
centYsort =zeros(289,2); %preallocate
for i = 1:17:289
    [sortedY,IY] = sortrows(centXsort(i:i+16,:),2); %sort by y within individual column
    centYsort(i:i+16,:) = sortedY; 
end
cent1mat = reshape(centYsort,17,17,2); %reshape into centre matrices

这不适用于第二张图片，因为一些圆圈在 x 或 y 方向重叠，但相邻的圆圈从不重叠。这意味着在第二组矩阵中，排序后相邻的圆圈不是相邻元素。

有没有办法将散点近似为矩阵？

Answer 1

这个答案并非在所有情况下都适用，但对于分数变化不大的情况似乎已经足够了。

我的想法是从网格角开始，沿着矩阵的外侧对角线工作，尝试 "grab" 最近的点，这些点看起来适合基于任何周围点的网格点我们已经捕获了。

您需要提供：

网格中的行数 (rows) 和列数 (cols)。
您的数据点 P 排列在 N x 2 数组中， 重新缩放 到 [0,1] x [0,1] 上的单位正方形。（我假设您可以通过目视检查原始数据的角点来做到这一点。）
一个权重参数edge_weight告诉算法应该将多少边界点吸引到网格边界。一些测试表明 3-5 左右是好的值。

包含测试用例的代码：

%// input parameters
rows = 11;     
cols = 11;     
edge_weight = 4;

%// function for getting squared errors between the points list P and a specific point pt
getErr =@(P,pt) sqrt(  sum( bsxfun(@minus,P,pt(:)').^2, 2 )  );    %'

output_grid = zeros(rows,cols,2);   %// output grid matrix
check_grid = zeros(rows,cols);      %// matrix flagging the gridpoints we have covered
[ROW,COL] = meshgrid(...            %// coordinate points of an "ideal grid"
    linspace(0,1,rows),...
    linspace(0,1,cols));


%// create a test case
G = [ROW(:),COL(:)];                       %// the actual grid-points
noise_factor = 0.35;                       %// noise radius allowed
rn = noise_factor/rows;
cn = noise_factor/cols;
row_noise = -rn + 2*rn*rand(numel(ROW),1);
col_noise = -cn + 2*cn*rand(numel(ROW),1);
P = G + [row_noise,col_noise];             %// add noise to get points


%// MAIN LOOP
d = 0;                
while ~isempty(P)                       %// while points remain...
    d = d+1;                            %// increase diagonal depth (d=1 are the outer corners)
    for ii = max(d-rows+1,1):min(d,rows)%// for every row number i...
        i = ii;
        j = d-i+1;                      %// on the dth diagonal, we have d=i+j-1          
        for c = 1:4                     %// repeat for all 4 corners
            if i<rows & j<cols & ~check_grid(i,j) %// check for out-of-bounds/repetitions

                check_grid(i,j) = true;         %// flag gridpoint
                current_gridpoint = [ROW(i,j),COL(i,j)];

                %// get error between all remaining points and the next gridpoint's neighbours
                if i>1
                    errI = getErr(P,output_grid(i-1,j,:));
                else
                    errI = edge_weight*getErr(P,current_gridpoint);
                end
                if check_grid(i+1,j)
                    errI = errI + edge_weight*getErr(P,current_gridpoint);
                end
                if j>1
                    errJ = getErr(P,output_grid(i,j-1,:));
                else
                    errJ = edge_weight*getErr(P,current_gridpoint);
                end
                if check_grid(i,j+1)
                    errJ = errJ + edge_weight*getErr(P,current_gridpoint);
                end

                err = errI.^2 + errJ.^2;


                %// find the point with minimal error, add it to the grid,
                %//     and delete it from the points list
                [~,idx] = min(err);                         
                output_grid(i,j,:) = permute( P(idx,:), [1 3 2] );
                P(idx,:) = [];

            end

            %// rotate the grid 90 degrees and repeat for next corner
            output_grid = cat(3, rot90(output_grid(:,:,1)), rot90(output_grid(:,:,2)));
            check_grid = rot90(check_grid);
            ROW = rot90(ROW);               
            COL = rot90(COL);

        end
    end
end

用边绘制结果点的代码：

%// plotting code
figure(1); clf; hold on;
axis([-0.1 1.1 -0.1 1.1])
for i = 1:size(output_grid,1)
    for j =  1:size(output_grid,2)
        scatter(output_grid(i,j,1),output_grid(i,j,2),'b')
        if i < size(output_grid,1)
            plot(   [output_grid(i,j,1),output_grid(i+1,j,1)],...
                [output_grid(i,j,2),output_grid(i+1,j,2)],...
                'r');
        end
        if j < size(output_grid,2)
            plot(   [output_grid(i,j,1),output_grid(i,j+1,1)],...
                [output_grid(i,j,2),output_grid(i,j+1,2)],...
                'r');
        end
    end
end

Answer 2

我开发了一个解决方案，它适用于我的情况，但可能不如某些人所需要的那样可靠。它需要 'square' 网格中的已知点数和点之间间距的粗略概念。我找到了点的 'AlphaShape' 和所有位于边缘的点。边缘向量被移动到从最小值开始，然后环绕一个矩阵，对应的点从顶点列表中被丢弃。对于大型点云来说可能不是最好的主意，但对我来说已经足够了。

R = 50; % search radius
xy = centres2;
x = centres2(:,1);
y = centres2(:,2);

for i = 1:8
    T = delaunay(xy); % delaunay
    [~,r] = circumcenter(triangulation(T,x,y)); % circumcenters
    T = T(r < R,:); % points within radius
    B = freeBoundary(triangulation(T,x,y)); % find edge vertices
    A = B(:,1);

    EdgeList = [x(A) y(A) x(A)+y(A)]; % find point closest to origin and rotate vector
    [~,I] = min(EdgeList);
    EdgeList = circshift(EdgeList,-I(3)+1);

    n = sqrt(length(xy)); % define zeros matrix
    matX = zeros(n); % wrap x vector around zeros matrix
    matX(1,1:n) = EdgeList(1:n,1);
    matX(2:n-1,n) = EdgeList(n+1:(2*n)-2,1);
    matX(n,n:-1:1) = EdgeList((2*n)-1:(3*n)-2,1);
    matX(n-1:-1:2,1) = EdgeList((3*n)-1:(4*n)-4,1);

    matY = zeros(n); % wrap y vector around zeros matrix
    matY(1,1:n) = EdgeList(1:n,2);
    matY(2:n-1,n) = EdgeList(n+1:(2*n)-2,2);
    matY(n,n:-1:1) = EdgeList((2*n)-1:(3*n)-2,2);
    matY(n-1:-1:2,1) = EdgeList((3*n)-1:(4*n)-4,2);

    centreMatX(i:n+i-1,i:n+i-1) = matX; % paste into main matrix
    centreMatY(i:n+i-1,i:n+i-1) = matY;

    xy(B(:,1),:) = 0; % discard values
    xy = xy(all(xy,2),:);   
    x = xy(:,1);
    y = xy(:,2);   

end

centreMatX(centreMatX == 0) = x;
centreMatY(centreMatY == 0) = y;

根据位置将两个列向量排序为 3D 矩阵

Sorting two column vectors into 3D matrix based on position

sorting

grid

matlab