Идентификация частотного диапазона: оценка моделей Используя данные о частотном диапазоне
В этом примере показано, как оценить модели с помощью данных о частотном диапазоне. Оценка и валидация моделей с помощью данных о частотном диапазоне работают одинаково, как они делают с данными об области времени. Это обеспечивает большое количество гибкости по оценке и анализу моделей, использующих время и частотный диапазон, а также спектральные данные (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.418 -1.785
B =
u1
x1 -8.3
x2 27.17
C =
x1 x2
y1 0.985 0.3947
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.311e-07
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
Используя данные частотного диапазона непрерывного времени, чтобы оценить модели непрерывного времени
Данные об области времени могут естественно только храниться и имели дело с как дискретное время, выборочные данные. Данные о частотном диапазоне имеют преимущество, что непрерывные данные времени могут быть представлены правильно. Предположим, что базовые непрерывные сигналы времени не имеют никакой информации о частоте выше частоты Найквиста, e.g. потому что они производятся быстро, или вход не имеет никакой частотной составляющей выше частоты Найквиста и что данные были собраны из установившегося эксперимента. Затем Дискретные преобразования Фурье (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.01831 s + 1.009
F(s) = s^2 + 1.001 s + 0.997
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')
Заключения
Мы видели, как время, частота и спектральные данные могут беспрепятственно использоваться, чтобы оценить множество линейных моделей и в областях непрерывного и в дискретного времени. Модели могут быть подтверждены и сравнены в областях, отличающихся от тех, они были оценены в. Форматы данных (время, частота и спектр) являются взаимозаменяемыми, с помощью методов, таких как fftifft, spafdr и spa. Кроме того, прямой, оценка непрерывного времени достижима при помощи tfest, ssest и procest стандартные программы оценки. Бесшовное использование данных в любой области для оценки и анализа является важной функцией System Identification Toolbox.