function [estParams] = ferraris_1994(trajData, accel_MPU)
%FERRARIS_1994_FCN Summary of this function goes here
%   Detailed explanation goes here

g = 9.81;

% Get accelerations for each position
av_acc = computeAvAccels(trajData, accel_MPU);

U_apos = av_acc(1:3,:); U_aneg = av_acc(4:6,:);

U_as = U_apos + U_aneg; B_a = 0.5 * U_as;
U_ad = U_apos - U_aneg;

k_as = sqrt((1/(4*g^2)) * diag(U_ad * U_ad')); % estimated scale factors
R_a = (1/(2*g)) * U_ad \ (diag(k_as)); % estimated misalignment matrix.

S_est = diag(k_as) * R_a; b_est = mean(B_a, 1);

estParams.Bias = b_est; % estimated accelerometer biases
estParams.AxesMisalignment = S_est; % estimated accelerometer axis sensitivities

% Return individual calibration components:
estParams.ScaleFactors = diag(S_est);
estParams.MisalignAngles = -rotm2eul(R_a', 'XYZ'); % in rad

end