M6-T6 — An Inertial Drift Budget

Learning objectives

After this lesson, you should be able to:

• Write the position error from an accelerometer bias, a gyro bias, angle random walk and velocity random walk as functions of time.
• Check each formula against a simulation.
• Combine independent terms by root sum square and identify the dominant one.
• Find the longest mission a grade can fly unaided within an accuracy and choose a grade.

The intuition

UAV M2-T8: Dead reckoning and inertial drift explains qualitatively why drift grows and sketches its dependence on time. This lesson makes it quantitative: four formulas, each with its power of time, each verified, and a way to read the dominant one from a plot. The result is the tool for sizing a sensor from a requirement.

Four terms, four powers of time

Each sensor error enters the position through integrations, and the number of integrations fixes the power of time. A constant accelerometer bias is integrated twice: error $\tfrac12 b t^2$. White noise on the accelerometer (velocity random walk) is integrated once to a velocity random walk and again to position: $V t^{3/2}/\sqrt3$. A constant gyro bias is integrated once to a tilt that grows as $\varepsilon t$, which turns gravity into a false acceleration, integrated twice more: $\tfrac16 g\varepsilon t^3$. White noise on the gyro (angle random walk) gives a random-walk tilt and the same two integrations: $gN t^{5/2}/\sqrt{20}$. The higher the power, the faster the term takes over as the mission gets longer.

Root sum square and the dominant term

The four terms come from independent sources, so their variances add and the total is the root of the sum of squares. At short times the accelerometer terms dominate, because their powers are lower; at long times the gyro terms win. For the industrial class on Hopper's hop, at 60 s the gyro bias gives 5.1 m, the angle random walk 1.8 m, the accelerometer bias 1.8 m and the velocity random walk 0.22 m.

What it means for Hopper

The error budget of the first module allots 0.60 m (1σ per axis) to navigation. At 60 s the total error is 35.8 m (consumer), 5.72 m (industrial), 0.91 m (tactical) and 0.097 m (navigation). The longest unaided time within 0.60 m is 14, 27, 51 and 140 s. A tactical gyro is therefore close, but short of the 60 s hop at the level the budget needs, so the inertial system has to be aided, which is the subject of module 7.

Position error from each term is a closed form with a distinct power of time. Add the four in root sum square, find the dominant one, and read off how long each grade can fly alone. For Hopper the numbers say tactical grade plus aiding.

Definitions

Accelerometer bias $b$ (m/s²)

Constant offset of the accelerometer, 1σ over the population of start-ups.

Gyro bias $\varepsilon$ (rad/s)

Constant offset of the gyro, 1σ.

Angle random walk $N$ (rad/√s)

White noise on the gyro output.

Velocity random walk $V$ (m/s/√s)

White noise on the accelerometer output.

Root sum square

$\sigma = \sqrt{\sigma_1^2 + \sigma_2^2 + \cdots}$: the combination of independent errors.

Derivation: the four formulas

All four follow from the flat-earth, short-time error model $\dot{x} = v,\ \dot{v} = -g\psi + a_{\mathrm{err}},\ \dot\psi = \omega_{\mathrm{err}}$, the version of the Schuler equations of the last lesson with the coupling dropped.

1. $$\delta x = \iint b\,dt\,dt = \tfrac12\,b\,t^2$$

   Constant accelerometer bias, twice integrated.

2. $$\psi = \varepsilon t, \quad \delta x = \iint -g\,\varepsilon\,t\,dt\,dt = -\tfrac16\,g\,\varepsilon\,t^3$$

   Constant gyro bias gives a linearly growing tilt, which gives a false acceleration, integrated twice. The magnitude is $g\varepsilon t^3/6$.

3. $$\psi(t) = N\,W(t), \quad \delta x = -g N\int_0^t \dfrac{(t - s)^2}{2}\,dW(s)$$

   White noise on the gyro: the tilt is a Wiener process of intensity N. Twice integrating a Wiener process is a single integral of $(t-s)^2/2$ against its increments.

4. $$\mathrm{Var}(\delta x) = g^2N^2\int_0^t \dfrac{(t-s)^4}{4}\,ds = \dfrac{g^2N^2t^5}{20}$$

   The variance of a stochastic integral is the integral of the squared kernel (Itô isometry). So $\sigma_x = gN\,t^{5/2}/\sqrt{20}$.

5. $$\mathrm{Var}(\delta x) = V^2\int_0^t (t - s)^2\,ds = \dfrac{V^2t^3}{3}$$

   White noise on the accelerometer: the position is $\int_0^t (t-s)\,dW$ with the same isometry, giving $\sigma_x = V\,t^{3/2}/\sqrt3$.

The power of time of each term counts the integrations: bias on the accelerometer gives t², on the gyro t³; the white-noise terms have half a power less (t³ᐟ² and t⁵ᐟ²). Independent terms add as squares.

A worked example

Hopper's hop with the four grade classes of the first lesson of this module and the 0.60 m navigation budget.

1. Check the formulas. For the industrial class at 60 s a simulation gives 1.766 m (accelerometer bias; formula 1.766), 5.143 m (gyro bias; 5.137) and, over 400 trials, 1.703 m (angle random walk; formula 1.779) and 0.2206 m (velocity random walk; 0.2236).
2. Combine. The root sum square at 60 s is 35.8 m for the consumer grade, 5.72 m for the industrial grade, 0.912 m for the tactical grade and 0.097 m for the navigation grade. At 20 s: 1.64, 0.30, 0.071 and 0.010 m. At 120 s: 278, 42.9, 5.5 and 0.43 m.
3. The dominant term. At 60 s it is the gyro bias for the consumer and industrial grades, the angle random walk for the tactical grade and the accelerometer bias for the navigation grade.
4. Unaided time. Within 0.60 m the consumer grade lasts 14 s, the industrial 27 s, the tactical 51 s and the navigation grade 140 s.
5. Choice. Only the navigation grade holds the budget for the whole hop. A tactical gyro is within a factor of 1.5, and in the module that follows a position fix at modest accuracy reduces the error far below the budget. The consumer and industrial grades are 60 and 10 times over at 60 s and need aiding as well.

import math

def lcg(seed):
    s = seed
    while True:
        s = (1103515245 * s + 12345) % 2147483648
        yield (s / 2147483648.0 - 0.5) * 3.4641016

G0 = 9.81
D2R = math.pi / 180
rt = math.sqrt

def terms(t, gb_dph, arw, ab_mg, vrw):
    """One-sigma horizontal position error (m) at time t (s) from each sensor error term, flat-earth, short-time model."""
    eps = gb_dph * D2R / 3600                  # rad/s
    N = arw * D2R / 60                         # rad/sqrt(s)
    b = ab_mg * 1e-3 * G0                      # m/s2
    V = vrw / 60                               # m/s/sqrt(s)
    return {'accel bias': 0.5 * b * t * t,
            'gyro bias': G0 * eps * t ** 3 / 6,
            'angle random walk': G0 * N * t ** 2.5 / rt(20),
            'velocity random walk': V * t ** 1.5 / rt(3)}

# Monte Carlo check of the four formulas for one grade at 60 s
rnd = lcg(11)
dt = 0.05
t_end = 60.0
n = int(t_end / dt)
gb, arw, ab, vrw = 3.0, 0.10, 0.10, 0.05
eps = gb * D2R / 3600; N = arw * D2R / 60; b = ab * 1e-3 * G0; V = vrw / 60
def sim(kind):
    th = v = x = 0.0
    for k in range(n):
        t = k * dt
        if kind == 'gyro bias': th += eps * dt
        elif kind == 'angle random walk': th += N * math.sqrt(dt) * next(rnd)
        acc = -G0 * th
        if kind == 'accel bias': acc += b
        elif kind == 'velocity random walk': acc += V / math.sqrt(dt) * next(rnd)
        x += v * dt + 0.5 * acc * dt * dt
        v += acc * dt
    return x
trials = 400
print(f'Monte Carlo against the formulas at {t_end:.0f} s (industrial grade, {trials} trials for the random terms):')
th_ = terms(t_end, gb, arw, ab, vrw)
for kind in ('accel bias', 'gyro bias'):
    print(f'  {kind:22s} simulated {abs(sim(kind)):7.3f} m, formula {th_[kind]:7.3f} m')
for kind in ('angle random walk', 'velocity random walk'):
    xs = [sim(kind) for _ in range(trials)]
    sd = rt(sum(x * x for x in xs) / trials)
    print(f'  {kind:22s} simulated {sd:7.4f} m, formula {th_[kind]:7.4f} m')

grades = [('consumer MEMS', 20.0, 0.30, 0.5, 0.10), ('industrial MEMS', 3.0, 0.10, 0.10, 0.05),
          ('tactical', 0.3, 0.03, 0.03, 0.01), ('navigation', 0.01, 0.002, 0.005, 0.002)]
print('Hopper hop, one-sigma position error per axis, root-sum-square of the four terms (budget 0.60 m):')
print('  grade            t = 20 s    t = 60 s    t = 120 s   dominant term at 60 s')
for name, gb, arw, ab, vrw in grades:
    row = []
    for t in (20.0, 60.0, 120.0):
        d = terms(t, gb, arw, ab, vrw)
        row.append(rt(sum(v * v for v in d.values())))
    d60 = terms(60.0, gb, arw, ab, vrw)
    dom = max(d60, key=d60.get)
    print(f'  {name:16s} {row[0]:8.3f} m  {row[1]:8.3f} m  {row[2]:8.2f} m   {dom}')
# longest time each grade meets 0.6 m
print('longest free-inertial time within 0.60 m:')
for name, gb, arw, ab, vrw in grades:
    t = 1.0
    while t < 3000 and rt(sum(v * v for v in terms(t, gb, arw, ab, vrw).values())) < 0.60:
        t += 1.0
    print(f'  {name:16s} {t:6.0f} s')

Monte Carlo against the formulas at 60 s (industrial grade, 400 trials for the random terms):
  accel bias             simulated   1.766 m, formula   1.766 m
  gyro bias              simulated   5.143 m, formula   5.137 m
  angle random walk      simulated  1.7029 m, formula  1.7793 m
  velocity random walk   simulated  0.2206 m, formula  0.2236 m
Hopper hop, one-sigma position error per axis, root-sum-square of the four terms (budget 0.60 m):
  grade            t = 20 s    t = 60 s    t = 120 s   dominant term at 60 s
  consumer MEMS       1.642 m    35.767 m    277.86 m   gyro bias
  industrial MEMS     0.299 m     5.720 m     42.90 m   gyro bias
  tactical            0.071 m     0.912 m      5.52 m   angle random walk
  navigation          0.010 m     0.097 m      0.43 m   accel bias
longest free-inertial time within 0.60 m:
  consumer MEMS        14 s
  industrial MEMS      27 s
  tactical             51 s
  navigation          140 s

In robotics

The budget turns requirements into specifications in both directions.

• From mission to sensor. Given an accuracy and a time, find the grades that meet it. For Hopper, 0.60 m at 60 s selects navigation grade for unaided flight, and a tactical grade with aiding.
• From sensor to mission. Given a sensor, the longest time it can go unaided. A tactical IMU is good for tens of seconds at the metre level, which is how guided munitions and short hops are planned.
• Where to spend. The dominant term shows where a better sensor or a better calibration pays. If the gyro bias dominates, improving the accelerometer buys nothing. For the tactical grade at 60 s the angle random walk is the largest, so a better bias calibration helps little.
• Allocation inside the error budget. The navigation share of the landing requirement (0.60 m) was an allocation, not a result. If a tactical IMU is all that fits the mass budget, the guidance and control shares have to shrink or the aiding has to be better, and this lesson says by how much.
• Limits. The formulas are short-time, flat-earth and one axis. Beyond about ten minutes the Schuler and vertical-channel effects of the last lesson change the picture, and the approximation of a level vehicle fails during manoeuvres, which add scale-factor and misalignment errors.

The lab M6-L2: Mechanize and Watch It Drift reproduces these growth laws in a planar mechanization. Module 7 begins with the aiding that makes a tactical grade enough, starting at M7-T1: The Position Fix: Pseudoranges, Least Squares and Dilution of Precision.

Common misconceptions

Key takeaways

1. Position error from each term: accelerometer bias $\tfrac12 bt^2$, gyro bias $\tfrac16 g\varepsilon t^3$, angle random walk $gNt^{5/2}/\sqrt{20}$, velocity random walk $Vt^{3/2}/\sqrt3$.
2. The powers of time count the integrations; independent terms combine by root sum square.
3. Monte Carlo simulation agrees with the formulas (1.766 m and 5.143 m against 1.766 and 5.137; the noise terms within 5%).
4. At 60 s the root sum square is 35.8, 5.72, 0.912 and 0.097 m for the consumer, industrial, tactical and navigation grades.
5. Within the 0.60 m budget the grades last 14, 27, 51 and 140 s unaided.
6. For Hopper the budget calls for a navigation-grade IMU, or a tactical one with the aiding of module 7.

References and further reading

Each source is marked free or paid. Free to read does not mean free to reuse: the course resources page lists every source with its licence.

• P. D. Groves, Principles of GNSS, Inertial, and Multisensor Integrated Navigation Systems, 2nd edition, Artech House, 2013, chapter 5, section on error growth with time. Paid.
• O. J. Woodman, "An introduction to inertial navigation", University of Cambridge Computer Laboratory Technical Report 696, 2007, with the same error growth laws for low-cost sensors. Free to read from the University of Cambridge repository.
• N. Barbour and G. Schmidt, "Inertial sensor technology trends", IEEE Sensors Journal 1(4), 2001, on sensor classes and error sizes. Paid (IEEE).
• MIT 16.322 Stochastic Estimation and Control (OpenCourseWare) , for random walks and the Itô isometry used in the derivation. Free; licence CC BY-NC-SA 4.0.

Knowledge check

Three questions. Write an answer, then reveal the explanation. Your answers stay in your browser.

What’s next

Related content

Nothing below is required to finish this lesson.

• Before
  M6-T5: Initial Alignment: Levelling and Gyrocompassing

  Find roll and pitch from the accelerometers and heading from the Earth rate at rest, and state the heading accuracy a given gyro bias allows.

• Lab
  M6-L2: Mechanize and Watch It Drift

  Write a planar strapdown mechanization, inject a known accelerometer bias and gyro bias, and match the growth of the position error to the error equations.

Other courses on this site treat neighbouring topics from their own angle:

• UAV M2-T7: Inertial Navigation Systems
• UAV M2-T8: Dead reckoning and inertial drift
• UAV M3-T7: IMU noise and bias modelling
• UAV M3-T8: IMU preintegration
• UAV M1-T8: IMU, barometer, magnetometer and GNSS
• Probability M2-T4: Distributions of Real Sensor Noise
• Probability M5-T4: Identifying Noise Parameters from Data
• Linear Algebra M10-T6: Small-Angle Approximations

This course leaves the following to the pages that already teach them:

• What an INS is and why it drifts, without the equations
• IMU preintegration for factor-graph estimators
• Wheel odometry and its motion models