Идентифицируйте стационарные параметры фильтра от данных о частотной характеристике
Преобразуйте простую передаточную функцию в данные частотной характеристики и затем обратно в исходные коэффициенты фильтра.
a = [1 2 3 2 1 4]; b = [1 2 3 2 3]; [h,w] = freqs(b,a,64); [bb,aa] = invfreqs(h,w,4,5)
bb = 1×5
1.0000 2.0000 3.0000 2.0000 3.0000
aa = 1×6
1.0000 2.0000 3.0000 2.0000 1.0000 4.0000
bb и aa эквивалентны b и a, соответственно. Однако система нестабильна потому что aa имеет полюса с положительной действительной частью. Просмотрите полюса bb и aa.
zplane(bb,aa)
Используйте итеративный алгоритм invfreqs найти устойчивое приближение к системе.
[bbb,aaa] = invfreqs(h,w,4,5,[],30)
bbb = 1×5
0.6816 2.1015 2.6694 0.9113 -0.1218
aaa = 1×6
1.0000 3.4676 7.4060 6.2102 2.5413 0.0001
Проверьте, что система устойчива путем графического вывода новых полюсов.
zplane(bbb,aaa)
Сгенерируйте два вектора, mag и phase, это симулирует величину и данные о фазе, собранные в лаборатории. Также сгенерируйте вектор, w, из частот.
rng('default')
fs = 1000;
t = 0:1/fs:2;
mag = periodogram(sin(2*pi*100*t)+randn(size(t))/10,[],[],fs);
phase = randn(size(mag))/10;
w = linspace(0,fs/2,length(mag))';Используйте invfreqs преобразовывать данные в передаточную функцию непрерывного времени. Постройте результат.
[b,a] = invfreqs(mag.*exp(1j*phase),w,2,2,[],4); freqs(b,a)
При создавании моделей высшего порядка с помощью высоких частот важно масштабировать частоты, делящиеся на фактор, такие как половина самой высокой частоты, существующей в w, чтобы получить хорошо подготовленные значения a и b. Это соответствует перемасштабированию времени.
По умолчанию, invfreqs использует ошибочный метод уравнения, чтобы идентифицировать лучшую модель из данных. Это находит b и a \in
путем создания системы линейных уравнений и решения их с MATLAB®
\ оператор. Здесь A (w (k)) и B (w (k)) является преобразованиями Фурье полиномов a и b, соответственно, на частоте w (k) и n является количеством точек частоты (длина h и w). Этот алгоритм основан на Леви [1]. Несколько вариантов были предложены в литературе, где функция взвешивания wt уделяет меньше внимания высоким частотам.
Начальник (“ошибка на выходе”) алгоритм использует ослабленный метод Ньютона Гаусса для итеративного поиска [2] с выходом первого алгоритма как первоначальная оценка. Это решает прямую задачу минимизации взвешенной суммы квадратичной невязки между фактическим и желаемыми точками частотной характеристики.
[1] Леви, E. C. “Комплексный Curve Fitting”. Сделка IRE на Автоматическом управлении. Издание AC-4, 1959, стр 37–44.
[2] Деннис, J. E. младший, и Р. Б. Шнабель. Численные методы для оптимизации без ограничений и нелинейных уравнений. Englewood Cliffs, NJ: Prentice Hall, 1983.