www.gusucode.com > wavelet工具箱matlab源码程序 > wavelet/wavelet/private/mwptdec.m
function [dec,coefs] = mwptdec(dirDec,x,lev,varargin) %MWPTDEC Multisignal wavelet packet 1-D decomposition. % M. Misiti, Y. Misiti, G. Oppenheim, J.M. Poggi 18-Feb-2011 % Last Revision: 28-Aug-2011. % Copyright 1995-2015 The MathWorks, Inc. % Check arguments. narginchk(4,11) nargoutchk(0,3); if errargt(mfilename,lev,'int') error(message('Wavelet:FunctionArgVal:Invalid_Input')); end if ischar(varargin{1}) wname = varargin{1}; [LoD,HiD,LoR,HiR] = wfilters(wname); varargin(1) = []; else wname = ''; [LoD,HiD,LoR,HiR] = deal(varargin{1:4}); varargin(1:4) = []; end % Type of decomposition. if ~isempty(dirDec) dirDec = dirDec(1); dirDec = lower(dirDec); else dirDec = 'r'; end dataSize = size(x); switch dirDec case 'c' , x = x'; dirCAT = 2; otherwise , dirDec = 'r'; dirCAT = 1; end x = double(x); % Default: Shift and Extension. dwtATTR = dwtmode('get'); shift = dwtATTR.shift1D; dwtEXTM = dwtATTR.extMode; % Check arguments for Extension and Shift. nbin = length(varargin); for k = 1:2:nbin-1 switch varargin{k} case 'mode' , dwtEXTM = varargin{k+1}; case 'shift' , shift = mod(varargin{k+1},2); end end perFLAG = isequal(dwtEXTM,'per'); % Initialization. sx = zeros(lev+2,2); sx(end,:) = size(x); first = 2 - rem(shift,2); % Full decomposition cfs = {x}; idxLST = 1; for j = 1:((2^lev)-1) [a,d] = dwtLOC(cfs{idxLST},LoD,HiD,dwtEXTM,perFLAG,first); cfs = [cfs,a,d]; cfs(idxLST) = []; k = log2(j+1); if isequal(k,fix(k)) , sx(lev+2-k,:) = size(a); end end sx(1,:) = sx(2,:); Filters = struct('LoD',LoD,'HiD',HiD,'LoR',LoR,'HiR',HiR); dec = struct(... 'dirDec',dirDec,'level',lev, ... 'wname',wname,'dwtFilters',Filters, ... 'dwtEXTM',dwtEXTM,'dwtShift',shift, ... 'dataSize',dataSize,'sx',sx); dec.cfs = cfs; coefs = cat(dirCAT,dec.cfs{:}); %------------------------------------------------------------------------ function [a,d] = dwtLOC(x,LoD,HiD,dwtEXTM,perFLAG,first) % Compute sizes. lf = length(LoD); lx = size(x,2); % Extend, Decompose & Extract coefficients. dCol = lf-1; if ~perFLAG lenEXT = lf-1; lenKEPT = lx+lf-1; else lenEXT = lf/2; lenKEPT = 2*ceil(lx/2); end idxCOL = (first + dCol : 2 : lenKEPT + dCol); y = wextend('addcol',dwtEXTM,x,lenEXT); a = conv2(y,LoD,'full'); a = a(:,idxCOL); d = conv2(y,HiD,'full'); d = d(:,idxCOL); %--------------------------------------------------------------------