空间域和频域的图像增强(2)--直方图均衡化与键索引计数法

接第一篇，第一篇整理了一些常用的灰度空间的函数映射变换，和一个基于比特位的Bit-plane slicing.继续往下整理；

仍然是基于point的灰度空间处理，直方图是非常常用而有效的技巧。

直方图(Histogram manipulation)可以非常有效的用于图像增强，此外，提供了对图像的统计信息，也可以用于其他一些图像处理，如图像压缩和图像分割。

较暗的图像的直方图统计成分集中在较低的灰度区间，  明亮图像的直方图统计成分集中在较高的灰度空间，低对比度的图像的直方图统计成分往往集中在灰度空间的中部，所占据的灰度区间也比较狭窄。而较好的具有高对比的图像直方图成分往往会散步整个灰度空间，而且接近均匀分布。直方图均衡化要做的事情，就是把集中与一定区间的直方图成分均匀散布到整个灰度空间，使暗图，亮图，低对比度图像成为高对比度，易于观测的图像。

直方图均衡化(Histogram Equalization)，打散密集于一定区间的灰度级至整个灰度空间，类似哈希过程：

首先考虑一个连续函数，r表示图像的灰度级，并被正则化到[0,1]区间

( s = T(r)  0<=r<=1)

对于原图的每一个灰度级r，产生一个新的灰度级s，函数(T(r))完成这样的过程，需满足的条件有：

1.(T(r))必须在定义域内是一个单值的单调递增函数;

2.(T(r))映射到的新灰度级不能跃迁到定义的灰度空间以外，即(T(r))的值域为[0,1];

r是[0,1]上的随机变量，可以用概率密度函数描述r，即(P_r(r)) 和(P_s(s))表示r ，s在[0,1]区间上的概率密度，如果已知(T(r))和(P_r(r)),由概率论知识，可以计算出(P_s(s)):

(P_s(s) = P_r(r)left|frac{dr}{ds} ight|)

因此可以得到转换公式，两边积分即可：

(s = T(r) = int_{0}^{r}P_r(w) dw)

(T(r))是随机变量r的CDF累积分布函数，并且满足作为转换函数的两个条件。

(注：此部分数学推导不是很明白，有一部分证明(P_s(s))概率密度与r无关，呈均匀分布的证明未理解);

此理论可用到离散图像的直方图均衡化，在离散图像中:

(P_r(r_k) = frac{n_k}{n}       k = 0,1,2, ... L-1 )

 n是point的总数,(n_k)是灰度级为K的point的数目。由此：

(s_k = T(r_k) = sum_{j=0}^{k} P_r(r_j) = sum_{j=0}^{k} frac{n_j}{n}   k = 0,1,2, ... L -1)

(r_k)是输入图像的灰度级， (s_k)是对应输出图像的灰度级，以上函数完成了一样一个映射。就是直方图均衡化过程。

MATLAB中实现用imhist,histeq函数

clc;
img = imread('D:MATLABimageDIP3E_Original_Images_CH03Fig0310(b)(washed_out_pollen_image).tif');
himg = imhist(img);
imshow(f);
himg1 = himg(1:256);
horz = 1:256;
figure,stem(horz,himg1,'fill')
axis([0 255 0 40000])
set(gca,'xtick',[0:10:256])
set(gca,'ytick',[0:2000:40000])
title('origin image histogram')
hqimg = histeq(img);
figure,imshow(hqimg);
hqhist = imhist(hqimg);
hqhist1=hqhist(1:256);
figure,stem(horz,hqhist1,'fill')
axis([0 255 0 40000])
set(gca,'xtick',[0:10:256])
set(gca,'ytick',[0:2000;40000]);

显示图：

 

 

 离散直方图均衡化的算法实现过程似曾相识，很容易想起一个比较轻量级但在特定情况很高效的排序算法--键索引计数法，键索引计数法是一些字符串排序算法的基础。因为是非基于比较的排序算法，所以打破了比较排序算法的下限(NlogN);

注：比较排序法的下限来自于香农的信息论;

 键索引计数法的java代码：

public class KeySort{
    public static sort(int[] a){
        int[] count = new int[a.length+1];
        int[] aux = new int[a.length];
        /*统计频率*/
        for(int i = 0;i < a.length;i++){
            count[a[i]+1]++;           
        }
        /*频率转索引*/
        for(int j = 0;j<count.length-1;j++){
            count[j+1]+=count[j];
        }
        /*回写*/
        for(int k = 0;k<a.length;k++){
            aux[count[k]++] = a[k];
        }
    }
}

 这是一个非原地排序程序，时间复杂度为(O(n))，空间复杂度为(O(N))，主要有三个部分构成，频率统计，统计过程相当于图像的直方图化，接着频率转索引，count数组相当于一个索引映射表，数组下标为原索引，值为映射索引。最后根据索引表依次安排原数组元素，相当于一个直方图均衡化过程。

由此，可以自己写出离散图像的直方图均衡化程序：

MATLAB代码：

clc;
img = imread('D:MATLABimageDIP3E_Original_Images_CH03Fig0310(b)(washed_out_pollen_image).tif');
[x,y] = size(img);
imshow(img);
Max = max(max(img));
Min = min(min(img));
count = zeros(1,256);
for i =1:x
    for j= 1:y
        %统计像素点
        count(img(i,j)+1)=count(img(i,j)+1)+1;
    end
end
figure,plot(count);
axis([0 255 0 40000]);
%直方图映射函数
p = count./(x*y);
%figure,plot(p)
%axis([0 255 0 1]);
for m=1:255
    p(m+1) = p(m)+p(m+1);
end
%离散图像映射
imghist = zeros(x,y);
for a = 1:x
    for b = 1:y
        imghist(a,b)=p(img(a,b)+1)*256;
    end
end
imghist = uint8(imghist);
figure,imshow(imghist);
count2 = zeros(x,y);
for i =1:x
    for j =1:y
        count2(imghist(i,j)+1)=count2(imghist(i,j)+1)+1;
    end
end
figure,plot(count2)
axis([0 255 0 40000]);

变换结果

 

排序的键索引技术法和图像的直方图均衡化在原理和实现上有很大的相似，只是目的不一样而已，第一步统计，第二部频率转索引都是相同的，最后一步不同，需要多加理解；