AHE_equalize图片_数字图像处理_

%this is the main function clear clc; img=((imread('24.jpg'))); figure; imshow(img); result=ahe(img,8,150); figure; imshow(result); disp 'done'; %% function img_new=ahe(img,grid,limit) img_new=img; %img_new is the img after AHE [m,n]=size(img); %ECR grid_cols=grid; grid_rows=grid; grid_width=int32(fix(m/grid_cols)); grid_height=int32(fix(n/grid_rows)); map=zeros(grid_cols,grid_rows,256); %for each grid,we create their mapping function for i=1:grid_cols for j=1:grid_rows map(i,j,:)=MakeHistogram(img,1+(i-1)*grid_width,1+(j-1)*grid_height,grid_width,grid_height,limit); end end %interpolate %boundary cases I followed the Karel Zuiderveld's implement(C version) xi = 1; for i = 1:grid_cols+1 if i == 1 subx = grid_width/2; xu = 1; xd = 1; elseif i == grid_cols+1 subx = grid_width/2; xu = grid_cols; xd = grid_cols; else subx = grid_width; xu = i - 1; xd = i; end yi = 1; for j = 1:grid_rows+1 if j == 1 suby = grid_height/2; yl = 1; yr = 1; elseif j == grid_rows+1 suby = grid_height/2; yl = grid_rows; yr = grid_rows; else suby = grid_height; yl = j - 1; yr = j; end UL = map(xu,yl,:); UR = map(xu,yr,:); DL = map(xd,yl,:); DR = map(xd,yr,:); subimg = img(xi:xi+subx-1,yi:yi+suby-1); subimg = Interpolate(subimg,UL,UR,DL,DR,subx,suby); img_new(xi:xi+subx-1,yi:yi+suby-1) = subimg; yi = yi + suby; end xi = xi + subx; end end %% function map=MakeHistogram(img,startx,starty,width,height,limit) %make histogram hist=zeros(1,256); for i=startx:startx+width-1 for j=starty:starty+height-1 value=img(i,j); hist(value+1)=hist(value+1)+1; end end %clip and allienate noise if (limit>0) excess = 0; %get the excess for degree = 1:256 excess=excess+max(0,hist(degree)-limit); end ts=limit-excess/256; avginc=excess/256; %cut them and compensate other gray level for degree=1:256 if (hist(degree)>limit) hist(degree)=limit; end %when compensating,should not larger than limit if (hist(degree)>ts) excess=excess-(limit-hist(degree)); hist(degree)=limit; else hist(degree)=hist(degree)+avginc; excess=excess-avginc; end end %if there are still some excess,equal divide them to all of the %gray levels avginc=excess/256; for degree=1:256 hist(degree)=hist(degree)+avginc; end end %create mapping function map=zeros(1,256); map(1)=hist(1); for i=2:256 map(i)=map(i-1)+hist(i); end for i=1:256 map(i) = map(i)*256/(width*height); if map(i) > 256 map(i) = 256; end end end %% function new=Interpolate(img,lutUL,lutUR,lutDL,lutDR,width,height) %interpolate points according to the following formula %f(x,y)=f(0,0)(1-x)(1-y) +f(0,1)(1-x)y+ f(1,1)xy+ f(1,0)x(1-y) new=zeros(width,height); for i=1:width for j=1:height value=img(i,j); a=double(width); c=double(i-1); b=double(height); d=double(j-1); x=c/a; y=d/b; new(i,j)=fix(lutUL(value+1)*(1-x)*(1-y)+lutDL(value+1)*x*(1-y)... +lutUR(value+1)*(1-x)*y+lutDR(value+1)*x*y); end end end % % %