688 lines
22 KiB
Matlab
688 lines
22 KiB
Matlab
%% Clear Workspace and Close figures
|
|
clear; close all; clc;
|
|
|
|
%% Intialize Laplace variable
|
|
s = zpk('s');
|
|
|
|
%% Path for functions, data and scripts
|
|
addpath('./mat/'); % Path for data
|
|
addpath('./src/'); % Path for functions
|
|
|
|
%% Colors for the figures
|
|
colors = colororder;
|
|
|
|
% FRF Identification - Setup
|
|
% Similarly to what was done for the identification of the APA, the identification is performed in three steps:
|
|
% 1. White noise excitation with small amplitude.
|
|
% This is used to determine the main resonance of the system.
|
|
% 2. Sweep sine excitation with the amplitude lowered around the resonance.
|
|
% The sweep sine is from 10Hz to 400Hz.
|
|
% 3. High frequency noise.
|
|
% The noise is band-passed between 300Hz and 2kHz.
|
|
|
|
% Then, the result of the second identification is used between 10Hz and 350Hz and the result of the third identification if used between 350Hz and 2kHz.
|
|
|
|
|
|
%% Sampling frequency/time
|
|
Ts = 1e-4; % Sampling Time [s]
|
|
Nfft = floor(1/Ts);
|
|
win = hanning(Nfft);
|
|
Noverlap = floor(Nfft/2);
|
|
|
|
%% Load Data
|
|
leg_sweep = load('frf_data_leg_1_sweep.mat', 'u', 'Vs', 'de', 'da');
|
|
leg_noise_hf = load('frf_data_leg_1_noise_hf.mat', 'u', 'Vs', 'de', 'da');
|
|
|
|
%% We get the frequency vector that will be the same for all the frequency domain analysis.
|
|
[~, f] = tfestimate(leg_sweep.u, leg_sweep.de, win, Noverlap, Nfft, 1/Ts);
|
|
i_lf = f <= 350; % Indices used for the low frequency
|
|
i_hf = f > 350; % Indices used for the low frequency
|
|
|
|
% FRF Identification - Interferometer
|
|
% In this section, the dynamics from the excitation voltage $u$ to the interferometer $d_a$ is identified.
|
|
% The transfer function from $u$ to the interferometer measured displacement $d_a$ is estimated and shown in Figure ref:fig:strut_1_frf_dvf_plant_tf.
|
|
|
|
%% Compute FRF function from u to da
|
|
[frf_sweep, ~] = tfestimate(leg_sweep.u, leg_sweep.da, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_noise_hf, ~] = tfestimate(leg_noise_hf.u, leg_noise_hf.da, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Combine the FRF
|
|
int_frf = [frf_sweep(i_lf); frf_noise_hf(i_hf)];
|
|
|
|
%% Plot the measured FRF
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(int_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-9, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(int_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% FRF Identification - IFF
|
|
% In this section, the dynamics from $u$ to $V_s$ is identified.
|
|
% Then the FRF are estimated and shown in Figure ref:fig:strut_1_frf_iff_plant_tf
|
|
|
|
%% Compute the FRF
|
|
[frf_sweep, ~] = tfestimate(leg_sweep.u, leg_sweep.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_noise_hf, ~] = tfestimate(leg_noise_hf.u, leg_noise_hf.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Combine the FRF
|
|
iff_frf = [frf_sweep(i_lf); frf_noise_hf(i_hf)];
|
|
|
|
%% Plot the measured FRF
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(iff_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $V_s/u$ [V/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-2, 1e2]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(iff_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% Measurement Data
|
|
% The measurements are loaded.
|
|
|
|
%% Load data
|
|
leg_enc_sweep = load('frf_data_leg_coder_1_noise.mat', 'u', 'Vs', 'de', 'da');
|
|
leg_enc_noise_hf = load('frf_data_leg_coder_1_noise_hf.mat', 'u', 'Vs', 'de', 'da');
|
|
|
|
% FRF Identification - Interferometer
|
|
% In this section, the dynamics from $u$ to $d_a$ is identified.
|
|
|
|
%% Compute FRF function from u to da
|
|
[frf_sweep, ~] = tfestimate(leg_enc_sweep.u, leg_enc_sweep.da, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_noise_hf, ~] = tfestimate(leg_enc_noise_hf.u, leg_enc_noise_hf.da, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Combine the FRF
|
|
int_with_enc_frf = [frf_sweep(i_lf); frf_noise_hf(i_hf)];
|
|
|
|
|
|
|
|
% The obtained FRF is very close to the one that was obtained when no encoder was fixed to the struts as shown in Figure ref:fig:strut_leg_compare_int_frf.
|
|
|
|
%% Plot the FRF from u to da with and without the encoder
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(int_with_enc_frf), '-', 'DisplayName', 'With encoder');
|
|
plot(f, abs(int_frf), '-', 'DisplayName', 'Without encoder');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_a/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-7, 1e-3]);
|
|
legend('location', 'northeast')
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(int_with_enc_frf), '-');
|
|
plot(f, 180/pi*angle(int_frf), '-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% FRF Identification - Encoder
|
|
% In this section, the dynamics from $u$ to $d_e$ (encoder) is identified.
|
|
|
|
% The FRF from $u$ to the encoder measured displacement $d_e$ is computed and shown in Figure ref:fig:strut_1_enc_frf_dvf_plant_tf.
|
|
|
|
%% Compute FRF function from u to da
|
|
[frf_sweep, ~] = tfestimate(leg_enc_sweep.u, leg_enc_sweep.de, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_noise_hf, ~] = tfestimate(leg_enc_noise_hf.u, leg_enc_noise_hf.de, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Combine the FRF
|
|
enc_frf = [frf_sweep(i_lf); frf_noise_hf(i_hf)];
|
|
|
|
%% Plot the FRF from u to de
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(enc_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-7, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(enc_frf), 'k-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
|
|
|
|
% #+name: fig:strut_1_enc_frf_dvf_plant_tf
|
|
% #+caption: Estimated FRF for the DVF plant (transfer function from $u$ to the encoder $d_e$)
|
|
% #+RESULTS:
|
|
% [[file:figs/strut_1_enc_frf_dvf_plant_tf.png]]
|
|
|
|
% The transfer functions from $u$ to $d_e$ (encoder) and to $d_a$ (interferometer) are compared in Figure ref:fig:strut_1_comp_enc_int.
|
|
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(enc_frf), 'DisplayName', 'Encoder');
|
|
plot(f, abs(int_with_enc_frf), 'DisplayName', 'Interferometer');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 2);
|
|
ylim([1e-8, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(enc_frf));
|
|
plot(f, 180/pi*angle(int_with_enc_frf));
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% APA Resonances Frequency
|
|
% As shown in Figure ref:fig:strut_1_spurious_resonances, we can clearly see three spurious resonances at 197Hz, 290Hz and 376Hz.
|
|
|
|
|
|
%% Transfer function from Vs to de with indicated resonances
|
|
figure;
|
|
hold on;
|
|
plot(f, abs(enc_frf), 'k-');
|
|
text(93, 4e-4, {'93Hz'}, 'VerticalAlignment','bottom','HorizontalAlignment','center')
|
|
text(200, 1.3e-4,{'197Hz'},'VerticalAlignment','bottom','HorizontalAlignment','center')
|
|
text(300, 4e-6, {'290Hz'},'VerticalAlignment','bottom','HorizontalAlignment','center')
|
|
text(400, 1.4e-6,{'376Hz'},'VerticalAlignment','bottom','HorizontalAlignment','center')
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
ylim([1e-7, 1e-3]); xlim([10, 2e3]);
|
|
|
|
% FRF Identification - Force Sensor
|
|
% In this section, the dynamics from $u$ to $V_s$ is identified.
|
|
|
|
|
|
%% Compute FRF function from u to da
|
|
[frf_sweep, ~] = tfestimate(leg_enc_sweep.u, leg_enc_sweep.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_noise_hf, ~] = tfestimate(leg_enc_noise_hf.u, leg_enc_noise_hf.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Combine the FRF
|
|
iff_with_enc_frf = [frf_sweep(i_lf); frf_noise_hf(i_hf)];
|
|
|
|
|
|
|
|
% Let's now compare the IFF plants whether the encoders are fixed to the APA or not (Figure ref:fig:strut_1_frf_iff_comp_enc).
|
|
|
|
%% Compare the IFF plant with and without the encoders
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(iff_with_enc_frf), 'DisplayName', 'With Encoder');
|
|
plot(f, abs(iff_frf), 'DisplayName', 'Without Encoder');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $V_s/u$ [V/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend('location', 'northeast', 'FontSize', 8);
|
|
ylim([1e-2, 1e2]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(iff_with_enc_frf));
|
|
plot(f, 180/pi*angle(iff_frf));
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% Non-Minimum phase zero?
|
|
|
|
% In order to determine if the complex conjugate zero of Figure ref:fig:strut_1_enc_frf_iff_plant_tf is minimum phase or non-minimum phase, longer measurements are performed.
|
|
|
|
|
|
long_noise = load('frf_struts_align_3_noise_long.mat', 't', 'u', 'Vs');
|
|
|
|
Ts = 1e-4; % Sampling Time [s]
|
|
Nfft = floor(10/Ts);
|
|
win = hanning(Nfft);
|
|
Noverlap = floor(Nfft/2);
|
|
|
|
%% Transfer function estimation
|
|
[frf_noise, f] = tfestimate(long_noise.u, long_noise.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
%% Bode plot of the FRF from u to de
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(frf_noise), '.-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(frf_noise), '.-');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 0]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([38, 45]);
|
|
|
|
% FRF Identification - Setup
|
|
% The identification of the struts dynamics is performed in two steps:
|
|
% 1. The excitation signal is a white noise with small amplitude.
|
|
% This is used to estimate the low frequency dynamics.
|
|
% 2. Then a high frequency noise band-passed between 300Hz and 2kHz is used to estimate the high frequency dynamics.
|
|
|
|
% Then, the result of the first identification is used between 10Hz and 350Hz and the result of the second identification if used between 350Hz and 2kHz.
|
|
|
|
% Here are the leg numbers that have been measured.
|
|
|
|
%% Numnbers of the measured legs
|
|
strut_nums = [1 2 3 4 5];
|
|
|
|
|
|
|
|
% The data are loaded for both the first and second identification:
|
|
|
|
%% First identification (low frequency noise)
|
|
leg_noise = {};
|
|
for i = 1:length(strut_nums)
|
|
leg_noise(i) = {load(sprintf('frf_data_leg_coder_%i_noise.mat', strut_nums(i)), 'u', 'Vs', 'de', 'da')};
|
|
end
|
|
|
|
%% Second identification (high frequency noise)
|
|
leg_noise_hf = {};
|
|
for i = 1:length(strut_nums)
|
|
leg_noise_hf(i) = {load(sprintf('frf_data_leg_coder_%i_noise_hf.mat', strut_nums(i)), 'u', 'Vs', 'de', 'da')};
|
|
end
|
|
|
|
Ts = 1e-4; % Sampling Time [s]
|
|
Nfft = floor(1/Ts);
|
|
win = hanning(Nfft);
|
|
Noverlap = floor(Nfft/2);
|
|
|
|
|
|
|
|
% We get the frequency vector that will be the same for all the frequency domain analysis.
|
|
|
|
% Only used to have the frequency vector "f"
|
|
[~, f] = tfestimate(leg_noise{1}.u, leg_noise{1}.de, win, Noverlap, Nfft, 1/Ts);
|
|
i_lf = f <= 350;
|
|
i_hf = f > 350;
|
|
|
|
% FRF Identification - Encoder
|
|
% In this section, the dynamics from $u$ to $d_e$ (encoder) is identified.
|
|
|
|
% Then, the transfer function from the DAC output voltage $u$ to the measured displacement by the encoder $d_e$ is computed:
|
|
|
|
%% Transfer function estimation
|
|
enc_frf = zeros(length(f), length(strut_nums));
|
|
|
|
for i = 1:length(strut_nums)
|
|
[frf_lf, ~] = tfestimate(leg_noise{i}.u, detrend(leg_noise{i}.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[frf_hf, ~] = tfestimate(leg_noise_hf{i}.u, detrend(leg_noise_hf{i}.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
enc_frf(:, i) = [frf_lf(i_lf); frf_hf(i_hf)];
|
|
end
|
|
|
|
|
|
|
|
% The obtained transfer functions are shown in Figure ref:fig:struts_frf_dvf_plant_tf.
|
|
|
|
%% Bode plot of the FRF from u to de
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, abs(enc_frf(:, i)), ...
|
|
'DisplayName', sprintf('Leg %i', strut_nums(i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 2);
|
|
ylim([1e-8, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, 180/pi*angle(enc_frf(:, i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% FRF Identification - Interferometer
|
|
% In this section, the dynamics from $u$ to $d_a$ (interferometer) is identified.
|
|
|
|
% Then, the transfer function from the DAC output voltage $u$ to the measured displacement by the Attocube is computed for all the struts and shown in Figure ref:fig:struts_frf_int_plant_tf.
|
|
% All the struts are giving very similar FRF.
|
|
|
|
%% Transfer function estimation
|
|
int_frf = zeros(length(f), length(strut_nums));
|
|
for i = 1:length(strut_nums)
|
|
[frf_lf, ~] = tfestimate(leg_noise{i}.u, leg_noise{i}.da, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_hf, ~] = tfestimate(leg_noise_hf{i}.u, leg_noise_hf{i}.da, win, Noverlap, Nfft, 1/Ts);
|
|
int_frf(:, i) = [frf_lf(i_lf); frf_hf(i_hf)];
|
|
end
|
|
|
|
%% Plot the FRF from u to de (interferometer)
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, abs(int_frf(:, i)), ...
|
|
'DisplayName', sprintf('Leg %i', strut_nums(i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_a/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 2);
|
|
ylim([1e-9, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, 180/pi*angle(int_frf(:, i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% FRF Identification - Force Sensor
|
|
% In this section, the dynamics from $u$ to $V_s$ is identified.
|
|
% Then the FRF are estimated and shown in Figure ref:fig:struts_frf_iff_plant_tf.
|
|
% They are also shown all to be very similar.
|
|
|
|
%% FRF estimation of the transfer function from u to Vs
|
|
iff_frf = zeros(length(f), length(strut_nums));
|
|
for i = 1:length(strut_nums)
|
|
[frf_lf, ~] = tfestimate(leg_noise{i}.u, leg_noise{i}.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
[frf_hf, ~] = tfestimate(leg_noise_hf{i}.u, leg_noise_hf{i}.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
iff_frf(:, i) = [frf_lf(i_lf); frf_hf(i_hf)];
|
|
end
|
|
|
|
%% Plot the FRF from u to Vs
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, abs(iff_frf(:, i)), ...
|
|
'DisplayName', sprintf('Leg %i', strut_nums(i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $V_s/u$ [V/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-2, 1e2]);
|
|
legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 2);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, 180/pi*angle(iff_frf(:, i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% Misalignment of the APA and flexible joints
|
|
|
|
% The misalignment between the two flexible joints and the APA has been measured for all the struts:
|
|
% - the strut is fixed to the mounting bench
|
|
% - using an indicator, the height difference from the flexible joints and the APA is measured both for the top and bottom joints and on both sides
|
|
% - then it is possible to obtain the misalignment for both flexible joints
|
|
|
|
% The raw measurements are shown in Table ref:tab:meas_misalignment_struts_raw.
|
|
|
|
% As the flexible joint's "thickness" is 1mm larger than the APA "thickness", ideally (i.e. if it were perfectly centered) we would measure =-0.50mm= each time.
|
|
|
|
|
|
strut_nums = [1, 2, 3, 4, 5];
|
|
|
|
% R Top B Top R Bot B Bot
|
|
strut_align = [[-0.40, -0.60, -0.16, -0.82] % Strut 1
|
|
[-0.67, -0.30, -0.34, -0.63] % Strut 2
|
|
[-0.07, -0.88, -0.16, -0.79] % Strut 3
|
|
[-0.48, -0.46, 0.07, -1.00] % Strut 4
|
|
[-0.33, -0.64, -0.48, -0.52]]; % Strut 5
|
|
|
|
%% Save the estimated FRF for further analysis
|
|
save('./mat/meas_struts_frf.mat', 'f', 'enc_frf', 'int_frf', 'iff_frf', 'strut_nums', 'strut_align');
|
|
|
|
% Measured misalignment of the APA and flexible joints
|
|
% The misalignment between the APA and the flexible joints are measured.
|
|
|
|
% The results are defined below and summarized in Table ref:tab:meas_misalignment_struts_new_raw.
|
|
|
|
|
|
% R Top B Top R Bot B Bot
|
|
strut_align = [[-0.54, -0.50, -0.50, -0.52] % strut 1
|
|
[-0.44, -0.55, -0.49, -0.49] % strut 2
|
|
[-0.48, -0.50, -0.50, -0.46] % strut 3
|
|
[-0.45, -0.51, -0.51, -0.45] % strut 4
|
|
[-0.50, -0.50, -0.50, -0.50] % strut 5
|
|
[-0.50, -0.49, -0.43, -0.54]]; % strut 6
|
|
|
|
% FRF Identification - Setup
|
|
% The excitation signal is a low pass filtered white noise.
|
|
% Both the encoder and the force sensor voltage are measured.
|
|
|
|
% Here are the leg numbers that have been measured.
|
|
|
|
%% Numnbers of the measured legs
|
|
strut_nums = [1 2 3 4 5 6];
|
|
|
|
%% First identification (low frequency noise)
|
|
leg_noise = {};
|
|
for i = 1:length(strut_nums)
|
|
leg_noise(i) = {load(sprintf('frf_struts_align_%i_noise.mat', strut_nums(i)), 'u', 'Vs', 'de')};
|
|
end
|
|
|
|
Ts = 1e-4; % Sampling Time [s]
|
|
Nfft = floor(1/Ts);
|
|
win = hanning(Nfft);
|
|
Noverlap = floor(Nfft/2);
|
|
|
|
|
|
|
|
% We get the frequency vector that will be the same for all the frequency domain analysis.
|
|
|
|
% Only used to have the frequency vector "f"
|
|
[~, f] = tfestimate(leg_noise{1}.u, leg_noise{1}.de, win, Noverlap, Nfft, 1/Ts);
|
|
|
|
% FRF Identification - Encoder
|
|
% In this section, the dynamics from $u$ to $d_e$ (encoder) is identified.
|
|
|
|
% Then, the transfer function from the DAC output voltage $u$ to the measured displacement by the encoder $d_e$ is computed:
|
|
|
|
%% Transfer function estimation
|
|
enc_frf = zeros(length(f), length(strut_nums));
|
|
|
|
for i = 1:length(strut_nums)
|
|
enc_frf(:, i) = tfestimate(leg_noise{i}.u, leg_noise{i}.de, win, Noverlap, Nfft, 1/Ts);
|
|
end
|
|
|
|
%% Transfer function estimation
|
|
iff_frf = zeros(length(f), length(strut_nums));
|
|
|
|
for i = 1:length(strut_nums)
|
|
iff_frf(:, i) = tfestimate(leg_noise{i}.u, leg_noise{i}.Vs, win, Noverlap, Nfft, 1/Ts);
|
|
end
|
|
|
|
|
|
|
|
% The obtained transfer functions are shown in Figure ref:fig:struts_align_frf_dvf_plant_tf.
|
|
|
|
%% Bode plot of the FRF from u to de
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, abs(enc_frf(:, i)), ...
|
|
'DisplayName', sprintf('Leg %i', strut_nums(i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 2);
|
|
ylim([1e-8, 1e-3]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(strut_nums)
|
|
plot(f, 180/pi*angle(enc_frf(:, i)));
|
|
end
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
yticks(-360:90:360); ylim([-180, 180]);
|
|
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([10, 2e3]);
|
|
|
|
% TODO Noise measurement :noexport:
|
|
|
|
%% Nothing connected to the actuator stacks
|
|
open_circuit = load('frf_struts_align_3_huddle_open_circuit.mat', 't', 'Vs', 'de');
|
|
|
|
%% PD200 connected but its input short-circuited
|
|
mid_voltage = load('frf_struts_align_3_huddle_mid_voltage_dac.mat', 't', 'Vs', 'de');
|
|
|
|
%% PD200 connected to the DAC that outputs 0V
|
|
zero_voltage = load('frf_struts_align_3_huddle_dac_zero.mat', 't', 'Vs', 'de');
|
|
|
|
%% PD200 connected to the DAC that outputs 3.25V
|
|
short_circuit = load('frf_struts_align_3_huddle_amp_short_circuit.mat', 't', 'Vs', 'de');
|
|
|
|
Ts = 1e-4; % Sampling Time [s]
|
|
Nfft = floor(2/Ts);
|
|
win = hanning(Nfft);
|
|
Noverlap = floor(Nfft/2);
|
|
|
|
[pxx_oc, f] = pwelch(detrend(open_circuit.Vs, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_mv, ~] = pwelch(detrend(mid_voltage.Vs, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_zv, ~] = pwelch(detrend(zero_voltage.Vs, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_sc, ~] = pwelch(detrend(short_circuit.Vs, 0), win, Noverlap, Nfft, 1/Ts);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, sqrt(pxx_oc), 'DisplayName', 'Open Circuit')
|
|
plot(f, sqrt(pxx_sc), 'DisplayName', 'Amp Short-Circuited')
|
|
plot(f, sqrt(pxx_zv), 'DisplayName', 'Zero Voltage (DAC)')
|
|
plot(f, sqrt(pxx_mv), 'DisplayName', 'Mid Voltage (DAC)')
|
|
hold off;
|
|
xlabel('Frequency [Hz]'); ylabel('ASD [$V/\sqrt{Hz}$]');
|
|
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
|
|
legend('location', 'northeast');
|
|
xlim([1, 5e3]);
|
|
|
|
[pxx_oc, f] = pwelch(detrend(open_circuit.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_mv, ~] = pwelch(detrend(mid_voltage.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_zv, ~] = pwelch(detrend(zero_voltage.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
[pxx_sc, ~] = pwelch(detrend(short_circuit.de, 0), win, Noverlap, Nfft, 1/Ts);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, sqrt(pxx_oc), 'DisplayName', 'Open Circuit')
|
|
plot(f, sqrt(pxx_sc), 'DisplayName', 'Amp Short-Circuited')
|
|
plot(f, sqrt(pxx_zv), 'DisplayName', 'Zero Voltage (DAC)')
|
|
plot(f, sqrt(pxx_mv), 'DisplayName', 'Mid Voltage (DAC)')
|
|
hold off;
|
|
xlabel('Frequency [Hz]'); ylabel('ASD [$m/\sqrt{Hz}$]');
|
|
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
|
|
legend('location', 'northeast');
|
|
xlim([1, 5e3])
|