WTMM基于WTMM的图像多重分形谱matlab仿真

Posted 2022-06-17 fpga和matlab

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了WTMM基于WTMM的图像多重分形谱matlab仿真相关的知识，希望对你有一定的参考价值。

1.软件版本

matlab2013b

2.系统程序

提出一种新的多重分形图象去噪算法,讨论了基于小波极大模的多重分形谱估计算法.在此基础上推导了图像取噪声的谱移位算子.该方法没有对噪声的类型提出任何假设条件,而是通过定义一个变换算子对每一点的 Hausdorf指数进行处理,使处理后的图象的Hausdorf指数接近于2,从而取得最佳效果.实验结果表明,该方法在去除噪声的同时可很好地保留了原始图像的纹理信息.

算法流程图如图所示：

部分程序如下：

clc;
clear;
close all;
warning off;
addpath 'FUNC\\'

sel = 0;%1：计算多重分形谱；0：分类


%%
%读取三个图片
I1s = imread('图像\\D10.gif');
I1  = imresize(I1s,[512,512]);
I2s = imread('图像\\D05.gif');
I2  = imresize(I2s,[512,512]);
figure;
subplot(121);imshow(I1);title('图片1');
subplot(122);imshow(I2);title('图片2');

%%
%分成多个子图片
[M,N]      = size(I1);
sub_length = 64;
n1         = floor(M/sub_length);
n2         = floor(N/sub_length);
DF         = zeros(1,n1*n2);
 
for num_image = 1:2
    k=1;
    for i=1:n1
        for j=1:n2
            x_start   = (i-1)*sub_length+1;
            x_end     =  i*sub_length;
            y_start   = (j-1)*sub_length+1;
            y_end     =  j*sub_length;
            if num_image == 1
               sub_image = I1(x_start:x_end,y_start:y_end);
               imwrite(sub_image,['save_images\\1\\I_',num2str(k),'.jpg']);
            end
            if num_image == 2
               sub_image = I2(x_start:x_end,y_start:y_end);
               imwrite(sub_image,['save_images\\2\\I_',num2str(k),'.jpg']);
            end                  
            %保存分割的多个子图片
            
            k = k + 1;
        end
    end
end

%%
%对保存的多张图片读取并调用WTMM方法求图像的多重分形谱，对得到的结果保存其特征值
if sel == 1
   k = 1;

    for i=1:2*n1*n2;
        
        if i<=n1*n2
           k      = i;
           folder = 'save_images\\1\\';
           lists  = dir('save_images\\1\\*.jpg');        
        end
        
        if i<=2*n1*n2 & i>n1*n2
           k      = i - n1*n2;
           folder = 'save_images\\2\\';
           lists  = dir('save_images\\2\\*.jpg');        
        end        
              
        
        i
        %read an image
        I                               = imread(fullfile(folder,lists(k).name));
        %调用分形函数
        [qt,rt,ft,fft,Dt,feature_data]  = func_Wavelet_multifractal(I);
        qi                            = qt;
        ri                            = rt;    
        fi                            = ft;   
        ffi                           = fft;   
        Di                            = Dt;
        Featurei                      = feature_data;
    end

    save result.mat q r f ff D Feature

    K = 120;
    figure;
    plot(rK,fK,'b-o'); 
    xlabel('奇异指数a');
    ylabel('多重分行谱f(a)') 
    grid on;
    
    figure;
    plot(qK+2,DK,'b-o'); 
    xlabel('q');
    ylabel('D(q)') 
    grid on;
    
    figure
    plot(qK,rK,'b-o');
    title('ｑ和阿尔法a'); 
    xlabel('权重因子q');
    ylabel('奇异指数a');
    grid on;
    
    figure;
    plot(qK,fK,'b-o');
    title('ｑ和f(a) '); 
    xlabel('权重因子q');
    ylabel('多重分行谱f(a)'); 
    grid on;
end 




%%
%调用分类器对特征参数进行分类
if sel == 0
    load result.mat %q r f ff Feature

    K = 120;
    figure;
    plot(rK,fK,'b-o'); 
    xlabel('奇异指数a');
    ylabel('多重分行谱f(a)') 
    grid on;
    
    figure;
    plot(qK+2,DK,'b-o'); 
    xlabel('q');
    ylabel('D(q)') 
    grid on;
    
    figure
    plot(qK,rK,'b-o');
    title('ｑ和阿尔法a'); 
    xlabel('权重因子q');
    ylabel('奇异指数a');
    grid on;
    
    figure;
    plot(qK,fK,'b-o');
    title('ｑ和f(a) '); 
    xlabel('权重因子q');
    ylabel('多重分行谱f(a)'); 
    grid on;   
    
    
    
    
    for i = 1:length(Feature)
        P(i) =  Featurei(3);
    end 

    T =  [1*ones(1,length(Feature)/2),2*ones(1,length(Feature)/2)];

    t1                      = clock;                              %计时开始
    net                     = fitnet(65);
    net.trainParam.epochs   = 1000;                               %设置训练次数
    net.trainParam.goal     = 0.0001;                             %设置性能函数
    net.trainParam.show     = 1;                                  %每10显示
    net.trainParam.Ir       = 0.005;                              %设置学习速率
    net                     = train(net,P,T);                     %训练BP网络
    datat                   = etime(clock,t1);
    Nets                    = net;
    view(Nets);
    figure;
    plot(P,'b-*');

    y = sim(net,P);  

    figure;
    stem(y,'r-*');hold on
    plot(T,'b-','LineWidth',2);hold on
    legend('预测数据','实际数据');
    title('输出1为第一类，输出2为第二类（即可对比实际的健康部分和肿瘤部分）');
    
    
    disp('预测正确率');
    error = 0;
    for i = 1:length(y)
        if i <= length(y)/2 
           if y(i) > 1.5
              error = error + 1;
           end
        else
           if y(i) < 1.5
              error = error + 1;
           end          
        end
    end
    1-error/length(y)
end