vsmp2.m

Sparse Signal Representation using Overlapping Frames (matlab toolbox)

function w=VSmp2(x,S,it)
% VSmp2     Matching Pursuit variant 2, Vector Selection algorithm 
% This is the same algorithm that sometimes is called Basic Matching Pursuit (BMP). 
% VSmp2 (this function) do 'it' iterations selecting new weights or improving
% already selected weights while 
%
% w=VSmp2(x,S,it);
% w=VSmp2(x,S);            % default value of 'it' used, it=2*S
%----------------------------------------------------------------------
% arguments:
%   w    - the weights where at most S of the elements are non-zero, size Kx1
%   x    - the vector to approximate, size Nx1
%   S    - number of vectors to select or non-zero weights 
%   it   - number of iterations to do
%-----------------------------------------------------------------------------------
% global variables:
%   F    - the dictionary (or F matrix), size NxK
%          the vectors of F must be normalized, F(:,k)'*F(:,k)==1 for all k
%          normalization may be done by: F=F*diag(1./sqrt(diag(F'*F)));
%-----------------------------------------------------------------------------------
% Note that we here use a simple complete F matrix, and do not make any assumtions
% regarding structure of F.


%----------------------------------------------------------------------
% Copyright (c) 1999-2001.  Karl Skretting.  All rights reserved.
% Hogskolen in Stavanger (Stavanger University), Signal Processing Group
% Mail:  [e-mail]   Homepage:  http://www.ux.his.no/~karlsk/
% 
% HISTORY:
% Ver. 1.0  19.04.2001  KS: function made based on VSmp
% Ver. 1.1  25.11.2002  KS: moved from ..\Frames to ..\FrameTools
%----------------------------------------------------------------------

global F

Mfile='VSmp2';
[N,K]=size(F);

if nargin<2
   error([Mfile,': wrong number of arguments, see help.']);
end
if (nargin<3); it=2*S; end;
if (it<S); it=S; end;

% start algorithm
w=zeros(K,1);      % the weights
wI=[];
wIinv=1:K;
r=x(:);
k=0;

while k<it
    c=(r'*F);                % the inner products
    if (length(wI)==S) | ((it>(4*S)) & ((k*S)<(length(wI)*(it-2*S))))
        c(wIinv)=0;           % we do not select any new vectors
    end
    [temp,i]=max(abs(c));i=i(1);
    if w(i)
        w(i)=w(i)+c(i);
    else
        w(i)=c(i);
        wI=find(w);
        wIinv=find(w==0);
    end   
    r=r-c(i)*F(:,i);
    k=k+1;                   % make sure the loop do not run forever
end

return;