diff solvers/SMALL_ompGabor/ompGabor.m @ 140:31d2864dfdd4 ivand_dev

Audio Impainting additional constraints with cvx added
author Ivan Damnjanovic lnx <ivan.damnjanovic@eecs.qmul.ac.uk>
date Mon, 25 Jul 2011 17:27:05 +0100
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/solvers/SMALL_ompGabor/ompGabor.m	Mon Jul 25 17:27:05 2011 +0100
@@ -0,0 +1,180 @@
+function gamma = omp(varargin)
+%OMP Sparsity-constrained Orthogonal Matching Pursuit.
+%  GAMMA = OMP(D,X,G,T) solves the optimization problem
+%
+%       min  |X - D*GAMMA|_2     s.t.  |GAMMA|_0 <= T
+%      gamma
+%
+%  for each of the signals in X, using Batch Orthogonal Matching Pursuit.
+%  Here, D is a dictionary with normalized columns, X is a matrix
+%  containing column signals, T is the # of non-zeros in each signal
+%  representation, and G is the Gramm matrix D'*D. The output GAMMA is a
+%  matrix containing the sparse representations as its columns. 
+%
+%  GAMMA = OMP(D,X,[],T) performs the same operation, but without the
+%  matrix G, using OMP-Cholesky. This call produces the same output as
+%  Batch-OMP, but is significantly slower. Using this syntax is only
+%  recommended when available memory is too small to store G.
+%
+%  GAMMA = OMP(DtX,G,T) is the fastest implementation of OMP, but also
+%  requires the most memory. Here, DtX stores the projections D'*X. In this
+%  case Batch-OMP is used, but without having to compute D'*X in advance,
+%  which slightly improves runtime. Note that in general, the call
+%
+%    GAMMA = OMP(D'*X,G,T);
+%
+%  will be faster than the call
+%
+%    GAMMA = OMP(D,X,G,T);
+%
+%  due to optimized matrix multiplications in Matlab. However, when the
+%  entire matrix D'*X cannot be stored in memory, one of the other two
+%  versions can be used. Both compute D'*X for just one signal at a time,
+%  and thus require much less memory.
+%
+%  GAMMA = OMP(...,PARAM1,VAL1,PARAM2,VAL2,...) specifies additional
+%  parameters for OMP. Available parameters are:
+%
+%    'gammamode' - Specifies the representation mode for GAMMA. Can be
+%                  either 'full' or 'sparse', corresponding to a full or
+%                  sparse matrix, respectively. By default, GAMMA is
+%                  returned as a sparse matrix.
+%    'messages'  - Specifies whether progress messages should be displayed.
+%                  When positive, this is the number of seconds between
+%                  status prints. When negative, indicates that no messages
+%                  should be displayed (this is the default).
+%    'checkdict' - Specifies whether dictionary normalization should be
+%                  verified. When set to 'on' (default) the dictionary
+%                  atoms are verified to be of unit L2-norm. Setting this
+%                  parameter to 'off' disables verification and accelerates
+%                  function performance. Note that an unnormalized
+%                  dictionary will produce invalid results.
+%    'profile'   - Can be either 'on' or 'off'. When 'on', profiling
+%                  information is displayed at the end of the funciton
+%                  execution.
+%
+%
+%  Summary of OMP versions:
+%
+%    version      |   speed     |   memory
+%  --------------------------------------------------
+%   OMP(DtX,G,T)  |  very fast  |  very large
+%   OMP(D,X,G,T)  |  fast       |  moderate
+%   OMP(D,X,[],T) |  very slow  |  small
+%  --------------------------------------------------
+%
+%
+%  References:
+%  [1] M. Elad, R. Rubinstein, and M. Zibulevsky, "Efficient Implementation
+%      of the K-SVD Algorithm using Batch Orthogonal Matching Pursuit",
+%      Technical Report - CS, Technion, April 2008.
+%
+%  See also OMP2.
+
+
+%  Ron Rubinstein
+%  Computer Science Department
+%  Technion, Haifa 32000 Israel
+%  ronrubin@cs
+%
+%  April 2009
+
+
+% default options
+
+sparse_gamma = 1;
+msgdelta = -1;
+checkdict = 1;
+profile = 0;
+
+
+% determine number of parameters
+
+paramnum = 1;
+while (paramnum<=nargin && ~ischar(varargin{paramnum}))
+  paramnum = paramnum+1;
+end
+paramnum = paramnum-1;
+
+
+% parse options
+
+for i = paramnum+1:2:length(varargin)
+  paramname = varargin{i};
+  paramval = varargin{i+1};
+
+  switch lower(paramname)
+
+    case 'gammamode'
+      if (strcmpi(paramval,'sparse'))
+        sparse_gamma = 1;
+      elseif (strcmpi(paramval,'full'))
+        sparse_gamma = 0;
+      else
+        error('Invalid GAMMA mode');
+      end
+      
+    case 'messages'
+      msgdelta = paramval;
+
+    case 'checkdict'
+      if (strcmpi(paramval,'on'))
+        checkdict = 1;
+      elseif (strcmpi(paramval,'off'))
+        checkdict = 0;
+      else
+        error('Invalid checkdict option');
+      end
+
+    case 'profile'
+      if (strcmpi(paramval,'on'))
+        profile = 1;
+      elseif (strcmpi(paramval,'off'))
+        profile = 0;
+      else
+        error('Invalid profile mode');
+      end
+
+    otherwise
+      error(['Unknown option: ' paramname]);
+  end
+  
+end
+
+
+% determine call type
+
+if (paramnum==3)
+  DtX = varargin{1};
+  G = varargin{2};
+  T = varargin{3};
+  D = [];
+  X = [];
+elseif (paramnum==4)
+  D = varargin{1};
+  X = varargin{2};
+  G = varargin{3};
+  T = varargin{4};
+  DtX = [];
+else
+  error('Invalid number of parameters');
+end
+
+
+% verify dictionary normalization
+
+if (checkdict)
+  if (isempty(G))
+    atomnorms = sum(D.*D);
+  else
+    atomnorms = diag(G);
+  end
+  if (any(abs(atomnorms-1) > 1e-2))
+    error('Dictionary columns must be normalized to unit length');
+  end
+end
+
+
+% omp
+
+gamma = ompmexGabor(D,X,DtX,G,T,sparse_gamma,msgdelta,profile);