-
-
Notifications
You must be signed in to change notification settings - Fork 2.9k
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Lower time required for attitude recovery after very hard (involuntary) landings #12937
base: master
Are you sure you want to change the base?
Conversation
Do you want to test this code? You can flash it directly from Betaflight Configurator:
WARNING: It may be unstable. Use only for testing! |
AUTOMERGE: (FAIL)
|
@tbolin please rebase on master. Please explain the tent function :) |
I'll rebase tomorrow. |
1bdecea
to
5b6372e
Compare
Rebased on master and fixed an issue with building the SITL target. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Great work, looks very promising.
Can you please add comments/(references) to concepts used in code, where relevant?
src/main/flight/imu.c
Outdated
@@ -290,10 +291,23 @@ STATIC_UNIT_TESTED void imuMahonyAHRSupdate(float dt, float gx, float gy, float | |||
ay *= recipAccNorm; | |||
az *= recipAccNorm; | |||
|
|||
fpVector3_t acc_upp = {.x = ax, .y = ay, .z = az}; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
fpVector3_t acc_upp = {.x = ax, .y = ay, .z = az}; | |
fpVector3_t acc_bf = {.x = ax, .y = ay, .z = az}; |
Not sure what upp
means. And better use _ef
and _bf
suffix, it makes the code slightly more readable
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'll change the suffixes. What does ef and bf stand for here?
And when I think about it it should probably be AccEF
and AccBF
to comply with Betaflight's coding standards.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
BodyFrame and EarthFrame. Rest of code uses underscore suffix, but using camel case ma be good idea. (personally, I slightly prefer underscore version in this case (reference frame is separated from variable name)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
upp should be up so acc_up
means the direction of up as measured by the accelerometer. In the same way est_up(p)
is the previous estimate of the upwards direction.
Unless I'm missing something both vectors are in the body frame here, and should probably have a bf suffix in addition to the direction.
I'll keep the snake case suffixes for now.
src/main/flight/imu.c
Outdated
@@ -290,10 +291,23 @@ STATIC_UNIT_TESTED void imuMahonyAHRSupdate(float dt, float gx, float gy, float | |||
ay *= recipAccNorm; | |||
az *= recipAccNorm; | |||
|
|||
fpVector3_t acc_upp = {.x = ax, .y = ay, .z = az}; | |||
fpVector3_t est_upp = {.x = rMat[2][0], .y = rMat[2][1], .z = rMat[2][2]}; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is up vector, but acceleration is pointing down.
The code needs cleanup, it was hacked together so that all signs ended up correctly in the end, but was wrong during calculation.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
If the accelerometer is at rest the measured acceleration will point up, in the opposite direction of gravitation.
I agree with the statement about the code being hacked together and needing a bit of cleanup though.
@@ -290,10 +291,23 @@ STATIC_UNIT_TESTED void imuMahonyAHRSupdate(float dt, float gx, float gy, float | |||
ay *= recipAccNorm; | |||
az *= recipAccNorm; | |||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
(You can rewrite code above to use new vector functions .. )
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm painfully aware. I've tried to refrain from vectorizing code that I haven't touched for other reasons. Otherwise that would be all I would be doing. I will definitely have another look at it later when the dust have cleared a bit.
src/main/flight/imu.c
Outdated
// To avoid the gain decreasing for angles > 90 degrees: | ||
// set magnitude of error vector to 1 + |cos| of angle between estimated | ||
// and measured downwards vectors if the absolute angle of the error is > 90 degrees | ||
acc_err = (dot > 0) ? acc_err : *vectorScale(&acc_err, vectorNormalize(&acc_err, &acc_err), 1.0f - dot); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
acc_err = (dot > 0) ? acc_err : *vectorScale(&acc_err, vectorNormalize(&acc_err, &acc_err), 1.0f - dot); | |
acc_err = (dot > 0) ? acc_err : *vectorNormalize(&acc_err, &acc_err); |
IMO keeping error vector constant (=1) after 90 deg is enough and may prevent some obscure stability problem.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Possibly, keeping it at 1 above 90 degrees would also simplify the code. However I like having the error increase if the angle gets larger, and I haven't seen any stability issues despite the current code sometimes using very large gains. The magnitude of acc_err is still bounded to 2 instead of 1, so it's not like it can run away to infinity.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'd prefer simpler version, but it's your call
src/main/flight/imu.c
Outdated
@@ -366,80 +380,36 @@ STATIC_UNIT_TESTED void imuUpdateEulerAngles(void) | |||
} | |||
} | |||
|
|||
static bool imuIsAccelerometerHealthy(float *accAverage) | |||
static float imuAccTrust(const float* accAverage, const float* gyroAverage, const float gRecip) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
gRecip
is really poor name. acc_1gRecip
?
Also, it may be better to convert ACC to g's outside of this function (and eventually move it to acc calibration code)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Agreed on both accounts. I'll rename the variable, and take a look at moving the scaling outside.
@@ -576,7 +564,7 @@ static int calculateThrottleAngleCorrection(void) | |||
|
|||
void imuUpdateAttitude(timeUs_t currentTimeUs) | |||
{ | |||
if (sensors(SENSOR_ACC) && acc.isAccelUpdatedAtLeastOnce) { | |||
if (sensors(SENSOR_ACC) && acc.isAccelUpdatedAtLeastOnce && gyroIsCalibrationComplete()) { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Gyro calibration takes considerable time and offset is not that big. It is IMO possible to converge IMU during this time
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I thought so too. However I had problems where the estimated attitude would end up way off if the quad was moved when connecting the power source and for perhaps 5 seconds after. It's almost as if there are some huge jump in the gyro data during or after calibration.
Adding this check removed the issue.
Might be worth looking into further, but currently this check is needed.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ok. It may also be possible to handle it similarly to gyro overflow case - set assumed gyro noise higher, so that attitude estimate is low weight and filter converges quickly.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The gyroADC and gyroADCf values won't be populated unless gyroCalibrationCompleted
is True, see
betaflight/src/main/sensors/gyro.c
Line 417 in d5d3ee8
if (isGyroSensorCalibrationComplete(&gyro.gyroSensor1)) { |
The estimate covariance is already initiated to the maximum value at start, and converges quickly once it's getting data. However, the covariance should probably reset each time
gyroCalibrationCompleted
is false, to avoid weirdness if another calibration is initiated somehow.
src/main/flight/imu.c
Outdated
@@ -489,10 +459,19 @@ static void imuComputeQuaternionFromRPY(quaternionProducts *quatProd, int16_t in | |||
} | |||
#endif | |||
|
|||
static float imuUpdateDeviation(const float deviation, const float imuDt, const float durationSaturated) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
static float imuUpdateDeviation(const float deviation, const float imuDt, const float durationSaturated) | |
static float imuUpdateDeviationGyro(const float deviation, const float imuDt, const float durationSaturated) |
@tbolin: Alternate approach may be to use principles from Kalman filter (yes, I know). Estimate (or more likely guess) noise value for accelerometer (expected noise in down vector) / gyro (saturated/calibrating/normal) and then keep track of attitude estimate noise, as single value. The math will be very similar to this PR, but with cleaner mathematical background. |
I've worked a little bit with Kalman filters before, and I will take a look at if anything in this PR can be brought in line with that. However, I think the biggest problem will be how to translate the various noises into something that can be used in the Mahoney filter to make it converge in at a rate appropriate for how fast the estimate noise converges. This is also a core problem with the current PR. There really isn't a well reasoned link between the convergence rates of the attitude estimate and the attitude deviation estimate. |
Mahony is in fact complementary filter (or close to one. It is described/implemented as PI regulator, but math behind it is equivalent). So Kp is (almost) equivalent to complementary filter gain. And if I understand Kalman correctly, it is calculating optimum gains based on measurement noise and state noise (and updating state noise). So this should provide gain quite easily. |
In the complementary filters I've seen the gain has been in the range [0, 1], and used for a linear interpolation between the state given by the two sensors e. g. |
You can rewrite it as (with Then, in Mahony,
Quaternion (rotation in general) and nonlinear, but most of Mahony (at least in BF) is using first-order approximations (or simply operations that point in correct direction), so setting Nice thing about Kalman [based derivation] is that tuning will work with noise magnitudes. The have some physical meaning, there are theories how to estimate them (blackbox logs are available). Actual gain calculated will be optimal - at least as far as Kalman theory assumptions are satisfies, which probably means 'not much'. |
04b098f
to
c459cb7
Compare
I have pushed a new version with a much more Kalman like dynamic gain. There is some trickery to avoid doing more than one integration, causing the order of updates to be a bit weird compared to a standard Kalman filter. I also added a check that prevents updates of the attitude if the time delta is greater than 0.1 sec, and that solved an issue with sudden jumps at start up. The roll/pitch estimate covariance should be reset if this happens, but I haven't implemented that yet. |
@tbolin : Great! I need to study the math, but code looks fine. |
I think I've found the source of the gyro wobble drift: #11894 I confirmed the drift when integrating from the BMI270 as in the issue. The scaling adjustment required also seems to be different for different axis. |
Added settings to adjust the gyro rate scaling. The settings are in the error per full rotation around that axis in decidegrees. I turned off the acc and did the wobble test with and without adjustments and the drift was was consistently less than half with the adjustments compared to without. I also made the gyro covariance scale with high gyro rates (still experimental, and should probably be turned down now), and fixed a bug that caused the error vector to always be normalized. Edit: In logs I have seen short bursts of accelerometer readings that are considered valid (the magnitude of the acc vector and gyro rates are within limits), but are off by more than 45 degrees. These readings are likely due to the quad accelerating in just the right way to spoof a valid gravity vector. A standard Kalman filter isn't great at handling outliers. Even if a sensor reading is extremely unlikely given the filters current estimate covariance, the reading can still have a large influence on the estimate. There are ways to handle outliers by e.g. ignoring readings beyond a certain Mahalanobis distance but that's another can of worms if the filter settles on an erroneous value. I think part of the problem in the current implementation is that the adjustments is proportional to the error. |
1d2ab2d
to
6ea3bd0
Compare
Did some flight testing with calibrated gyro rate gains, and I think that helped a lot with the drift during flight. I also had some hard bumps that saturated the gyro and the attitude seems to have recovered quickly, and to a correct value. |
src/main/flight/imu.c
Outdated
|
||
#ifdef USE_MAG | ||
// Use measured magnetic field vector | ||
fpVector3_t mag_bf = {{mag.magADC[X], mag.magADC[Y], mag.magADC[Z]}}; | ||
float recipMagNorm = vectorNormSquared(&mag_bf); | ||
if (useMag && recipMagNorm > 0.01f) { | ||
if (useMag > 0.0f && recipMagNorm > 0.01f) { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Should probably not have been changed, marking so I can get back to it
Finally got the tests to compile locally. Now the old tests compile and passes as they should. Tests for the new code is still a todo. |
348efc7
to
1bf9283
Compare
The error between the estimated up direction and the upward direction measured by the accelerometer is calculated by a cross product. Both vectors have magnitude 1, so the norm of the error is just the absolute value of sin of the angle between the vectors. For angles with an absolute value of more than 90 degrees the magnitude of the error would decrease, slowing down convergence. This commit introduces a check for those large angles and ensure that the magnitude of the error always increases when the angle increases.
- Added tent function to math.h - Decoupled cog calculation from dcmKp gain in imuMahonyAHRSupdate (cog gain shoul be the same as before, but will not be affected by the new dynamicer acceleroemter gain) - Add IMU_GAIN debug mode - Add function (gyroGetDurationSpentSaturated) to gyro.c that returns the time the gyro has been saturated since last heading update - Add estimated attitude deviation to imu.c - Replace imuIsAccelerometerHealthy function with imuAccTrust that calculates a 1.0 - 0.0 scalar instead of just a yes/no. Also checks for high gyroa rates. - Replace old kp gain scheduling function with one that's based on the trust level and assumed deviation instead of armed/disarmed etc. - Only update attitude if the gyro calibration is completed
- Calc duration saturated from gyro looptime instead of PID looptime - Add check that the duration spent saturated is equal to or smaller than the total duration when updating the acc deviation - Lower limit for acceptable gyro rate norm to 100 degrees (was 500 degrees). - Increas limit for acceptable accelerometer norm to 0.2 g (was 0.1 g). - Implement experimental 4th order gyro integration, turned off by define by default - Change normal gyro integration accumulation rate to 1.0 deg/s (was 5.0) - Use base kpgain multiplied by acc trust to update attitude deviation instead of calculated kp
- Implement more kalman like dynamic gain for roll pitch updates from accelerometer. - Add comments to tent function in math.h and remove center parameter - Add guard to prevent IMU heading from being updated if the delta T is greater than 0.1s. Prevents huge timestep after boot. - Remove experimental 4th order quaternion integration - Minor refactoring
- Use durationSaturated instead of durationSaturated^2 when updating the RP estimate covariance in imuUpdateRPEstimateCovariance. The squared duration was a missguided attempt at modeling the gyro drift when saturated as a static bias. The a bias should probably be modeled by having the standard deviation increase linearly, instead of the covariance as it is now. However, that would require taking the square root of the total estimate variance, adding the saturated duration times the assumed drift and then squaring the result. The rate of drift is a crude guess and the extra code complication and perforrmance required are most likely not worth it. - Reset roll/pitch estimate covariance if delta T is greater than 0.1s
- Noise is set in decidegrees, e.g. 5 will lead to a value of 0.5 degrees
- fix bug causing the normalized accDiff vector to always be used - Make it possible to not use acc in attitude update by setting imu_acc_noise_std = 0 - set debug mode IMU_GAIN 2 to display acc angle error in decidegrees - decreas acceleration and gyro rate norm limits to 0.1g and 50 deg/s (from 0.2 g and 100 deg/s) - refactoring
- Turn some magic config variables in functions into static constants on top of file - Remove unused defines related to configuring previous acc gain calculation - Change IMU_ESTIMATE_COVARIANCE_MAXIMUM to 180 (from 360) - Change IMU_ACC_COVARIANCE_CALC_ACC_NORM_LIMIT to 0.2 (from 0.1) - Change IMU_GYRO_COVARIANCE_CALC_RATE_SCALING to 1/10 (from 1/5)
1bf9283
to
04ebc9e
Compare
Work in progress for several small improvements to speed up the recovery after crashing.
Improvements so far
Avoid decreasing accelerometer gain at large errors
The error between the estimated up direction and the upward direction measured by the accelerometer is calculated by a cross product. Both vectors have magnitude 1, so the norm of the error is just the absolute value of sin of the angle between the vectors. For angles with an absolute value of more than 90 degrees the magnitude of the error would decrease, slowing down convergence. This commit introduces a check for those large angles and ensure that the magnitude of the error always increases when the angle increases.
Similar to #12792 but for the accelerometer.
Reliability based accelerometer gain
Basically adds a variable keeps track of how off the current attitude is likely to be. The gain is then calculated based on that variable.
The reliability is implemented as a Kalman filter for a linear system in a single variable, and the reliability is measured as the covariance of the estimated attitude.
The prediction step integrates the gyro rates, and the update step uses the error from the accelerometer.
To save calculations the update step is effectively delayed one cycle. The update and prediction steps can then be carried out in the same quaternion multiplication (this is essentially the same as what was previously done in the Mahoney filter, only now with a different gain for the accelerometer error). The output from the filter is therefore effectively the predicted value.
Settings
imu_gyro_noise_std
standard deviation of the gyro noise in [deci-degrees] . Default: 5.imu_acc_noise_std
standard deviation of the accelerometer noise in [deci-degrees] . Default: 50.Trust based accelerometer covariance
The covariance of the accelerometer noise is increase when the magnitude of the acceleration vector differs from 1g.
When the difference is greater than 0.2g the acceleration vector is not used to update the roll/pitch estimate.
The reasoning behind using this scaling is pretty obvious (and a simpler version is already in BF): if the acceleration is not close to 1g it's probably due to the quad accelerating and the vector is not pointing in the correct direction.
The covariance is also increased at high gyro rates. At rates over 50 deg/s the accelerometer is not used to update the estimate.
The reasons for including the gyroscope rates are that the accelerometer vector is filtered and will lag behind slightly at high rotational rates, and that the centripetal force will create biased readings unless the accelerometer is mounted exactly at the center of gravity.
TODO: add formulas
Trust based gyro covariance
The roll/pitch estimate seems to drift much faster at high gyro rates than when the quad is not moving.
The covariance of the gyro is therefore increased with the gyro rates.
The current formula is
gyro_covariance = imu_gyro_noise_std^2 + gyro_rate_magnitude * imu_gyro_noise_std^2 / 10
Gyro rate scaling adjustments
The sensitivity of a gyro (especially the BMI270s) can differ slightly from what is printed in the data sheet. Even a small error in the sensitivity will give large errors when integrated. e.g. on of my BMI270s have an error of about 2.5% on the yaw axis, which will result in a 9 degree error in heading if turned 360 degrees.
The scaling errors also seems to be the main driving factor behind the weird drift of the roll/pitch estimate when either wobbling the quad by hand or when wobbling quickly in flight.
The PR adds a
gyro_scaling_adjustment
setting that can be used to adjust the sensitivity scaling.gyro_scaling_adjustment
is a 3 value vector, where the values represent the drift in roll, pitch, and yaw resulting from turning the model 360 degrees counter clockwise around the respective axis.TODO: Add detailed calibration guide.
The drunk quad syndrome
Sometimes a quad will no longer be able to fly straight in angle mode after a crash, making it appear drunk.
The exact reason for this is not entirely known.
One factor could be if the gyro is saturated during the crash. An incorrect value will then be integrated and the attitude will then have to be corrected by the accelerometer.
This hypothesis can be tested by quickly flicking an FC (or a whole quad) and checking for drift in the attitude.
However, in the current code the accelerometer gain is constant when the quad is armed, but uses much higher gains when disarmed, so problems that appear during flight can be hard to observe in the configurator.
By commenting out the gain scheduling I could however make the quad model in the configurator diverge when turning the FC quick enough, and due to the low accelerometer gain it takes a long time to get back to the correct value when the FC is placed back on the table.
I also found that I could make the attitude diverge by just wobbling the quad around for a few seconds at far bellow the saturation rate, and I'm not sure why this is.
The attitude should drift a little bit due to integration of the gyro data.
I would expect about 1 degree/sec at worst just from the static offset from the gyro measurements, but considering how steady the heading is (and that bf does a short calibration of the gyro each power up) the drift should be less. By wobbling I could make the attitude deviate by over 20 degrees in a few seconds.
I increased the attitude update rate to 1kHz (from 100Hz) to see if there was something with the integration, but I couldn't see any notable difference.
I also experimented with changing the filtering for the accelerometer and the gyro. Removing all the filtering from the gyro (including the hard-coded filter for data sent to the attitude calculation), or setting the accelerometer filter to any extreme didn't make a difference as far as I could see. Heavy filtering of the gyro (2 pt3 at 10Hz) made the deviation worse, especially the heading deviated much more than before, but not by much.
My best guess is that the wobble deviation come from accelerometer readings that are included despite being biased.