wavelet transform and MATLAB source __matlab

When the call to the Wavefast function found that MATLAB does not have this function, by looking up the wavelet transformation of the M file, saved into a function file can be directly invoked.

Matlab Source code:
(1) wave2gray.m

function w = Wave2gray (c, S, scale, border)%wave2gray Display wavelet decomposition.  
% 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 ' Bor                   Der 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, vertic  

Al, and diagonal detail% coefficients, respectively.  

% Check input arguments for reasonableness.  

Error (Nargchk (2, 4, nargin));   if (Ndims (c) ~= 2) | (Size (c, 1) ~= 1) error (' C must be a row vector. '); End If (ndims (s) ~= 2) | ~isreal (s) |   ~isnumeric (s) | (Size (S, 2) ~= 2) error (' s must being a real, numeric two-column array. ');  

End elements = Prod (s, 2);   if (length (c) < elements (end)) | ... ~ (elements (1) + 3 * SUM (Elements (2:end-1)) >= elements (end)) error ([' [C S] must is a standard WAV  
                Elet ' ...  
' decomposition structure. ']); End If (Nargin > 2) & (~isreal (scale) | ~isnumeRic (Scale)) error (' Scale must being a real, numeric scalar. ');  
End If (Nargin > 3) & (~ischar (border)) error (' Border must is character string. ');    End If Nargin = 2 scale = 1;  
% Default scale. End If Nargin < 4 border = ' absorb ';  
% Default border.  
End% Scale coefficients and determine pad fill.  
Absflag = Scale < 0;  
Scale = ABS (scale);  
If scale = = 0 scale = 1;  End [CD, W] = Wavecut (' A ', C, s);  
W = Mat2gray (w);  
CDX = MAX (ABS (CD (:)))/scale; If Absflag cd = Mat2gray (ABS (CD), [0, CDX]);  
Fill = 0; else CD = Mat2gray (CD, [-CDX, CDX]);  
Fill = 0.5;  
End% builds gray image one decomposition at a time.  
    For i = size (s, 1)-2:-1:1 ws = Size (w);  
    h = wavecopy (' h ', CD, S, i);    Pad = ws-size (h);  
    Frontporch = round (PAD/2);  
    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/2);  
    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/2);  
    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, [1 1], 1, ' post ');  
            h = Padarray (h, [1 0], 1, ' post ');  
        v = Padarray (V, [0 1], 1, ' post ');         Case ' Absorb ' W (:, end) = 1;  
            W (end,:) = 1;          H (end,:) = 1;  
        V (:, end) = 1;  
        Otherwise error (' unrecognized BORDER parameter. ');                                          End w = [w h; v d];  
% concatenate coefs. End If nargout = 0 
    Imshow (w);  
% Display result.  
 End

(2) WAVEBACK.M

 function [Varargout] = Waveback (c, S, Varargin)%waveback performs a multi-level two-dimensional inverse. % [Varargout] = Waveback (C, S, Varargin) computes a 2D n-level% partial or complete wavelet reconstruction of decom   
Position% 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  
Reconstruction.  

% of also wavefast and Wavefilter. % Copyright 2002-2004 R. Gonzalez, R. E Woods, & S. L.  Eddins% Digital Image processing Using MATLAB, Prentice-hall,% $Revision: 1.5 $ $Date: 2003/10/13 01:29:36  
$% Check the input and output arguments for reasonableness.  
Error (Nargchk (3, 5, nargin));  

Error (NARGCHK (1, 2, nargout)); if (Ndims (c) ~= 2) |     
(Size (c, 1) ~= 1) error (' C must be a row vector. '); End If (ndims (s) ~= 2) | ~isreal (s) | ~isnumeric (s) |     
(Size (s,2) ~= 2) error (' s must is a real, numeric two-column array. ');  
End elements = Prod (s, 2); if (length (c) < elements (end)) |  ... ~ (elements (1) + 3 * SUM (Elements (2:end-1)) >= elements (end)) error ([' [C S] must is a standard wavelet  
         ' ...   
' decomposition structure. ']);  
End% Maximum levels in [C, S].          
Nmax = Size (s, 1)-2;  
% get third input parameter and init check flags.  Wname = Varargin{1};   Filterchk = 0;      

Nchk = 0;   Switch Nargin Case 3 if Ischar (Wname) [LP, HP] = Wavefilter (wname, ' R ');  
n = nmax;   else error (' Undefined filter. ');    
   End If nargout ~= 1 error (' Wrong number of output arguments. ');     
      End Case 4 if Ischar (Wname) [LP, HP] = Wavefilter (wname, ' R ');    n = varargin{2};  
   Nchk = 1;   Else LP = Varargin{1};     
      HP = varargin{2};   Filterchk = 1;  
      n = nmax;    
      If nargout ~= 1 error (' Wrong number of output arguments. ');   End of CASE 5 LP = Varargin{1};   HP = varargin{2};  
    Filterchk = 1;    n = varargin{3};  
Nchk = 1;       
Otherwise error (' Improper number of input arguments. ');  
End fl = length (LP);  
  If Filterchk% Check filters. if (NDIMS (LP) ~= 2) | ~ISREAL (LP) | ~ISNUMERIC (LP) ... | (Ndims (HP) ~= 2) | ~isreal (HP) | ~isnumeric (HP) ... | (fl ~= Length (HP)) |  
             REM (fl, 2) ~= 0 error ([' LP and HP must is 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 < 1) error (' Invalid number (n) of reconstructions requested. ');   
End if (n ~= nmax) & (nargout ~= 2) error (' Not enough output arguments. ');    End NC = C;    NS = s;             Nnmax = Nmax;  
% Init decomposition.  
   For i = 1:n% Compute a new approximation.  
                 A = Symconvup (Wavecopy (' A ', NC, NS), LP, LP, FL, NS (3,:)) + ... symconvup (wavecopy (' H ', NC, NS, Nnmax), ...  HP, LP, FL, NS (3,:)) + ... symconvup (wavecopy (' V ', NC, NS, Nnmax), ... LP, HP,  

    FL, NS (3,:)) + ... symconvup (wavecopy (' d ', NC, NS, Nnmax), ... hp, HP, FL, NS (3,:));  
    % Update decomposition.     NC = NC (4 * PROD (ns (1,:)) + 1: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 ([1 2]. * Size (y));  
Z (:, 1:2:end) = y;  
z = Conv2 (z, F2);   Z = Z (FLN-1:FLN + Keep (1)-2, FLN-1:FLN + Keep (2)-2);

(3) WAVECOPY.M

Function y = wavecopy (type, C, S, N)%wavecopy fetches coefficients of a wavelet decomposition.     
% 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 ' d   
Iagonal details% [C, S] is a wavelet data structure.   
% N Specifies a decomposition level (ignored if TYPE = ' a ').   

% of also wavework, Wavecut, and Wavepaste. % Copyright 2002-2004 R. Gonzalez, R. Woods, & L. Eddins% Digital Image processing Using MATLAB, Pren   
Tice-hall,% $Revision: 1.4 $ $Date: 2003/10/13 01:20:41 $ error (NARGCHK (3, 4, nargin));   
if Nargin = = 4 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. % [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 ' App  Roximation 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 ').  

% of also wavework, Wavecopy, and Wavepaste. % Copyright 2002-2004 R. Gonzalez, R. Woods, & L. Eddins% Digital Image processing Using MATLAB, prent  
Ice-hall,% $Revision: 1.4 $ $Date: 2003/10/13 01:20:09 $ error (NARGCHK (3, 4, nargin)); If Nargin = = 4 [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 2-dimensional fast wavelet. % [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 is 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. This is,% if X = [C1 c2 C3 ... cn] (in 1D), then it symmetric extension% would be [... c3 C2 ...   
CN cn cn-1 Cn-2 ...].  % outputs:% Matrix C is a coefficient decomposition vector:% c = [A (n) h (n) v (n) d (n) H (n-1) ... v (1) d (1)]% where a, H, V, and D are columnwise vectors containing% approximation, horizontal, vertical, and diagonal coefficient% matrices, respectively.    
C has 3n + 1 sections where n is the% of wavelet decompositions. % Matrix S is a (n+2) x 2 bookkeeping Matrix:% s = [SA (n,:); SD (n,:); SD (n-1,:);.;; SD (1,: );   
SX]% where SA and SD are approximation and detail size entries.   

% of also waveback and Wavefilter. % Copyright 2002-2004 R. Gonzalez, R. Woods, & L. Eddins% Digital Image processing Using MATLAB, Pren   
Tice-hall% $Revision: 1.5 $ $Date: 2003/10/13 01:14:17 $% Check The input arguments for reasonableness.   

Error (NARGCHK (3, 4, nargin));   
   if Nargin = = 3 if Ischar (Varargin{1}) [LP, HP] = Wavefilter (Varargin{1}, ' d ');      
   else error (' Missing wavelet name. ');     END Else LP = Varargin{1};      
HP = varargin{2};      End fl = length (LP);   

SX = size (x); if (Ndims (x) ~= 2) | (min (SX) < 2) | ~isreal (x) |        
~isnumeric (x) error (' x must being a real, numeric matrix. '); End If (NDIMS (LP) ~= 2) | ~ISREAL (LP) | ~ISNUMERIC (LP) ... | (Ndims (HP) ~= 2) | ~isreal (HP) | ~isnumeric (HP) ... | (fl ~= Length (HP)) |   
          REM (fl, 2) ~= 0 error ([' LP and HP must is even and equal length real, ' ...    
' Numeric filter vectors. '); End If ~isreal (n) | ~isnumeric (n) | (N < 1) |   
          (N > log2 (max SX)) error ([' n must be a real scalar between 1 and ' ...)       
' Log2 (Max (((X))). "   
End% Init The starting output data structures and initial approximation.       c = [];     s = SX;   

App = double (x);   
   % for each decomposition ... for i = 1: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)% Comput E The number of coefficients to keep after convolution% and downsampling.   

Then extend x in both dimensions.   
Keep = Floor (fl + size (x)-1)/2);   

y = Padarray (x, [(FL-1) (FL-1)], ' symmetric ', ' both '); %-------------------------------------------------------------------% function y = symconv (x, H, type, FL, keep)% convolve the rows or columns of x with H, Downsamp   

Le,% and extract the center section since symmetrically extended.   
   If strcmp (type, ' row ') y = conv2 (x, h);   
   y = y (:, 1:2:end); y = y (:, FL/2 + 1:FL/