scikit-image直方图匹配方法源代码分析

源代码分析

直方图匹配又称为直方图规定化，是指将一幅图像的直方图变成规定形状的直方图而进行的图像增强方法，即将某幅影像或某一区域的直方图匹配到另一幅影像上，使两幅影像的色调保持一致。该方法可以在单波段影像直方图之间进行匹配，也可以对多波段影像进行同时匹配。
scikit-image库中有一个直方图匹配方法，官方描述为：对一个图像进行调整，使其累积直方图与另一个图像的直方图相匹配。通过from skimage.exposure import match_histograms调用，该方法有一个参数channel_axis可以指定通道位置，但打开其源代码发现该参数并未使用，源代码始终使用-1作为通道参数，即指定最后一个维度为通道数。

• line 67-76

if channel_axis is not None:
    if image.shape[-1] != reference.shape[-1]:
        raise ValueError('Number of channels in the input image and '
                         'reference image must match!')

    matched = np.empty(image.shape, dtype=image.dtype)
    for channel in range(image.shape[-1]):
        matched_channel = _match_cumulative_cdf(image[..., channel],
                                                reference[..., channel])
        matched[..., channel] = matched_channel

尝试修改

我在对源代码进行分析后，认为这可能是一个bug，于是尝试对该源代码进行修改，对于传入的不同channel_axis值，可以进行不同的操作，修改内容如下图所示：

修改之后的代码可以实现下述功能：

1. 判断待修改图像与参考图像指定维度的通道数是否一致
2. 判断用户输入channels_axis是否超出维度
3. 对不同的通道位置进行不同的处理

深度分析

在进行修改后我给官方提了PR，发现很多测试都没有通过check，查看错误细节，首先可知官方测试的channel_axis有-1, 0和1，而我提交的代码并未考虑-1和1的情况（最早的提交考虑过为1的情况，但因为没有见过通道位置在第二维的图像便删去了），其次官方代码中考虑到了不同channel_axis的情况，但是histogram_matching.py中确实并未使用channel_axis指定通道位置，为了验证设置不同channel_axis是否可以得到不一样的结果，以及探究如何实现这一功能，我写了一段测试代码如下

from skimage.exposure import match_histograms
import numpy as np

image = np.array([[[143, 120, 104],
                   [143, 120, 104],
                   [141, 118, 102],
                   [45, 27, 13],
                   [125, 86, 53],
                   [161, 137, 127],
                   [161, 137, 127],
                   [162, 138, 128]]])
reference = np.array([[[154, 147, 151],
                       [109, 103, 124],
                       [63, 58, 102],
                       [127, 120, 115],
                       [180, 162, 171],
                       [0, 0, 0],
                       [1, 1, 1],
                       [0, 0, 0]]])

print('channel_axis is 0, result is\n', match_histograms(image, reference, channel_axis=0))
print('channel_axis is 1, result is\n', match_histograms(image, reference, channel_axis=1))
print('channel_axis is 2, result is\n', match_histograms(image, reference, channel_axis=2))

并且在histogram_matching.py中添加了以下语句打印channel_axis

运行test.py结果如下图所示，通过结果可知

1. 指定不同的channel_axis会有不同的结果
2. 在histogram_matching.py中channel_axis进行了处理，始终为-1

接下来进一步探究如何实现这一功能，将channel_axis设置为3运行test.py，程序报错如下

在划线部分代码可知程序通过np.moveaxis进行了移轴操作，np.moveaxis()方法介绍如下图所示，用于将数组的指定轴移动到新位置，其他轴保持原来的顺序，传入参数为：

• a：要操作的数组
• source：要移动的轴的初始位置
• destination：要移动的轴的目标位置

在histogram_matching.py中通过from .._shared import utils导入了这一方法，并通过@utils.channel_as_last_axis(channel_arg_positions=(0, 1))实现了通道转换。
打开utils.py可知官方通过划线部分代码先将原始图片的通道轴转到-1，并将传入的channel_axis参数修改为-1（此时只修改了histogram_matching.py中channel_axis的值，utils.py中变量channel_axis的值仍为用户最初传入的值，并未被修改，否则后续操作无法实现），接着执行直方图匹配操作，得到matched图像，最后再将结果中通道轴的位置从-1转到channel_axis，从而实现了对不同通道轴图片的直方图匹配。

总结

在对代码进行修改前，没有仔细检查完整代码，分析代码之间的调用关系，更为不应该的是没有使用测试代码先行进行测试，所以出现了上述的问题，特此记录下来为了提醒自己，也希望可以对在看到官方源代码后有同样困惑的用户有所帮助，不过channel_axis在histogram_matching.py中只在下面一行代码使用确实有点奇怪。·

if channel_axis is not None:

可能官方觉得如果在下述代码使用channel_axis，而在matched[..., channel]处始终将最后一个轴作为通道轴也不容易解释，所以直接使用-1了。

for channel in range(image.shape[-1]):
    matched_channel = _match_cumulative_cdf(image[..., channel],
                                            reference[..., channel])
    matched[..., channel] = matched_channel