2012-08-31 41 views
7

我有一個閉合的非自相交多邊形。它的頂點保存在兩個向量X和Y中。最後,X和Y的值被綁定在0和22之間。什麼是計算圖像與多邊形重疊的簡單方法?

我想構造一個大小爲22x22的矩陣,並設置每個bin的值等於如果部分多邊形與該部分重疊,則爲true,否則爲false。

我最初的想法是生成一個由[a, b] = meshgrid(1:22)定義的點的網格,然後使用inpolygon來確定網格中的哪些點在多邊形中。

[a b] = meshgrid(1:22); 
inPoly1 = inpolygon(a,b,X,Y); 

然而,這僅返回true,如果如果倉的中心包含在多邊形,也就是說,它下面的圖片返回紅色形狀。然而,什麼需要更多地沿着綠色形狀的線條(儘管它仍然是一個不完整的解決方案)。

爲了獲得綠色斑點,我對inpolygon執行了四次調用。對於每次比較,我將NE,NW,SE或SW的網格點移動1/2。這相當於測試垃圾箱的角落是否在多邊形中。

inPoly2 = inpolygon(a-.5,b-.5,X,Y) | inpolygon(a+.5,b-.5,X,Y) | inpolygon(a-.5,b+5,X,Y) | inpolygon(a+.5,b+.5,X,Y); 

雖然這並提供給我的情況下,失敗的部分解決方案時,一個頂點是包含一個箱子,但沒有垃圾桶的角落都是。

有沒有更直接的方法來攻擊這個問題,最好是一個能產生更多可讀代碼的解決方案?

enter image description here

這個陰謀被繪製:

imagesc(inPoly1 + inPoly2); hold on; 
line(a, b, 'w.'); 
line(X, Y, 'y); 
+0

必須離開計算機,但認爲我會提供一個可能有所幫助的通用解決方案。首先將網格網格縮放至22的倍數,以便以等於或大於頂點使用的密度來定義面積 - 這將消除拐角問題。然後再回到22乘22的網格,您可以簡單地按照您放大的相同因子進行分割,將頂部/左側的點鋪在地面上,並將底部/右側的點鋪在天花板上。希望有幫助 – Salain

回答

5

一個建議是使用polybool功能(在2008b或更早版本中不可用)。它找到兩個多邊形的相交點並返回結果頂點(如果沒有頂點存在,則返回空向量)。要在這裏使用它,我們迭代(使用arrayfun)在網格檢查中的所有方塊中,查看polybool的輸出參數是否爲空(例如,不重疊)。

N=22; 
sqX = repmat([1:N]',1,N); 
sqX = sqX(:); 
sqY = repmat(1:N,N,1); 
sqY = sqY(:); 

intersects = arrayfun((@(xs,ys) ... 
     (~isempty(polybool('intersection',X,Y,[xs-1 xs-1 xs xs],[ys-1 ys ys ys-1])))),... 
     sqX,sqY); 

intersects = reshape(intersects,22,22); 

這裏是產生的圖像:

enter image description here

代碼繪製:

imagesc(.5:1:N-.5,.5:1:N-.5,intersects'); 
hold on; 
plot(X,Y,'w'); 
for x = 1:N 
    plot([0 N],[x x],'-k'); 
    plot([x x],[0 N],'-k'); 
end 
hold off; 
+0

不是我在找什麼,但它的作品! – slayton

1

這個怎麼樣僞代碼算法:

For each pair of points p1=p(i), p2=p(i+1), i = 1..n-1 
    Find the line passing through p1 and p2 
    Find every tile this line intersects // See note 
    Add intersecting tiles to the list of contained tiles 

Find the red area using the centers of each tile, and add these to the list of contained tiles 

注意:此線將採取的努力,實現一點點,但我認爲它有一個相當直接的,衆所周知的算法。另外,如果我使用.NET,我會簡單地定義一個對應於每個網格拼塊的矩形,然後查看哪些拼塊與多邊形相交。但是,我不知道在Matlab中檢查交叉點是否容易。

0

首先我定義了一個低分辨率圓這個例子

X=11+cos(linspace(0,2*pi,10))*5; 
Y=11+sin(linspace(0,2.01*pi,10))*5; 

就像你的例子是,在中〜22個單元網格適合。然後,跟隨你的主角,我們聲明一個meshgrid並檢查點是否在多邊形中。

stepSize=0.1; 
[a b] = meshgrid(1:stepSize:22); 
inPoly1 = inpolygon(a,b,X,Y); 

唯一的區別是,如果您的原始解決方案採取的步驟之一,此網格可以採取較小的步驟。最後,以包括正方形

inPolyFull=unique(round([a(inPoly1) b(inPoly1)]) ,'rows'); 

round只是需要我們的高分辨率電網和舍入點適當地離他們最近的分辨率低當量「邊緣」中的任何東西。然後我們使用'rows'限定符呼叫unique來刪除矢量樣式或成對樣式中的所有副本。這就是它

要查看結果,

[aOrig bOrig] = meshgrid(1:22); 
imagesc(1:stepSize:22,1:stepSize:22,inPoly1); hold on; 
plot(X,Y,'y'); 
plot(aOrig,bOrig,'k.'); 
plot(inPolyFull(:,1),inPolyFull(:,2),'w.'); hold off; 

Example polygon

更改stepSize具有改善的速度和內存的成本結果的預期效果。

如果您需要的結果是在相同的格式,在您的示例inPoly2,您可以使用

inPoly2=zeros(22); 
inPoly2(inPolyFull(:,1),inPolyFull(:,2))=1 

希望有所幫助。我可以想到其他一些方法,但這似乎是最直接的。

1

我建議在圖像處理工具箱使用poly2mask,它或多或少什麼你想要的,我想,也或多或少都是你自己和Salain建議的。

+1

是的,我已經嘗試過'poly2mask'。它給出了與我在我的問題中描述的「inpolygon」調用相同的結果。 – slayton

1

略有改善

首先,簡化你的「部分解決」 - 你在做什麼,只是看着在角落。如果不考慮22x22點的網格,您可以考慮23x23網格的角點(它將從較小的網格偏移(-0.5,-0.5)。一旦你有了,你可以在22x22網格上標記點有至少一個角在多邊形

完整的解決方案:。

但是,你真正尋找的是多邊形是否與周圍的每一個像素的1x1盒相交,這並不一定包括任何角落,但它確實要求多邊形與盒子的四個邊中的一個相交

一個你能找到的像素而多邊形與收納箱相交的方法是使用下面的算法:

For each pair of adjacent points in the polygon, calling them pA and pB: 
    Calculate rounded Y-values: Round(pA.y) and Round(pB.y) 
    For each horizontal pixel edge between these two values: 
     * Solve the simple linear equation to find out at what X-coordinate 
      the line between pA and pB crosses this edge 
     * Round the X-coordinate 
     * Use the rounded X-coordinate to mark the pixels above and below 
      where it crosses the edge 
    Do a similar thing for the other axis 

因此,舉例來說,假設我們正在尋找pA = (1, 1)pB = (2, 3)

  • 首先,我們計算了四捨五入的Y值:1和3。
  • 然後,我們看一下這些值之間的像素邊緣:y = 1.5y = 2.5(像素邊緣是從像素的半偏移
  • 對於這些,我們求解線性方程以找到pA - >pB與相交我們的。邊緣這使我們:x = 1.25, y = 1.5,和x = 1.75, y = 2.5
  • 對於每個這些交叉點,我們採取圓形X值,並用它來標記的邊緣的像素兩側
    • x = 1.25被四捨五入爲1。 (邊緣y = 1.5)。因此,我們可以標記像素爲(1, 1)(1, 2)作爲我們的一部分。
    • x = 1.75被舍入爲2(邊緣y = 2.5)。因此,我們可以將像素標記爲(2, 2)(2, 3)

所以這是的照顧水平邊緣。接下來,讓我們看看那些垂直:

  • 首先我們計算圓形的X值:1和2
  • 然後,我們來看看像素邊緣。這裏只有一個:x = 1.5
  • 對於這個邊緣,我們找到它符合線條的地方pA - >pB。這給了我們x = 1.5, y = 2
  • 對於此交匯,我們採取圓形Y值,並用它來標記像素任一邊緣的側:
    • y = 2舍入成2因此,我們可以在(1, 2)(2, 2)標記像素。

完成!

那麼,有點。這會給你的邊緣,但它不會填充多邊形的主體。但是,您可以將這些與您以前的(紅色)結果結合起來以獲得完整的集合。

0

嗯,我想我遲到了,儘管嚴格來說,賞金時間一直持續到明天;)。但是,這是我的嘗試。首先,標記包含/觸摸點的單元格的函數。考慮到與間隔LX,LY,和具有座標(XP,YP)的一組點的結構化格柵,設置含有細胞:

function cells = mark_cells(lx, ly, Xp, Yp, cells) 

% Find cell numbers to which points belong. 
% Search by subtracting point coordinates from 
% grid coordinates and observing the sign of the result. 
% Points lying on edges/grid points are assumed 
% to belong to all surrounding cells. 

sx=sign(bsxfun(@minus, lx, Xp')); 
sy=sign(bsxfun(@minus, ly, Yp')); 
cx=diff(sx, 1, 2); 
cy=diff(sy, 1, 2); 

% for every point, mark the surrounding cells 
for i=1:size(cy, 1) 
    cells(find(cx(i,:)), find(cy(i,:)))=1; 
end 
end 

現在,代碼的其餘部分。對於多邊形中的每個分段(您必須逐個穿過分段),將分段與網格線相交。使用給定的網格點座標分別進行水平線和垂直線的交叉操作,以避免數值不準確。對於發現交點我打電話mark_cells到周圍的細胞標記爲1:

% example grid 
nx=21; 
ny=51; 
lx = linspace(0, 1, nx); 
ly = linspace(0, 1, ny); 
dx=1/(nx-1); 
dy=1/(ny-1); 
cells = zeros(nx-1, ny-1); 

% for every line in the polygon... 
% Xp and Yp contain start-end points of a single segment 
Xp = [0.15 0.61]; 
Yp = [0.1 0.78]; 

% line equation 
slope = diff(Yp)/diff(Xp); 
inter = Yp(1) - (slope*Xp(1)); 

if isinf(slope) 
    % SPECIAL CASE: vertical polygon segments 
    % intersect horizontal grid lines 
    ymax = 1+floor(max(Yp)/dy); 
    ymin = 1+ceil(min(Yp)/dy); 
    x=repmat(Xp(1), 1, ymax-ymin+1); 
    y=ly(ymin:ymax); 
    cells = mark_cells(lx, ly, x, y, cells); 
else 
    % SPECIAL CASE: not horizontal polygon segments 
    if slope ~= 0 
     % intersect horizontal grid lines 
     ymax = 1+floor(max(Yp)/dy); 
     ymin = 1+ceil(min(Yp)/dy); 
     xmax = (ly(ymax)-inter)/slope; 
     xmin = (ly(ymin)-inter)/slope; 
     % interpolate in x... 
     x=linspace(xmin, xmax, ymax-ymin+1); 
     % use exact grid point y-coordinates! 
     y=ly(ymin:ymax); 
     cells = mark_cells(lx, ly, x, y, cells); 
    end 

    % intersect vertical grid lines 
    xmax = 1+floor(max(Xp)/dx); 
    xmin = 1+ceil(min(Xp)/dx); 
    % interpolate in y... 
    ymax = inter+slope*lx(xmax); 
    ymin = inter+slope*lx(xmin); 
    % use exact grid point x-coordinates! 
    x=lx(xmin:xmax); 
    y=linspace(ymin, ymax, xmax-xmin+1); 
    cells = mark_cells(lx, ly, x, y, cells); 
end 

輸出爲例如網孔/區段: output

通過所有多邊形段步行給你多邊形「暈輪」。最後,使用標準的inpolygon函數獲得多邊形的內部。讓我知道你是否需要更多關於代碼的細節。

相關問題