[c3bd51]: frames / franaiter.m Maximize Restore History

Download this file

franaiter.m    106 lines (89 with data), 3.0 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
function [c,relres,iter]=franaiter(F,f,varargin)
%FRANAITER Iterative analysis
% Usage: c=franaiter(F,f);
% [c,relres,iter]=franaiter(F,f,...);
%
% Input parameters:
% F : Frame.
% f : Signal.
% Ls : Length of signal.
% Output parameters:
% c : Array of coefficients.
% relres : Vector of residuals.
% iter : Number of iterations done.
%
% `c=franaiter(F,f)` computes the frame coefficients *c* of the signal *f*
% using an iterative method such that perfect reconstruction can be
% obtained using |frsyn|. `franaiter` always works, even when |frana|
% cannot generate perfect reconstruction coefficients.
%
% `[c,relres,iter]=franaiter(...)` additionally returns the relative
% residuals in a vector *relres* and the number of iteration steps *iter*.
%
% **Note:** If it is possible to explicitly calculate the canonical dual
% frame then this is usually a much faster method than invoking
% `franaiter`.
%
% `franaiter` takes the following parameters at the end of the line of
% input arguments:
%
% 'tol',t Stop if relative residual error is less than the
% specified tolerance. Default is 1e-9
%
% 'maxit',n Do at most n iterations.
%
% 'pg' Solve the problem using the Conjugate Gradient
% algorithm. This is the default.
%
% 'pcg' Solve the problem using the Preconditioned Conjugate Gradient
% algorithm.
%
% 'print' Display the progress.
%
% 'quiet' Don't print anything, this is the default.
%
% Examples
% --------
%
% The following example shows how to rectruct a signal without ever
% using the dual frame:::
%
% f=greasy;
% F=frame('dgtreal','gauss',40,60);
% [c,relres,iter]=franaiter(F,f,'tol',1e-14);
% r=frsyn(F,c);
% norm(f-r)/norm(f)
% semilogy(relres);
% title('Conversion rate of the CG algorithm');
% xlabel('No. of iterations');
% ylabel('Relative residual');
%
% See also: frame, frana, frsyn, frsyniter
% AUTHORS: Peter L. S��ndergaard
if nargin<2
error('%s: Too few input parameters.',upper(mfilename));
end;
definput.keyvals.Ls=[];
definput.keyvals.tol=1e-9;
definput.keyvals.maxit=100;
definput.flags.alg={'cg','pcg'};
definput.keyvals.printstep=10;
definput.flags.print={'quiet','print'};
[flags,kv,Ls]=ltfatarghelper({'Ls'},definput,varargin);
% Determine L from the first vector, it must match for all of them.
L=framelength(F,size(f,1));
A=@(x) frsyn(F,frana(F,x));
% An explicit postpad is needed for the pcg algorithm to not fail
f=postpad(f,L);
if flags.do_pcg
d=framediag(F,L);
M=spdiags(d,0,L,L);
[fout,flag,~,iter,relres]=pcg(A,f,kv.tol,kv.maxit,M);
else
[fout,flag,~,iter,relres]=pcg(A,f,kv.tol,kv.maxit);
end;
c=frana(F,fout);
if nargout>1
relres=relres/norm(fout(:));
end;
end