THESE/codes/wave/allcode/dualtree2D.m

   1 function w = dualtree2D(x, J, Faf, af)
   2
   3 % 2D Dual-Tree Discrete Wavelet Transform
   4 %
   5 % USAGE:
   6 %   w = dualtree2D(x, J, Faf, af)
   7 % INPUT:
   8 %   x - M by N array
   9 %   J - number of stages
  10 %   Faf - first stage filters
  11 %   af - filters for remaining stages
  12 % OUPUT:
  13 %   w{j}{d1}{d2} - DWT coefficients
  14 %       j = 1..J, k = 1..2, d = 1..3
  15 %   w{J+1}{k} - lowpass coefficients
  16 %       k = 1..2
  17 % % EXAMPLE:
  18 %   x = rand(256);
  19 %   J = 3;
  20 %   [Faf, Fsf] = FSfarras;
  21 %   [af, sf] = dualfilt1;
  22 %   w = dualtree2D(x, J, Faf, af);
  23 %   y = idualtree2D(w, J, Fsf, sf);
  24 %   err = x - y;
  25 %   max(max(abs(err)))
  26 %
  27 % WAVELET SOFTWARE AT POLYTECHNIC UNIVERSITY, BROOKLYN, NY
  28 % http://taco.poly.edu/WaveletSoftware/
  29
  30 % normalization
  31 x = x/sqrt(2);
  32
  33 % Tree 1
  34 [x1 w{1}{1}] = afb2D(x, Faf{1});      % stage 1
  35 for j = 2:J
  36     [x1 w{j}{1}] = afb2D(x1, af{1});  % remaining stages
  37 end
  38 w{J+1}{1} = x1;                       % lowpass subband
  39
  40 % Tree 2
  41 [x2 w{1}{2}] = afb2D(x, Faf{2});      % stage 1
  42 for j = 2:J
  43     [x2 w{j}{2}] = afb2D(x2, af{2});  % remaining stages
  44 end
  45 w{J+1}{2} = x2;                       % lowpass subband
  46
  47 % sum and difference
  48 for j = 1:J
  49     for m = 1:3
  50         A = w{j}{1}{m};
  51         B = w{j}{2}{m};
  52         w{j}{1}{m} = (A+B)/sqrt(2);
  53         w{j}{2}{m} = (A-B)/sqrt(2);
  54     end
  55 end
  56