МЭМС акселерометры, магнитометры и углы ориентации



При разработке устройств, в которых необходимо оценивать углы ориентации, может встать вопрос — какой МЭМС датчик выбрать. Производители датчиков приводят множество параметров, из которых нам требуется получить полезную информацию о качестве датчика.

Чтобы понять, на какие точности углов мы можем рассчитывать, нужно приложить некоторое количество усилий.

TL;DR: Описан небольшой скрипт для Octave/MATLAB, позволяющий оценить ошибки расчёта углов ориентации по измерениям МЭМС акселерометров и магнитометров. На входе скрипта — параметры датчиков из даташитов (и/или погрешности калибровки). Статья может быть полезна тем, кто начинает использовать инерциальные датчики в своих устройствах. Небольшой ликбез по датчикам прилагается. Ссылка на гитхаб тоже.

Сходу примем такие условия:

  • Мы хотим оценить углы ориентации неподвижного устройства.

    Почему?
    Ориентацию подвижного устройства просто по формулам не посчитать, нужно использовать хитрые алгоритмы.
  • Для оценки углов мы будем использовать измерения МЭМС акселерометров и магнитометров.

1. Краткий ликбез


Углы ориентации




Будем понимать под углами ориентации объекта углы Эйлера — крен (roll), тангаж (pitch), рыскание (yaw), связывающие собственную систему координат XYZ объекта и локальную систему координат восток-север-верх (ENU — East North Up). Углы roll, pitch, yaw обозначают поворот, который нужно совершить осям XYZ чтобы перейти в оси ENU. Соответственно, нулевые углы означают, что ось X объекта смотрит на восток, ось Y объекта смотрит на север, ось Z — вверх.

Порядок поворота осей — начиная с последнего угла: сначала на yaw (вокруг оси Z), потом на pitch (вокруг оси Y), потом на roll (вокруг оси X).

Акселерометр


Это датчик, измеряющий проекцию кажущегося ускорения на ось чувствительности. Кажущегося — потому что измеряет и силу тяжести тоже, даже в то время как акселерометр неподвижен. Проще всего представить акселерометр как грузик на пружинке, его выдаваемые измерения пропорциональны степени растяжения пружины. Если акселерометр покоится — пружина растянута лишь силой тяжести. Если ускоряется — то будет сумма сил: инерции грузика $\left (F=m\overrightarrow{a}\right )$ и силы тяжести $\left (F_{g}=m\overrightarrow{g}\right )$

Примем следующую модель измерений триады ортогональных (взаимно перпендикулярных) акселерометров:

$a_{XYZ}= m_{a} \cdot A_{XYZ}+b_{a}+n_{a},$

где $a_{XYZ}$ – измеряемое ускорение в ССК (собственной системе координат) XYZ, $m_{a}$ – матрица перекоса осей и масштабных коэффициентов акселерометра, $A_{XYZ}$ – вектор истинного ускорения в ССК XYZ, $b_{a}$ – вектор смещения нуля акселерометра, $n_{a}$ – шум измерений.

Матрица перекоса осей и масштабных коэффициентов выглядит следующим образом:

$m_{a}=\begin{bmatrix}1+m_{a,1,1} & m_{a,1,2} & m_{a,1,3}\\ m_{a,2,1}& 1+m_{a,2,2} & m_{a,2,3}\\ m_{a,3,1} & m_{a,3,2} & 1+m_{a,3,3} \end{bmatrix},$

где элементы, расположенные по главной диагонали ($ 1+m_{a,1,1}, 1+m_{a,2,2}, 1+m_{a,3,3}$) — это масштабные коэффициенты и их погрешности по трём осям акселерометра, а остальные элементы матрицы — перекосы осей акселерометра.

Выбор параметров акселерометра из даташита

Акселерометр MPU-9250


  • Смещение нуля акселерометра — Zero-G Initial Calibration Tolerance ($\pm60mg$ для компонент $X, Y$, $\pm80mg$ для компоненты $Z$) — для расчётов переводим в единицы $g$ домножив на $10^{-3};$
  • Погрешность масштабного коэффициента — Initial Tolerance ($\pm3\%$) — выражается в процентах. Для расчётов надо перевести в разы, домножив на $10^{-2};$
  • Перекосы осей — Cross Axis Sensitivity ($\pm2\%$) — также умножаем на $10^{-2};$
  • Спектральная плотность мощности шума акселерометра — Noise Power Spectral Density $\left (300\frac{\mu g}{\sqrt{Hz}} \right )$ — переводим числитель в $g$ домножая все на $10^{-6};$
  • Полоса пропускания — Low Pass Filter Response $\left (5-260 Hz\right )$ — приведены границы, в пределах которых её можно изменять. Установим максимальную полосу. Все равно ошибки будут определяться не шумами;

Зная спектральную плотность мощности шума и полосу пропускания датчика можно рассчитать СКО шума на выходе датчика:

$\sigma _{noise}=G_{0}\cdot\sqrt{\Pi_{noise}};$


Акселерометр ADIS16488A:


  • Смещение нуля — Bias Repeatability ($\pm 16mg$) — переводим в $g$ домножая на $10^{-3};$
  • Погрешность масштабного коэффициента — (Sensitivity) Repeatability ($\pm 0.5\%$) — переводим из процентов в разы;
  • Перекосы осей — Misalignment Axis to frame ($\pm1^{\circ}$) — в градусах, переводим в разы (радианы, поскольку величины малые);
  • Спектральная плотность мощности шума — Noise Density $\left ( 0.063\frac{mg}{\sqrt{Hz}}rms \right )$ — переводим числитель в $g;$
  • Полоса пропускания — $\left (3dB Bandwidth\right)$ — выберем такой же, как у MPU-9250;


Магнитометр


Датчик, который измеряет проекцию индукции магнитного поля на ось чувствительности. Магнитометру свойственны искажения hard-iron и soft-iron. Hard-iron искажение — это аддитивный эффект, когда к измеряемому полю добавляется постоянная составляющая. Причиной может быть, например, действие постоянного магнита или собственное смещение нуля датчика. Искажение soft-iron — мультипликативный эффект, отражающий изменение направления и/или ослабление вектора магнитной индукции. Этот эффект может быть вызван наличием металлического предмета в непосредственной близости от магнитометра или же собственными искажениями датчика — погрешностью масштабного коэффициента или перекосом его оси чувствительности.
Примем модель измерений триады магнитометров:

$m_{XYZ}=S_{m}\cdot M_{XYZ}+b_{m}+n_{m},$

где $m_{XYZ}$ – измерения магнитометра в ССК XYZ, $S_{m}$ – диагональная матрица перекоса осей и масштабных коэффициентов (которая описывает эффект soft–iron), $M_{XYZ}$ – вектор истинной магнитной индукции в ССК, $b_{m}$ – смещение нулей магнитометра (описывает эффект hard–iron), $n_{m}$ – шум измерений.
Матрица перекоса осей и масштабных коэффициентов магнитометра:

$S_{m}=\begin{bmatrix}1+S_{m,1,1} & S_{m,1,2} &S_{m,1,3}\\ S_{m,2,1}& 1+S_{m,2,2} & S_{m,2,3}\\ S_{m,3,1} & S_{m,3,2}&1+ S_{m,3,3} \end{bmatrix},$

элементы, расположенные на главной диагонали ($S_{m,1,1}, S_{m,2,2}, S_{m,3,3}$) — масштабные коэффициенты и их погрешности по трём осям магнитометра, остальные элементы матрицы — перекосы осей магнитометра. Все элементы матрицы также учитывают эффект soft-iron.

Выбор параметров магнитометра из даташита

Магнитометр MPU-9250


В даташите нужных нам параметров нет, поэтому предположим, что магнитометр откалиброван и возьмем следующие числа:

  • смещение нуля — $(1\mu T);$
  • погрешность масштабных коэффициентов — $(5\%);$
  • перекосы осей — предположим, что они такие же, как у акселерометров — $(\pm2\%);$
  • шум на выходе — $(0.6\mu T);$

Магнитометр ADIS16488A


  • Смещение нуля — Initial Bias Error $\left (\pm15 mgauss = 1.5 \mu T \right )$ — будем считать, что мы его откалибровали до $(0.5 \mu T)$;
  • Погрешность масштабного коэффициента — Initial Sensitivity Tolerance $\left (2\% \right );$
  • Перекосы осей — Misalignment Axis to axis $\left (0.35^{\circ} \right )$ — в градусах, переводим в разы (радианы, так как величина маленькая);
  • Спектральная плотность мощности шума — Noise Density $\left (0.042\frac{mgauss}{\sqrt{Hz}} \right )$ — переводим в $\left (\frac{Tesla}{\sqrt{Hz}} \right );$
  • Полоса пропускания — возьмем для модели значение $260 Hz;$


Расчет углов ориентации


Благодаря наличию на Земле силы тяжести, акселерометры «чувствуют» направление вниз. Их измерения используются для расчета углов крена и тангажа. Формулы для расчёта можно найти тут. Третий — угол рыскания (а в данном случае — магнитного азимута), может быть определен благодаря наличию у Земли магнитного поля. Вектор индукции магнитного поля измеряется магнитометрами и их измерения участвуют в расчете угла рыскания. Нужно отметить, что в расчёте магнитного азимута используются измерения магнитометра, пересчитанные в плоскость. Здесь можно найти формулу для расчёта магнитного азимута.

$roll=atan\left ( \frac{a_{Y}}{a_{Z}} \right ),$

$pitch=atan\left ( \frac{-a_{X}}{\sqrt{a_{Y}^{2}+a_{Z}^{2}}} \right ),$

$yaw=atan2\left ( \frac{m_{E}}{m_{N}} \right ),$

где $atan2$ — функция полного арктангенса, $a_{X}$, $a_{Y}$, $a_{Z}$ — измерения акселерометра по трём осям в ССК, $m_{E}$, $m_{N}$ — измерения магнитометра по осям X', Y' (измерения магнитометров пересчитаны в плоскость).

2. Ошибки оценивания углов ориентации


Описание алгоритма




  • Сформируем массивы случайных углов Эйлера roll, pitch, yaw. Они будут задавать наборы вариантов истинной ориентации объекта в модели.
    Зачем много углов?
    Потому что ошибки зависят от значения углов ориентации, и если мы хотим получить представление об их величине во всем диапазоне изменения — то это самый простой способ.
  • Из случайных углов roll, pitch, yaw формируется матрица преобразования из ССК XYZ в ЛСК ENU:

    $C_{XYZ}^{ENU}=\begin{vmatrix} cy\cdot cp &-cr\cdot sy+sr\cdot cy\cdot sp &sr\cdot sy+cr\cdot cy\cdot sp \\ sy\cdot cp &cr\cdot cy+sr\cdot sy\cdot sp & -sr\cdot cy+cr\cdot sy\cdot sp \\ -sp &sr\cdot cp &cr\cdot cp \end{vmatrix},$

    где $cr=\cos (roll)$, $sr=\sin(roll)$, $cp=\cos(pitch)$, $sp=\sin(pitch)$, $cy=\cos(yaw)$, $sy=\sin(yaw)$.
  • Используя данную матрицу можно получить выражение для истинных ускорений в ССК:

    $A_{XYZ}=\left ( C_{XYZ}^{ENU} \right )^{T}\cdot \begin{vmatrix} 0\\ 0\\ -1\\ \end{vmatrix},$


    $\begin{vmatrix} 0\\ 0\\ -1\\ \end{vmatrix}$ — вектор, определяющий направление гравитационного ускорения, выраженный в единицах g, ${(C_{XYZ}^{ENU} )}^{T}$ — матрица преобразования координат из ЛСК в ССК (обратная матрице преобразования из ССК в ЛСК).
  • Применяем модель измерения акселерометра:

    $a_{XYZ}=\left ( I+m_{a} \right )\cdot A_{XYZ}+b_{a}+n_{a},$

  • По измерениям акселерометра рассчитываются новые углы крена и тангажа (оценки) по формулам:

    $roll'=atan \left ( \frac{a_{Y}}{a_{Z}} \right ),$

    $pitch'=atan\left ( \frac{-a_{X}}{\sqrt{a_{Y}^{2}+a_{Z}^{2}}} \right ).$

  • Также необходимо сформировать матрицу пересчета в «горизонт» из этих углов, для этого воспользуемся функцией rpy2mat:

    $C_{XYZ}^{XYZ'}=rpy2mat\left ( \begin{bmatrix} roll'\\ pitch'\\ 0 \end{bmatrix}^{T} \right ),$

    где углы roll' и pitch' — это углы, рассчитанные по измерениям акселерометра, а третий угол — нулевой.
  • Возьмем вектор истинных магнитных плотностей в ЛСК ENU и пересчитаем его в ССК XYZ:

    $M_{XYZ}= {(C_{XYZ}^{ENU} )}^{T}\cdot M_{ENU}.$

  • Применяем модель измерений магнитометра:

    $m_{XYZ}=S_{m}\cdot M_{XYZ}+b_{m}+n_{m}.$

  • Осталось пересчитать измерения магнитометра из ССК в «горизонт»:

    $m_{XYZ'}=C_{XYZ}^{XYZ'}\cdot m_{XYZ}.$

  • По «горизонтированным» измерениям магнитометра расcчитывается угол магнитного азимута:

    $yaw'=atan2\left ( \frac{m_{Y'}}{m_{X'}} \right ).$

  • Ошибки оценивания углов ориентации рассчитываются как разность между истинными углами roll, pitch, yaw и рассчитанными по измерениям датчиков — roll', pitch', yaw'.

3. Результаты — расчет погрешностей оценки углов


Для двух датчиков, которые мы взяли для примера — ADIS16488A и MPU-9250, получены предельные ошибки оценивания углов ориентации при совместном влиянии погрешностей акселерометра и магнитометра.

В таблице ниже — максимальные значения полученных ошибок:
Угол MPU-9250 ADIS16488A
Крен

$30^{\circ}$

$8^{\circ}$

Тангаж

$10^{\circ}$

$2^{\circ}$

Магнитный азимут

$30^{\circ}$

$20^{\circ}$


Совместное влияние погрешностей акселерометра и магнитометра на ошибки оценивания углов ориентации:

  • Вот так выглядят ошибки оценивания крена в зависимости от значений крена и тангажа:

  • Ошибки оценивания тангажа в зависимости от значений крена и тангажа:

  • Ошибки оценивания магнитного азимута от углов крена и тангажа:

  • Ошибки оценивания магнитного азимута от углов крена и магнитного азимута:

  • Ошибки оценивания магнитного азимута от углов тангажа и магнитного азимута:


Как можно заметить - величина ошибок растёт с приближением к границе диапазонов измерений углов. Почему?
Это можно понять из рисунка ниже.



Допустим мы поворачиваем ось чувствительности Z ($z1\rightarrow z2$) акселерометра так, чтобы проекция силы тяжести на эту ось стала меньше ($g'\rightarrow g"$). Значение проекции силы тяжести плюс погрешность акселерометра дадут область возможных значений измерения (розовая область). Погрешность оценки угла при этом возрастает ($\Delta _{1}\rightarrow \Delta _{2}$). Таким образом, при уменьшении проекции вектора силы тяжести на ось чувствительности ошибка акселерометра начинает вносить бОльшую ошибку в оценку угла.

Поиграем с входными параметрами чтобы понять чем определяются ошибки
Влияние только погрешностей акселерометра (магнитометр считаем идеальным) на ошибки оценивания углов ориентации:

  • Влияние погрешностей акселерометра на ошибки оценивания крена от крена и тангажа
  • Влияние погрешностей акселерометра на ошибки оценивания тангажа от крена и тангажа
  • Влияние погрешностей акселерометра на ошибки оценивания магнитного азимута от крена и тангажа

  • Влияние погрешностей акселерометра на ошибки оценивания магнитного азимута от крена и магнитного азимута
  • Влияние погрешностей акселерометра на ошибки оценивания магнитного азимута от тангажа и магнитного азимута


Влияние погрешностей только магнитометра (акселерометр считаем идеальным) на ошибки оценивания углов ориентации:

  • Влияние погрешностей магнитометра на ошибки оценивания крена от крена и тангажа
  • Влияние погрешностей магнитометра на ошибки оценивания тангажа от крена и тангажа
  • Влияние погрешностей магнитометра на ошибки оценивания магнитного азимута от крена и тангажа
  • Влияние погрешностей магнитометра на ошибки оценивания магнитного азимута от крена и магнитного азимута
  • Влияние погрешностей магнитометра на ошибки оценивания магнитного азимута от тангажа магнитного азимута


Итог


  • Разработана модель, позволяющая оценивать погрешности расчёта углов ориентации. Углы ориентации рассчитаны по измерениям акселерометров и магнитометров. В расчете погрешности углов учитываются модели погрешностей этих датчиков.
  • Ошибки оценивания углов ориентации являются функциями одновременно всех углов ориентации. При этом максимальные значения ошибок наблюдаются на границах диапазонов измерений углов.

Спойлер - дисклеймер:
  • Не рекламируем и не рекомендуем покупать эти датчики (они уже довольно старые).
  • Не учитываем влияние нелинейности измерений, влияние вибрации, нестабильности и т.д. — используем лишь первое приближение к модели ошибок датчиков.

Литература



Авторы

  • Илья Нагин — Stormgazer
  • Дарья Малафеева

Московский Энергетический Институт, Кафедра Радиотехнических Систем
AdBlock похитил этот баннер, но баннеры не зубы — отрастут

Подробнее
Реклама

Комментарии 23

    0
    Круто
      0
      Полезно, но ещё бы гироскоп добавить.
        0
        но ещё бы гироскоп добавить.

        Для неподвижного обьекта гироскоп будет показывать omx=0, omy=0, omz=0. Какой смысл его добавлять?

          0
          Для непожвижного — да. Но на картинке — самолет.
            0

            и? в статье же написано в самом начале:


            Сходу примем такие условия:
            Мы хотим оценить углы ориентации неподвижного устройства.
            ....

              0
              Тогда снимается. Пропустил.
            +4
            Это не так, Земля вращается, поэтому если у Вас есть гироскоп, то Вы можете определить направление на истинный север(а не магнитное склонение) и широту на которой Вы находитесь
              +1
              Полагаю вы имеете в виду нормальный такой гирокомпас весом так это под 20-50 кг. Современные т.н. «гироскопы» не имеют к нему никакого отношения. В последнее время гироскопами принято обзывать акселерометры углового вращения. Ну это то г.но, которое встраивается в смартфоны.
                +1
                Ну не 20-50. Мы, например, делаем гироскоп с электростатическим подвесом ротора (кардана нет). Банка с ротором весит 1 кг. Ну и ещё 1 кг платы питания, обработки информации и подвеса ротора. Есть, конечно, модель сильно покрупнее и уже с карданом (как датчик ноля). Там да, вес приличный.
                0

                ок. Давайте прикинем угловую скорость вращения: 360/24/3600= 1/240= 0.0041(6) гр/с и откроем дашит на хотя бы тот же mpu9250: 131 lsb/гр./с или один мл. бит соотвествует 1/131 гр./с=0.00763 гр./с.

                  0
                  Вы правы только в том, что тема данной статьи МЭМС и ему подобные по массе и размерам датчики, и действительно сейчас нет столь малых и дешевых гирокомпасов.
                  А так, гирокомпасы высокой точности есть, но их масса начинается от 2 а то 3 кг и совершенно несопоставимая цена, выше даже не на порядок а скорее на 2 а то и 3 порядка
                +1
                Ноль будет постоянно ползти.
                К сожалению, для высокоточной техники, у которой нельзя регулярно делать коррекцию, лазерные гироскопы и ВОГ незаменимы.
                  +1
                  лазерные гироскопы и ВОГ незаменимы.


                  Они довольно прилично уплывают. А вот электростатический гироскоп в кардане с обкруткой как датчик ноля (чтобы электроды всегда были одинаково расположены относительно ротора — это сильно улучшает уходы), вот это другое дело.
                    0
                    Есть еще и с магнитным подвесом, правда они тяжелее конечно но зато кардан не нужен
                      0
                      Тут вот в чём дело: если у вас нет кардана, то списывание информации возможно либо оптически, либо ёмкостным методом (по биениям ротора). Эти варианты приводят к искажению формы ротора (и нанесённый рисунок и внесённый дисбаланс портят форму ротора). Есть ещё сквид магнитометр для списывания, но там уже отдельное веселье. Если кардан есть, то положение кардана можно списать да хоть высокоточным энкодером, а форма ротора относительно оси кинетического момента может быть очень точной. Кроме того, одинаковость условий для ротора, бесконтактно подвешенного в кардане, очень способствует уменьшению уходов. Прибор, о котором я говорю, имеет уход лучше чем минус пятая градуса в час.
                        0
                        Вы правы для электростатического ротора, магнитный работает немного по другому(Вообще говоря он хуже почти во всем кроме ухода и снятия информации) — у магнитного ротора есть полюса и поэтому его положение определяется очень легко, второе — он очень жесткий и массивный поэтому его форма не искажается и держит уход он лучше, но это же делает его хуже из-за того что он не выносит сильных изменений ускорения. очень долго запускается и долго входит в меридиан да и борьба с прецессией это даже не песня а скорее опера, причем трагическая. Реально его можно использовать только на больших подлодках да и то в погруженном состоянии когда нет волнения
                          0
                          Вы правы для электростатического ротора,


                          Для магнитного так же используются эти методы (См. Малеев, «Новые типы гироскопов».)

                          Вообще говоря он хуже почти во всем кроме ухода и снятия информации


                          И по этим параметрам он тоже хуже (за исключением криогенного подвеса).

                          у магнитного ротора есть полюса и поэтому его положение определяется очень легко


                          Это смотря с какой точностью. Обычный индукционный метод считывания создаёт дополнительные уводящие моменты. Сквид магнетометр же весьма непростая штука.
              0
              А какие сейчас нужно покупать датчики? Например для динамической ориентации объекта?
                0
                Обычно требования к датчикам определяются проектом (задачи, время работы, точность, цена итд). А уже потом выбираются датчики.
                  0
                  Да, это я бодро спросил, не поставив граничные условия.
                  Если говорить, например, о датчиках стоимостью до 10 USD основной задачей которых является курсовая устойчивость в условиях сильной вибрации, сможете что-то подсказать?
                    0
                    для 10 USD? да любой берите, подороже, они все примерно одинаковые.
                    0

                    А что на счёт БПЛА? Есть возможность коррекции по ГНСС. Цена не определяющий фактор, главное: серийное производство и свободная продажа. Что бы в госдеп не отчитываться.

                      0
                      А что на счёт БПЛА?

                      Все то же самое. Есть условия применения, есть ТЗ они и определяют требования к датчикам ИНС.


                      Есть возможность коррекции по ГНСС.

                      Смотря где, попробуйте откорректироваться по гпс в помещении, или в подвале.


                      Цена не определяющий фактор, главное: серийное производство и свободная продажа.

                      Да ладно? Тоже смотря какой БПЛА, если это кусок пенопласта на него что то точное==дорогое ставить смысла нет. Если это серьезный агрегат, то частному лицу фиг такие датчики продадут да и не каждому предприятию продадут, и ценник у них совсем не гуманный.

                Только полноправные пользователи могут оставлять комментарии. Войдите, пожалуйста.

                Самое читаемое