hsss1.m

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

function sy = hsss1(an,sl,sh)
%sy = hsss1(an,lpf,hpf) performs a one level reconstruction of a 
%1 level DWT, an, using half sample symmetric (HSS) synthesis filters 
%lpf and hpf. The resultant synthesized 2-D function is returned to sy.
%
%This routine is used by the routine synthh2d.
%
%Author: Ajit S. Bopardikar
%Copyright (c) 1998 by Addison Wesley Longman, Inc.
%

  [m,n] = size(an); %size of the original image

  ll = length(sl);
  lh = length(sh);
  le = ll/2;
  he = lh/2;

%first reconstruct the columns

  if(m/2 == round(m/2)) %even number of rows
    
    if(le/2 ~=round(le/2)) %unequal subbands
      lo = an(1:m/2+1,:)';
      hi = an(m/2+2:m,:)'; %isolate the subbands

      %upsample
      u0 = zeros([n m+1]);
      u1 = zeros([n m-1]);

      u0(:,1:2:m+1) = lo;
      u1(:,2:2:m-1) = hi;

      %extensions
      extl = [u0(:,le+1:-1:2) u0 u0(:,m:-1:max(1,m-le+1))];
      exth = [-u1(:,he:-1:1) zeros([n 1]) u1 zeros([n 1]) -u1(:,m-1:-1:max(1,m-he))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      syc = r0(:,ll+1:ll+m)' + r1(:,lh+1:lh+m)';
    else %equal subbands
      %isolate the subbands
      lo = an(1:m/2,:)';
      hi = an(m/2+1:m,:)';

      %upsample
      u0 = zeros([n m]);
      u1 = zeros([n m+1]);

      u0(:,1:2:m) = lo;
      u1(:,2:2:m+1) = hi; 

      %extensions
      extl=[u0(:,le-1:-1:1) zeros([n 1]) u0 u0(:,m-1:-1:max(1,m-le-1))];
      exth=[-u1(:,he+1:-1:2) u1 -u1(:,m:-1:max(1,m-he-1))]; 

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      syc = r0(:,ll:ll+m-1)' + r1(:,lh+1:lh+m)';
    end %end inner if
  else %odd number of rows

    %isolate the subbands
    lo = an(1:ceil(m/2),:)';
    hi = an(ceil(m/2)+1:m,:)';

    u0 = zeros([n m]);
    u1 = zeros([n m]);

    if (le/2 ~= round(le/2)) 
      %upsample     
      u0(:,1:2:m) = lo;
      u1(:,2:2:m) = hi;

      extl=[u0(:,le+1:-1:2) u0 zeros([n 1]) u0(:,m:-1:max(1,m-le-1))];
      exth=[-u1(:,he:-1:1) zeros([n 1]) u1 -u1(:,m-1:-1:max(1,m-he-1))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      syc = r0(:,ll+1:ll+m)'+r1(:,lh+1:lh+m)';
    else
      %upsample
      u0(:,1:2:m) = lo;
      u1(:,2:2:m) = hi;
      
      extl=[u0(:,le-1:-1:1) zeros([n 1]) u0 u0(:,m-1:-1:max(1,m-le-1))];
      exth=[-u1(:,he+1:-1:2) u1 zeros([n 1]) -u1(:,m:-1:max(1,m-he-2))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      syc = r0(:,ll:ll+m-1)'+r1(:,lh+1:lh+m)';
    end %end inner if
%this finishes reconstructing the columns-now for the rows
  end %end if

if(n/2 == round(n/2)) %even number of columns
    
    if(le/2 ~=round(le/2)) %unequal subbands
      lo = syc(:,1:n/2+1);
      hi = syc(:,n/2+2:n); %isolate the subbands

      %upsample
      u0 = zeros([m n+1]);
      u1 = zeros([m n-1]);

      u0(:,1:2:n+1) = lo;
      u1(:,2:2:n-1) = hi;

      %extensions
      extl = [u0(:,le+1:-1:2) u0 u0(:,n:-1:max(1,n-le+1))];
      exth = [-u1(:,he:-1:1) zeros([m 1]) u1 zeros([m 1]) -u1(:,n-1:-1:max(1,n-he))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      sy = r0(:,ll+1:ll+n) + r1(:,lh+1:lh+n);
    else %equal subbands
      %isolate the subbands
      lo = syc(:,1:n/2);
      hi = syc(:,n/2+1:n);

      %upsample
      u0 = zeros([m n]);
      u1 = zeros([m n+1]);

      u0(:,1:2:n) = lo;
      u1(:,2:2:n+1) = hi; 

      %extensions
      extl=[u0(:,le-1:-1:1) zeros([m 1]) u0 u0(:,n-1:-1:max(1,n-le-1))];
      exth=[-u1(:,he+1:-1:2) u1 -u1(:,n:-1:max(1,n-he-1))]; 

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      sy = r0(:,ll:ll+n-1) + r1(:,lh+1:lh+n);
    end %end inner if
  else %odd number of rows

    %isolate the subbands
    lo = syc(:,1:ceil(n/2));
    hi = syc(:,ceil(n/2)+1:n);

    u0 = zeros([m n]);
    u1 = zeros([m n]);

    if (le/2 ~= round(le/2)) 
      %upsample     
      u0(:,1:2:n) = lo;
      u1(:,2:2:n) = hi;

      extl=[u0(:,le+1:-1:2) u0 zeros([m 1]) u0(:,n:-1:max(1,n-le-1))];
      exth=[-u1(:,he:-1:1) zeros([m 1]) u1 -u1(:,n-1:-1:max(1,n-he-1))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      sy = r0(:,ll+1:ll+n) + r1(:,lh+1:lh+n);
    else
      %upsample
      u0(:,1:2:n) = lo;
      u1(:,2:2:n) = hi;
      
      extl=[u0(:,le-1:-1:1) zeros([m 1]) u0 u0(:,n-1:-1:max(1,n-le-1))];
      exth=[-u1(:,he+1:-1:2) u1 zeros([m 1]) -u1(:,n:-1:max(1,n-he-2))];

      r0 = conv2(extl,sl);
      r1 = conv2(exth,sh);
      sy = r0(:,ll:ll+n-1)+r1(:,lh+1:lh+n);
    end %end inner if                
  end %endif