[数字图像处理]灰度变换——直⽅图处理
直⽅图均衡
⼀副图像的直⽅图，表⽰了其灰度分布的特性。对于数字图像来说，假设灰度值k出现了次，那么其概率密度函数如下所⽰。
这个式⼦，表⽰了像素的灰度值为k概率。其中，M与N为图像的尺⼨。
对于⼀幅动态范围较窄的图像，其归⼀化灰度直⽅图如下所⽰。
对于此类图像，使⽤灰度拉伸，也能使其的动态范围得到改善，增强图像对⽐度。但是，灰度拉伸从本质上来讲，是很暧昧的，灰度拉伸没有什么明确的⽬的，只是依靠某个函数的做了灰度变换，使得图像的在灰度直⽅图的分布范围扩⼤。灰度拉伸对于结果没有严格的要求，因此根据变换函数的不同，灰度
拉伸有⽆数种结果。
与灰度拉伸不同的是，直⽅图均衡的⽬的就⽐较明确，使⽤某个特殊的函数，使原图像的灰度分布平均化。为了将其平均化，这⾥需要⽤到累积分布函数（Probability Density Function）的概念。
所谓的平均化，直⽅图均衡的⽬的是将⼀幅图的累积分布的曲线（下图左），变为下图右的分布形式。
（左）原图的累积分布曲线                        （右）理想化直⽅图均衡后的累积分布曲线
从右图中可以看出，为了保证直⽅图均衡后的累积分布曲线是⼀条倾⾓为45°的直线，那么，变换后的函数的概率密度函数应与输⼊值⽆关，且保持⼀个常量。输⼊r与输出s有如下关系，
将变换关系
带⼊上式⼦，可得
将其带⼊输⼊r与输出s关系，可得到
到此，已经可以看出，输出的概率密度函数，与s的取值⽆关。很明显，这个分布是绝对平均的。
但是，事实上，上⾯的式⼦均是在连续的情况下的。在数字图像的领域，我们必须将其离散化，才能使⽤。离散的情况下，其累计分布函数如下所⽰。
⾄此，直⽅图均衡就完成了。
为了验证结果，其直⽅图均衡后的图像与灰度直⽅图如下所⽰。
结果就是这样了。但是，呃问题来了，这是直⽅图均衡么？这跟灰度拉升有何区别，⽽且其均衡后的图像的直⽅图也不是⼀个很平均的值，和我们的构想有很⼤区别。这个结果真的没问么？为了回答这个问题，将其原图像的累积分布函数与直⽅均衡后的累积分布函数画出来，如下。
（左）原图的累积分布曲线                                            （右）直⽅图均衡后的累积分布曲线
从累积分布曲线中，我们可以看到，直⽅图均衡后的曲线，确实跟我们设想的理想曲线很接近！由于我们证明的时候，是在连续的情况下去思考的，所以，将其离散之后，多少是有不理想的情况发⽣的。均衡后的累积分布曲线确实是⼀个不错的结果。这样看来，我们才能看清楚直⽅图均衡的本质。
若不是这样，还真不知道直⽅均衡与灰度拉伸的区别。
为此，直⽅图均衡成功，Matlab代码如下。
close all;
clear all;
%% -------------Histogram Equalization-----------------
close all;
clear all;
f = imread('washed_out_pollen_image.tif');
f = mat2gray(f,[0 255]);
[M,N] = size(f);
g = zeros(M,N);
r = imhist(f)/(M*N);  %(500*500 is size of image)
s = zeros(256,1);
for k = 1:1:256
for j = 1:1:k
s(k,1) = s(k,1) + r(j);
end
end
end
for x = 1:1:M        %(500*500 is size of image)
for y = 1:1:N
g(x,y) = s(uint8(f(x,y)*255)+1,1);
end
end
figure();
subplot(2,2,1);
imshow(f);
xlabel('a).Original Image');
subplot(2,2,2);
h = imhist(f)/(M*N);
matlab直方图bar(0:1/255:1,h);
axis([0 1 0 0.1]),grid;
%axis square;
xlabel('b).The Histogram of a');
ylabel('Number of pixels');
subplot(2,2,3);
imshow(g);
xlabel('c).Histogram Equalization');
subplot(2,2,4);
h = imhist(g)/(M*N);
bar(0:1/255:1,h);
axis([0 1 0 0.1]),grid;
%axis square;
xlabel('d).The Histogram of c');
ylabel('Number of pixels');
figure();
x=0:1/255:1;
plot(x,s);
axis([0,1,0,1]),grid;
axis square;
xlabel('intensity level of r');
ylabel('P_{r}(r)');
r = imhist(g)/(M*N);
s = zeros(256,1);
for k = 1:1:256
for j = 1:1:k
s(k,1) = s(k,1) + r(j);
end
end
figure();
x=0:1/255:1;
plot(x,s);
axis([0,1,0,1]),grid;
axis square;
xlabel('intensity level of s');
ylabel('P_{s}(s)');
（作者注：其实，直⽅图均衡有直接的函数，这⾥只是想按照⾃⼰理解的⽅法去做个均衡出来看看，所以⾼级⼀点的库函数基本没有，都⽤循环来的）
直⽅图匹配
前⾯，我们说了直⽅图的均衡，直⽅均衡可以得到⼀个灰度直⽅图分布平均的图像。但是有⼀些问题，直⽅图均衡的结果是唯⼀的。这是其优点，也是其缺点。优点在于，⽆需多余的参数，就可以实现对⽐度的增强。缺点在于，没有参数，对于结果⽆法调整。举个栗⼦（输
⼊法君卖萌了╭(╯3╰)╮），⼀幅很⿊的图像，使⽤直⽅图均衡的话，这幅图像的直⽅均衡结果，将使得图像偏⽩。使⽤ 《Digital Image Processing》 Rafael C. Gonzalez / Richard E. Woods 的例⼦来说，看下图。
如上图所⽰，直⽅均衡的结果不是很理想。由直⽅均衡的结果看来，我们能看到，其实对⽐原图，原图⼀些细节已经显露出来了，只要稍加修改，就可以得到很理想的增强结果。为此，直觉告诉我们（不要吐槽这句话，有时候直觉就是那么简单。\(≧▽≦)/），只需要将原图的灰度直⽅图稍微左移⼀些即可。但是，直⽅均衡却是不可调整的，为此，就有了直⽅图匹配。
对于原图像，其概率分布函数是
其输出也就是直⽅均衡的结果。为了实现直⽅图匹配，我们还需要⼀个函数，也就是我们所期待的概率分布函数。这个概率分布函数的概率分布函数如下所⽰。
其中，z代表了期望的分布。我们假设如下关系成⽴。
也就是，我们所期待的直⽅图分布的概率分布函数与原图像的直⽅图概率分布函数相等。那么，我们可以利⽤G()的反函数，我们也就可以得到我们所期待的直⽅图。
上⾯的描述可能太过拖沓，我们采⽤简洁⼀点的描述。
⼀. 将原图进⾏直⽅图均衡。
⼆. 假设我们期待的分布，直⽅均衡的结果与原图直⽅图均衡的结果⼀致。利⽤反函数，还原成我们期待的分布。
所实现的直⽅图均衡结果如下所⽰。
根据直⽅图匹配的结果，我们可以看出，所得出的结果和我们所期待的分布很接近。其结果也和我们的初衷⼀样，使得原图的直⽅图稍微左移了⼀些。所得到的结果⽐起直⽅均衡要好得多。
其Matlab代码如下。
%% -------------Histogram Matching-----------------
close all;
clear all;
f = imread('mars_moon_phobos.tif');
f = mat2gray(f,[0 255]);
[M,N] = size(f);
g = zeros(M,N);
r = imhist(f)/(M*N);
s = zeros(256,1);
for k = 1:1:256
for j = 1:1:k
s(k,1) = s(k,1) + r(j);
end
end
for x = 1:1:M
for y = 1:1:N
g(x,y) = s(uint8(f(x,y)*255)+1,1);
end
end
figure();
subplot(2,2,1);
imshow(f);
xlabel('a).Original Image');
subplot(2,2,2);
h = imhist(f)/(M*N);
bar(0:1/255:1,h);
axis([0 1 0 0.1]),grid;
xlabel('b).The Histogram of a');
subplot(2,2,3);
imshow(g);
xlabel('c).Histogram Equalization');
subplot(2,2,4);
h = imhist(g)/(M*N);
bar(0:1/255:1,h);
axis([0 1 0 0.1]),grid;
xlabel('d).The Histogram of c');
%---------------pz-----------------------
r = r';
p = [0:r(1,1)/15:r(1,1) 17*r(1,2:241)];
%p = [0:r(1,1)/15:r(1,1) r(1,2:241)];
%p = [0:55:255 , 255:-25:24 , 24:-0.14:0 , 0:1:15 , 15:-0.285:0];
p = (p/sum(p));
p = p';
r = r';
G = zeros(256,1);
for k = 1:1:256
for j = 1:1:k
G(k,1) = G(k,1) + p(j);
end

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。