3
$\begingroup$

I have a digital notch filter (band-stop) implemented in C. It works very well. The issue is that varying the filter coefficients in order to track a noise signal that varies in frequency results in a transient response that can be more noisy than the incoming signal. We refer to this as "shot" noise and it appears this is a well-known problem with time-varying notch filters. I have looked at a couple of papers that appear to damp the changes to the coefficients in order to avoid this effect, but the maths is beyond me and I am really just looking for a C example that does something like this.

EDIT: this problem appears worse in our case because we have cascading notches - at least 24 all being updated at the same time - it seems that "shot" noise gets amplified in this case if you get the wrong combination of inputs.

Picture of the frequency changes

EDIT2: An additional problem we are seeing is that the filter cascade becomes unstable at certain points - so we get much larger output than is possible on the input. This only happens when we are linking several notch filters together. This is all in floating point. Is there a way of avoiding this instability?

$\endgroup$
7
  • 1
    $\begingroup$ I'm assuming this is an IIR? The best solution I've found is making a second filter, feeding it with the same input, initialize the internal state with the same as the currently active filter, and sloooowly "blend" over to the new filter. $\endgroup$ Aug 4, 2022 at 16:47
  • $\begingroup$ the alternative of course is not having a hard "new" set of filter coefficients, but "training" a filter live, an adaptive filter. Not sure of any keywords that you could search for in terms of adaptive notch filters, though. $\endgroup$ Aug 4, 2022 at 16:53
  • $\begingroup$ Thanks, helpful - although timely updates are important here so not sure slowly blending is going to work $\endgroup$
    – Andy Piper
    Aug 4, 2022 at 19:14
  • 1
    $\begingroup$ The real issue is how rapidly are the coefficients changing and, more specifically, how they are changing. If you have unlimited computational bandwidth, you can recompute the two key coefficients each and every sample and update them. Then you want to make sure that the notch frequency $f_0$ and $Q$ are not allowed to change too rapidly. This means actually low-pass filtering either the $f_0$ and $Q$ parameters or the coefficients themselves. $\endgroup$ Aug 5, 2022 at 18:39
  • 1
    $\begingroup$ Is it really called shot noise in the larger DSP community, or is that just your in-house name for it? Because in electronics "shot noise" has a specific meaning that's been around for over 100 years and it doesn't describe what you're using the term for! $\endgroup$
    – TimWescott
    Aug 5, 2022 at 18:56

3 Answers 3

7
$\begingroup$

A few ideas:

If you are implementing this as an IIR filter, use Direct Form I. This minimizes the discontinuity of the state variables (which are simply input and output).

A biquad notch filter can be implemented as the sum of unity and a second order allpass, i.e.

$$H(z) = \frac{1}{2}\left(1+ \frac{a_2+a_1z^{-1} + z^{-2}}{1+a_1z^{-1} + a_2 z^{-2}} \right) $$

In this structure only two coefficients need to be updated.

Don't update too fast. Big jumps in frequency will create big jumps in the frequency response. To slow things down you can put a first order lowpass on the notch frequency and adjust the corner frequency to dial in the trade off between responsiveness and artifacts.

implemented in C.

I would first develop, verify and test the algorithm in something that can easily be instrumented like Matlab, Octave or Python. Once you have a stable fully unit tested reference, you can port to C .

$\endgroup$
3
  • $\begingroup$ Thanks, that's helpful. Unfortunately the frequency needs to be updated fast because the noise peak it is tracking changes quickly. To put some flesh on the bone, the sampling rate is 4Khz and this is the rate at which the notch filters are running. The center frequency update rate is 800Hz, and really it needs to be this quick. One solution someone suggested was to modify the lagged signal state to produce the same output after the parameter update as before. It did not seem to help but wondered if this was a valid approach. $\endgroup$
    – Andy Piper
    Aug 6, 2022 at 13:32
  • $\begingroup$ If you use Direct Form I there is NO lagged signal state other than the input and the output (which should not be updated). What's the Q of your notch filter and how big are the frequency jumps, i.e. what is the spectrum of your control signal? $\endgroup$
    – Hilmar
    Aug 6, 2022 at 13:47
  • $\begingroup$ The notches are relatively narrow, so 3db bandwidth about f/15. I'll try and post a picture of the frequency jumps. Generally it works ok but there seems to be some kind of feedback loop that causes it all to go awray. $\endgroup$
    – Andy Piper
    Aug 6, 2022 at 13:56
2
$\begingroup$

After quite a lot of analysis we concluded that the main culprit causing the instability is the rate of change of the center frequency driving the coefficients. We found that a simple slew of:

$$ f'= \begin{cases} f(1-s) & \text{if $f_\lambda < f(1-s)$}\\ f\frac{1}{(1-s)} & \text{if $f_\lambda > f\frac{1}{(1-s)}$}\\ f_\lambda & \text{otherwise} \end{cases} $$

where $f_\lambda$ is the requested frequency, $f$ is the current frequency and $f'$ is the new frequency, worked quite well for reasonable update rates.

$\endgroup$
0
$\begingroup$

In the most simple cases in discrete time, it is enough to assume your noise is white, for designing a linear filter.

When the noise is not white, such as what you describe, an additional term could be useful.

Linear parameter estimation, such as the ARMAX type models, are based on a specific noise probabilistic distribution:

$$ A(z)y(t)=B(z)u(t)+C(z)e(t) $$

Without more information about these spots, it is impossible to say anything else. But with data, you can estimate $C(z)$ and eliminate those noise effects straightforwardly without trouble.

Most discrete time filter design packages usually integrate the 'Moving Average' component, which you can use specifically in this case for this effect. But if that is not the case, or it results too cumbersome to obtain from your package, perhaps you must guide the design manually to address this, or even, to calculate manually a estimate of the noise response $F(z)$, and including that estimate explicitly as a first stage compensator $C(z)=1/F(z)$ in a proper smooth way for this compensator to work, to make the noise $C(z)e(t)$ enough white for the subsequent filter to operate well.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.