synthw.m

这是伯克里wavelet transforms一书中的例子的代码

function sy = synthw(anal,sl,sh,levs)
%Routine to reconstruct the  signal from its multilevel DWT using an
%odd length or whole sample symmetric (WSS)  biorthogonal filter bank. 
%Symmetric extension used to extend the subbands. This is complement to
%the corresponding analysis function 'analw'.
%
%sy = synthw(an,lpf,hpf,L) reconstructs the original signal from 
%its DWT, 'an' (which could be the output of the 'analw' routine).
%The image is assumed to be decomposed to L levels. 
%The WSS filter bank used is comprised of lpf and hpf. If this 
%corresponds to the analysis filter bank used for decomposition 
%then you get perfect reconstruction. The reconstructed sequence is 
%returned in the array 'sy'.
%
%sy = synthw(an,fopt,L) reconstructs the original signal from 
%its DWT, 'an' (which could be the output of the 'analw' routine). 
%The image is assumed to be decomposed to L levels. 
%fopt chooses a WSS synthesis filter bank provided in the routine. 
%If this corresponds to the analysis filter bank used for the 
%decomposition then you get perfect reconstruction. The reconstructed 
%sequence is returned in the array 'sy'. 
%
%The input variables are
%an      : L-level DWT 
%lpf, hpf: these two filters define the synthesis filter bank
%OR
%fopt    :Takes the values of 1 and 2 and correspond to the 7/9 and 
%         9/7 filter banks used by FBI for their fingerprint compression
%         scheme. They correspond to options 1 and 2 in the analysis 
%         routine respectively.
%L       : Number of levels of decomposition and hence synthesis.
%
%The final result of the synthesis procedure is plotted in the end.
%
%Refer to Chapter 4 for information on biorthogonal wavelet 
%decomposition and reconstruction. 
%
%Author: Ajit S. Bopardikar
%Copyright (c) 1998 by Addison Wesley Longman, Inc.
%

num = nargin; %number of input arguments

  if (num == 3) %if 3 input arguments then choose your filters 
    levs = sh;
   if(sl == 1)
     sl =[-0.06453888262894 -0.04068941760956 0.41809227322221 0.78848561640566 0.41809227322221 -0.04068941760956 -0.06453888262894];
     sh =[0.03782845550700 0.02384946501938 -0.11062440441840 -0.37740285561265 0.85269867900940 -0.37740285561265 -0.11062440441840 0.02384946501938 0.03782845550700];
   elseif (sl >=2)
      if (sl > 2)
      fprintf('fopt chosen to be greater than 2. Using fopt=2 instead\n');
      end; 
     sl =[0.03782845550700 -0.02384946501938 -0.11062440441842 0.37740285561265 0.85269867900940 0.37740285561265 -0.11062440441842 -0.02384946501938 0.03782845550700];
     sh =[0.06453888262894 -0.04068941760956 -0.41809227322221 0.78848561640566 -0.41809227322221 -0.04068941760956 0.06453888262894];
   end %end inner if
  end %end if

 l = length(anal); %calculate the length...

 lens = [l]; %initialize a lengths array
 l1 =l;

 for i=1:levs-1
  if (l1/2 == round(l1/2))
    l1=l1/2;
  elseif (length(sl) > length(sh))  %analysis LPF smaller than HPF
      l1 = floor(l1/2);
   else                    %analysis LPF larger than LPF
      l1 = ceil(l1/2); 
   end %endif
   
   lens = [l1 lens];
 end;

    for i=1:levs
    sy = anal(1:lens(i));

   if (lens(i)/2 == round(lens(i)/2)) %even length
    sy = bisynt1(sy,sl,sh);%one level decomposition for even length signal
   else
    sy = bisynt3(sy,sl,sh);%one level decomposition for odd length signal
  end; %endif
  
  anal(1:lens(i)) = sy;%update the an array with the reconstructed signal
  end; %endfor 

  %plot the reconstructed signal
  figure;plot(sy);title('IDWT of Input DWT Array');