裂缝类型选择

裂缝类型有很多种，这里我们判断类型的目的是要搞明白是否有必要检测裂缝的长度。在本文中，需要判断的裂缝类型共有四种：横向裂缝、纵向裂缝、斜裂缝、网状裂缝。

环境搭建

上一节骨架图提取部分，我们已经安装了skimage，并且在这里最好还是要有torch环境。接下来，如果你已经有的就没必要下了。

opencv-python:

​pip install opencv-python -i https://pypi.tuna.tsinghua.edu.cn/simple

opencv-contrib-python:

​pip install opencv-contrib-python -i https://pypi.tuna.tsinghua.edu.cn/simple

matplotlib:

pip install -i https://pypi.tuna.tsinghua.edu.cn/simple matplotlib==3.5.2

pyzjr:

pip install pyzjr -i https://pypi.tuna.tsinghua.edu.cn/simple

最好按照顺序来，因为pyzjr需要的依赖比较的多。

裂缝类型判断

获取图片文件

首先，我们肯定要有目标文件夹和结果保存文件夹。从文件夹中读取下面的图片路径，这里有一个bug，具体可以看这里。

masks_dir = r"E:\pythonconda2\dimension2_data\num" 
results_save_dir = "A_results"
os.makedirs(results_save_dir, exist_ok=True)
imgfile = pz.getPhotopath(masks_dir)

在我们的imgfile里面就储存了该文件夹下图片的路径。

获取最小外接矩形信息

接下来，get_minAreaRect_information函数会从二值化掩膜图像中提取最小外接矩形的相关信息，包括中心点坐标、宽高和旋转角度。inference_minAreaRect函数用于计算最小外接矩形框的宽、高和角度信息，并将角度转换为相对于图像水平方向的夹角。

def inference_minAreaRect(minAreaRect):
    w, h = minAreaRect[1]
    if w > h:
        angle = int(minAreaRect[2])
    else:
        angle = -(90 - int(minAreaRect[2]))
    return w, h, angle

def get_minAreaRect_information(mask):
    mask = pz.BinaryImg(mask)
    contours, hierarchy = cv2.findContours(mask, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    contour_merge = np.vstack(contours)
    minAreaRect = cv2.minAreaRect(contour_merge)
    return minAreaRect

pz.BinaryImg获取图像二值图，请确保读取时候为BGR的图片。 

初始化分类裂缝的参数

创建了一个ClassificationCrack类，并且对裂缝的分类参数进行初始化，有分类裂缝的阈值threshold，分类裂缝的高宽比阈值HWration，用于分类裂缝的直方图比例阈值Histration。

class ClassificationCrack():
    def __init__(self,threshold=3,HWratio=10,Histratio=0.5):
        self.threshold=threshold
        self.HWratio=HWratio
        self.Histratio=Histratio
        self.types = {0: 'Horizontal',
                      1: 'Vertical',
                      2: 'Oblique',
                      3: 'Mesh'}

这里我们使用字典self.types，这样就可以通过键值对判断裂缝的类型了。

骨骼点投影直方图

在ClassificationCrack类下，我们再定义了一个hist_judge的方法，less_than_T统计直方图中大于 0 且小于等于阈值 self.threshold 的像素数量，more_than_T统计直方图中大于阈值 self.threshold 的像素数量。通过more_than_T / （less_than_T + 1e-5）来比较是否超过了直方图比例阈值。

    def hist_judge(self, hist_v):
        less_than_T = np.count_nonzero((hist_v > 0) & (hist_v <= self.threshold))
        more_than_T = np.count_nonzero(hist_v > self.threshold)
        return more_than_T / (less_than_T + 1e-5) > self.Histratio

裂缝分类

classify 方法是 ClassificationCrack 类中的另一个成员方法，它接收三个值，minAreaRect 是一个元组，表示最小外接矩形框的信息，包括中心点坐标、宽高和旋转角度；skeleton_pts是一个数组，表示骨骼点的坐标；HW是当前 patch 的高和宽。

    def classify(self, minAreaRect, skeleton_pts, HW):
        H, W = HW
        w, h, angle = inference_minAreaRect(minAreaRect)
        if w / h < self.HWratio or h / w < self.HWratio:
            pts_y, pts_x = skeleton_pts[:, 0], skeleton_pts[:, 1]
            hist_x = np.histogram(pts_x, W)
            hist_y = np.histogram(pts_y, H)
            if self.hist_judge(hist_x[0]) and self.hist_judge(hist_y[0]):
                return 3

        return self.angle2cls(angle)

    @staticmethod
    def angle2cls(angle):
        angle = abs(angle)
        assert 0 <= angle <= 90, "ERROR: The angle value exceeds the limit and should be between 0 and 90 degrees!"
        if angle < 35:
            return 0
        elif 35 <= angle <= 55:
            return 2
        elif angle > 55:
            return 1
        else:
            return None

利用 inference_minAreaRect 函数从 minAreaRect 中获取旋转矩形框的宽度 w、高度 h 和角度 angle。接下来，通过判断 w / h 和 h / w 是否小于 self.HWratio 来判断旋转矩形框的长宽比是否满足分类条件。

如果长宽比满足条件，则将 skeleton_pts 按照 x 和 y 方向投影到直方图 hist_x 和 hist_y，然后通过 self.hist_judge 方法判断这两个直方图是否满足分类条件。以上条件均满足，则会认为是网状裂缝，否则就使用angle2cls来进行角度分类。

根据角度的大小将裂缝分为以下三类：

• 如果角度小于 35 度，则返回 0，表示水平裂缝。
• 如果角度在 35 到 55 度之间，则返回 2，表示倾斜裂缝。
• 如果角度大于 55 度，则返回 1，表示垂直裂缝。
• 如果角度不在上述范围内，则返回 None 。

主文件

if __name__ == '__main__':
    masks_dir = r"E:\pythonconda2\dimension2_data\num"  # 这里改为存放上面图片的路径
    results_save_dir = "A_results"
    os.makedirs(results_save_dir, exist_ok=True)
    classifier = ClassificationCrack()
    imgfile = pz.getPhotopath(masks_dir)
    for path in imgfile:
        mask = cv2.imread(path)
        H, W = mask.shape[:2]
        mask_copy = mask.copy()
        skeimage, skepoints = pz.ske_information(mask_copy)

        minAreaRect=get_minAreaRect_information(mask)
        timer=pz.Timer()
        w, h, angle = inference_minAreaRect(minAreaRect)
        print(w, h, angle)
        pts_y, pts_x = skepoints[:, 0], skepoints[:, 1]
        hist_x = np.histogram(pts_x, W)
        hist_y = np.histogram(pts_y, H)
        timer.stop()
        print("time cost: ", timer.total())

        result = classifier.classify(minAreaRect, skepoints, HW=(H, W))
        crack_type = classifier.types[result]
        print(crack_type)

        T = classifier.threshold

        plt.figure()
        plt.subplot(221)
        plt.imshow(mask_copy)
        plt.title("crack type: {}".format(crack_type))
        plt.subplot(222)
        plt.plot(hist_x[1][:-1], [T] * len(hist_x[0]), 'r')
        plt.bar(hist_x[1][:-1], hist_x[0])
        plt.subplot(224)
        plt.plot(hist_y[1][:-1], [T] * len(hist_y[0]), 'r')
        plt.bar(hist_y[1][:-1], hist_y[0])
        plt.subplot(223)
        plt.imshow(skeimage)
        # plt.savefig(os.path.join(results_save_dir, mask_file))
        plt.show()

检测的效果图：

图1 横向裂缝

图2 纵向裂缝

图3 斜裂缝

图4 网状裂缝

检测效果均不错，threshold，HWratio，Histratio这三个初始值均为经验所得，大家在用的时候还是要依照自己的数据来设定，可以看见左下的图，此图为没有经过消除毛刺的骨架图，但并不影响，我们采用的是骨骼点投影直方图的方法，些许毛刺影响不了我们的判断。

全文代码如下：

"""
裂缝分类如何判断
横向、纵向、网状、斜裂缝
"""
import os
import time
import matplotlib.pyplot as plt
import numpy as np
import cv2
import pyzjr as pz

def inference_minAreaRect(minAreaRect):
    w, h = minAreaRect[1]
    if w > h:
        angle = int(minAreaRect[2])
    else:
        angle = -(90 - int(minAreaRect[2]))
    return w, h, angle


class ClassificationCrack():
    def __init__(self,threshold=3,HWratio=10,Histratio=0.5):
        self.threshold=threshold
        self.HWratio=HWratio
        self.Histratio=Histratio
        self.types = {0: 'Horizontal',
                      1: 'Vertical',
                      2: 'Oblique',
                      3: 'Mesh'}

    def hist_judge(self, hist_v):
        less_than_T = np.count_nonzero((hist_v > 0) & (hist_v <= self.threshold))
        more_than_T = np.count_nonzero(hist_v > self.threshold)
        return more_than_T / (less_than_T + 1e-5) > self.Histratio

    def classify(self, minAreaRect, skeleton_pts, HW):
        H, W = HW
        w, h, angle = inference_minAreaRect(minAreaRect)
        if w / h < self.HWratio or h / w < self.HWratio:
            pts_y, pts_x = skeleton_pts[:, 0], skeleton_pts[:, 1]
            hist_x = np.histogram(pts_x, W)
            hist_y = np.histogram(pts_y, H)
            if self.hist_judge(hist_x[0]) and self.hist_judge(hist_y[0]):
                return 3

        return self.angle2cls(angle)

    @staticmethod
    def angle2cls(angle):
        angle = abs(angle)
        assert 0 <= angle <= 90, "ERROR: The angle value exceeds the limit and should be between 0 and 90 degrees!"
        if angle < 35:
            return 0
        elif 35 <= angle <= 55:
            return 2
        elif angle > 55:
            return 1
        else:
            return None

def get_minAreaRect_information(mask):
    mask = pz.BinaryImg(mask)
    contours, hierarchy = cv2.findContours(mask, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    contour_merge = np.vstack(contours)
    minAreaRect = cv2.minAreaRect(contour_merge)
    return minAreaRect


if __name__ == '__main__':
    masks_dir = r"E:\pythonconda2\dimension2_data\num"  
    results_save_dir = "A_results"
    os.makedirs(results_save_dir, exist_ok=True)
    classifier = ClassificationCrack()
    imgfile = pz.getPhotopath(masks_dir)
    for path in imgfile:
        mask = cv2.imread(path)
        H, W = mask.shape[:2]
        mask_copy = mask.copy()
        skeimage, skepoints = pz.ske_information(mask_copy)

        minAreaRect=get_minAreaRect_information(mask)
        timer=pz.Timer()
        w, h, angle = inference_minAreaRect(minAreaRect)
        print(w, h, angle)
        pts_y, pts_x = skepoints[:, 0], skepoints[:, 1]
        hist_x = np.histogram(pts_x, W)
        hist_y = np.histogram(pts_y, H)
        timer.stop()
        print("time cost: ", timer.total())

        result = classifier.classify(minAreaRect, skepoints, HW=(H, W))
        crack_type = classifier.types[result]
        print(crack_type)

        T = classifier.threshold

        plt.figure()
        plt.subplot(221)
        plt.imshow(mask_copy)
        plt.title("crack type: {}".format(crack_type))
        plt.subplot(222)
        plt.plot(hist_x[1][:-1], [T] * len(hist_x[0]), 'r')
        plt.bar(hist_x[1][:-1], hist_x[0])
        plt.subplot(224)
        plt.plot(hist_y[1][:-1], [T] * len(hist_y[0]), 'r')
        plt.bar(hist_y[1][:-1], hist_y[0])
        plt.subplot(223)
        plt.imshow(skeimage)
        # plt.savefig(os.path.join(results_save_dir, mask_file))
        plt.show()