Идентификация частотного диапазона: оценка моделей Используя данные о частотном диапазоне
Этот пример показывает, как оценить модели с помощью данных о частотном диапазоне. Оценка и валидация моделей с помощью данных о частотном диапазоне работают одинаково, как они делают с данными об области времени. Это обеспечивает большое количество гибкости по оценке и анализу моделей, использующих время и частотный диапазон, а также спектральные данные (FRF). Можно одновременно оценить модели с помощью данных в обеих областях, сравнить и объединить эти модели. Модель, оцененная с помощью данных об области времени, может быть подтверждена с помощью спектральных данных или наоборот.
Данные о частотном диапазоне не могут использоваться для оценки или валидации нелинейных моделей.
Введение
Экспериментальные данные частотного диапазона распространены во многих приложениях. Могло случиться так, что данные были собраны как данные о частотной характеристике (функции частоты: FRF) от процесса с помощью частоты анализатор. Могло также случиться так, что это более практично, чтобы работать с преобразованиями Фурье ввода и вывода (БПФ данных временного интервала), например, обработать периодические или ограниченные данные полосы. (Ограниченный полосой непрерывный сигнал времени не имеет никаких частотных составляющих выше частоты Найквиста). В System Identification Toolbox данные о вводе-выводе частотного диапазона представлены тот же путь как данные временного интервала, т.е. использование объекты iddata. Свойство 'Domain' объекта должно быть установлено в 'Частоту'. Данные о частотной характеристике представлены как комплексные векторы или как векторы значения/фазы как функция частоты. Объекты IDFRD в тулбоксе используются, чтобы инкапсулировать FRFs, где пользователь задает комплексные данные об ответе и вектор частоты. Такой IDDATA или объекты IDFRD (и также объекты FRD Control System Toolbox) могут использоваться беспрепятственно с любой стандартной программой оценки (такой как procest, tfest и т.д.).
Осмотр данных о частотном диапазоне
Давайте начнем путем загрузки некоторых данных о частотном диапазоне:
load demofr
Этот MAT-файл содержит данные о частотной характеристике на частотах W с амплитудным ответом AMP и фазовый отклик PHA. Давайте сначала взглянем на данные:
subplot(211), loglog(W,AMP),title('Amplitude Response') subplot(212), semilogx(W,PHA),title('Phase Response')
Это экспериментальные данные будет теперь храниться как объект IDFRD. Сначала преобразуйте амплитуду, и фаза к комплексу оценила ответ:
zfr = AMP.*exp(1i*PHA*pi/180); Ts = 0.1; gfr = idfrd(zfr,W,Ts);
Ts является шагом расчета базовых данных. Если данные соответствуют непрерывному времени, например, поскольку вход был ограничен полосой, используйте Ts = 0.
Примечание: Если у вас есть Control System Toolbox™, вы могли бы использовать объект FRD вместо объекта IDFRD. IDFRD имеет опции для получения дополнительной информации, как спектры воздействия и меры по неуверенности, которые не доступны в объектах FRD.
Объект IDFRD gfr теперь содержит данные, и это может строиться и анализироваться по-разному. Чтобы просмотреть данные, мы можем использовать plot или bode:
clf
bode(gfr), legend('gfr')
Оценка моделей Используя частотную характеристику (FRF) данные
Чтобы оценить модели, можно теперь использовать gfr в качестве набора данных со всеми командами тулбокса прозрачным способом. Единственное ограничение - то, что шумовые модели не могут быть созданы. Это означает, что для полиномиальных моделей только OE (модели ошибки на выходе) применяются, и для моделей в пространстве состояний, необходимо зафиксировать K = 0.
m1 = oe(gfr,[2 2 1]) % Discrete-time Output error (transfer function) model ms = ssest(gfr) % Continuous-time state-space model with default choice of order mproc = procest(gfr,'P2UDZ') % 2nd-order, continuous-time model with underdamped poles compare(gfr,m1,ms,mproc) L = findobj(gcf,'type','legend'); L.Location = 'southwest'; % move legend to non-overlapping location
m1 =
Discrete-time OE model: y(t) = [B(z)/F(z)]u(t) + e(t)
B(z) = 0.9986 z^-1 + 0.4968 z^-2
F(z) = 1 - 1.499 z^-1 + 0.6998 z^-2
Sample time: 0.1 seconds
Parameterization:
Polynomial orders: nb=2 nf=2 nk=1
Number of free coefficients: 4
Use "polydata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using OE on frequency response data "gfr".
Fit to estimation data: 88.04%
FPE: 0.2512, MSE: 0.2492
ms =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2
x1 -1.785 6.193
x2 -3.417 -1.785
B =
u1
x1 -8.3
x2 27.17
C =
x1 x2
y1 0.9848 0.3948
D =
u1
y1 0
K =
y1
x1 0
x2 0
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: none
Disturbance component: none
Number of free coefficients: 8
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on frequency response data "gfr".
Fit to estimation data: 88.04%
FPE: 0.2512, MSE: 0.2492
mproc =
Process model with transfer function:
1+Tz*s
G(s) = Kp * ---------------------- * exp(-Td*s)
1+2*Zeta*Tw*s+(Tw*s)^2
Kp = 7.4619
Tw = 0.20245
Zeta = 0.36242
Td = 0
Tz = 0.013617
Parameterization:
'P2DUZ'
Number of free coefficients: 5
Use "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using PROCEST on frequency response data "gfr".
Fit to estimation data: 88.03%
FPE: 0.2517, MSE: 0.2492
Как показано выше множества линейных типов модели может быть оценен и в областях непрерывного и в дискретного времени, с помощью спектральных данных. Эти модели могут быть подтверждены с помощью, данные временного интервала. Набор данных ввода-вывода временного интервала ztime, например, собран из той же системы и может использоваться для валидации m1, ms и mproc:
compare(ztime,m1,ms,mproc) %validation in a different domain
Мы можем также посмотреть на невязки, чтобы подтвердить качество модели с помощью данных о валидации ztime. Как наблюдается, невязки являются почти белыми:
resid(ztime,mproc) % Residuals plot
Сжатие данных Используя SPAFDR
Важная причина работать с данными о частотной характеристике состоит в том, что легко уплотнить информацию с небольшой потерей. Команда SPAFDR позволяет вам вычислять сглаживавшие данные об ответе по ограниченным частотам, например, с логарифмическим интервалом. Вот пример, где данные gfr сжаты к 100 логарифмически расположенным с интервалами значениям частоты. С подобным методом также могут быть сжаты исходные данные об области времени:
sgfr = spafdr(gfr) % spectral estimation with frequency-dependent resolution sz = spafdr(ztime); % spectral estimation using time-domain data clf bode(gfr,sgfr,sz) axis([pi/100 10*pi, -272 105]) legend('gfr (raw data)','sgfr','sz','location','southwest')
sgfr = IDFRD model. Contains Frequency Response Data for 1 output(s) and 1 input(s), and the spectra for disturbances at the outputs. Response data and disturbance spectra are available at 100 frequency points, ranging from 0.03142 rad/s to 31.42 rad/s. Sample time: 0.1 seconds Output channels: 'y1' Input channels: 'u1' Status: Estimated using SPAFDR on frequency response data "gfr".
Диаграммы Боде показывают, что информация в сглаживавших данных была взята хорошо забота за. Теперь, эти записи данных с 100 точками могут очень хорошо использоваться для образцовой оценки. Например:
msm = oe(sgfr,[2 2 1]);
compare(ztime,msm,m1) % msm has the same accuracy as M1 (based on 1000 points)
Оценка Используя данные о вводе-выводе частотного диапазона
Может случиться так, что измерения доступны как преобразования Фурье входных параметров и вывода. Такие данные о частотном диапазоне из системы даны как сигналы Y и U. В графиках loglog они похожи
Wfd = (0:500)'*10*pi/500; subplot(211),loglog(Wfd,abs(Y)),title('The amplitude of the output') subplot(212),loglog(Wfd,abs(U)),title('The amplitude of the input')
Данные о частотной характеристике являются по существу отношением между Y и U. Чтобы собрать данные о частотном диапазоне как объект IDDATA, сделайте можно следующим образом:
ZFD = iddata(Y, U, 'Ts', 0.1, 'Frequency', Wfd)
ZFD =
Frequency domain data set with responses at 501 frequencies.
Frequency range: 0 to 31.416 rad/seconds
Sample time: 0.1 seconds
Outputs Unit (if specified)
y1
Inputs Unit (if specified)
u1
Теперь, снова набор данных частотного диапазона ZFD может использоваться в качестве данных во всех стандартных программах оценки, в качестве данных об области времени и данных о частотной характеристике:
mf = ssest(ZFD) % SSEST picks best order in 1:10 range when called this way mfr = ssregest(ZFD) % an alternative regularized reduction based state-space estimator clf compare(ztime,mf,mfr,m1)
mf =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2
x1 -1.78 6.189
x2 -3.406 -1.78
B =
u1
x1 1.32
x2 14.31
C =
x1 x2
y1 2 1.522e-05
D =
u1
y1 0
K =
y1
x1 0
x2 0
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: none
Disturbance component: none
Number of free coefficients: 8
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on frequency domain data "ZFD".
Fit to estimation data: 97.21%
FPE: 0.04288, MSE: 0.04186
mfr =
Discrete-time identified state-space model:
x(t+Ts) = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2 x3 x4 x5 x6
x1 0.7607 0.3671 -0.3322 0.08998 0.01548 0.09474
x2 -0.227 -0.4068 -0.3153 -0.233 -0.1859 0.2795
x3 0.262 -0.3038 0.6404 -0.1926 -0.02235 -0.2056
x4 0.0205 0.1497 -0.03978 -0.1501 0.353 0.3459
x5 0.03394 0.01352 0.2328 0.4938 -0.2778 -0.09414
x6 0.01766 0.3918 0.1908 -0.2817 -0.0382 0.141
x7 -0.03115 -0.3107 -0.02915 0.2754 0.3491 0.5479
x8 0.01892 0.06892 0.09164 0.1246 -0.4316 0.01104
x9 -0.003668 -0.1277 -0.1057 -0.03904 0.1206 -0.423
x10 0.0174 0.04442 -0.02417 -0.1074 0.02164 -0.2264
x7 x8 x9 x10
x1 -0.1456 0.007799 0.05833 -0.1396
x2 -0.336 0.009541 0.2668 -0.1218
x3 0.233 -0.1177 0.03473 -0.06575
x4 0.1299 -0.166 -0.1195 0.2518
x5 -0.6382 -0.2387 0.04806 0.04836
x6 -0.2343 -0.4258 0.6103 -0.1953
x7 -0.2887 -0.6828 0.136 -0.4042
x8 0.5226 -0.2932 -0.007668 0.09615
x9 -0.1939 -0.3768 -0.09409 0.507
x10 -0.02294 -0.6167 -0.4447 -0.5943
B =
u1
x1 -1.588
x2 0.1338
x3 0.09661
x4 -0.02391
x5 0.02324
x6 -0.03246
x7 -0.0002968
x8 0.03677
x9 -0.06961
x10 0.006027
C =
x1 x2 x3 x4 x5 x6 x7
y1 -0.7727 0.5047 -2.644 1.195 0.5607 -1.542 1.08
x8 x9 x10
y1 -0.06959 0.9512 -0.2948
D =
u1
y1 0
K =
y1
x1 0.0306
x2 0.01498
x3 -0.08894
x4 -0.04447
x5 -0.04233
x6 0.01548
x7 -0.01024
x8 0.0004959
x9 0.003125
x10 0.001307
Sample time: 0.1 seconds
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: none
Disturbance component: estimate
Number of free coefficients: 130
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSREGEST on frequency domain data "ZFD".
Fit to estimation data: 76.61% (prediction focus)
FPE: 3.448, MSE: 2.938
Преобразования между представлениями данных (время - частота)
Время и наборы данных ввода - вывода частотного диапазона могут быть преобразованы к любой области при помощи БПФ и ОБПФ. Эти команды адаптируются к объектам IDDATA:
dataf = fft(ztime) datat = ifft(dataf)
dataf =
Frequency domain data set with responses at 501 frequencies.
Frequency range: 0 to 31.416 rad/seconds
Sample time: 0.1 seconds
Outputs Unit (if specified)
y1
Inputs Unit (if specified)
u1
datat =
Time domain data set with 1000 samples.
Sample time: 0.1 seconds
Outputs Unit (if specified)
y1
Inputs Unit (if specified)
u1
Время и данные ввода - вывода частотного диапазона могут быть преобразованы к данным о частотной характеристике SPAFDR, SPA и ETFE:
g1 = spafdr(ztime) g2 = spafdr(ZFD); clf; bode(g1,g2)
g1 = IDFRD model. Contains Frequency Response Data for 1 output(s) and 1 input(s), and the spectra for disturbances at the outputs. Response data and disturbance spectra are available at 100 frequency points, ranging from 0.06283 rad/s to 31.42 rad/s. Sample time: 0.1 seconds Output channels: 'y1' Input channels: 'u1' Status: Estimated using SPAFDR on time domain data "ztime".
Данные о частотной характеристике могут также быть преобразованы к более сглаживавшим данным (меньше разрешения и меньше данных) SPAFDR и SPA;
g3 = spafdr(gfr);
Данные о частотной характеристике могут быть преобразованы к сигналам ввода - вывода частотного диапазона командой IDDATA:
gfd = iddata(g3) plot(gfd)
gfd =
Frequency domain data set with responses at 100 frequencies.
Frequency range: 0.031416 to 31.416 rad/seconds
Sample time: 0.1 seconds
Outputs Unit (if specified)
y1
Inputs Unit (if specified)
u1
Используя непрерывно-разовые данные частотного диапазона, чтобы оценить непрерывно-разовые модели
Данные об области времени могут естественно только храниться и имели дело с как дискретное время, выборочные данные. Данные о частотном диапазоне имеют преимущество, что непрерывные данные времени могут быть представлены правильно. Предположим, что базовые непрерывные сигналы времени не имеют никакой информации о частоте выше частоты Найквиста, например, потому что они выбираются быстро, или вход не имеет никакой частотной составляющей выше частоты Найквиста и что данные были собраны из установившегося эксперимента. Затем Дискретные преобразования Фурье (DFT) данных также являются преобразованиями Фурье непрерывных сигналов времени на выбранных частотах. Они могут поэтому использоваться к непосредственно подходящим непрерывным моделям времени.
Это будет проиллюстрировано следующим примером.
Рассмотрите непрерывную систему времени:
m0 = idpoly(1,1,1,1,[1 1 1],'ts',0)
m0 =
Continuous-time OE model: y(t) = [B(s)/F(s)]u(t) + e(t)
B(s) = 1
F(s) = s^2 + s + 1
Parameterization:
Polynomial orders: nb=1 nf=2 nk=0
Number of free coefficients: 3
Use "polydata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Created by direct construction or transformation. Not estimated.
Загрузите данные, которые прибывают из установившейся симуляции этой системы с помощью периодических входных параметров. Собранные данные были преобразованы в частотный диапазон и сохраненные в файле CTFDData.mat.
load CTFDData.mat dataf % load continuous-time frequency-domain data.
Взгляд на данные:
plot(dataf)
set(gca,'XLim',[0.1 10])
Используя dataf для оценки по умолчанию даст непрерывные модели времени: пространство состояний:
m4 = ssest(dataf,2); %Second order continuous-time model
Поскольку полиномиальная модель с коэффициентом числителя nb = 2 и nf = 2 оценила содействующее использование знаменателя:
nb = 2; nf = 2; m5 = oe(dataf,[nb nf])
m5 =
Continuous-time OE model: y(t) = [B(s)/F(s)]u(t) + e(t)
B(s) = -0.01814 s + 1.008
F(s) = s^2 + 1.001 s + 0.9967
Parameterization:
Polynomial orders: nb=2 nf=2 nk=0
Number of free coefficients: 4
Use "polydata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using OE on frequency domain data "dataf".
Fit to estimation data: 70.15%
FPE: 0.00491, MSE: 0.00468
Сравните переходные процессы с неуверенностью в истинной системе m0 и модели m4 и m5. Доверительные интервалы показывают с закрашенными фигурами в графике.
clf h = stepplot(m0,m4,m5); showConfidence(h,1) legend('show','location','southeast')
Несмотря на то, что это не было необходимо в этом случае, обычно рекомендуется фокусировать подгонку к ограниченному диапазону частот (фильтр нижних частот данные) при оценке использования непрерывных данных времени. Система имеет пропускную способность приблизительно 3 рад/с и была взволнована синусоидами до 6,2 рад/с. Разумный частотный диапазон, чтобы фокусировать подгонку к затем [0 7] rad/s:
m6 = ssest(dataf,2,ssestOptions('WeightingFilter',[0 7])) % state space model
m6 =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2
x1 -0.5011 1.001
x2 -0.7446 -0.5011
B =
u1
x1 -0.01706
x2 1.016
C =
x1 x2
y1 1.001 -0.0005347
D =
u1
y1 0
K =
y1
x1 0
x2 0
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: none
Disturbance component: none
Number of free coefficients: 8
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on frequency domain data "dataf".
Fit to estimation data: 87.03% (data prefiltered)
FPE: 0.004832, MSE: 0.003631
m7 = oe(dataf,[1 2],oeOptions('WeightingFilter',[0 7])) % polynomial model of Output Error structure
m7 =
Continuous-time OE model: y(t) = [B(s)/F(s)]u(t) + e(t)
B(s) = 0.9861
F(s) = s^2 + 0.9498 s + 0.9704
Parameterization:
Polynomial orders: nb=1 nf=2 nk=0
Number of free coefficients: 3
Use "polydata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using OE on frequency domain data "dataf".
Fit to estimation data: 86.81% (data prefiltered)
FPE: 0.004902, MSE: 0.003752
opt = procestOptions('SearchMethod','lsqnonlin',... 'WeightingFilter',[0 7]); % Requires Optimization Toolbox(TM) m8 = procest(dataf,'P2UZ',opt) % process model with underdamped poles
m8 =
Process model with transfer function:
1+Tz*s
G(s) = Kp * ----------------------
1+2*Zeta*Tw*s+(Tw*s)^2
Kp = 1.0124
Tw = 1.0019
Zeta = 0.5021
Tz = -0.017474
Parameterization:
'P2UZ'
Number of free coefficients: 4
Use "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using PROCEST on frequency domain data "dataf".
Fit to estimation data: 87.03% (data prefiltered)
FPE: 0.004832, MSE: 0.003631
opt = tfestOptions('SearchMethod','lsqnonlin',... 'WeightingFilter',[0 7]); % Requires Optimization Toolbox m9 = tfest(dataf,2,opt) % transfer function with 2 poles
m9 = From input "u1" to output "y1": -0.01662 s + 1.007 ------------------ s^2 + s + 0.995 Continuous-time identified transfer function. Parameterization: Number of poles: 2 Number of zeros: 1 Number of free coefficients: 4 Use "tfdata", "getpvec", "getcov" for parameters and their uncertainties. Status: Estimated using TFEST on frequency domain data "dataf". Fit to estimation data: 87.03% (data prefiltered) FPE: 0.00491, MSE: 0.003629
h = stepplot(m0,m6,m7,m8,m9);
showConfidence(h,1)
legend('show')
Заключения
Мы видели, как время, частота и спектральные данные могут беспрепятственно использоваться, чтобы оценить множество линейных моделей и в областях непрерывного и в дискретного времени. Модели могут быть подтверждены и сравнены в областях, отличающихся от тех, они были оценены в. Форматы данных (время, частота и спектр) являются взаимозаменяемыми, с помощью методов, таких как fft, ifft, spafdr и spa. Кроме того, прямая, непрерывно-разовая оценка достижима при помощи tfest, ssest и стандартных программ оценки procest. Бесшовное использование данных в любой области для оценки и анализа является важной функцией System Identification Toolbox.