Как извлечь полуточные частоты из файла WAV с помощью преобразований Фурье

Файлы WAV содержат данные с линейной импульсной кодовой модуляцией (LPCM) . Это просто означает, что это последовательность значений амплитуды с фиксированной частотой дискретизации. Заголовок RIFF содержится в начале файла для передачи такой информации, как частота дискретизации и количество бит на выборку (например, 8 кГц, 16 бит со знаком).

Формат очень простой, и вы легко можете создать свой собственный. Однако для ускорения процесса доступно несколько библиотек, например libsndfile . Simple Direct-media Layer (SDL) / SDL_mixer и PortAudio - две хорошие библиотеки для воспроизведения.

Что касается подачи данных в FFTW, вам необходимо буферизовать 1-секундные фрагменты (размер определяется частотой дискретизации и битами на выборку). Затем преобразуйте все образцы в числа с плавающей запятой IEEE (т.е. float или double в зависимости от конфигурации FFTW - libsndfile может сделать это за вас). Затем создайте еще один массив для хранения выходных данных частотной области. Наконец, создайте и выполните план FFTW, передав оба буфера в fftw_plan_dft_r2c_1d и вызвав fftw_execute с возвращенным дескриптором fftw_plan .

Для получения спектра мощности участка вашего файла:

• соберите N выборок, где N - степень 2 - если ваша частота дискретизации 44,1 кГц, и вы хотите делать выборку примерно каждую секунду, тогда возьмите N = 32768 выборок.

• примените подходящую оконную функцию к выборкам, например. Хэннинга

• передайте образцы с окном в FFT - в идеале вам нужно вещественно-комплексное БПФ, но если у вас есть только комплексно-комплексное БПФ, то передайте 0 для всех образцов. тогда передайте 0 для всех мнимых входных частей

• вычислите квадрат величины выходных бинов БПФ (re * re + im * im)

• (необязательно) вычислите 10 * log10 каждого квадрата величины выходного бина, чтобы получить значение величины в дБ

Теперь, когда у вас есть спектр мощности, вам просто нужно определить пик(ы), что должно быть довольно просто, если у вас разумное отношение S/N. Обратите внимание, что частотное разрешение улучшается с увеличением N. Для примера 44,1 кГц и N = 32768 частотное разрешение каждого бина составляет 44100 / 32768 = 1,35 Гц.

Вас в основном интересует оценка спектра - предполагая, что вы уже прошли этап чтения WAV и преобразования его в сигнал дискретного времени.

Среди различных методов самым простым является метод периодограммы, который сводится к получению оконного дискретного преобразования Фурье (с помощью БПФ) и сохранению его квадрата величины. Это соответствует ответу Пола. Вам нужно окно, которое охватывает несколько периодов самой низкой частоты, которую вы хотите обнаружить. Пример: если ваши синусоиды могут иметь частоту до 10 Гц (период = 100 мс), вы должны взять окно 200 мс или 300 мс или около того (или больше). Однако у периодограммы есть некоторые недостатки, хотя она проста в вычислении и ее более чем достаточно, если не требуется высокая точность:

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

Периодограмма может работать лучше, если усреднить несколько окон с разумным выбором ширины (метод Бартлета). И есть много других методов для оценки спектра (AR моделирование).

На самом деле, вы не совсем заинтересованы в оценке полного спектра, а только местоположения одной частоты. Это можно сделать, ища пик оценочного спектра (делается, как объяснялось), но также и более специфическими и мощными (и сложными) методами (Писаренко, алгоритм MUSIC). В вашем случае они, вероятно, будут излишни.