Преобразования Фурье

Преобразование Фурье является математической формулой, которая связывает сигнал, произведенный вовремя или пробел к тому же сигналу, произведенному в частоте. В обработке сигналов преобразование Фурье может показать важные характеристики сигнала, а именно, его частотных составляющих.

Преобразование Фурье задано для вектора x с n однородно произведенные точки

yk+1=j=0n-1ωjkxj+1.

ω=e-2πi/n один из n комплексные корни из единицы, где i мнимая единица. Для x и y, индексы j и k расположитесь от 0 к n-1.

fft функция в MATLAB® использует алгоритм быстрого преобразования Фурье, чтобы вычислить преобразование Фурье данных. Рассмотрите синусоидальный x сигнала это - функция времени t с частотными составляющими 15 Гц и 20 Гц. Используйте временной вектор, произведенный с шагом 150 из секунды в течение 10 секунд.

Ts = 1/50;
t = 0:Ts:10-Ts;                     
x = sin(2*pi*15*t) + sin(2*pi*20*t);
plot(t,x)
xlabel('Time (seconds)')
ylabel('Amplitude')

Figure contains an axes object. The axes object contains an object of type line.

Вычислите преобразование Фурье сигнала и создайте векторный f это соответствует выборке сигнала в пространстве частоты.

y = fft(x);   
fs = 1/Ts;
f = (0:length(y)-1)*fs/length(y);

Когда вы строите величину сигнала в зависимости от частоты, скачки в величине соответствуют частотным составляющим сигнала 15 Гц и 20 Гц.

plot(f,abs(y))
xlabel('Frequency (Hz)')
ylabel('Magnitude')
title('Magnitude')

Figure contains an axes object. The axes object with title Magnitude contains an object of type line.

Преобразование также производит зеркальную копию скачков, которые соответствуют отрицательным частотам сигнала. Чтобы лучше визуализировать эту периодичность, можно использовать fftshift функция, которая выполняет сосредоточенный нулем, циклический сдвиг на преобразовании.

n = length(x);                         
fshift = (-n/2:n/2-1)*(fs/n);
yshift = fftshift(y);
plot(fshift,abs(yshift))
xlabel('Frequency (Hz)')
ylabel('Magnitude')

Figure contains an axes object. The axes object contains an object of type line.

Сигналы с шумом

В научных приложениях сигналы часто повреждаются со случайным шумом, маскируя их частотные составляющие. Преобразование Фурье может обработать случайный шум и показать частоты. Например, создайте новый сигнал, xnoise, путем введения Гауссова шума в исходный сигнал, x.

rng('default')
xnoise = x + 2.5*randn(size(t));

Степень сигнала в зависимости от частоты является общей метрикой, используемой в обработке сигналов. Степень является величиной в квадрате преобразования Фурье сигнала, нормированного на количество выборок частоты. Вычислите и постройте спектр мощности сигнала с шумом, сосредоточенного на нулевой частоте. Несмотря на шум, можно все еще разобрать частоты сигнала из-за скачков в степени.

ynoise = fft(xnoise);
ynoiseshift = fftshift(ynoise);    
power = abs(ynoiseshift).^2/n; 
plot(fshift,power)
title('Power')
xlabel('Frequency (Hz)')
ylabel('Power')

Figure contains an axes object. The axes object with title Power contains an object of type line.

Вычислительная эффективность

Используя формулу преобразования Фурье непосредственно, чтобы вычислить каждый из n элементы y требует порядка n2 операции с плавающей точкой. Алгоритм быстрого преобразования Фурье требует только порядка nlogn операции, чтобы вычислить. Этот вычислительный КПД является большим преимуществом при обработке данных, которые имеют миллионы точек данных. Много специализированных реализаций алгоритма быстрого преобразования Фурье являются четными более эффективные когда n степень 2.

Считайте аудиоданные собранными из подводных микрофонов недалеко от берегов Калифорнии. Эти данные могут быть найдены в библиотеке, обеспеченной Программой исследований Биоакустики Корнелльского университета. Загрузите и отформатируйте подмножество данных в bluewhale.au, который содержит Тихоокеанскую вокализацию голубого кита. Поскольку вызовы голубого кита являются низкочастотными звуками, они являются едва слышимыми людям. Масштаб времени в данных сжат на коэффициент 10, чтобы повысить тангаж и выполнить более явно слышимый вызов. Можно использовать команду sound(x,fs) слушать целый звуковой файл.

whaleFile = 'bluewhale.au';
[x,fs] = audioread(whaleFile);
whaleMoan = x(2.45e4:3.10e4);
t = 10*(0:1/fs:(length(whaleMoan)-1)/fs);

plot(t,whaleMoan)
xlabel('Time (seconds)')
ylabel('Amplitude')
xlim([0 t(end)])

Figure contains an axes object. The axes object contains an object of type line.

Задайте новую длину сигнала, которая является следующей степенью 2 больших, чем исходная длина. Затем используйте fft вычислить преобразование Фурье с помощью новой длины сигнала. fft автоматически заполняет данные нулями, чтобы увеличить объем выборки. Это дополнение может сделать расчет преобразования значительно быстрее, особенно для объемов выборки с большими простыми множителями.

m = length(whaleMoan); 
n = pow2(nextpow2(m));
y = fft(whaleMoan,n);        

Постройте спектр мощности сигнала. График показывает, что стон состоит из основной частоты приблизительно 17 Гц и последовательности гармоник, где вторая гармоника подчеркнута.

f = (0:n-1)*(fs/n)/10; % frequency vector
power = abs(y).^2/n;   % power spectrum      
plot(f(1:floor(n/2)),power(1:floor(n/2)))
xlabel('Frequency (Hz)')
ylabel('Power')

Figure contains an axes object. The axes object contains an object of type line.

Смотрите также

| | | | | |

Похожие темы