當調用wavefast函數發現matlab沒有該函數，通過查閱找到了小波變換的m檔案，儲存成函數檔案後就可以直接調用。
matlab源碼：
(1) wave2gray.m
function w = wave2gray(c, s, scale, border)
%WAVE2GRAY Display wavelet decomposition coefficients.
% W = WAVE2GRAY(C, S, SCALE, BORDER) displays and returns a
% wavelet coefficient image.
% EXAMPLES:
% wave2gray(c, s); Display w/defaults.
% foo = wave2gray(c, s); Display and return.
% foo = wave2gray(c, s, 4); Magnify the details.
% foo = wave2gray(c, s, -4); Magnify absolute values.
% foo = wave2gray(c, s, 1, 'append'); Keep border values.
%
% INPUTS/OUTPUTS:
% [C, S] is a wavelet decomposition vector and bookkeeping
% matrix.
%
% SCALE Detail coefficient scaling
%--------------------------------------------------------------------
% 0 or 1 Maximum range (default)
% 2, 3... Magnify default by the scale factor
% -1, -2... Magnify absolute values by abs(scale)
%
% BORDER Border between wavelet decompositions
%--------------------------------------------------------------------
% 'absorb' Border replaces image (default)
% 'append' Border increases width of image
%
% Image W: ------- ------ -------------- -------------------
% | | | |
% | a(n) | h(n) | |
% | | | |
% ------- ------ h(n-1) |
% | | | |
% | v(n) | d(n) | | h(n-2)
% | | | |
% ------- ------ --------------
% | | |
% | v(n-1) | d(n-1) |
% | | |
% -------------- -------------- -------------------
% | |
% | v(n-2) | d(n-2)
% | |
%
% Here, n denotes the decomposition step scale and a, h, v, d are
% approximation, horizontal, vertical, and diagonal detail
% coefficients, respectively.
% Check input arguments for reasonableness.
error(nargchk(, , nargin));
if (ndims(c) ~= ) | (size(c, ) ~= )
error('C must be a row vector.'); end
if (ndims(s) ~= ) | ~isreal(s) | ~isnumeric(s) | (size(s, ) ~= )
error('S must be a real, numeric two-column array.'); end
elements = prod(s, );
if (length(c) < elements(end)) | ...
~(elements() + * sum(elements(:end - )) >= elements(end))
error(['[C S] must be a standard wavelet ' ...
'decomposition structure.']);
end
if (nargin > ) & (~isreal(scale) | ~isnumeric(scale))
error('SCALE must be a real, numeric scalar.');
end
if (nargin > ) & (~ischar(border))
error('BORDER must be character string.');
end
if nargin ==
scale =; % Default scale.
end
if nargin <
border = 'absorb'; % Default border.
end
% Scale coefficients and determine pad fill.
absflag = scale < ;
scale = abs(scale);
if scale ==
scale = ;
end
[cd, w] = wavecut('a', c, s); w = mat2gray(w);
cdx = max(abs(cd(:))) / scale;
if absflag
cd = mat2gray(abs(cd), [, cdx]); fill = ;
else
cd = mat2gray(cd, [-cdx, cdx]); fill = ;
end
% Build gray image one decomposition at a time.
for i = size(s, ) - :-:
ws = size(w);
h = wavecopy('h', cd, s, i);
pad = ws - size(h); frontporch = round(pad / );
h = padarray(h, frontporch, fill, 'pre');
h = padarray(h, pad - frontporch, fill, 'post');
v = wavecopy('v', cd, s, i);
pad = ws - size(v); frontporch = round(pad / );
v = padarray(v, frontporch, fill, 'pre');
v = padarray(v, pad - frontporch, fill, 'post');
d = wavecopy('d', cd, s, i);
pad = ws - size(d); frontporch = round(pad / );
d = padarray(d, frontporch, fill, 'pre');
d = padarray(d, pad - frontporch, fill, 'post');
% Add 1 pixel white border.
switch lower(border)
case 'append'
w = padarray(w, [ ], , 'post');
h = padarray(h, [ ], , 'post');
v = padarray(v, [ ], , 'post');
case 'absorb'
w(:, end) = ; w(end, :) = ;
h(end, :) = ; v(:, end) = ;
otherwise
error('Unrecognized BORDER parameter.');
end
w = [w h; v d]; % Concatenate coefs.
end
if nargout ==
imshow(w); % Display result.
end
(2) waveback.m
function [varargout] = waveback(c, s, varargin)
%WAVEBACK Performs a multi-level two-dimensional inverse FWT.
% [VARARGOUT] = WAVEBACK(C, S, VARARGIN) computes a 2D N-level
% partial or complete wavelet reconstruction of decomposition
% structure [C, S].
%
% SYNTAX:
% Y = WAVEBACK(C, S, 'WNAME'); Output inverse FWT matrix Y
% Y = WAVEBACK(C, S, LR, HR); using lowpass and highpass
% reconstruction filters (LR and
% HR) or filters obtained by
% calling WAVEFILTER with 'WNAME'.
%
% [NC, NS] = WAVEBACK(C, S, 'WNAME', N); Output new wavelet
% [NC, NS] = WAVEBACK(C, S, LR, HR, N); decomposition structure
% [NC, NS] after N step
% reconstruction.
%
% See also WAVEFAST and WAVEFILTER.
% Copyright - R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall,
% $Revision: $ $Date: // :: $
% Check the input and output arguments for reasonableness.
error(nargchk(, , nargin));
error(nargchk(, , nargout));
if (ndims(c) ~= ) | (size(c, ) ~= )
error('C must be a row vector.');
end
if (ndims(s) ~= ) | ~isreal(s) | ~isnumeric(s) | (size(s,) ~= )
error('S must be a real, numeric two-column array.');
end
elements = prod(s, );
if (length(c) < elements(end)) | ...
~(elements() + * sum(elements(:end - )) >= elements(end))
error(['[C S] must be a standard wavelet ' ...
'decomposition structure.']);
end
% Maximum levels in [C, S].
nmax = size(s, ) - ;
% Get third input parameter and init check flags.
wname = varargin{}; filterchk = ; nchk = ;
switch nargin
case
if ischar(wname)
[lp, hp] = wavefilter(wname, 'r'); n = nmax;
else
error('Undefined filter.');
end
if nargout ~=
error('Wrong number of output arguments.');
end
case
if ischar(wname)
[lp, hp] = wavefilter(wname, 'r');
n = varargin{}; nchk = ;
else
lp = varargin{}; hp = varargin{};
filterchk = ; n = nmax;
if nargout ~=
error('Wrong number of output arguments.');
end
end
case
lp = varargin{}; hp = varargin{}; filterchk = ;
n = varargin{}; nchk = ;
otherwise
error('Improper number of input arguments.');
end
fl = length(lp);
if filterchk % Check filters.
if (ndims(lp) ~= ) | ~isreal(lp) | ~isnumeric(lp) ...
| (ndims(hp) ~= ) | ~isreal(hp) | ~isnumeric(hp) ...
| (fl ~= length(hp)) | rem(fl, ) ~=
error(['LP and HP must be even and equal length real, ' ...
'numeric filter vectors.']);
end
end
if nchk & (~isnumeric(n) | ~isreal(n)) % Check scale N.
error('N must be a real numeric.');
end
if (n > nmax) | (n < )
error('Invalid number (N) of reconstructions requested.');
end
if (n ~= nmax) & (nargout ~= )
error('Not enough output arguments.');
end
nc = c; ns = s; nnmax = nmax; % Init decomposition.
for i = :n
% Compute a new approximation.
a = symconvup(wavecopy('a', nc, ns), lp, lp, fl, ns(, :)) + ...
symconvup(wavecopy('h', nc, ns, nnmax), ...
hp, lp, fl, ns(, :)) + ...
symconvup(wavecopy('v', nc, ns, nnmax), ...
lp, hp, fl, ns(, :)) + ...
symconvup(wavecopy('d', nc, ns, nnmax), ...
hp, hp, fl, ns(, :));
% Update decomposition.
nc = nc( * prod(ns(, :)) + :end); nc = [a(:)' nc];
ns = ns(3:end, :); ns = [ns(1, :); ns];
nnmax = size(ns, 1) - 2;
end
% For complete reconstructions, reformat output as 2-D.
if nargout == 1
a = nc; nc = repmat(0, ns(1, :)); nc(:) = a;
end
varargout{1} = nc;
if nargout == 2
varargout{2} = ns;
end
%-------------------------------------------------------------------%
function z = symconvup(x, f1, f2, fln, keep)
% Upsample rows and convolve columns with f1; upsample columns and
% convolve rows with f2; then extract center assuming symmetrical
% extension.
y = zeros([2 1] .* size(x)); y(1:2:end, :) = x;
y = conv2(y, f1');
z = zeros([ ] .* size(y)); z(:, ::end) = y;
z = conv2(z, f2);
z = z(fln - :fln + keep() - , fln - :fln + keep() - );
(3) wavecopy.m
function y = wavecopy(type, c, s, n)
%WAVECOPY Fetches coefficients of a wavelet decomposition structure.
% Y = WAVECOPY(TYPE, C, S, N) returns a coefficient array based on
% TYPE and N.
%
% INPUTS:
% TYPE Coefficient category
% -------------------------------------
% 'a' Approximation coefficients
% 'h' Horizontal details
% 'v' Vertical details
% 'd' Diagonal details
%
% [C, S] is a wavelet data structure.
% N specifies a decomposition level (ignored if TYPE = 'a').
%
% See also WAVEWORK, WAVECUT, and WAVEPASTE.
% Copyright 2002-2004 R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall, 2004
% $Revision: 1.4 $ $Date: 2003/10/13 01:20:41 $
error(nargchk(, , nargin));
if nargin ==
y = wavework('copy', type, c, s, n);
else
y = wavework('copy', type, c, s);
end
(4) wavecut.m
function [nc, y] = wavecut(type, c, s, n)
%WAVECUT Zeroes coefficients in a wavelet decomposition structure.
% [NC, Y] = WAVECUT(TYPE, C, S, N) returns a new decomposition
% vector whose detail or approximation coefficients (based on TYPE
% and N) have been zeroed. The coefficients that were zeroed are
% returned in Y.
%
% INPUTS:
% TYPE Coefficient category
% -------------------------------------
% 'a' Approximation coefficients
% 'h' Horizontal details
% 'v' Vertical details
% 'd' Diagonal details
%
% [C, S] is a wavelet data structure.
% N specifies a decomposition level (ignored if TYPE = 'a').
%
% See also WAVEWORK, WAVECOPY, and WAVEPASTE.
% Copyright 2002-2004 R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall, 2004
% $Revision: 1.4 $ $Date: 2003/10/13 01:20:09 $
error(nargchk(, , nargin));
if nargin ==
[nc, y] = wavework('cut', type, c, s, n);
else
[nc, y] = wavework('cut', type, c, s);
end
(5) wavefast.m
function [c, s] = wavefast(x, n, varargin)
%WAVEFAST Perform multi-level -dimensional fast wavelet transform.
% [C, L] = WAVEFAST(X, N, LP, HP) performs a 2D N-level FWT of
% image (or matrix) X with respect to decomposition filters LP and
% HP.
%
% [C, L] = WAVEFAST(X, N, WNAME) performs the same operation but
% fetches filters LP and HP for wavelet WNAME using WAVEFILTER.
%
% Scale parameter N must be less than or equal to log2 of the
% maximum image dimension. Filters LP and HP must be even. To
% reduce border distortion, X is symmetrically extended. That is,
% if X = [c1 c2 c3 ... cn] (in 1D), then its symmetric extension
% would be [... c3 c2 c1 c1 c2 c3 ... cn cn cn- cn- ...].
%
% OUTPUTS:
% Matrix C is a coefficient decomposition vector:
%
% C = [ a(n) h(n) v(n) d(n) h(n-) ... v() d() ]
%
% where a, h, v, and d are columnwise vectors containing
% approximation, horizontal, vertical, and diagonal coefficient
% matrices, respectively. C has 3n + sections where n is the
% number of wavelet decompositions.
%
% Matrix S is an (n+) x bookkeeping matrix:
%
% S = [ sa(n, :); sd(n, :); sd(n-, :); ... ; sd(, :); sx ]
%
% where sa and sd are approximation and detail size entries.
%
% See also WAVEBACK and WAVEFILTER.
% Copyright - R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall,
% $Revision: $ $Date: // :: $
% Check the input arguments for reasonableness.
error(nargchk(, , nargin));
if nargin ==
if ischar(varargin{})
[lp, hp] = wavefilter(varargin{}, 'd');
else
error('Missing wavelet name.');
end
else
lp = varargin{}; hp = varargin{};
end
fl = length(lp); sx = size(x);
if (ndims(x) ~= ) | (min(sx) < ) | ~isreal(x) | ~isnumeric(x)
error('X must be a real, numeric matrix.');
end
if (ndims(lp) ~= ) | ~isreal(lp) | ~isnumeric(lp) ...
| (ndims(hp) ~= ) | ~isreal(hp) | ~isnumeric(hp) ...
| (fl ~= length(hp)) | rem(fl, ) ~=
error(['LP and HP must be even and equal length real, ' ...
'numeric filter vectors.']);
end
if ~isreal(n) | ~isnumeric(n) | (n < ) | (n > log2(max(sx)))
error(['N must be a real scalar between 1 and ' ...
'log2(max(size((X))).']);
end
% Init the starting output data structures and initial approximation.
c = []; s = sx; app = double(x);
% For each decomposition ...
for i = :n
% Extend the approximation symmetrically.
[app, keep] = symextend(app, fl);
% Convolve rows with HP and downsample. Then convolve columns
% with HP and LP to get the diagonal and vertical coefficients.
rows = symconv(app, hp, 'row', fl, keep);
coefs = symconv(rows, hp, 'col', fl, keep);
c = [coefs(:)' c]; s = [size(coefs); s];
coefs = symconv(rows, lp, 'col', fl, keep);
c = [coefs(:)' c];
% Convolve rows with LP and downsample. Then convolve columns
% with HP and LP to get the horizontal and next approximation
% coeffcients.
rows = symconv(app, lp, 'row', fl, keep);
coefs = symconv(rows, hp, 'col', fl, keep);
c = [coefs(:)' c];
app = symconv(rows, lp, 'col', fl, keep);
end
% Append final approximation structures.
c = [app(:)' c]; s = [size(app); s];
%-------------------------------------------------------------------%
function [y, keep] = symextend(x, fl)
% Compute the number of coefficients to keep after convolution
% and downsampling. Then extend x in both dimensions.
keep = floor((fl + size(x) - ) / );
y = padarray(x, [(fl - ) (fl - )], 'symmetric', 'both');
%-------------------------------------------------------------------%
function y = symconv(x, h, type, fl, keep)
% Convolve the rows or columns of x with h, downsample,
% and extract the center section since symmetrically extended.
if strcmp(type, 'row')
y = conv2(x, h);
y = y(:, ::end);
y = y(:, fl / + :fl / + keep());
else
y = conv2(x, h');
y = y(1:2:end, :);
y = y(fl / 2 + 1:fl / 2 + keep(1), :);
end
(6) wavefilter.m
function [varargout]=wavefilter(wname,type)
error(nargchk(,,nargin));
if(nargin== & nargout ~=) | (nargin== & nargout~=)
error('Invalid number of output arguments.');
end
if nargin== & ~ischar(wname)
error('WNAME must be a string.');
end
if nargin== & ~ischar(type)
error('TYPE must be a string.');
end
switch lower(wname)
case {'haar','db1'}
ld=[ ]/sqrt(); hd=[- ]/sqrt();
lr=ld; hr=-hd;
case 'db4'
ld=[- - ...
- ];
t=[:];
hd=ld; hd(end:-:)=cos(pi*t).*ld;
lr=ld; lr(end:-:)=ld;
hr=cos(pi*t).*ld;
case 'sym4'
ld=[- - ...
- - ];
t=[:];
hd=ld; hd(end:-:)=cos(pi*t).*ld;
lr=ld; lr(end:-:)=ld;
hr=cos(pi*t).*ld;
case 'bior6.8'
ld=[ - - ...
- - ...
- - ...
- - ...
];
hd=[ - - ...
- ...
- - ];
t=[:];
lr=cos(pi*(t+)).*hd;
hr=cos(pi*t).*ld;
case 'jpeg9.7'
ld=[ - - ...
- - ...
];
hd=[ - - ...
- ];
t=[:];
lr=cos(pi*(t+)).*hd;
hr=cos(pi*t).*ld;
otherwise
error('Unrecongizable wavelet name (WNAME).');
end
if(nargin==)
varargout(:)={ld,hd,lr,hr};
else
switch lower(type())
case 'd'
varargout={ld,hd};
case 'r'
varargout={lr,hr};
otherwise
error('Unrecongizable filter TYPE.');
end
end
(7) wavepaste.m
function nc = wavepaste(type, c, s, n, x)
%WAVEPASTE Puts coefficients in a wavelet decomposition structure.
% NC = WAVEPASTE(TYPE, C, S, N, X) returns the new decomposition
% structure after pasting X into it based on TYPE and N.
%
% INPUTS:
% TYPE Coefficient category
% -------------------------------------
% 'a' Approximation coefficients
% 'h' Horizontal details
% 'v' Vertical details
% 'd' Diagonal details
%
% [C, S] is a wavelet data structure.
% N specifies a decomposition level (Ignored if TYPE = 'a').
% X is a two-dimensional approximation or detail coefficient
% matrix whose dimensions are appropriate for decomposition
% level N.
%
% See also WAVEWORK, WAVECUT, and WAVECOPY.
% Copyright - R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall,
% $Revision: $ $Date: // :: $
error(nargchk(, , nargin))
nc = wavework('paste', type, c, s, n, x);
(8) wavework.m
function [varargout] = wavework(opcode, type, c, s, n, x)
%WAVEWORK is used to edit wavelet decomposition structures.
% [VARARGOUT] = WAVEWORK(OPCODE, TYPE, C, S, N, X) gets the
% coefficients specified by TYPE and N for access or modification
% based on OPCODE.
%
% INPUTS:
% OPCODE Operation to perform
% --------------------------------------------------------------
% 'copy' [varargout] = Y = requested (via TYPE and N)
% coefficient matrix
% 'cut' [varargout] = [NC, Y] = New decomposition vector
% (with requested coefficient matrix zeroed) AND
% requested coefficient matrix
% 'paste' [varargout] = [NC] = new decomposition vector with
% coefficient matrix replaced by X
%
% TYPE Coefficient category
% --------------------------------------------
% 'a' Approximation coefficients
% 'h' Horizontal details
% 'v' Vertical details
% 'd' Diagonal details
%
% [C, S] is a wavelet toolbox decomposition structure.
% N is a decomposition level (Ignored if TYPE = 'a').
% X is a two-dimensional coefficient matrix for pasting.
%
% See also WAVECUT, WAVECOPY, and WAVEPASTE.
% Copyright - R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall,
% $Revision: $ $Date: // :: $
error(nargchk(, , nargin));
if (ndims(c) ~= ) | (size(c, ) ~= )
error('C must be a row vector.');
end
if (ndims(s) ~= ) | ~isreal(s) | ~isnumeric(s) | (size(s, ) ~= )
error('S must be a real, numeric two-column array.');
end
elements = prod(s, ); % Coefficient matrix elements.
if (length(c) < elements(end)) | ...
~(elements() + * sum(elements(:end - 1)) >= elements(end))
error(['[C S] must form a standard wavelet decomposition ' ...
'structure.']);
end
if strcmp(lower(opcode(:)), 'pas') & nargin <
error('Not enough input arguments.');
end
if nargin <
n = ; % Default level is
end
nmax = size(s, ) - ; % Maximum levels in [C, S].
aflag = (lower(type()) == 'a');
if ~aflag & (n > nmax)
error('N exceeds the decompositions in [C, S].');
end
switch lower(type()) % Make pointers into C.
case 'a'
nindex = ;
start = ; stop = elements(); ntst = nmax;
case {'h', 'v', 'd'}
switch type
case 'h', offset = ; % Offset to details.
case 'v', offset = ;
case 'd', offset = ;
end
nindex = size(s, ) - n; % Index to detail info.
start = elements() + * sum(elements(:nmax - n + )) + ...
offset * elements(nindex) + ;
stop = start + elements(nindex) - ;
ntst = n;
otherwise
error('TYPE must begin with "a", "h", "v", or "d".');
end
switch lower(opcode) % Do requested action.
case {'copy', 'cut'}
y = repmat(, s(nindex, :));
y(:) = c(start:stop); nc = c;
if strcmp(lower(opcode(:)), 'cut')
nc(start: stop) = ; varargout = {nc, y};
else
varargout = {y};
end
case 'paste'
if prod(size(x)) ~= elements(end - ntst)
error('X is not sized for the requested paste.');
else
nc = c; nc(start:stop) = x(:); varargout = {nc};
end
otherwise
error('Unrecognized OPCODE.');
end
(9) wavezero.m
function [nc, g8] = wavezero(c, s, l, wname)
%WAVEZERO Zeroes wavelet transform detail coefficients.
% [NC, G8] = WAVEZERO(C, S, L, WNAME) zeroes the level L detail
% coefficients in wavelet decomposition structure [C, S] and
% computes the resulting inverse transform with respect to WNAME
% wavelets.
% Copyright 2002-2004 R. C. Gonzalez, R. E. Woods, & S. L. Eddins
% Digital Image Processing Using MATLAB, Prentice-Hall, 2004
% $Revision: 1.4 $ $Date: 2003/10/13 01:31:35 $
[nc, foo] = wavecut('h', c, s, l);
[nc, foo] = wavecut('v', nc, s, l);
[nc, foo] = wavecut('d', nc, s, l);
i = waveback(nc, s, wname);
g8 = im2uint8(mat2gray(i));
figure; imshow(g8);
測試程式：
clear
f=imread('cameraman.tif');
[c,s]=wavefast(f,,'sym4');%快速小波變換
figure;
wave2gray(c,s,-);
[nc,y]=wavecut('a',c,s);%近似系數置
figure;
wave2gray(nc,s,-);
edges=abs(waveback(nc,s,'sym4'));%邊緣重構
figure;
imshow(mat2gray(edges));
結果：
![POPUP_TO_DECIDE_LIST 彈出選擇框_SAP劉夢_新浪部落格[圖]](data:image/gif;base64,R0lGODlhAQABAIAAAP///wAAACwAAAAAAQABAAACAkQBADs=)