Search This Blog

Everyday DSP for Programmers: DC and Impulsive Noise Removal

For the last installment of Everyday DSP for Programmers, we'll take a look at how to remove unwanted cruft from a signal. Real-world signals can come with all kinds of unwanted noise and biases. A DC bias can wreak havoc on many DSP algorithms because the bias accumulates in the algorithm and saturates the output. The bias doesn't even have to be pure DC to cause problems. The low frequency hum of power lines is pervasive, and most electronic signals ride on top of a 60 Hz sine wave in the US and 50 Hz in Europe. Other data sources can have their own biases that need to be removed.

Another insidious type of noise is impulsive noise. Sometimes called spikes, sparkles, fliers, or outliers, these high-intensity deviations from normal input values can be hard to remove from a signal without affecting the values around them. Averaging tends to spread them out when they should be simply removed from the input signal altogether.

We'll explore how to remove both of these undesirable signal components, starting with DC removal.

DC Removal

The obvious way to remove the DC component of a signal is to take a long moving average of the signal, and subtract the average out of the signal. This works okay, but due to the frequency response of a moving average, some of the low frequencies just above DC will get attenuated as well. Designing a high-pass FIR filter will result in similar issues with the added cost of it being fairly computationally expensive. We can do better.

We can use a simple IIR filter, similar to the complex resonator used for frequency detection to remove the DC bias. Whereas the complex resonator only responded to a narrow frequency band, this filter will only attenuate a narrow frequency band at DC. Here's the formula for it:

w(t) = x(t) + α·w(t-1)
y(t) = w(t) - w(t-1)

Where y(t) is the output of the filter, x(t) is the current input sample, w(t) is an intermediate value that acts like a history of the DC value of the signal, and α is a scale factor that widens or narrows the filter. If α is 1, the filter will pass everything through, and if it's 0, nothing gets through. Neither of those values are very useful, but if α is close to 1, it creates a narrow stop band at the DC frequency.

To implement this formula in JavaScript is fairly straightforward:
function dcRemoval(x, w, alpha) {
  var w_n = x + alpha*w;
  return { w: w_n, y: w_n - w };
}
prev = dcRemoval(signal, prev.w, 0.9);
The function needs to return both the output value and the intermediate value w so that it can be provided to the next function call. In this code the stream of prev.y values constitute the output of the filter, and it looks like the following graph for an input signal that changes its DC value over time (click to start and stop the animation):


The blue signal is a sine wave that rides on top of a DC signal that shifts up or down periodically. The red signal shows the output of the filter, and you can see that the DC signal transitions are slow enough that the filtered signal barely moves before the filter compensates for the change in the DC signal. If the DC transitions were sharper, the filter wouldn't be able to compensate immediately, but it would still return the signal to zero bias pretty quickly. The filter also doesn't attenuate the sine wave at all, even though it's a fairly low frequency sine wave.

DC Removal Frequency Response

We can look at the impulse response of this IIR filter to get a sense of how it looks in the frequency domain. Here is what we get for two different values of α:


The blue signal shows an α of 0.9, and the red signal shows an α of 0.99. We can see that as the impulse runs through the filters, both versions start out passing DC, but fairly quickly start attenuating it. The blue filter stops the DC faster, but it also attenuates higher frequencies more than the red filter. The red filter takes longer to stop the DC values, but has a much narrower stop band. Both filters end up with a very smooth response after running through enough samples. It's a nice, simple filter that can come in handy when you need to strip out DC bias from a high frequency signal.

Impulsive Noise

Impulsive noise can be particularly nasty because once these spikes get into your signal, they are hard to remove by conventional methods without messing up the underlying signal. As an example, say we have a sine wave input, but it has some white noise in it that we want to remove with a moving average filter. Unfortunately, we also have some random periodic spikes in the signal that are coming in for some reason out of our control. Here's what happens to the filtered output:


Not good. The signal gets highly distorted because the moving average spreads each spike out over the signal. The filter is even a smooth FIR filter with a cut-off frequency set at four times the signal frequency. The problem is the spikes have infinite frequency content, so some of their energy passes right through the filter. They are very hard to remove completely.

Impulsive Noise Removal

What we need is a way to detect these impulse spikes and adjust them to be of similar value to the signal around them. Conceptually, this could be done by putting bounds on how much the signal could change in a single sample and detecting when a sample exceeds those bounds with normal samples on either side of it. Here's a neat little algorithm that does essentially that, but with bounds that are dynamically determined from the local behavior of the signal:

y(t) = mean(t..t-N+1) + (Npos - Nneg)*Dpos / N2

Where y(t) is the output, mean(t..t-N+1) is the average of N samples, Npos and Nneg are the number of samples above and below the average, and Dpos is the sum of the differences between the Npos samples and the average (Dneg could also be used, as it should be the same value). If you think about it, an impulse spike will skew the average so that there are fewer samples on its side of the average, and more samples on the opposite side. The algorithm adjusts the average in the opposite direction of the impulse by an amount proportional to the size of the impulse to cancel it out.

This algorithm does require a bit more bookkeeping, but the code is still pretty reasonable:
function impulseRemoval(samples) {
  var N = samples.length;
  var mean = samples.reduce(function(sum, x) { return sum + x }) / N;
  var tallies = samples.reduce(function(tallies, x) {
      if (x > mean) {
        tallies.pos++;
        tallies.diff += x - mean;
      }
      if (x < mean) tallies.neg++;
      return tallies;
  }, { pos: 0, neg: 0, diff: 0 });
  return mean + (tallies.pos - tallies.neg)*tallies.diff / (N*N);
} 
The tallies object is the result of a reduce function that runs through the samples tallying up the samples greater than and less than the mean and accumulating the positive differences from the mean. It may seem like this task could be accomplished in less operations, but after all is said and done, detecting the impulses and adjusting the output will still require a similar amount of work if done a different way. Plus, this method gives us some extra smoothing with the averaging.

Here's what the impulsive noise removal function does to our noisy, sparkly input signal when added just in front of the FIR filter:


As you can see, it doesn't perfectly recover the desired signal, but it does a much better job than the FIR filter alone. If there's no other way to remove impulsive noise upstream, this impulsive noise removal algorithm can be a good option for significantly cleaning up a signal.

So there you have it, two more simple algorithms for removing undesirable bias and noise from your signals. That about wraps up this series on DSP programming. We've covered a lot of ground, but of course, there's plenty more to digital signal processing. If you're interested in learning some of these concepts in more depth or learning some more advanced stuff, I'll give another plug for Richard Lyons' Understanding Digital Signal Processing. Otherwise, I hope you enjoyed learning the basics of DSP, and that you appreciate how much you can do with some simple addition and multiplication. It's hard to believe, but that's most of what DSP algorithms are about.


Other DSP posts:
Basic Signals
Transforms
Sampling Theory
Averaging
Step Response of Averaging
Signal Variance
Edge Detection
Signal Envelopes
Frequency Measurement
The Discrete Fourier Transform
The DFT in Use
Spectral Peak Detection
FIR Filters
Frequency Detection
DC and Impulsive Noise Removal

2 comments:

  1. Hi sam, I am really not sure if you are gonna ever read this, but first I found your posts really interesting, they have help me a lot. Second I have a question, I have been stocked in a project for weeks, and the thing is that I am doing a Digital PI control, and I have to kill DC coming from the ADC signal (reading through SPI). And I have used this filter that you are describing. Problem is that the Phase and attenuation are to big! this is because I am using an Fs of 180KHz and my signal to process has 50Hz. So no metter if I use Alpha=0.995 for the filter, attenuation and Phase will distort my signal. do you have any recommendations for filtering DC in a 50 - 60 Hz signals that are being sampled with a high frequency (180KHz in my case).
    Best regards!
    and again, really good posts!
    Max

    ReplyDelete
    Replies
    1. If all you need to do is remove DC from your signal and the DC component is fixed, (instead of being a very low-frequency signal like in my example) then the method of subtracting the average may work better. If you can tolerate some distortion at the beginning during a warm-up period, then using an accumulating average may work best. The average will be way off when it starts, but it will approach the correct DC value fairly quickly. One issue with this accumulating average is overflow, so you may want to use an exponential average with a really long time constant instead.

      Delete