在 Matlab 中的 pcolor 中查找 contour/edge
Find contour/edge in pcolor in Matlab
我正在尝试在 Matlab 中的 pcolor 图中绘制一个轮廓，该轮廓遵循 'pixels' 的边缘。这可能最好用图片来解释。这是我的数据图。黄色数据（data==1）和蓝色数据（data==0）之间有明显的界限：
请注意，这是一个 pcolor 图，因此每个 'square' 本质上都是一个像素。我想要 return 一个跟随黄色数据 faces 的轮廓 pixels, not 只是黄色数据的边缘。
因此输出轮廓（绿线）通过像素的面（红点）的中点。
请注意，我不希望轮廓跟随数据的中心点（黑点），这会像这条绿线那样。这可以通过 contour.
轻松实现
此外，如果有任何帮助，我有一些可能有用的网格。我在像素中间有点（很明显，这就是我在这里绘制的），我在角落也有点，我在 west/east 面和 north/south 面孔。如果您熟悉 Arakawa grids，这是一个 Arakawa-C 网格，所以我有 rho-、u-、v- 和 psi- 点。
我试过插值、交织网格和其他一些东西，但我没有任何运气。任何帮助将不胜感激，并会阻止我发疯。
干杯，戴夫
编辑：
抱歉，我简化了图像以使我试图解释的内容更加明显，但这是我试图分离的区域的更大（缩小）图像：
如您所见，它是一个复杂的轮廓，在环绕并返回 "northeast" 之前朝 "southwest" 方向前进。这是我想通过黑点画的红线：
看看下面的代码：
% plotting some data:
data = [0 0 0 0 0 0 1 1
0 0 0 0 0 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 0 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 1 1 1 1];
p = pcolor(data);
axis ij
% compute the contour
x = size(data,2)-cumsum(data,2)+1;
x = x(:,end);
y = (1:size(data,1));
% compute the edges shift
Y = get(gca,'YTick');
y_shift = (Y(2)-Y(1))/2;
% plot it:
hold on
plot(x,y+y_shift,'g','LineWidth',3,'Marker','o',...
'MarkerFaceColor','r','MarkerEdgeColor','none')
它产生这个：
这是你要找的吗？
上面最重要的几行是：
x = size(data,2)-cumsum(data,2)+1;
x = x(:,end);
它为每一行找到在 0 到 1 之间移动的位置（假设一行中只有一个）。
然后，在 plot 内，我将 y 移动两个相邻 y 轴刻度之间距离的一半，因此它们将被放置在边缘的中心。
编辑：
经过对这类数据的一些试验，我得到了这样的结果：
imagesc(data);
axis ij
b = bwboundaries(data.','noholes');
x = b{1}(:,1);
y = b{1}(:,2);
X = reshape(bsxfun(@plus,x,[0 -0.5 0.5]),[],1);
Y = reshape(bsxfun(@plus,y,[0 0.5 -0.5]),[],1);
k = boundary(X,Y,1);
hold on
plot(X(k),Y(k),'g','LineWidth',3,'Marker','o',...
'MarkerFaceColor','r','MarkerEdgeColor','none')
它并不完美，但可以通过更简单的方法让您更接近您想要的：
好的，我想我已经解决了...很接近了，很高兴。
首先我获取原始数据（我称之为 mask_rho 并使用它来制作掩码 mask_u、mask_v，这与 mask_rho 类似但被移动了分别在水平和垂直方向略微。
%make mask_u and mask_v
for i = 2:size(mask_rho,2)
for j = 1:size(mask_rho,1)
mask_u(j, i-1) = mask_rho(j, i) * mask_rho(j, i-1);
end
end
for i = 1:size(mask_rho,2)
for j = 2:size(mask_rho,1)
mask_v(j-1, i) = mask_rho(j, i) * mask_rho(j-1, i);
end
end
然后我制作修改后的蒙版 mask_u1 和 mask_v1，它们与 mask_rho 相同，但分别在水平和垂直方向上与相邻点进行平均。
%make mask which is shifted E/W (u) and N/S (v)
mask_u1 = (mask_rho(1:end-1,:)+mask_rho(2:end,:))/2;
mask_v1 = (mask_rho(:,1:end-1)+mask_rho(:,2:end))/2;
然后我利用掩码之间的差异来定位掩码在水平方向（在u掩码中）和垂直方向（在v掩码中）从0到1和1到0变化的地方。
% mask_u-mask_u1 gives the NEXT row with a change from 0-1.
diff_mask_u=logical(mask_u-mask_u1);
lon_u_bnds=lon_u.*double(diff_mask_u);
lon_u_bnds(lon_u_bnds==0)=NaN;
lat_u_bnds=lat_u.*double(diff_mask_u);
lat_u_bnds(lat_u_bnds==0)=NaN;
lon_u_bnds(isnan(lon_u_bnds))=[];
lat_u_bnds(isnan(lat_u_bnds))=[];
%now same for changes in mask_v
diff_mask_v=logical(mask_v-mask_v1);
lon_v_bnds=lon_v.*double(diff_mask_v);
lon_v_bnds(lon_v_bnds==0)=NaN;
lat_v_bnds=lat_v.*double(diff_mask_v);
lat_v_bnds(lat_v_bnds==0)=NaN;
lon_v_bnds(isnan(lon_v_bnds))=[];
lat_v_bnds(isnan(lat_v_bnds))=[];
bnd_coords_cat = [lon_u_bnds,lon_v_bnds;lat_u_bnds,lat_v_bnds]'; %make into 2 cols, many rows
结果抓取边界边缘的所有坐标：
现在我的回答有点不对劲。如果我将上面的向量绘制为点 plot(bnd_coords_cat(:,1),bnd_coords_cat(:,2),'kx'，我会得到上面的图像，这很好。但是，如果我加入该行，如： plot(bnd_coords_cat(:,1),bnd_coords_cat(:,2),'-' 那么该行会跳来跳去，因为这些点没有排序。当我进行排序（使用 sort 和 pdist2）以按最接近的点排序时，Matlab 有时会选择奇数点......不过我想我会将这段代码作为附录和可选的额外内容包含在内。有人可能知道更好的排序输出向量的方法bnds_coords_cat:
% now attempt to sort
[~,I]=sort([lon_u_bnds,lon_v_bnds]);
bnd_coords_inc1 = bnd_coords_cat(I,1);
bnd_coords_inc2 = bnd_coords_cat(I,2);
bnd_coords = [bnd_coords_inc1,bnd_coords_inc2];
bnd_coords_dist = pdist2(bnd_coords,bnd_coords);
bnd_coords_sort = nan(1,size(bnd_coords,1));
bnd_coords_sort(1)=1;
for ii=2:size(bnd_coords,1)
bnd_coords_dist(:,bnd_coords_sort(ii-1)) = Inf; %don't go backwards?
[~,closest_idx] = min(bnd_coords_dist(bnd_coords_sort(ii-1),:));
bnd_coords_sort(ii)=closest_idx;
end
bnd_coords_final(:,1)=bnd_coords(bnd_coords_sort,1);
bnd_coords_final(:,2)=bnd_coords(bnd_coords_sort,2);
请注意，pdist2 方法是由同事建议的，也是来自这个 SO 答案 Sort coordinates points in matlab。这是最终结果：
老实说，没有台词的剧情还是可以的。因此，就我而言，这已经足够接近可以回答了！
您可以通过对 to a . I used a section of the sample image mask in the question for data. First, you will need to fill the holes in the mask, which you can do using imfill from the the Image Processing Toolbox:
进行一些修改来解决此问题
x = 1:15; % X coordinates for pixels
y = 1:17; % Y coordinates for pixels
mask = imfill(data, 'holes');
接下来，应用我的另一个答案中的方法来计算一组有序的轮廓坐标（位于像素角上）：
% Create raw triangulation data:
[cx, cy] = meshgrid(x, y);
xTri = bsxfun(@plus, [0; 1; 1; 0], cx(mask).');
yTri = bsxfun(@plus, [0; 0; 1; 1], cy(mask).');
V = [xTri(:) yTri(:)];
F = reshape(bsxfun(@plus, [1; 2; 3; 1; 3; 4], 0:4:(4*nnz(mask)-4)), 3, []).';
% Trim triangulation data:
[V, ~, Vindex] = unique(V, 'rows');
V = V-0.5;
F = Vindex(F);
% Create triangulation and find free edge coordinates:
TR = triangulation(F, V);
freeEdges = freeBoundary(TR).';
xOutline = V(freeEdges(1, [1:end 1]), 1); % Ordered edge x coordinates
yOutline = V(freeEdges(1, [1:end 1]), 2); % Ordered edge y coordinates
最后，您可以像这样在像素边缘的中心获得所需的坐标：
ex = xOutline(1:(end-1))+diff(xOutline)./2;
ey = yOutline(1:(end-1))+diff(yOutline)./2;
这是显示结果的图表：
imagesc(x, y, data);
axis equal
set(gca, 'XLim', [0.5 0.5+size(mask, 2)], 'YLim', [0.5 0.5+size(mask, 1)]);
hold on;
plot(ex([1:end 1]), ey([1:end 1]), 'r', 'LineWidth', 2);
plot(ex, ey, 'k.', 'LineWidth', 2);
我正在尝试在 Matlab 中的 pcolor 图中绘制一个轮廓，该轮廓遵循 'pixels' 的边缘。这可能最好用图片来解释。这是我的数据图。黄色数据（data==1）和蓝色数据（data==0）之间有明显的界限：
请注意，这是一个 pcolor 图，因此每个 'square' 本质上都是一个像素。我想要 return 一个跟随黄色数据 faces 的轮廓 pixels, not 只是黄色数据的边缘。
因此输出轮廓（绿线）通过像素的面（红点）的中点。
请注意，我不希望轮廓跟随数据的中心点（黑点），这会像这条绿线那样。这可以通过 contour.
此外，如果有任何帮助，我有一些可能有用的网格。我在像素中间有点（很明显，这就是我在这里绘制的），我在角落也有点，我在 west/east 面和 north/south 面孔。如果您熟悉 Arakawa grids，这是一个 Arakawa-C 网格，所以我有 rho-、u-、v- 和 psi- 点。
我试过插值、交织网格和其他一些东西，但我没有任何运气。任何帮助将不胜感激，并会阻止我发疯。
干杯，戴夫
编辑：
抱歉，我简化了图像以使我试图解释的内容更加明显，但这是我试图分离的区域的更大（缩小）图像：
如您所见，它是一个复杂的轮廓，在环绕并返回 "northeast" 之前朝 "southwest" 方向前进。这是我想通过黑点画的红线：
看看下面的代码：
% plotting some data:
data = [0 0 0 0 0 0 1 1
0 0 0 0 0 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 0 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 1 1 1 1
0 0 0 0 1 1 1 1];
p = pcolor(data);
axis ij
% compute the contour
x = size(data,2)-cumsum(data,2)+1;
x = x(:,end);
y = (1:size(data,1));
% compute the edges shift
Y = get(gca,'YTick');
y_shift = (Y(2)-Y(1))/2;
% plot it:
hold on
plot(x,y+y_shift,'g','LineWidth',3,'Marker','o',...
'MarkerFaceColor','r','MarkerEdgeColor','none')
它产生这个：
这是你要找的吗？
上面最重要的几行是：
x = size(data,2)-cumsum(data,2)+1;
x = x(:,end);
它为每一行找到在 0 到 1 之间移动的位置（假设一行中只有一个）。
然后，在 plot 内，我将 y 移动两个相邻 y 轴刻度之间距离的一半，因此它们将被放置在边缘的中心。
编辑：
经过对这类数据的一些试验，我得到了这样的结果：
imagesc(data);
axis ij
b = bwboundaries(data.','noholes');
x = b{1}(:,1);
y = b{1}(:,2);
X = reshape(bsxfun(@plus,x,[0 -0.5 0.5]),[],1);
Y = reshape(bsxfun(@plus,y,[0 0.5 -0.5]),[],1);
k = boundary(X,Y,1);
hold on
plot(X(k),Y(k),'g','LineWidth',3,'Marker','o',...
'MarkerFaceColor','r','MarkerEdgeColor','none')
它并不完美，但可以通过更简单的方法让您更接近您想要的：
好的，我想我已经解决了...很接近了，很高兴。
首先我获取原始数据（我称之为 mask_rho 并使用它来制作掩码 mask_u、mask_v，这与 mask_rho 类似但被移动了分别在水平和垂直方向略微。
%make mask_u and mask_v
for i = 2:size(mask_rho,2)
for j = 1:size(mask_rho,1)
mask_u(j, i-1) = mask_rho(j, i) * mask_rho(j, i-1);
end
end
for i = 1:size(mask_rho,2)
for j = 2:size(mask_rho,1)
mask_v(j-1, i) = mask_rho(j, i) * mask_rho(j-1, i);
end
end
然后我制作修改后的蒙版 mask_u1 和 mask_v1，它们与 mask_rho 相同，但分别在水平和垂直方向上与相邻点进行平均。
%make mask which is shifted E/W (u) and N/S (v)
mask_u1 = (mask_rho(1:end-1,:)+mask_rho(2:end,:))/2;
mask_v1 = (mask_rho(:,1:end-1)+mask_rho(:,2:end))/2;
然后我利用掩码之间的差异来定位掩码在水平方向（在u掩码中）和垂直方向（在v掩码中）从0到1和1到0变化的地方。
% mask_u-mask_u1 gives the NEXT row with a change from 0-1.
diff_mask_u=logical(mask_u-mask_u1);
lon_u_bnds=lon_u.*double(diff_mask_u);
lon_u_bnds(lon_u_bnds==0)=NaN;
lat_u_bnds=lat_u.*double(diff_mask_u);
lat_u_bnds(lat_u_bnds==0)=NaN;
lon_u_bnds(isnan(lon_u_bnds))=[];
lat_u_bnds(isnan(lat_u_bnds))=[];
%now same for changes in mask_v
diff_mask_v=logical(mask_v-mask_v1);
lon_v_bnds=lon_v.*double(diff_mask_v);
lon_v_bnds(lon_v_bnds==0)=NaN;
lat_v_bnds=lat_v.*double(diff_mask_v);
lat_v_bnds(lat_v_bnds==0)=NaN;
lon_v_bnds(isnan(lon_v_bnds))=[];
lat_v_bnds(isnan(lat_v_bnds))=[];
bnd_coords_cat = [lon_u_bnds,lon_v_bnds;lat_u_bnds,lat_v_bnds]'; %make into 2 cols, many rows
结果抓取边界边缘的所有坐标：
现在我的回答有点不对劲。如果我将上面的向量绘制为点 plot(bnd_coords_cat(:,1),bnd_coords_cat(:,2),'kx'，我会得到上面的图像，这很好。但是，如果我加入该行，如： plot(bnd_coords_cat(:,1),bnd_coords_cat(:,2),'-' 那么该行会跳来跳去，因为这些点没有排序。当我进行排序（使用 sort 和 pdist2）以按最接近的点排序时，Matlab 有时会选择奇数点......不过我想我会将这段代码作为附录和可选的额外内容包含在内。有人可能知道更好的排序输出向量的方法bnds_coords_cat:
% now attempt to sort
[~,I]=sort([lon_u_bnds,lon_v_bnds]);
bnd_coords_inc1 = bnd_coords_cat(I,1);
bnd_coords_inc2 = bnd_coords_cat(I,2);
bnd_coords = [bnd_coords_inc1,bnd_coords_inc2];
bnd_coords_dist = pdist2(bnd_coords,bnd_coords);
bnd_coords_sort = nan(1,size(bnd_coords,1));
bnd_coords_sort(1)=1;
for ii=2:size(bnd_coords,1)
bnd_coords_dist(:,bnd_coords_sort(ii-1)) = Inf; %don't go backwards?
[~,closest_idx] = min(bnd_coords_dist(bnd_coords_sort(ii-1),:));
bnd_coords_sort(ii)=closest_idx;
end
bnd_coords_final(:,1)=bnd_coords(bnd_coords_sort,1);
bnd_coords_final(:,2)=bnd_coords(bnd_coords_sort,2);
请注意，pdist2 方法是由同事建议的，也是来自这个 SO 答案 Sort coordinates points in matlab。这是最终结果：
老实说，没有台词的剧情还是可以的。因此，就我而言，这已经足够接近可以回答了！
您可以通过对 data. First, you will need to fill the holes in the mask, which you can do using imfill from the the Image Processing Toolbox:
x = 1:15; % X coordinates for pixels
y = 1:17; % Y coordinates for pixels
mask = imfill(data, 'holes');
接下来，应用我的另一个答案中的方法来计算一组有序的轮廓坐标（位于像素角上）：
% Create raw triangulation data:
[cx, cy] = meshgrid(x, y);
xTri = bsxfun(@plus, [0; 1; 1; 0], cx(mask).');
yTri = bsxfun(@plus, [0; 0; 1; 1], cy(mask).');
V = [xTri(:) yTri(:)];
F = reshape(bsxfun(@plus, [1; 2; 3; 1; 3; 4], 0:4:(4*nnz(mask)-4)), 3, []).';
% Trim triangulation data:
[V, ~, Vindex] = unique(V, 'rows');
V = V-0.5;
F = Vindex(F);
% Create triangulation and find free edge coordinates:
TR = triangulation(F, V);
freeEdges = freeBoundary(TR).';
xOutline = V(freeEdges(1, [1:end 1]), 1); % Ordered edge x coordinates
yOutline = V(freeEdges(1, [1:end 1]), 2); % Ordered edge y coordinates
最后，您可以像这样在像素边缘的中心获得所需的坐标：
ex = xOutline(1:(end-1))+diff(xOutline)./2;
ey = yOutline(1:(end-1))+diff(yOutline)./2;
这是显示结果的图表：
imagesc(x, y, data);
axis equal
set(gca, 'XLim', [0.5 0.5+size(mask, 2)], 'YLim', [0.5 0.5+size(mask, 1)]);
hold on;
plot(ex([1:end 1]), ey([1:end 1]), 'r', 'LineWidth', 2);
plot(ex, ey, 'k.', 'LineWidth', 2);