Куйте возбуждение

Загрузите файл, содержащий three-input/three-output данные о возбуждении молотка, произведенные на уровне 4 кГц. Используйте первое ${10}^4$ выборки для оценки и выборки $2\times {10}^4$ к $5\times {10}^4$ для качественной валидации модели. Задайте шаг расчета как инверсию частоты дискретизации. Храните данные как @iddata объекты.

load modaldata XhammerMISO1 YhammerMISO1 fs

rest = 1:1e4;
rval = 2e4:5e4;
Ts = 1/fs;

Estimation = iddata(YhammerMISO1(rest,:),XhammerMISO1(rest,:),Ts);
Validation = iddata(YhammerMISO1(rval,:),XhammerMISO1(rval,:),Ts,'Tstart',rval(1)*Ts);

Отобразите на графике данные об оценке и данные о валидации.

plot(Estimation,Validation)
legend(gca,'show')

Используйте ssest функционируйте, чтобы оценить модель в пространстве состояний 7-го порядка системы, которая минимизирует ошибку симуляции между измеренными выходными параметрами и выходными параметрами модели. Укажите, что модель в пространстве состояний имеет сквозное соединение.

Orders = 7;
opt = ssestOptions('Focus','simulation');

sys = ssest(Estimation,Orders,'Feedthrough',true,'Ts',Ts,opt);

(Чтобы найти порядок модели, который дает лучший компромисс между точностью и сложностью, устанавливает Orders к 1:15 в предыдущем коде. ssest выводит логарифмический график сингулярных значений, который позволяет вам задать порядок в интерактивном режиме. Функция также рекомендует порядок модели 7.)

Подтвердите качество модели на наборе данных валидации. Постройте меру качества подгонки нормированной среднеквадратичной ошибки (NRMSE). Модель описывает точно выходные сигналы данных о валидации.

compare(Validation,sys)

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

[frf,f] = modalfrf(sys);
modalfrf(sys)

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

Modes = 3;
[fn,dr,ms] = modalfit(sys,f,Modes)

fn = 3×1
103 ×

    0.3727
    0.8525
    1.3706

dr = 3×1

    0.0008
    0.0018
    0.0029

ms = 3×3 complex

   0.0036 - 0.0019i   0.0039 - 0.0005i   0.0021 + 0.0006i
   0.0043 - 0.0023i   0.0010 - 0.0001i  -0.0033 - 0.0010i
   0.0040 - 0.0021i  -0.0031 + 0.0004i   0.0011 + 0.0003i

Вычислите и отобразите восстановленную функцию частотной характеристики. Опишите величины в децибелах.

[~,~,~,ofrf] = modalfit(sys,f,Modes);

clf
for ij = 1:3
    for ji = 1:3
        subplot(3,3,3*(ij-1)+ji)
        plot(f/1000,20*log10(abs(ofrf(:,ji,ij))))
        axis tight
        title(sprintf('In%d -> Out%d',ij,ji))
        if ij==3
            xlabel('Frequency (kHz)')
        end
    end
end

Управляемый нестабильный процесс

Загрузите файл, содержащий высокое модальное измерение частотной характеристики плотности. Данные соответствуют нестабильному процессу, обеспеченному в равновесии с помощью управления с обратной связью. Храните данные как idfrd объект для идентификации. Постройте диаграмму Боде.

load HighModalDensData FRF f

G = idfrd(permute(FRF,[2 3 1]),f,0,'FrequencyUnit','Hz');
figure
bodemag(G)
xlim([0.01,2e3])

Идентифицируйте передаточную функцию с 32 полюсами и 32 нулями.

sys = tfest(G,32,32);

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

bodemag(G,sys)
xlim([0.01,2e3])
legend(gca,'show')

Извлеките собственные частоты и коэффициенты затухания первых 10 наименее ослабленных колебательных режимов. Сохраните результаты в таблице.

[fn,dr] = modalfit(sys,[],10);
T = table((1:10)',fn,dr,'VariableNames',{'Mode','Frequency','Damping'})

T=10×3 table
    Mode    Frequency     Damping 
    ____    _________    _________

      1      82.764       0.011304
      2      85.013       0.015632
      3      124.04       0.025252
      4      142.04       0.017687
      5      251.46      0.0062182
      6      332.79      0.0058266
      7      401.21      0.0043645
      8      625.14      0.0039247
      9      770.49       0.002795
     10      943.64      0.0019943