linbwt.m

matlab程序

function X = linbwt(X, Level, Dim)
%LINBWT  Linear binary wavelet transform.
%   Y = LINBWT(X,L) computes the L level decomposition of a signal X
%   using a linear binary wavelet transform [1].  The signal length
%   must be divisible by 2^L.  If X is a matrix, the transform is
%   applied to each column.
%
%   LINBWT(X,-L) is the inverse transform, reversing L levels.
%
%   LINBWT(X,L,DIM) applies the transform across the dimension DIM.
%
%   Example:
%   Y = linbwt(X,3);   % Perform 3 levels of decomposition on X
%   R = linbwt(Y,-3);  % Recover X from Y
%
%   Reference:
%   [1] M. Swanson and A. Tewfik.  ``A Binary Wavelet Decomposition of
%       Binary Images.''  IEEE Transactions on Image Processing, Vol. 5,
%       No. 12, Dec. 1996.
%
%   See also LINBWT2, MORPHBWT2.

% Pascal Getreuer 2005

if nargin < 2, error('Not enough input arguments.'); end
if nargin < 3, Dim = min(find(size(X) ~= 1)); end

XSize = size(X);
N = XSize(Dim);
Perm = [Dim:max(length(XSize),Dim) 1:Dim-1];
X = reshape(permute(X,Perm),N,prod(XSize)/N);

if rem(N,pow2(abs(Level))), error('Invalid input size.'); end

if Level > 0
   for k = 1:Level
      N = size(X,1).*pow2(1-k);
      s = filter([1,1,1],1,X([2:N,1],:),[1,1;1,0]*X([1,N],:),1);
      d = filter([1,1],1,X([2:N,1],:),X(1,:),1);
      X(1:N,:) = bitand([s(1:2:N,:);d(1:2:N,:)],1);
   end
elseif Level < 0
   for k = Level:-1
      N = size(X,1)*pow2(k+1);
      s(2:2:N,:) = X(1:N/2,:);
      d(2:2:N,:) = X(N/2+1:N,:);
      X(1:N,:) = bitand(filter([1,1],1,s([min(4,N):N,1:min(3,N-1)],:),s(min(3,N-1),:),1) ...
         + filter([1,1,1],1,d([min(4,N):N,1:min(3,N-1)],:),[1,1;1,0]*d([min(3,N-1),2],:),1),1);
   end
end

X = ipermute(reshape(X,XSize(Perm)),Perm);