Added identification data and scripts from students

data/README.md

@@ -57,4 +57,9 @@ You will get something like the figures below
 ![figures/simulated_pos_ol_response.png](figures/simulated_pos_ol_response.png)
 ![figures/measured_acc_ol_response.png](figures/measured_acc_ol_response.png)
 
-Apparently, the simulation responses do not quite fit the experimental ones. Some optimization-based fitting might be needed here.
+Apparently, the simulation responses do not quite fit the experimental ones. Some optimization-based fitting might be needed here.
\ No newline at end of file
+
+## Additional student-provided experiments
+
+The `students` subdirectory contains additional experiments and associated
+scripts provided by past students of B3M35ORR.

data/students/README.md

+# Identification experiments by students
+
+This directory is split in two parts:
+
+ * Identification data itself in `data` subdirectory. The `.mat` files
+   contain recorded responses of the physical system to different input
+   signals.
+ * Identification scripts in the `scripts` subdirectory. These scripts
+   can be used to fit a preexisting linear/nonlinear model into the
+   measured data (-> greybox identification).

data/students/data/README.md

+# Measured data
+
+This dataset contains different types of experiments:
+ * responses to steps in road height,
+ * responses to steps in actuator force,
+ * responses to sine sweep in road height,
+ * responses to sine sweep in actuator force.
+
+The data format is similar to the format used for the original identification
+dataset. The data can be loaded into Matlab workspace using
+
+```matlab
+load road_step_slow.mat
+```
+
+Up to three new data structs will then appear in the workspace:
+
+- `measured_pos_ol_response`
+  - This struct contains the measured road height and the heights of lower and upper masses.
+  - The columns correspond to `[zr zus zs]` -- road height, lower mass height and upper mass height.
+  - This struct is always present.
+- `measured_acc_ol_response`
+  - This struct contains a **low-pass-filetered** signal from the accelerometer mounted on the upper mass.
+  - The columns correspond to `[acc]` -- just the accelerometer signal.
+  - The low-pass filter is a second-order low pass filter with ωₙ = 70 rad/sec and ζ = 1.
+  - This struct is sometimes missing.
+- `control_signal`
+  - This struct contains the signals that were sent into the physical model.
+  - The columns correspond to `[zr_ref actuator_out]` -- road height reference and actuator output signal.
+  - This struct is sometimes missing. However, it is never missing when the
+    file contains a response to a change in actuator force -- this struct
+    contains which exact signal was applied.
+
+This can be demonstrated by plotting the data:
+```matlab
+figure(1) 
+plot(measured_pos_ol_response.time, measured_pos_ol_response.signals.values)
+legend('zr','zus','zs')
+xlabel('Time [s]')
+ylabel('Measured positions [m]')
+grid on
+
+figure(2)
+plot(measured_acc_ol_response.time, measured_acc_ol_response.signals.values)
+legend('d2/dt2 zs')
+xlabel('Time [s]')
+ylabel('Measured acceleration [m/s2]')
+grid on
+
+figure(3)
+plot(control_signal.time, control_signal.signals.values)
+legend('target road height', 'actuator output signal')
+xlabel('Time [s]')
+ylabel('Control signal [-]')
+grid on
+```
+
+## List of experiments:
+
+### Actuator step
+ - `actuator_step.mat` - basically that.
+
+### Actuator sweep
+ - `actuator_sweep.mat` - sweep from "middle" to "high" frequencies.
+ - `actuator_sweep_lf.mat` - sweep covering "low" frequencies.
+ - `actuator_sweep_lf_longer.mat` - sweep trying to cover the resonant peak of the model.
+ - `actuator_sweep_lf_longer_v2.mat` - sweep trying to cover the resonant peak of the model.
+ - `actuator_sweep_reversed.mat` - sweep from going from "high" to "middle" frequencies. The goal of this experiment to reduce the effect of static friction at low frequencies.
+
+### Road step
+ - `road_step_slow.mat` - response of the system when the road platform is subjected to ramped steps in the reference signal.
+ - `road_step_fast.mat` - response of the system when the road platform is subjected to immediate jumps in the reference signal.
+   This goal of this experiment was to excite the system more than in `road_step_slow.mat`.
+
+### Road sweep
+ - `road_sweep_normal.mat` - sweep from "middle" to "semi-high" frequencies.
+ - `road_sweep_normal_2.mat` - sweep from "middle" to "semi-high" frequencies.
+ - `road_sweep_long.mat` - sweep from "middle" to "high" frequencies.
+ - `road_sweep_normal_reversed.mat` - sweep from "high" to "middle" frequencies. The goal was again to reduce the effect of static friction.
+ - `road_sweep_lf.mat` - sweep covering "low" frequencies.

data/students/scripts/README.md

+# Scripts for system identification
+
+These scripts can be used to determine values of constants in predetermined
+models of the active suspension system.
+
+The format of the model corresponds to the "physical" system state described
+in the `<repo>/software/students/README.md` file.
+
+There are two main files:
+ - `ident_linear_model.m` will try to fit constants in the fully linear model.
+   This model may be most suitable for linear controller design.
+   The downside of this model is that it does not model static friction and
+   so simulations with this model sometimes does not match the reality very well.
+ - `ident_nonlinear_model.m` will try to fit constants in the nonlinear model.
+   This model may be most suitable for simulation because it at least
+   somehow models static friction.
+
+The meaning of the fitted parameters is briefly explained in the
+`as_model_lin.m` and `as_model_nl.m` files where the greybox models are
+implemented.

data/students/scripts/as_model_lin.m

+function [A, B, C, D] = as_model_lin(ks, ku, r, bs, bu, gf, Ts)
+    % AS_MODEL_LIN(...) Generate linear model matrices for greybox identification
+    %
+    % This function is intended to be given to greyest(). Please see
+    % the source code of as_model_lin() for explanation of what 
+    % its parameters mean.
+
+    ms = 2.5;
+    mu = ms / r;
+
+    % explanation of parameters:
+    % - ks = stiffness of the spring between the upper and lower masses
+    % - ku = stiffness of the spring between the lower mass and the road platform
+    % - ms = mass of the upper mass
+    % - mu = mass of the lower mass
+    % - r  = ratio of ms/mu. This is what is actually optimized - it seems
+    %        that this removes an unconstrained degree of freedom in the
+    %        parameters.
+    % - bs = coefficient of viscous friction between the upper mass and the guiding rail connected to the base
+    % - bu = coefficient of viscous friction between the lower mass and the guiding rail connected to the base
+    % - gf = relative strength of the actuator between lower and upper masses
+    % !! no static friction here, as it is nonlinear
+
+    % explanation of states/inputs/outputs:
+    % - x(1) = height of the upper mass above the base, minus its idle position
+    % - x(2) = speed of the upper mass wrt. the base
+    % - x(3) = height of the lower mass above the base, minus its idle position
+    % - x(4) = speed of the lower mass wrt. the base
+    % - u(1) = road height above the base, minus its idle position
+    % - u(2) = actuator control signal
+    % - y(1) = x(1)
+    % - y(2) = x(3)
+    
+    A = [
+             0,      1,           0,      0;
+        -ks/ms, -bs/ms,      +ks/ms,      0;
+             0,      0,           0,      1;
+        +ks/mu,      0, (-ks-ku)/mu, -bu/mu;
+    ];
+    B = [
+             0,      0;
+             0, +gf/ms;
+             0,      0;
+        +ku/mu, -gf/mu;
+    ];
+    C = [
+        1, 0, 0, 0;
+        0, 0, 1, 0;
+    ];
+    D = [
+        0, 0;
+        0, 0;
+    ];
+end

data/students/scripts/as_model_nl.m

+function [dx, y] = as_model_nl(t, x, u, ks, ku, r, bs, bu, fs, fu, gf, FileArgument)
+    % AS_MODEL_NL(t, x, u, ...) Simulate a system step for nonlinear greybox identification
+    % 
+    % This function is intended to be given to nlgreyest(). Please see
+    % the source code of as_model_nl() for explanation of what 
+    % its parameters mean.
+
+    ms = 2.5;
+    mu = ms / r;
+    
+    % explanation of parameters:
+    % - ks = stiffness of the spring between the upper and lower masses
+    % - ku = stiffness of the spring between the lower mass and the road platform
+    % - ms = mass of the upper mass
+    % - mu = mass of the lower mass
+    % - r  = ratio of ms/mu. This is what is actually optimized - it seems
+    %        that this removes an unconstrained degree of freedom in the
+    %        parameters.
+    % - bs = coefficient of viscous friction between the upper mass and the guiding rail connected to the base
+    % - bu = coefficient of viscous friction between the lower mass and the guiding rail connected to the base
+    % - fs = coefficient of static friction between the upper mass and the guiding rail connected to the base
+    % - fu = coefficient of static friction between the lower mass and the guiding rail connected to the base
+    % - gf = relative strength of the actuator between lower and upper masses
+
+    % explanation of states/inputs/outputs:
+    % - x(1) = height of the upper mass above the base, minus its idle position
+    % - x(2) = speed of the upper mass wrt. the base
+    % - x(3) = height of the lower mass above the base, minus its idle position
+    % - x(4) = speed of the lower mass wrt. the base
+    % - u(1) = road height above the base, minus its idle position
+    % - u(2) = actuator control signal
+    % - y(1) = x(1)
+    % - y(2) = x(3)
+
+
+    zs = x(1);
+    dzs = x(2);
+    zus = x(3);
+    dzus = x(4);
+    zr = u(1);
+    act = u(2);
+
+    sprung_Fext   = +gf*act -(zs-zus)*ks;
+    unsprung_Fext = -gf*act +(zs-zus)*ks + (zr-zus)*ku;
+    sprung_Ffri   = -dzs *bs - sign(dzs) *fs;
+    unsprung_Ffri = -dzus*bu - sign(dzus)*fu;
+
+    dx = [
+        dzs;
+        (sprung_Fext + sprung_Ffri) / ms;
+        dzus;
+        (unsprung_Fext + unsprung_Ffri) / mu
+    ];
+    y = [zs; zus];
+end

data/students/scripts/ident_linear_model.m

+% Run linear model identification
+
+clear variables;
+
+load_ident_data;
+
+% see the as_model_lin.m file for explanation of parameters
+i_vars = {
+    'ks', 1000,
+    'ku', 2500,
+    'r',  2.45/1.3,
+    'bs', 1,
+    'bu', 1,
+    'gf', 1
+};
+model = idgrey(@as_model_lin, i_vars, 'c');
+model.InputName = {'road', 'actuator'};
+model.InputUnit = {'m', 'N'};
+model.OutputName = {'zs', 'zus'};
+model.OutputUnit = {'m', 'm'};
+
+merged_idd = merge(idd_hurak, idd_faststep, idd_forcesweep2, idd_forcestep);
+
+% sometimes it might be better to preprocess the merged data via spafdr()
+% to convert the time-domain data to a frequency-domain representation
+
+opts = greyestOptions('Display', 'on');
+opts.SearchOptions.MaxIterations = 500;
+opts.Focus = 'simulation';
+model_est = greyest(merged_idd, model, opts);
+model_est