modalfrf

Функция частотной характеристики для модального анализа

Описание

пример

frf = modalfrf(x,y,fs,window) оценивает матрицу функций частотной характеристики, frf, от сигналов возбуждения, x, и сигналы ответа, y, все произведенные на уровне fs. Выход, frf, H 1 оценка, вычисленная с помощью метода валлийцев с window к окну сигналы. x и y должен иметь одинаковое число строк. Если x или y матрица, каждый столбец представляет сигнал.

Отклик системы, y, принят, чтобы содержать ускоряющие измерения. Чтобы вычислить функцию частотной характеристики, начинающую со смещения или скоростных измерений, используйте 'Sensor' аргумент. modalfrf всегда выводит функцию частотной характеристики в динамической гибкости (receptance) формат независимо от типа датчика.

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

пример

frf = modalfrf(___,Name,Value) задает опции с помощью аргументов name-value, с помощью любой комбинации входных параметров от предыдущих синтаксисов. Опции включают средство оценки, настройку измерения и тип датчика, измеряющего отклик системы.

пример

[frf,f,coh] = modalfrf(___) также возвращает вектор частоты, соответствующий каждой функции частотной характеристики, а также нескольким матрица когерентности.

[frf,f] = modalfrf(sys) вычисляет функцию частотной характеристики идентифицированной модели sys. Используйте команды оценки как ssest (System Identification Toolbox), n4sid (System Identification Toolbox), или tfest (System Identification Toolbox), чтобы создать sys от сигналов ввода и вывода временного интервала. Этот синтаксис позволяет использование только 'Sensor' аргумент значения имени. У вас должна быть лицензия System Identification Toolbox™, чтобы использовать этот синтаксис.

frf = modalfrf(sys,f) задает частоты, на которых можно вычислить frf. Этот синтаксис позволяет использование только 'Sensor' аргумент значения имени. У вас должна быть лицензия System Identification Toolbox, чтобы использовать этот синтаксис.

пример

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

Примеры

свернуть все

Визуализируйте функцию частотной характеристики возбуждения молотка single-input/single-output.

Загрузите файл данных, который содержит:

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

  • Yhammer Ответ системы к входу. Yhammer измеряется как смещение.

Сигналы производятся на уровне 4 кГц. Постройте возбуждение и выходные сигналы.

load modaldata

subplot(2,1,1)
plot(thammer,Xhammer(:))
ylabel('Force (N)')
subplot(2,1,2)
plot(thammer,Yhammer(:))
ylabel('Displacement (m)')
xlabel('Time (s)')

Figure contains 2 axes objects. Axes object 1 contains an object of type line. Axes object 2 contains an object of type line.

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

clf
winlen = size(Xhammer,1);
modalfrf(Xhammer(:),Yhammer(:),fs,winlen,'Sensor','dis')

Figure contains 2 axes objects. Axes object 1 with title FRF11 contains an object of type line. Axes object 2 contains an object of type line.

Вычислите функцию частотной характеристики для two-input/two-output системы, взволнованной случайным шумом.

Загрузите файл данных, который содержит Xrand, входной сигнал возбуждения и Yrand, отклик системы. Вычислите функцию частотной характеристики с помощью окна Hann с 5000 выборками и 50%-го перекрытия между смежными сегментами данных. Укажите, что выходные измерения являются смещениями.

load modaldata
winlen = 5000;

frf = modalfrf(Xrand,Yrand,fs,hann(winlen),0.5*winlen,'Sensor','dis');

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

modalfrf(Xrand,Yrand,fs,hann(winlen),0.5*winlen,'Sensor','dis')

Figure contains 8 axes objects. Axes object 1 with title FRF11 contains an object of type line. Axes object 2 contains an object of type line. Axes object 3 with title FRF12 contains an object of type line. Axes object 4 contains an object of type line. Axes object 5 with title FRF21 contains an object of type line. Axes object 6 contains an object of type line. Axes object 7 with title FRF22 contains an object of type line. Axes object 8 contains an object of type line.

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

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

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

Fs = 1;
dt = 1/Fs;
N = 3000;
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=r измеренный выход. Матрицы пространства состояний

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

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(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

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

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

Сгенерируйте в 30,000 раз выборки, эквивалентные 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=[r1r2]T выходной вектор. Матрицы пространства состояний

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

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

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

Набор 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')

Figure contains 8 axes objects. Axes object 1 with title FRF11 contains an object of type line. Axes object 2 contains an object of type line. Axes object 3 with title FRF12 contains an object of type line. Axes object 4 contains an object of type line. Axes object 5 with title FRF21 contains an object of type line. Axes object 6 contains an object of type line. Axes object 7 with title FRF22 contains an object of type line. Axes object 8 contains an object of type line.

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

[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

Вычислите функцию частотной характеристики two-input/six-output набора данных, соответствующего стальной конструкции.

Загрузите структуру, содержащую входные возбуждения и выходные измерения акселерометра. Система производится на уровне 1 024 Гц в течение приблизительно 3,9 секунд.

load modaldata SteelFrame
X = SteelFrame.Input;
Y = SteelFrame.Output;
fs = SteelFrame.Fs;

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

[frf,f] = modalfrf(X,Y,fs,1000,'Estimator','subspace','Order',36);

Визуализируйте диаграмму стабилизации для системы. Идентифицируйте до 15 физических режимов.

modalsd(frf,f,fs,'MaxModes',15)

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.

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

свернуть все

Возбуждение сигнализирует в виде вектора или матрицы.

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

Ответ сигнализирует в виде вектора или матрицы.

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

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

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

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

  • Если window целое число, затем modalfrf делит x и y в сегменты длины window и окна каждый сегмент с прямоугольным окном той длины.

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

  • Если 'Estimator' задан как 'subspace'то modalfrf игнорирует форму window и использует его длину, чтобы определить количество точек частоты в возвращенной функции частотной характеристики.

Если длина 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.

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

Идентифицированная система в виде модели идентифицированными параметрами. Используйте команды оценки как ssest (System Identification Toolbox), n4sid (System Identification Toolbox), или tfest (System Identification Toolbox), чтобы создать sys от сигналов ввода и вывода временного интервала. Смотрите Модальный Анализ Идентифицированных Моделей для примера. Синтаксисы то использование 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 Гц.

Частоты в виде вектора описываются в Гц.

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

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

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

Пример: 'Sensor','vel','Est','H1' указывает, что входной сигнал состоит из скоростных измерений и что предпочтительное средство оценки является H1.

Средство оценки в виде 'H1', 'H2', 'Hv', или 'subspace'. Смотрите Передаточную функцию для получения дополнительной информации о H 1 и H 2 средства оценки.

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

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

  • Используйте 'Hv' минимизировать несоответствие между смоделированными и предполагаемыми данными об ответе путем минимизации трассировки ошибочной матрицы. H v является геометрическим средним значением H 1 и H 2: H v = (H 1H2)1/2

    Измерение должно быть single-input/single-output (SISO).

  • Используйте 'subspace' вычислить функцию частотной характеристики с помощью модели в пространстве состояний. В этом случае, noverlap аргумент проигнорирован. Этот метод обычно требует меньшего количества данных, чем непараметрические подходы. Смотрите n4sid (System Identification Toolbox) для получения дополнительной информации.

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

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

Настройка измерения для равных количеств возбуждения и ответа образовывает канал в виде 'fixed', 'rovinginput', или 'rovingoutput'.

  • Используйте 'fixed' когда существуют источники возбуждения и датчики в фиксированных местоположениях системы. Каждое возбуждение способствует каждому ответу.

  • Используйте 'rovinginput' когда измерения следуют из мобильного возбуждения (или roving hammer) тест. Один датчик сохранен в фиксированном местоположении системы. Один источник возбуждения помещается в несколько местоположений и производит один ответ датчика на местоположение. Функциональный выход frf(:,:,i) = modalfrf(x(:,i),y(:,i)).

  • Используйте 'rovingoutput' когда измерения следуют из теста roving sensor. Один источник возбуждения сохранен в фиксированном местоположении системы. Один датчик помещается в несколько местоположений и отвечает на одно возбуждение на местоположение. Функциональный выход frf(:,i) = modalfrf(x(:,i),y(:,i)).

Порядок модели в пространстве состояний в виде целочисленного или вектора-строки из целых чисел. Если вы задаете вектор из целых чисел, то функция выбирает оптимальную стоимость заказа из заданной области. Этот аргумент доступен только если 'Estimator' задан как 'subspace'.

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

Тип датчика в виде 'acc', 'vel', или 'dis'.

  • 'acc' — Указывает, что сигнал ответа системы пропорционален ускорению.

  • 'vel' — Указывает, что сигнал ответа системы пропорционален скорости.

  • 'dis' — Указывает, что сигнал ответа системы пропорционален смещению.

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

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

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

H(z)=NSensor(z)1-2z-1cosΔt+z-2,

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

  • Смещение: Ndis(z)=(z-1+z-2)(1-cosΔt)

  • Скорость: Nvel(z)=(z-1-z-2)sinΔt

  • Ускорение: Nacc(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.

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

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

свернуть все

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

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

Частоты, возвращенные как вектор.

Несколько матрица когерентности, возвращенная как матрица. coh имеет один столбец для каждого сигнала ответа.

Ссылки

[1] "Динамическая Жесткость, Податливость, Мобильность и другие..." Siemens, в последний раз измененный 2019, https://community.sw.siemens.com/s/article/dynamic-stiffness-compliance-mobility-and-more.

[2] Брандт, Андерс. Шум и анализ вибрации: анализ сигнала и экспериментальные процедуры. Чичестер, Великобритания: John Wiley & Sons, 2011.

[3] Ирвин, Том. "Введение в Функции Частотной характеристики", Vibrationdata, 2000, https://vibrationdata.com/tutorials2/frf.pdf.

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

Введенный в R2017a