5.0 KiB
+++ title = "System Identification" author = ["Dehaeze Thomas"] draft = false +++
- Tags
- [Modal Analysis]({{< relref "modal_analysis.md" >}})
SISO Identification
Problem Description
{{< figure src="/ox-hugo/siso_identification_schematic.png" caption="<span class="figure-number">Figure 1: Block diagram of the SISO system identification" >}}
{{< figure src="/ox-hugo/siso_identification_schematic_simplier.png" caption="<span class="figure-number">Figure 2: Simpler Block diagram of the SISO system identification" >}}
If the open-loop system is unstable, first design a simple controller that stabilizes the system and then identify the closed-loop system.
Design of the Excitation Signal
Introduction
There are several choices for excitation signals:
- Impulse, Steps
- Sweep Sinus
- Random noise, Periodic signals
Random noise with specific ASD
The ASD of the measured output is:
\begin{equation} \Gamma_{y_m}(\omega) = \Gamma_d(\omega) + \Gamma_u(\omega) \cdot |G(j\omega)| \end{equation}
And we want the effect of the excitation signal to be much higher than the effect of the exogenous signals (measurement noise, input noise, disturbances).
\begin{equation} \Gamma_u(\omega) \gg \Gamma_d(\omega) \cdot |G(j\omega)|^{-1} \end{equation}
Note that \(\Gamma_d(\omega)\) can be estimated by measuring the system output in the absence of any excitation signal. The plant magnitude \(|G(j\omega)|\) can be roughly estimated from a first identification with bad coherence.
In order to design a random excitation signal with specific spectral characteristics, first a signal with an ASD equal to one is generated (i.e. white noise with unity ASD):
Ts = 1e-4; % Sampling Time [s]
t = 0:Ts:10; % Time Vector [s]
%% Signal with an ASD equal to one
u_norm = sqrt(1/2/Ts)*randn(length(t), 1);
Then, a transfer function whose magnitude \(|G_u(j\omega)|\) has the same shape as the wanted excitation ASD \(\Gamma_u(\omega)\) is designed:
%% Transfer function representing the wanted ASD
G_u = tf([1], [1/2/pi/100 1]);
Finally lsim is used to compute the shaped excitation signal.
%% Shape the ASD of the excitation signal
u = lsim(G_u, u_norm, t);
Choose Sampling Frequency and Duration of Excitation
The sampling frequency \(F_s\) will determine the maximum frequency \(F_{\text{max}}\) that can be estimated (see Nyquist theorem):
\begin{equation} F_{\text{max}} = \frac{1}{2} F_s \end{equation}
The duration of excitation \(T_{\text{exc}}\) will determine the minimum frequency \(F_{\text{min}}\) that can be estimated:
\begin{equation} F_{\text{min}} = \frac{1}{T_{\text{exc}}} \end{equation}
It will also corresponds to the frequency resolution \(\Delta f\):
\begin{equation} \Delta f = \frac{1}{T_{\text{exc}}} \end{equation}
In order to increase the estimation quality, averaging can be use with a longer excitation duration. A factor 10 is usually good enough, therefore the excitation time can be taken as:
\begin{equation} T_{\text{exc}} \approx \frac{10}{F_{\text{min}}} \end{equation}
Therefore, if the system has to be identified from 1Hz up to 500Hz, the sampling frequency should be:
\begin{equation} F_s = 2 F_{\text{max}} = 1\,\text{kHz} \end{equation}
Then, the excitation duration should be (10 averaging):
\begin{equation} T_{\text{exc}} = \frac{10}{1} = 10\,s \end{equation}
Computation of the Frequency Response Function
Windowing Function
Example
tfestimate
[G, f] = tfestimate(u, y, win, [], [], 1/Ts);
Verification of the Identification Quality
mscohere
[coh, f] = mscohere(u, y, win, [], [], 1/Ts);