数字图像处理个人练习06--图像分割

实验环境配置

本次实验使用python语言，版本为3.9.1
使用如下库

import numpy as np# 图形处理
import cv2# 读入读出和一些特殊函数使用
import matplotlib.pyplot as plt#绘画灰度直方图

1. 实现2种基于阈值的分割算法

1. 确定一个合适的阈值T
2. 将大于等于阈值的像素作为物体或背景，生成一个二值图像
3. 在四邻域中有背景的像素，既是边界像素。

阈值处理操作
$$T=T[x,y,p(x,y),f(x,y)]$$
这里f是对应点的灰度，p表示局部性质：中心为邻域的平均灰度级

阈值处理后
$$g(x,y)\{\begin{array}{lc}1& f(x,y)>T\\ 0& f(x,y)\le T\end{array}$$

4. 标记为1的像素对应于前景，标记为0的像素对应于背景
5. 当T仅取决于f(x,y)，阈值称为全局的
6. 当T取决于f(x,y)和p(x,y)，阈值是局部的
7. 当T取决于空间坐标x和y，阈值就是动态的或自适应的

适用于物体与背景有较强对比的情
况，重要的是背景或物体的灰度比
较单一。（可通过先求背景，然后
求反得到物体）

这种方法总可以得到封闭且连通区
域的边界。

基本全局阈值分割以及局部阈值分割

基本思路

寻找灰度直方图的第一个最高点和第二个最高点，然后根据这两个灰度的值来找到谷底，对此本次实验中使用了以下四种思路来进行操作

1. 直接找最小值
2. 对每个点进行局部处理（平均灰度）再选择最小值（这里为局部阈值）
3. 对找到的两个最高点去除边缘再对剩余的点平均灰度并选择最小值
4. 通过固定一个点分开两个区域，计算两个区域的平均灰度是否近似，一直迭代直到近似为止
5. 求出两个峰值后直接取中点作为threshold

核心函数

class segmentation():
    def basic_glo_threshold(self, image):
#这里在类下进行操作，只需传入image的灰度图像即可

后续将通过这四个函数

来获取“基本思路”中四个不同阈值的threshold值
而第五个方法直接用公式$mean(peak)$得到即可

peak1	peak2	gray_hist	img	T
第一个峰值	第二个峰值	灰度直方分布图	灰度图	初始化的threshold点

实现效果

图像灰度直方图

原图像

处理后的图像

左上 法1

右上 法2

左下 法3

右下 法4

下图为方法五的效果

自适应阈值

基本思路

➢ 将图像进一步细分为子图像，

➢ 并对不同的子图像使用不同的阈值处理

➢ 解决的关键问题：如何将图像进行细分和如何为得到的子图像估计
阈值

➢ 自适应阈值：取决于像素在子图像中的位置

◆ 算法的实现：

1. 对图像进行梯度计算，得到梯度图像。
2. 得到梯度值最大的那一部分（比如10%）的像素直方图
3. 通过直方图的谷底，得到阈值T

◆ 用拉普拉斯算子，不通过直方图，可直接得到阈值。

➢ 方法：使用拉普拉斯算子过滤图像，将0跨越点对应的灰度值为阈值T

核心函数

def adaptive(self, image):

只需要传入image的灰度图即可

实现效果

左边为原图，中间为直接调库，右边为利用滤波函数求出threshold的关联参数再处理结果

2. 实现1种基于种子点生长的分割算法

基本思路

① 根据图像的不同应用选择一个或一组种子，它或者是最亮或最暗的点，
或者是位于点簇中心的点

② 选择一个描述符（条件）

③ 从该种子开始向外扩张，首先把种子像素加入结果集合，然后不断将与
集合中各个像素连通、且满足描述符的像素加入集合

④ 上一过程进行到不再有满足条件的新结点加入集合为止

综上我们直接使用BFS即可，在函数中随意选取一个点（这里选择[259,297]这个点

核心函数

while(len(queue_point) > 0):
            print("len = ", len(queue_point))
            if count == 0:
                point = queue_point.pop(0)
                x = point[0]
                y = point[1]
            
            loc = image[y][x]
            less = loc - 8
            more = loc + 8

            for i in range(8):
                if x + dx[i] < 0 or x + dx[i] >= width \
                    or y + dy[i] < 0 or y + dy[i] >= height:
                    continue

                elif(seed_mark[y + dy[i]][x + dx[i]] != 1):

                    if less < image[y + dy[i]][x + dx[i]] < more:
                        seed_mark[y + dy[i]][x + dx[i]] = 1
                        p = [x + dx[i], y + dy[i]]
                        if p not in queue_point:
                            if  0 < p[0] < width and 0 < p[1] < height:
                                queue_point.append(p)
                    else:
                        seed_mark[y + dy[i]][x + dx[i]] = 0
            point = queue_point.pop(0)
            x = point[0]
            y = point[1]            
            count += 1

实现效果

3. 实现1种分裂合并分割算法

基本思路

1）先确定一个分裂合并的准则，即区域特征一致性的测度

2）当图像中某个区域的特征不一致时就将该区域分裂成4个相等的子区域；

3）当相邻的子区域满足一致性特征时，则将它们合成一个大区域；

4）重复进行步骤(2)和(3)直至所有区域不再满足分裂合并的条件为止

核心函数

if flag and (min(w1, h1) > 7):
    self.seperate_merge(image, w0, h0, int(w1 / 2), int(h1 / 2))
    self.seperate_merge(image, w0 + int(w1 / 2), h0, int(w1 / 2), int(h1 / 2))
    self.seperate_merge(image, w0, h0 + int(h1 / 2), int(w1 / 2), int(h1 / 2))
    self.seperate_merge(image, w0 + int(w1 / 2), h0 + int(h1 / 2), int(w1 / 2), int(h1 / 2))

实现效果

4. 附录

import numpy as np
import cv2
import matplotlib.pyplot as plt


class segmentation():
    
    def __init__(self) -> None:
        pass
    def th1(self, gray_hist, peak1, peak2):
        temp = gray_hist[int(peak1) : int(peak2)]
        valley = np.where(temp == np.min(temp))
        threshold = peak1 + valley[0][0] + 1
        return threshold
    def th2(self, gray_hist, peak1, peak2):
        temp = gray_hist[int(peak1) : int(peak2)].copy()
        if(peak2 - peak1 < 7):
            valley = np.where(temp == np.min(temp))
            threshold1 = peak1 + valley[0][0] + 1
        else:
            for i in range(peak2 - peak1):
                if(i - peak1 == 0):
                    temp[i] = temp[i:i+4].mean()
                elif(i - peak1 == 1):
                    temp[i] = temp[i-1:i+4].mean()
                elif(i - peak1 == 2):
                    temp[i] = temp[i-2:i+4].mean()

                elif(peak2 - i == 1):
                    temp[i] = temp[i-3:i+1].mean()
                elif(peak2 - i == 2):
                    temp[i] = temp[i-3:i+2].mean()
                elif(peak2 - i == 3):
                    temp[i] = temp[i-3:i+3].mean()
                else:
                    temp[i] = temp[i-3:i+4].mean()
            
            valley = np.where(temp == np.min(temp))
        return peak1 + valley[0][0] + 1
    def th3(self, gray_hist, peak1, peak2):
        edge = (peak2 - peak1)//8
        temp = gray_hist[int(peak1) + edge : int(peak2) - edge].copy()
        
        gray_temp = gray_hist.copy()
        for i in range(peak2 - peak1 - 2*edge):
            temp[i] = gray_temp[int(peak1) + edge - 5 : int(peak2) - edge + 6].mean() *\
                gray_temp[int(peak1) + edge - 5 : int(peak2) - edge + 6].std()
        valley = np.where(temp == np.min(temp))
        return peak1 + valley[0][0] + 1 + edge
    def th4_func(self,img,T):
        h=img.shape[0]
        w=img.shape[1]
        G1=G2=0
        g1=g2=0
        for i in range (h):
            for j in range (w):
                if img[i,j]>T:
                    G1+=img[i,j]
                    g1+=1
                else:
                    G2+=img[i,j]
                    g2+=1
        m1=int(G1/g1)
        m2=int(G2/g2)   # m1，m2计算两组像素均值
        T0=int((m1+m2)/2)   # 据公式计算新的阈值
        return T0
    def th4(self, img, T):
        height, width = img.shape[0:2]
        img1=np.zeros((height,width),np.uint8)
        T0=T
        T1=self.th4_func(img,T0)
        for k in range (100):   # 迭代次数为经验值，可据实际情况选定
            if abs(T1-T0)!=0:   # 若新阈值减旧阈值差值为零，则为二值图最佳阈值
                T2=self.th4_func(img,T1)
                T0=T1
                T1=T2   # 变量转换，保证if条件为新阈值减旧阈值
        return T1
        
    def basic_glo_threshold(self, image):
        
        #得到灰度直方图
        height, width = image.shape[:2]
        gray_hist = np.zeros([256], np.uint8)
        for y in range(height):
            for x in range(width):
                gray_hist[image[y][x]] += 1
        #print(gray_hist)

        

        # ① 规定一个阈值T，逐行扫描图像。
        kurtosis1 = np.where(gray_hist == np.max(gray_hist))
        #print(kurtosis1)
        peak1 = kurtosis1[0][0]; 
        
        search_peak2 = np.zeros([256], np.float32)
        for i in range(256):
            search_peak2[i] = pow(i - peak1, 2) * gray_hist[i]
            #利用公式表明远近并乘上峰值 d(i,max)^2 + hist(i)
        kurtosis2 = np.where(search_peak2 == np.max(search_peak2))
        peak2 = kurtosis2[0][0]
        
        #接下来找阈值T，两者之间的最小值
        if(peak1 > peak2):
            peak1,peak2 = peak2, peak1
        
        #方法1
        threshold = self.th1(gray_hist, peak1, peak2)
        #实际上该计算有点粗糙
        #应该是找相对谷底，而不是单纯的一点最低
        #不然会导致局部点过低而影响整个算法的图片显示
        
        #法2 用平均值
        
        threshold1 = self.th2(gray_hist, peak1, peak2)
        # print(valley)
        # print(threshold, threshold1)
        # plt.plot(np.arange(256), gray_hist)
        # plt.show()


        #法3 平均值并去边缘
        threshold2 = self.th3(gray_hist, peak1, peak2)

        #法4 直接找中位数
        threshold3 = self.th4(image, 127)

        # ② 凡灰度级大于T的，颜色置为255；凡灰度级小于T的，颜色置为0
        img1 = image.copy()
        img1[image > threshold] = 255
        img1[image <= threshold] = 0
        plt.plot(np.arange(256), gray_hist)
        plt.show()
        print(peak1, threshold, peak2)
        #cv2.imshow('1', img1)
        
        img2 = image.copy()
        img2[image > threshold1] = 255
        img2[image <= threshold1] = 0
        print(peak1, threshold1, peak2)
        #cv2.imshow('2', img)

        img3 = image.copy()
        img3[image > threshold2] = 255
        img3[image <= threshold2] = 0
        print(peak1, threshold2, peak2)
        #cv2.imshow('3', img3)

        img4 = image.copy()
        img4[image > threshold3] = 255
        img4[image <= threshold3] = 0
        print(peak1, threshold3, peak2)
        img1 = np.vstack((img1, img2))
        img3 = np.vstack((img3, img4))
        img1 = np.hstack((img1,img3))
        cv2.imshow('1 2 3 4', img1)

        threshold4 = peak1+peak2/2
        img5 = image.copy()
        img5[image > threshold4] = 255
        img5[image <= threshold4] = 0
        cv2.imshow('5', img5)
        cv2.waitKey(0)


    def adaptive(self, image):
        image_mean = cv2.boxFilter(image, cv2.CV_32FC1,(10,11))

        after = image - \
            0.87 * \
            image_mean
        print(after)
        #cv2.imshow('1', np.uint8(after))
        after[after >= 0] = 255
        after[after < 0] = 0

        dst = cv2.adaptiveThreshold(image, 255, cv2.ADAPTIVE_THRESH_GAUSSIAN_C, cv2.THRESH_BINARY, 3, 2)  # 自适应阈值
        result = np.hstack((image, dst, np.uint8(after)))
        cv2.imshow('library', result)
        cv2.waitKey(0)

    def seed(self, image):
        height, width = image.shape[:2]

        #init_point = [56,297]
        init_point = [259,297]

        seed_mark = np.zeros((height, width))
        seed_image = np.zeros((height, width))
        queue_point = [init_point]

        count = 0
        dx = [-1, 0, 1, -1, 1, -1, 0, 1]
        dy = [-1, -1, -1, 0, 0, 1, 1, 1]
        while(len(queue_point) > 0):
            if count == 0:
                point = queue_point.pop(0)
                x = point[0]
                y = point[1]
            
            loc = image[y][x]
            less = loc - 8
            more = loc + 8

            for i in range(8):
                if x + dx[i] < 0 or x + dx[i] >= width \
                    or y + dy[i] < 0 or y + dy[i] >= height:
                    continue

                elif(seed_mark[y + dy[i]][x + dx[i]] != 1):

                    if less < image[y + dy[i]][x + dx[i]] < more:
                        seed_mark[y + dy[i]][x + dx[i]] = 1
                        p = [x + dx[i], y + dy[i]]
                        if p not in queue_point:
                            if  0 < p[0] < width and 0 < p[1] < height:
                                queue_point.append(p)
                    else:
                        seed_mark[y + dy[i]][x + dx[i]] = 0
            point = queue_point.pop(0)
            x = point[0]
            y = point[1]            
            count += 1
        
        cv2.imshow('1',seed_mark)
        cv2.waitKey(0)



    def seperate_merge(self, image, w0, h0, w1, h1):
        a = image[h0: h0 + h1, w0: w0 + w1]
        #递归的时候还能改对应的numpy
        ave = np.mean(a)
        sigma = np.std(a, ddof=1)
        count = 0
        total = 0
        for x in range(w0, w0 + w1):
            for y in range(h0, h0 + h1):
                if abs(image[y, x] - ave) < 1 * sigma:
                    #灰度小于均方差
                    count += 1
                total += 1
        flag = 1 if(count / total) < 0.68 else 0
        #小于均方差的点比较多
        if flag and (min(w1, h1) > 7):
            self.seperate_merge(image, w0, h0, int(w1 / 2), int(h1 / 2))
            self.seperate_merge(image, w0 + int(w1 / 2), h0, int(w1 / 2), int(h1 / 2))
            self.seperate_merge(image, w0, h0 + int(h1 / 2), int(w1 / 2), int(h1 / 2))
            self.seperate_merge(image, w0 + int(w1 / 2), h0 + int(h1 / 2), int(w1 / 2), int(h1 / 2))
        else:
            for x in range(w0, w0 + w1):
                for y in range(h0, h0 + h1):
                    if image[y, x] > 127:
                        image[y, x] = 255
                    else:
                        image[y, x] = 0


def test():
    image = cv2.imread("work07\\pic4.png", cv2.IMREAD_GRAYSCALE)
    #image = cv2.imread("work07\\pic2.png", cv2.IMREAD_GRAYSCALE)
    #segmentation().basic_glo_threshold(image)
    #segmentation().adaptive(image)
    #segmentation().seed(image)

    copy1 = image.copy()
    h,w = copy1.shape[:2]
    segmentation().seperate_merge(copy1,0,0,w,h)
    copy1 = np.hstack((image, copy1))
    cv2.imshow('1',copy1)
    cv2.waitKey(0)
if __name__ == '__main__':
    test()