This is the implementation of FCND Estimation project.
-
Capture
config/log/Graph1.txt
(GPS X data) andconfig/log/Graph2.txt
(Accelerometer X data) -
Put results to
config/6_Sensornoise.txt
- In
QuadEstimatorEKF.cpp
, update functionUpdateFromIMU()
as described in 7.1.2 Nonlinear Complementary Filter:
- create a quaternion
qt
from Euler angles usingFromEuler123_RPY
function:
Quaternion<float> qt = Quaternion<float>::FromEuler123_RPY(rollEst, pitchEst, ekfState(6))
- use
IntegrateBodyRate()
function (which is just multiplication of 2 quaternions) to integrate gyro rate - use
ToEulerRPY()
to convert resulting quaternion back to Euler angles Roll, Pitch and Yaw
V3D predictedRPY = qt.IntegrateBodyRate(gyro, dtIMU).ToEulerRPY();
float predictedPitch = predictedRPY.y; // Pitch
float predictedRoll = predictedRPY.x; // Roll
ekfState(6) = predictedRPY.z; // Yaw
- normalize yaw to -pi .. pi
// normalize yaw to -pi .. pi
if (ekfState(6) > F_PI) ekfState(6) -= 2.f*F_PI;
if (ekfState(6) < -F_PI) ekfState(6) += 2.f*F_PI;
- Implement the state prediction step in the
PredictState()
functon as described in 7.2 Transition Model: taking into account that control vectorut
is x,y,z accelerations in body frame + Yaw rate: Note, we don't really have to calculate the rotation matrixRbg
, just rotate the accelerations part of control vectorut
usingRotate_BtoI()
function of the attitude quaternion.
// X,Y,Z
predictedState[0] = predictedState[0] + predictedState[3] * dt;
predictedState[1] = predictedState[1] + predictedState[4] * dt;
predictedState[2] = predictedState[2] + predictedState[5] * dt;
// velocities X, Y, Z
V3F accelXYZ = attitude.Rotate_BtoI(accel);
predictedState[3] = predictedState[3] + accelXYZ.x * dt;
predictedState[4] = predictedState[4] + accelXYZ.y * dt;
predictedState[5] = predictedState[5] - CONST_GRAVITY * dt + accelXYZ.z * dt;
Note, there is no need to update predictedState[6]
as it was updated at the previous step.
- Implement function
GetRbgPrime()
, just fill coefficients as following: - Implement
Predict()
function:
gPrime(0,3) = dt;
gPrime(1,4) = dt;
gPrime(2,5) = dt;
VectorXf ut(3);
ut << accel[0], accel[1], accel[2];
gPrime.block<3,1>(3,6) = RbgPrime * ut * dt;
ekfCov = gPrime * ekfCov * gPrime.transpose() + Q;
- Tune the process parameters in
QuadEstimatorEKF.txt
hPrime(0,6) = 1;
zFromX(0) = ekfState(6);
then normalize the difference between measured and estimated yaw
//Normalize
if (z(0)-ekfState(6) > F_PI) z(0) -= 2.f*F_PI;
if (z(0)-ekfState(6) < -F_PI) z(0) += 2.f*F_PI;
hPrime.block<6,6>(0,0).setIdentity();
zFromX = ekfState.head(6);