Более точная частота от FFT с чисто синусоидальными тонами

Question

Я в настоящее время использую код FFT здесь: https://github.com/syedhali/EZAudio/tree/master/EZAudioExamples/iOS/EZAudioFFTExample

Вот код из 2 соответствующих методов:

-(void)createFFTWithBufferSize:(float)bufferSize withAudioData:(float*)data { 

    // Setup the length 
    _log2n = log2f(bufferSize); 

    // Calculate the weights array. This is a one-off operation. 
    _FFTSetup = vDSP_create_fftsetup(_log2n, FFT_RADIX2); 

    // For an FFT, numSamples must be a power of 2, i.e. is always even 
    int nOver2 = bufferSize/2; 

    // Populate *window with the values for a hamming window function 
    float *window = (float *)malloc(sizeof(float)*bufferSize); 
    vDSP_hamm_window(window, bufferSize, 0); 
    // Window the samples 
    vDSP_vmul(data, 1, window, 1, data, 1, bufferSize); 
    free(window); 

    // Define complex buffer 
_A.realp = (float *) malloc(nOver2*sizeof(float)); 
_A.imagp = (float *) malloc(nOver2*sizeof(float)); 

} 

-(void)updateFFTWithBufferSize:(float)bufferSize withAudioData:(float*)data { 

    // For an FFT, numSamples must be a power of 2, i.e. is always even 
    int nOver2 = bufferSize/2; 

    // Pack samples: 
    // C(re) -> A[n], C(im) -> A[n+1] 
    vDSP_ctoz((COMPLEX*)data, 2, &_A, 1, nOver2); 

    // Perform a forward FFT using fftSetup and A 
    // Results are returned in A 
    vDSP_fft_zrip(_FFTSetup, &_A, 1, _log2n, FFT_FORWARD); 

    // Convert COMPLEX_SPLIT A result to magnitudes 
    float amp[nOver2]; 
    float maxMag = 0; 

    for(int i=0; i<nOver2; i++) { 
    // Calculate the magnitude 
    float mag = _A.realp[i]*_A.realp[i]+_A.imagp[i]*_A.imagp[i]; 
    maxMag = mag > maxMag ? mag : maxMag; 
    } 
    for(int i=0; i<nOver2; i++) { 
    // Calculate the magnitude 
    float mag = _A.realp[i]*_A.realp[i]+_A.imagp[i]*_A.imagp[i]; 
    // Bind the value to be less than 1.0 to fit in the graph 
    amp[i] = [EZAudio MAP:mag leftMin:0.0 leftMax:maxMag rightMin:0.0 rightMax:1.0]; 
} 

Я модифицировал updateFFTWithBufferSize Описанный выше метод, так что я мог бы получить частоту в Гц, как это:

for(int i=0; i<nOver2; i++) { 
    // Calculate the magnitude 
    float mag = _A.realp[i]*_A.realp[i]+_A.imagp[i]*_A.imagp[i]; 
    if(maxMag < mag) { 
     _i_max = i; 
    } 
    maxMag = mag > maxMag ? mag : maxMag; 
} 

float frequency = _i_max/bufferSize * 44100; 
NSLog(@"FREQUENCY: %f", frequency); 

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

Например: Синусоидальный сигнал, генерируемый на 19255 Гц, будет отображаться с FFT как 19293.750000Hz. Таким образом, синусовый тон генерируется на частоте 19330 Гц. В расчетах должно быть что-то не так.

Признается, что любая помощь в том, как я могу изменить приведенный выше код, чтобы получить более точное считывание частоты FFT для чистых синусоидов. Спасибо!

Answer 1

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

Лучшая оценка может быть создана с помощью преобразования вашего окна FFT в качестве ядра интерполяции и выполнения последовательной аппроксимации для уточнения оценки максимумов интерполированных точек. (Нулевое заполнение и использование гораздо более длинного БПФ даст вам аналогичный тип интерполированной оценки.)

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

Answer 2

Возможно, ваша частота дискретизации не достаточно высокая?

Answer 3

У вас есть целый ряд проблем, происходящих здесь:

1) Ваша частота разнос оси е _макс/N, или около 80Гц, так что вы не собираетесь, чтобы получить разрешение гораздо лучше чем это.

2) Вы сигнал очень близко к Nyquist frequency (т.е. 20 кГц/44,1 кГц почти 0,5), и когда вы приближаетесь к пределу Найквиста, вам нужно быть очень осторожным, если вы хотите получить точные результаты , (Т. Е. На частоте 20 кГц вы записываете только две точки данных для каждого полного цикла колебаний.)

3) Поскольку 20 кГц находится на грани человеческого слуха (и выше для большинства людей), многие микрофоны надеются, действительно беспокоиться об этом. Here's измерение для iPhone.

Answer 4

FFT - очень хороший способ получить спектр, если вы ничего не знаете о вводе. Если вы знаете, что вход представляет собой чистую синусоидальную волну, вы можете сделать гораздо лучше. Начните с вычисления FFT, чтобы получить приблизительную идею, где синус. Получите минимум и максимум, чтобы оценить амплитуду [или получить от FFT - квадрат все входы, добавить их, взять квадратный корень], получить фазу в начале и конце с учетом предполагаемой частоты и амплитуды.

В общем, вы обнаружите, что фаза не соответствует. Это потому, что фаза на конце отключена на 2 * Δ f * N.f - Δ f - лучшая оценка частоты. Имейте в виду, что такой метод является сверхчувствительным. Метод работает, потому что вход представляет собой чистую синусоидальную волну, а шум - это все, кроме этого. Использование этого метода итеративно взрывается быстро; вы даже сталкиваетесь с ошибками округления (не синусоидальными)

Другой подобный трюк - это вычитание оцененной волны. Разница между двумя синусами - это результат двух синусов: один с добавленными частотами (в вашем случае ± 38,5 кГц) и один с вычитаемыми частотами (Δ_f_, менее 100 Гц). См. Также Heterodyne detection