update scripts
This commit is contained in:
		
							
								
								
									
										
											BIN
										
									
								
								matlab/frf_data.mat
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										
											BIN
										
									
								
								matlab/frf_data.mat
									
									
									
									
									
										Normal file
									
								
							
										
											Binary file not shown.
										
									
								
							
							
								
								
									
										94
									
								
								matlab/identif_analyze_all.m
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										94
									
								
								matlab/identif_analyze_all.m
									
									
									
									
									
										Normal file
									
								
							@@ -0,0 +1,94 @@
 | 
			
		||||
%% Clear Workspace and Close figures
 | 
			
		||||
clear; close all; clc;
 | 
			
		||||
 | 
			
		||||
%% Intialize Laplace variable
 | 
			
		||||
s = zpk('s');
 | 
			
		||||
 | 
			
		||||
addpath('./src/');
 | 
			
		||||
 | 
			
		||||
% Test with one APA
 | 
			
		||||
 | 
			
		||||
%% Load measurement data for APA number 1
 | 
			
		||||
meas_data = {};
 | 
			
		||||
for i = 1:6
 | 
			
		||||
    meas_data(i) = {load(sprintf('mat/frf_data_exc_strut_%i_noise_lf.mat', i), 't', 'Va', 'Vs', 'de')};
 | 
			
		||||
%     load(sprintf('mat/frf_data_exc_strut_%i_noise_hf.mat', i), 't', 'Va', 'Vs', 'de');
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
%%
 | 
			
		||||
Ts = (meas_data{1}.t(end) - (meas_data{1}.t(1)))/(length(meas_data{1}.t)-1);
 | 
			
		||||
Fs = 1/Ts;
 | 
			
		||||
 | 
			
		||||
win = hanning(ceil(1*Fs)); % Hannning Windows
 | 
			
		||||
 | 
			
		||||
%% DVF
 | 
			
		||||
[~, f] = tfestimate(meas_data{1}.Va, meas_data{1}.de, win, [], [], 1/Ts);
 | 
			
		||||
 | 
			
		||||
G_dvf = zeros(length(f), 6, 6);
 | 
			
		||||
for i = 1:6
 | 
			
		||||
    G_dvf(:, :, i) = tfestimate(meas_data{i}.Va, meas_data{i}.de, win, [], [], 1/Ts);
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
%% IFF
 | 
			
		||||
[~, f] = tfestimate(meas_data{1}.Va, meas_data{1}.de, win, [], [], 1/Ts);
 | 
			
		||||
 | 
			
		||||
G_iff = zeros(length(f), 6, 6);
 | 
			
		||||
for i = 1:6
 | 
			
		||||
    G_iff(:, :, i) = tfestimate(meas_data{i}.Va, meas_data{i}.Vs, win, [], [], 1/Ts);
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
%%
 | 
			
		||||
figure;
 | 
			
		||||
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
 | 
			
		||||
 | 
			
		||||
ax1 = nexttile([2,1]);
 | 
			
		||||
hold on;
 | 
			
		||||
for i =1:6
 | 
			
		||||
    plot(f, abs(G_dvf(:,i, i)));
 | 
			
		||||
end
 | 
			
		||||
hold off;
 | 
			
		||||
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 | 
			
		||||
ylabel('Amplitude $V_{out}/V_{in}$ [V/V]'); set(gca, 'XTickLabel',[]);
 | 
			
		||||
hold off;
 | 
			
		||||
 | 
			
		||||
ax2 = nexttile;
 | 
			
		||||
hold on;
 | 
			
		||||
for i =1:6
 | 
			
		||||
    plot(f, 180/pi*angle(G_dvf(:,i, i)));
 | 
			
		||||
end
 | 
			
		||||
hold off;
 | 
			
		||||
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 | 
			
		||||
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
 | 
			
		||||
hold off;
 | 
			
		||||
yticks(-360:90:360);
 | 
			
		||||
 | 
			
		||||
linkaxes([ax1,ax2],'x');
 | 
			
		||||
xlim([5, 5e3]);
 | 
			
		||||
 | 
			
		||||
%%
 | 
			
		||||
figure;
 | 
			
		||||
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
 | 
			
		||||
 | 
			
		||||
ax1 = nexttile([2,1]);
 | 
			
		||||
hold on;
 | 
			
		||||
for i =1:6
 | 
			
		||||
    plot(f, abs(G_iff(:,i, i)));
 | 
			
		||||
end
 | 
			
		||||
hold off;
 | 
			
		||||
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 | 
			
		||||
ylabel('Amplitude $V_{out}/V_{in}$ [V/V]'); set(gca, 'XTickLabel',[]);
 | 
			
		||||
hold off;
 | 
			
		||||
 | 
			
		||||
ax2 = nexttile;
 | 
			
		||||
hold on;
 | 
			
		||||
for i =1:6
 | 
			
		||||
    plot(f, 180/pi*angle(G_iff(:,i, i)));
 | 
			
		||||
end
 | 
			
		||||
hold off;
 | 
			
		||||
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 | 
			
		||||
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
 | 
			
		||||
hold off;
 | 
			
		||||
yticks(-360:90:360);
 | 
			
		||||
 | 
			
		||||
linkaxes([ax1,ax2],'x');
 | 
			
		||||
xlim([5, 5e3]);
 | 
			
		||||
@@ -24,7 +24,7 @@ t  = data(:, end); % Time [s]
 | 
			
		||||
 | 
			
		||||
 | 
			
		||||
% And we save this to a =mat= file:
 | 
			
		||||
strut_number = 1;
 | 
			
		||||
strut_number = 6;
 | 
			
		||||
 | 
			
		||||
% save(sprintf('mat/frf_data_exc_strut_%i_noise.mat', strut_number), 't', 'Va', 'Vs', 'de');
 | 
			
		||||
% save(sprintf('mat/frf_data_exc_strut_%i_noise_lf.mat', strut_number), 't', 'Va', 'Vs', 'de');
 | 
			
		||||
 
 | 
			
		||||
							
								
								
									
										53
									
								
								matlab/src/generateShapedNoise.m
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										53
									
								
								matlab/src/generateShapedNoise.m
									
									
									
									
									
										Normal file
									
								
							@@ -0,0 +1,53 @@
 | 
			
		||||
function [U_exc] = generateShapedNoise(args)
 | 
			
		||||
% generateShapedNoise - Generate a Shaped Noise excitation signal
 | 
			
		||||
%
 | 
			
		||||
% Syntax: [U_exc] = generateShapedNoise(args)
 | 
			
		||||
%
 | 
			
		||||
% Inputs:
 | 
			
		||||
%    - args - Optinal arguments:
 | 
			
		||||
%        - Ts              - Sampling Time                                              - [s]
 | 
			
		||||
%        - V_mean          - Mean value of the excitation voltage                       - [V]
 | 
			
		||||
%        - V_exc           - Excitation Amplitude, could be numeric or TF               - [V rms]
 | 
			
		||||
%        - t_start         - Time at which the noise begins                             - [s]
 | 
			
		||||
%        - exc_duration    - Duration of the noise                                      - [s]
 | 
			
		||||
%        - smooth_ends     - 'true' or 'false': smooth transition between 0 and V_mean  - [-]
 | 
			
		||||
 | 
			
		||||
arguments
 | 
			
		||||
    args.Ts              (1,1) double  {mustBeNumeric, mustBePositive} = 1e-4
 | 
			
		||||
    args.V_mean          (1,1) double  {mustBeNumeric} = 0
 | 
			
		||||
    args.V_exc                                         = 1
 | 
			
		||||
    args.t_start         (1,1) double  {mustBeNumeric, mustBePositive} = 5
 | 
			
		||||
    args.exc_duration    (1,1) double  {mustBeNumeric, mustBePositive} = 10
 | 
			
		||||
    args.smooth_ends           logical {mustBeNumericOrLogical} = true
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
t_noise = 0:args.Ts:args.exc_duration;
 | 
			
		||||
 | 
			
		||||
if isnumeric(args.V_exc)
 | 
			
		||||
    V_noise = args.V_mean + args.V_exc*sqrt(1/args.Ts/2)*randn(length(t_noise), 1)';
 | 
			
		||||
elseif isct(args.V_exc)
 | 
			
		||||
    V_noise = args.V_mean + lsim(args.V_exc, sqrt(1/args.Ts/2)*randn(length(t_noise), 1), t_noise)';
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
t_smooth_start = args.Ts:args.Ts:args.t_start;
 | 
			
		||||
 | 
			
		||||
V_smooth_start = zeros(size(t_smooth_start));
 | 
			
		||||
V_smooth_end   = zeros(size(t_smooth_start));
 | 
			
		||||
 | 
			
		||||
if args.smooth_ends
 | 
			
		||||
    Vd_max = args.V_mean/(0.7*args.t_start);
 | 
			
		||||
 | 
			
		||||
    V_d = zeros(size(t_smooth_start));
 | 
			
		||||
    V_d(t_smooth_start < 0.2*args.t_start) = t_smooth_start(t_smooth_start < 0.2*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
    V_d(t_smooth_start > 0.2*args.t_start & t_smooth_start < 0.7*args.t_start) = Vd_max;
 | 
			
		||||
    V_d(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) = Vd_max - (t_smooth_start(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) - 0.7*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
 | 
			
		||||
    V_smooth_start = cumtrapz(V_d)*args.Ts;
 | 
			
		||||
 | 
			
		||||
    V_smooth_end = args.V_mean - V_smooth_start;
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
V_exc = [V_smooth_start, V_noise, V_smooth_end];
 | 
			
		||||
t_exc = args.Ts*[0:1:length(V_exc)-1];
 | 
			
		||||
 | 
			
		||||
U_exc = [t_exc; V_exc];
 | 
			
		||||
							
								
								
									
										54
									
								
								matlab/src/generateSinIncreasingAmpl.m
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										54
									
								
								matlab/src/generateSinIncreasingAmpl.m
									
									
									
									
									
										Normal file
									
								
							@@ -0,0 +1,54 @@
 | 
			
		||||
function [U_exc] = generateSinIncreasingAmpl(args)
 | 
			
		||||
% generateSinIncreasingAmpl - Generate Sinus with increasing amplitude
 | 
			
		||||
%
 | 
			
		||||
% Syntax: [U_exc] = generateSinIncreasingAmpl(args)
 | 
			
		||||
%
 | 
			
		||||
% Inputs:
 | 
			
		||||
%    - args - Optinal arguments:
 | 
			
		||||
%        - Ts              - Sampling Time                                              - [s]
 | 
			
		||||
%        - V_mean          - Mean value of the excitation voltage                       - [V]
 | 
			
		||||
%        - sin_ampls       - Excitation Amplitudes                                      - [V]
 | 
			
		||||
%        - sin_freq        - Excitation Frequency                                       - [Hz]
 | 
			
		||||
%        - sin_num         - Number of period for each amplitude                        - [-]
 | 
			
		||||
%        - t_start         - Time at which the excitation begins                        - [s]
 | 
			
		||||
%        - smooth_ends     - 'true' or 'false': smooth transition between 0 and V_mean  - [-]
 | 
			
		||||
 | 
			
		||||
arguments
 | 
			
		||||
    args.Ts              (1,1) double  {mustBeNumeric, mustBePositive} = 1e-4
 | 
			
		||||
    args.V_mean          (1,1) double  {mustBeNumeric} = 0
 | 
			
		||||
    args.sin_ampls             double  {mustBeNumeric, mustBePositive} = [0.1, 0.2, 0.3]
 | 
			
		||||
    args.sin_period      (1,1) double  {mustBeNumeric, mustBePositive} = 1
 | 
			
		||||
    args.sin_num         (1,1) double  {mustBeNumeric, mustBePositive, mustBeInteger} = 3
 | 
			
		||||
    args.t_start         (1,1) double  {mustBeNumeric, mustBePositive} = 5
 | 
			
		||||
    args.smooth_ends           logical {mustBeNumericOrLogical} = true
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
t_noise = 0:args.Ts:args.sin_period*args.sin_num;
 | 
			
		||||
sin_exc = [];
 | 
			
		||||
 | 
			
		||||
for sin_ampl = args.sin_ampls
 | 
			
		||||
    sin_exc = [sin_exc, args.V_mean + sin_ampl*sin(2*pi/args.sin_period*t_noise)];
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
t_smooth_start = args.Ts:args.Ts:args.t_start;
 | 
			
		||||
 | 
			
		||||
V_smooth_start = zeros(size(t_smooth_start));
 | 
			
		||||
V_smooth_end   = zeros(size(t_smooth_start));
 | 
			
		||||
 | 
			
		||||
if args.smooth_ends
 | 
			
		||||
    Vd_max = args.V_mean/(0.7*args.t_start);
 | 
			
		||||
 | 
			
		||||
    V_d = zeros(size(t_smooth_start));
 | 
			
		||||
    V_d(t_smooth_start < 0.2*args.t_start) = t_smooth_start(t_smooth_start < 0.2*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
    V_d(t_smooth_start > 0.2*args.t_start & t_smooth_start < 0.7*args.t_start) = Vd_max;
 | 
			
		||||
    V_d(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) = Vd_max - (t_smooth_start(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) - 0.7*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
 | 
			
		||||
    V_smooth_start = cumtrapz(V_d)*args.Ts;
 | 
			
		||||
 | 
			
		||||
    V_smooth_end = args.V_mean - V_smooth_start;
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
V_exc = [V_smooth_start, sin_exc, V_smooth_end];
 | 
			
		||||
t_exc = args.Ts*[0:1:length(V_exc)-1];
 | 
			
		||||
 | 
			
		||||
U_exc = [t_exc; V_exc];
 | 
			
		||||
							
								
								
									
										76
									
								
								matlab/src/generateSweepExc.m
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										76
									
								
								matlab/src/generateSweepExc.m
									
									
									
									
									
										Normal file
									
								
							@@ -0,0 +1,76 @@
 | 
			
		||||
function [U_exc] = generateSweepExc(args)
 | 
			
		||||
% generateSweepExc - Generate a Sweep Sine excitation signal
 | 
			
		||||
%
 | 
			
		||||
% Syntax: [U_exc] = generateSweepExc(args)
 | 
			
		||||
%
 | 
			
		||||
% Inputs:
 | 
			
		||||
%    - args - Optinal arguments:
 | 
			
		||||
%        - Ts              - Sampling Time                                              - [s]
 | 
			
		||||
%        - f_start         - Start frequency of the sweep                               - [Hz]
 | 
			
		||||
%        - f_end           - End frequency of the sweep                                 - [Hz]
 | 
			
		||||
%        - V_mean          - Mean value of the excitation voltage                       - [V]
 | 
			
		||||
%        - V_exc           - Excitation Amplitude for the Sweep, could be numeric or TF - [V]
 | 
			
		||||
%        - t_start         - Time at which the sweep begins                             - [s]
 | 
			
		||||
%        - exc_duration    - Duration of the sweep                                      - [s]
 | 
			
		||||
%        - sweep_type      - 'logarithmic' or 'linear'                                  - [-]
 | 
			
		||||
%        - smooth_ends     - 'true' or 'false': smooth transition between 0 and V_mean  - [-]
 | 
			
		||||
 | 
			
		||||
arguments
 | 
			
		||||
    args.Ts              (1,1) double  {mustBeNumeric, mustBePositive} = 1e-4
 | 
			
		||||
    args.f_start         (1,1) double  {mustBeNumeric, mustBePositive} = 1
 | 
			
		||||
    args.f_end           (1,1) double  {mustBeNumeric, mustBePositive} = 1e3
 | 
			
		||||
    args.V_mean          (1,1) double  {mustBeNumeric} = 0
 | 
			
		||||
    args.V_exc                                         = 1
 | 
			
		||||
    args.t_start         (1,1) double  {mustBeNumeric, mustBeNonnegative} = 5
 | 
			
		||||
    args.exc_duration    (1,1) double  {mustBeNumeric, mustBePositive} = 10
 | 
			
		||||
    args.sweep_type            char    {mustBeMember(args.sweep_type,{'log', 'lin'})} = 'lin'
 | 
			
		||||
    args.smooth_ends           logical {mustBeNumericOrLogical} = true
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
t_sweep = 0:args.Ts:args.exc_duration;
 | 
			
		||||
 | 
			
		||||
if strcmp(args.sweep_type, 'log')
 | 
			
		||||
    V_exc = sin(2*pi*args.f_start * args.exc_duration/log(args.f_end/args.f_start) * (exp(log(args.f_end/args.f_start)*t_sweep/args.exc_duration) - 1));
 | 
			
		||||
elseif strcmp(args.sweep_type, 'lin')
 | 
			
		||||
    V_exc = sin(2*pi*(args.f_start + (args.f_end - args.f_start)/2/args.exc_duration*t_sweep).*t_sweep);
 | 
			
		||||
else
 | 
			
		||||
    error('sweep_type should either be equal to "log" or to "lin"');
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
if isnumeric(args.V_exc)
 | 
			
		||||
    V_sweep = args.V_mean + args.V_exc*V_exc;
 | 
			
		||||
elseif isct(args.V_exc)
 | 
			
		||||
    if strcmp(args.sweep_type, 'log')
 | 
			
		||||
        V_sweep = args.V_mean + abs(squeeze(freqresp(args.V_exc, args.f_start*(args.f_end/args.f_start).^(t_sweep/args.exc_duration), 'Hz')))'.*V_exc;
 | 
			
		||||
    elseif strcmp(args.sweep_type, 'lin')
 | 
			
		||||
        V_sweep = args.V_mean + abs(squeeze(freqresp(args.V_exc, args.f_start+(args.f_end-args.f_start)/args.exc_duration*t_sweep, 'Hz')))'.*V_exc;
 | 
			
		||||
    end
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
if args.t_start > 0
 | 
			
		||||
    t_smooth_start = args.Ts:args.Ts:args.t_start;
 | 
			
		||||
 | 
			
		||||
    V_smooth_start = zeros(size(t_smooth_start));
 | 
			
		||||
    V_smooth_end   = zeros(size(t_smooth_start));
 | 
			
		||||
 | 
			
		||||
    if args.smooth_ends
 | 
			
		||||
        Vd_max = args.V_mean/(0.7*args.t_start);
 | 
			
		||||
 | 
			
		||||
        V_d = zeros(size(t_smooth_start));
 | 
			
		||||
        V_d(t_smooth_start < 0.2*args.t_start) = t_smooth_start(t_smooth_start < 0.2*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
        V_d(t_smooth_start > 0.2*args.t_start & t_smooth_start < 0.7*args.t_start) = Vd_max;
 | 
			
		||||
        V_d(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) = Vd_max - (t_smooth_start(t_smooth_start > 0.7*args.t_start & t_smooth_start < 0.9*args.t_start) - 0.7*args.t_start)*Vd_max/(0.2*args.t_start);
 | 
			
		||||
 | 
			
		||||
        V_smooth_start = cumtrapz(V_d)*args.Ts;
 | 
			
		||||
 | 
			
		||||
        V_smooth_end = args.V_mean - V_smooth_start;
 | 
			
		||||
    end
 | 
			
		||||
else
 | 
			
		||||
    V_smooth_start = [];
 | 
			
		||||
    V_smooth_end = [];
 | 
			
		||||
end
 | 
			
		||||
 | 
			
		||||
V_exc = [V_smooth_start, V_sweep, V_smooth_end];
 | 
			
		||||
t_exc = args.Ts*[0:1:length(V_exc)-1];
 | 
			
		||||
 | 
			
		||||
U_exc = [t_exc; V_exc];
 | 
			
		||||
		Reference in New Issue
	
	Block a user