goertzel
Дискретное преобразование Фурье с алгоритмом Goertzel второго порядка
Описание
Примеры
Оцените телефонные частоты клавиатуры
Оцените частоты тона, сгенерированного путем нажатия 1 кнопки на телефонной клавиатуре.
Номер 1 производит тон с частотами 697 и 1 209 Гц. Сгенерируйте 205 выборок тона с частотой дискретизации 8 кГц.
Fs = 8000; N = 205; lo = sin(2*pi*697*(0:N-1)/Fs); hi = sin(2*pi*1209*(0:N-1)/Fs); data = lo + hi;
Используйте алгоритм Goertzel, чтобы вычислить дискретное преобразование Фурье (DFT) тона. Выберите индексы, соответствующие частотам раньше, генерировали числа 0 до 9.
f = [697 770 852 941 1209 1336 1477]; freq_indices = round(f/Fs*N) + 1; dft_data = goertzel(data,freq_indices);
Постройте величину ДПФ.
stem(f,abs(dft_data)) ax = gca; ax.XTick = f; xlabel('Frequency (Hz)') ylabel('DFT Magnitude')
Разрешите частотные составляющие шумного тона
Сгенерируйте шумный косинус с частотными составляющими на уровне 1,24 кГц, 1,26 кГц и 10 кГц. Задайте частоту дискретизации 32 кГц. Сбросьте генератор случайных чисел для восстанавливаемых результатов.
rng default Fs = 32e3; t = 0:1/Fs:2.96; x = cos(2*pi*t*10e3) + cos(2*pi*t*1.24e3) + cos(2*pi*t*1.26e3) ... + randn(size(t));
Сгенерируйте вектор частоты. Использование алгоритм Goertzel, чтобы вычислить ДПФ ограничивает область значений частот к между 1.2 и 1,3 кГц.
N = (length(x)+1)/2; f = (Fs/2)/N*(0:N-1); indxs = find(f>1.2e3 & f<1.3e3); X = goertzel(x,indxs);
Постройте среднеквадратический спектр, выраженный в децибелах.
plot(f(indxs)/1e3,mag2db(abs(X)/length(X))) title('Mean Squared Spectrum') xlabel('Frequency (kHz)') ylabel('Power (dB)') grid
Дискретное преобразование Фурье массива N-D
Сгенерируйте двухканальный сигнал, произведенный на уровне 3,2 кГц в течение 10 секунд и встроенный в белый Гауссов шум. Первый канал сигнала является синусоидой на 124 Гц. Второй канал является комплексной экпонентой с частотой 126 Гц. Измените сигнал в 3D массив, таким образом, что ось времени запускается по третьему измерению.
fs = 3.2e3; t = 0:1/fs:10-1/fs; x = [cos(2*pi*t*124);exp(2j*pi*t*126)] + randn(2,length(t))/100; x = permute(x,[3 1 2]); size(x)
ans = 1×3
1 2 32000
Вычислите дискретное преобразование Фурье сигнала с помощью алгоритма Goertzel. Ограничьте область значений частот к между 120 Гц и 130 Гц.
N = (length(x)+1)/2; f = (fs/2)/N*(0:N-1); indxs = find(f>=120 & f<=130); X = goertzel(x,indxs,3);
Постройте величину дискретного преобразования Фурье, выраженного в децибелах.
plot(f(indxs),mag2db(abs(squeeze(X)))) xlabel('Frequency (Hz)') ylabel('DFT Magnitude (dB)') grid
Входные параметры
Выходные аргументы
Алгоритмы
Алгоритм Goertzel реализует дискретное преобразование Фурье X (k) как свертка N - входа x точки (n), n = 0, 1, …, N – 1, с импульсной характеристикой
где u (n), модульная последовательность шага, 1 для n ≥ 0 и 0 в противном случае. k не должен быть целым числом. На частоте f = k f s/N, где f s является частотой дискретизации, преобразование имеет значение
где
и x (N) = 0. Z-преобразование импульсной характеристики
с этой прямой формой II реализаций:
Сравните выход goertzel к результату прямой реализации алгоритма Goertzel. Для входного сигнала используйте щебет, произведенный на уровне 50 Гц в течение 10 секунд и встроенный в белый Гауссов шум. Частота щебета увеличивается линейно с 15 Гц до 20 Гц во время измерения. Вычислите дискретное преобразование Фурье на частоте, которая не является целочисленным кратным f s/N. При вызове goertzel, имейте в виду что векторы MATLAB®, запущенные от 1 до N вместо от 0 до N – 1. Результаты соглашаются на высокую точность.
fs = 50; t = 0:1/fs:10-1/fs; N = length(t); xn = chirp(t,15,t(end),20)+randn(1,N)/100; f0 = 17.36; k = N*f0/fs; ykn = filter([1 -exp(-2j*pi*k/N)],[1 -2*cos(2*pi*k/N) 1],[xn 0]); Xk = exp(-2j*pi*k)*ykn(end) dft = goertzel(xn,k+1) df = abs(Xk-dft)
df = 4.3634e-12
Альтернативы
Можно также вычислить ДПФ с:
fft: менее эффективный, чем алгоритм Goertzel, когда вам только нужен ДПФ на нескольких частотах.fftболее эффективно, чемgoertzelкогда необходимо оценить преобразование на больше, чем log2N частотах, где N является длиной входного сигнала.czt:cztвычисляет Z-преобразование ЛЧМ входного сигнала на проспекте, или спираль очерчивают, и включает ДПФ как особый случай.
Ссылки
[1] Burrus, К. Сидни и Томас В. Парки. DFT/FFT и Алгоритмы Свертки: Теория и Реализация. Нью-Йорк: John Wiley & Sons, 1985.
[2] Proakis, Джон Г. и Димитрис Г. Манолакис. Цифровая обработка сигналов: принципы, алгоритмы и приложения. 3-й выпуск. Верхний Сэддл-Ривер, NJ: Prentice Hall, 1996.
[3] Sysel, Петр и Павел Раймик. “Алгоритм Goertzel, Обобщенный ко Множителям Нецелого числа Основной Частоты”. Журнал EURASIP на Усовершенствованиях в Обработке сигналов. Издание 2012, Номер 1, декабрь 2012, стр 56-1–56-8. https://doi.org/10.1186/1687-6180-2012-56.