diff solvers/SMALL_MP.m @ 1:7750624e0c73 version0.5

(none)

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/solvers/SMALL_MP.m	Thu Nov 05 16:36:01 2009 +0000
@@ -0,0 +1,85 @@
+function [A]=SMALL_MP(Dict,X, m, maxNumCoef, errorGoal, varargin) 
+%%
+%=============================================
+% Sparse coding of a group of signals based on a given 
+% dictionary and specified number of atoms to use. 
+%  input arguments: Dict - the dictionary
+%                   X - the signals to represent
+%                   m - number of atoms in Dictionary
+%                   errorGoal - the maximal allowed representation error for
+%                   each signal.
+%
+% optional:         if Dict is function handle then Transpose Dictionary
+%                   handle needs to be specified.
+%
+% based on KSVD toolbox solver found on Miki Elad webpage (finding inverse
+% with pinv() is changed with matching pursuit)
+% Ivan Damnjanovic 2009
+%=============================================
+%%
+%%
+% This Dictionary check is based on Thomas Blumensath work in sparsify 0_4 greedy solvers 
+explicitD=0;
+if     isa(Dict,'float')      
+            D =@(z) Dict*z;  
+            Dt =@(z) Dict'*z;
+            explicitD=1;
+elseif isobject(Dict) 
+            D =@(z) Dict*z;  
+            Dt =@(z) Dict'*z;
+elseif isa(Dict,'function_handle') 
+    try
+        DictT=varargin{1};
+        if isa(DictT,'function_handle'); 
+                D=Dict; 
+                Dt=DictT;
+        else
+                error('If Dictionary is a function handle,Transpose Dictionary also needs to be a function handle. ');
+        end
+    catch
+        error('If Dictionary is a function handle, Transpose Dictionary  needs to be specified. Exiting.');
+    end
+else
+    error('Dictionary is of unsupported type. Use explicit matrix, function_handle or object. Exiting.');
+end
+%%
+[n,P]=size(X);
+
+E2 = errorGoal^2;
+
+
+A = sparse(m,size(X,2));
+for k=1:1:P,
+   
+    x = X(:,k);
+    residual=x;
+	indx = [];
+	j=0;
+    
+	currResNorm = norm(residual);
+    errorGoal=errorGoal*currResNorm;
+	
+    a = zeros(m,1);
+
+
+    while currResNorm>errorGoal & j < maxNumCoef,
+		
+        j = j+1;
+        dir=Dt(residual);
+        
+        [tmp__, pos]=max(abs(dir));
+        indx(j)=pos;
+        a(pos)=dir(pos);
+        e = zeros(m,1);
+        e(pos) = 1;
+        residual=residual-D(e)*a(pos);
+    
+		currResNorm = norm(residual);
+   end;
+   if (~isempty(indx))
+       A(indx,k)=a(indx);
+   end
+end;
+return;
+
+