is there a bandwidth limit on the bandstop filter

I have the following code which is stopping a super narrow band on a signal. It runs but won't complete and locks up matlab. I assume becasue I'm pushing the filter past it's limit. Can anyone confirm this?

y = bandstop(signal, [4500 950000], 5e9);

Thanks

[Walter Roberson]

What is size(signal) ?

[STEPHEN BARRETT]

10M samples.

[Walter Roberson]

On my system, with 4 CPUs, that runs 4 threads (automatically) doing the calculation. The amount of memory involved is not much comparatively -- the entire MATLAB process is sitting at less than 1.5 gigabytes including the desktop.

I ran 1 million earlier and it took a bit but finished. I am timing 10M now.

[Walter Roberson]

970.3 seconds for 10M on my system.

Answer 1

[Star Strider]

The easiest way to see what’s going on with it is to use fvtool to look at it, and to give it a short random signal to look at the output:

signal = rand(1, 1E+4);
[y,df] = bandstop(signal, [4500 950000], 5e9);

fvtool(df)

figure
plot(y)

It appears to be a well-designed filter, and works with a random signal. I suspect that something about your ‘signal’ vector could be part of the problem, possibly its amplitude, since the output of the random signal (with an amplitude between 0 and 1) produces an output on the order of . The normalised frequency of the stopband center frequency is about . (In my experience, this is a bit extreme.)

In situations like this, it’s sometimes worthwhile to design both a lowpass filter and a highpass filter, and then add the output vectors. (There is no way to cascade them, since the output of one filter will not overlap with the output of the other filter.)

The problem does not appear to be with the filter itself, since the fvtool analysis suggests that it’s stable and works well.

[STEPHEN BARRETT]

It might be worthwhile to explain more about what I'm actually doing. I need to remove the noise from a signal. There is frequency content in just about everywhere in the spectrum so filtering has been a challenge. A very very small slice is below. The code above attemps to remove the noise in between just two of the peaks which you can see is extremely narrow. i'd need to do this a few thousand times to make it to the end of the spectrum (2.5 GHz). And i can't even get the first to work. The time series vector is basically a series of 200 sets of 5 pulses spaced at 100kHz apart. Each pulse set has a slight decay to it. I'm not even really sure this is the right approach.

[Star Strider]

That very much looks like the Fourier transform of a square pulse. If you want to see if the signal actually represents noise, plot it using semilogy (or dB). My guess is that the information you want is actually the envelope of that Fourier transform. The actual ‘noiise’ may simply be the quantization effects typically due to sampling a continuous-time signal, and the computation of the fft.

[STEPHEN BARRETT]

Here's a slice of the time series from the scope. There's 200 of these peaks in the 10M samples i have. Again, just trying to filter out as much noise everywhere to get a more accurate height measure of each peak. But sounds to me by your post that this noise is unavoidable becasue of quantization effects. Is it possible to filter out those types of effects?

And now zoomed into just one of the peak:

Spectrum of 10M sample signal in db:

Amplitude spectrum

[Star Strider]

I can’t contribute more than I already have. The problem isn’t the filter. The ‘noise’ such as it is, is low-level quantization noise present near the baseline. The signal you want is the upper envelope of the ‘spikes’.

[Walter Roberson]

Try using tic/toc to filter for 1000 samples, 10000 samples, 100000 samples, a million samples... if 1000 works then you know the bandwidth is not the problem but the length of the signal might be. WIth that range of sizes you will be able to interpolate how long 10 million will take on your system

[Star Strider]

If you absolutely must eliminate the ‘spikes’ (that I believe are not noise and are actually what you want), you could always try a FIR ’comb’ filter. See How do I eliminate multiple 60Hz harmonics within 1 set of data in an ECG (link) for an example. It will be a very long filter, so it may not be possible to use it with your data, and it may not be an improvement in terms of speed or efficiency. I mention it simply as an option.

[STEPHEN BARRETT]

Yes, i thought i clarified this. The peaks are certainly what i want to keep. It's all the noise in between the several thousands of peaks i have that i want to get rid of. Does that change anything you said?

[Star Strider]

I wasn’t certain what you intended by ‘noise’. In any event, nothing else changes. The ‘comb’ filter may do what you want. I would also experiment with wavelet denoising if you have the Wavelet Toolbox. See Denoising Signals and Images (link).

[STEPHEN BARRETT]

The denoising works quite well, although I have to admit I don't really know what it's actually doing. It also decreases the peaks height much more than the amplitude of the actual noise in the flat sections. So I'm not really sure how much of my actual signal is getting filtered out too.

[Star Strider]

Part of the problem is that wavelets will find it difficult to differentiate completely the spikes from the noise unless you tell them to. That’s simply the nature of wavelets. I¹ve not worked with wavelets in a while (I would have to back and review my wavelet references), however concentrating on the higher-frequency wavelets (that would eliminate the noise) might be worthwhile, and would not attenuate the low-frequency parts of your signal (spikes). For this, you would use wavelet filters.