0

I'm trying to implement compressed sensing of an image using the cvx library for MATLAB. This is the same library used by Steve Brunton in his example here. My test image is Lenna.

Here's my MATLAB script:

close all; clear; clc

% read in an image
lenna = imread('lenna.png');
X = lenna(:, :, 2); % greem channel only
X = X(200:249, 100:149); % just a chunk
M = size(X, 1); N = size(X, 2);
x = vectorify(X);

figure, subplot(1, 3, 1), imshow(X), title(strcat("Original Image: ", num2str(M), " by ", num2str(N), " Pixels"))

% Sample randomly
K = 80; %number of random sampled pixels
c = randsample(numel(X), K); %locations of pixels as a list
y = x(c); %values
C = zeros(K, numel(X)); %C(c) = y;%(1:K);
Y = zeros(M, N); % there must be a more elegant way to do this...
for k=1:K
    C(k, c(k)) = y(k);
    Y2 = rectanglefy(C(k, :), M, N);
    Y=Y+Y2;
end
% Y=rectanglefy(C, M, N)
C = C>0; %convert C to a binary matrix
C = double(C); %cast
subplot(1, 3, 2), imshow(Y/255), title(strcat(num2str(K), ' Sampled Pixels'))


% Solve for sparse representation
psi = dftmtx(M*N); psi=real(psi);
Theta = C*psi; n = M*N;

Theta = double(Theta); y = double(y); %cast required
cvx_begin;
    variable s_L1(n)
    minimize( norm(s_L1, 1) );
    subject to
        Theta * s_L1 ==y; 
cvx_end

xr = psi*s_L1;
Xr = rectanglefy(xr, M, N);
subplot(1, 3, 3), imshow(Xr/255), title('Reconstructed Image')

function x = vectorify(X)
    x = reshape(X, [numel(X), 1]);
end

function X = rectanglefy(x, M, N)
    X = reshape(x, [M, N]);
end

When I run it with, say, K=80, it's able to reconstruct an image, although of course it's not very close to the real thing: Image reconstruction of 50x50 pixels, with 80 randomly sampled pixels.

If I up the number of samples to 800, however - which I would expect would make the result closer - the solver tells me it's infeasible: Attempted image reconstruction with 800 randomly sampled pixels, and text from cvx solver in command window stating that solving is infeasible.

In between, if I set K to 180 random samples, the solver yields "Failed," which is evidently a different result from "Infeasible": Attempted image reconstruction with 180 randomly sampled pixels, and text from cvx solver in command window stating that solving failed.

I'm not sure if the problem is in my code or my (rather poor) understanding of convex optimization, but my question is: (a) why is this happening? and (b) what should I be doing differently to reconstruct this image?

Iotatron
  • 65
  • 7

1 Answers1

0

so I'm a bit of an ape without much experience on this stuff but I've tried to replicate your code and I got 2 tips as of now:

1: if you define the cvx part like this you'll get better results (which get better by increasing the sampling points:

epsilon=0.01;
%Using the cvx
cvx_begin;
    variable s_L1(n) complex
    minimize( norm(s_L1, 1) );
    subject to
        norm(y - Theta*s_L1,2) <= epsilon; 
cvx_end

2: You're turning a 2D image into a 1D vector, column by column, which means that your first element and 51st element have a very high correlation. I guess you'd be better off using a 2d fft, like this, but I'm not sure how one would implement that as of now. I'd be happy to share whatever I find if you're still interested in the topic. :)

hope that has helped.

sucksatmath
  • 1