16

I agree that the windowing filter design method is not one of the most important design methods anymore, and it might indeed be the case that it is overrepresented in traditional textbooks, probably due to historical reasons. However, I think that its use can be justified in certain situations. I do not agree that computational complexity is no issue ...


7

You can think of linear-least squares in single dimension. The cost function is something like $a^{2}$. The first derivative (Jacobian) is then $2a$, hence linear in $a$. The second derivative (Hessian) is $2$ - a constant. Since the second derivative is positive, you are dealing with convex cost function. This is eqivalent to positive definite Hessian ...


7

The LMS algorithm is based on the idea of gradient descent to search for the optimal (minimum error) condition, with a cost function equal to the mean squared error at the filter output. However, it doesn't actually calculate the gradient directly, as that would require knowing: $$ E(\mathbf{x}[n]e^*[n]) $$ where $\mathbf{x}[n]$ is the input vector and $e[...


7

Windowed Sinc filters can be adaptively generated on the fly on processors barely powerful enough to run the associated FIR filter. Windowed Sinc filters can be generated in finite bounded time. The generation of some simple windowed Sinc filters can be completely described (and inspected for malware, etc.) in a few lines of code, versus blind use of some ...


6

That really depends on context, but generally adaptive implies that the calculations are done on-line / on the fly. In some applications, the filter is updated for a while, then the adaptation is turned off and the last lot of weights are used.


6

If you want to implement the "standard" NLMS algorithm without cutting any corners, then you're probably not going to find a structure that is significantly more efficient. Block forms of LMS filtering aim to use fast convolution techniques (like overlap-save or overlap-add) to speed that part of the process. However, as you noted, the filter coefficients ...


6

To expand on what Peter K. has said, if the signals being used by the filter are stationary, then the filter weights or coefficients can be determined and the filter operates as it was designed without further updates to the filter weights. However, if the signals change, or become quasi-stationary, the filter will adapt continuously.


6

I'll show here one benefit of a windowed design and a trick to get the same benefit from Parks–McClellan. For half-band, quarter-band etc. filters windowing retains the time-domain zeros of the scaled sinc function, which is the prototypical ideal low-pass filter. The zeros end up in the coefficients, reducing the computational cost of the filters. For a ...


6

Slope from all samples obtained To summarize the question's problem, you want to calculate the slope based on all samples obtained thus far, and as new samples are obtained, update the slope without going through all the samples again. On the page you cite is the equation for calculation of the slope $m_n$ that together with $b_n$ minimizes the sum of ...


5

While I completely agree with Jason R's answer, I would like to add a few things that I consider important. First of all, it is a misunderstanding to believe that for least squares designs the transition band width is smaller than for equi-ripple designs. The width of the transition band depends on many design parameters but it is independent of the ...


5

I'm unsure if this is the answer you are looking for, but why not save and share $|H_i|^2$ in addition to the least squares estimates? If you have $w_1^*=\frac 1 {|H_1|^2} H_1^t d_1$, $w_2^*=\frac 1 {|H_2|^2} H_2^t d_2$, and also know $|H_1|^2$ and $|H_2|^2$, you get the total least squares estimate with: $w^*=\frac 1 {|H_1|^2+|H_2|^2} (H_1^t d_1 + H_2^t ...


5

For anyone who is interested, i coincidentally found a paper describing the method implemented in matlab's filtfilt.m. A link to the paper is attached. At least to my understanding matlab's filtfilt.m doesn't implement the Gustafson algorithm. Sadovsky, P.; Bartusek, K: Optimisation of the Transient Response of a Digital Filter, Radioengineering Vol. 9, No. ...


5

It's a bandpass filter, which is fairly obvious. FIR or IIR is hard to tell without phase. A FIR design would likely to have assumed linear phase, and omitted the phase response. and a package probably wouldn't calculate it, like in an IIR design, so I would tend towards a FIR. The band goes to .5 which implies whoever designed the filter was taught unit ...


5

The Jacobian is not computed numerically but analytically and then just evaluated. The frequency response of the IIR filter is $$H(e^{j\omega})=\frac{b_0+b_1e^{-j\omega}+\ldots+b_Me^{-jM\omega}}{1+a_1e^{-j\omega}+\ldots+a_Ne^{-jN\omega}}=\frac{B(e^{j\omega})}{A(e^{j\omega})}\tag{1}$$ Now you need the derivative with respect to the filter coefficients: $$\...


4

First of all, the minimum norm least square solution is $A^+b$, where $A^+$ is the pseudoinverse. Only when the left inverse $A_L^{-1}$ (or right inverse $A_R^{-1}$) exists, you have $A^+=A_L^{-1}=(A^TA)^{-1}A^T$ (or $A^+=A_R^{-1}=A^T(AA^T)^{-1}$). Likewise, the projection matrix onto the column space is $P=AA^+$, or $P=AA_R^{-1}=AA^T(AA^T)^{-1}=I$ if $A_R^{...


4

Here's the way I think about a discrete Wiener Filter Consider a sequence of observations $\mathbf{y} \in \Re^n $ Form a matrix from the input $\mathbf{x} \in \Re^{n+r-1}$ by shifting columns one sample each: $$ X= \begin{bmatrix} x_1 & x_2 & ... & x_r \\ x_2 & x_3 & & x_{r+1} \\ x_3 & x_4 & & x_{r+2} \\ ... & &...


4

Kalman filter is the best linear estimator regardless of stationarity or Gaussianity. Also in the Gaussian case it does not require stationarity (unlike Wiener filter). In the linear Gaussian case Kalman filter is also a MMSE estimator or the conditional mean.


3

I don't understand the subscript $n$ notation, however, in the least squares problem that is given by: \begin{equation} {\bf{y}}={\bf{H}}{\theta}+\bf{n}, \end{equation} where ${\bf{n}}\sim\mathcal{N}(\bf{0}, \sigma^2I_N)$ is a zero mean additive white Gaussian noise and $I_N$ is the $N \times N$ identity matrix, the maximum likelihood and the least squares ...


3

One main difference is the cost function used in the two design methods: Equiripple filters seek to minimize the maximum error between the desired filter response and the designed approximation. Least-squares filters seek to minimize the total squared error betwen the desired filter response and the designed approximation. These different strategies lead ...


3

From here: $$ (AB)(B^{-1}A^{-1}) = A (BB^{-1}) A^{-1} = A A^{-1} = I $$ So $(AB)^{-1} = B^{-1} A^{-1}$ provided $A$ and $B$ are invertible.


3

Our system has the impulse response (why we changed to imaginary instead of real? :) ): $$h(t)=\mathbb{Im}(Ae(iwt)e(-bt))=Ae^{-bt}sin(wt)$$ With the following structure (ref.) (discretized): $$H(s)=\frac{Kw}{(s+b)^2+w^2}$$ Under matlab, you only need to use procest: procest(data,'P2U'); Which models the same structure under the form: $$H(s)=\frac{K}{(1+2\...


3

The equation you're trying to solve is $$ \mathbf{y}=\mathbf{X}\mathbf{h}, $$ where $\mathbf{h}$ is your unknown. The matrix $\mathbf{X}$ is going to have a time-shifted structure that reflects the convolution operator. If we assume that the $\mathbf{y}$ vector starts with y(3) i.e. ignores the first two zeroed out elements of y, then the corresponding $\...


3

I think that the term "constraints" is not a very fortunate choice in this context, but what is meant is the number of frequency points that are specified: $$H(\omega_i)\stackrel{!}{=}D(\omega_i),\qquad i=1,2,\ldots,K\tag{1}$$ where $H(\omega_i)$ is the actual frequency response evaluated at frequency $\omega_i$, $D(\omega_i)$ is the desired response at $\...


3

Constant zero sequence is the only finite sequence for which $a_i$ are not uniquely defined. Indeed, $a_i$ are a solution to the equation $$ \begin{bmatrix} R_0 && R_1 && \ldots && R_{k-1} \\ R_1 && R_0 && \ldots && R_{k-2} \\ \vdots && \vdots && \ddots && \vdots \\ R_{k-1} && ...


3

I had a look at your specs, and I designed a few filters to see what is going on. First of all, we shouldn't expect that an IIR filter should perform much better than an FIR filter for the given specifications, because poles are mainly useful if sharp transitions from pass bands to stop bands must be realized, as is the case for frequency selective filters. ...


3

There are really great answers. I will try to give the Sequential Least Squares approach which generalizes to any Linear Model. Sequential Least Squares Model We're after solving the Linear Least Squares model: $$ \arg \min_{\boldsymbol{\theta}} {\left\| H \boldsymbol{\theta} - \boldsymbol{x} \right\|}_{2}^{2} $$ Now imagine that we have new measurement ...


2

You observe the vector $\mathbf{b}$ which is given by a linear transformation of the unknown data $\mathbf{x}$ plus noise (measurement error) $\mathbf{r}$: $$\mathbf{b}=\mathbf{A}\mathbf{x}+\mathbf{r}\tag{1}$$ The error $\mathbf{r}$ is modeled as a random vector with known mean and known covariance. Note that the transformation matrix $\mathbf{A}$ is ...


2

The linear model is $b=Ax+r$, so you can see the errors $r$ as the noise added to the measurements $b$. The covariance $\Sigma_b=E\{(b-\bar{b})(b-\bar{b})^T\}=E\{(b-Ax)(b-Ax)^T\}=E\{rr^T\}$ $E\{rr^T\} \neq E\{bb^T\}$ unless $Ax=0$


2

Their are some problems in the code and same algorithmic issues: Code: You are estimating the wrong frequency: f2=1/((365*24*60)/2); should be f2=1/((365*24*60))/2; Algorithm: Basically you are using a fourier series (http://en.wikipedia.org/wiki/Fourier_series) with those terms corresponding to the known frequencies. But you should be aware of the ...


2

One approach you can take is to try to fit your data to an ARMA model. There are several implementations (as that link suggests). Also a good reference (if you're mathematically minded) is Lennart Ljung's book, System Identification: Theory for the User Most algorithms mentioned in the links use a mean square error criterion. Many algorithms work well (...


Only top voted, non community-wiki answers of a minimum length are eligible