Audioread в matlab из wav файла и FFT

Я работаю над Matlab, я хочу выполнить FFT на wav-файле, который я ранее записал на Matlab.

fs = 44100; % Hz
t = 0:1/fs:1; % seconds
f = 600; % Hz

y = sin(2.*pi.*f.*t);

audiowrite('600freq.wav',y,fs)

Именно так я записываю в wav-файл. Теперь к чтению и части FFT: 

[y,Fs] = audioread('600freq.wav');
sound(y)
plot(fft(y))

Это сюжет БПФ, который я получаю:

Введите описание изображения здесь

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

Заранее благодарю за любую помощь.

1 3
matlab audio signal-processing fft

1 ответ:

Это потому, что вы не строите график величины. То, что вы строите, - это коэффициенты, но они имеют комплексное значение. Из-за этого горизонтальная ось является действительной компонентой, а вертикальная ось - мнимой компонентой. Кроме того, когда вы используете sound сама по себе частота дискретизации по умолчанию составляет 8 кГц (8192 Гц, если быть точным), что объясняет, почему ваш звук имеет более низкий шаг. Вы должны использовать частоту дискретизации в качестве второго аргумента в sound, и это дается вам вторым выходом audioread.

Итак, попробуйте поместить abs после вызова fft, а также использовать Fs в sound: 

[y,Fs] = audioread('600freq.wav');
sound(y, Fs);
plot(abs(fft(y)))

Кроме того, приведенный выше код не строит горизонтальную ось должным образом. Если вы хотите сделать это, убедитесь, что Вы fftshift ваши спектры после того, как вы возьмете преобразование Фурье, затем правильно обозначьте свою ось. Если вы хотите определить, что такое каждое горизонтальное значение с точки зрения частоты, этот удивительный пост от Paul R делает фокус: Как получить частоты каждого значения в БПФ?

В принципе, каждое горизонтальное значение в вашем БПФ таково, что: 

F = i * Fs / N

i это номер ячейки, Fs - частота дискретизации и N - Количество точек, которые вы используете для БПФ. F - это интерпретируемая частота компонента, на который вы смотрите.

По умолчанию fft предполагает, что N - это общее число точек в массиве. Для одностороннего FFT, i идет от 0, 1, 2, до floor((N-1)/2) из-за теоремы выборки Найквиста.

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

Мы можем включить это в ваш код здесь: 

[y,Fs] = audioread('600freq.wav');
sound(y, Fs);
F = fftshift(abs(fft(y)));
f = linspace(-Fs/2, Fs/2, numel(y)+1);
f(end) = [];    
plot(f, F);

Горизонтальная ось теперь отражает правильную частоту. каждого компонента, а также вертикальной оси, отражающей величину каждого компонента.

Запустив код генерации звука, который генерирует синусоидальный сигнал на частоте 600 Гц, а затем приведенный выше код для построения спектров, я получаю следующее:

Введите описание изображения здесь

Обратите внимание, что я вставил наконечник инструмента прямо в положительную сторону спектров... и это около 600 Гц!
7