101 lines
2.2 KiB
Mathematica
101 lines
2.2 KiB
Mathematica
|
%% Clear Workspace and Close figures
|
||
|
clear; close all; clc;
|
||
|
|
||
|
%% Intialize Laplace variable
|
||
|
s = zpk('s');
|
||
|
|
||
|
addpath('./mat/');
|
||
|
|
||
|
% Load Data
|
||
|
% The data from the "noise test" and the identification test are loaded.
|
||
|
|
||
|
ht = load('huddle_test.mat', 't', 'y');
|
||
|
load('apa95ml_5kg_Amp_E505.mat', 't', 'um', 'y');
|
||
|
|
||
|
|
||
|
|
||
|
% Any offset value is removed:
|
||
|
|
||
|
um = detrend(um, 0); % Generated DAC Voltage [V]
|
||
|
y = detrend(y , 0); % Mass displacement [m]
|
||
|
|
||
|
ht.y = detrend(ht.y, 0);
|
||
|
|
||
|
|
||
|
|
||
|
% Now we add a factor 10 to take into account the gain of the voltage amplifier.
|
||
|
|
||
|
Va = 10*um;
|
||
|
|
||
|
% Comparison of the PSD with Huddle Test
|
||
|
|
||
|
Ts = t(end)/(length(t)-1);
|
||
|
Fs = 1/Ts;
|
||
|
|
||
|
win = hanning(ceil(1*Fs));
|
||
|
|
||
|
[pxx, f] = pwelch(y, win, [], [], Fs);
|
||
|
[pht, ~] = pwelch(ht.y, win, [], [], Fs);
|
||
|
|
||
|
figure;
|
||
|
hold on;
|
||
|
plot(f, sqrt(pxx), 'DisplayName', '5kg');
|
||
|
plot(f, sqrt(pht), 'DisplayName', 'Huddle Test');
|
||
|
hold off;
|
||
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
||
|
xlabel('Frequency [Hz]'); ylabel('ASD [$m/\sqrt{Hz}$]');
|
||
|
legend('location', 'northeast');
|
||
|
xlim([1, Fs/2]);
|
||
|
|
||
|
% Compute TF estimate and Coherence
|
||
|
|
||
|
[tf_est, f] = tfestimate(Va, -y, win, [], [], 1/Ts);
|
||
|
[co_est, ~] = mscohere( Va, -y, win, [], [], 1/Ts);
|
||
|
|
||
|
figure;
|
||
|
|
||
|
hold on;
|
||
|
plot(f, co_est, 'k-')
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
||
|
ylabel('Coherence'); xlabel('Frequency [Hz]');
|
||
|
hold off;
|
||
|
xlim([10, 5e3]);
|
||
|
|
||
|
|
||
|
|
||
|
% #+name: fig:apa95ml_5kg_PI_coh
|
||
|
% #+caption: Coherence
|
||
|
% #+RESULTS:
|
||
|
% [[file:figs/apa95ml_5kg_PI_coh.png]]
|
||
|
|
||
|
% Comparison with the FEM model
|
||
|
|
||
|
load('mat/fem_simscape_models.mat', 'Ghm');
|
||
|
|
||
|
freqs = logspace(0, 4, 1000);
|
||
|
|
||
|
figure;
|
||
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
||
|
|
||
|
ax1 = nexttile([2,1]);
|
||
|
hold on;
|
||
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
||
|
plot(freqs, abs(squeeze(freqresp(Ghm, freqs, 'Hz'))), 'DisplayName', 'FEM')
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
||
|
ylabel('Amplitude [m/V]'); set(gca, 'XTickLabel', []);
|
||
|
legend('location', 'northeast')
|
||
|
hold off;
|
||
|
|
||
|
ax2 = nexttile;
|
||
|
hold on;
|
||
|
plot(f, 180/pi*angle(tf_est))
|
||
|
plot(freqs, 180/pi*angle(squeeze(freqresp(Ghm, freqs, 'Hz'))))
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
||
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
||
|
hold off;
|
||
|
ylim([-180, 180]);
|
||
|
yticks(-180:90:180);
|
||
|
|
||
|
linkaxes([ax1,ax2], 'x');
|
||
|
xlim([10, 5e3]);
|