1 內(nèi)容介紹

本文介紹了一種基于和聲搜索算法（HSA）的多級(jí)閾值（MT）算法。HSA 是一種進(jìn)化方法，其靈感來自音樂家在演奏時(shí)即興創(chuàng)作新的和聲。與其他進(jìn)化算法不同，HSA 展示了有趣的搜索能力，仍然保持較低的計(jì)算開銷。所提出的算法將來自圖像直方圖中可行搜索空間的隨機(jī)樣本編碼為候選解決方案，而考慮 Otsu 或 Kapur 方法所采用的目標(biāo)函數(shù)來評(píng)估它們的質(zhì)量。在這些目標(biāo)值的指導(dǎo)下，候選解集通過 HSA 算子進(jìn)行演進(jìn)，直到找到最優(yōu)解。實(shí)驗(yàn)結(jié)果證明了所提出的數(shù)字圖像分割方法的高性能。

2 仿真代碼

%Intructions:

% I -> Original Image, could be RGB or Gray Scale

% level -> Number of threshold to find

% This version works with KAPUR as fitness function.

close all

clear all

% Se carga la imagen RGB o escala de grises

I1 = imread(’Picture 148710088.jpg’);

I=rgb2gray(I1);

level = 3;

% Se obtienen los histogramas si la imagen es RGB uno por cada canal si es

% en escala de grises solamente un historgrama.

if size(I,3) == 1 %grayscale image

    [n_countR, x_valueR] = imhist(I(:,:,1));

elseif size(I,3) == 3 %RGB image

    %histograma para cada canal RGB

    [n_countR, x_valueR] = imhist(I(:,:,1));

    [n_countG, x_valueG] = imhist(I(:,:,2));

    [n_countB, x_valueB] = imhist(I(:,:,3));

end

Nt = size(I,1) * size(I,2); %Cantidad total de pixeles en la imagen RENG X COL

%Lmax niveles de color a segmentar 0 - 256

Lmax = 256;   %256 different maximum levels are considered in an image (i.e., 0 to 255)

% Distribucion de probabilidades de cada nivel de intensidad del histograma 0 - 256 

for i = 1:Lmax

    if size(I,3) == 1 

        %grayscale image

        probR(i) = n_countR(i) / Nt;

    elseif size(I,3) == 3 

        %RGB image    

        probR(i) = n_countR(i) / Nt;

        probG(i) = n_countG(i) / Nt;

        probB(i) = n_countB(i) / Nt;

    end

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Parametros del problema de segmentacion

N_PAR = level; %number of thresholds (number of levels-1) (dimensiones)

ndim = N_PAR;  

%Parametros Harmony Search %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

MaxAttempt = 25000;  % Max number of Attempt

% Initial parameter setting

HS_size = 50;        %Length of solution vector

HMacceptRate = 0.95; %HM Accepting Rate

PArate = 0.5;        %Pitch Adjusting rate

if size(I,3) == 1 

    %Imagen escala de grises

    range = ones(ndim,2);

    range(:,2) = range(:,2) * Lmax;

    

    %initializa harmony memory

    HM = zeros(HS_size,ndim);

    

    % Pitch range for pitch adjusting

    pa_range = ones(ndim);

    pa_range = pa_range * 100;

elseif size(I,3) == 3

    %Imagen RGB

    range = ones(ndim,2);

    range(:,2) = range(:,2) * Lmax;

    %IR

    xR = zeros(HS_size,ndim);

    %IG

    xG = zeros(HS_size,ndim);

    %IB

    xB = zeros(HS_size,ndim);

end

C_Func = 0;

tic

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Generating Initial Solution Vector

for i = 1:HS_size,

    for j = 1:ndim,

        x(j) = range(j,1) + (range(j,2) - range(j,1)) * rand;

    end

    x = fix(sort(x));

    HM(i, :) = x;

end %% for i

    %evalua x en la funcion objetivo

     %[HMbest, fitBestR] = fitnessIMG(I, HS_size, Lmax, level, HM, probR);

%     C_Func = length(HMbest);

    

    %evalua x en la funcion objetivo

    x = fix(sort(x));

    %evalua x en la funcion objetivo

    %[fbest, fitBestR] = fitnessIMG(I, 1, Lmax, level, x, probR);

    fbest = Kapur(1,level,x,probR);

    C_Func = C_Func + 1;

    

    

    % Find the best in the HS solution vector   

    [HStemp, ii] = sort(HMbest, ’descend’); %Maximiza

    HMbest = HMbest(ii);

    HM = HM(ii,:);

    

    % Updating the current solution if better

    if fbest > HMbest(HS_size), %maximiza

        HM(HS_size, :) = x;

        HMbest(HS_size) = fbest;

    end

    solution = x;   % Record the solution

    %Obtiene los mejores valores de cada attempt y los alamacena

    [mm,ii] = max(HMbest); %maximiza

    Fit_bests(count) = mm; %Mejores Fitness

    HS_elem(count,:) = HM(ii,1:ndim-1); %Mejores Elementos de HM

    HS_bestit = HM(ii,1:ndim-1); %Guarda el mejor HS

    HS_bestF = mm; %Guarda el mejor fitness

    

   % Output the results  to screen

   str=strcat(’Best estimates: =’,num2str(HS_bestit));

   str=strcat(str,’  fmin=’); str=strcat(str,num2str(HS_bestF));

   str=strcat(str,’  Iteration=’); str=strcat(str,num2str(count));

   disp(str);   

   

   %Save the best values that will be chek in the stop criterion

   if count == 1  || HS_bestF > HS_ant

        HS_ant = HS_bestF;

        cc = 0;

   elseif HS_bestF == HS_ant

        cc = cc + 1;

   end

   

   if cc > (MaxAttempt * 0.10)

       break;

   end

   

end %% for count (harmony search)

toc

%plot fitness

plot(Fit_bests)

%Prepare results to be show

 gBestR = sort(HS_bestit);

    Iout = imageGRAY(I,gBestR);

    Iout2 = mat2gray(Iout);

    %Show results

    intensity = gBestR(1:ndim-1)     

    STDR =  std(Fit_bests)      %Standar deviation of fitness       

    MEANR = mean(Fit_bests)     %Mean of fitness

    PSNRV = PSNR(I, Iout)       %PSNR between original image I and the segmented image Iout

    Fit_bests(count)            %Best fitness

    %Show results on images

    figure

    subplot(121)

       imshow(I);title(’原圖’)

    subplot(122)

    imshow(Iout);title(’分割圖’)

 

    

    %Plot the threshold values over the histogram

    figure 

    plot(probR)

    hold on

    vmax = max(probR);

    for i = 1:ndim-1

        line([intensity(i), intensity(i)],[0 vmax],[1 1],’Color’,’r’,’Marker’,’.’,’LineStyle’,’-’)

        %plot(lineas(i,:))

        hold on

    end

     hold off

3 運(yùn)行結(jié)果

4 參考文獻(xiàn)

[1] Oliva D ,  Cuevas E ,  Pajares G , et al. Multilevel Thresholding Segmentation Based on Harmony Search Optimization[J]. Journal of Applied Mathematics, 2013, 2013:1-12.

博主簡介：擅長智能優(yōu)化算法、神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)、信號(hào)處理、元胞自動(dòng)機(jī)、圖像處理、路徑規(guī)劃、無人機(jī)等多種領(lǐng)域的Matlab仿真，相關(guān)matlab代碼問題可私信交流。

部分理論引用網(wǎng)絡(luò)文獻(xiàn)，若有侵權(quán)聯(lián)系博主刪除。