View Raw SPL
/*****************************************************************************
* *
* NOTCHFILT.SPL Copyright (C) 2008 DSP Development Corporation *
* All Rights Reserved *
* *
* Author: Randy Race *
* *
* Synopsis: Filters a series with a notch comb filter *
* *
* Revisions: 14 Nov 2008 RRR Creation *
* *
*****************************************************************************/
#if @HELP_NOTCHFILT
NOTCHFILT
Purpose: Filters a series with a notch comb filter.
Format: NOTCHFILT(s, k, bw, Ab)
s - A series, the input to filter.
k - An integer, the notch filter order. The order is
equal to Fs/Fo where Fs is the sample frequency and
Fo is the fundamental notch filter frequency.
Defaults to 10.
bw - A real, the bandwidth of each notch. The bandwith
is equal to (Fo/(Fs/2))/Q where Q is the Q factor of
the filter. Defaults to K * 2 / 35.
Ab - An optional real, the dB level for bandwidth BW.
Defaults to 10*log10(.5) or -3 dB.
Returns:
The filtered output series.
Example:
L = 10000; Fs = 1800; Fo = 60; Q = 10;
W1: gsin(L, 1/Fs, 12) + gsin(L, 1/Fs, 60) + gsin(L, 1/Fs, 180)
W2: notchfilt(w1, Fs/Fo, (Fo/(Fs/2))/Q)
W3: spectrum(w1, 32*1024)
W4: spectrum(w2, 32*1024)
The series in W1 is the sum of 3 sinusoids with frequencies
of 12, 60 and 180 Hz.
W2 uses a notch filter to attenuate the 60 and 180 Hz components
of W1.
W3 shows the original spectrum and W4 displays the filtered
spectrum.
Remarks:
The filter order K specifies the number of notch frequencies
over the range 0 to Fs. The specific notch frequncies are
located at F(n) = (n * Fs) / K.
The transfer function if the notch comb filter is:
1 - z^-k
H(z) = gain * ----------
1 - a*z^-k
where the gain is determined to normalize the filter output
to unity gain.
See NOTCHCOMB to design a notch comb filter in direct form.
See Also:
FILTEQ
NOTCHCOMB
PEAKCOMB
PEAKFILT
References:
[1] Sophocles J. Orfanidis
Introduction To Signal Processing
Prentice-Hall, 1996
#endif
/* filter a series with a notch comb filter */
notchfilt(s, K, BW, Ab)
{
local b, a;
if (argc < 4)
{
if (argc < 3)
{
if (argc < 2)
{
if (argc < 1)
{
error(sprintf("%s - input series required", __FUNC__));
}
K = 10;
}
BW = K * 2 / 35;
}
Ab = 10 * log10(0.5);
}
(b, a) = notchcomb(K, BW, Ab, rate(s));
s = filteq(b, a, s);
return(s);
}