FFT-based FIR filtering
FFT-based FIR filtering using the overlap-add method.
fftfilt(b, x, n = NULL) ## Default S3 method: fftfilt(b, x, n = NULL) ## S3 method for class 'Ma' fftfilt(b, x, n = NULL)
b |
moving average (Ma) coefficients of a FIR filter, specified as a vector. |
x |
the input signal to be filtered. If x is a matrix, its columns are filtered. |
n |
FFT length, specified as a positive integer. The FFT size must be an
even power of 2 and must be greater than or equal to the length of
|
This function combines two important techniques to speed up filtering of long signals, the overlap-add method, and FFT convolution. The overlap-add method is used to break long signals into smaller segments for easier processing or preventing memory problems. FFT convolution uses the overlap-add method together with the Fast Fourier Transform, allowing signals to be convolved by multiplying their frequency spectra. For filter kernels longer than about 64 points, FFT convolution is faster than standard convolution, while producing exactly the same result.
The overlap-add technique works as follows. When an N length signal is
convolved with a filter kernel of length M, the output signal is
N + M - 1 samples long, i.e., the signal is expanded 'to the right'.
The signal is then broken into k smaller segments, and the convolution
of each segment with the f kernel will have a result of length N / k +
M -1. The individual segments are then added together. The rightmost M
- 1 samples overlap with the leftmost M - 1 samples of the next
segment. The overlap-add method produces exactly the same output signal as
direct convolution.
FFT convolution uses the principle that multiplication in the frequency domain corresponds to convolution in the time domain. The input signal is transformed into the frequency domain using the FFT, multiplied by the frequency response of the filter, and then transformed back into the time domain using the inverse FFT. With FFT convolution, the filter kernel can be made very long, with very little penalty in execution time.
The filtered signal, returned as a vector or matrix with the same
dimensions as x.
Kurt Hornik, Kurt.Hornik@wu-wien.ac.at,
adapted by John W. Eaton.
Conversion to R by Tom Short,
adapted by Geert van Boxtel, G.J.M.vanBoxtel@gmail.com.
t <- seq(0, 1, len = 10000) # 1 second sample x <- sin(2* pi * t * 2.3) + 0.25 * rnorm(length(t)) # 2.3 Hz sinusoid+noise filt <- rep(0.1, 10) # filter kernel y1 <- filter(filt, 1, x) # use normal convolution y2 <- fftfilt(filt, x) # FFT convolution plot(t, x, type = "l") lines(t, y1, col = "red") lines(t, y2, col = "blue") ## use 'filter' with different classes t <- seq(0, 1, len = 10000) # 1 second sample x <- sin(2* pi * t * 2.3) + 0.25 * rnorm(length(t)) # 2.3 Hz sinusoid+noise ma <- Ma(rep(0.1, 10)) # filter kernel y1 <- filter(ma, x) # convulution filter y2 <- fftfilt(ma, x) # FFT filter all.equal(y1, y2) # same result
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.