python实现灰度共生矩阵

import math
import itertools
import numpy as np
import cv2


def calculate_glcm(patch, ranges, dires, nb_bins, normalized=False, symmetric=False, SUM=False):
    """计算灰度共生矩阵

    Args:
        patch ([ndarray]): 输入的图像，可分成小块
        ranges ([list]]): 像素相邻间隔的大小
        dirs ([list]): 遍历的方向。改用1,2,...，7，8对应0,45,...，270,315度，
        nb_bins ([int]): 灰度级
        normalized (bool, optional): 是否归一化. Defaults to False.
        symmetric (bool, optional): 在对称模式下运行将删除相反的方向，
            并将用对应的方向替换5到8 !为了避免混淆输出，请确保只在对称模式下使用方向1 - 4.
             Defaults to False.
        SUM (bool, optional): 是否将所有需要的方向被加在一起. Defaults to False.

    Returns:
        [ndarray]: 返回灰度共生矩阵.该glcm数组的形状为:[距离][方向][nb_bins][nb_bins]
    """    
    nb_xdims = patch.shape[0]
    nb_ydims = patch.shape[1]
    
    coords = [[0, 0], [0, 1], [1, 1], [1, 0], [1, -1], [0, -1], [-1, -1], [-1, 0], [-1, 1]]
    # 所有需要的方向被加在一起，所以每距离只生成一个glcm。
    if (SUM):
        glcm = np.zeros((len(ranges), nb_bins, nb_bins)).astype(np.int32)

        for dist, dire, px, py in itertools.product(range(len(ranges)),
                                                   range(len(dires)),
                                                   range(nb_xdims),
                                                   range(nb_ydims)):
            direction = coords[dires[dire]] * ranges[dist]

            if nb_xdims > px + direction[0] >= 0 and nb_ydims > py + direction[1] >= 0:
                src = patch[px, py]
                dst = patch[px + direction[0], py + direction[1]]

                if (symmetric):
                    glcm[dist, dst, src] += 1
                glcm[dist, src, dst] += 1

        if (normalized):
            for dist in itertools.product(range(len(ranges))):
                glcm = glcm.astype(np.float64)
                glcm[dist, :, :] /= np.sum(glcm[dist, :, :])

        return glcm


    else:
        glcm = np.zeros((len(ranges), len(dires), nb_bins, nb_bins)).astype(np.int32)

        for dist, dire, px, py in itertools.product(range(len(ranges)),
                                                   range(len(dires)),
                                                   range(nb_xdims),
                                                   range(nb_ydims), ):
            direction = coords[dires[dire]] * ranges[dist]

            if nb_xdims > px + direction[0] >= 0 and nb_ydims > py + direction[1] >= 0:
                src = patch[px, py]
                dst = patch[px + direction[0], py + direction[1]]

                if (symmetric):
                    glcm[dist, dire, dst, src] += 1
                glcm[dist, dire, src, dst] += 1

        if (normalized):
            for dist, dire in itertools.product(range(len(ranges)),
                                               range(len(dires))):
                glcm = glcm.astype(np.float64)
                glcm[dist, dire, :, :] /= np.sum(glcm[dist, dire, :, :])

        return glcm