ECE 265 Assignment 3
Please work alone on this assignment.
Question 1: Slow, periodic sampling
In this and the next problem, we’re going to investigate compressed sensing using a system where we know a priori that the signal we’re trying to measure has only a few periodic components.
First, let’s generate a signal we know fully – read and run the code “genp1.m” to make a signal with two different frequencies in it.
Under one sampling scheme, uniform temporal sampling, assume you have access to the signal 20% of the time, at times given by “sigMask1”. The available signal is called “available1”.
Q 1a: What is the Nyquist frequency of the system?
Q 1b: Given the convolution theorem, what can you say about the signal you should observe if you multiply the native signal with sigMask1? Plot the Fourier transform of sigMask1: what do you see?
Q 1c: Plot the Fourier transforms of the true signal and the available signal on the same axis: figure; plot(abs(fft(theSig))); hold on; plot(abs(fft(available1)),’r’). What do you observe, and why?
Note that in general, the signal might have complex Fourier amplitudes: a magnitude and a phase. In some senses, this means that a single Fourier component really can have two free parameters. In this assignment, let’s simplify things somewhat: assume the (potentially complex) amplitudes of any signal components have zero imaginary part – this simplifies the bookkeeping and makes it easier to understand the compressed sensing involved, which is the point of the exercise.
Question 2: Slow, randomized sampling
Under a second sampling criterion, the signal samples available are spaced at times given by “sigMask2”. There are the same number of available samples as in Question 1, but these samples are now distributed pseudorandomly on the interval.
Q 2a: Plot the Fourier transform of sigMask2: what do you see? How does it compare with the Fourier transform of sigMask1?
Q 2b: Given the convolution theorem, what can you say about the expected effect of sampling at the times given by sigMask2?
Q 2c: Plot the Fourier transforms of the true signal and the available signal on the same axis: figure; plot(abs(fft(theSig))); hold on; plot(abs(fft(available2)),’r’). What do you observe, and why?
Question 3: Implementing the first iteration of CoSaMP
Now, use an iterative procedure to recover a few of the Fourier components. First, load the A matrix in “A.mat” relating each possible frequency to the observations at the times in sigMask2. It’s normalized by the L2 norm of the columns of A as follows: A = A ./repmat(sqrt(sum(A.^2)),size(A,1),1);
Plot the A matrix, and then confirm that the A matrix does an equivalent operation to the code in “genp1.m” by generating a test signal with your matrix and running (for instance):
testVect = zeros(129,1);
testVect(114) = 4;
testVect(79) = 6;
testSig = A * testVect;
figure; plot(testSig,available2(sigMask2),’.’). Normalizing the A matrix leads to a scaling factor in this xy plot.
Note the first entry I used here is a DC term, followed by 128 potentially independent frequencies, which is the cause of an off-by-one relationship between the entries of testVect and the frequencies (113) and (78) used in genp1.m.
Now, plot the adjoint product u of A’ and the test signal, and note if any peaks are observable. Do a Tikhonov or pseudoinverse reconstruction of the signal; compare the Tikhonov reconstruction with the adjoint product. Why isn’t Tikhonov reconstruction working?
Now, implement one implementation of CoSaMP with k = 5. To do this, find the 5 components of u (the product of A’ and the observations) with the greatest magnitude, and construct the 52 x 5 matrix Ak on the subspace spanned by the 5 chosen components with the largest adjoint product. Use this matrix to do a pseudoinverse; try Tikhonov regularization parameter gamma = 0. What are you able to recover?
Question 4: Finding an unknown signal with two iterations of CoSaMP
Now, there’s another signal with a small number of nonzero unknows, generated by available3 = A*hiddenSig using the same A matrix as in question 3; where “available3.mat” contains the 52 available measurements.
This time, use k = 7 and apply one iteration of cosamp. The first iteration will not return the correct result; a second iteration is needed. To see that there is a residual that needs to be addressed, generate the following figure. Assuming “recon1” is Ak*xk, the data you’d get from the solution after the first CoSaMP iteration, use figure; plot(available3,’r’); hold on; plot(recon1,’b’).
There is a residual, which will need to be fit with components not currently in the active set of k components. To find it, define residual = available3′ – recon1′;, find a new adjoint product with the residual, and plot it. Plot the new adjoint product together with the first adjoint product. What do you notice?
Next, form a new Ak matrix using both the components from the first iteration and the second, and check the reconstruction accuracy now. Which nonzero components are present in the signal?