一文弄懂图像的矩和相关应用

InfoRich 2021-11-01

展开全文

引言

在数字图像处理、计算机视觉与相关领域中，图像矩(Image moments)是指图像的某些特定像素灰度的加权平均值（矩），或者是图像具有类似功能或意义的属性。

图像矩通常用来描述分割后的图像对象。可以通过图像的矩来获得图像的部分性质，包括面积(或总体亮度)，以及有关几何中心和方向的信息。

原始矩

对于灰度图像f,其(p+q)阶矩(有时称为'原始矩'),定义如下:

原始矩中包含以下有关原始图像属性的信息:

二值图像或者或度图像的像素总和,可以表示为: $m_{00}$ ,即

图像的几何中心也叫质心 $(c_x,c_y)$ ,计算公式如下:

中心矩

中心矩是以质心为中心的矩,相比原始矩,只要添加一个平移即可,中心矩的定义如下:

上述公式中的 $(c_x,c_y)$ 为图像的质心.

有了中心矩的定义,我们可以得到图像f的协方差矩阵,其定义为:

接着我们可以计算该协方差矩阵的特征向量,最大特征值对应的特征向量就是物体长轴的角度,我们可以将其定义为物体的方向.

该角度计算公式如下:

中心矩的性质

接下来我们推导一下阶数较低的中心矩的性质:

应用

推导了一堆公式,不知道小伙伴有没有觉得很枯燥。

好的吧,接下来我们来举个栗子,进入轻松的写代码实战环节.

本节的目的在于利用矩求出前景物体的中心点以及长轴的角度.

1）读取图像

我们这里以一群飞翔的大雁组成的图像为例来进行讲解,首先我们读入图像并执行灰度化,代码如下:

img_file = './images/bird_swarm.jpg'img = cv2.imread(img_file,cv2.IMREAD_COLOR)gray_img = cv2.cvtColor(img,cv2.COLOR_RGB2GRAY)

效果如下:

上图中,左侧为原始输入图,右侧为灰度化后的图像.

2）二值化

通过第二章节我们知道,二值图像的原始矩就是所有像素值之和,所以我们可以先进行二值化来提取我们的前景物体,二值化代码如下:

retval, bin_img = cv2.threshold(gray_img, 50, 1, cv2.THRESH_BINARY)show_img = bin_img * 255

结果如下:

可以看出,我们使用较低的阈值,将前景大雁(白色部分)通过二值化提取出来.

3）图像矩的计算

接下来我们来计算图像矩,我们得到的二值化图中,对于前景大雁来说 $(x,y)=1$ ,对于背景黑色像素值来说 $f(x,y)=0$ .

根据上文矩的定义,编写代码如下:

m00 = m01 = m10 = m11 = m20 = m02 = m21 = m12 = 0height, width = bin_img.shapefor y in range(height): for x in range(width): m00 += bin_img[y, x] m10 += x * bin_img[y, x] m01 += y * bin_img[y, x] m11 += x * y * bin_img[y, x] m20 += x * x * bin_img[y, x] m02 += y * y * bin_img[y, x] m21 += x * x * y * bin_img[y, x] m12 += x * y * y * bin_img[y, x]

进而利用以下代码可以得到中心点的位置坐标:

cx = m10/m00cy = m01/m00print('Centriod: ({0:.2f}, {1:.2f})'.format(cx, cy))

输出如下:

Centriod: (1039.64, 660.24)

画到图像上,结果如下:

上图中的红色小圆圈即为计算出来的大雁群的中心点.

4）计算角度

根据上述中心矩的定义和物体长轴角度的定义,我们可以计算出大雁群长轴和水平线的夹角,代码如下:

mu00 = m00mu11 = m11 - cx*m01mu20 = m20 - cx*m10mu02 = m02 - cy*m01theta = 1/2*np.arctan2(2*mu11/mu00, (mu20 - mu02)/mu00)print('Angle {0:.2f}'.format(theta*180/np.pi))

输出如下:

Angle 37.96

5）结果可视化

虽然通过上述代码,我们可以计算出大雁群的质心和长轴角度,但是只有数值不太直观,那么我们接下来来可视化上述输出.

我们现在在原始图像上绘制长轴和短轴（这里设置默认长度）以可视化说明质心和角度的位置。

代码如下:

# visualrho = 800dx_major = rho * np.cos(theta)dy_major = rho * np.sin(theta)dx_minor = 0.3 * rho * np.cos(theta - np.pi / 2)dy_minor = 0.3 * rho * np.sin(theta - np.pi / 2)# shortshort_axis=[(int(cx-dx_minor),int(cy-dy_minor)),(int(cx),int(cy)),(int(cx+dx_minor),int(cy+dy_minor))]for i in range(len(short_axis)-1):    cv2.line(img,short_axis[i],short_axis[i+1],color=(255,0,0),thickness=2)for pt in short_axis:    cv2.circle(img,pt,radius=5,color=(255,0,0),thickness=3,lineType=-1)# longlong_axis = [(int(cx - dx_major), int(cy - dy_major)), (int(cx), int(cy)), (int(cx + dx_major), int(cy + dy_major))]for i in range(len(long_axis) - 1):    cv2.line(img, long_axis[i], long_axis[i + 1], color=(0, 0, 255), thickness=2)for pt in long_axis:    cv2.circle(img,pt,radius=5,color=(0,0,255),thickness=3,lineType=-1)# centercv2.circle(img, (int(cx),int(cy)), radius=5, color=(0, 255, 0), thickness=3, lineType=-1)# showcv2.imshow('img',img)cv2.waitKey(0)