**4.1 Signal composed of a single sinusoid**

Consider first a signal of the form given by eq. (4) with *K* = 1, *f*1 = 0.44681715350529533 and *a***1** = 0.8018794857541801. Fig. 1 shows a plot of *G***(***f***,** *y***)** for *N* = 1024, *M* = 4. Finding the argmax of *G***(***f***,** *y***)** evaluated for 32,768 frequencies between *f***min** and *f***max** (*N***f** = 32) yields an initial value for the frequency and amplitude of *f*1 = 0.446808 and *a***1** = 0.801070+0.026421i. Minimization of *R*[{)*x(f*1)}] starting at *f***1** = 0.446808 yields a final value *f*1 equal to the input frequency *f*1 with error less than 1x10-16 (machine precision) and the amplitude *a*1 through *A*[{) *x(f*1)}] equal to the input amplitude with error less than 4x10-16.

Fig. 1. *G*(*f*, *y*) as a function of frequency for a signal composed of a single sinusoid mixed with *N* = 1024x4. Note the appearance of a single strong peak in the estimator that serves as an excellent starting value for minimizing the functional *R*(*f*) given in eq. (13).

### **4.2 Signal composed of a 20 sinusoids**

The algorithm also works for multiple frequencies. More than 20 independent tests were performed for an input signal composed of 20 independent frequencies randomly chosen between 0 and 1; all frequency components have amplitude of 1. In all tests our algorithm recovered the 20 frequencies to machine precision with a 128x1024 mixing matrix. For test 1,

Applications of the Orthogonal Matching Pursuit/ Nonlinear

Baraniuk, 2010).

less than 1% failure.

**4.3 Signal composed of 2 sinusoids with large dynamic range** 

performance is not found in the presence of noise as discussed below.

**4.5 Dependence on dimensions of the mixing matrix** 

7 are equispaced and correspond to *C* ~ 1.5.

**4.4 Signal composed of 2 sinusoids with closely spaced frequencies** 

Least Squares Algorithm to Compressive Sensing Recovery 179

For signals composed of 2 well separated frequencies but widely different amplitudes in the absence of noise, we recover the amplitude and frequency of the 2 sinusoids when *a*1 = 1 and *a*2 is as small as 10-14 with an 8x1024 mixing matrix. For this case the amplitude and frequency of the large signal are recovered to machine precision while the frequency and amplitude error of the weak signal are 3x10-4 and 1%, respectively. Naturally, such

We have also input a signal with 2 very closely spaced frequencies and unity amplitudes. For frequencies {0.3389, 0.3390} we recover the frequencies to machine precision with a 16x1024 mixing matrix. Smaller values of *M* for the mixing matrix yield one root half way between the two frequencies. For frequencies {0.33899,0.33900} mixed with 16x1024 and 32x1024 matrices the OMP part of our algorithm yields a signal with one frequency at 0.338995 and an amplitude of 1.9996. Attempts to find a second frequency yield a badly conditioned matrix for *U***H***WU* and the inversion required to find the 2nd amplitude in eq. (11) fails. For a 64x1024 mixing matrix OMP finds two separated estimates of the frequencies and this allows NLS determination of both frequencies to an accuracy of a few parts in 105. These results are in contrast to those obtained using the "spectral compressive sensing" algorithms that use "a signal model that inhibits closely spaced sinusoids" (Duarte and

We have investigated the requirements on *M*, the small dimension of the measurement matrix, to recover a signal composed of a small number of sinusoids using the OMP-NLS algorithm. Fig. 4 shows the fraction of failed recoveries as a function of M for a problem in which the signal is composed of 1,3,5, or 7 sinusoids and *N* = 128. For each value of *K* we performed 1000 trials so a failure fraction of 0.1 corresponds to 100 failures. The conventional relation between *K*, *M*, and *N* for recovery is given by *M* **=** *C K* **log(***N/K***)** (Baraniuk, 2007; Candes and Wakin, 2008). From Fig. 4 we see that the curves for *K* = 3,5 and

We have also investigated several different types of the measurement matrix as displayed in Fig. 5. The three curves correspond to three different measurement matrices. For the blue curve the mixing matrix is generated from the sum of random integers drawn from {-1,0,1} plus i times different random integers drawn from {-1,0,1}; for the red curve, complex numbers with the real and imaginary parts given by reals uniformly distributed between -1 and 1 and i times uniformly distributed reals; for the magenta curve, the mixing matrix is generated from randomly chosen -1's and 1's. The magenta curve for a real mixing matrix made from 1's and -1's is inferior to the blue and red curves for the two complex mixing matrices. The differences between the red and blue curves in Fig. 5 appear to be random fluctuations and are in agreement with other CS results that Gaussian and Bernoulli measurement matrices perform equally well (Baraniuk, 2007; Candes and Wakin, 2008). Fig. 6 compares calculations with the weighting matrix given by eq. (7) to calculations with the weighting matrix set to the identity matrix. One can see that the green curve with the weighting matrix set to the identity matrix is significantly worse in the important region of

shown in detail here, the closest frequency pairs in the signal are {0.2663, 0.2689} and {0.7715, 0.7736}, but while signals with nearly the same frequency are difficult cases, here the combined OMP/NLS recovers all the sinusoids to machine precision. Fig. 2 shows the initial calculation of *G***(***f***,** *y***)** for a 128x1024 mixing matrix and 8192 frequency points (*N***<sup>f</sup>** = 8). Note that most, but not all of the frequencies have peaks in the initial scan of *G***(***f***,** *y***)**. Fig. 3 shows *G***(***f***,** *r***19)** during the 20th iteration of the Do loop in the algorithm shown in Table 1. After refining the frequencies by finding the minimum of *R***(***f***)** in (10), the frequency errors are reduced to less than 10-16 and the amplitude errors are reduced to 4x10-14. Our results compare favorably to those obtained using other recovery methods for a test problem with 20 arbitrary frequency complex sinusoids, *N* = 1024, and variable numbers of measurements *M* (Duarte and Baraniuk, 2010).

Fig. 2. The initial calculation of *G*(*f*, *y*) for a signal with 20 input frequencies mixed with a 128 x 1024 matrix. The red dots indicate the input frequencies.

Fig. 3. The next to last calculation of *G*(*f*, *r*19) for a signal with 20 input frequencies mixed with a 128x1024 matrix showing a large peak near the frequency of the only remaining sinusoid.

178 Applications of Digital Signal Processing

shown in detail here, the closest frequency pairs in the signal are {0.2663, 0.2689} and {0.7715, 0.7736}, but while signals with nearly the same frequency are difficult cases, here the combined OMP/NLS recovers all the sinusoids to machine precision. Fig. 2 shows the initial calculation of *G***(***f***,** *y***)** for a 128x1024 mixing matrix and 8192 frequency points (*N***<sup>f</sup>** = 8). Note that most, but not all of the frequencies have peaks in the initial scan of *G***(***f***,** *y***)**. Fig. 3 shows *G***(***f***,** *r***19)** during the 20th iteration of the Do loop in the algorithm shown in Table 1. After refining the frequencies by finding the minimum of *R***(***f***)** in (10), the frequency errors are reduced to less than 10-16 and the amplitude errors are reduced to 4x10-14. Our results compare favorably to those obtained using other recovery methods for a test problem with 20 arbitrary frequency complex sinusoids, *N* = 1024, and variable numbers of measurements

Fig. 2. The initial calculation of *G*(*f*, *y*) for a signal with 20 input frequencies mixed with a

Fig. 3. The next to last calculation of *G*(*f*, *r*19) for a signal with 20 input frequencies mixed with a 128x1024 matrix showing a large peak near the frequency of the only remaining

128 x 1024 matrix. The red dots indicate the input frequencies.

*M* (Duarte and Baraniuk, 2010).

sinusoid.
