令F表示标记图像,令G表示模板图像。在讨论中,我们假设两幅图像都是二值图像,且F包含于G。标记图像相对于模板大小为1的测地膨胀定义为
F相对于G的大小为n的测地膨胀定义为:
式中,n≥1是整数,D0=F。在这个递归公式中,每一步都是膨胀后取交集。交集运算可以保证模板G限制标记F的生长(膨胀)。下图是大小为1的测地膨胀。
- import cv2
- import numpy as np
- import matplotlib.pyplot as plt
-
- origin = cv2.imread("a.jpg")
- imgGray = cv2.imread("a.jpg", flags=0) # flags=0 灰度图像
- ret, imgBin = cv2.threshold(imgGray, 205, 255,
- cv2.THRESH_BINARY) # 二值化处理 (黑色0/白色1) ret:阈值 205:阈值 255:最大值 cv2.THRESH_BINARY:二值化
- imgMask = cv2.bitwise_not(imgBin) # 二值图像的补集 (黑色1/白色0)
-
- # 构造标记图像
- marker = np.zeros_like(imgBin, dtype=np.uint8) # 与imgBin同尺寸的全0图像
- marker[0, :] = imgBin[0, :] # 第一行
- marker[-1, :] = imgBin[-1, :] # 最后一行
- marker[:, 0] = imgBin[:, 0] # 第一列
- marker[:, -1] = imgBin[:, -1] # 最后一列
- markerIni = marker.copy() # 标记图像: 边框 f(x,y)=I(x,y),其它为 0
-
- # 形态学重建
- element = cv2.getStructuringElement(cv2.MORPH_CROSS, (3, 3))
- # 创建数组画出迭代图
- diffList = []
- while True: # 进行测地膨胀 是有条件的膨胀
- marker_pre = marker # 保存 F(n-1)
- dilation = cv2.dilate(marker, kernel=element) # 膨胀重建 目的是扩大 F(n-1) 的边界
- marker = cv2.bitwise_and(dilation, imgMask) # 原图像作为模板用来约束重建,按位与,有 0 得 0 目的是保持 F(n-1) 的内部
- # 计算重建后的标记图像与重建前的标记图像的差异
- diff = cv2.subtract(marker, marker_pre)
- print(diff.sum())
- diffList.append(diff.sum())
- if (marker_pre == marker).all(): # F(n)=F(n-1)?,判断是否达到稳定收敛状态
- break
- imgRebuild = cv2.bitwise_or(imgBin, marker) # 按位或,有 1 得 1
- # 绘制迭代图,命名为"图1:多次迭代的图像差异"
- plt.figure("图1:多次迭代的图像差异")
- plt.plot(diffList)
- plt.show()
-
- # 显示
- # cv2.imshow("imgGray", imgGray)
- # cv2.imshow("imgBin", imgBin)
- cv2.imshow("figure2:after Expansion", imgMask)
- # cv2.imshow("markerIni", markerIni)
- # cv2.imshow("imgRebuild", imgRebuild)
- # cv2.waitKey(0)
- # cv2.destroyAllWindows()
-
-
- # 统计重建后的图像中的连通域个数
- a = 1 # 腐蚀膨胀次数
- img = imgMask
- gray = img
- kernel = np.ones((2, 2), np.uint8) # 进行腐蚀膨胀操作 (2,2)的矩阵 1代表白色 0代表黑色 uint8代表8位无符号整数
- erosion = cv2.erode(gray, kernel, iterations=a) # 膨胀 iterations代表迭代次数
- dilation = cv2.dilate(erosion, kernel, iterations=a) # 腐蚀 先腐蚀再膨胀的目的是去除噪声
- ret, thresh = cv2.threshold(dilation, 150, 255, cv2.THRESH_BINARY) # 阈值处理 二值法 150代表阈值 255代表最大值 cv2.THRESH_BINARY代表二值法
- thresh1 = cv2.GaussianBlur(thresh, (15, 15), 0) # 高斯滤波 3,3代表高斯核的大小 0代表标准差
- contours, hirearchy = cv2.findContours(thresh1, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE) # 找出连通域
- # 对连通域面积进行比较
- area = [] # 建立空数组,放连通域面积
- contours1 = [] # 建立空数组,放减去最小面积的数
-
- for i in contours:
- # area.append(cv2.contourArea(i))
- # print(area)
- if cv2.contourArea(i) > 10: # 计算面积 去除面积小的 连通域
- contours1.append(i)
- # 输出细胞数目
- print(len(contours1))
- draw = cv2.drawContours(origin, contours1, -1, (0, 255, 0), 1) # 描绘连通域
- # 求连通域重心 以及 在重心坐标点描绘数字
- for i, j in zip(contours1, range(len(contours1))): # zip函数将两个数组合并
- M = cv2.moments(i) # 计算矩
- cX = int(M["m10"] / M["m00"]) # 计算重心横坐标
- cY = int(M["m01"] / M["m00"]) # 计算重心纵坐标
- draw1 = cv2.putText(draw, str(j), (cX, cY), 1, 1, (255, 0, 255), 1) # 在中心坐标点上描绘数字
- # 展示图片
- cv2.imshow("figure3:img_after_operation", thresh1)
- cv2.imshow("figure4:result", draw1)
- cv2.waitKey()
测地膨胀和腐蚀的概念:【图像处理笔记】形态学重建 - 湾仔码农 - 博客园
例程 10.24:基于形态学重建的细胞计数 第二部分的代码来源于此:【youcans 的 OpenCV 例程200篇】134. 形态学重建之细胞计数_YouCans的博客-CSDN博客
第二第三步的代码来源于此
opencv——python细胞计数实现_microspore的博客-CSDN博客_opencv细胞计数
如果不进行测地膨胀直接进行第二第三步会怎样
效果不理想的代码:
- import cv2
- import numpy as np
- a=13
- img=cv2.imread(r'C:\Users\Administrator\PycharmProjects\pythonProject2\a.jpg',1) #读取图片
- #去除斑点
-
- gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY) #将图片变为灰度图片
- kernel=np.ones((2,2),np.uint8) #进行腐蚀膨胀操作 (2,2)的矩阵 1代表白色 0代表黑色 uint8代表8位无符号整数
- erosion=cv2.erode(gray,kernel,iterations=a) #膨胀 iterations代表迭代次数
- dilation=cv2.dilate(erosion,kernel,iterations=a) #腐蚀
- ret, thresh = cv2.threshold(dilation, 150, 255, cv2.THRESH_BINARY) # 阈值处理 二值法 150代表阈值 255代表最大值 cv2.THRESH_BINARY代表二值法
- thresh1 = cv2.GaussianBlur(thresh,(3,3),0)# 高斯滤波 3,3代表高斯核的大小 0代表标准差
- contours,hirearchy=cv2.findContours(thresh1, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)# 找出连通域
- #对连通域面积进行比较
- area=[] #建立空数组,放连通域面积
- contours1=[] #建立空数组,放减去最小面积的数
- count=0 #计数
- total_area=0 #总面积
- for i in contours:
- total_area+=cv2.contourArea(i) #计算连通域面积
- avg_area=total_area/len(contours) #计算平均面积
- for i in contours:
- # area.append(cv2.contourArea(i))
- # print(area)
- if cv2.contourArea(i)>10: # 计算面积 去除面积小的 连通域
- contours1.append(i)
- count+=cv2.contourArea(i)/avg_area #计算面积比
- print(int(count)) #计算连通域个数
- draw=cv2.drawContours(img,contours1,-1,(0,255,0),1) #描绘连通域
- #求连通域重心 以及 在重心坐标点描绘数字
- for i,j in zip(contours1,range(len(contours1))): #zip函数将两个数组合并
- M = cv2.moments(i) # 计算矩
- cX=int(M["m10"]/M["m00"]) #计算重心横坐标
- cY=int(M["m01"]/M["m00"]) #计算重心纵坐标
- draw1=cv2.putText(draw, str(j), (cX, cY), 1,1, (255, 0, 255), 1) #在中心坐标点上描绘数字
- #展示图片
- cv2.imshow("draw",draw1)
- cv2.imshow("thresh1",thresh1)
- cv2.waitKey()
- cv2.destroyWindow()