annotate toolboxes/AudioInpaintingToolbox/Solvers/inpaintFrame_OMP_Gabor.m @ 233:88a5c02d20d3

Added tag ver_2.0_alpha1 for changeset 15fd4c5be71a
author luisf <luis.figueira@eecs.qmul.ac.uk>
date Thu, 19 Apr 2012 15:55:59 +0100
parents 56d719a5fd31
children
rev   line source
ivan@138 1 function [y , Coeff]= inpaintFrame_OMP_Gabor(problemData,param)
ivan@138 2 % Inpainting method based on OMP using the Gabor dictionary
ivan@138 3 % generated by Gabor_Dictionary.m. The method jointly selects
ivan@138 4 % cosine and sine atoms at the same frequency
ivan@138 5 %
ivan@138 6 % Usage: y = inpaintFrame_OMP_Gabor(problemData,param)
ivan@138 7 %
ivan@138 8 %
ivan@138 9 % Inputs:
ivan@138 10 % - problemData.x: observed signal to be inpainted
ivan@138 11 % - problemData.Imiss: Indices of clean samples
ivan@138 12 % - param.D - the dictionary matrix (optional if param.D_fun is set)
ivan@138 13 % - param.D_fun - a function handle that generates the dictionary
ivan@138 14 % matrix param.D if param.D is not given. See Gabor_Dictionary.m
ivan@138 15 % - param.wa - Analysis window
ivan@138 16 %
ivan@138 17 % Outputs:
ivan@138 18 % - y: estimated frame
ivan@138 19 %
ivan@138 20 %
ivan@138 21 % -------------------
ivan@138 22 %
ivan@138 23 % Audio Inpainting toolbox
ivan@138 24 % Date: June 28, 2011
ivan@138 25 % By Valentin Emiya, Amir Adler, Michael Elad, Maria Jafari
ivan@138 26 % This code is distributed under the terms of the GNU Public License version 3 (http://www.gnu.org/licenses/gpl.txt).
ivan@138 27 % ========================================================
ivan@138 28 % To do next: use a faster implementation of OMP
ivan@138 29
ivan@138 30
ivan@138 31 %% Load data and parameters
ivan@138 32
ivan@138 33
ivan@138 34 x = problemData.x;
ivan@138 35 IObs = find(~problemData.IMiss);
ivan@138 36 p.N = length(x);
ivan@138 37 E2 = param.OMPerr^2;
ivan@138 38 E2M=E2*length(IObs);
ivan@138 39 wa = param.wa(param.N);
ivan@138 40
ivan@138 41 %% Build and normalized dictionary
ivan@138 42 % build the dictionary matrix if only the dictionary generation function is given
ivan@138 43 if ~isfield(param,'D')
ivan@138 44 param.D = param.D_fun(param);
ivan@138 45 end
ivan@138 46 Dict=param.D(IObs,:);
ivan@138 47 W=1./sqrt(diag(Dict'*Dict));
ivan@138 48 Dict=Dict*diag(W);
ivan@138 49 Dict1 = Dict(:,1:end/2);
ivan@138 50 Dict2 = Dict(:,end/2+1:end);
ivan@138 51 Dict1Dict2 = sum(Dict1.*Dict2);
ivan@138 52 n12 = 1./(1-Dict1Dict2.^2);
ivan@138 53 xObs=x(IObs);
ivan@138 54
ivan@138 55 %% OMP iterations
ivan@138 56 residual=xObs;
ivan@138 57 maxNumCoef = param.sparsityDegree;
ivan@138 58 indx = [];
ivan@138 59 % currResNorm2 = sum(residual.^2);
ivan@138 60 currResNorm2 = E2M*2; % set a value above the threshold in order to have/force at least one loop executed
ivan@138 61 j = 0;
ivan@138 62 while currResNorm2>E2M && j < maxNumCoef,
ivan@138 63 j = j+1;
ivan@138 64 proj=residual'*Dict;
ivan@138 65 proj1 = proj(1:end/2);
ivan@138 66 proj2 = proj(end/2+1:end);
ivan@138 67
ivan@138 68 alpha_j = (proj1-Dict1Dict2.*proj2).*n12;
ivan@138 69 beta_j = (proj2-Dict1Dict2.*proj1).*n12;
ivan@138 70
ivan@138 71 err_j = sum(abs(repmat(residual,1,size(Dict1,2))-Dict1*sparse(diag(alpha_j))-Dict2*sparse(diag(beta_j))).^2);
ivan@138 72 [dum pos] = min(err_j);
ivan@138 73
ivan@138 74 indx(end+1)=pos;
ivan@138 75 indx(end+1)=pos+size(Dict1,2);
ivan@138 76 a=pinv(Dict(:,indx(1:2*j)))*xObs;
ivan@138 77 residual=xObs-Dict(:,indx(1:2*j))*a;
ivan@138 78 currResNorm2=sum(residual.^2);
ivan@138 79 end;
ivan@138 80
ivan@138 81 %% Frame Reconstruction
ivan@138 82 indx(length(a)+1:end) = [];
ivan@138 83
ivan@138 84 Coeff = sparse(size(param.D,2),1);
ivan@138 85 if (~isempty(indx))
ivan@138 86 Coeff(indx) = a;
ivan@138 87 Coeff = W.*Coeff;
ivan@138 88 end
ivan@138 89 y = param.D*Coeff;
ivan@138 90
ivan@138 91 return
ivan@138 92
ivan@138 93