goertzel

Дискретное преобразование Фурье с алгоритмом Goertzel второго порядка

Описание

пример

dft = goertzel(data) возвращает дискретное преобразование Фурье (ДПФ) входного массива data использование алгоритма Goertzel второго порядка. Если data имеет больше чем одну размерность, затем goertzel действует вдоль первого измерения массива с размером, больше, чем 1.

пример

dft = goertzel(data,findx) возвращает ДПФ для индексов частоты, заданных в findx.

пример

dft = goertzel(data,findx,dim) вычисляет ДПФ по измерению dim. Вводить размерность и использовать значение по умолчанию findx, задайте второй аргумент как пустой, [].

Примеры

свернуть все

Оцените частоты тона, сгенерированного путем нажатия 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, чтобы вычислить дискретное преобразование Фурье (ДПФ) тона. Выберите индексы, соответствующие частотам раньше, генерировали числа 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')

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

Сгенерируйте шумный косинус с частотными составляющими на уровне 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

Figure contains an axes object. The axes object with title Mean Squared Spectrum contains an object of type line.

Сгенерируйте двухканальный сигнал, произведенный на уровне 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

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

Входные параметры

свернуть все

Входной массив в виде вектора, матрицы или N-D массив.

Пример: sin(2*pi*(0:255)/4) задает синусоиду как вектор-строку.

Пример: sin(2*pi*[0.1;0.3]*(0:39))' задает двухканальную синусоиду.

Типы данных: single | double
Поддержка комплексного числа: Да

Индексы частоты в виде вектора. Индексы могут соответствовать целому числу или множителям нецелого числа f s/N, где f s является частотой дискретизации, и N является длиной сигнала.

Типы данных: single | double

Величина для работы, заданная как положительный целый скаляр.

Типы данных: single | double

Выходные аргументы

свернуть все

Дискретное преобразование Фурье, возвращенное как вектор, матрица или N-D массив.

Алгоритмы

Алгоритм Goertzel реализует дискретное преобразование Фурье X (k) как свертка N - входа x точки (n), n = 0, 1, …, N – 1, с импульсной характеристикой

hk(n)=ej2πkej2πkn/Nu(n)ej2πkWNknu(n),

где u (n), модульная последовательность шага, 1 для n ≥ 0 и 0 в противном случае. k не должен быть целым числом. На частоте f = k f s/N, где f s является частотой дискретизации, преобразование имеет значение

X(k)=yk(n)|n=N,

где

yk(n)=m=0Nx(m)hk(nm)

и x (N) = 0. Z-преобразование импульсной характеристики

Hk(z)=(1WNkz1)ej2πk12cos(2πkN)z1+z2,

с этой прямой формой 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.

Расширенные возможности

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

|

Представлено до R2006a