54 lines
1.6 KiB
Mathematica
54 lines
1.6 KiB
Mathematica
|
%% Load Data
|
||
|
iff_g10 = load('../mat/apa95ml_iff_g10_res.mat', 'u', 't', 'y', 'v');
|
||
|
iff_g100 = load('../mat/apa95ml_iff_g100_res.mat', 'u', 't', 'y', 'v');
|
||
|
iff_of = load('../mat/apa95ml_iff_off_res.mat', 'u', 't', 'y', 'v');
|
||
|
|
||
|
%% Compute TF
|
||
|
Ts = 1e-4;
|
||
|
win = hann(ceil(10/Ts));
|
||
|
|
||
|
[tf_iff_g10, f] = tfestimate(iff_g10.u, iff_g10.y, win, [], [], 1/Ts);
|
||
|
[co_iff_g10, ~] = mscohere(iff_g10.u, iff_g10.y, win, [], [], 1/Ts);
|
||
|
|
||
|
[tf_iff_g100, f] = tfestimate(iff_g100.u, iff_g100.y, win, [], [], 1/Ts);
|
||
|
[co_iff_g100, ~] = mscohere(iff_g100.u, iff_g100.y, win, [], [], 1/Ts);
|
||
|
|
||
|
[tf_iff_of, ~] = tfestimate(iff_of.u, iff_of.y, win, [], [], 1/Ts);
|
||
|
[co_iff_of, ~] = mscohere(iff_of.u, iff_of.y, win, [], [], 1/Ts);
|
||
|
|
||
|
%% Coherence
|
||
|
figure;
|
||
|
|
||
|
hold on;
|
||
|
plot(f, co_iff_of, '-', 'DisplayName', 'g=0')
|
||
|
plot(f, co_iff_g10, '-', 'DisplayName', 'g=10')
|
||
|
plot(f, co_iff_g100, '-', 'DisplayName', 'g=100')
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
||
|
ylabel('Coherence'); xlabel('Frequency [Hz]');
|
||
|
hold off;
|
||
|
legend();
|
||
|
xlim([60, 600])
|
||
|
|
||
|
%% Compare Dynamics IFF ON/OFF
|
||
|
figure;
|
||
|
ax1 = subplot(2, 1, 1);
|
||
|
hold on;
|
||
|
plot(f, abs(tf_iff_of), '-', 'DisplayName', 'g=0')
|
||
|
plot(f, abs(tf_iff_g10), '-', 'DisplayName', 'g=10')
|
||
|
plot(f, abs(tf_iff_g100), '-', 'DisplayName', 'g=100')
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
||
|
ylabel('Amplitude'); xlabel('Frequency [Hz]');
|
||
|
hold off;
|
||
|
legend();
|
||
|
|
||
|
ax2 = subplot(2, 1, 2);
|
||
|
hold on;
|
||
|
plot(f, 180/pi*angle(-tf_iff_of), '-')
|
||
|
plot(f, 180/pi*angle(-tf_iff_g10), '-')
|
||
|
plot(f, 180/pi*angle(-tf_iff_g100), '-')
|
||
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
||
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
||
|
hold off;
|
||
|
|
||
|
linkaxes([ax1,ax2], 'x');
|
||
|
xlim([60, 600]);
|