CODEWITHSUNDEEP
rfft
user3910073
user3910073

Reputation: 531

Bandpassfilter R using fft

I have a time series z with sampling frequeny fs = 12(monthly data) and I would like to perform a bandpass filter using the fftat 10 months and 15 months. This is how I would proceed:

y <- as.data.frame(fft(z))
y$freq <- ..
y$y <-  ifelse(y$freq>= 1/10 & y$freq<= 1/15,y$y,0)
zz <- fft(y$y, inverse = TRUE)/length(z)
plot zz in the time domain...

However, I don't know how to derive the frequencies of the fft and I don't know how to plot zz in the time domain. Can someone help me?

Upvotes: 3

Views: 2344

Answers (1)

Andrey Kolyadin
Andrey Kolyadin

Reputation: 1321

I have a function, that wraps fft() a bit:

    function(y, samp.freq, ...){
      N <- length(y)
      fk <- fft(y)
      fk <- fk[2:length(fk)/2+1]
      fk <- 2*fk[seq(1, length(fk), by = 2)]/N
      freq <- (1:(length(fk)))* samp.freq/(2*length(fk))
      return(data.frame(fur = fk, freq = freq))
    }

y is values of your signal, and samp.freq is it's sample frequency. It's output is data.frame with two columns - fur is complex numbers we get after fast fourier transform (Mod(fur) will be an amplitude, Arg(fur) - a phase) and freq is vector of corresponding frequencies.

But for frequency filtering I highly reccomend using signal package.

For example using Butterworth filter:

     library('signal')
     bf <- butter(2, c(low, high), type = "pass")
     signal.filtered <- filtfilt(bf, signal.noisy)

In this case interval should be defined as c(Low.freq, High.freq) * (2/samp.freq), where Low.freq and High.freq - borders of frequency intervals. More information can be found in package documentation and octave reference guide.

Also, notice that with fft you can get only frequencies up to (sample frequency)/2.

Upvotes: 5

Related Questions