mapcd.m

Sparse Signal Representation using Overlapping Frames (matlab toolbox)

function [Ctab,Dtab]=MapCD(G)
% MapCD     Makes the tables that map from WWT (or W*W') to C 
%           and from XWT (or X*W') to D.
%
% [Ctab,Dtab]=MapCD(G);
% ----------------------------------------------
% arguments:
%  G         - The structure of the frame, only used if Type='g', size is NxKxP
%              Note that Q is given by max(G(:))
%  Ctab      - a map from WWT (size KxKxIc) to C (size QxQ)          
%              Now we assume 1D signal, then Ic=P. 
%              The structure is more complex when we have 2D or multi-D signals.
%              Note that we use single indexing into both C and A
%              C(q1,q2)=C(ci)  where ci=(q2-1)*Q+q1 and
%              WWT(k1,k2,p)=WWT(ai)  where ai=((p-1)*K+(k2-1))*K+k1
%              We may then build C from WWT with the following statements:
%                  C=zeros(Q,Q);
%                  for i=1:size(Ctab,1)
%                     ci=Ctab(i,1);ai=Ctab(i,2)
%                     C(ci)=C(ci)+sign(ai)*WWT(abs(ai));
%                  end
%                  C=C+C'-diag(diag(C));   % IMPORTANT!
%  Dtab      - a map from XWT (size NxKxP) to D (size Qx1)          
%               We also use single indexing into XWT
%               XWT(n,k,p)=XWT(bi)  where bi=((p-1)*K+(k-1))*N+n
%               We may then build D from XWT with the following statements:
%                  D=zeros(Q,1);
%                  for i=1:size(Dtab,1)
%                     di=Dtab(i,1);bi=Dtab(i,2)
%                     D(di)=D(di)+sign(bi)*XWT(abs(bi));
%                  end
%----------------------------------------------------------------------

%----------------------------------------------------------------------
% Copyright (c) 2000.  Karl Skretting.  All rights reserved.
% Hogskolen in Stavanger (Stavanger University), Signal Processing Group
% Mail:  [e-mail]   Homepage:  http://www.ux.his.no/~karlsk/
% 
% HISTORY:  dd.mm.yyyy
% Ver. 1.0  15.03.2000  KS  made function
% Ver. 1.1  09.01.2001  KS  function renamed MapCD (from MapACBd)
% Ver. 1.2  04.12.2002  KS: moved from ..\Frames\ to ..\FrameTools
%----------------------------------------------------------------------

Mfile='MapCD';

[N,K,P]=size(G);
Q=max(G(:));

% try to find lengths of these tables
Gtemp=abs(sign(G));
Gtemp=reshape(Gtemp,N,K*P);
temp=sum(Gtemp');
DtabL=sum(temp);
CtabL=temp*temp';  % this is size if we include upper right triangular of C
                   % (which is unneccessary since C=C')
% find elements on diagonal
Gtemp=reshape(abs(G),N,K*P)';
temp=zeros(1,N);
for q=1:Q
   temp=temp+(sum(Gtemp==q)).^2;
end
CtabL=(CtabL+sum(temp))/2;   % and now we have correct length

Ctab=zeros(CtabL,2);
Dtab=zeros(DtabL,2);
DtabI=0;CtabI=0;         % the indexes

for n=1:N 
   for p1=1:P
      for k1=1:K
         q1=G(n,k1,p1);
         if q1
            for p2=1:P 
               for k2=1:K
                  q2=G(n,k2,p2);
                  if q2
                     if (abs(q1) >= abs(q2))
                        ci=(abs(q2)-1)*Q+abs(q1);
                        if (p1>p2)
                           ai=((p1-p2)*K+k1-1)*K+k2;
                        else
                           ai=((p2-p1)*K+k2-1)*K+k1;
                        end
                        CtabI=CtabI+1;
                        Ctab(CtabI,:)=[ci,sign(q1)*sign(q2)*ai];
                     end
                  end
               end 
            end
            DtabI=DtabI+1;
            Dtab(DtabI,:)=[abs(q1),sign(q1)*(((p1-1)*K+k1-1)*N+n)];
         end
      end
   end
end

return