# Filter design by convex iteration

### From Wikimization

(Difference between revisions)

Line 77: | Line 77: | ||

<pre> | <pre> | ||

N = 32; | N = 32; | ||

- | delta1 = | + | delta1 = 1.01; |

delta2 = 0.01; | delta2 = 0.01; | ||

+ | w = linspace(0,pi,N); | ||

- | cvx_begin | + | cvx_begin |

- | variable h(N,1); | + | variable h(N,1); |

- | + | variable r(2*N-1,1); | |

- | G = | + | variable g(N^2,1); |

- | minimize( | + | variable G(N^2,N^2); |

- | r = | + | variable R(2*N-1,1); |

- | + | ||

- | 1/delta1^2 < | + | for I=1:N |

+ | g((I-1)*N+1:I*N,1) = circshift(h,[0,I-1]); | ||

+ | end | ||

+ | G = g*g'; | ||

+ | % G >= 0; | ||

+ | minimize(max(diag(G))); | ||

+ | for I=1:N | ||

+ | r(I) = 1/I*trace(diag(ones(1,N),I-N)*G); | ||

+ | end | ||

+ | for I=N+1:2*N-1 | ||

+ | r(I) = 1/(I-N)*trace(diag(ones(1,N),I-N)*G); | ||

+ | end | ||

+ | for I=1:N | ||

+ | temp = 0; | ||

+ | for J=1:N-1 | ||

+ | temp = temp + 2*r(J)*cos(w(J)*J); | ||

+ | end | ||

+ | R(I) = r(1) + temp; | ||

+ | end | ||

+ | for I=1:N/2 | ||

+ | 1/delta1^2 < R(I); | ||

+ | R(I) < delta1^2; | ||

+ | end | ||

+ | for I=N/2+1:N | ||

+ | R(I) < delta2^2; | ||

+ | end | ||

+ | |||

+ | R > 0; | ||

+ | |||

+ | |||

cvx_end | cvx_end | ||

</pre> | </pre> |

## Revision as of 05:46, 24 August 2010

where denotes impulse response.

For a low pass filter, frequency domain specifications are:

To minimize peak magnitude of , the problem becomes

But this problem statement is nonconvex.

So instead, a new vector

is defined by concatenation of time-shifted versions of .

Then

is a positive semidefinite rank 1 matrix.

Summing along each of subdiagonals gives entries of the autocorrelation function of .

In particular, the main diagonal of holds squared entries of .

Minimizing is therefore equivalent to minimizing .

Define ...

By spectral factorization, , an equivalent problem is:

Excepting the rank constraint, this problem is convex.

N = 32; delta1 = 1.01; delta2 = 0.01; w = linspace(0,pi,N); cvx_begin variable h(N,1); variable r(2*N-1,1); variable g(N^2,1); variable G(N^2,N^2); variable R(2*N-1,1); for I=1:N g((I-1)*N+1:I*N,1) = circshift(h,[0,I-1]); end G = g*g'; % G >= 0; minimize(max(diag(G))); for I=1:N r(I) = 1/I*trace(diag(ones(1,N),I-N)*G); end for I=N+1:2*N-1 r(I) = 1/(I-N)*trace(diag(ones(1,N),I-N)*G); end for I=1:N temp = 0; for J=1:N-1 temp = temp + 2*r(J)*cos(w(J)*J); end R(I) = r(1) + temp; end for I=1:N/2 1/delta1^2 < R(I); R(I) < delta1^2; end for I=N/2+1:N R(I) < delta2^2; end R > 0; cvx_end