ivan@137
|
1 function gamma = omp2(varargin)
|
ivan@137
|
2 %OMP2 Error-constrained Orthogonal Matching Pursuit.
|
ivan@137
|
3 % GAMMA = OMP2(D,X,G,EPSILON) solves the optimization problem
|
ivan@137
|
4 %
|
ivan@137
|
5 % min |GAMMA|_0 s.t. |X - D*GAMMA|_2 <= EPSILON
|
ivan@137
|
6 % gamma
|
ivan@137
|
7 %
|
ivan@137
|
8 % for each of the signals in X, using Batch Orthogonal Matching Pursuit.
|
ivan@137
|
9 % Here, D is a dictionary with normalized columns, X is a matrix
|
ivan@137
|
10 % containing column signals, EPSILON is the error target for each signal,
|
ivan@137
|
11 % and G is the Gramm matrix D'*D. The output GAMMA is a matrix containing
|
ivan@137
|
12 % the sparse representations as its columns.
|
ivan@137
|
13 %
|
ivan@137
|
14 % GAMMA = OMP2(D,X,[],EPSILON) performs the same operation, but without
|
ivan@137
|
15 % the matrix G, using OMP-Cholesky. This call produces the same output as
|
ivan@137
|
16 % Batch-OMP, but is significantly slower. Using this syntax is only
|
ivan@137
|
17 % recommended when available memory is too small to store G.
|
ivan@137
|
18 %
|
ivan@137
|
19 % GAMMA = OMP2(DtX,XtX,G,EPSILON) is the fastest implementation of OMP2,
|
ivan@137
|
20 % but also requires the most memory. Here, DtX stores the projections
|
ivan@137
|
21 % D'*X, and XtX is a row vector containing the squared norms of the
|
ivan@137
|
22 % signals, sum(X.*X). In this case Batch-OMP is used, but without having
|
ivan@137
|
23 % to compute D'*X and XtX in advance, which slightly improves runtime.
|
ivan@137
|
24 % Note that in general, the call
|
ivan@137
|
25 %
|
ivan@137
|
26 % GAMMA = OMP2(D'*X, sum(X.*X), G, EPSILON);
|
ivan@137
|
27 %
|
ivan@137
|
28 % will be faster than the call
|
ivan@137
|
29 %
|
ivan@137
|
30 % GAMMA = OMP2(D,X,G,EPSILON);
|
ivan@137
|
31 %
|
ivan@137
|
32 % due to optimized matrix multiplications in Matlab. However, when the
|
ivan@137
|
33 % entire matrix D'*X cannot be stored in memory, one of the other two
|
ivan@137
|
34 % versions can be used. Both compute D'*X for just one signal at a time,
|
ivan@137
|
35 % and thus require much less memory.
|
ivan@137
|
36 %
|
ivan@137
|
37 % GAMMA = OMP2(...,PARAM1,VAL1,PARAM2,VAL2,...) specifies additional
|
ivan@137
|
38 % parameters for OMP2. Available parameters are:
|
ivan@137
|
39 %
|
ivan@137
|
40 % 'gammamode' - Specifies the representation mode for GAMMA. Can be
|
ivan@137
|
41 % either 'full' or 'sparse', corresponding to a full or
|
ivan@137
|
42 % sparse matrix, respectively. By default, GAMMA is
|
ivan@137
|
43 % returned as a sparse matrix.
|
ivan@137
|
44 % 'maxatoms' - Limits the number of atoms in the representation of each
|
ivan@137
|
45 % signal. If specified, the number of atoms in each
|
ivan@137
|
46 % representation does not exceed this number, even if the
|
ivan@137
|
47 % error target is not met. Specifying maxatoms<0 implies
|
ivan@137
|
48 % no limit (default).
|
ivan@137
|
49 % 'messages' - Specifies whether progress messages should be displayed.
|
ivan@137
|
50 % When positive, this is the number of seconds between
|
ivan@137
|
51 % status prints. When negative, indicates that no messages
|
ivan@137
|
52 % should be displayed (this is the default).
|
ivan@137
|
53 % 'checkdict' - Specifies whether dictionary normalization should be
|
ivan@137
|
54 % verified. When set to 'on' (default) the dictionary
|
ivan@137
|
55 % atoms are verified to be of unit L2-norm. Setting this
|
ivan@137
|
56 % parameter to 'off' disables verification and accelerates
|
ivan@137
|
57 % function performance. Note that an unnormalized
|
ivan@137
|
58 % dictionary will produce invalid results.
|
ivan@137
|
59 % 'profile' - Can be either 'on' or 'off'. When 'on', profiling
|
ivan@137
|
60 % information is displayed at the end of the funciton
|
ivan@137
|
61 % execution.
|
ivan@137
|
62 %
|
ivan@137
|
63 %
|
ivan@137
|
64 % Summary of OMP2 versions:
|
ivan@137
|
65 %
|
ivan@137
|
66 % version | speed | memory
|
ivan@137
|
67 % -------------------------------------------------------------
|
ivan@137
|
68 % OMP2(DtX,XtX,G,EPSILON) | very fast | very large
|
ivan@137
|
69 % OMP2(D,X,G,EPSILON) | fast | moderate
|
ivan@137
|
70 % OMP2(D,X,[],EPSILON) | very slow | small
|
ivan@137
|
71 % -------------------------------------------------------------
|
ivan@137
|
72 %
|
ivan@137
|
73 %
|
ivan@137
|
74 % References:
|
ivan@137
|
75 % [1] M. Elad, R. Rubinstein, and M. Zibulevsky, "Efficient Implementation
|
ivan@137
|
76 % of the K-SVD Algorithm using Batch Orthogonal Matching Pursuit",
|
ivan@137
|
77 % Technical Report - CS, Technion, April 2008.
|
ivan@137
|
78 %
|
ivan@137
|
79 % See also OMP.
|
ivan@137
|
80
|
ivan@137
|
81
|
ivan@137
|
82 % Ron Rubinstein
|
ivan@137
|
83 % Computer Science Department
|
ivan@137
|
84 % Technion, Haifa 32000 Israel
|
ivan@137
|
85 % ronrubin@cs
|
ivan@137
|
86 %
|
ivan@137
|
87 % April 2009
|
ivan@137
|
88
|
ivan@137
|
89
|
ivan@137
|
90 % default options
|
ivan@137
|
91
|
ivan@137
|
92 sparse_gamma = 1;
|
ivan@137
|
93 msgdelta = -1;
|
ivan@137
|
94 maxatoms = -1;
|
ivan@137
|
95 checkdict = 1;
|
ivan@137
|
96 profile = 0;
|
ivan@137
|
97
|
ivan@137
|
98
|
ivan@137
|
99 % determine number of parameters
|
ivan@137
|
100
|
ivan@137
|
101 paramnum = 1;
|
ivan@137
|
102 while (paramnum<=nargin && ~ischar(varargin{paramnum}))
|
ivan@137
|
103 paramnum = paramnum+1;
|
ivan@137
|
104 end
|
ivan@137
|
105 paramnum = paramnum-1;
|
ivan@137
|
106
|
ivan@137
|
107
|
ivan@137
|
108 % parse options
|
ivan@137
|
109
|
ivan@137
|
110 for i = paramnum+1:2:length(varargin)
|
ivan@137
|
111 paramname = varargin{i};
|
ivan@137
|
112 paramval = varargin{i+1};
|
ivan@137
|
113
|
ivan@137
|
114 switch lower(paramname)
|
ivan@137
|
115
|
ivan@137
|
116 case 'gammamode'
|
ivan@137
|
117 if (strcmpi(paramval,'sparse'))
|
ivan@137
|
118 sparse_gamma = 1;
|
ivan@137
|
119 elseif (strcmpi(paramval,'full'))
|
ivan@137
|
120 sparse_gamma = 0;
|
ivan@137
|
121 else
|
ivan@137
|
122 error('Invalid GAMMA mode');
|
ivan@137
|
123 end
|
ivan@137
|
124
|
ivan@137
|
125 case 'maxatoms'
|
ivan@137
|
126 maxatoms = paramval;
|
ivan@137
|
127
|
ivan@137
|
128 case 'messages'
|
ivan@137
|
129 msgdelta = paramval;
|
ivan@137
|
130
|
ivan@137
|
131 case 'checkdict'
|
ivan@137
|
132 if (strcmpi(paramval,'on'))
|
ivan@137
|
133 checkdict = 1;
|
ivan@137
|
134 elseif (strcmpi(paramval,'off'))
|
ivan@137
|
135 checkdict = 0;
|
ivan@137
|
136 else
|
ivan@137
|
137 error('Invalid checkdict option');
|
ivan@137
|
138 end
|
ivan@137
|
139
|
ivan@137
|
140 case 'profile'
|
ivan@137
|
141 if (strcmpi(paramval,'on'))
|
ivan@137
|
142 profile = 1;
|
ivan@137
|
143 elseif (strcmpi(paramval,'off'))
|
ivan@137
|
144 profile = 0;
|
ivan@137
|
145 else
|
ivan@137
|
146 error('Invalid profile mode');
|
ivan@137
|
147 end
|
ivan@137
|
148
|
ivan@137
|
149 otherwise
|
ivan@137
|
150 error(['Unknown option: ' paramname]);
|
ivan@137
|
151 end
|
ivan@137
|
152
|
ivan@137
|
153 end
|
ivan@137
|
154
|
ivan@137
|
155
|
ivan@137
|
156 % determine call type
|
ivan@137
|
157
|
ivan@137
|
158 if (paramnum==4)
|
ivan@137
|
159
|
ivan@137
|
160 n1 = size(varargin{1},1);
|
ivan@137
|
161 n2 = size(varargin{2},1);
|
ivan@137
|
162 n3 = size(varargin{3},1);
|
ivan@137
|
163
|
ivan@137
|
164 if ( (n1>1 && n2==1) || (n1==1 && n2==1 && n3==1) ) % DtX,XtX,G,EPSILON
|
ivan@137
|
165
|
ivan@137
|
166 DtX = varargin{1};
|
ivan@137
|
167 XtX = varargin{2};
|
ivan@137
|
168 G = varargin{3};
|
ivan@137
|
169 epsilon = varargin{4};
|
ivan@137
|
170 D = [];
|
ivan@137
|
171 X = [];
|
ivan@137
|
172
|
ivan@137
|
173 else % D,X,G,EPSILON
|
ivan@137
|
174
|
ivan@137
|
175 D = varargin{1};
|
ivan@137
|
176 X = varargin{2};
|
ivan@137
|
177 G = varargin{3};
|
ivan@137
|
178 epsilon = varargin{4};
|
ivan@137
|
179 DtX = [];
|
ivan@137
|
180 XtX = [];
|
ivan@137
|
181
|
ivan@137
|
182 end
|
ivan@137
|
183
|
ivan@137
|
184 else
|
ivan@137
|
185 error('Invalid number of parameters');
|
ivan@137
|
186 end
|
ivan@137
|
187
|
ivan@137
|
188
|
ivan@137
|
189 % verify dictionary normalization
|
ivan@137
|
190
|
ivan@137
|
191 if (checkdict)
|
ivan@137
|
192 if (isempty(G))
|
ivan@137
|
193 atomnorms = sum(D.*D);
|
ivan@137
|
194 else
|
ivan@137
|
195 atomnorms = diag(G);
|
ivan@137
|
196 end
|
ivan@137
|
197 if (any(abs(atomnorms-1) > 1e-2))
|
ivan@137
|
198 error('Dictionary columns must be normalized to unit length');
|
ivan@137
|
199 end
|
ivan@137
|
200 end
|
ivan@137
|
201
|
ivan@137
|
202
|
ivan@137
|
203 % omp
|
ivan@137
|
204
|
ivan@137
|
205 gamma = omp2mexGabor(D,X,DtX,XtX,G,epsilon,sparse_gamma,msgdelta,maxatoms,profile);
|