Bandpassfilter R using fft

Question

У меня есть временные ряды z с частотой выборки fs = 12 (ежемесячные данные), и я бы хотел использовать полосовой фильтр, используя fft через 10 месяцев и 15 месяцев. Это, как я бы продолжить:

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... 

Однако, я не знаю, как получить частоты БПФ, и я не знаю, как построить ZZ во временной области. Кто-нибудь может мне помочь?

Answer 1

У меня есть функция, которая надстроена fft() немного:

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 является значением вашего сигнала, и samp.freq является это образец частота. Его вывод data.frame с двумя столбцами - fur - это комплексные числа, которые мы получаем после быстрого преобразования Фурье (Mod(fur) будет амплитудой, Arg(fur) - фазой) и freq является вектором соответствующих частот.

Но для частотной фильтрации я настоятельно рекомендую использовать пакет сигналов.

Например, используя фильтр Баттерворта:

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

В этом случае интервал должен быть определен как с (Low.freq, High.freq) * (2/samp.freq), где Low.freq и высокий .freq - границы частотных интервалов. Более подробную информацию можно найти в документации пакета и octave reference guide.

Также обратите внимание, что с fft вы можете получить только частоты до (частота выборки)/2.