数字图像处理——直方图匹配

步骤：

1.求原图像和匹配图像的灰度分布直方图

2.根据灰度变换函数分别构建原图像和匹配图像的映射表

3.根据单映射规则（最近映射）构建原图像到目标图像的映射表

4.根据映射表生成匹配图像

Matlab代码：

clear,clc;
f = imread('LENA.png');
f_ref = imread('EightAM.png');
[h1,w1] = size(f);
[h2,w2] = size(f_ref);

% 求图像f的灰度分布直方图
hist1 = zeros(1,256);   % hist1为灰度分布向量，0无法作为索引，索引1对应灰度值0
for row = 1:h1
    for col = 1:w1
        hist1(f(row,col)+1) = hist1(f(row,col)+1)+1;    
    end
end
% 求直方图均衡的变换函数
f_table = zeros(1,256);   % f_table为映射表，索引i代表灰度值+1，f_table[i+1]代表原来的灰度值i经过变换后的灰度值
cum_sum = 0;
for index = 1:256
    cum_sum = cum_sum + hist1(index);
    f_table(index) = (255/(h1*w1))*cum_sum;
end

% 求图像f_ref的灰度分布直方图
hist1_ref = zeros(1,256);   % hist1_ref为灰度分布向量，0无法作为索引，索引0对应灰度值1
for row = 1:h2
    for col = 1:w2
        hist1_ref(f_ref(row,col)+1) = hist1_ref(f_ref(row,col)+1)+1;    % 0无法作为索引，灰度值0对应下标1
    end
end
% 求直方图均衡的变换函数
f_ref_table = zeros(1,256);   % f_ref_table为映射表，索引i代表灰度值+1，f_ref_table[i+1]代表原来的灰度值i经过变换后的灰度值
cum_sum = 0;
for index = 1:256
    cum_sum = cum_sum + hist1_ref(index);
    f_ref_table(index) = (255/(h2*w2))*cum_sum;
end

% 以下四条语句等效于上述求两幅图像的累计分布函数
% hist1 = imhist(f);
% hist1_ref = imhist(f_ref);
% f_table = cumsum(hist1)/numel(f);
% f_ref_table = cumsum(hist1_ref)/numel(f_ref);

% 根据单映射规则构造映射表
map_table = zeros(1,256);
for index = 1:256
    [temp,ind] = min(abs(f_table(index)-f_ref_table));
    map_table(index) = ind-1;
end

% 根据映射表生成匹配后的图像
f_match = zeros(h1,w1);  % 先初始化f_match
for row = 1:h1
    for col = 1:w1
        f_match(row,col) = map_table(double(f(row,col)+1));  % 先进行类型转换，防止溢出
    end
end
% 上述循环等效于f_match = map_table(double(f)+1);
f_match = uint8(f_match);

% 显示原图像，匹配图像，匹配后的图像
figure;
subplot(1,3,1),imshow(f),title('原图像');
subplot(1,3,2),imshow(f_ref),title('匹配图像');
subplot(1,3,3),imshow(f_match),title('匹配后的图像');

% 显示原图像，匹配图像，匹配后的图像的直方图
figure;
subplot(1,3,1),imhist(f),title('原直方图');
subplot(1,3,2),imhist(f_ref),title('匹配图像的直方图');
subplot(1,3,3),imhist(f_match),title('匹配后图像的直方图');

测试图像：左边是待处理图像，右边是匹配图像

运行结果：