Skip to content

FIR filters

Linear-phase signal filters

Command Description
FILTER Apply a FIR filter to one or more signals
FILTER-DESIGN Display filter properties

FILTER

Applies a linear-phase FIR filter to a signal

This command modifies the in-memory signal, by applying a finite impulse response (FIR) filter, which can be either a low-pass, high-pass, band-pass or band-stop filter.

Parameters

Core parameters are as follows:

Parameter Example Description
sig sig=C3,F3 Restrict analysis to these channels
bandpass bandpass=0.3,35 Band-pass filter between 0.3 and 35 Hz
lowpass lowpass=35 Low-pass filter with cutoff of 35 Hz
highpass highpass=0.3 High-pass filter between 0.3 and 35 Hz
bandstop bandstop=55,65 Band-stop filter between 0.3 and 35 Hz
fft Use FFT to implement the filter rather than time-domain convolution

By default, the FILTER command uses the Kaiser window approach to define the filter, which requires the two following parameters:

Parameter Example Description
ripple ripple=0.02 Ripple (as a proportion)
tw tw=1 Transition width (in Hz)

Alternatively, it is possible to a) read the FIR coefficients in from a file, or b) use the window method to design a FIR with a fixed filter order, using either a Bartlett, Hann, Blackman or rectangular window:

Parameter Example Description
file file=fir1.txt Read FIR coefficients from a file
order order=30 Fix the FIR filter order
rectangular Specify a rectangular window
bartlett Specify a Bartlett window
hann Specify a Hann window
blackman Specify a Blackman window

Output

After running FILTER, the in-memory signal for a filtered channel will represent the filtered signal. No explicit output is generated by the FILTER command.

Examples

As a first example, here we apply a bandpass filter to the EEG channel C3, with transition frequencies of 0.3 and 35 Hz, using a FIR filter designed via the Kaiser window method with a ripple parameter (ripple) of 0.01 and a transition width parameter (tw) of 0.5 Hz:

FILTER sig=C3 bandpass=0.3,35 ripple=0.01 tw=0.5

In this second example, we use lunaR to copy a signal, apply different bandpass filters to each duplicate (corresponding to delta, theta, alpha, sigma and beta frequency bands) and then plot the results. Here we use the individual nsrr02 from the tutorial dataset. In R, where we assume that the current directory is the one in which we unzipped the tutorial dataset (i.e. and so the s.lst sample-list is present), we'd enter the following:

library(luna)
lattach( lsl( "s.lst" ) , "nsrr02" )

nsrr02 : 14 signals, 10 annotations, 09:57:30 duration
For simplicity, we first drop all signals except EEG using SIGNALS:

leval( "SIGNALS keep=EEG" )

nsrr02 : 1 signals, 10 annotations, 09:57:30 duration

We then duplicate this signal using COPY five times, giving new labels to each:

leval( "COPY sig=EEG tag=DELTA" )
leval( "COPY sig=EEG tag=THETA" )
leval( "COPY sig=EEG tag=ALPHA" )
leval( "COPY sig=EEG tag=SIGMA" )
leval( "COPY sig=EEG tag=BETA" )
By the end, the message in the console reads: 
nsrr02 : 6 signals, 10 annotations, 09:57:30 duration
Using lchs(), we confirm that all channels have been appropriately added: 
lchs()

[1] "EEG"       "EEG_DELTA" "EEG_THETA" "EEG_ALPHA" "EEG_SIGMA" "EEG_BETA"
We then apply the FILTER command to apply a different band-pass filter to each of the channels: 
leval( "FILTER sig=EEG_DELTA bandpass=0.5,4 tw=1 ripple=0.02" )
leval( "FILTER sig=EEG_THETA bandpass=4,8   tw=1 ripple=0.02" )
leval( "FILTER sig=EEG_ALPHA bandpass=8,12  tw=1 ripple=0.02" )
leval( "FILTER sig=EEG_SIGMA bandpass=12,15 tw=1 ripple=0.02" )
leval( "FILTER sig=EEG_BETA  bandpass=15,30 tw=1 ripple=0.02" )

We then epoch the data using lepoch(): 

ne <- lepoch()

nsrr02 : 6 signals, 10 annotations, 09:57:30 duration, 1195 unmasked 30-sec epochs, and 0 masked

We can use the ldata() function to pull out epochs of data for all 6 signals (i.e. as specified by giving lchs() as the second argument), e.g. here for epoch 135: 

d <- ldata( 135 , lchs() )
Looking at the first few rows of the returned data frame d: 
> head(d)
               INT   E      SEC        EEG EEG_DELTA   EEG_THETA  EEG_ALPHA
1 4020.00->4050.00 135 4020.000 -0.4901961 -3.473108 -0.06260304  3.4154174
2 4020.00->4050.00 135 4020.008 -2.4509804 -3.452390 -0.53536473  2.2660953
3 4020.00->4050.00 135 4020.016 -4.4117647 -3.307357 -1.13706143  0.4709635
4 4020.00->4050.00 135 4020.024 -5.3921569 -3.051824 -1.76561958 -1.5923434
5 4020.00->4050.00 135 4020.032 -6.3725490 -2.678883 -2.30284877 -3.4422047
6 4020.00->4050.00 135 4020.040 -7.3529412 -2.188535 -2.65204775 -4.5915269
   EEG_SIGMA   EEG_BETA
1 -1.9357408 -0.2801375
2 -2.7737269  1.0227155
3 -2.3793805  1.6093081
4 -0.9104402  0.8806983
5  0.9774931 -0.6320931
6  2.4513626 -1.3668775

We can plot a data frame in this format with a R's internal functions, such as:

f1 <- function(d) {
par(mfcol=c(6,1),mar=c(0,5,0,0))
plot( d$EEG , type="l" , col=rgb(50,50,50,100,max=255) , axes=F , ylab="Raw" )
plot( d$EEG_DELTA , type="l" , col=rgb(200,100,00,150,max=255) , axes=F , ylab="Delta" )
plot( d$EEG_THETA , type="l" , col=rgb(100,200,0,150,max=255) , axes=F , ylab="Theta" )
plot( d$EEG_ALPHA , type="l" , col=rgb(0,200,100,150,max=255) , axes=F , ylab="Alpha" )
plot( d$EEG_SIGMA , type="l" , col=rgb(0,0,100,150,max=255) , axes=F , ylab="Sigma" )
plot( d$EEG_BETA ,  type="l" , col=rgb(100,0,100,150,max=255) , axes=F , ylab="Beta" )
}

Then, running f1( d ) will produce this plot (although, note that different signals are not uniformly scaled in this particular representation):

img

FILTER-DESIGN

Frequency and impulse responses for FIR filters designed via the Kaiser window method

Filters specified by the FILTER command use the Kaiser window method to design the filter. The FILTER-DESIGN command (or Luna --fir option) can be used to show the properties of these filters.

This command does not depend on any EDFs to be present, and so can be run without a sample-list or EDF (see the example below).

Parameters
Parameter Example Description
fs fs=256 Sampling rate
bandpass bandpass=0.3,35 Band-pass filter between 0.3 and 35 Hz
lowpass lowpass=35 Low-pass filter with cutoff of 35 Hz
highpass highpass=0.3 High-pass filter with cutoff of 0.3 Hz
bandstop bandstop=55,65 Band-stop filter between 0.3 and 35 Hz

The FIR design design approaches are as for the FILTER command: either through the window method (with either a Kaiser window -- tw and ripple -- or fixing the FIR order) or reading from a file:

Parameter Description
ripple Ripple (as a proportion)
tw Transition width (in Hz)
file Read FIR coefficients from a file
order Fix FIR order
rectangular Specify a rectangular window
bartlett Specify a Bartlett window
hann Specify a Hann window
blackman Specify a Blackman window
Output

Per-filter basics (strata: none)

Variable Description
FIR Label for FIR filter (constructed from input parameters)
FS Sampling rate (from fs input parameter)
NTAPS Filter order (number of taps)

Filter coefficients (strata: F x TAP)

Variable Description
TAP Tap
W Filter coefficient

Frequency response characteristics (strata: F x FIR)

Variable Description
F Frequency (Hz)
FIR FIR filter label
MAG Magnitude
MAG_DB Magnitude (dB)
PHASE Phase

Impulse response (strata: FIR x SEC)

Variable Description
F Time (seconds)
FIR FIR filter label
IR Impulse response
Example

Consider a band-pass filter in the sigma band, 11 to 15Hz, applied to a signal with 200 Hz sampling rate. Transition frequencies 11Hz and 15Hz have -6dB attenuation (i.e. 50% amplitude ratio).

If specifying the --fir option, Luna expects the parameters above via standard input, e.g. piped from the echo shell command, with the output still sent to a database out.db:

echo "fs=200 bandpass=11,15 ripple=0.02 tw=1" | luna --fir -o out.db

Examining the out.db, we see that this command produces the following output:

destrat out.db

--------------------------------------------------------------------------------
out.db: 1 command(s), 1 individual(s), 7 variable(s), 4207 values
--------------------------------------------------------------------------------
  command #1:   c1  Mon Feb 25 16:11:40 2019    FIR-DESIGN  
--------------------------------------------------------------------------------
distinct strata group(s):
  commands      : factors           : levels        : variables 
----------------:-------------------:---------------:---------------------------
  [FIR-DESIGN]  : FIR               : 1 level(s)    : FS NTAPS
                :                   :               : 
  [FIR-DESIGN]  : F FIR             : 1025 level(s) : MAG MAG_DB PHASE
                :                   :               : 
  [FIR-DESIGN]  : FIR TAP           : 365 level(s)  : W
                :                   :               : 
  [FIR-DESIGN]  : FIR SEC           : 765 level(s)  : IR
                :                   :               : 
----------------:-------------------:---------------:---------------------------

The first stratum (FIR) gives the sampling rate (FS) as input, and the derived order of the filter (NTAPS):

destrat out.db +FIR-DESIGN -r FIR 

ID    FIR                        FS      NTAPS
.     BANDPASS_11..15_0.02_1     200     365

The actual filter coefficients are in the stratum defined by FIR x TAP: 

destrat out.db +FIR-DESIGN -r FIR TAP > fir.txt

In R, we can plot the filter coefficients: 

w <- read.table("fir.txt",header=T)
plot( w$TAP/200 , w$W , type="l", lwd=2, col="blue", ylab="Coef", xlab="Time" )

img

The amplitude (magnitude) response of the filter is in the F x FIR stratum: 

destrat out.db +FIR-DESIGN -r F FIR > fir2.txt
In R: 
d <- read.table("fir2.txt", header=T)
d <- d[ d$F < 20 , ]
par(mfcol=c(1,2))
plot(d$F, d$MAG   , type="l",lwd=2,col="blue", xlab="Hz",ylab="Magnitude")
plot(d$F, d$MAG_DB, type="l",lwd=2,col="blue", xlab="Hz",ylab="Magnitude (dB)")

img

Back to top