基于Mass-Radius法的分形维度计算代码优化需求

需求背景 #

通过Mass-Radius法计算图像（示例为Ottawa.png）的分形维度，核心步骤为：

1. 寻找图像中非黑色像素的质心
2. 统计不同半径圆内的非黑色像素数量
3. 将半径与像素数量取对数后存入列表
4. 对对数数据做线性回归，回归系数a即为分形维度

原代码运行速度极慢，以下是针对性优化方案：

原代码性能瓶颈分析 #

• 嵌套循环逐像素遍历，时间复杂度随图像尺寸呈指数增长
• 逐像素判断是否为黑色，未利用NumPy向量化运算的高效性
• D_R函数中重复遍历环形区域，存在大量冗余计算

优化方案 #

1. 向量化处理像素判断 #

用NumPy批量判断所有像素是否为黑色，替代逐像素循环判断。

2. 质心计算优化 #

直接利用NumPy数组运算计算质心，避免嵌套循环累加坐标。

3. 预计算像素到质心的距离 #

一次性计算所有非黑色像素到质心的欧氏距离，后续仅需对距离数据统计，无需重复遍历图像。

4. 直方图统计半径区间像素数 #

用NumPy的histogram函数高效统计不同半径区间内的像素数量，替代环形区域的嵌套遍历。

优化后的代码 #

from scipy.stats import linregress
from skimage import io
import numpy as np

def main():
    # 读取图像
    img = io.imread('Ottawa.png')
    height, width, _ = img.shape

    # 1. 生成非黑色像素的掩码（向量化判断）
    is_non_black = ~np.all(img == [0, 0, 0], axis=2)
    # 获取所有非黑色像素的坐标
    y_coords, x_coords = np.where(is_non_black)
    total_non_black = len(x_coords)

    if total_non_black == 0:
        print("图像中无非黑色像素")
        return

    # 2. 计算质心（向量化计算）
    centroid_x = int(np.mean(x_coords))
    centroid_y = int(np.mean(y_coords))

    # 3. 计算所有非黑色像素到质心的距离
    dx = x_coords - centroid_x
    dy = y_coords - centroid_y
    distances = np.sqrt(dx**2 + dy**2)

    # 确定最大半径（质心到图像边缘的最大距离）
    max_radius = int(max(centroid_x, width - centroid_x, centroid_y, height - centroid_y))
    # 生成半径区间：1到max_radius的整数
    radius_bins = np.arange(1, max_radius + 1)

    # 4. 统计每个半径内的累计像素数（直方图累加）
    counts, _ = np.histogram(distances, bins=np.append(radius_bins, max_radius + 1))
    cumulative_counts = np.cumsum(counts)

    # 过滤掉像素数为0的半径（避免log(0)错误）
    valid_mask = cumulative_counts > 0
    log_radii = np.log(radius_bins[valid_mask])
    log_counts = np.log(cumulative_counts[valid_mask])

    # 线性回归计算分形维度
    if len(log_radii) < 2:
        print("有效数据点不足，无法进行回归")
        return
    a, b, r_value, p_value, std_err = linregress(log_radii, log_counts)

    print(f"分形维度a: {a:.4f}")
    print(f"回归截距b: {b:.4f}")
    print(f"R²值: {r_value**2:.4f}")

if __name__ == "__main__":
    main()

优化效果说明 #

• 完全摒弃嵌套循环，所有核心计算均采用NumPy向量化操作，速度提升数十倍甚至上百倍
• 预计算距离后仅需一次统计，避免了原代码中多次遍历图像的冗余操作
• 代码结构更简洁，逻辑清晰，可维护性更强

内容的提问来源于stack exchange，提问作者Jotadiolyne Dicci