高效地找到标记图像区域的质心

问题描述：

我有一个分割图像作为独特标签1 ... k的2维矩阵。例如：高效地找到标记图像区域的质心

img = 
    [1 1 2 2 2 2 2 3 3] 
    [1 1 1 2 2 2 2 3 3] 
    [1 1 2 2 2 2 3 3 3] 
    [1 4 4 4 2 2 2 2 3] 
    [4 4 4 5 5 5 2 3 3] 
    [4 4 4 5 5 6 6 6 6] 
    [4 4 5 5 5 5 6 6 6]

我想确定区域的质心。也就是说，根据标签，质心的X，Y坐标是什么？例如，标签1的质心是（1.25，0.625）。我只知道行数（(0 + 0 + 1 + 1 + 1 + 2 + 2 + 3)/8 = 1.25）和列数（(0 + 0 + 0 + 0 + 1 + 1 + 1 + 2)/8 = 0.625）

我知道如何做到这一点的唯一方法是使用for循环从1到k（或在我的例子中，1到6），找到每个标签点的索引，并通过索引图像的网格来平均它们的坐标。

但是，我正在寻找以GPU计算优化的方式做到这一点。因此，使用for循环并不理想（在一个漂亮的GPU上，每个图像需要花费大约1秒时间才能完成几百个标签）。我正在使用PyTorch，但真正的任何numpy解决方案都应该足够了。

是否有GPU高效的解决方案来完成此任务？

答

这个计算需要积累，我不知道它在GPU上有多高效。这是伪代码的顺序算法：

int n[k] = 0 
int sx[k] = 0 
int sy[k] = 0 
loop over y: 
    loop over x: 
     i = img[x,y] 
     ++n[i] 
     sx[i] += x 
     sy[i] += y 
for i = 1 to k 
    sx[i] /= n[i] 
    sy[i] /= n[i]

然后当然，(sx[i],sy[i])是对象i的重心。

它在CPU上真的很快，除非它已经存在，否则不值得将数据发送到GPU。

是的这是我可以提出的最佳解决方案。数据已经在GPU上，因此我为什么要这么做： - /。我会看看CPU如何去。这当然是高效的 - 每一个像素都会触及一次，但我想知道是否有办法在GPU上实现并行化... – marcman

有一件事你可以做，如果上述算法的天真GPU实现是太慢了，就是为每个图像列创建数组'n'，'sx'和'sy'，然后将它们添加到一起。这可能会减少等待对数组值进行原子更新的内核数量。 –

我不确定，以原子方式更新全局数组中的单个值是否有效？你在GPU上获得虚假分享问题吗？ –

答

考虑使用scikit-image或重新使用它们的code（基于numpy/scipy）。

这里是一个演示：

import numpy as np 
from skimage import measure 
from time import perf_counter as pc 

img = np.array([[1, 1, 2, 2, 2, 2, 2, 3, 3], 
       [1, 1, 1, 2, 2, 2, 2, 3, 3], 
       [1, 1, 2, 2, 2, 2, 3, 3, 3], 
       [1, 4, 4, 4, 2, 2, 2, 2, 3], 
       [4, 4, 4, 5, 5, 5, 2, 3, 3], 
       [4, 4, 4, 5, 5, 6, 6, 6, 6], 
       [4, 4, 5, 5, 5, 5, 6, 6, 6]]) 

# assuming already labels of 1, 2, ... n 
times = [pc()] 
props = measure.regionprops(img) 
times.append(pc()) 
for i in range(np.unique(img).shape[0]): 
    print(props[i].centroid) 
    times.append(pc()) 

print(np.diff(times))

输出：

(1.25, 0.625) 
(1.5555555555555556, 4.4444444444444446) 
(1.8999999999999999, 7.4000000000000004) 
(4.3636363636363633, 1.1818181818181819) 
(5.1111111111111107, 3.6666666666666665) 
(5.4285714285714288, 6.7142857142857144) 
[ 9.05569615e-05 8.51235438e-04 2.48126075e-04 2.59294767e-04 
    2.42692657e-04 2.00734598e-04 2.34542530e-04]

答

一种想法是使用bincount累积的行索引和列索引为使用数字输入阵列中的每个区域该箱，因此有一个矢量化的解决方案，像这样 -

m,n = a.shape 
r,c = np.mgrid[:m,:n] 
count = np.bincount(a.ravel()) 
centroid_row = np.bincount(a.ravel(),r.ravel())/count 
centroid_col = np.bincount(a.ravel(),c.ravel())/count

Sa mple run -

In [77]: a 
Out[77]: 
array([[1, 1, 2, 2, 2, 2, 2, 3, 3], 
     [1, 1, 1, 2, 2, 2, 2, 3, 3], 
     [1, 1, 2, 2, 2, 2, 3, 3, 3], 
     [1, 4, 4, 4, 2, 2, 2, 2, 3], 
     [4, 4, 4, 5, 5, 5, 2, 3, 3], 
     [4, 4, 4, 5, 5, 6, 6, 6, 6], 
     [4, 4, 5, 5, 5, 5, 6, 6, 6]]) 

In [78]: np.c_[centroid_row, centroid_col] 
Out[78]: 
array([[ nan, nan], 
     [ 1.25, 0.62], # centroid for region-1 
     [ 1.56, 4.44], # centroid for region-2 
     [ 1.9 , 7.4 ], # centroid for region-3 and so on. 
     [ 4.36, 1.18], 
     [ 5.11, 3.67], 
     [ 5.43, 6.71]])

高效地找到标记图像区域的质心

相关推荐