Документация

fn = modalfit(frf,f,fs,mnum) оценивает собственные частоты mnum режимы системы с измеренной функцией частотной характеристики frf заданный на частотах f и для частоты дискретизации fsИспользование modalfrf сгенерировать матрицу функции частотной характеристики от результатов измерений. frf принят, чтобы быть в динамической гибкости (receptance) формат.

fn = modalfit(frf,f,fs,mnum,Name,Value) задает дополнительные опции с помощью аргументов name-value.

[fn,dr,ms] = modalfit(___) также возвращает коэффициенты затухания и векторы формы режима, соответствующие каждой собственной частоте в fn, использование любой комбинации входных параметров от предыдущих синтаксисов.

[fn,dr,ms,ofrf] = modalfit(___) также возвращает восстановленный массив функции частотной характеристики на основе предполагаемых модальных параметров.

[___] = modalfit(sys,f,mnum,Name,Value) оценивает модальные параметры идентифицированной модели sys. Используйте команды оценки как ssest (System Identification Toolbox) или tfest (System Identification Toolbox), чтобы создать sys запуск с измеренной функции частотной характеристики или с сигналов ввода и вывода временного интервала. Этот синтаксис позволяет использование 'DriveIndex', 'FreqRange', и 'PhysFreq' аргументы name-value. Обычно требуется меньше данных, чем синтаксисы, которые используют непараметрические методы. У вас должна быть лицензия System Identification Toolbox™, чтобы использовать этот синтаксис.

Примеры

Функция частотной характеристики системы SISO

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

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

Сгенерируйте в 3000 раз выборки. Задайте интервал выборки $Δ t = 1 / F_{s}$ .

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

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

$\begin{array}{c} x (k + 1) = A x (k) + B u (k), \\ y (k) = C x (k) + D u (k), \end{array}$

где $x = {[\begin{matrix} r & v \end{matrix}]}^{T}$ вектор состояния, $r$ и $v$ соответственно смещение и скорость массы, $u$ движущая сила, и $y = r$ измеренный выход. Матрицы пространства состояний

$A = \exp (A_{c} Δ t), B = A_{c}^{- 1} (A - I) B_{c}, C = [\begin{array}{c} 1 & 0 \end{array}], D = 0,$

$I$ $2 \times 2$ идентичность и матрицы пространства состояний непрерывного времени

$A_{c} = [\begin{array}{c} 0 & 1 \\ - 1 & - b \end{array}], B_{c} = [\begin{array}{c} 0 \\ 1 \end{array}] .$

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

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

C = [1 0];
D = 0;

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

rng default
u = randn(1,N)/2;
u(2001:end) = 0;

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)

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

Оцените модальную функцию частотной характеристики системы. Используйте окно Hann, вдвое менее длинное, чем измеренные сигналы. Укажите, что выход является смещением массы.

wind = hann(N/2);

[frf,f] = modalfrf(u',y',Fs,wind,'Sensor','dis');

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

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

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

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

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

Оцените собственную частоту и коэффициент затухания для режима вибрации.

[fn,dr] = modalfit(frf,f,Fs,1,'FitMethod','PP')

fn = 0.1593

dr = 0.0043

Сравните собственную частоту с $1 / 2 π$ , который является теоретическим значением для незатухающей системы.

theo = 1/(2*pi)

theo = 0.1592

Модальные параметры Используя метод рациональной функции наименьших квадратов

Вычислите модальные параметры модуля Космической станции, начинающего с его массива функции частотной характеристики (FRF).

Загрузите структуру, содержащую three-input/three-output FRF массив. Система производится на уровне 320 Гц.

load modaldata SpaceStationFRF

frf = SpaceStationFRF.FRF; 
f = SpaceStationFRF.f; 
fs = SpaceStationFRF.Fs;

Извлеките модальные параметры самых низких 24 режимов с помощью метода рациональной функции наименьших квадратов.

[fn,dr,ms,ofrf] = modalfit(frf,f,fs,24,'FitMethod','lsrf');

Сравните восстановленный массив FRF с измеренным.

for ij = 1:3
 for ji = 1:3
    subplot(3,3,3*(ij-1)+ji)
    loglog(f,abs(frf(:,ji,ij)))
    hold on
    loglog(f,abs(ofrf(:,ji,ij)))
    hold off
    axis tight
    title(sprintf('In%d -> Out%d',ij,ji))
    if ij==3
        xlabel('Frequency (Hz)')
    end
 end
end

Модальные параметры генератора 2D тела

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

Идеальная одномерная колеблющаяся система состоит из двух масс, $m_{1}$ и $m_{2}$ , ограниченный между двумя стенами. Модули таковы что $m_{1} = 1$ и $m_{2} = μ$ . Каждая масса присоединена к самой близкой стене к пружине с эластичной константой $k$ . Идентичная пружина соединяет эти две массы. Три демпфера препятствуют движению масс путем проявления на них, обеспечивает пропорциональный, чтобы ускориться, с затуханием постоянного $b$ . Выборка датчиков $r_{1}$ и $r_{2}$ , смещения масс, в $F_{s} = 50$ Гц.

Сгенерируйте в 30,000 раз выборки, эквивалентные 600 секундам. Задайте интервал выборки $Δ t = 1 / F_{s}$ .

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

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

$\begin{array}{c} x (k + 1) = A x (k) + B u (k), \\ y (k) = C x (k) + D u (k), \end{array}$

где $x = {[\begin{array}{c} r_{1} & v_{1} & r_{2} & v_{2} \end{array}]}^{T}$ вектор состояния, $r_{i}$ и $v_{i}$ соответственно местоположение и скорость $i$ масса th, $u = {[\begin{array}{c} u_{1} & u_{2} \end{array}]}^{T}$ вектор из входных движущих сил, и $y = {[\begin{array}{c} r_{1} & r_{2} \end{array}]}^{T}$ выходной вектор. Матрицы пространства состояний

$A = \exp (A_{c} Δ t), B = A_{c}^{- 1} (A - I) B_{c}, C = [\begin{array}{c} 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 \end{array}], D = [\begin{array}{c} 0 & 0 \\ 0 & 0 \end{array}],$

$I$ $4 \times 4$ идентичность и матрицы пространства состояний непрерывного времени

$A_{c} = [\begin{array}{c} 0 & 1 & 0 & 0 \\ - 2 k & - 2 b & k & b \\ 0 & 0 & 0 & 1 \\ k / μ & b / μ & - 2 k / μ & - 2 b / μ \end{array}], B_{c} = [\begin{array}{c} 0 & 0 \\ 1 & 0 \\ 0 & 0 \\ 0 & 1 / μ \end{array}] .$

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

k = 400;
b = 0.1;
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 = [1 0 0 0;0 0 1 0];
D = zeros(2);

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

rng default
u = randn(2,N);

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

Используйте входные и выходные данные, чтобы оценить передаточную функцию системы в зависимости от частоты. Используйте окно Hann с 15000 выборками с 9 000 выборок перекрытия между смежными сегментами. Укажите, что измеренные выходные параметры являются смещениями.

wind = hann(15000);
nove = 9000;
[FRF,f] = modalfrf(u',y',Fs,wind,nove,'Sensor','dis');

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

nfs = 2048;
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(f,20*log10(abs(FRF(:,jk,kj))))
        hold on
        plot(fz*Fs,20*log10(abs(frf(jk,:,kj))))
        hold off
        axis([0 Fs/2 -100 0])
        title(sprintf('Input %d, Output %d',jk,kj))
    end
end

Figure contains 4 axes objects. Axes object 1 with title Input 1, Output 1 contains 2 objects of type line. Axes object 2 with title Input 1, Output 2 contains 2 objects of type line. Axes object 3 with title Input 2, Output 1 contains 2 objects of type line. Axes object 4 with title Input 2, Output 2 contains 2 objects of type line.

Постройте оценки при помощи синтаксиса modalfrf без выходных аргументов.

figure
modalfrf(u',y',Fs,wind,nove,'Sensor','dis')

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

[fn,dr,ms] = modalfit(FRF,f,Fs,2,'FitMethod','pp');
fn

fn = 
fn(:,:,1) =

    3.8466    3.8466
    3.8495    3.8495


fn(:,:,2) =

    3.8492    3.8490
    3.8552   14.4684

Сравните собственные частоты с теоретическими предсказаниями для незатухающей системы.

undamped = sqrt(eig([2*k -k;-k/m 2*k/m]))/2/pi

undamped = 2×1

    3.8470
   14.4259

Модальные параметры системы MIMO

Вычислите собственные частоты, коэффициенты затухания и формы режима для two-input/three-output системы, взволнованной несколькими пакетами случайного шума. Каждый пакет служит в течение 1 секунды, и существует 2 секунды между концом каждого пакета и запуском следующего. Данные производятся на уровне 4 кГц.

Загрузите файл данных. Постройте входные сигналы и выходные сигналы.

load modaldata

subplot(2,1,1)
plot(Xburst)
title('Input Signals')
subplot(2,1,2)
plot(Yburst)
title('Output Signals')

Figure contains 2 axes objects. Axes object 1 with title Input Signals contains 2 objects of type line. Axes object 2 with title Output Signals contains 3 objects of type line.

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

burstLen = 12000;
[frf,f] = modalfrf(Xburst,Yburst,fs,burstLen);

Визуализируйте диаграмму стабилизации и возвратите устойчивые собственные частоты. Задайте максимальный порядок модели 30 режимов.

figure
modalsd(frf,f,fs,'MaxModes',30);

Figure contains an axes object. The axes object with title Stabilization Diagram contains 4 objects of type line. These objects represent Stable in frequency, Stable in frequency and damping, Not stable in frequency, Averaged response function.

Увеличьте масштаб графика. Усредненная функция отклика имеет максимумы на уровне 373 Гц, 852 Гц и 1 371 Гц, который соответствует физическим частотам системы. Сохраните максимумы в переменную.

phfr = [373 852 1371];

Вычислите модальные параметры с помощью алгоритма комплексной экпоненты наименьших квадратов (LSCE). Задайте порядок модели 6 режимов и задайте физические частоты для этих 3 режимов, определенных из диаграммы стабилизации. Функция генерирует один набор собственных частот и коэффициентов затухания для каждой входной ссылки.

[fn,dr,ms,ofrf] = modalfit(frf,f,fs,6,'PhysFreq',phfr);

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

for k = 1:2
    for m = 1:3
        subplot(2,3,m+3*(k-1))
        plot(f/1000,10*log10(abs(frf(:,m,k))))
        hold on
        plot(f/1000,10*log10(abs(ofrf(:,m,k))))
        hold off
        text(1,-50,[['Output ';' Input '] num2str([m k]')])
        ylim([-100 -40])
    end
end
subplot(2,3,2)
title('Frequency-Response Functions')

Figure contains 6 axes objects. Axes object 1 contains 3 objects of type line, text. Axes object 2 with title Frequency-Response Functions contains 3 objects of type line, text. Axes object 3 contains 3 objects of type line, text. Axes object 4 contains 3 objects of type line, text. Axes object 5 contains 3 objects of type line, text. Axes object 6 contains 3 objects of type line, text.

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

`frf` — Функция частотной характеристики
вектор | матрица | трехмерный массив

Функция частотной характеристики в виде вектора, матрицы или трехмерного массива. frf имеет размер p-by-m-by-n, где p является количеством интервалов частоты, m является количеством сигналов ответа, и n является количеством сигналов возбуждения, используемых, чтобы оценить передаточную функцию. frf принят, чтобы быть в динамической гибкости (receptance) формат.

Использование modalfrf сгенерировать матрицу функции частотной характеристики от результатов измерений.

Пример: незатухающий гармонический генератор

Движение простого незатухающего гармонического генератора модульной массы и эластичной константы, произведенной на уровне $f_{s} = 1 / Δ t$ описан передаточной функцией

$H (z) = \frac{N_{S e n s o r} (z)}{1 - 2 z^{- 1} \cos Δ t + z^{- 2}}$ ,

где числитель зависит от измеряемой величины:

Смещение: $N_{d i s} (z) = (z^{- 1} + z^{- 2}) (1 - \cos Δ t)$
Скорость: $N_{v e l} (z) = (z^{- 1} - z^{- 2}) \sin Δ t$
Ускорение: $N_{a c c} (z) = (1 - z^{- 1}) - (z^{- 1} - z^{- 2}) \cos Δ t$

Вычислите функцию частотной характеристики для трех возможных типов датчика отклика системы. Используйте частоту дискретизации 2 Гц и 30 000 выборок белого шума, как введено.

fs = 2;
dt = 1/fs;
N = 30000;

u = randn(N,1);

ydis = filter((1-cos(dt))*[0 1 1],[1 -2*cos(dt) 1],u);
[frfd,fd] = modalfrf(u,ydis,fs,hann(N/2),Sensor="dis");

yvel = filter(sin(dt)*[0 1 -1],[1 -2*cos(dt) 1],u);
[frfv,fv] = modalfrf(u,yvel,fs,hann(N/2),Sensor="vel");

yacc = filter([1 -(1+cos(dt)) cos(dt)],[1 -2*cos(dt) 1],u);
[frfa,fa] = modalfrf(u,yacc,fs,hann(N/2),Sensor="acc");

loglog(fd,abs(frfd),fv,abs(frfv),fa,abs(frfa))
grid
legend(["dis" "vel" "acc"],Location="best")

Figure contains an axes object. The axes object contains 3 objects of type line. These objects represent dis, vel, acc.

Во всех случаях сгенерированная функция частотной характеристики находится в формате, соответствующем смещению. Скорость и ускоряющие измерения являются производными первого и второго раза, соответственно, измерений смещения. Функция частотной характеристики эквивалентна в области значений вокруг собственной частоты системы. Далеко от собственной частоты, функция частотной характеристики отличается.

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

`f` — Частоты
вектор

Частоты в виде вектора. Число элементов f должен равняться количеству строк frf.

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

`fs` — Частота дискретизации данных об измерении
положительная скалярная величина

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

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

`mnum` — Количество режимов
положительное целое число

Количество режимов в виде положительного целого числа.

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

`sys` — Идентифицированная система
модель идентифицированными параметрами

Идентифицированная система в виде модели идентифицированными параметрами. Используйте команды оценки как ssest (System Identification Toolbox), n4sid (System Identification Toolbox), или tfest (System Identification Toolbox), чтобы создать sys запуск с измеренной функции частотной характеристики или с сигналов ввода и вывода временного интервала. Смотрите Модальный Анализ Идентифицированных Моделей для примера. У вас должна быть лицензия System Identification Toolbox, чтобы использовать этот входной параметр.

Пример: idss([0.5418 0.8373;-0.8373 0.5334],[0.4852;0.8373],[1 0],0,[0;0],[0;0],1) генерирует идентифицированную модель в пространстве состояний, соответствующую модульной массе, присоединенной к стене к пружине модуля эластичная константа и демпфер с постоянными 0.01. Смещение массы производится на уровне 1 Гц.

Пример: idtf([0 0.4582 0.4566],[1 -1.0752 0.99],1) генерирует идентифицированное соответствие модели передаточной функции модульной массе, присоединенной к стене к пружине модуля эластичная константа и демпфер с постоянными 0.01. Смещение массы производится на уровне 1 Гц.

Аргументы name-value

Задайте дополнительные разделенные запятой пары Name,Value аргументы. Name имя аргумента и Value соответствующее значение. Name должен появиться в кавычках. Вы можете задать несколько аргументов в виде пар имен и значений в любом порядке, например: Name1, Value1, ..., NameN, ValueN.

Пример: 'FitMethod','pp','FreqRange',[0 500] использует выбирающий пик метод, чтобы выполнить подгонку и ограничивает частотный диапазон между 0 и 500 Гц.

`Feedthrough` — Присутствие сквозного соединения в предполагаемой передаточной функции
`false` (значение по умолчанию) | `true`

Присутствие сквозного соединения в предполагаемой передаточной функции в виде логического значения. Этот аргумент доступен только если 'FitMethod' задан как 'lsrf'.

Типы данных: логический

`FitMethod` — Алгоритм подбора
`'lsce'` (значение по умолчанию) | `'lsrf'` | `'pp'`

Алгоритм подбора в виде 'lsce', 'lsrf', или 'pp'.

'lsce' — Метод Комплексной экпоненты наименьших квадратов. Если вы задаете 'lsce', затем fn вектор с mnum элементы, независимые от размера frf.
'lsrf' — Метод оценки рациональной функции наименьших квадратов. Если вы задаете 'lsrf', затем fn вектор с mnum элементы, независимые от размера frf. Метод описан в [3]. Смотрите, что Оценка Передаточной функции Непрерывного времени Использует Данные частотной области Непрерывного времени (System Identification Toolbox) для получения дополнительной информации. Этот алгоритм обычно требует меньшего количества данных, чем непараметрические подходы и является единственным, который работает на неоднородный f.
'pp' — Выбирающий пик метод. Для frf вычисленный из сигналов возбуждения n и сигналов ответа m, fn mnum- m n массивом с одной оценкой fn и одна оценка dr на frf.

`FreqRange` — Частотный диапазон
двухэлементный вектор из увеличения положительных значений

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

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

`PhysFreq` — Собственные частоты для физических режимов
вектор

Собственные частоты для физических режимов, чтобы включать в анализ в виде вектора из значений частоты в области значений, заполненной f. Функция включает в анализ те режимы с собственными частотами, самыми близкими к значениям, заданным в векторе. Если вектор содержит значения частоты m, то fn и dr имейте строки m каждый и ms имеет столбцы m. Если вы не задаете этот аргумент, то функция использует целый частотный диапазон в f.

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

`DriveIndex` — Индексы функции частотной характеристики ведущей точки
[1 1] (значение по умолчанию) | двухэлементный вектор из положительных целых чисел

Индексы функции частотной характеристики ведущей точки в виде двухэлементного вектора из положительных целых чисел. Первый элемент вектора должен быть меньше чем или равен количеству откликов системы. Второй элемент вектора должен быть меньше чем или равен количеству системных возбуждений. Формы режима нормированы к единице, модальной на основе ведущей точки.

Пример: 'DriveIndex',[2 3] указывает, что функцией частотной характеристики ведущей точки является frf(:,2,3).

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

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

`fn` — Собственные частоты
матрица | трехмерный массив

Собственные частоты, возвращенные как матрица или трехмерный массив. Размер fn зависит от выбора алгоритма подбора, заданного с 'FitMethod':

Если вы задаете 'lsce' или 'lsrf', затем fn вектор с mnum элементы, независимые от размера frf. Если система имеет больше, чем mnum колебательные режимы, затем 'lsrf' метод возвращает первый mnum наименее ослабленные режимы отсортированы в порядке увеличения собственной частоты.
Если вы задаете 'pp', затем fn массив размера mnum- m n с одной оценкой fn и одна оценка dr на frf.

`dr` — Коэффициенты затухания
матрица | трехмерный массив

Коэффициенты затухания для собственных частот в fn, возвращенный как матрица или трехмерный массив одного размера с fn.

`ms` — Векторы формы режима
матрица

Векторы формы режима, возвращенные как матрица. ms имеет mnum столбцы, каждый содержащий вектор формы режима из длины q, где q является большим из количества каналов возбуждения и количества каналов ответа.

`ofrf` — Восстановленная функция частотной характеристики
вектор | матрица | трехмерный массив

Восстановленная функция частотной характеристики, возвращенная как вектор, матрица или трехмерный массив с тем же размером как frf.

Алгоритмы