aboutsummaryrefslogtreecommitdiff
path: root/tvdantzig_example.m
diff options
context:
space:
mode:
authorEugeniy E. Mikhailov <evgmik@gmail.com>2021-01-29 16:23:05 -0500
committerEugeniy E. Mikhailov <evgmik@gmail.com>2021-01-29 16:23:05 -0500
commit3983eb46023c1edd00617729ba929057fda8d0bd (patch)
tree816ad084f355000656c43da9160f1c257bbb1ddc /tvdantzig_example.m
downloadl1magic-3983eb46023c1edd00617729ba929057fda8d0bd.tar.gz
l1magic-3983eb46023c1edd00617729ba929057fda8d0bd.zip
Initial import from https://statweb.stanford.edu/~candes/software/l1magic/v1.11
Additional Clean up of Mac dirs and tex generated files
Diffstat (limited to 'tvdantzig_example.m')
-rw-r--r--tvdantzig_example.m72
1 files changed, 72 insertions, 0 deletions
diff --git a/tvdantzig_example.m b/tvdantzig_example.m
new file mode 100644
index 0000000..f169594
--- /dev/null
+++ b/tvdantzig_example.m
@@ -0,0 +1,72 @@
+% tvdantzig_example.m
+%
+% Test out tvdantzig code (TV minimization with bounded residual correlation).
+%
+% Written by: Justin Romberg, Caltech
+% Email: jrom@acm.caltech.edu
+% Created: October 2005
+%
+
+% use implicit, matrix-free algorithms ?
+largescale = 0;
+
+path(path, './Optimization');
+path(path, './Measurements');
+path(path, './Data');
+
+
+% test image = 32x32 piece of cameraman's arm
+load camera
+I = camera(81:112,37:68);
+n = 32;
+N = n*n;
+I = I/norm(I(:));
+I = I - mean(I(:));
+x = reshape(I,N,1);
+
+
+% num obs
+K = 300;
+
+% permutation P and observation set OMEGA
+P = randperm(N)';
+q = randperm(N/2-1)+1;
+OMEGA = q(1:K/2)';
+
+
+
+% measurement matrix
+if (largescale)
+ A = @(z) A_f(z, OMEGA, P);
+ At = @(z) At_f(z, N, OMEGA, P);
+ % obsevations
+ b = A(x);
+ % initial point
+ x0 = At(b);
+else
+ FT = 1/sqrt(N)*fft(eye(N));
+ A = sqrt(2)*[real(FT(OMEGA,:)); imag(FT(OMEGA,:))];
+ A = [1/sqrt(N)*ones(1,N); A];
+ At = [];
+ % observations
+ b = A*x;
+ % initial point
+ x0 = A'*b;
+end
+
+
+epsilon = 5e-3;
+
+
+
+tvI = sum(sum(sqrt([diff(I,1,2) zeros(n,1)].^2 + [diff(I,1,1); zeros(1,n)].^2 )));
+disp(sprintf('Original TV = %.3f', tvI));
+
+time0 = clock;
+xp = tvdantzig_logbarrier(x0, A, At, b, epsilon, 1e-3, 5, 1e-8, 1500);
+Ip = reshape(xp, n, n);
+disp(sprintf('Total elapsed time = %f secs\n', etime(clock,time0)));
+
+
+
+