關于圖像增強必須清楚的基本概念
1.圖像增強的目的:
1)改善圖像的視覺效果,
2)轉換為更适合于人或機器分析處理的形式
3)突出對人或機器分析有意義的資訊
4)抑制無用資訊,提高圖像的使用價值
5)增強後的圖像并不一定保真
2,圖像增強的方法分類:
1)從處理對象分類:灰階圖像,(僞)彩色圖像
2)從處理政策分類:全局處理,局部處理(ROI ROI,Region of Interest Interest)
3)從處理方法分類:空間域(點域運算,即灰階變換;鄰域方法,即空域濾波),頻域方法
4)從處理目的分類:圖像銳化,平滑去噪,灰階調整(對比度增強)
3,圖像增強的方法之對比度增強
1)灰階變換法
線性變換(已實作)
對數變換(已實作)
指數變換(已實作)
2)直方圖調整法
直方圖均衡化(已實作)
直方圖比對(未實作)
一,直方圖均衡化
直方圖均衡化的英文名稱是Histogram Equalization.
圖像對比度增強的方法可以分成兩類:一類是直接對比度增強方法;另一類是間接對比度增強方法。直方圖拉伸和直方圖均衡化是兩種最常見的間接對比度增強方法。直方圖拉伸是通過對比度拉伸對直方圖進行調整,進而“擴大”前景和背景灰階的差别,以達到增強對比度的目的,這種方法可以利用線性或非線性的方法來實作;直方圖均衡化則通過使用累積函數對灰階值進行“調整”以實作對比度的增強。
直方圖均衡化處理的“中心思想”是把原始圖像的灰階直方圖從比較集中的某個灰階區間變成在全部灰階範圍内的均勻分布。直方圖均衡化就是對圖像進行非線性拉伸,重新配置設定圖像像素值,使一定灰階範圍内的像素數量大緻相同。直方圖均衡化就是把給定圖像的直方圖分布改變成“均勻”分布直方圖分布。
缺點:
1)變換後圖像的灰階級減少,某些細節消失;
2)某些圖像,如直方圖有高峰,經處理後對比度不自然的過分增強。
直方圖均衡化是圖像處理領域中利用圖像直方圖對對比度進行調整的方法。
這種方法通常用來增加許多圖像的局部對比度,尤其是當圖像的有用資料的對比度相當接近的時候。通過這種方法,亮度可以更好地在直方圖上分布。這樣就可以用于增強局部的對比度而不影響整體的對比度,直方圖均衡化通過有效地擴充常用的亮度來實作這種功能。
這種方法對于背景和前景都太亮或者太暗的圖像非常有用,這種方法尤其是可以帶來X光圖像中更好的骨骼結構顯示以及曝光過度或者曝光不足照片中更好的細節。這種方法的一個主要優勢是它是一個相當直覺的技術并且是可逆操作,如果已知均衡化函數,那麼就可以恢複原始的直方圖,并且計算量也不大。這種方法的一個缺點是它對處理的資料不加選擇,它可能會增加背景雜訊的對比度并且降低有用信号的對比度。
關于程式設計實作,同樣是不調用matlab庫函數,自己程式設計實作。這樣可以更深刻地了解直方圖均衡化技術,提高程式設計能力。
實作代碼(matlab):
clc;
close all;
clear all;
src_img = imread('flyman_gray.bmp');
figure (1)
subplot(321),imshow(src_img),title('原圖像');%顯示原始圖像
subplot(322),imhist(src_img),title('原圖像直方圖');%顯示原始圖像直方圖
matlab_eq=histeq(src_img); %利用matlab的函數直方圖均衡化
subplot(323),imshow(matlab_eq),title('matlab直方圖均衡化原圖像');%顯示原始圖像
subplot(324),imhist(matlab_eq),title('matlab均衡化後的直方圖');%顯示原始圖像直方圖
dst_img=myHE(src_img); %利用自己寫的函數直方圖均衡化
subplot(325),imshow(dst_img),title('手寫均衡化效果');%顯示原始圖像
subplot(326),imhist(dst_img),title('手寫均衡化直方圖');%顯示原始圖像直方圖
直方圖均衡化函數的實作:
function dst_img=myHE(src_img)
[height,width] = size(src_img);
dst_img=uint8(zeros(height,width));
%進行像素灰階統計;
NumPixel = zeros(1,256);%統計各灰階數目,共256個灰階級
for i = 1:height
for j = 1: width
NumPixel(src_img(i,j) + 1) = NumPixel(src_img(i,j) + 1) + 1;%對應灰階值像素點數量增加一
end
end
%計算灰階分布密度
ProbPixel = zeros(1,256);
for i = 1:256
ProbPixel(i) = NumPixel(i) / (height * width * 1.0);
end
%計算累計直方圖分布
CumuPixel = zeros(1,256);
for i = 1:256
if i == 1
CumuPixel(i) = ProbPixel(i);
else
CumuPixel(i) = CumuPixel(i - 1) + ProbPixel(i);
end
end
% 指定範圍進行均衡化
% pixel_max=max(max(I));
% pixel_min=min(min(I));
pixel_max=255;
pixel_min=0;
%對灰階值進行映射(均衡化)
for i = 1:height
for j = 1: width
dst_img(i,j) = CumuPixel(src_img(i,j)+1)*(pixel_max-pixel_min)+pixel_min;
end
end
return;
為什們和matlab的直方圖不一樣呢???
二,指數變換
指數變換(Power-Law )的公式:S=c*R^r, 通過合理的選擇c和r可以壓縮灰階範圍,算法以c=1.0/255.0, r=2實作。
要做該圖像增強變換需要先做歸一化,再指數變換,最後反歸一化 增強效果展示:可以看見,改增強算法并不能很好的将像素盡可能的碾平。
指數增強參考程式為:
clc;
close all;
clear all;
% -------------Gamma Transformations-----------------
%f = imread('Fig0316(4)(bottom_left).tif');
f = imread('seed.tif');
Gamma = 0.4;
g2 = myExpEnhance(f,Gamma);
figure();
subplot(221); imshow(f); xlabel('a).Original Image');
subplot(222),imhist(f),title('原圖像直方圖');%顯示原始圖像直方圖
subplot(223); imshow(g2); xlabel('b).Gamma Transformations \gamma = 0.4');
subplot(224),imhist(g2),title('增強圖像直方圖');%顯示原始圖像直方圖
指數增強核心函數為:
function dst_img=myExpEnhance(src_img,Gamma)
src_img = mat2gray(src_img,[0 255]);%将圖像矩陣A中介于amin和amax的資料歸一化處理, 其餘小于amin的元素都變為0, 大于amax的元素都變為1。
C = 1;
g2 = C*(src_img.^Gamma);
%反歸一化
max=255;
min=0;
dst_img=uint8(g2*(max-min)+min);
三,對數變換
對數變換主要用于将圖像的低灰階值部分擴充,将其高灰階值部分壓縮,以達到強調圖像低灰階部分的目的。變換方法由下式給出。
這裡的對數變換,底數為(v+1),實際計算的時候,需要用換底公式。其輸入範圍為歸一化的【0-1】,其輸出也為【0-1】。對于不同的底數,其對應的變換曲線如下圖所示。
底數越大,對低灰階部分的強調就越強,對高灰階部分的壓縮也就越強。相反的,如果想強調高灰階部分,則用反對數函數就可以了。看下面的實驗就可以很直覺的了解,下圖是某圖像的二維傅裡葉變換圖像,其為了使其灰階部分較為明顯,一般都會使用灰階變換處理一下。
效果圖:
參考代碼:
clc;
close all;
clear all;
%-------------Log Transformations-----------------
f = imread('seed.tif');
g_1 = myLogEnhance(f,10);
g_2 = myLogEnhance(f,100);
g_3 = myLogEnhance(f,200);
figure();
subplot(2,2,1);
imshow(f);xlabel('a).Original Image');
subplot(2,2,2);
imshow(g_1);xlabel('b).Log Transformations v=10');
subplot(2,2,3);
imshow(g_2);xlabel('c).Log Transformations v=100');
subplot(2,2,4);
imshow(g_3);
xlabel('d).Log Transformations v=200');
對數變換核心函數
function dst_img=myLogEnhance(src_img,v)
c=1.0;
src_img = mat2gray(src_img,[0 255]);
g =c*log2(1 + v*src_img)/log2(v+1);
%反歸一化
max=255;
min=0;
dst_img=uint8(g*(max-min)+min);
四,灰階拉伸
灰階拉伸也用于強調圖像的某個部分,與伽馬變換與對數變換不同的是,灰階拉升可以改善圖像的動态範圍。可以将原來低對比度的圖像拉伸為高對比度圖像。實作灰階拉升的方法很多,其中最簡單的一種就是線性拉伸。而這裡介紹的方法稍微複雜一些。灰階拉伸所用數學式如下所示。
同樣的,其輸入r為【0-1】,其輸出s也為【0-1】。這個式子再熟悉不過了,跟巴特沃斯高通濾波器像極了,其輸入輸出關系也大緻能猜到是個什麼形狀的。但是,這裡就出現一個問題了,輸入為0時候,式子無意義了。是以,在用Matlab計算的時候,将其變為如下形式。
這裡的eps,就是Matlab裡面,一個很小數。如此做的話,式子變得有意義了。但是,其輸入範圍為【0-1】的時候,其輸出範圍變為了
。輸出範圍大緻為【0-1】,為了精确起見,使用mat2gray函數将其歸一化到精确的[0-1]。調用格式如下。
五,線性拉伸
為了突出感興趣的目标或者灰階區間,相對抑制那些不感興趣的灰階區域,可采用分段線性法,常用的是三段線性變換
參考程式:
clc;
close all;
clear all;
I=imread('seed.tif');
[m,n,k]=size(I);
figure (1)
imshow('seed.tif');title(' 原圖像');
mid=mean(mean(I));
%橫軸
fa=20; fb=80;
%縱軸
ga=50; gb=230;
J=myLinearEnhance(I,fa,fb,ga,gb);
figure (2)
imshow(J);title(' 線性拉伸圖像');
pixel_f=1:256;
pixel_g=zeros(1,256);
%三段斜率,小于1表示該段将會被收縮
k1=double(ga/fa);
k2=(gb- ga)/(fb- fa);
k3=(256- gb)/(256- fb);
for i=1:256
if i <= fa
pixel_g(i)= k1*i;
elseif fa < i && i <= fb
pixel_g(i)= k2*( i- fa)+ ga;
else
pixel_g(i)= k3*( i - fb)+ gb;
end
end
figure (3)
plot(pixel_f,pixel_g);
核心函數:
function dst_img=myLinearEnhance(src_img,fa,fb,ga,gb)
[height,width] = size(src_img);
dst_img=uint8(zeros(height,width));
src_img=double(src_img);
%三段斜率
k1=ga/fa;
k2=(gb- ga)/(fb- fa);
k3=(255- gb)/(255- fb);
for i=1:height
for j=1:width
if src_img(i,j) <= fa
dst_img(i,j)= k1*src_img(i,j);
elseif fa < src_img(i,j) && src_img(i,j) <= fb
dst_img(i,j)= k2*( src_img(i,j)- fa)+ ga;
else
dst_img(i,j)= k3*( src_img(i,j)- fb)+ gb;
end
end
end
dst_img=uint8(dst_img);
附錄:
附錄網上的另一份講解:
直方圖均衡化算法分為三個步驟,第一步是統計直方圖每個灰階級出現的次數,第二步是累計歸一化的直方圖,第三步是計算新的像素值。
第一步:
for(i=0;i<height;i++)
for(j=0;j<width;j++)
n[s[i][j]]++;
for(i=0;i<L;i++)
p[i]=n[i]/(width*height);
這裡,n[i]表示的是灰階級為i的像素的個數,L表示的是最大灰階級,width和height分别表示的是原始圖像的寬度和高度,是以,p[i]表示的就是灰階級為i的像素在整幅圖像中出現的機率(其實就是p[]這個數組存儲的就是這幅圖像的歸一化之後的直方圖)。
第二步:
for(i=0;i<=L;i++)
for(j=0;j<=i;j++)
c[i]+=p[j];
c[]這個數組存儲的就是累計的歸一化直方圖。
第三步:
max=min=s[0][0];
for(i=0;i<height;i++)
for(j=0;j<width;j++)
if(max<s[i][j]){
max=s[i][j];
}else if(min>s[i][j]){
min=s[i][j];
}
找出像素的最大值和最小值。
for(i=0;i<height;i++)
for(j=0;j<width;j++)
t[i][j]=c[s[i][j]]*(max-min)+min;
t[][]就是最終直方圖均衡化之後的結果。
收錄優秀代碼:
這份代碼寫得不錯,學習了,原部落格位址見參考資源【3】!
#include <stdio.h>
#include <iostream>
#include "fftw3.h"
#include "string"
#include "vector"
#include <windows.h>
#include <opencv2/legacy/legacy.hpp>
#include <opencv2/nonfree/nonfree.hpp>//opencv_nonfree子產品:包含一些擁有專利的算法,如SIFT、SURF函數源碼。
#include "opencv2/core/core.hpp"
#include "opencv2/features2d/features2d.hpp"
#include "opencv2/highgui/highgui.hpp"
#include <opencv2/nonfree/features2d.hpp>
using namespace cv;
using namespace std;
class hisEqt
{
public:
hisEqt::hisEqt();
hisEqt::~hisEqt();
public:
int w;
int h;
int nlen;
int *pHis;
float *pdf;
//=====求像素分布機率密度====
void getPdf();
//======統計像素個數=======
void getHis(unsigned char*imgdata);
//==========畫統計分布直方圖===============
void drawHistogram(const float*pdf,Mat &hist1);
//===========直方圖均衡化==========
void hisBal();
//====直方圖均衡化後的圖像===
void imgBal(unsigned char* img);
};
hisEqt::hisEqt() :nlen(0){
pHis = new int[256 * sizeof(int)];
memset(pHis, 0, 256 * sizeof(int));
pdf = new float[255 * sizeof(float)];
memset(pdf, 0, 255 * sizeof(float));
}
hisEqt::~hisEqt(){
delete[]pHis;
delete[]pdf;
}
//======統計像素個數=======
void hisEqt::getHis(unsigned char*imgdata){
for (int i = 0; i<nlen; i++)
{
pHis[imgdata[i]]++;
}
}
//=====求像素分布機率密度====
void hisEqt::getPdf(){
for (int k = 0; k<256; k++)
{
pdf[k] = pHis[k] / float(nlen);
}
}
//===========直方圖均衡化==========
void hisEqt::hisBal(){
for (int k = 1; k<256; k++)
{
pdf[k] += pdf[k - 1];
}
for (int k = 0; k<256; k++)
{
pHis[k] = 255 * pdf[k];
}
}
//====直方圖均衡化
void hisEqt::imgBal(unsigned char* img){
for (int i = 0; i<nlen; i++)
{
img[i] = pHis[img[i]];
}
}
void hisEqt::drawHistogram(const float *pdf, Mat& hist1){
for (int k = 0; k<256; k++)
{
if (k % 2 == 0)
{
Point a(k, 255), b(k, 255 - pdf[k] * 2550);
line(hist1,
a,
b,
Scalar(0, 0, 255),
1);
}
else
{
Point a(k, 255), b(k, 255 - pdf[k] * 2550);
line(hist1,
a,
b,
Scalar(0, 255, 0),
1);
}
}
}
int main()
{
Mat image = imread("Fig0651(a)(flower_no_compression).tif");
if (!image.data)
return -1;
Mat hist2(256, 256, CV_8UC3, Scalar(0, 0, 0));
Mat hist1(256, 256, CV_8UC3, Scalar(0, 0, 0));
Mat imgOut = Mat(image.rows, image.cols, CV_8UC3, Scalar(0, 0, 0));
vector<Mat> planes;
int chn = image.channels();
if (chn == 3)
{
split(image, planes);
}
while (chn)
{
chn--;
unsigned char* imageData = new unsigned char[sizeof(unsigned char)*(image.cols*image.rows)];
memcpy(imageData, planes[chn].data, planes[chn].cols*planes[chn].rows);
hisEqt his;//自定義的類
his.nlen = image.rows*image.cols;
his.getHis(imageData);
his.getPdf();
// //======畫原圖直方圖并儲存============
his.drawHistogram(his.pdf, hist1);
string pic_name = "hisline";
pic_name = pic_name + to_string(chn);
pic_name=pic_name+ ".jpg";
imwrite(pic_name, hist1);
his.hisBal();
his.getPdf();
// //======畫均衡化後直方圖并儲存============
his.drawHistogram(his.pdf, hist2);
string pic_name0 = "his_balanceline";
pic_name0 = pic_name0 + to_string(chn);
pic_name0 = pic_name0 + ".jpg";
imwrite(pic_name0, hist2);
// //=====圖像均衡化===
his.imgBal(imageData);
memcpy(planes[chn].data, imageData, planes[chn].cols*planes[chn].rows);
delete[] imageData;
imageData = NULL;
}
merge(planes, imgOut);//單通道合并
imwrite("result.jpg", imgOut);
return 0;
}
參考資源
【1】http://blog.csdn.net/xiajun07061225/article/details/6910129
【2】數字圖像處理,岡薩雷斯著
【3】http://blog.csdn.net/bettyshasha/article/details/46940805
【4】http://blog.csdn.net/terryzero/article/details/6043821
【5】http://www.myexception.cn/image/1450848.html