tfestimate

Оценка передаточной функции

Синтаксис

txy = tfestimate(x,y)
txy = tfestimate(x,y,window)
txy = tfestimate(x,y,window,noverlap)
txy = tfestimate(x,y,window,noverlap,nfft)
txy = tfestimate(___,'mimo')
[txy,w] = tfestimate(___)
[txy,f] = tfestimate(___,fs)
[txy,w] = tfestimate(x,y,window,noverlap,w)
[txy,f] = tfestimate(x,y,window,noverlap,f,fs)
[___] = tfestimate(x,y,___,freqrange)
[___] = tfestimate(___,'Estimator',est)
tfestimate(___)

Описание

txy = tfestimate(x,y) находит оценку передаточной функции, txy, учитывая входной сигнал, x, и выходной сигнал, y.

  • Если x и y являются оба векторами, у них должна быть та же длина.

  • Если один из сигналов является матрицей, и другой вектор, то длина вектора должна равняться количеству строк в матрице. Функция расширяет вектор и возвращает матрицу оценок передаточной функции столбца столбцом.

  • Если x и y являются матрицами с одинаковым числом строк, но различные количества столбцов, то txy является мультивходом/мультивыводом (MIMO) передаточная функция, которая комбинирует все сигналы ввода и вывода. txy является 3D массивом. Если x имеет столбцы m, и y имеет столбцы n, то txy имеет столбцы n и страницы m. Смотрите Передаточную функцию для получения дополнительной информации.

  • Если x и y являются матрицами равного размера, то tfestimate действует по столбцам: txy(:,n) = tfestimate(x(:,n),y(:,n)). Чтобы получить оценку MIMO, добавьте 'mimo' к списку аргументов.

пример

txy = tfestimate(x,y,window) использование window, чтобы разделить x и y в сегменты и выполнить работу с окнами.

txy = tfestimate(x,y,window,noverlap) использование выборки noverlap перекрытия между смежными сегментами.

txy = tfestimate(x,y,window,noverlap,nfft) использование выборка nfft указывает, чтобы вычислить дискретное преобразование Фурье.

txy = tfestimate(___,'mimo') вычисляет передаточную функцию MIMO для матричных входных параметров. Этот синтаксис может включать любую комбинацию входных параметров от предыдущих синтаксисов.

[txy,w] = tfestimate(___) возвращает вектор нормированных частот, w, в котором оценивается передаточная функция.

пример

[txy,f] = tfestimate(___,fs) возвращает вектор частот, f, выраженного с точки зрения частоты дискретизации, fs, в котором оценивается передаточная функция. fs должен быть шестым числовым входом к tfestimate. Чтобы ввести частоту дискретизации и все еще использовать значения по умолчанию предыдущих дополнительных аргументов, задайте эти аргументы как пустые, [].

[txy,w] = tfestimate(x,y,window,noverlap,w) возвращает оценку передаточной функции на нормированных частотах, заданных в w.

[txy,f] = tfestimate(x,y,window,noverlap,f,fs) возвращает оценку передаточной функции на частотах, заданных в f.

[___] = tfestimate(x,y,___,freqrange) возвращает оценку передаточной функции по частотному диапазону, заданному freqrange. Допустимыми опциями для freqrange является 'onesided', 'twosided' и 'centered'.

пример

[___] = tfestimate(___,'Estimator',est) оценочные передаточные функции с помощью средства оценки est. Допустимыми опциями для est является 'H1' и 'H2'.

tfestimate(___) без выходных аргументов строит оценку передаточной функции в окне текущей фигуры.

Примеры

свернуть все

Вычислите и постройте оценку передаточной функции между двумя последовательностями, x и y. Последовательность x состоит из белого Гауссова шума. результаты y фильтрации x с 30-м порядком lowpass фильтруют с нормированной частотой среза 0.2π рад/выборка. Используйте прямоугольное окно, чтобы разработать фильтр. Задайте частоту дискретизации 500 Гц и Окно Хэмминга длины 1024 для оценки передаточной функции.

h = fir1(30,0.2,rectwin(31));
x = randn(16384,1);
y = filter(h,1,x);

fs = 500;
tfestimate(x,y,1024,[],[],fs)

Используйте fvtool, чтобы проверить, что передаточная функция аппроксимирует частотную характеристику фильтра.

fvtool(h,1,'Fs',fs)

Получите тот же результат путем возврата оценки передаточной функции в переменной и графического вывода ее абсолютного значения в децибелах.

[Txy,f] = tfestimate(x,y,1024,[],[],fs);

plot(f,mag2db(abs(Txy)))

Оцените передаточную функцию для простой single-input/single-output системы и сравните ее с определением.

Одномерное дискретное время колеблющаяся система состоит из модульной массы, m, присоединенный к стене к пружине модуля эластичная константа. Датчик выбирает ускорение, a, из массы в Fs=1 Гц. Демпфер препятствует движению массы путем проявления на нем силы, пропорциональной, чтобы ускориться с затуханием постоянного b=0.01.

Сгенерируйте в 2000 раз выборки. Задайте интервал выборки Δt=1/Fs.

Fs = 1;
dt = 1/Fs;
N = 2000;
t = dt*(0:N-1);
b = 0.01;

Система может быть описана моделью в пространстве состояний

x(k+1)=Ax(k)+Bu(k),y(k)=Cx(k)+Du(k),

где x=[rv]T вектор состояния, r и v соответственно положение и скорость массы, u движущая сила, и y=a измеренный вывод. Матрицы пространства состояний

A=exp(AcΔt),B=Ac-1(A-I)Bc,C=[-1-b],D=1,

I 2×2 идентичность и непрерывно-разовые матрицы пространства состояний

Ac=[01-1-b],Bc=[01].

Ac = [0 1;-1 -b];
A = expm(Ac*dt);

Bc = [0;1];
B = Ac\(A-eye(size(A)))*Bc;

C = [-1 -b];
D = 1;

Масса управляется случайным входом для половины интервала измерения. Используйте модель в пространстве состояний, чтобы вычислить эволюцию времени системы, начинающей со все-нулевого начального состояния. Постройте ускорение массы как функция времени.

rng default

u = zeros(1,N);
u(1:N/2) = randn(1,N/2);

y = 0;
x = [0;0];
for k = 1:N
    y(k) = C*x + D*u(k);
    x = A*x + B*u(k);
end

plot(t,y)

Оцените передаточную функцию системы как функция частоты. Используйте точки ДПФ 20:48 и задайте окно Kaiser с форм-фактором 15. Используйте значение по умолчанию перекрытия между смежными сегментами.

nfs = 2048;
wind = kaiser(N,15);

[txy,ft] = tfestimate(u,y,wind,[],nfs,Fs);

Функция частотной характеристики системы дискретного времени может быть выражена как Z-преобразование передаточной функции временного интервала системы, оцененной в модульном кругу. Проверьте, что оценка, вычисленная tfestimate, совпадает с этим определением.

[b,a] = ss2tf(A,B,C,D);

fz = 0:1/nfs:1/2-1/nfs;
z = exp(2j*pi*fz);
frf = polyval(b,z)./polyval(a,z);

plot(ft,20*log10(abs(txy)))
hold on
plot(fz,20*log10(abs(frf)))
hold off
grid
ylim([-60 40])

Постройте оценку с помощью встроенной функциональности tfestimate.

tfestimate(u,y,wind,[],nfs,Fs)

Оцените передаточную функцию для простой системы мультивхода/мультивывода.

Идеальная одномерная колеблющаяся система состоит из двух масс, m1 и m2, ограниченный между двумя стенами. Модули таковы что m1=1 и m2=μ. Каждая масса присоединена к самой близкой стене к пружине с эластичной константой k. Идентичная пружина соединяет эти две массы. Три демпфера препятствуют движению масс путем проявления на них, обеспечивает пропорциональный, чтобы ускориться, с затуханием постоянного b. Выборка датчиков a1 и a2, ускорения масс, в Fs=50 Гц.

Сгенерируйте в 30000 раз выборки, эквивалентные 600 секундам. Задайте интервал выборки Δt=1/Fs.

Fs = 50;
dt = 1/Fs;
N = 30000;
t = dt*(0:N-1);

Система может быть описана моделью в пространстве состояний

x(k+1)=Ax(k)+Bu(k),y(k)=Cx(k)+Du(k),

где x=[r1v1r2v2]T вектор состояния, ri и vi соответственно местоположение и скорость iмасса th, u=[u1u2]T вектор входных движущих сил, и y=[a1a2]T выходной вектор. Матрицы пространства состояний

A=exp(AcΔt),B=Ac-1(A-I)Bc,C=[-2k-2bkbk/μb/μ-2k/μ-2b/μ],D=[1001/μ],

I 4×4 идентичность и непрерывно-разовые матрицы пространства состояний

Ac=[0100-2k-2bkb0001k/μb/μ-2k/μ-2b/μ],Bc=[00100001/μ].

Набор k=400, b=0, и μ=1/10.

k = 400;
b = 0;
m = 1/10;

Ac = [0 1 0 0;-2*k -2*b k b;0 0 0 1;k/m b/m -2*k/m -2*b/m];
A = expm(Ac*dt);
Bc = [0 0;1 0;0 0;0 1/m];
B = Ac\(A-eye(4))*Bc;
C = [-2*k -2*b k b;k/m b/m -2*k/m -2*b/m];
D = [1 0;0 1/m];

Массы управляются случайным входом в течение измерения. Используйте модель в пространстве состояний, чтобы вычислить эволюцию времени системы, начинающей со все-нулевого начального состояния.

rng default
u = randn(2,N);

x = [0;0;0;0];
for kk = 1:N
    y(:,kk) = C*x + D*u(:,kk);
    x = A*x + B*u(:,kk);
end

Используйте входные и выходные данные, чтобы оценить передаточную функцию системы как функция частоты. Задайте опцию 'mimo', чтобы произвести все четыре передаточных функции. Используйте окно Hann с 5000 выборками, чтобы разделить сигналы на сегменты. Задайте 2 500 выборок перекрытия между смежными сегментами и 214 Точки ДПФ. Постройте оценки.

wind = hann(5000);
nov = 2500;

[q,fq] = tfestimate(u',y',wind,nov,2^14,Fs,'mimo');

Вычислите теоретическую передаточную функцию как Z-преобразование передаточной функции временного интервала, оцененной в модульном кругу.

nfs = 2^14;

fz = 0:1/nfs:1/2-1/nfs;
z = exp(2j*pi*fz);

[b1,a1] = ss2tf(A,B,C,D,1);
[b2,a2] = ss2tf(A,B,C,D,2);

frf(1,:,1) = polyval(b1(1,:),z)./polyval(a1,z);
frf(1,:,2) = polyval(b1(2,:),z)./polyval(a1,z);

frf(2,:,1) = polyval(b2(1,:),z)./polyval(a2,z);
frf(2,:,2) = polyval(b2(2,:),z)./polyval(a2,z);

Постройте теоретические передаточные функции и их соответствующие оценки.

for jk = 1:2
    for kj = 1:2
        subplot(2,2,2*(jk-1)+kj)
        plot(fq,20*log10(abs(q(:,jk,kj))))
        hold on
        plot(fz*Fs,20*log10(abs(frf(jk,:,kj))))
        hold off
        grid
        title(['Input ' int2str(kj) ', Output ' int2str(jk)])
        axis([0 Fs/2 -50 100])
    end
end

Передаточные функции имеют максимумы в ожидаемых значениях, ω1,2/2π, где ω собственные значения модальной матрицы.

sqrt(eig(k*[2 -1;-1/m 2/m]))/(2*pi)
ans = 2×1

    3.8470
   14.4259

Добавьте затухание в систему путем установки b=0.1. Вычислите эволюцию времени ослабленной системы с теми же движущими силами. Вычислите H2 оценка передаточной функции MIMO с помощью того же окна и перекрытия. Постройте оценки с помощью функциональности tfestimate.

b = 0.1;

Ac = [0 1 0 0;-2*k -2*b k b;0 0 0 1;k/m b/m -2*k/m -2*b/m];
A = expm(Ac*dt);
B = Ac\(A-eye(4))*Bc;
C = [-2*k -2*b k b;k/m b/m -2*k/m -2*b/m];

x = [0;0;0;0];
for kk = 1:N
    y(:,kk) = C*x + D*u(:,kk);
    x = A*x + B*u(:,kk);
end

clf
tfestimate(u',y',wind,nov,[],Fs,'mimo','Estimator','H2')
legend('I1, O1','I1, O2','I2, O1','I2, O2')

yl = ylim;

Сравните оценки с теоретическими прогнозами.

[b1,a1] = ss2tf(A,B,C,D,1);
[b2,a2] = ss2tf(A,B,C,D,2);

frf(1,:,1) = polyval(b1(1,:),z)./polyval(a1,z);
frf(1,:,2) = polyval(b1(2,:),z)./polyval(a1,z);

frf(2,:,1) = polyval(b2(1,:),z)./polyval(a2,z);
frf(2,:,2) = polyval(b2(2,:),z)./polyval(a2,z);

plot(fz*Fs,20*log10(abs(reshape(permute(frf,[2 1 3]),[nfs/2 4]))))
legend('I1, O1','I1, O2','I2, O1','I2, O2')
ylim(yl)
grid

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

свернуть все

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

Пример: cos(pi/4*(0:159))+randn(1,160) задает синусоиду, встроенную в белый Гауссов шум.

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

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

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

Окно, заданное как целое число или как строка или вектор-столбец. Используйте window, чтобы разделить сигнал на сегменты.

  • Если window является целым числом, то tfestimate делит x и y в сегменты длины window и окна каждый сегмент с Окном Хэмминга той длины.

  • Если window является вектором, то tfestimate делит x и y в сегменты той же длины как вектор и окна каждый сегмент с помощью window.

Если длина x и y не может быть разделена точно в целое число сегментов с noverlap перекрывающиеся выборки, то сигналы являются усеченными соответственно.

Если вы задаете window как пустой, то tfestimate использует Окно Хэмминга, таким образом, что x и y разделены на восемь сегментов с noverlap перекрывающиеся выборки.

Для списка доступных окон смотрите Windows.

Пример: hann(N+1) и (1-cos(2*pi*(0:N)'/N))/2 оба задают окно Hann длины N + 1.

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

Количество перекрытых выборок, заданных как положительное целое число.

  • Если window является скаляром, то noverlap должен быть меньшим, чем window.

  • Если window является вектором, то noverlap должен быть меньшим, чем длина window.

Если вы задаете noverlap как пустой, то tfestimate использует номер, который производит 50%-е перекрытие между сегментами. Если длина сегмента не задана, функция устанавливает noverlap на ⌊N/4.5 ⌋, где N является длиной сигналов ввода и вывода.

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

Количество точек ДПФ, заданных как положительное целое число. Если вы задаете nfft как пустой, то tfestimate устанавливает этот аргумент на макс. (256,2p), где p = ⌈log2  N  для входных сигналов длины N.

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

Частота дискретизации, заданная как положительная скалярная величина. Частота дискретизации является количеством выборок в единицу времени. Если модуль времени является секундами, то частота дискретизации имеет модули Гц.

Нормированные частоты, заданные как строка или вектор-столбец по крайней мере с двумя элементами. Нормированные частоты находятся в раде/выборке.

Пример: w = [pi/4 pi/2]

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

Частоты, заданные как строка или вектор-столбец по крайней мере с двумя элементами. Частоты находятся в циклах в единицу времени. Единица времени задана частотой дискретизации, fs. Если fs имеет модули выборок/секунда, то f имеет модули Гц.

Пример: fs = 1000; f = [100 200]

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

Частотный диапазон для оценки передаточной функции, заданной как та 'onesided', 'twosided' или 'centered'. Значением по умолчанию является 'onesided' для сигналов с действительным знаком и 'twosided' для сигналов с комплексным знаком.

  • 'onesided' — Возвращает одностороннюю оценку передаточной функции между двумя входными сигналами с действительным знаком, x и y. Если nfft даже, txy имеет nfft/2 + 1 строка и вычисляется на интервале [0, π] рад/выборка. Если nfft нечетен, txy имеет (nfft + 1)/2 строки, и интервал [0, π), рад/выборка. Если вы задаете fs, соответствующие интервалы являются [0, fs/2] циклами/единицей времени для даже nfft и [0, fs/2), циклы/единица времени для нечетного nfft.

  • 'twosided' — Возвращает двухстороннюю оценку передаточной функции между двумя входными сигналами с комплексным знаком или с действительным знаком, x и y. В этом случае txy имеет строки nfft и вычисляется на интервале [0,2π), рад/выборка. Если вы задаете fs, интервал [0, fs), циклы/единица времени.

  • 'centered' — Возвращает двухстороннюю оценку в центре передаточной функции между двумя входными сигналами с комплексным знаком или с действительным знаком, x и y. В этом случае txy имеет строки nfft и вычисляется на интервале (–π, π] рад/выборка для даже nfft и (–π, π) рад/выборка для нечетного nfft. Если вы задаете fs, соответствующие интервалы (–fs/2, fs/2] циклы/единица времени для даже nfft и (–fs/2, fs/2) циклы/единица времени для нечетного nfft.

Средство оценки передаточной функции, заданное как 'H1' или 'H2'.

  • Используйте 'H1', когда шум будет некоррелированым с входными сигналами.

  • Используйте 'H2', когда шум будет некоррелированым с выходными сигналами. В этом случае количество входных сигналов должно равняться количеству выходных сигналов.

Смотрите Передаточную функцию для получения дополнительной информации.

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

свернуть все

Оценка передаточной функции, возвращенная как вектор, матрица или 3D массив.

Нормированные частоты, возвращенные как вектор-столбец с действительным знаком.

Циклические частоты, возвращенные как вектор-столбец с действительным знаком.

Больше о

свернуть все

Передаточная функция

Отношение между входом x и выводом y моделируется линейной, независимой от времени передаточной функцией txy. В частотном диапазоне, Y (f) = H (f) X (f).

  • Для single-input/single-output системы H 1 оценкой передаточной функции дают

    H1(f)=Pyx(f)Pxx(f),

    где Pyx является перекрестной степенью спектральная плотность x и y, и Pxx является степенью спектральная плотность x. Эта оценка принимает, что шум не коррелируется с системным входом.

    Для мультивхода/мультивывода (MIMO) системы H 1 средство оценки становится

    H1(f)=PYX(f)PXX1(f)=[Py1x1(f)Py1x2(f)Py1xm(f)Py2x1(f)Py2x2(f)Py2xm(f)Pynx1(f)Pynx2(f)Pynxm(f)][Px1x1(f)Px1x2(f)Px1xm(f)Px2x1(f)Px2x2(f)Px2xm(f)Pxmx1(f)Pxmx2(f)Pxmxm(f)]1

    для входных параметров m и n выходные параметры, где:

    • Pyixk является перекрестной степенью спектральная плотность k th вход и i th вывод.

    • Pxixk является перекрестной степенью спектральная плотность k th и i th входные параметры.

    Для двух входных параметров и двух выходных параметров, средство оценки является матрицей

    H1(f)=[Py1x1(f)Px2x2(f)Py1x2(f)Px2x1(f)Py1x2(f)Px1x1(f)Py1x1(f)Px1x2(f)Py2x1(f)Px2x2(f)Py2x2(f)Px2x1(f)Py2x2(f)Px1x1(f)Py2x1(f)Px1x2(f)]Px1x1(f)Px2x2(f)Px1x2(f)Px2x1(f).

  • Для single-input/single-output системы H 2 оценками передаточной функции дают

    H2(f)=Pyy(f)Pxy(f),

    где Pyy является степенью спектральная плотность y и Pxy = , P *yx является сопряженным комплексным числом перекрестной степени спектральная плотность x и y. Эта оценка принимает, что шум не коррелируется с системой вывод.

    Для систем MIMO H 2 средства оценки четко определены только для равных количеств вводов и выводов: n = m. Средство оценки становится

    H2(f)=PYY(f)PXY1(f)=[Py1y1(f)Py1y2(f)Py1yn(f)Py2y1(f)Py2y2(f)Py2yn(f)Pyny1(f)Pyny2(f)Pynyn(f)][Px1y1(f)Px1y2(f)Px1yn(f)Px2y1(f)Px2y2(f)Px2yn(f)Pxny1(f)Pxny2(f)Pxnyn(f)]1,

    где:

    • Pyiyk является перекрестной степенью спектральная плотность k th и i th выходные параметры.

    • Pxiyk является сопряженным комплексным числом перекрестной степени спектральная плотность i th вход и k th вывод.

Алгоритмы

tfestimate использует усредненный метод периодограммы валлийцев. Смотрите pwelch для деталей.

Ссылки

[1] Vold, Håvard, Джон Кроули и Г. Томас Роклин. “Новые Способы Оценить Функции Частотной характеристики”. Звук и Вибрация. Издание 18, ноябрь 1984, стр 34–38.

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

| | |

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

Для просмотра документации необходимо авторизоваться на сайте