JP4093861B2 - 電子コンパス及び全オリエンテーション動作に対する大磁気誤差の補償 - Google Patents

Description

本発明は契約番号DAAB07-95-C-M028につき米国政府支援の下になされた。政府は本発明について所定の権利を有する。
本発明は電子コンパスに関する。特に、本発明は磁気誤差について電子コンパスを補償する技術に関する。

電子コンパスは当技術分野において周知である。代表的に、この種の装置は磁束ゲート等の磁気センサを用いて磁界を測定し、地球磁界すなわち地球の磁界に関する方位を決定する。しかしながら、磁針や方向目盛板の付いたコンパス（羅針盤）と同様に、電子コンパスを鉄金属や摂動界の環境下で使用すると、磁力計センサの検出する磁界はかき乱され、地球磁界の測定値及びコンパスの磁針方位角に誤差が生じる。

正確な測定値を得るには、これらの磁気摂動について補償する必要がある。磁針や方向目盛板の付いたコンパスを補償するには棒磁石や軟鉄の塊を用いて摂動界を物理的に中和又は相殺しようとする。これらの磁石や軟鉄塊は、磁気摂動を相殺し、偏差誤差が例えば３〜５°に低減されるように、コンパスのまわりに注意深く配置しなければならない。さらに、このような低減を行った後でも、残留偏差誤差を真方位に対してプロットしてユーザが方位角の値を訂正できるようにする必要がある。

電子コンパスは、磁力計センサから渡されたデータをマイクロプロセッサで処理する構成であり、数値解法を用いて補償を行うことができる。数値解法を使用する補償の一具体例としてレベルコンパスに適用される伝統的な５項補償公式がある。この補償公式では、真方位に対するコンパス測定方位のずれ（偏差）を真方位角θの関数として度（°）で示す。式で示すと、
偏差＝A＋B・sin(θ)＋C・cos(θ)＋D・sin(2θ)＋E・cos(2θ)
ここに、A、B、C、D、Eは較正手続きを通じて値が求められる係数である。この補償技術には一定の限界がある。例えば、上記補償公式は近似にすぎず、偏差が小さくA、B、C、D、Eの値が小さい場合のみ有効である。したがって、この補償技術は、磁石や軟鉄塊を加える等の物理的補償により偏差誤差を低減した後でしか利用されない。さらに、補償公式はレベルコンパスに対してのみ有効であり、例えば傾いた船内のようにコンパスのオリエンテーションが三次元空間内で大きく変動するような環境下には適さない。また緯度が変わると補償精度に影響が出る。この補償技術は航空機の磁束ゲートコンパスにも適用することができるが、係数がチルトオリエンテーションおよび磁気緯度に依存するので普通の垂下式磁束ゲートコンパスに対してのみ有効である。

一般に補償技術はこれまで実施されてきたが、その多くは２軸コンパスに限定され、コンパスが三次元空間で様々なオリエンテーションを取るような場合に十分な成果を出すことができない。従来の補償技術のなかには３軸コンパスに適用されるものもある。しかしながら、この種の従来技術はの対称係数行列に依拠しているため、ある種の誤差について十分な補償を行うことができない。

本発明は上記及びその他の問題に取り組むために、３軸物理モデルを使用して電子コンパスにおける磁界測定値の誤差について数値補償する。このモデルは近似ではなく物理の原理に基づくもので、精度の高い補償を行う。さらに、このモデルは補償の適用及び補償に必要なパラメータ決定を促進する線形代数法を利用する。３軸物理モデルを使用することにより、本発明は、ほぼ水平なオリエンテーションでの使用に限定されず、任意のオリエンテーションに対して電子コンパスを補償する。

一実施形態によれば、本発明は、摂動界の影響下でも正確な方位データが得られるように電子コンパスを補償する電子コンパス方位補償方法に係る。オリエンテーションと方位の組合せ毎に、３軸について磁界強度を表す測定磁界ベクトルH_MEAS、及び３軸について重力場強度を表す測定重力ベクトルG_MEASを取得する。これらの測定ベクトルに基づいてシステム方程式を用いることにより、行列補償係数L_E及びベクトル補償係数H_PEを計算する。計算した行列補償係数L_E及びベクトル補償係数H_PEは後に測定した磁力計出力データを訂正するのに利用し、これにより正確な地球磁界データを得る。

他の実施形態は本発明の方法を実施する電子コンパス装置及びマイクロプロセッサ読み取り可能な媒体に係る。この欄の説明は図示した全ての実施の形態あるいは本発明の全ての具体例を記載したものではない。以下説明する発明の実施の形態及び図面においてこれらの実施の形態をより詳細に説明する。

本発明の上記及びその他の利点は図面を参照してなされる以下の発明の実施の形態において明らかにされる。
本発明は種々の変形及び代替形態が可能である。特定の形態について図面に例示し、以下で説明する。しかしながら、その意図するところは、本発明を記載される特定形態に限定することではない。むしろ、本発明は特許請求の範囲に規定されるように、本発明の範囲内でなされる全ての変形、均等及び代替形態を包含するものである。

本発明は電子コンパスに適用可能である。本発明は、大きな磁気摂動がある不良な磁界環境下で様々なオリエンテーションで動作するコンパスに対して特に有効である。様々な本発明の実施形態に関する説明を通じて本発明の種々の特徴を理解することができる。

本発明はコンパスの取り得る全オリエンテーションにおいて正確な補償及び動作を効果的に行うもので、これによりコンパスが設置されるビークル又は装置は任意のオリエンテーションを取ることが可能になる。「鋼鉄」及び「軟鉄」による誤差がともに補償される。地球磁界（地磁気）より大きな固定磁気摂動であっても補償される。磁界信号による方位の計算に先立って磁界信号に補償を施しているので、従来の技術に比べ高精度が実現される。本発明において使用する補償式は磁界の物理学に基づくものであり、数学的な近似を行う必要がないため、精度はさらに改善される。

他の利点として、同一の磁界補償式及び較正手続きにより磁界によらない様々なセンサ誤差、例えばゼロオフセット、不等目盛係数、軸の不一致、磁力計とオリエンテーションセンサ間の軸の不一致等についても補償することができる。

本発明の具体的な一実施形態によれば、電子コンパスはそれが設置されるビークル又は他の装置に対して固定された３軸電子磁力計に基づいている。３軸加速度計又はオリエンテーションセンサは重力の方向を測定し、コンパスの方位が定義される局所的な水平面がこれから求められる。当業者には明らかなように、重力方向は他の構成、例えば航空機のフライトコントロールセンサ系を使用して検出してもよい。コンピュータのマイクロプロセッサは測定磁界及び測定重力場について取得した数値を使用することにより、コンパスのオリエンテーションに関わらずその軸の磁気方位を計算する。結果として、電子コンパスは任意のオリエンテーションで動作することができ、水平を維持する必要がない。

電子コンパスを較正するために、磁界と重力場の測定はピッチ、ロール及び方位のオリエンテーションの組合せについて様々な、偏りのない分布で行う。次に各オリエンテーションの組合せ毎に、決定すべき、行列係数、ベクトル係数及び測定磁界の関数として地球磁界の推定値を与える方程式を立てる。このようにして一連の測定に基づいてシステム方程式（連立方程式）が作成される。これらの方程式の一部は一次（線形）であり、他は二次方程式である。二次方程式は近似化した一次方程式として再構成される。

次に、一次方程式の部分集合を先ず解いて行列係数及びベクトル係数について予備（初期）推定値を得る。次に全システム方程式を解いて予備の推定値を細分化する。全システム方程式を必要に応じ反復して解き、さらに推定値を細分化して補償精度を改善する。次にこのようにして得られた行列係数及びベクトル係数を記憶、保存し、電子コンパスの通常動作において磁力計の誤差を補償するのに使用する。

本開示の残りでは、先ず本発明に内在する原理について説明する。次に電子コンパスを較正する際に使用する測定値の取得処理について説明する。次に電子コンパスを補償するのに使用する補償式について説明する。最後にシステム方程式を解いて所要の数値係数について所要値を求める処理について説明する。
基本原理
上述したように、測定磁界は動作環境に存在する、例えば鉄金属により生じる誤差の影響を受ける。これらの誤差は加法的な磁気摂動に起因し、これにより方位に偏差誤差（真方位からのズレ）が生じる。これらの摂動のため、測定磁界H_MEASは地球磁界H_EARTHに一致せず、次のような関係になる。

ここにH_INDUCEDは誘導磁界摂動を表し、H_PERMは永久磁界摂動を表す。
摂動が測定磁界に及ぼす影響を図１Ａ及び１Ｂに２軸の非常に簡単なケースについて示す。図１Ａは、レベルコンパスと３６０°のコンパスの振れを想定し、摂動なしの場合の地球磁界の軌跡を示す。図１Ａから明らかなように、軌跡は円であり、方位角は磁北から時計回りに定義される。H_EARTHはコンパスの振れとともに反時計回りに回転し、垂直（鉛直）磁界成分H_Zは定数である。

図１ＢはH_PERM及びH_INDUCEDによる合成影響を示したもので、前者により軌跡の中心がずれ、後者により楕円形の軌跡になる。したがって、軌跡は中心のずれた楕円となり、任意のオリエンテーション（方向）を取り得る。さらにポイントは楕円に沿って回転する。H_PERMがH_EARTHより大きいならばH_MEASは振り切れず、ポイントによっては１８０°の誤差が発生する。さらに、H_INDUCED摂動によりH_MEASは全方位で一定の偏差誤差を含むものとなり、全てのオリエンテーションで不正確になる。

永久磁界摂動H_PERMはコンパスの任意部分及びその設置環境における永久磁化及び直流電流に基づくものである。したがって、コンパスの場合、H_PERMはコンパスとその設置に対しては一定で変化しない。一方、H_EARTHはオリエンテーションに依存して変わるので、H_EARTHに対するH_PERMのオリエンテーションはコンパスのオリエンテーションにより変化する。H_PERMは３軸全てについて磁界のゼロオフセットを発生させ、方位の偏差誤差となる。これらの誤差は方位及び傾きオリエンテーションに応じて変動する。

誘導磁界摂動H_INDUCEDはコンパスが動作する環境にある磁化可能な軟鉄金属に起因する。H_INDUCEDは大きさと方向がともに可変量でH_EARTHに比例する。

ここにKは磁化率に似た行列である。
誘導磁界摂動ベクトルH_INDUCEDは常に３つの線形独立な寄与に分解することができる。式で示すと、

ここにKはH_EARTHの関数としてH_INDUCEDを特徴づける３×３の磁化率行列である。磁気センサにおける３つの磁界成分H₁、H₂、及びH₃に沿う単位ベクトルh₁、h₂、h₃の方向は固定され、かつ互いに直交し、地球磁界H_EARTHの方向に依存しない。H_EARTHのオリエンテーションの変化によりこれらの成分のいずれかが最大又は最小値になると残る２成分はゼロになる。

地球磁界H_EARTH自体の向きは単位ベクトルe₁、e₂、e₃で示される３つの対応する直交方向を有し、各単位ベクトルに沿うH_EARTHの各成分は正又は負方向で最大値となる。これらの最大値はk₁、k₂、k₃で示される別個の値を有する。H_INDUCEDとH_EARTHに係る、２つの直交三組は互いに独立している。三組と最大値は、いずれも特異値分解法又はSVD法と呼ばれる技法により磁化率行列Kを用いて推定され、表現される。

すなわち、Kは次のように３つのSVD行列因子に分解される。

ここにh₁、h₂、h₃はセンサにおける誘導磁界方向に対応する単位列ベクトル、k₁、k₂、k₃は各成分のピーク値、e₁、e₂、e₃は地球磁界H_EARTHに関する直交成分方向に対応する単位行ベクトルである。H_EARTHの左側から作用する、列ベクトルh〉と行ベクトル〈eから得られる３個の行列の和としてKの３項ダイアディック表現が導かれる。

KがH_EARTHに作用するとき、上記方程式は次のようになる。

ここに、

である。
任意の値のKについて、Kを特異値分解するとこれらの表現に関する数値、及びH_INDUCEDとH_EARTHに対するそれぞれの幾何学的寄与が得られる。

非電子コンパスを物理的に補償する伝統技術でよく知られているように、誘導磁化摂動H_INDUCEDに関して独立した９種類の誤差がある。さらに個別の種類毎に特定の強度値がある。この伝統技術では、さらにこれらの各種誤差を、それぞれコンパスのまわりに特定の態様で配置される９種類の透磁性材料の「棒」によって物理的に記述する。これらの仮想的な棒はコンパス設置下における鉄金属の分布と結果的に等価であり、この種の分布を表している。実際には、９種類の棒は磁化率行列Kを構成する９つの要素に対応する。なんらかの特定の性質を磁化率行列Kに帰属させるような、磁気学上の物理的法則は全くない。例えば、磁化率行列Kは、６個の要素値だけで行列を特定することができるような対称行列ではない。

図２は、３軸空間磁力計により測定される、永久磁界及び誘導磁界の双方による磁気摂動が地球磁界の軌跡に対して与える合成影響を示したものである。線形代数学でよく知られるように、３×３の対称行列は、図２に示すように、センサの測定した磁界の大きさに関して三次元の楕円面を定義する。しかしながら、一般的な行列における非対称部は、形状を変えることなく、楕円面に沿う軌跡点を様々に回転させ、角度変位を生じさせる。したがって、９種類の誤差のうち３つは、楕円面が表す磁界強度における摂動とは別に、測定した磁界方向の回転として記述することができる。コンパスはこの回転のみによる誤差によって相当変動するが、他の点に関し誤差は非常によく似たものになる。

図２は様々な固定チルトオリエンテーション及び磁気緯度について軌跡をシミュレーションした各種例を示したものである。チルトオリエンテーションと磁気緯度の組合せ毎に、軌跡は図２の楕円体の一断面に相当する、閉じた平面楕円を形成し、コンパスとその設置に特徴的なものになる。楕円体の中心はH_PERMベクトルの分だけ原点から変位する。個々の軌跡において、＊印はコンパスの方位が北となるときの測定磁界H_MEASを示している。図２は磁気摂動の補償方法を決定する際の複雑さと難しさを表している。コンパスのこれらの特性に関わる１２個の独立定数がある。

個々の種類の誤差は、コンパスの設置環境において任意の度合いの正あるいは負の強度で、独立して存在し得る。H_INDUCEDに起因するコンパス誤差を完全に補償するには、独立した９種類の誤差の補償を正確に行わなければならない。永久磁化摂動H_PERMにより、補償すべき誤差がさらに３種類加わり、全部で１２種類になる。本発明はH_INDUCEDとH_PERMに起因する全種類の誤差を十分かつ正確に補償するものである。

磁化率行列Kは通常、小さい。測定磁界H_MEASにおいて、H_EARTH及び

は常に加法関係にあるので、計算を容易にするために透磁率行列Mを用いることができる。ここにMはKと次の関係にある。

したがって、

誘導磁界摂動H_INDUCEDと永久磁界摂動H_PERMの両方に配慮すると、磁力計の測定する磁界H_MEASは次式で表される。

ここに、Mは上述したように透磁率行列である。行列M及びベクトルH_PERMの値を決定することにより、補償前のコンパスを完全に特徴づけることができる。例えば温度及びヒステリシス条件、負荷シフトや装置誤差等のパラメータの不安定性に関係する因子を除けば、この式はコンパスに関係するほとんど全ての誤差を考慮している。
実施例
図３は本発明の一実施形態に基づいた、電子コンパス補償構成例３００を示すブロック図である。電子コンパス３００はそれが設置されるビークル又は装置に対して固定された３軸磁力計３０２を有する。磁力計３０２は図４に示すように、３軸の内部座標系において磁界を測定する。図４において、電子コンパスの内部座標系はX軸、Y軸、及びZ軸で表され、重力ベクトルGと直交する地平面と比較して図示してある。必然性はないが通常、内部座標系のZ軸はコンパスが設置される装置に関して下方向にとられる。Y軸はZ軸に直交し、右方向に伸び、X軸はY軸及びZ軸に直交し、コンパス及び装置の概ね「前方」の方向にとられる。X軸は特に意味があり、地球の磁北に対するX軸のオリエンテーションにより方位角が定められる。具体的には、北方向からX軸方向まで測定される角度が方位を定める。この角度は水平面で定義される。コンパスが水平でなく、ある程度のピッチやロール角に置かれる場合、方位角はX軸を水平面に投影することによって定義される。このように、磁気方位角は地球磁界の水平成分からコンパスに係る内部座標系のX軸の水平成分までを測った角度である。この角度は磁北から時計回りに測った場合に正となる。

磁力計３０２は全計算を司るマイクロプロセッサ３０４に磁気センサデータを供給する。マイクロプロセッサ３０４は生成した補償データをメモリ３０６にストアし、後でメモリ３０６からデータを検索し、測定した磁力計データを磁気誤差について補償するのに利用する。重力センサ３０８はマイクロプロセッサ３０４に電子コンパス３００のオリエンテーションに関する情報を提供する。このオリエンテーションデータは磁力計データとともに、電子コンパス３００を較正するのに使用される。重力センサ３０８は、例えば３軸加速度計又はオリエンテーションセンサにより実現される。この代わりに、外部情報を利用してオリエンテーションデータを取得してもよい。例えば、電子コンパス３００が航空機に設置されるのであれば、航空機のフライトコントロール又はナビゲーションシステムを利用してオリエンテーションデータを取得することができ、専用の重力センサは不要となる。電子コンパス３００を較正した後、補償データを利用して磁力計測定データの磁気誤差を補償し、これにより磁気方位の正確な計算が可能となる。

図５は図３に概要を示した電子コンパス３００の具体的構成例を示したものである。磁気抵抗（MR）ブリッジで構成される磁気センサ５０２は、図４に例示するような内部座標系を定める３軸のそれぞれについて磁界測定値を取得する。MRブリッジは小さなサイズで電力消費も少なく、磁気センサ５０２を実現するのに特に適している。磁気センサ５０２は電子コンパスが設置されるビークル又は装置に対して固定される。測定信号はセンサ・エレクトロニクス・モジュール５０４に供給される。センサ・エレクトロニクス・モジュール５０４も、例えば、周知である普通の種類のセンサ処理回路で実現できる。センサ・エレクトロニクス・モジュール５０４の出力は、順番に各モジュール出力を選択するマルチプレクサ５０８に供給される。

図５に示す具体的構成例において、マルチプレクサ５０８は加速度計５１４からの加速度情報も受け取る。これらの加速度計５１４は電子コンパスの内部座標系を構成する３軸のそれぞれに沿った加速度及び重力を測定する。

マルチプレクサ５０８は磁気センサ及び加速度計の信号をアナログ／デジタル（A/D）変換器５１８に供給し、ここで信号はデジタルデータに変換され、変換データはマイクロプロセッサ５２２とメモリ５２３を有するマイクロプロセッサ装置５２０に供給される。マイクロプロセッサ装置５２０はこの情報を利用して方位角を計算する。ただし、磁気誤差があるので補償を行わなければ正確な方位は得られない。

測定磁界を補正するために、マイクロプロセッサ５２２は関連するメモリ５２３にロードされたソフトウェア５２４に従ってデータを処理する。このソフトウェア５２４は、電子コンパス３００を較正し、電子コンパス３００の較正後に磁気補償を行い、磁気補償実行後に磁気方位の計算を行う各種ルーチンを実現したものである。補償処理の結果として、マイクロプロセッサ５２２は方位角を計算し、計算された方位角は、例えば普通のRS-232データインタフェース等を介してユーザに出力される。

測定磁界H_MEASから地球磁界H_EARTHを決定するために、物理モデル式

の代数的逆数を用いて磁気誤差を補償する。この代数的逆数は次のようになる。

係数について次の定義により若干形式を変える。

したがって、上記式は次のようになる。

この式において、＿_EARTHは求めるべきゴールであり、磁力計出力H_MEASから実時間で計算される、地球磁界（真地磁気）に関する正確な推定値（推定値）を表す。行列L_EとベクトルH_PEは補償係数である。これらの補償係数の数値を決定することにより、原モデルのパラメータM及びH_PERMを計算する必要はなくなる。

いったんL_E（行列補償係数）とH_PE（ベクトル補償係数）が定まれば、上記式は３軸空間における任意のオリエンテーションに対して有効であり、電子コンパスを完全な３軸コンパスとして使用することができる。この式（補償式）は物理モデルに基づくもので、この補償式を導出する過程で近似は一切使用していないため、任意の大きさの磁気摂動及び任意の大きさの地球磁界に対して有効である。さらに補償式はH_MEASにのみ依存し、地球磁界、重力ベクトル、局所的水平面等に関する情報等の外部情報を一切必要としない。いったんコンパスを較正したならば、磁気緯度に関わりなく世界のどこでも使用することができる。

本発明に係る補償式において、行列L_Eを構成する９要素は上述したように従来技術で知られている、H_INDUCED及び磁化率Kに関わる９種類の誤差を十分に補償するのに必要な数値である。さらにベクトルH_PEを構成する三組はH_PERMについて補償するのに必要な値である。したがって、この１２要素について適切な値を設定することにより、補償式は摂動界に起因する種類の誤差を全て補償することができる。

行列Kについて述べたように、行列L_Eには特徴的といえる特別な性質はない。すなわち、行列L_Eは対称でなく、その他、全９種類の誤差を完全に補償するのに必要な独立変数の数を低減するような性質を持ち合わせていない。電子コンパスを補償する従来の技術において、しばしば近似が行われるがこの種の近似手法は行列L_Eが対称であって６要素の値だけで行列が特定されると仮定しているのと等価である。その他のケースでは、仮定をさらに付加し、決定すべき係数の数を例えば４個のようにさらに減らそうとした。その目的は較正手続きに係る代数演算を簡単にすることにあった。今日のマイクロプロセッサ及びメモリによれば、かつては意味のあった簡略化の動機は失われている。L_Eについて対称行列を想定するのに等しいかあるいは、さらに簡略化を行うような補償方法はH_INDUCED摂動に関係する全９種類の誤差を補償できない。この種の補償方法は不完全である。本発明は補償係数行列L_Eを構成する９個の数値を決定し、全誤差を十分に補償する方法を提供するものである。

L_Eは３×３行列で、H_PEは３軸ベクトルであるので、この補償式を用いてコンパスを補償するには補償係数を構成する１２個の実数値要素を決定する必要がある。補償精度は係数の精度に依存する。したがって、正確な係数値を求めることが特に有益である。

補償係数値は様々な較正手続きを任意に用いて決定される。定式化した較正においてコンパスの内部座標系を基準座標系として使用し、ユーザの地球座標系よりセンサの視点を優先するのが好ましい。コンパスの内部座標系は図４において、X、Y、Z方向の軸で示される。この座標系において、一定の大きさである地球磁界ベクトルH_EARTHと重力単位ベクトルＧの双方はその間に一定の角度を形成する剛性な対を構成する。これらのベクトルは使用する唯一の座標系を構成するX、Y、Z軸に対して任意のオリエンテーションを取り得る。このように地球座標系は較正手続きにおいて必要でない。

図６Ａは本発明のさらに別の実施形態に基づいた、電子コンパスを補償する補償方法例を示すフローチャートである。補償は電子コンパスの通常動作中に行われる。先ずブロック６０２から始めて、３軸磁界ベクトルH_MEAS及び３軸重力ベクトルG_MEASを測定する。次に、ブロック６０４において、補償係数すなわち、行列係数L_E及びベクトル係数H_PEを電子コンパス３００のメモリから検索する。これらの補償係数は、例えば下方の図６Ｂに説明するような方法で取得される。ブロック６０６において、行列係数L_E及びベクトル係数H_PEの検索値を用いて補償式

を適用する。この補償式の適用により磁力計データの磁気誤差を補正する。補正した磁力計データに基づいてブロック６０８で磁気方位を計算する。ついでフローはブロック６０２に戻り処理を繰り返す。

このようにして電子コンパス３００を補償することによりコンパスをどのオリエンテーションでも使用することができる。さらに、磁気摂動が地球磁界より大きくても全磁気誤差を補償することができる。図１０は非常に大きな磁気摂動を補償する際の本発明の有効性を示している。

しかしながら、電子コンパス３００の補償を行うには前もってコンパスを較正しておく必要がある。較正は、補償式

で使用する行列係数L_Eとベクトル係数H_PEの値を決定する処理を含むもので、時々行えば十分である。図６Ｂは電子コンパスを較正する較正方法例を示したものである。先ず、ブロック６５２から開始し、選択したオリエンテーションと方位の組合せについて較正測定値すなわち３軸磁界ベクトルH_MEAS及び３軸重力ベクトルG_MEASを取得する。測定の各セットにおいてベクトルH_MEAS及びG_MEASはコンパスの座標系を構成するX、Y、Z軸にそれぞれ対応する単位ベクトルx、y、zにより表現される。具体的構成例において、測定は、例えば図７に示す北、東、南及び西のような４つの方位のそれぞれについて４つのオリエンテーション（ピッチとロールの組合せ）のセットについて行い、測定は全部で１６種類となる。なお、特定の方位角を選択する必要はなく、方位間が概ね９０°程度離れていればよい。例えば、北西、北東、南東、及び南西に対応する方位を選択しても等しく動作する。各方位において同一のチルトオリエンテーションを使用してよい。例えば、オリエンテーションのセットとして、ピッチ角とロール角が全て約３０°の条件下で、上／下ピッチ角と左／右ロール角の組合せ、すなわち上／左、上／右、下／左、下／右を使用することができる。なお、方位角と同様に、ピッチ角やロール角の大きさは問題にならない。さらに、較正手続きは測定の順序とは全く無関係である。ただし、測定中はコンパスを安定した状態に保ち、特にGについて、正確な測定及び較正が行われるようにしなければならない。方位、ピッチ及びロール角について測定を偏りのない分布で（ほぼ均等な分布で）行うことにより較正精度を改善することができる。他の実施例において、２４あるいは３２種類の測定セットが使用される。ついでブロック６５４において較正用測定値を電子コンパス３００のメモリにストアする。

次に、このようにして得られた１６種類の測定セット、すなわち１６組あるピッチとロール角の組合せ毎に測定した三組ベクトルH_MEAS及び三組ベクトルG_MEASを使用して、連立方程式を立てそれを解いて図６Ｂのブロック６５６に示すように、L_E及びH_PEを決定する。較正用測定セットにおける各測定値H_MEASについて、

が与えられる。さらに全測定を通じて、＿_EARTHは一定であるがその大きさと傾斜は未知である。これらの性質から係数を計算するのに使用するシステム方程式が導かれる。
先ず、＿_EARTHは一定の大きさ（絶対値）H₀を有する。したがって大きさの平方も定数である。

この方程式を測定セット毎の上記式と組み合わせることにより、記録した較正用測定セットの各値について次式が得られる。

ここにH₀ ²は未知の定数で地球磁界の大きさの平方（２乗）を表す。全部で１６種類の測定を行ったとすると、上記の形式を有する方程式が全部で１６個得られる。この連立方程式に、G_MEASの値は使用されない。これらの方程式は未知の係数L_E及びH_PEについて二次式である。

第二に、傾斜角が一定であることから全測定のベクトル内積について次式が成り立つ。

ここにH_GRAVは全測定を通じて値が変わらない未知の定数で地球磁界の垂直成分を表す。したがって、各測定について次式が成立する。

同じく、全部で（測定セットとして）１６種類の測定を行ったとすると、上記の形式を有するベクトル方程式が１６個得られる。これらの方程式は未知の係数について一次である。まとめると、１６個の二次方程式及び１６個の一次方程式

は、３２個の方程式から成るシステム方程式を構成し、これには１４個の未知数として９要素のL_E、三組のH_PEとともにH₀及びH_GRAVが含まれ、９６個の既知の測定量として１６回の測定毎に三組のH_MEASと三組のG_MEASが含まれる。これらの方程式は未知数について同次であり、H₀ ²とH_GRAVに代入可能な、誤差のない定数ソース値はない。H₀対H_GRAVの比は重要であるがそれぞれの大きさはそうでない。したがって量H₀ ²、G_MEAS、及び未知係数を適当に正規化することにより、未知数を１４個から１３個に減らすことができる。

システム方程式の解法に進む前に注意事項を述べると、H_MEASとG_MEASについて１６回の測定を行うことは必然的でなく例示にすぎない。より少数又は多数の測定を行っても本発明に相反することにはならないが較正の精度に影響は出る。例えば、２４又は３２種の測定をセットとして行えば精度が相当改善される。

システム方程式の解法について説明すると、先ずシステム方程式を一次方程式のセット（連立一次方程式）に再構成する。すなわち、各測定から見積もった H_EARTHを北／東／下の基準系に回転する。この回転は係数値の初期推定値及び３個の新パラメータを使用し、北、東、下の各成分について１６個のシステム方程式を生成し、方程式の数は全部で４８個となる。

詳細に述べると、推定値H_EARTHを北／東／下の基準座標系に回転するために、下記の方程式を使用する。

ここにR1(θ)及びR2(φ)は各G_MEAS値から導出され、H_MEASは対応する磁力計出力値である。ψは各G_MEAS値から導出され、H_MEASは対応する磁力計出力値である。ψは各測定においてH_EASTがゼロになるようにして計算される量である。上記方程式は未知数H_EARTH、L_E及びH_PEについて一次である。H_DOWNに関する１６個の方程式から成る部分集合は前にH_GRAV＝G_MEAS'・＿_EARTHについて使用した部分集合と同じである。

これらの方程式のうち、１６個の「ダウン（down）」方程式は未知の係数H_DOWN、 L_E、及びH_PEに関して真に一次（線形）であり、既存の標準アルゴリズムによりこれらの値について解くことができる。便宜上これらの１６個の方程式のことをまとめてA方程式と呼ぶことにする。A方程式を含む、４８個の方程式から成る全セット（全システム方程式）も同様にまとめてB方程式と呼ぶことにする。L_E及びH_PEの値は先ずA方程式を解いて求め、次に大きなB方程式を解くことによりL_E及びH_PEの値を細分化する。

コンパスの較正を行う上で、H_MEAS及びG_MEASの測定値は測定により既知である一方、補償係数L_E及びH_PEは未知である。１６種の測定の各々についてH_MEAS及びG_MEASの対を記録する。

上述したように、これらの方程式が満足すべき補償式は次の通りである。

ここにH_Eは未知の地球磁界を表し、全測定を通じて同一である。地球座標系においてH_Eは固定ベクトルであるが、内部座標系から見ると、H_Eの数値は、較正中に様々なオリエンテーションを取るので、測定毎に変わる。

未知の補償係数について解くべきシステム方程式を構成するために、この方程式はH_MEASとG_MEASの測定毎に１つずつ、計１６個立てられる。しかしながら、標準の数値解法に合わせるため、これらの方程式は適当に再構成する必要がある。これにより全ての未知数は列ベクトルの要素となり、これに既知の値の大きなシステム行列が掛けられる。再構成は以下のステップで行われる。

例えば、G_MEASとH_Eのベクトル内積を取ることにより１６個の一次方程式が得られる。ベクトル内積の結果は地球磁界の垂直下成分、すなわちH_GRAVであり、これは未知であるが全測定を通じて同一である。この方程式は次のように書かれる。

次に、行列L_Eは３つの列ベクトルにより次のように表される。

そしてH_GRAVに関する上記システム方程式は、測定したH_MEASとG_MEASの任意の対について次のようになる。

これは、右辺のH_GRAVが単に数値であってベクトルではないので、スカラー方程式である。左辺のG_MEASに関する最初の積（product）はG_MEASと３列のL_Eの各々との積であるので行ベクトルとなり、これを測定したH_Mである列ベクトルに掛ける。

全ての未知数を１つの長い列ベクトルに集めるために、この方程式を再構成すると次のような等価な形式となる。

左辺の行ベクトルは全測定値を含んでおり、測定毎に異なる。Pで記す右辺の列ベクトルは全ての未知補償係数から構成され、全測定に共通である。いずれも１三組のベクトルである。右辺のゼロは単一要素でありベクトルではない。

この方程式は、Pを決定すべき係数ベクトル、Aを測定値から構成される行ベクトルとして、次のように簡単に書き表せる。

この方程式を、較正手続きのために測定毎にオリエンテーションが異なる、全部で１６種類の測定の各々について書き直すと、Aベクトルは測定毎に１行増える、１６×１３の大きなシステム行列になる。このシステムはAシステム方程式と呼ばれ、未知数について一次（線形）である。このAシステムを係数について解く。ただし、Aシステムは全システムではないので、完全な解を与えない。Aシステム単独では一次的な近似しか得られない。

全システム方程式を構成するには、上述した方法と同様に、ベクトルPにH_Eの大きさ（絶対値）の２乗あるいはH_Eの水平・北成分を付加する必要がある。目的の大きなシステムのために、Aシステムで述べたのと同じステップを使用することにより、Pベクトルに１４番目の要素が追加された、全システム方程式を定式化する。これにより、システムは１４個の未知数を含む、３２個の方程式を構成する。

また、上述したように、各H_Eを北、東、下の成分に回転させてもよく、この回転によりこれらの成分はPベクトルの最後の三組となり、全部で１５となる。これにより全システムは１５個の未知数を含む４８個の方程式から構成される。これらの代替システム方程式はいずれも補償係数の決定を行う較正計算法として実証された。

ベクトルPについて解を求める場合、その要素値を補償係数及び地球磁界の適当な成分に書き改める。すなわち、ベクトルPの最初の９値を行列L_Eの３列に、次の三組をベクトルH_PEに書き改める。残りの要素は地球磁界成分に関する適当な値に書き改められるが、具体的にどの成分が該当するかについては使用した計算方法及びシステム方程式に依存する。

さらに、測定結果をシステムに組み込む（又はシステムに使用する）順序は得られる解になんら違いをもたらさない。これは非常に有意義な実用的利点であり、較正手続きに取り込む測定結果の順序は問題にならない。

この較正手続きの更なる利点として、対を成す測定値H_MEAS及びG_MEAS以外に情報を付加する必要がない。すなわち、測定を行ったときのピッチやロールの傾斜角度や方位角に関する数値情報をシステム方程式に与える必要はない。また、較正に係るコンパスの試験範囲も必要なく、既知の方位における１つの照準すら不要である。

本発明の具体的構成例において、「A」及び「B」のシステム方程式を解くために、例えばマスワークス社（Math Works Inc）から入手可能なMATLABのようなコンピュータ言語を使用する。この構成例によれば、両システムを解くために特異値分解（SVD）法を以下のように使用する。本開示に添付するのは更なる実施可能資料を含むアペンディックスであり、この資料に較正計算を行う方法の詳細が示される。ここに、これらのアルゴリズムのプロトタイプを作製する際に使用できるMATLAB（登録商標）プログラムのリストが与えられる。

A・P＝0＋δ
M個の未知数を含むN個の連立一次方程式を行列形式で書き表すと、システム行列AはN×M行列となり、PはM個の未知数を表すM要素の列ベクトルになる。不能形のシステム、すなわちN＞Mであるシステムでは、行列Aは正方行列でなく、その逆が存在せず、Pについて一意的な解が存在しない。したがって、行列積A・Pを取ると最良でも残余としてN要素の列ベクトルが発生する。

ベクトルPに関する最適解（best solution）を求めるには、SVD解析に特有な、｜δ｜²を最小にする最小２乗法を用いて上記システム方程式（連立方程式）を解く。
非正方行列Aの特異値分解（SVD）は３つの行列U、S、及びVに因数分解することであり、式で示すと次のようになる。

ここにSは特異値行列であり、左上のM×Mの正方ブロックが対角行列になっている。これらの特異値の一部はゼロであり得る。UとVは正方かつ直交行列、すなわちユニタリ行列である。すなわち、U'・U≡I_Nで、V'・V≡I_Mである。これらの３行列を計算するために既存の各種数値ソフトウェア解析法を任意に使用できる。いったんこれらの行列を計算すると、最小の特異値、行列Sの（通常の順序における）最後の要素S_MM、及び行列U及びVのM番目の列ベクトルから、未知量Pに関する完全な最小２乗解が得られる。すなわち、最小の２乗和平方根（RSS）による残差値σは次式で与えられる。

σ≡S_MM
未知量に対する解ベクトルPは次のように定義される。

すなわち、SVD行列因子VのM番目の列である。解ベクトルPはM×1の列ベクトルであり、適宜正規化し直すことができる。最後に、方程式の残差δはσに比例するN×1のベクトルであり、次式で定義される。

U、S、及びV行列の他の要素値は意味がなく、捨てられる。
解ベクトルを解析して行列係数L_E及びベクトル係数H_PEの値を求める。すなわち、解ベクトルPの成分をL_E、H_PE、H_NORTH、H_EAST、及びH_DOWNの要素値に読み込む。このようにして行列係数L_Eとベクトル係数H_PEが求められると、電子コンパスの較正は完了する。具体的構成例では、「B」システムを反復して解いてより正確な係数L_E及びH_PEを求める。図８は最小特異値とともに反復回数に伴う、その収束を示したものである。

「東」方程式は解としてゼロを含み、残差が小さいことが観察された。したがって、別の具体的構成例では、これらの「東」方程式を「B」システムから消去することにより収束が効果的に改善される。この構成例において、「Ｂ」システムは、解セットの一部を構成するH_NORTH及びH_DOWNが未知量である、３２×１４システムとして再定義される。１３個の未知量を含む１６個の方程式から成る「A」システムを先ず解いて係数L_E及びH_PEの初期推定を出す。次にこれらの結果を、１４個の未知量を含む３２個の方程式から成る「B」システムに適用してL_E及びH_PEについて細分化した値を求める。次に、「B」システムを再度解いてL_E及びH_PEをさらに細分化する。この構成例では、図９に示すように、反復を２回行うだけで十分な収束を得ることができる。

本発明の他の具体例では、測定結果から代替方程式セットを導いてL_E及びH_PEを決定する。上述したように、個々の測定について、地球磁界ベクトルの推定値＿_EARTHは次式で与えられる。

地球磁界に関する別の等価な基本性質から代替方程式が立てられる。すなわち、全測定を通じてH_EARTHの垂直成分と水平（北）成分はそれぞれ一定である。垂直成分は一定であるので、ベクトル内積

において、H_GRAVは未知の定数となる。１６種類の測定の各々についてこれらの方程式をまとめると次のシステム方程式（連立方程式）が得られる。

一定の水平（北）成分から１６種類の測定の各々について次の方程式が導かれる。

ここにH_NORTHは未知の定数である。
以上により、３２個の方程式から成るシステムが定義される。これらの方程式をL_E、H_PE、H_GRAV、及び、H_NORTHの値について解くために、例えばコンピュータで実行される数値解法を利用する。方程式

は未知の係数について一次（線形）であり、先ずこれについて解く。方程式

は水平成分を計算するのにL_E及びH_PEの初期推定値を必要とする。見積もったH_EARTHを北／東／下の基準座標系に回転するときの計算について上述したのと同様の仕方でG_MEASに関して回転の計算を行う。

次に、全方程式を、３２個の一次（線形）方程式から成る単一のシステムとみなして解く。方程式の解き方は「A」及び「B」システムに関して述べたのとほぼ同様の方法で行われる。すなわち、

の形式を有する最初の１６個の方程式を先ず解いてL_E及びH_PEの予備推定値を求め、ついでこれを用いて大きなシステム方程式を解くことでL_E及びH_PEに対する新たな推定値を得る。そして再度、大きなシステム方程式を解いてL_E及びH_PEの推定値をさらに細分化する。僅か２回の反復計算でH_GRAV及びH_NORTHについて十分に収束した値が得られることが観察された。図６Ｂのブロック６５８に示すように、較正計算完了後、得られた補償係数をメモリ５３２に記憶する。

図６Ａに戻り、コンパスの通常動作において、ブロック６０２に示すように、H_MEAS及びG_MEASについて連続的かつ反復して測定を行う。ブロック６０４において補償係数L_E及びH_PEの値をメモリから検索した後、ブロック６０６に進み、これらの補償係数値により磁力計の測定出力データを補正して正確な方位データを求める。すなわち、

の式に行列補償係数L_E及びベクトル補償係数H_PEを適用して、磁力計の測定出力データH_MEASからターゲットである、真地球磁界（真地磁気）の推定値（推定値）＿_EARTHを取得する。

上記のようにして電子コンパスを補償することにより、任意のオリエンテーションで電子コンパスを使用できる。さらに、地球磁界より大きな磁気摂動を含む、全ての磁気誤差が補償される。

上記較正法（補償法）にはその他の利点も確実に認められる。すなわち、磁気誤差を補償すること以外に、上記補償法は磁力計センサに関する様々な種類のセンサ誤差も正確に補償する。例えば、センサのバイアスオフセットは永久磁界摂動と全く同様であり、H_PE係数に計上される（織り込まれる）。その他の磁力計センサ誤差、例えば不等スケールファクタ、軸の一致や非直交性、及び磁力計と加速度計の軸の不一致等は誘導磁界摂動と全く同様であり、L_E係数に計上される（織り込まれる）。この種の誤差は全て本発明の補償法及び較正手続きにより十分かつ正確に補償される。

上述した各種実施の形態は単に例示を目的とするものであり、本発明を制限するものではない。当業者には明らかなように、ここで図示し説明した例示の実施形態及び応用形態に厳密に従うことなく、また、特許請求の範囲に記載される、本発明の精神及び範囲から逸脱することなく、種々の変形及び変更をこれらの実施形態に施すことができる。例えば、補償係数の値はコンパスの使用中に連続的に更新することができる。さらに、各種フィルタ、例えばカルマンフィルタを較正処理で使用することができる。

図１Ａは摂動のない水平面における地球磁界の軌跡を示す。 図１Ｂは地球磁界の軌跡に対して、磁力計により測定される、永久磁界及び誘導磁界の双方による磁気摂動が与える合成影響を示す。 ３軸磁力計により測定される、永久磁界及び誘導磁界の双方による磁気摂動が地球磁界の軌跡に対して与える合成影響を示す。 本発明の一実施形態に基づいた、電子コンパス補償構成を示すブロック図である。 図３の電子コンパスの座標系とともに局所的な地球座標系を示す。 本発明の別の実施形態に基づいた、図３の電子コンパスの具体的構成例を示すブロック図である。 図６Ａは本発明のさらに別の実施形態に基づいた、電子コンパスを補償する補償方法例を示すフローチャートである。 図６Ｂは本発明のさらに別の実施形態に基づいた、電子コンパスを較正する較正方法例を示すフローチャートである。 本発明の特定の実施形態に基づいた、電子コンパスを較正する際に使用するオリエンテーションセットの具体例を示す。 本発明の特定の実施形態に基づいた、較正方法を用いて得られた数値補償係数の解セットに係る最小特異値及び収束を示す２つのグラフである。 本発明の別の特定実施形態に基づいた、較正方法を用いて得られた数値補償係数の解セットに係る最小特異値及び収束を示す２つのグラフである。 本発明の実施形態に基づいた、補償方法を用いて行われる数値補償による電子コンパス精度の効果を図示する。

Claims (10)

1. ３軸について磁界強度を表す測定磁界ベクトルH_MEASを、オリエンテーションと方位の複数の各組合せに対して与える３軸磁力計と、
   ３軸について重力場強度を表す測定重力ベクトルG_MEASを、オリエンテーションと方位の複数の各組合せに対して与える３軸重力センサと、
   前記３軸磁力計及び前記３軸重力センサの装置と通信可能に結合されたプロセッサと、を備え、
   前記プロセッサは、システム方程式を解くことにより、複数の測定磁界ベクトル H _MEAS 、及び複数の測定重力ベクトル G _MEAS の関数として、行列補償係数 L _E 及びベクトル補償係数 H _PE を計算し、
   前記プロセッサは行列補償係数L_Eとベクトル補償係数H_PEとを使用し、次式、
   

   

   を適用することにより測定磁界データを補正し、そして前記補正磁界データの関数として方位データを計算する、ここに＿_EARTHは真の地球磁界のための補償値である、
   前記行列補償係数 L _E は誘導磁界摂動および磁化率行列に関連する誤差を補償する要素を含み、
   前記ベクトル補償係数 H _PE は永久磁界摂動を補償する要素を含む、
   電子コンパス装置。
2. 前記プロセッサと結合され、前記行列補償係数L_E及び前記ベクトル補償係数H_PEを記憶するメモリをさらに備える、請求項１に記載の電子コンパス装置。
3. 前記複数の測定磁界ベクトルH_MEAS、及び複数の測定重力ベクトルG_MEASは４つの方位のそれぞれについて４つのオリエンテーションで測定され、前記４つの方位は、９０°離れている、請求項１に記載の電子コンパス装置。
4. 前記４つの各方位に対する前記４つのオリエンテーションは、（ｉ）第１のオリエンテーションは上ピッチ角と左ロール角により特徴付けられ、（ｉｉ）第２のオリエンテーションは上ピッチ角と右ロール角により特徴付けられ、（ｉｉｉ）第３のオリエンテーションは下ピッチ角と左ロール角により特徴付けられ、（ｉｖ）第４のオリエンテーションは下ピッチ角と右ロール角により特徴付けられる、請求項３に記載の電子コンパス装置。
5. 各前記ピッチ角と各前記ロール角は３０°である、請求項４に記載の電子コンパス装置。
6. 前記システム方程式は、次の形式を有する複数の方程式を含み、
   

   

   ここにH₀は＿_EARTHの大きさである、及び次の形式を有する複数の方程式を含む、
   

   

   ここにH_GRAVは＿_EARTHの未知な垂直成分を表す定数である、
   請求項１に記載の電子コンパス装置。
7. さらに前記プロセッサは、前記複数の測定磁界ベクトルH_MEASのそれぞれについて＿_EARTHを基準の北／東／下座標系に回転させることにより全システム方程式を導き出す、請求項１に記載の電子コンパス装置。
8. 前記システム方程式は、前記測定磁界ベクトル H _MEAS 及び前記測定重力ベ クトル G _MEAS の対に対して解かれる一次方程式の第１の部分集合を有する、請求項１に記載の電子コンパス装置。
9. さらに前記プロセッサは、
   前記一次方程式の前記第１の部分集合を解いて前記行列補償係数L_Eの推定値及び前記ベクトル補償係数H_PEの推定値を計算し、
   前記行列補償係数L_Eの前記推定値及び前記ベクトル補償係数H_PEの前記推定値を用いて前記全システム方程式を解くことにより前記行列補償係数L_E及び前記ベクトル補償係数H_PEを計算する、
   請求項８記載の電子コンパス装置。
10. さらに前記プロセッサは、前記システム方程式を反復法で解くことにより前記行列補償係数L_E及び前記ベクトル補償係数H_PEの計算値を改善する、請求項９記載の電子コンパス装置。

Images