IMUによるモーションセンシング #4 〜加速度信号のキャリブレーション〜
MEMS型IMUは出荷された状態では,キャリブレーションは行われていないに等しい.そのまま使用すると恐ろしく誤差があるということに気がついている人はあまりいないかもしれない.こここでは,みなさんも自分でできる加速度センサのキャリブレーション方法について述べていく.加速度センサのキャリブレーションに関する方法は数多く発表されているが,ここでは比較的簡便で正確な方法の一つを示す.
今後,このキャリブレーション方法を実装したセンサをリリースする予定で,ご要望があれば,旧タイプのモーションセンサもキャリブレーションを行うサービスも検討している.
加速度信号に含まれる誤差
そもそも加速度センサには多くの種類があるが,身近なセンサはやはりMEMS型のIMUに含まれる加速度センサであろう.MEMSセンサについては
などを参照されたい.チップに組み込んだセンサであるので超小型であるが,そのためノイズの割合が高く,精度も総合的にひずみゲージ型加速度センサなどと比べて若干低い.ただし日進月歩で,限界はあるだろうが精度も様々な面で改良が進められている.ひずみゲージ型加速度センサは大きさも含めて数十年同じスペックのまま販売し続け,ほぼ変わらないことには色々な意味で驚かされる.
さて,この記事でも対象としているモーションセンサは,MEMS型のセンサを使用しやはり誤差を含んでいる.加速度センサで計測される信号には多くの種類の誤差が含まれるが,ここでは,そのうちいくつか実際に使用していく上で大切なものだけをピックアップして述べていく.
1.非線形性(nonlinearity):
運動している加速度が増加すれば,当然比例した加速度を出力すべきで,これらの関係には直線性があるべきだが,加速度の増加に比例して,少しずつ減少するなど,線形性が保たれないことがある.これはセンサのスペックを見ると記述が書かれていることが多いので,よく見てみるとよいだろう.0.1%とか1%などと表示がある.これはMEMSでもひずみゲージ型でも違いはなさそうだが,ユーザレベルでは検証しようがなく,ある程度直線性が保たれていると信じるしかなく,この記事では補正の対象外である.
2.ノイズ(SN比, signal-noise ratio):
図1は200Gの加速度センサを机の上において静止した状態で計測した信号の例である.x軸だけ重力加速度を含んでいる.信号に対するノイズの大きさが一目瞭然だろう.200Gなどの容量の大きな加速度はノイズの割合が増える.16Gなどになれば,加速度としてみた絶対的なノイズの大きさは相対的に小さくなるが,それでもこのようなノイズはMEMS型センサでは避けられない.もちろん運動を行い,100Gもの加速度信号を出力すれば気にもならないだろうが,このようなノイズを含むことは念頭に置くべきだろう.
一般的にこのようなノイズ特性はSN比(signal-noise ratio)と呼ばれる単位で測られ,信号の分散をノイズの分散で割った比で,単位はdBである.これも,補正は難しく,センサの性能そのものである.
3.感度(sensitivity):
ゲイン(gain),スケールファクター(scale factor)などとも呼ばれる.センサの入力を加速度とし,センサの出力を加速度信号とすると,その比を指す.いわゆる加速度の校正係数(キャリブレーション・パラメータ)である.通常,ひずみゲージ型加速度センサには,個々のセンサに対して,この感度のキャリブレーションを行っており校正係数表が一つ一つのセンサに対して添付されている.
しかし,MEMS型センサは,センサの型番が決まればスペックシートに書いてあるので,それに従い計測された信号を加速度に変換する.この違いを見れば,MEMSセンサの精度がいかにいい加減かということがわかる.
4.オフセット(offset):
静止時に加速度センサの水平方向には,一切加速度は検出されないので,センサの水平成分は静止していれば信号は零となるはずだが,実際には0からの偏差が平均的に発生する.これをオフセットと呼ぶ.またこの誤差をバイアス誤差とも呼ぶ.
図1の破線は信号の平均値を示している.この場合センサのx軸を鉛直方向に向けているので,センサのyとz成分の平均は0となるはずだが,オフセットによるバイアス誤差が発生している.z軸も1G(約9.79 m/s2)だけ発生すべきだが,若干小さい.
オフセットは様々要因で混入するが電気的な問題に依存しやすく,温度でも変化する.スポーツセンシング社のモーションセンサも含めて多くのセンセでは,電源を投入してから約20分ほど経過して,このオフセットの変化が安定してから使用することを推奨する.温度補償のないセンサ(ほとんどのMEMS型センサにはない)には,かなりの変動が発生する.最も変化しやすく,最も精度に影響を与える.また,重力加速度は場所によって異なる.このことも影響するが,この違いは通常のバイオメカニクスの利用ではそれほど神経質になる必要ないだろう.9.8でも9.79 m/s$${{}^2}$$でも,スイング運動のような短時間な運動であれば二階積分を行ってもその違いの影響は少ない.それよりは,オフセットを適切に取り除くことだ.
5.非直交性誤差(non-orthogonal error):
センサの3軸(3つの検出方向)は本来直交すべきだが,MEMSであろうが,どのようなセンサでもその直交性には取り付け誤差が発生する(図2).
ここでは,3つの誤差パラメータを図2のように割り当てている.Z軸には誤差がなく,微小な3つの軸まわりの角度誤差があると仮定する.この誤差が,どの程度,精度に影響を及ぼすかの検討はできていないが,ここでは補正の対象とする.
1〜4をまとめると,図3のように示すことができる.
ノイズや誤差については,IMU Errors and Their Effects も参照するとよいだろう.
キャリブレーションモデル
前述の誤差モデルを踏まえて,キャリブレーションは以下のモデルを検討する(文献1).
ここで$${\bm{a}}$$は計測する加速度信号で,キャリブレーション前の信号で,$${\hat{\bm{a}}}$$はキャリブレーション後の加速度である.
$${\bm{T}}$$は非直交性誤差でこのような行列で表現し,$${\bm{S}}$$は感度(ゲイン)で対角行列で表し,$${\bm{o}}$$はオフセットベクトルである.これらのキャリブレーションパラメータはベクトル$${\bm{\theta}}$$にまとめた.
この式で補正される予定の加速度$${\hat{\bm{a}}}$$の大きさは,センサをどのような姿勢で傾けても大きさが常に1Gであるはずだから,
$$
|\hat{\bm{a}}|^2 = \hat{\bm{a}}^T \hat{\bm{a}} =
\begin{bmatrix} \hat{a}_x & \hat{a}_y & \hat{a}_z \end{bmatrix}
\begin{bmatrix} \hat{a}_x \\ \hat{a}_y \\ \hat{a}_z \end{bmatrix}
= \hat{a}_x^2+\hat{a}_y^2 +\hat{a}_z^2
=1
$$
となるはずである.そこで,いろいろな姿勢でIMUを静止させて計測し,その加速度信号を計測し,各平均値を計測する.そして目的関数
$$
L(\bm{\theta}) = \sum_{i=1}^n (\bm{g}^T \bm{g} -\hat{\bm{a}}_i^T \hat{\bm{a}}_i)^2 = \sum_{i=1}^n (1 -\hat{\bm{a}}_i^T \hat{\bm{a}}_i)^2
$$
を最小化する非線形最小問題として解けば,各パラメータ$${\bm{\theta}}$$を計算できる.ここで,$${i}$$は各姿勢を表す.ここでは30の姿勢を与えて計算することになる.
キャリブレーション方法
センサは電源を投入後の約20分後にキャリブレーションデータの取得を開始した.
正確な水平面でなくて良いので,センサを各x,y,z軸の正負の向きをかえて6姿勢計測し,さらに図4の治具に様々な方向を向くように24姿勢固定し,合計30姿勢の静止した加速度信号を計測した.図4にモーションセンサを固定する治具の例を示した.治具は45度傾けて設置しているが,角度は正確でなくても良い.
センサのZ軸の精度が低いことを考えると,多少z軸を鉛直方向に近い姿勢を増やしてもよいだろう.
図5にキャリブレーション前後の加速度信号のノルム(大きさ)を示したが,キャリブレーションを行う前の大きさは1Gから大きくずれている.1Gが1.8G近くになるのだから,かなりの誤差である.ただし,200Gに対してみれば小さい誤差で,これは致し方ないとも言える.
また,これを見ていると,静止時のオフセットが大きいこともうかがえる.
30姿勢分の,この各xyz軸の加速度の平均値
array([[-0.41110849, 0.46267138, 1.25347828],
[ 0.43639215, 1.2984439 , -0.03676222],
[ 0.47153497, -0.26463761, -0.08583459],
[-0.40883585, 0.40048813, -0.06362025],
[ 0.35639765, 1.24939969, 1.34995804],
[-0.28648554, 0.57271477, 1.38092257],
[-0.23495765, -0.29553266, 0.6699422 ],
[ 1.06672854, -0.25744006, 0.49939472],
[-0.29976642, 1.26680148, 0.54795955],
[ 1.07543632, 0.39547795, 1.3763359 ],
[ 0.40222426, -0.19154498, -0.16494066],
[ 0.34101382, 1.12673915, -0.22250941],
[-0.37772971, -0.16851161, 0.59824566],
[ 1.2372339 , -0.08304045, 0.60397591],
[ 1.0069298 , 1.32238043, 0.66020388],
[ 0.40202616, -0.16485611, 1.38717541],
[ 0.99981618, 0.56644808, -0.25670663],
[ 0.44625796, -0.49818785, 0.60046152],
[-0.43364714, 1.11314785, 0.68226333],
[ 0.38219684, 0.52263983, 1.62641298],
[-0.29419147, 0.45656039, -0.20297686],
[ 0.40712141, 1.18934618, 1.39325812],
[ 1.40245984, 0.57672268, 0.57882571],
[ 1.17624931, 0.54677384, -0.07407806],
[ 1.19097004, 0.51384712, 1.28159228],
[ 1.16995439, 1.15485314, 0.55484705],
[-0.6474664 , 0.47789972, 0.55704857],
[ 0.40455249, -0.24458296, 1.32125122],
[ 0.38942328, 0.51774223, -0.44362555],
[ 0.36300253, 1.52517707, 0.58936809]])
を利用し,目的関数$${L(\bm{\theta})}$$を最小化するように,非線形最小二乗法でキャリブレーションパラメータ$${\bm{\theta}=[\alpha_{xy}~\alpha_{xz}~\alpha_{yz}~s_x~s_y~s_z~o_x~o_y ~ o_z ]^T}$$を求めることができる.
この最適化計算では,適切な初期値$${\bm{\theta}_0}$$を与える必要があり,下記のようにパラメータ$${\bm{\theta}}$$を
array([-0.00593346, -0.00974296, 0.01212713, 0.96939586, 0.98263737,
0.95195915, 0.38859128, 0.50761141, 0.60035874])
得ることができる.ここではPythonの Levenberg-Marquardt法(Pythonのscipy.optimize.leastsq)を使用した.
このパラメータ$${\bm{\theta}}$$で補正した$${\hat{\bm{a}}}$$の各姿勢のノルム$${|\hat{\bm{a}}|}$$を図5のオレンジ色の線に示した.このオレンジ色の線の軸を40倍拡大し,図6に示した.
キャリブレーション後は,おおよそ,2%ほどの誤差に収まっている.これはキャリブレーションに使用したデータでキャリブレーションを行った再構成誤差である.
そこで,このキャリブレーションパラメータを実装した加速度センサで,同様に30姿勢を別途計測し,その加速度のノルムを計算すると,おおよそ3%に誤差が収まる.
おわりに
加速度センサのキャリブレーションに関する方法は,数多く発表されている.ここでは,比較的簡便で正確な方法を示した.
ただし,1Gを基準に校正しているので,おのずと200Gの加速度計に対して補正するには限界がある.キャリブレーションを行っても,誤差はそれなりに含んでいることは留意すべきである.
キャリブレーションしていないセンサからキャリブレーションパラメータを求めて,センサから計測した値を補正することは可能だが,オフセットが容易に変化することを考えると,キャリブレーション済みのセンサを利用する方がもちろん断然楽であるし,精度も若干向上する.
Pyhtonコードはどこかで,別途示す予定であるが,$${L(\theta)}$$の目的関数を$${\hat{\bm{a}}}$$を使って展開して記述し,それをPythonのscipyの最適化計算の関数scipy.optimize.leastsq()の中で使用するだけで,非常に短いコードですむので,練習問題としてご自身で試していただければと思う.
今後,このキャリブレーションを実装したセンサをリリースする予定で,ご要望があれば,旧タイプのご利用の方も,キャリブレーションを行うサービスも検討している.
また,下記の学会でもこのキャリブレーション技術をされに応用し,加速度信号を積分し軌道推定する内容で発表予定(文献2)である.
参考文献
D. Tedaldi, A. Pretto and E. Menegatti, “A Robust and Easy to Implement Method for IMU Calibration without External Equipments”, IEEE International Conference on Robotics & Automation (ICRA), 2014, pp.3042–3049.
太田,徳永,大平,センサ座標系におけるバイアス誤差を考慮したIMUによるスイング運動の軌道推定,日本機械学会 シンポジウム:スポーツ工学・ヒューマンダイナミクス2023(SHD2023)@名城大学,2023.11.10-12
【著作権・転載・免責について】
権利の帰属
本ホームページで提示しているソフトウェアならびにプログラムリストは,スポーツセンシング社の著作物であり,スポーツセンシング社に知的所有権がありますが,自由にご利用いただいて構いません.
本ページに掲載されている記事,ソフトウェア,プログラムなどに関する著作権および工業所有権については,株式会社スポーツセンシングに帰属するものです.非営利目的で行う研究用途に限り,無償での使用を許可します.
転載
本ページの内容の転載については非営利目的に限り,本ページの引用であることを明記したうえで,自由に行えるものとします.
免責
本ページで掲載されている内容は,特定の条件下についての内容である場合があります. ソフトウェアやプログラム等,本ページの内容を参照して研究などを行う場合には,その点を十分に踏まえた上で,自己責任でご利用ください.また,本ページの掲載内容によって生じた一切の損害については,株式会社スポーツセンシングおよび著者はその責を負わないものとします.
【解析・受託開発について】
スポーツセンシングでは,豊富な知見を持つ,研究者や各種エンジニアが研究・開発のお手伝いをしております.研究・開発でお困りの方は,ぜひスポーツセンシングにご相談ください.
【例】
・データ解析の代行
・受託開発
(ハードウェア、組込みソフトウェア、PC/モバイルアプリ)
・測定システム構築に関するコンサルティング など
その他,幅広い分野をカバーしておりますので,まずはお気軽にお問い合わせください.
【データの計測について】
スポーツセンシング社のスタジオで,フォースプレートやモーションキャプチャを利用した計測も行えます.出力されるデータと,ここで示したプログラム(入力データの取り込み関数を少々改変する必要があるが)で,同様な解析を行えますので,まずはお気軽にお問い合わせください.