gsignal (version 0.3-1)

fftfilt: FFT-based FIR filtering

Description

FFT-based FIR filtering using the overlap-add method.

Usage

fftfilt(b, x, n = NULL)
# S3 method for default
fftfilt(b, x, n = NULL)
# S3 method for Ma
fftfilt(b, x, n = NULL)

Arguments

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 filt. If the specified n does not meet these criteria, it is automatically adjusted to the nearest value that does. If n = NULL (default), then the overlap-add method is not used.

Value

The filtered signal, returned as a vector or matrix with the same dimensions as x.

Details

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.

References

https://en.wikipedia.org/wiki/Overlap-add_method.

See Also

filter

Examples

# NOT RUN {
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

# }