From 5cfb7b8de17b99733f188fd6d9c69cc0b2bac537 Mon Sep 17 00:00:00 2001 From: Thomas Dehaeze Date: Wed, 21 Aug 2019 16:34:47 +0200 Subject: [PATCH] Update matlab file --- matlab/index.html | 301 +++++++++++++++++++++++----------------------- matlab/index.org | 20 +-- 2 files changed, 163 insertions(+), 158 deletions(-) diff --git a/matlab/index.html b/matlab/index.html index 0079520..98cd2c2 100644 --- a/matlab/index.html +++ b/matlab/index.html @@ -3,7 +3,7 @@ "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"> - + On the Design of Complementary Filters for Control - Computation with Matlab @@ -279,52 +279,52 @@ for the JavaScript code in this tag.

Table of Contents

@@ -344,31 +344,31 @@ To achieve this, the sensors included in the filter should complement one anothe

When blending two sensors using complementary filters with unknown dynamics, phase lag may be introduced that renders the close-loop system unstable.

Then, three design methods for generating two complementary filters are proposed:

-
-

1 Optimal Sensor Fusion for noise characteristics

+
+

1 Optimal Sensor Fusion for noise characteristics

- +

The idea is to combine sensors that works in different frequency range using complementary filters. @@ -389,11 +389,11 @@ All the files (data and Matlab scripts) are accessible -

1.1 Architecture

+
+

1.1 Architecture

-Let's consider the sensor fusion architecture shown on figure 1 where two sensors 1 and 2 are measuring the same quantity \(x\) with different noise characteristics determined by \(W_1\) and \(W_2\). +Let's consider the sensor fusion architecture shown on figure 1 where two sensors 1 and 2 are measuring the same quantity \(x\) with different noise characteristics determined by \(W_1\) and \(W_2\).

@@ -401,19 +401,20 @@ Let's consider the sensor fusion architecture shown on figure -

fusion_two_noisy_sensors_with_dyn.png +

-We consider that the two sensor dynamics \(G_1\) and \(G_2\) are ideal (\(G_1 = G_2 = 1\)). We obtain the architecture of figure 2. +We consider that the two sensor dynamics \(G_1\) and \(G_2\) are ideal (\(G_1 = G_2 = 1\)). We obtain the architecture of figure 2.

-
-

fusion_two_noisy_sensors.png +

+

fusion_two_noisy_sensors.png

Figure 2: Fusion of two sensors with ideal dynamics

@@ -447,8 +448,8 @@ For that, we will use the \(\mathcal{H}_2\) Synthesis.
-
-

1.2 Noise of the sensors

+
+

1.2 Noise of the sensors

Let's define the noise characteristics of the two sensors by choosing \(W_1\) and \(W_2\): @@ -468,7 +469,7 @@ W2 = ( +

nosie_characteristics_sensors.png

Figure 3: Noise Characteristics of the two sensors (png, pdf)

@@ -476,16 +477,16 @@ W2 = ( -

1.3 H-Two Synthesis

+
+

1.3 H-Two Synthesis

-We use the generalized plant architecture shown on figure 4. +We use the generalized plant architecture shown on figure 4.

-
-

h_infinity_optimal_comp_filters.png +

+

h_infinity_optimal_comp_filters.png

Figure 4: \(\mathcal{H}_2\) Synthesis - Generalized plant used for the optimal generation of complementary filters

@@ -502,7 +503,7 @@ Thus, if we minimize the \(\mathcal{H}_2\) norm of this transfer function, we mi

-We define the generalized plant \(P\) on matlab as shown on figure 4. +We define the generalized plant \(P\) on matlab as shown on figure 4.

P = [0   W2  1;
@@ -532,30 +533,30 @@ Finally, we define \(H_2 = 1 - H_1\).
 
-
-

1.4 Analysis

+
+

1.4 Analysis

-The complementary filters obtained are shown on figure 5. The PSD of the 6. -Finally, the RMS value of \(\hat{x}\) is shown on table 1. +The complementary filters obtained are shown on figure 5. The PSD of the 6. +Finally, the RMS value of \(\hat{x}\) is shown on table 1. The optimal sensor fusion has permitted to reduced the RMS value of the estimation error by a factor 8 compare to when using only one sensor.

-
+

htwo_comp_filters.png

Figure 5: Obtained complementary filters using the \(\mathcal{H}_2\) Synthesis (png, pdf)

-
+

psd_sensors_htwo_synthesis.png

Figure 6: Power Spectral Density of the estimated \(\hat{x}\) using the two sensors alone and using the optimally fused signal (png, pdf)

- +
@@ -589,18 +590,19 @@ The optimal sensor fusion has permitted to reduced the RMS value of the estimati -
-

2 Robustness to sensor dynamics uncertainty

+ +
+

2 Robustness to sensor dynamics uncertainty

- +

-Let's first consider ideal sensors where \(G_1 = 1\) and \(G_2 = 1\) (figure 7). +Let's first consider ideal sensors where \(G_1 = 1\) and \(G_2 = 1\) (figure 7).

-
+

fusion_two_noisy_sensors_with_dyn_bis.png

Figure 7: Fusion of two sensors

@@ -630,17 +632,17 @@ All the files (data and Matlab scripts) are accessible -

2.1 Unknown sensor dynamics dynamics

+
+

2.1 Unknown sensor dynamics dynamics

In practical systems, the sensor dynamics has always some level of uncertainty. -Let's represent that with multiplicative input uncertainty as shown on figure 8. +Let's represent that with multiplicative input uncertainty as shown on figure 8.

-
-

fusion_gain_mismatch.png +

+

fusion_gain_mismatch.png

Figure 8: Fusion of two sensors with input multiplicative uncertainty

@@ -684,7 +686,7 @@ Which is approximately the same as requiring

-The uncertainty set of the transfer function from \(\hat{x}\) to \(x\) is bounded in the complex plane by a circle centered on 1 and with a radius equal to \(\epsilon\) (figure 9). +The uncertainty set of the transfer function from \(\hat{x}\) to \(x\) is bounded in the complex plane by a circle centered on 1 and with a radius equal to \(\epsilon\) (figure 9).

@@ -693,8 +695,8 @@ We then have that the angle introduced by the super sensor is bounded by \(\arcs

-
-

uncertainty_gain_phase_variation.png +

+

uncertainty_gain_phase_variation.png

Figure 9: Maximum phase variation

@@ -705,8 +707,8 @@ Thus, we choose should choose \(\epsilon\) so that the maximum phase uncertainty
-
-

2.2 Design the complementary filters in order to limit the phase and gain uncertainty of the super sensor

+
+

2.2 Design the complementary filters in order to limit the phase and gain uncertainty of the super sensor

Let's say the two sensors dynamics \(H_1\) and \(H_2\) have been identified with the associated uncertainty weights \(W_1\) and \(W_2\). @@ -739,8 +741,8 @@ This is of primary importance in order to ensure the stability of the feedback l

-
-

2.3 First Basic Example with gain mismatch

+
+

2.3 First Basic Example with gain mismatch

Let's consider two ideal sensors except one sensor has not an expected gain of one but a gain of \(0.6\). @@ -752,19 +754,19 @@ G2 = 0.10. +Let's design two complementary filters as shown on figure 10. The complementary filters shown in blue does not present a bump as the red ones but provides less sensor separation at high and low frequencies.

-
+

comp_filters_robustness_test.png

Figure 10: The two complementary filters designed for the robustness test (png, pdf)

-We then compute the bode plot of the super sensor transfer function \(H_1*G_1 + H_2*G_2\) for both complementary filters pair (figure 11). +We then compute the bode plot of the super sensor transfer function \(H_1*G_1 + H_2*G_2\) for both complementary filters pair (figure 11).

@@ -772,7 +774,7 @@ We see that the blue complementary filters with a lower maximum norm permits to

-
+

tf_super_sensor_comp.png

Figure 11: Comparison of the obtained super sensor transfer functions (png, pdf)

@@ -780,15 +782,15 @@ We see that the blue complementary filters with a lower maximum norm permits to
-
-

2.4 TODO More Complete example with model uncertainty

+
+

2.4 TODO More Complete example with model uncertainty

-
-

3 Complementary filters using analytical formula

+
+

3 Complementary filters using analytical formula

- +

@@ -798,8 +800,8 @@ All the files (data and Matlab scripts) are accessible -

3.1 Analytical 1st order complementary filters

+
+

3.1 Analytical 1st order complementary filters

First order complementary filters are defined with following equations: @@ -810,7 +812,7 @@ First order complementary filters are defined with following equations: \end{align}

-Their bode plot is shown figure 12. +Their bode plot is shown figure 12.

@@ -822,7 +824,7 @@ Hl1 = 1
-
+

comp_filter_1st_order.png

Figure 12: Bode plot of first order complementary filter (png, pdf)

@@ -830,8 +832,8 @@ Hl1 = 1
-
-

3.2 Second Order Complementary Filters

+
+

3.2 Second Order Complementary Filters

We here use analytical formula for the complementary filters \(H_L\) and \(H_H\). @@ -857,12 +859,12 @@ where:

-This is illustrated on figure 13. +This is illustrated on figure 13. The slope of those filters at high and low frequencies is \(-2\) and \(2\) respectively for \(H_L\) and \(H_H\).

-
+

comp_filters_param_alpha.png

Figure 13: Effect of the parameter \(\alpha\) on the shape of the generated second order complementary filters (png, pdf)

@@ -880,7 +882,7 @@ xlabel('$ +

param_alpha_hinf_norm.png

Figure 14: Evolution of the H-Infinity norm of the complementary filters with the parameter \(\alpha\) (png, pdf)

@@ -888,8 +890,8 @@ xlabel('$ -

3.3 Third Order Complementary Filters

+
+

3.3 Third Order Complementary Filters

The following formula gives complementary filters with slopes of \(-3\) and \(3\): @@ -908,7 +910,7 @@ The parameters are:

-The filters are defined below and the result is shown on figure 15. +The filters are defined below and the result is shown on figure 15.

@@ -922,7 +924,7 @@ Hl3_ana = ( +

complementary_filters_third_order.png

Figure 15: Third order complementary filters using the analytical formula (png, pdf)

@@ -931,11 +933,11 @@ Hl3_ana = ( -

4 H-Infinity synthesis of complementary filters

+
+

4 H-Infinity synthesis of complementary filters

- +

@@ -945,8 +947,8 @@ All the files (data and Matlab scripts) are accessible -

4.1 Synthesis Architecture

+
+

4.1 Synthesis Architecture

We here synthesize the complementary filters using the \(\mathcal{H}_\infty\) synthesis. @@ -954,18 +956,18 @@ The goal is to specify upper bounds on the norms of \(H_L\) and \(H_H\) while en

-In order to do so, we use the generalized plant shown on figure 16 where \(w_L\) and \(w_H\) weighting transfer functions that will be used to shape \(H_L\) and \(H_H\) respectively. +In order to do so, we use the generalized plant shown on figure 16 where \(w_L\) and \(w_H\) weighting transfer functions that will be used to shape \(H_L\) and \(H_H\) respectively.

-
-

sf_hinf_filters_plant_b.png +

+

sf_hinf_filters_plant_b.png

Figure 16: Generalized plant used for the \(\mathcal{H}_\infty\) synthesis of the complementary filters

-The \(\mathcal{H}_\infty\) synthesis applied on this generalized plant will give a transfer function \(H_L\) (figure 17) such that the \(\mathcal{H}_\infty\) norm of the transfer function from \(w\) to \([z_H,\ z_L]\) is less than one: +The \(\mathcal{H}_\infty\) synthesis applied on this generalized plant will give a transfer function \(H_L\) (figure 17) such that the \(\mathcal{H}_\infty\) norm of the transfer function from \(w\) to \([z_H,\ z_L]\) is less than one: \[ \left\| \begin{array}{c} H_L w_L \\ (1 - H_L) w_H \end{array} \right\|_\infty < 1 \]

@@ -984,16 +986,16 @@ We then see that \(w_L\) and \(w_H\) can be used to shape both \(H_L\) and \(H_H

-
-

sf_hinf_filters_b.png +

+

sf_hinf_filters_b.png

Figure 17: \(\mathcal{H}_\infty\) synthesis of the complementary filters

-
-

4.2 Weights

+
+

4.2 Weights

omegab = 2*pi*9;
@@ -1004,7 +1006,7 @@ wL = (s 
+

weights_wl_wh.png

Figure 18: Weights on the complementary filters \(w_L\) and \(w_H\) and the associated performance weights (png, pdf)

@@ -1012,8 +1014,8 @@ wL = (s -

4.3 H-Infinity Synthesis

+
+

4.3 H-Infinity Synthesis

We define the generalized plant \(P\) on matlab. @@ -1055,7 +1057,7 @@ Test bounds: 0.0000 < gamma <= 1.7285

-We then define the high pass filter \(H_H = 1 - H_L\). The bode plot of both \(H_L\) and \(H_H\) is shown on figure 19. +We then define the high pass filter \(H_H = 1 - H_L\). The bode plot of both \(H_L\) and \(H_H\) is shown on figure 19.

Hh_hinf = 1 - Hl_hinf;
@@ -1064,15 +1066,15 @@ We then define the high pass filter \(H_H = 1 - H_L\). The bode plot of both \(H
 
-
-

4.4 Obtained Complementary Filters

+
+

4.4 Obtained Complementary Filters

-The obtained complementary filters are shown on figure 19. +The obtained complementary filters are shown on figure 19.

-
+

hinf_filters_results.png

Figure 19: Obtained complementary filters using \(\mathcal{H}_\infty\) synthesis (png, pdf)

@@ -1081,11 +1083,11 @@ The obtained complementary filters are shown on figure 19<
-
-

5 Feedback Control Architecture to generate Complementary Filters

+
+

5 Feedback Control Architecture to generate Complementary Filters

- +

The idea is here to use the fact that in a classical feedback architecture, \(S + T = 1\), in order to design complementary filters. @@ -1102,12 +1104,12 @@ All the files (data and Matlab scripts) are accessible -

5.1 Architecture

+
+

5.1 Architecture

-
-

complementary_filters_feedback_architecture.png +

+

complementary_filters_feedback_architecture.png

Figure 20: Architecture used to generate the complementary filters

@@ -1133,8 +1135,8 @@ Which contains two integrator and a lead. \(\omega_c\) is used to tune the cross
-
-

5.2 Loop Gain Design

+
+

5.2 Loop Gain Design

Let's first define the loop gain \(L\). @@ -1148,7 +1150,7 @@ L = (wc +

loop_gain_bode_plot.png

Figure 21: Bode plot of the loop gain \(L\) (png, pdf)

@@ -1156,8 +1158,8 @@ L = (wc -

5.3 Complementary Filters Obtained

+
+

5.3 Complementary Filters Obtained

We then compute the resulting low pass and high pass filters. @@ -1169,7 +1171,7 @@ Hh = 1/

-
+

low_pass_high_pass_filters.png

Figure 22: Low pass and High pass filters \(H_L\) and \(H_H\) for different values of \(\alpha\) (png, pdf)

@@ -1178,16 +1180,16 @@ Hh = 1/
-
-

6 Analytical Formula found in the literature

+
+

6 Analytical Formula found in the literature

- +

-
-

6.1 Analytical Formula

+
+

6.1 Analytical Formula

min15_compl_filter_desig_angle_estim @@ -1236,8 +1238,8 @@ Hh = 1/

-
-

6.2 Matlab

+
+

6.2 Matlab

omega0 = 1*2*pi; % [rad/s]
@@ -1255,7 +1257,7 @@ HL3 = (
+

comp_filters_literature.png

Figure 23: Comparison of some complementary filters found in the literature (png, pdf)

@@ -1263,8 +1265,8 @@ HL3 = ( -

6.3 Discussion

+
+

6.3 Discussion

Analytical Formula found in the literature provides either no parameter for tuning the robustness / performance trade-off. @@ -1273,11 +1275,11 @@ Analytical Formula found in the literature provides either no parameter for tuni

-
-

7 Comparison of the different methods of synthesis

+
+

7 Comparison of the different methods of synthesis

- + The generated complementary filters using \(\mathcal{H}_\infty\) and the analytical formulas are very close to each other. However there is some difference to note here:

    @@ -1286,6 +1288,7 @@ The generated complementary filters using \(\mathcal{H}_\infty\) and the analyti
+

Bibliography

@@ -1300,7 +1303,7 @@ The generated complementary filters using \(\mathcal{H}_\infty\) and the analyti

Author: Thomas Dehaeze

-

Created: 2019-08-14 mer. 12:13

+

Created: 2019-08-21 mer. 13:17

Validate

diff --git a/matlab/index.org b/matlab/index.org index 1419b8d..964afc1 100644 --- a/matlab/index.org +++ b/matlab/index.org @@ -1,4 +1,4 @@ -#+TITLE: On the Design of Complementary Filters for Control - Computation with Matlab +#+TITLE: Robust and Optimal Sensor Fusion - Matlab Computation :DRAWER: #+HTML_LINK_HOME: ../index.html #+HTML_LINK_UP: ../index.html @@ -92,13 +92,13 @@ $n_1$ and $n_2$ are white noise (constant power spectral density over all freque #+name: fig:fusion_two_noisy_sensors_with_dyn #+caption: Fusion of two sensors -[[file:figs/fusion_two_noisy_sensors_with_dyn.png]] +[[file:figs-tikz/fusion_two_noisy_sensors_with_dyn.png]] We consider that the two sensor dynamics $G_1$ and $G_2$ are ideal ($G_1 = G_2 = 1$). We obtain the architecture of figure [[fig:fusion_two_noisy_sensors]]. #+name: fig:fusion_two_noisy_sensors #+caption: Fusion of two sensors with ideal dynamics -[[file:figs/fusion_two_noisy_sensors.png]] +[[file:figs-tikz/fusion_two_noisy_sensors.png]] $H_1$ and $H_2$ are complementary filters ($H_1 + H_2 = 1$). The goal is to design $H_1$ and $H_2$ such that the effect of the noise sources $n_1$ and $n_2$ has the smallest possible effect on the estimation $\hat{x}$. @@ -162,7 +162,7 @@ We use the generalized plant architecture shown on figure [[fig:h_infinity_optim #+name: fig:h_infinity_optimal_comp_filters #+caption: $\mathcal{H}_2$ Synthesis - Generalized plant used for the optimal generation of complementary filters -[[file:figs/h_infinity_optimal_comp_filters.png]] +[[file:figs-tikz/h_infinity_optimal_comp_filters.png]] The transfer function from $[n_1, n_2]$ to $\hat{x}$ is: \[ \begin{bmatrix} W_1 H_1 \\ W_2 (1 - H_1) \end{bmatrix} \] @@ -249,6 +249,7 @@ The optimal sensor fusion has permitted to reduced the RMS value of the estimati | Sensor 1 | 1.1e-02 | | Sensor 2 | 1.3e-03 | | Optimal Sensor Fusion | 1.5e-04 | + * Robustness to sensor dynamics uncertainty :PROPERTIES: :header-args:matlab+: :tangle matlab/comp_filter_robustness.m @@ -327,7 +328,7 @@ Let's represent that with multiplicative input uncertainty as shown on figure [[ #+name: fig:fusion_gain_mismatch #+caption: Fusion of two sensors with input multiplicative uncertainty -[[file:figs/fusion_gain_mismatch.png]] +[[file:figs-tikz/fusion_gain_mismatch.png]] We have: \begin{align*} @@ -360,7 +361,7 @@ We then have that the angle introduced by the super sensor is bounded by $\arcsi #+name: fig:uncertainty_gain_phase_variation #+caption: Maximum phase variation -[[file:figs/uncertainty_gain_phase_variation.png]] +[[file:figs-tikz/uncertainty_gain_phase_variation.png]] Thus, we choose should choose $\epsilon$ so that the maximum phase uncertainty introduced by the sensors is of an acceptable value. @@ -765,7 +766,7 @@ In order to do so, we use the generalized plant shown on figure [[fig:sf_hinf_fi #+name: fig:sf_hinf_filters_plant_b #+caption: Generalized plant used for the $\mathcal{H}_\infty$ synthesis of the complementary filters -[[file:figs/sf_hinf_filters_plant_b.png]] +[[file:figs-tikz/sf_hinf_filters_plant_b.png]] The $\mathcal{H}_\infty$ synthesis applied on this generalized plant will give a transfer function $H_L$ (figure [[fig:sf_hinf_filters_b]]) such that the $\mathcal{H}_\infty$ norm of the transfer function from $w$ to $[z_H,\ z_L]$ is less than one: \[ \left\| \begin{array}{c} H_L w_L \\ (1 - H_L) w_H \end{array} \right\|_\infty < 1 \] @@ -782,7 +783,7 @@ We then see that $w_L$ and $w_H$ can be used to shape both $H_L$ and $H_H$ while #+name: fig:sf_hinf_filters_b #+caption: $\mathcal{H}_\infty$ synthesis of the complementary filters -[[file:figs/sf_hinf_filters_b.png]] +[[file:figs-tikz/sf_hinf_filters_b.png]] ** Weights @@ -934,7 +935,7 @@ Thus, all the tools that has been developed for classical feedback control can b ** Architecture #+name: fig:complementary_filters_feedback_architecture #+caption: Architecture used to generate the complementary filters -[[file:figs/complementary_filters_feedback_architecture.png]] +[[file:figs-tikz/complementary_filters_feedback_architecture.png]] We have: \[ y = \underbrace{\frac{L}{L + 1}}_{H_L} y_1 + \underbrace{\frac{1}{L + 1}}_{H_H} y_2 \] @@ -1138,6 +1139,7 @@ Analytical Formula found in the literature provides either no parameter for tuni The generated complementary filters using $\mathcal{H}_\infty$ and the analytical formulas are very close to each other. However there is some difference to note here: - the analytical formula provides a very simple way to generate the complementary filters (and thus the controller), they could even be used to tune the controller online using the parameters $\alpha$ and $\omega_0$. However, these formula have the property that $|H_H|$ and $|H_L|$ are symmetrical with the frequency $\omega_0$ which may not be desirable. - while the $\mathcal{H}_\infty$ synthesis of the complementary filters is not as straightforward as using the analytical formula, it provides a more optimized procedure to obtain the complementary filters + * Bibliography :ignore: bibliographystyle:unsrt bibliography:ref.bib
Table 1: RMS value of the estimation error when using the sensor individually and when using the two sensor merged using the optimal complementary filters