swt1d_decomp.m

编程实现基于小波变换尺度相关性去噪算法

function [a, d, hj, gj] = swt1d_decomp(s,n,h,g)
%UNTITLED1 Summary of this function goes here
%  Detailed explanation goes here

% Preserve initial size.
% figure(1)    %显示原始图像
% imshow(x)
l = length(s);
pow = 2^n;  
if any(rem(l,pow))
    sOK = ceil(l/pow)*pow;
    oriStr = [int2str(l)];
    sugStr = [int2str(sOK)];
    msg = strvcat(...
        ['The level of decomposition ' int2str(n)],...
        ['and the length of the signal ' oriStr],...
        'are not compatible.',...
        ['Suggested length: ' sugStr],...
        ' ', ...
        ['2^Level has to divide the length of the signal.'] ...
            );
    errargt(mfilename,msg,'msg');
    varargout = {[] };
   return
end        

% 设置 DWT_Mode 为 'per',即周期延拓
modeDWT = 'per';

% 一系列参数的初始化操作.
evenoddVal = 0;  %用于控制滤波器的上抽样
evenLEN    = 2;  %用于控制滤波器的上抽样
a = zeros(l,n); %初始化逼近矩阵
d = zeros(l,n); %初始化细节信息矩阵 
x = s;
for k=1:n %开始逐级计算分解
    
    % 延拓
    lf = length(h);
    x  = wextend('1D',modeDWT,x,round(lf/2));
    
    if(mod(lf,2)==0)
        if(k==1)
            shift = -1;
        else
            shift = ceil(2^(k-2)-0.5)-1; %shift用于抵消 平移的误差
        end
    else
        shift = 0;
    end

    
%     size(x)
    % 分解. 得逼近矩阵a
    a(:,k) = wkeep(wconv('1D',x,h),l,floor(lf-shift)); % x*h  
      
    % 分解. 得细节信息矩阵d
    d(:,k)=wkeep(wconv('1D',x,g),l,floor(lf-shift));  % x*g
    
    hj{k} = h;
    gj{k} = g;

    % 对滤波器进行上抽样
    h = dyadup(h,evenoddVal,evenLEN);
    g = dyadup(g,evenoddVal,evenLEN);

    % 更新x的值 进行下一级分解.
    x = a(:,k);
end