diff --git a/C1-test-bench-apa/test_apa_1_basic_meas.m b/C1-test-bench-apa/test_apa_1_basic_meas.m
index 4bb5805..64ecd44 100644
--- a/C1-test-bench-apa/test_apa_1_basic_meas.m
+++ b/C1-test-bench-apa/test_apa_1_basic_meas.m
@@ -63,8 +63,9 @@ for i = 1:7
 end
 hold off;
 xlabel('Voltage [V]'); ylabel('Displacement [$\mu m$]')
-legend('location', 'southwest', 'FontSize', 8)
 xlim([-20, 150]); ylim([-250, 0]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 %% X-Bending Identification
 % Load Data
@@ -96,4 +97,5 @@ hold off;
 set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
 xlabel('Frequency [Hz]'); ylabel('Amplitude');
 xlim([100, 1e3]); ylim([5e-5, 2e-1]);
-legend('location', 'northeast', 'FontSize', 8)
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
diff --git a/C1-test-bench-apa/test_apa_2_dynamics.m b/C1-test-bench-apa/test_apa_2_dynamics.m
index 6f492ba..155a6ea 100644
--- a/C1-test-bench-apa/test_apa_2_dynamics.m
+++ b/C1-test-bench-apa/test_apa_2_dynamics.m
@@ -30,7 +30,8 @@ for i = [6,5,4,2]
 end
 hold off;
 xlabel('Stack Voltage $V_a$ [V]'); ylabel('Displacement $d_e$ [$\mu$m]');
-legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 xlim([-20, 150]);
 ylim([-120, 120]);
 
@@ -132,8 +133,9 @@ hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
 hold off;
-legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 2);
 ylim([1e-8, 1e-3]);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 2);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
@@ -164,7 +166,8 @@ 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);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 2);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
@@ -209,9 +212,9 @@ ylim([0, 1]);
 
 %% Bode plot of the FRF from u to de
 figure;
-tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
+tiledlayout(2, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
-ax1 = nexttile([2,1]);
+ax1 = nexttile;
 hold on;
 plot(f, abs(frf_noise), '.-');
 hold off;
@@ -258,22 +261,23 @@ tiledlayout(2, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile();
 hold on;
-plot(f, abs(frf_wo_k), 'DisplayName', 'Without $R$');
-plot(f, abs(frf_wi_k), 'DisplayName', 'With $R$');
-plot(f, abs(squeeze(freqresp(G_hpf, f, 'Hz'))), 'k--', 'DisplayName', 'RC model');
+plot(f, abs(frf_wo_k), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'Without $R$');
+plot(f, abs(frf_wi_k), '-', 'color', [colors(2,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'With $R$');
+plot(f, abs(squeeze(freqresp(G_hpf, f, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', 'RC model');
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 ylabel('Amplitude [V/V]'); set(gca, 'XTickLabel',[]);
 hold off;
 ylim([2e-1, 1e0]);
 yticks([0.2, 0.5, 1]);
-legend('location', 'southeast', 'FontSize', 8);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
-plot(f, 180/pi*angle(frf_wo_k));
-plot(f, 180/pi*angle(frf_wi_k));
-plot(f, 180/pi*angle(squeeze(freqresp(G_hpf, f, 'Hz'))), 'k--', 'DisplayName', 'RC');
+plot(f, 180/pi*angle(frf_wo_k), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5);
+plot(f, 180/pi*angle(frf_wi_k), '-', 'color', [colors(2,:), 0.5], 'linewidth', 2.5);
+plot(f, 180/pi*angle(squeeze(freqresp(G_hpf, f, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', 'RC');
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
@@ -312,23 +316,25 @@ tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile([2,1]);
 hold on;
-plot(f, abs(G_iff), 'color', colors(2,:), 'DisplayName', 'Identified plant')
-plot(f, abs(squeeze(freqresp(G_iff_model, f, 'Hz'))), 'k--', 'DisplayName', 'Manual fit')
+plot(f, abs(G_iff), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5, 'DisplayName', 'Measured FRF');
+plot(f, abs(squeeze(freqresp(G_iff_model, f, 'Hz'))), '--', 'color', colors(1,:), 'DisplayName', 'Manual fit')
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 ylabel('Amplitude $V_s/u$ [V/V]'); set(gca, 'XTickLabel',[]);
 legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+ylim([1e-2, 1e2]);
+yticks([1e-2, 1e-1, 1, 1e1, 1e2]);
 
 ax2 = nexttile;
 hold on;
-plot(f, 180/pi*angle(G_iff), 'color', colors(2,:));
-plot(f, 180/pi*angle(squeeze(freqresp(G_iff_model, f, 'Hz'))), 'k--')
+plot(f, 180/pi*angle(G_iff), '-', 'color', [colors(1,:), 0.5], 'linewidth', 2.5);
+plot(f, 180/pi*angle(squeeze(freqresp(G_iff_model, f, 'Hz'))), '--', 'color', colors(1,:))
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
 hold off;
-yticks(-360:45:360);
 ylim([-90, 180])
+yticks(-360:90:360);
 
 linkaxes([ax1,ax2],'x');
 xlim([0.2, 1e3]);
@@ -398,14 +404,15 @@ tiledlayout(1, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 ax1 = nexttile();
 hold on;
 for i = 1:length(i_kept)
-    plot(f, abs(G_dL_frf{i}), 'color', [colors(i,:), 1], 'DisplayName', sprintf('g = %.0f', data.gains(i_kept(i))))
+    plot(f, abs(G_dL_frf{i}), '-', 'color', [colors(i,:), 0.5], 'linewidth', 2.5, 'DisplayName', sprintf('g = %.0f', data.gains(i_kept(i))))
     plot(f, abs(squeeze(freqresp(G_dL_id{i}, f, 'Hz'))), '--', 'color', [colors(i,:), 1], 'HandleVisibility', 'off')
 end
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 xlabel('Frequency [Hz]'); ylabel('Amplitude $d_e/u^\prime$ [m/V]');
-xlim([10, 1e3]);
-legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+xlim([10, 800]);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 %% Root Locus of the APA300ML with Integral Force Feedback
 % Comparison between the computed root locus from the plant model and the root locus estimated from the damped plant pole identification
@@ -416,7 +423,7 @@ hold on;
 G_iff_poles = pole(pade(G_iff_model));
 i = imag(G_iff_poles) > 100; % Only keep relevant poles
 plot(real(G_iff_poles(i)),  imag(G_iff_poles(i)),  'kx', ...
-    'DisplayName', '$g = 0$');
+    'HandleVisibility', 'off');
 G_iff_zeros = tzero(G_iff_model);
 i = imag(G_iff_zeros) > 100; % Only keep relevant zeros
 plot(real(G_iff_zeros(i)), imag(G_iff_zeros(i)), 'ko', ...
@@ -437,4 +444,5 @@ ylabel('Imaginary Part')
 axis equal
 ylim([0, 610]);
 xlim([-300,0]);
-legend('location', 'southwest', 'FontSize', 8);
+leg = legend('location', 'southwest', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
diff --git a/C1-test-bench-apa/test_apa_3_model_2dof.m b/C1-test-bench-apa/test_apa_3_model_2dof.m
index a523de7..7539fe6 100644
--- a/C1-test-bench-apa/test_apa_3_model_2dof.m
+++ b/C1-test-bench-apa/test_apa_3_model_2dof.m
@@ -105,9 +105,9 @@ tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile([2,1]);
 hold on;
-plot(f, abs(enc_frf(:, 1)), 'color', [0,0,0,0.2], 'DisplayName', 'Identified');
+plot(f, abs(enc_frf(:, 1)), '-', 'color', [colors(2,:), 0.1], 'linewidth', 2.5, 'DisplayName', 'Identified');
 for i = 1:length(apa_nums)
-    plot(f, abs(enc_frf(:, i)), 'color', [0,0,0,0.2], 'HandleVisibility', 'off');
+    plot(f, abs(enc_frf(:, i)), '-', 'color', [colors(2,:), 0.1], 'linewidth', 2.5, 'HandleVisibility', 'off');
 end
 plot(freqs, abs(squeeze(freqresp(G_2dof('de', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', '2DoF Model')
 hold off;
@@ -115,12 +115,13 @@ set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
 hold off;
 ylim([1e-8, 1e-3]);
-legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
 for i = 1:length(apa_nums)
-    plot(f, 180/pi*angle(enc_frf(:, i)), 'color', [0,0,0,0.2]);
+    plot(f, 180/pi*angle(enc_frf(:, i)), '-', 'color', [colors(2,:), 0.1], 'linewidth', 2.5);
 end
 plot(freqs, 180/pi*angle(squeeze(freqresp(G_2dof('de', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:))
 hold off;
@@ -138,24 +139,25 @@ tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile([2,1]);
 hold on;
-plot(f, abs(iff_frf(:, 1)), 'color', [0,0,0,0.2], 'DisplayName', 'Identified');
+plot(f, abs(iff_frf(:, 1)), '-', 'color', [colors(1,:), 0.1], 'linewidth', 2.5, 'DisplayName', 'Identified');
 for i = 2:length(apa_nums)
-    plot(f, abs(iff_frf(:, i)), 'color', [0,0,0,0.2], 'HandleVisibility', 'off');
+    plot(f, abs(iff_frf(:, i)), '-', 'color', [colors(1,:), 0.1], 'linewidth', 2.5, 'HandleVisibility', 'off');
 end
-plot(freqs, abs(squeeze(freqresp(G_2dof('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', '2DoF Model')
+plot(freqs, abs(squeeze(freqresp(G_2dof('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(1,:), 'DisplayName', '2DoF Model')
 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', 1);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
 for i = 1:length(apa_nums)
-    plot(f, 180/pi*angle(iff_frf(:, i)), 'color', [0,0,0,0.2]);
+    plot(f, 180/pi*angle(iff_frf(:, i)), '-', 'color', [colors(1,:), 0.2], 'linewidth', 2.5);
 end
-plot(freqs, 180/pi*angle(squeeze(freqresp(G_2dof('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:))
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_2dof('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(1,:))
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
diff --git a/C1-test-bench-apa/test_apa_4_model_flexible.m b/C1-test-bench-apa/test_apa_4_model_flexible.m
index 15d0295..65a452d 100644
--- a/C1-test-bench-apa/test_apa_4_model_flexible.m
+++ b/C1-test-bench-apa/test_apa_4_model_flexible.m
@@ -71,9 +71,9 @@ tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile([2,1]);
 hold on;
-plot(f, abs(enc_frf(:, 1)), 'color', [0,0,0,0.2], 'DisplayName', 'Identified');
+plot(f, abs(enc_frf(:, 1)), 'color', [colors(2,:), 0.1], 'linewidth', 2.5, 'DisplayName', 'Identified');
 for i = 1:length(apa_nums)
-    plot(f, abs(enc_frf(:, i)), 'color', [0,0,0,0.2], 'HandleVisibility', 'off');
+    plot(f, abs(enc_frf(:, i)), 'color', [colors(2,:), 0.1], 'linewidth', 2.5, 'HandleVisibility', 'off');
 end
 plot(freqs, abs(squeeze(freqresp(G_flex('de', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', '"Flexible" Model')
 hold off;
@@ -81,12 +81,13 @@ set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
 ylabel('Amplitude $d_e/u$ [m/V]'); set(gca, 'XTickLabel',[]);
 hold off;
 ylim([1e-8, 1e-3]);
-legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg = legend('location', 'northeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
 for i = 1:length(apa_nums)
-    plot(f, 180/pi*angle(enc_frf(:, i)), 'color', [0,0,0,0.2]);
+    plot(f, 180/pi*angle(enc_frf(:, i)), 'color', [colors(2,:), 0.1], 'linewidth', 2.5);
 end
 plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex('de', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:))
 hold off;
@@ -104,24 +105,25 @@ tiledlayout(3, 1, 'TileSpacing', 'Compact', 'Padding', 'None');
 
 ax1 = nexttile([2,1]);
 hold on;
-plot(f, abs(iff_frf(:, 1)), 'color', [0,0,0,0.2], 'DisplayName', 'Identified');
+plot(f, abs(iff_frf(:, 1)), '-', 'color', [colors(1,:), 0.1], 'linewidth', 2.5, 'DisplayName', 'Identified');
 for i = 2:length(apa_nums)
-    plot(f, abs(iff_frf(:, i)), 'color', [0,0,0,0.2], 'HandleVisibility', 'off');
+    plot(f, abs(iff_frf(:, i)), '-', 'color', [colors(1,:), 0.1], 'linewidth', 2.5, 'HandleVisibility', 'off');
 end
-plot(freqs, abs(squeeze(freqresp(G_flex('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:), 'DisplayName', '"Flexible" Model')
+plot(freqs, abs(squeeze(freqresp(G_flex('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(1,:), 'DisplayName', '"Flexible" Model')
 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', 1);
+leg = legend('location', 'southeast', 'FontSize', 8, 'NumColumns', 1);
+leg.ItemTokenSize(1) = 15;
 
 ax2 = nexttile;
 hold on;
 for i = 1:length(apa_nums)
-    plot(f, 180/pi*angle(iff_frf(:, i)), 'color', [0,0,0,0.2]);
+    plot(f, 180/pi*angle(iff_frf(:, i)), '-', 'color', [colors(1,:), 0.1], 'linewidth', 2.5);
 end
-plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(2,:))
+plot(freqs, 180/pi*angle(squeeze(freqresp(G_flex('Vs', 'u'), freqs, 'Hz'))), '--', 'color', colors(1,:))
 hold off;
 set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
 xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
diff --git a/README.md b/README.md
index 4d1c91a..647b3df 100644
--- a/README.md
+++ b/README.md
@@ -74,6 +74,7 @@ The repository is organized as follows:
 ├── C3-test-bench-struts
 ├── C4-test-bench-nano-hexapod
 ├── C5-test-bench-id31
+├── dehaeze-nass-phd-thesis.pdf     # PhD Thesis Manuscript
 ├── init.m                          # Optional: Some Matlab configuration to have same figure display
 ├── LICENSE                         # MIT License is used for all the code in this repository
 └── README.md                       # The present file
diff --git a/dehaeze-nass-phd-thesis.pdf b/dehaeze-nass-phd-thesis.pdf
new file mode 100644
index 0000000..54ece65
Binary files /dev/null and b/dehaeze-nass-phd-thesis.pdf differ