(57)【要約】【課題】計測対象とする微弱な磁場に対して大きな外部磁場が存在する場合でも、計測対象とする磁場を精度良く計測することができる磁気計測システムを提供する。【解決手段】磁気計測システム１は、第１磁場と第２磁場とを計測する心磁センサー１０と、第２磁場を計測するノイズ磁気センサー３０と、ノイズ磁気センサー３０の計測値と多変数多項式とを用いて心磁センサー１０における第２磁場の近似値を計算する磁気計測装置２とを備える。磁気計測装置２は、心磁センサー１０の計測値から第２磁場の近似値を引き算する。【選択図】図１
2016217930.tif 000002
【特許請求の範囲】
【請求項１】
第１磁場と第２磁場とを計測する第１磁気センサーと、
前記第２磁場を計測する第２磁気センサーと、
前記第２磁気センサーの計測値と多変数多項式とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を備えたことを特徴とする磁気計測システム。
【請求項２】
第１磁場と第２磁場とを計測する第１磁気センサーと、
前記第２磁場を計測する第２磁気センサーと、
前記第２磁気センサーの計測値と非線形多項式とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を備えたことを特徴とする磁気計測システム。
【請求項３】
請求項１または２に記載の磁気計測システムであって、
前記処理装置は、前記第１磁気センサーの計測値から前記第２磁場の近似値を引き算することを特徴とする磁気計測システム。
【請求項４】
請求項１に記載の磁気計測システムであって、
前記多変数多項式は数式１で表されることを特徴とする磁気計測システム。
【数１】
2016217930.tif 000071
（なお、数式１で、ａij（ｉは１から３の整数、ｊは１から７の整数）は係数であり、ｘ、ｙ、ｚは磁場の近似値Ｂの空間座標であり、Ｂiは磁場の近似値Ｂの第ｉ成分である。）
【請求項５】
請求項２に記載の磁気計測システムであって、
前記非線形多項式は前記数式１で表されることを特徴とする磁気計測システム。
【請求項６】
請求項１から５のいずれか一項に記載の磁気計測システムであって、
前記第２磁気センサーの計測値から最小二乗法を用いて前記多項式の解を求めることを特徴とする磁気計測システム。
【請求項７】
請求項４または５に記載の磁気計測システムであって、
前記第２磁気センサーは、前記第２磁場の磁場ベクトル成分の２１個以上を計測することを特徴とする磁気計測システム。
【請求項８】
請求項７に記載の磁気計測システムであって、
前記数式１の未知数から作られる第１行列を数式２で表されるａとし、前記第２磁気センサーの計測値から作られる第２行列を数式３で表されるＭとし、前記第２磁気センサーの位置から作られる第３行列を数式４で表されるＰとし、前記第１行列ａを数式５または数式６にて求めることを特徴とする磁気計測システム。
【数２】
2016217930.tif 000072
【数３】
2016217930.tif 000073
【数４】
2016217930.tif 000074
【数５】
2016217930.tif 000075
【数６】
2016217930.tif 000076
（なお、数式５でＰ-1は第３行列Ｐの逆行列であり、数式６でＰ+1は第３行列Ｐの疑似逆行列である。）
【請求項９】
請求項７に記載の磁気計測システムであって、
前記数式１の未知数から作られる第１ベクトルを数式７で表されるｂとし、前記第２磁気センサーの計測値から作られる第２ベクトルを数式８で表されるＮとし、前記第２磁気センサーの位置から作られる第４行列を数式９で表されるＱとし、前記第１ベクトルｂを数式１０または数式１１にて求めることを特徴とする磁気計測システム。
【数７】
2016217930.tif 000077
【数８】
2016217930.tif 000078
【数９】
2016217930.tif 000079
【数１０】
2016217930.tif 000080
【数１１】
2016217930.tif 000081
【請求項１０】
請求項４に記載の磁気計測システムであって、
前記多変数多項式は数式１２を考慮した前記数式１で表されることを特徴とする磁気計測システム。
【数１２】
2016217930.tif 000082
【請求項１１】
請求項５に記載の磁気計測システムであって、
前記非線形多項式は前記数式１２を考慮した前記数式１で表されることを特徴とする磁気計測システム。
【請求項１２】
請求項１０または１１に記載の磁気計測システムであって、
前記第２磁気センサーは、前記第２磁場の磁場ベクトルの１７成分以上を計測することを特徴とする磁気計測システム。
【請求項１３】
請求項１２に記載の磁気計測システムであって、
前記数式１の未知数から作られる第３ベクトルを数式１３で表されるｃとし、前記第２磁気センサーの計測値から作られる第２ベクトルを前記数式８で表されるＮとし、前記第２磁気センサーの位置から作られる第５行列を数式１４および数式１５、または数式１６で表されるＳとし、前記第３ベクトルｃを数式１７にて求めることを特徴とする磁気計測システム。
【数１３】
2016217930.tif 000083
【数１４】
2016217930.tif 000084
【数１５】
2016217930.tif 000085
【数１６】
2016217930.tif 000086
【数１７】
2016217930.tif 000087
【請求項１４】
第１磁場と第２磁場とを計測する第１磁気センサーと、
前記第１磁気センサーの周囲に分布する第１の第２磁気センサーと第２の第２磁気センサーと、
前記第１の第２磁気センサーの計測値と前記第２の第２磁気センサーの計測値とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を含み、
前記第１磁気センサーは、前記第１の第２磁気センサーと前記第２の第２磁気センサーとの重心を含む位置に配置されたことを特徴とする磁気計測システム。
【請求項１５】
請求項１４に記載の磁気計測システムであって、
第３の第２磁気センサーと第４の第２磁気センサーとをさらに含み、
前記第３の第２磁気センサーと前記第４の第２磁気センサーとは前記重心に対して対称となる位置に配置されており、
前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、は交差することを特徴とする磁気計測システム。
【請求項１６】
第１磁場と第２磁場とを計測する第１磁気センサーと、
前記第１磁気センサーの周囲に分布する第１の第２磁気センサーと第２の第２磁気センサーと第３の第２磁気センサーと第４の第２磁気センサーと、
前記第１の第２磁気センサーの計測値と前記第２の第２磁気センサーの計測値と前記第３の第２磁気センサーの計測値と前記第４の第２磁気センサーの計測値とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を含み、
前記第１磁気センサーは、前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、の交差部を含む位置に配置されたことを特徴とする磁気計測システム。
【請求項１７】
請求項１５または１６に記載の磁気計測システムであって、
前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、は直交することを特徴とする磁気計測システム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、磁気計測システムに関する。
【背景技術】
【０００２】
地磁気に比べて微弱な心臓の磁場や脳の磁場等を測定するための磁気計測装置が提案されている（例えば、特許文献１参照）。磁気計測装置は非侵襲であり、被検体（生体）に負荷をかけずに臓器の状態を計測できる。特許文献１には、生体が発する磁場を計測するセンサー（ピックアップコイル）と、ノイズとなる環境磁場を計測するセンサー（リファレンスコイル）とを有し、ピックアップコイルが計測する磁場に含まれる環境磁場（ノイズ）をリファレンスコイルが計測する環境磁場に基づいて除去する構成の生体磁気計測装置が開示されている。
【０００３】
特許文献１に記載の生体磁気計測装置では、リファレンスコイルの数がピックアップコイルよりも少なく、一部のピックアップコイルに対応する位置に配置されたリファレンスコイルの計測データから各ピックアップコイルの位置における環境磁場が求められる。また、リファレンスコイルが対応して配置されない他のピックアップコイルの位置における環境磁場は、リファレンスコイルの計測データから線形補間して推定することにより求められる。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特開平５−２９７０８７号公報
【発明の概要】
【発明が解決しようとする課題】
【０００５】
特許文献１に記載の生体磁気計測装置では、リファレンスコイルがピックアップコイルとは異なる（離れた）位置に配置されている。そのため、対応するリファレンスコイルの計測データから求められた環境磁場は、ピックアップコイルの位置における環境磁場とは必ずしも一致しない。また、このようなリファレンスコイルの計測データから、対応するリファレンスコイルが配置されていないピックアップコイルの位置における環境磁場が推定される。そして、これらの環境磁場の計測データおよび推定データをピックアップコイルの位置における計測データから差し引いて、生体が発する磁場が求められる。したがって、ピックアップコイルの位置における環境磁場として求められたデータに誤差が生じ易く、生体が発する磁場を精度良く計測できないおそれがある。
【課題を解決するための手段】
【０００６】
本発明は、上述の課題を解決するためになされたものであり、以下の形態または適用例として実現することが可能である。
【０００７】
［適用例１］
本適用例にかかる磁気計測システムは、第１磁場と第２磁場とを計測する第１磁気センサーと、前記第２磁場を計測する第２磁気センサーと、前記第２磁気センサーの計測値と多変数多項式とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を備えたことを特徴とする。
【０００８】
本適用例の構成によれば、第２磁気センサーが計測する第２磁場の計測値と多変数多項式とを用いて、処理装置が第１磁気センサーにおける第２磁場の近似値を計算する。そのため、第１磁場と第２磁場とを計測する第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。
【０００９】
［適用例２］
本適用例にかかる磁気計測システムは、第１磁場と第２磁場とを計測する第１磁気センサーと、前記第２磁場を計測する第２磁気センサーと、前記第２磁気センサーの計測値と非線形多項式とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を備えたことを特徴とする。
【００１０】
本適用例の構成によれば、第２磁気センサーが計測する第２磁場の計測値と非線形多項式とを用いて、処理装置が第１磁気センサーにおける第２磁場の近似値を計算する。そのため、第１磁場と第２磁場とを計測する第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。
【００１１】
［適用例３］
上記適用例に係る磁気計測システムであって、前記処理装置は、前記第１磁気センサーの計測値から前記第２磁場の近似値を引き算することが好ましい。
【００１２】
本適用例の構成によれば、処理装置が、第１磁場と第２磁場とを計測する第１磁気センサーの計測値から、第１磁気センサーにおける第２磁場の近似値を引き算するので、これにより第１磁気センサーが計測する第１磁場が求められる。第２磁場の近似値は第２磁気センサーの計測値に基づいて高精度に計算されているので、第１磁場を高精度に算出することができる。
【００１３】
［適用例４］
上記適用例に係る磁気計測システムであって、前記多変数多項式は数式１で表されることが好ましい。
【００１４】
【数１】
2016217930.tif 000003
（なお、数式１で、ａij（ｉは１から３の整数、ｊは１から７の整数）は係数であり、ｘ、ｙ、ｚは磁場の近似値Ｂの空間座標であり、Ｂiは磁場の近似値Ｂの第ｉ成分である。）
【００１５】
本願発明者が鋭意研究したところ、数式１に示される多変数多項式により、計測対象空間の任意の位置における磁場ベクトルを精度良く近似できることが判明した。数式１の右辺の第１項ａi1は全体での平行磁場を表し、第２項ａi2ｘから第４項ａi4ｙは線形磁場を表し、第５項ａi5ｘｙから第７項ａi7ｚｘは交互磁場（磁場のねじれ）を表す。磁場には、ビオ・サバールの法則にて電流要素ベクトルと位置ベクトルとの外積項に比例する成分が存在する。そこで、本願発明者は、磁場の近似式にねじれ成分を表す第５項から第７項のｘｙ項とｙｚ項とｚｘ項とを導入した。本適用例の構成によれば、数式１に示される多変数多項式を用いて第２磁場の近似値を計算するので、第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。
【００１６】
［適用例５］
上記適用例に係る磁気計測システムであって、前記非線形多項式は前記数式１で表されることが好ましい。
【００１７】
本適用例の構成によれば、数式１に示される非線形多項式を用いて第２磁場の近似値を計算するので、第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。
【００１８】
［適用例６］
上記適用例に係る磁気計測システムであって、前記第２磁気センサーの計測値から最小二乗法を用いて前記多項式の解を求めることが好ましい。
【００１９】
本適用例の構成によれば、第２磁気センサーの計測値から最小二乗法を用いて数式１の解を求めて第２磁場の近似値を計算するので、精度良く第２磁場の近似値を算出することができる。
【００２０】
［適用例７］
上記適用例に係る磁気計測システムであって、前記第２磁気センサーは、前記第２磁場の磁場ベクトル成分の２１個以上を計測することが好ましい。
【００２１】
本適用例の構成によれば、数式１に示される多変数多項式または非線形多項式において、ＸＹＺの各成分に対して７つの未知数が存在するので、求めるべき未知数は３×７＝２１個となる。したがって、第２磁場の磁場ベクトル成分の２１個以上を計測することにより、数式１を用いて第２磁場の近似値を高精度に計算することができる。
【００２２】
［適用例８］
上記適用例に係る磁気計測システムであって、前記数式１の未知数から作られる第１行列を数式２で表されるａとし、前記第２磁気センサーの計測値から作られる第２行列を数式３で表されるＭとし、前記第２磁気センサーの位置から作られる第３行列を数式４で表されるＰとし、前記第１行列ａを数式５または数式６にて求めることが好ましい。
【００２３】
【数２】
2016217930.tif 000004
【００２４】
【数３】
2016217930.tif 000005
【００２５】
【数４】
2016217930.tif 000006
【００２６】
【数５】
2016217930.tif 000007
【００２７】
【数６】
2016217930.tif 000008
（なお、数式５でＰ-1は第３行列Ｐの逆行列であり、数式６でＰ+1は第３行列Ｐの疑似逆行列である。）
【００２８】
本適用例の構成によれば、第１磁気センサーと第２磁気センサーとが配置された空間の任意の位置を位置ベクトルｒで表し、数式１の未知数から作られる３行７列の第１行列をａとすると、任意の位置における第２磁場の磁場ベクトルＢは数式２で表される。第２磁気センサーの計測値（α個の検出磁場ベクトルＢk）から、数式２で表される３行α列の第２行列Ｍが作られる。そして、第２磁気センサーの位置（α個の磁気センサー項ベクトルＲk）から、数式４で表される７行α列の第３行列Ｐが作られる。第２磁気センサーの数αが７の場合は、第３行列Ｐの逆行列が存在するので、第１行列ａは、数式５のように第２行列Ｍに第３行列Ｐの逆行列（Ｐ-1）を掛けて求められる。また、第２磁気センサーの数αが８以上の場合は、第３行列Ｐの逆行列が存在しないので、数式６のように第２行列Ｍに第３行列Ｐの疑似逆行列（Ｐ+）を掛けて求められる。このように、数式１の未知数である第１行列ａが数式５または数式６にて求められるので、精度良く第２磁場の近似値を算出することができる。
【００２９】
［適用例９］
上記適用例に係る磁気計測システムであって、前記数式１の未知数から作られる第１ベクトルを数式７で表されるｂとし、前記第２磁気センサーの計測値から作られる第２ベクトルを数式８で表されるＮとし、前記第２磁気センサーの位置から作られる第４行列を数式９で表されるＱとし、前記第１ベクトルｂを数式１０または数式１１にて求めることが好ましい。
【００３０】
【数７】
2016217930.tif 000009
【００３１】
【数８】
2016217930.tif 000010
【００３２】
【数９】
2016217930.tif 000011
【００３３】
【数１０】
2016217930.tif 000012
【００３４】
【数１１】
2016217930.tif 000013
【００３５】
本適用例の構成によれば、数式１の未知数（３×７個）を１列に並べて、数式７で表される２１行１列の第１ベクトルｂが作られる。第２磁気センサーの計測値（α個の検出磁場ベクトルＢk）から数式８で表される（３×α）行１列の第２ベクトルＮが作られる。そして、第２磁気センサーの位置（α個の磁気センサー項ベクトルＲk）から、数式９で表される第４行列Ｑが作られる。数式９の中式で、行ベクトルＲkTは磁気センサー項ベクトルＲkの転置行列で１行７列の行ベクトルであり、ゼロベクトル０は行列要素が全てゼロの１行７列の行ベクトルである。第２磁気センサーの数αが７の場合は、第４行列Ｑの逆行列が存在するので、未知数である第１ベクトルｂは、数式１０のように第４行列Ｑの逆行列（Ｑ-1）を第２ベクトルＮに掛けて求められる。また、第２磁気センサーの数αが８以上の場合は、第４行列Ｑの逆行列が存在しないので、数式１１のように第４行列Ｑの疑似逆行列（Ｑ+）を第２ベクトルＮに掛けて求められる。このように、数式１の未知数である第１ベクトルｂが数式１０または数式１１にて求められるので、精度良く第２磁場の近似値を算出することができる。
【００３６】
［適用例１０］
上記適用例に係る磁気計測システムであって、前記多変数多項式は数式１２を考慮した前記数式１で表されることが好ましい。
【００３７】
【数１２】
2016217930.tif 000014
【００３８】
本適用例の構成によれば、数式１を解く際に、磁場に関するガウスの法則を適用することで、２１個の未知数のうち４個は数式１２で表される。したがって、数式１２の左辺の４つの未知数を解く必要がなくなるので、求めるべき未知数を１７個に減少させることができる。
【００３９】
［適用例１１］
上記適用例に係る磁気計測システムであって、前記非線形多項式は前記数式１２を考慮した前記数式１で表されることが好ましい。
【００４０】
本適用例の構成によれば、数式１を解く際に、磁場に関するガウスの法則を適用することで、２１個の未知数のうち４個は数式１２で表される。したがって、数式１２の左辺の４つの未知数を解く必要がなくなるので、求めるべき未知数を１７個に減少させることができる。
【００４１】
［適用例１２］
上記適用例に係る磁気計測システムであって、前記第２磁気センサーは、前記第２磁場の磁場ベクトルの１７成分以上を計測することが好ましい。
【００４２】
本適用例の構成によれば、求めるべき未知数が１７個となるので、第２磁場の磁場ベクトルの１７成分以上を計測すれば、第２磁場の近似値を高精度に計算することができる。
【００４３】
［適用例１３］
上記適用例に係る磁気計測システムであって、前記数式１の未知数から作られる第３ベクトルを数式１３で表されるｃとし、前記第２磁気センサーの計測値から作られる第２ベクトルを前記数式８で表されるＮとし、前記第２磁気センサーの位置から作られる第５行列を数式１４および数式１５、または数式１６で表されるＳとし、前記第３ベクトルｃを数式１７にて求めることを特徴とする磁気計測システム。
【００４４】
【数１３】
2016217930.tif 000015
【００４５】
【数１４】
2016217930.tif 000016
【００４６】
【数１５】
2016217930.tif 000017
【００４７】
【数１６】
2016217930.tif 000018
【００４８】
【数１７】
2016217930.tif 000019
【００４９】
本適用例の構成によれば、２１個の未知数から数式１２の左辺の４つの未知数を除いた１７個の未知数を１列に並べて、数式１３で表される１７行１列の第３ベクトルｃが作られる。第２磁気センサーの位置（α個の磁気センサー項ベクトルＲk）から、数式１６で表される３α行１７列の第５行列Ｓが作られる。ここで、数式１６において、３行毎に一つの磁気センサー項ベクトルＲkが対応しており、３ｋ−２行目と３ｋ−１行目と３ｋ行目とを３行１７列の部分行列Ｔkにて表すと、部分行列Ｔkは数式１５で表される。第５行列Ｓは、部分行列Ｔｋを用いると、ｋ＝１の部分行列Ｔ1からｋ＝αの部分行列Ｔαまでのα個をα行１列に並べた行列となり、数式１４にて表される。１７個の未知数に対して第２磁気センサーの数αが６以上であれば、数式８で表される第２ベクトルＮは１８行以上となるので、最小二乗法を適応して１７個の未知数を特定できる。この場合、第５行列Ｓの逆行列は存在しないので、未知数である第３ベクトルｃは、数式１７のように第５行列Ｓの疑似逆行列（Ｓ+）を第２ベクトルＮに掛けて求められる。
【００５０】
［適用例１４］
本適用例にかかる磁気計測システムは、第１磁場と第２磁場とを計測する第１磁気センサーと、前記第１磁気センサーの周囲に分布する第１の第２磁気センサーと第２の第２磁気センサーと、前記第１の第２磁気センサーの計測値と前記第２の第２磁気センサーの計測値とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を含み、前記第１磁気センサーは、前記第１の第２磁気センサーと前記第２の第２磁気センサーとの重心を含む位置に配置されたことを特徴とする。
【００５１】
本適用例の構成によれば、第１磁場と第２磁場とを計測する第１磁気センサーが、第１の第２磁気センサーと第２の第２磁気センサーとの重心を含む位置に配置されている。すなわち、第１の第２磁気センサーと第２の第２磁気センサーとは、第１磁気センサーに関して対称的に配置されている。そのため、処理装置は、第１の第２磁気センサーの計測値と第２の第２磁気センサーの計測値とを同じ重要度で用いるので、第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。
【００５２】
［適用例１５］
上記適用例に係る磁気計測システムであって、第３の第２磁気センサーと第４の第２磁気センサーとをさらに含み、前記第３の第２磁気センサーと前記第４の第２磁気センサーとは前記重心に対して対称となる位置に配置されており、前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、は交差することが好ましい。
【００５３】
本適用例の構成によれば、第３の第２磁気センサーと第４の第２磁気センサーとは、第１の第２磁気センサーと第２の第２磁気センサーとの重心に対して対称となる位置に配置されているので、第１磁気センサーに関して対称的に配置されている。そのため、処理装置は、第３の第２磁気センサーの計測値と第４の第２磁気センサーの計測値とを同じ重要度で用い、第１磁気センサーにおける第２磁場の近似値を高精度に計算することができる。さらに、第１の第２磁気センサーと第２の第２磁気センサーと第３の第２磁気センサーと第４の第２磁気センサーとが直線上に位置することなく、平面的（二次元的）に配置されるので、第１磁気センサーにおける第２磁場の近似値を二次元的に高精度に計算することができる。
【００５４】
［適用例１６］
本適用例にかかる磁気計測システムは、第１磁場と第２磁場とを計測する第１磁気センサーと、前記第１磁気センサーの周囲に分布する第１の第２磁気センサーと第２の第２磁気センサーと第３の第２磁気センサーと第４の第２磁気センサーと、前記第１の第２磁気センサーの計測値と前記第２の第２磁気センサーの計測値と前記第３の第２磁気センサーの計測値と前記第４の第２磁気センサーの計測値とを用いて前記第１磁気センサーにおける前記第２磁場の近似値を計算する処理装置と、を含み、前記第１磁気センサーは、前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、の交差部を含む位置に配置されたことを特徴とする。
【００５５】
本適用例の構成によれば、第１磁場と第２磁場とを計測する第１磁気センサーが、第１の第２磁気センサーと第２の第２磁気センサーとを結ぶ線分と、第３の第２磁気センサーと第４の第２磁気センサーとを結ぶ線分との交差部を含む位置に配置されている。すなわち、第１の第２磁気センサーと第２の第２磁気センサーと、および、第３の第２磁気センサーと第４の第２磁気センサーとは、直線上に位置することなく平面的（二次元的）に配置されている。そのため、処理装置が第１磁気センサーにおける第２磁場の近似値を二次元的に高精度に計算することができる。
【００５６】
［適用例１７］
上記適用例に係る磁気計測システムであって、前記第１の第２磁気センサーと前記第２の第２磁気センサーとを結ぶ線分と、前記第３の第２磁気センサーと前記第４の第２磁気センサーとを結ぶ線分と、は直交することが好ましい。
【００５７】
本適用例の構成によれば、第１の第２磁気センサーと第２の第２磁気センサーとを結ぶ線分に沿った方向の第１磁気センサーにおける第２磁場の近似値は、第１の第２磁気センサーの計測値と第２の第２磁気センサーの計測値とから高精度に計算することができる。また、第３の第２磁気センサーと第２の第４磁気センサーとを結ぶ線分に沿った方向の第１磁気センサーにおける第２磁場の近似値は、第３の第２磁気センサーの計測値と第４の第２磁気センサーの計測値とから高精度に計算することができる。さらに、第１の第２磁気センサーと第２の第２磁気センサーとを結ぶ線分と、第３の第２磁気センサーと第４の第２磁気センサーとを結ぶ線分とが直交しているので、第１磁気センサーにおける第２磁場の近似値を二次元的により高精度に計算することができる。
【図面の簡単な説明】
【００５８】
【図１】第１の実施形態に係る磁気計測システムの構成の一例を示す概略側面図。
【図２】第１の実施形態に係る心磁センサーの構造を示す模式図。
【図３】第１の実施形態に係る心磁場計測処理の概略を示すフローチャート。
【図４】実施例１−１に係るノイズ磁気センサーの配置を説明する図。
【図５】実施例１−１に係るノイズ磁気センサーの配置を説明する図。
【図６】実施例１−３に係るノイズ磁気センサーの配置を説明する図。
【図７】実施例１−３に係るノイズ磁気センサーの配置を説明する図。
【図８】実施例３−１に係るノイズ磁気センサーの配置を説明する図。
【図９】実施例３−１に係るノイズ磁気センサーの配置を説明する図。
【発明を実施するための形態】
【００５９】
以下、実施形態について図面に従って説明する。
尚、各図面における各部材は、各図面上で認識可能な程度の大きさとするため、各部材毎に縮尺を異ならせて図示している。
【００６０】
（第１の実施形態）
＜磁気計測システム＞
まず、第１の実施形態に係る磁気計測システムの構成例を説明する。図１は、第１の実施形態に係る磁気計測システムの構成の一例を示す概略側面図である。図１に示す磁気計測システム１は、計測対象物としての被検体（生体）９の心臓から発せられる心磁場を計測するシステムである。図１に示すように、磁気計測システム１は、第１磁気センサーとしての心磁センサー１０と、第２磁気センサーとしてのノイズ磁気センサー３０と、処理装置としての磁気計測装置２とを備えている。
【００６１】
心磁センサー１０は、心磁や脳磁等の微弱な第１磁場と、外部磁場（磁気ノイズ）等の第２磁場とを計測するセンサーであり、心磁計や脳磁計等として使用される。ノイズ磁気センサー３０は、外部磁場（磁気ノイズ）等の第２磁場を計測するセンサーである。心磁センサー１０およびノイズ磁気センサー３０としては、光ポンピング式磁気センサー、ＳＱＵＩＤ式磁気センサー、フラックスゲート磁気センサー、ＭＩセンサー、ホール素子等を用いることができる。
【００６２】
磁気計測装置２は、土台３と、テーブル４と、磁気シールド装置６とを備えている。磁気計測装置２の高さ方向（図１における上下方向）をＺ方向とする。Ｚ方向は鉛直方向である。土台３、テーブル４の上面が延在する方向をＸ方向およびＹ方向とする。Ｘ方向およびＹ方向は水平方向であり、Ｘ方向とＹ方向とは直交する方向である。横たわった状態の被検体９の身長方向（図１における左右方向）をＸ方向とする。
【００６３】
土台３は磁気シールド装置６（本体部６ａ）の内側の底面上に配置され、本体部６ａの外側にまで、Ｘ方向（被検体９の移動可能方向）に沿って延在している。テーブル４は、Ｘ方向テーブル４ａと、Ｚ方向テーブル４ｂと、Ｙ方向テーブル４ｃとを有している。土台３上には、Ｘ方向直動機構３ａによりＸ方向に沿って移動するＸ方向テーブル４ａが設置されている。Ｘ方向テーブル４ａの上には、図示しない昇降装置によりＺ方向に沿って昇降するＺ方向テーブル４ｂが設置されている。Ｚ方向テーブル４ｂの上には、図示しないＹ方向直動機構によりレール上をＹ方向に沿って移動するＹ方向テーブル４ｃが設置されている。
【００６４】
磁気シールド装置６は、開口部６ｃを有する角筒状の本体部６ａを備えている。本体部６ａの内部は空洞となっており、Ｙ方向およびＺ方向を通る面（Ｙ−Ｚ断面でＸ方向に直交した平面）の断面形状は概ね四角形になっている。心磁場を計測する際は、本体部６ａの内部に被検体９がテーブル４上に横たわった状態で収容される。本体部６ａはＸ方向に延在しており、これ自体でパッシブ磁気シールドとして機能する。
【００６５】
心磁センサー１０およびノイズ磁気センサー３０は、磁気シールド装置６の本体部６ａの内部に配置されている。磁気シールド装置６は、地磁気等の外部磁場が、心磁センサー１０が配置された空間へ流入する事態を抑制している。すなわち、磁気シールド装置６により、心磁センサー１０が配置された空間は外部磁場に比べて著しく低磁場とされ、外部磁場の心磁センサー１０への影響が抑制されている。
【００６６】
本体部６ａの開口部６ｃから＋Ｘ方向に土台３が突出している。磁気シールド装置６の大きさは、例えば、Ｘ方向の長さが約２００ｃｍ程度であり、開口部６ｃの一辺が９０ｃｍ程度である。そして、開口部６ｃから、磁気シールド装置６内に、テーブル４に横たわった被検体９がテーブル４と共に土台３上をＸ方向に沿って移動して出入することができる。
【００６７】
図示しないが、磁気計測装置２は、電気信号を用いて磁気計測装置２を制御する制御部を備えている。この電気信号により磁場や残留磁場が発生して心磁センサー１０に検出されるときノイズとなる。制御部は、発生される磁場や残留する磁場が心磁センサー１０に到達し難くなるように、磁気シールド装置６の開口部６ｃから離れた場所に設置されている。
【００６８】
制御部は、表示装置および入力装置を有している。表示装置は、ＬＣＤやＯＬＥＤ等の表示装置で構成され、測定の状況や測定結果等を表示する。入力装置は、キーボードや回転つまみ等から構成され、操作者が操作することにより磁気計測装置２の測定開始指示や測定条件等の各種指示が入力される。
【００６９】
磁気シールド装置６の本体部６ａは、比透磁率が例えば数千以上の強磁性体、または、高伝導率の導体によって形成される。強磁性体にはパーマロイ、フェライト、または鉄、クロムもしくはコバルト系のアモルファス等を用いることができる。高伝導率の導体には、例えば、アルミニウム等で、渦電流効果によって磁場低減効果を有するものを用いることができる。なお、強磁性体と高伝導率の導体とを交互に積層して本体部６ａを形成することも可能である。
【００７０】
本体部６ａおよび土台３の＋Ｘ方向側および−Ｘ方向側の端には補正コイル（ヘルムホルツコイル）６ｂが設置されている。補正コイル６ｂの形状は枠状であり、本体部６ａを囲むように配置されている。補正コイル６ｂは、本体部６ａの内部空間へ流入する流入磁場を補正するためのコイルである。流入磁場は、外部磁場が開口部６ｃを通過して内部空間に入り込む磁場を指す。流入磁場は開口部６ｃに対してＸ方向で最も強くなる。補正コイル６ｂは、制御部から供給される電流により流入磁場をキャンセルするように磁界を発生させる。
【００７１】
心磁センサー１０は、本体部６ａの天井に支持部材７を介して固定されている。心磁センサー１０は、Ｚ方向における磁場の強度成分を計測する。被検体９の心磁場を計測する際は、被検体９における計測位置である胸部９ａが心磁センサー１０と対向する位置になるようにＸ方向テーブル４ａおよびＹ方向テーブル４ｃを移動させ、胸部９ａが心磁センサー１０に接近するようにＺ方向テーブル４ｂを上昇させる。
【００７２】
ノイズ磁気センサー３０は、心磁センサー１０の周囲に複数個（６個以上）配置されている。詳細は後述するが、磁気計測システム１は、ノイズ磁気センサー３０として、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２と第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４と第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６とを含む（図８（ａ）参照）。
【００７３】
ノイズ磁気センサー３０（３１，３２，３３，３４，３５，３６）のそれぞれは、Ｘ方向とＹ方向とＺ方向とにおける磁場の３成分を計測する。これにより、ノイズ磁気センサー３０周辺の磁場分布を特定することができる。ノイズ磁気センサー３０は、磁場分布を特定したい空間（以下では、計測対象空間と称する）、すなわち、被検体９が配置される空間を取り囲むように立体的に配置されることが好ましい。
【００７４】
磁気計測装置２の制御部は、心磁センサー１０で計測された第１磁場および第２磁場の計測値と、ノイズ磁気センサー３０で計測された第２磁場の計測値とから、第１磁場である被検体９の心磁場を算出する機能を有している。より具体的には、制御部は、ノイズ磁気センサー３０（３１，３２，３３，３４，３５，３６）の計測値を用いて、心磁センサー１０における第２磁場の近似値を計算し、心磁センサー１０の計測値から差し引いて第１磁場を算出する。
【００７５】
心磁センサー１０は、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２とを結ぶ線分と、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４とを結ぶ線分と、第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６とを結ぶ線分と、の交差部を含む位置に配置されることが好ましい。
【００７６】
このように、ノイズ磁気センサー３０の数が２ｎ個、または２ｎ＋１個（ｎは３以上の整数）の場合、ノイズ磁気センサー３０をｎ個の対とし、各対に関して２つのノイズ磁気センサー３０の重心を含む位置に心磁センサー１０が配置されることが、高精度で計測する観点から好ましい。言い換えれば、各対に関して２つのノイズ磁気センサー３０は心磁センサー１０に関して対称的に配置されることが好ましい。
【００７７】
さらに、心磁センサー１０は、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２とを結ぶ線分と、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４とを結ぶ線分と、第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６とを結ぶ線分と、のうちの少なくとも２本の線分は直交し、この２本の線分に平行な平面に対して、残りの１本の線分が交差する配置であることが好ましい。要するに、第１のノイズ磁気センサー３１から第６のノイズ磁気センサー３６までのノイズ磁気センサー３０が立体的に配置されることが好ましい。ノイズ磁気センサー３０の配置の詳細については、後述する。
【００７８】
＜心磁センサー＞
次に、心磁センサー１０の概略構造を説明する。図２は、第１の実施形態に係る心磁センサーの構造を示す模式図である。詳しくは、図２（ａ）は心磁センサーの模式側面図であり、図２（ｂ）は心磁センサーの模式平面図である。
【００７９】
図２（ｂ）に示すように、心磁センサー１０には、レーザー光源１８からレーザー光１８ａが供給される。レーザー光源１８は制御部に設置され、レーザー光源１８から発せられたレーザー光１８ａは光ファイバー１９を通って心磁センサー１０に供給される。心磁センサー１０と光ファイバー１９とは、光コネクター２０を介して接続されている。
【００８０】
レーザー光源１８は、セシウムの吸収線に応じた波長のレーザー光１８ａを出力する。レーザー光１８ａの波長は特に限定されないが、本実施形態では、例えば、Ｄ１線に相当する８９４ｎｍの波長に設定している。レーザー光源１８はチューナブルレーザーであり、レーザー光源１８から出力されるレーザー光１８ａは一定の光量を有する連続光である。
【００８１】
光コネクター２０を介して供給されたレーザー光１８ａは、−Ｙ方向に進行して偏光板２１に入射する。偏光板２１を通過したレーザー光１８ａは、直線偏光になっている。そして、レーザー光１８ａは、第１ハーフミラー２２、第２ハーフミラー２３、第３ハーフミラー２４、第１反射ミラー２５に順次入射する。
【００８２】
第１ハーフミラー２２、第２ハーフミラー２３および第３ハーフミラー２４は、レーザー光１８ａの一部を反射して＋Ｘ方向に進行させ、一部のレーザー光１８ａを通過させて−Ｙ方向に進行させる。第１反射ミラー２５は、入射したレーザー光１８ａを全て＋Ｘ方向に反射する。第１ハーフミラー２２、第２ハーフミラー２３、第３ハーフミラー２４、第１反射ミラー２５により、レーザー光１８ａは４つの光路に分割される。各光路のレーザー光１８ａの光強度が同じ光強度になるように、各ミラーの反射率が設定されている。
【００８３】
次に、図２（ａ）に示すように、レーザー光１８ａは第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９に順次照射入射する。第４ハーフミラー２６、第５ハーフミラー２７および第６ハーフミラー２８は、レーザー光１８ａの一部を反射して＋Ｚ方向に進行させ、一部のレーザー光１８ａを通過させて＋Ｘ方向に進行させる。第２反射ミラー２９は、入射したレーザー光１８ａを全て＋Ｚ方向に反射する。
【００８４】
第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９により、１つの光路のレーザー光１８ａは４つの光路に分割される。各光路のレーザー光１８ａの光強度が同じ光強度になるように、各ミラーの反射率が設定されている。したがって、レーザー光１８ａは１６個の光路に分離される。そして、各光路のレーザー光１８ａの光強度が同じ強度になるように、各ミラーの反射率が設定されている。
【００８５】
第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９の＋Ｚ方向側には、レーザー光１８ａの各光路に、４行４列の１６個のガスセル１２が設置されている。そして、第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９にて反射したレーザー光１８ａは、ガスセル１２を通過する。ガスセル１２は、内部に空隙を有する箱であり、この空隙にはアルカリ金属のガスが封入されている。アルカリ金属は特に限定されず、カリウム、ルビジウムまたはセシウムを用いることができる。本実施形態では、例えばアルカリ金属にセシウムを用いている。
【００８６】
各ガスセル１２の＋Ｚ方向側には、偏光分離器１３が設置されている。偏光分離器１３は、入射したレーザー光１８ａを、互いに直交する２つの偏光成分のレーザー光１８ａに分離する素子である。偏光分離器１３には、例えば、ウォラストンプリズムまたは偏光ビームスプリッターを用いることができる。
【００８７】
偏光分離器１３の＋Ｚ方向側には第１光検出器１４が設置され、偏光分離器１３の＋Ｘ方向側には第２光検出器１５が設置されている。偏光分離器１３を通過したレーザー光１８ａは第１光検出器１４に入射し、偏光分離器１３にて反射したレーザー光１８ａは第２光検出器１５に入射する。第１光検出器１４および第２光検出器１５は、入射したレーザー光１８ａの光量に応じた電流を制御部に出力する。
【００８８】
第１光検出器１４および第２光検出器１５が磁場を発生すると測定に影響を与える可能性があるので、第１光検出器１４および第２光検出器１５は非磁性の材料で構成されることが望ましい。心磁センサー１０は、Ｘ方向の両面およびＹ方向の両面に設置されたヒーター１６を有している。ヒーター１６は磁界を発生しない構造であることが好ましく、例えば、流路中に蒸気や熱風を通過させて加熱する方式のヒーターを用いることができる。ヒーターの代わりに、高周波電圧によりガスセル１２を誘電加熱してもよい。
【００８９】
心磁センサー１０は、被検体９（図１参照）の＋Ｚ側に配置される。被検体９が発する磁気ベクトルＢは、−Ｚ方向側から心磁センサー１０に入る。磁気ベクトルＢは、第４ハーフミラー２６〜第２反射ミラー２９を通過し、ガスセル１２を通過した後、偏光分離器１３を通過して心磁センサー１０から出る。
【００９０】
心磁センサー１０は、光ポンピング式磁気センサーや光ポンピング原子磁気センサーと称されるセンサーである。ガスセル１２内のセシウムは、加熱されてガス状態になっている。そして、直線偏光になったレーザー光１８ａをセシウムガスに照射することにより、セシウム原子が励起され磁気モーメントの向きが揃えられる。この状態でガスセル１２に磁気ベクトルＢが通過するとき、セシウム原子の磁気モーメントが磁気ベクトルＢの磁場により歳差運動する。この歳差運動をラーモア歳差運動と称する。
【００９１】
ラーモア歳差運動の大きさは、磁気ベクトルＢの強さと正の相関を有している。ラーモア歳差運動は、レーザー光１８ａの偏向面を回転させる。ラーモア歳差運動の大きさとレーザー光１８ａの偏向面の回転角の変化量とは、正の相関を有する。したがって、磁気ベクトルＢの強さとレーザー光１８ａの偏向面の回転角の変化量とは、正の相関を有している。心磁センサー１０の感度は、磁気ベクトルＢのＺ方向において高く、Ｚ方向と直交する方向において低くなっている。
【００９２】
偏光分離器１３は、レーザー光１８ａを直交する２成分の直線偏光に分離する。そして、第１光検出器１４および第２光検出器１５は、直交する２成分の直線偏光の強さを検出する。これにより、第１光検出器１４および第２光検出器１５は、レーザー光１８ａの偏向面の回転角を検出することができる。そして、レーザー光１８ａの偏向面の回転角の変化から、心磁センサー１０は磁気ベクトルＢの強さを検出することができる。
【００９３】
ガスセル１２、偏光分離器１３、第１光検出器１４、および第２光検出器１５からなる素子をセンサー素子１１と称する。心磁センサー１０には、センサー素子１１が４行４列の１６個配置されている。心磁センサー１０におけるセンサー素子１１の個数および配置は特に限定されない。センサー素子１１は、３行以下でもよく５行以上でもよい。同様にセンサー素子１１は、３列以下でもよく５列以上でもよい。センサー素子１１の個数が多い程空間分解能を高くすることができる。
【００９４】
心磁センサー１０が配置される計測対象空間は、磁気シールド装置６（図１参照）により外部磁場の流入が抑制されているが、外部磁場の流入を皆無とすることは困難である。換言すれば、心磁センサー１０には、心磁場と外部磁場（磁気ノイズ）とが印加される。そのため、心磁センサー１０で計測して得られる計測値には、心磁場に基づく信号成分と外部磁場に基づくノイズ成分とが含まれる。したがって、被検体９の心磁場を正確に取得するためには、心磁センサー１０で得られる計測値からノイズ成分を高精度に除去する必要がある。
【００９５】
＜ノイズ磁気センサー＞
図１に戻り、ノイズ磁気センサー３０は、心磁センサー１０が配置される計測対象空間における外部磁場（磁気ノイズ）を計測するためのものである。ノイズ磁気センサー３０で得られる計測値から計測対象空間における外部磁場を特定することにより、心磁センサー１０で得られる計測値から外部磁場（磁気ノイズ）成分を除去することができる。ノイズ磁気センサー３０は、外部磁場（磁気ノイズ）等の第２磁場を検知し、第１磁場を検知しないものとする。なお、ノイズ磁気センサー３０が高感度であれば、ノイズ磁気センサー３０の計測値から第１磁場と第２磁場との合成磁場を計測することもできる。
【００９６】
ノイズ磁気センサー３０として用いるセンサーの種類は限定されないが、例えば、上述の心磁センサー１０と同様の光ポンピング式磁気センサーを用いることができる。光ポンピング式磁気センサーを用いる場合、例えば、１箇所のノイズ磁気センサー３０として、図２（ａ），（ｂ）に示すセンサー素子１１を３個組み合わせて用いてもよい。この場合、３個のセンサー素子１１で、Ｘ方向、Ｙ方向、およびＺ方向の各方向に沿った磁気ベクトルを計測する。また、１箇所のノイズ磁気センサー３０として１個のセンサー素子１１を用い、レーザー光１８ａをＸ方向、Ｙ方向、およびＺ方向の各方向から順次照射して、各方向に沿った磁気ベクトルを時系列的に計測することとしてもよい。
【００９７】
＜心磁場計測処理の概略＞
磁気計測装置２の制御部による、心磁場計測処理の概略を説明する。図３は、第１の実施形態に係る心磁場計測処理の概略を示すフローチャートである。図３の左側に心磁センサー１０の計測値に係る処理フローを示し、図３の右側にノイズ磁気センサー３０の計測値に係る処理フローを示している。
【００９８】
ステップＳ１１では、制御部は、心磁センサー１０の計測値を取得する。心磁センサー１０の計測値は、微弱な心磁場（第１磁場）と磁気ノイズ等の外部磁場（第２磁場）とを含んでいる。また、ステップＳ２１では、制御部は、ノイズ磁気センサー３０の計測値を取得する。ノイズ磁気センサー３０の計測値は、心磁センサー１０の位置を含む計測対象空間における磁気ノイズ等の外部磁場（第２磁場）を含んでいる。
【００９９】
ステップＳ１１の心磁センサー１０の計測値取得と、ステップＳ２１のノイズ磁気センサー３０の計測値取得とは、並行して行ってもよいし、別々に行ってもよい。なお、ノイズ磁気センサー３０が外部磁場だけでなく心磁場も検出する場合、すなわち、ノイズ磁気センサー３０の計測値に心磁場と外部磁場とが含まれる場合は、ステップＳ１１の心磁センサー１０の計測値取得を行う前に、被検体９がいない状態で、ステップＳ２１のノイズ磁気センサー３０の計測値取得を行う必要がある。
【０１００】
ステップＳ２２では、制御部は、ノイズ磁気センサー３０の計測値を、磁場を表す関数に適応させる。ステップＳ２２では、計測対象空間における外部磁場の分布を高精度に近似できる関数を用いることが好ましい。続くステップＳ２３では、制御部は、計測対象空間における外部磁場の近似値を計算する。そして、ステップＳ２４では、制御部は、心磁センサー１０の位置（計測位置とも称する）における外部磁場の近似値（近似値ベクトルＡとも称する）を算出する。ステップＳ２２〜ステップＳ２４において、ステップＳ２１で取得したノイズ磁気センサー３０の計測値を関数に適応して外部磁場の近似値を計算する方法については後述する。
【０１０１】
次に、ステップＳ１２では、制御部は、ステップＳ１１で取得した心磁センサー１０の計測値から、ステップＳ２４で算出した心磁センサー１０の位置における外部磁場の近似値を引き算する。具体的には、心磁センサー１０が磁場ベクトルを計測する場合、近似値のベクトルをベクトル的に減じてもよいし、特定成分として減じてもよい。本実施形態に係る心磁センサー１０は、直線偏光の進行方向の成分（本実施形態の例ではＺ成分）を高感度に計測する。こうした場合、近似値の当該成分を減じてもよい。
【０１０２】
一般的には、心磁センサー１０の計測方向を第４ベクトルｄとした場合、心磁センサー１０は、心磁センサー１０が位置する場所（計測位置）での磁場Ｂ0と第４ベクトルｄとの内積値Ｂ0・ｄを計測するので、計測位置における近似値ベクトルＡと第４ベクトルｄとの内積値Ａ・ｄを心磁センサー１０の計測値（内積値Ｂ0・ｄ）から減じてもよい。この結果、心磁センサー１０の計測値に含まれるノイズ成分である外部磁場（第２磁場）が除去されるので、信号成分である心磁場（第１磁場）の計測値が得られる。
【０１０３】
なお、本実施形態に係る磁気計測システム１では、図３に示す処理フローを繰り返すことにより、時間的に変動する磁場を計測することが可能である。したがって、高速で演算処理することが可能な（例えば、１００Ｈｚの時間分解能を有する）構成であれば、心磁波形をほぼリアルタイムで取得することができる。高速で演算処理することが困難な場合は、取得した計測値を時系列で記憶しておき、計測が終了した後で記憶しておいた計測値の演算処理を行うこととしてもよい。
【０１０４】
＜外部磁場の近似値の計算方法＞
次に、外部磁場の近似値を計算する方法について説明する。磁気計測装置２の制御部は、上述のステップＳ２１で取得したノイズ磁気センサー３０の計測値を多変数多項式に適応して、計測対象空間における外部磁場の近似値を計算する。本実施形態では、第１磁場である心磁場が、第２磁場である外部磁場に対して微弱である。
【０１０５】
計測対象空間の任意の位置を、数式１８に示すように、位置ベクトルｒ（以下では、位置ｒとも称する）にて表す。
【０１０６】
【数１８】
2016217930.tif 000020
【０１０７】
また、位置ｒにおける磁場を、数式１９に示すように、磁場ベクトルＢ（以下では、磁場Ｂとも称する）にて表す。
【０１０８】
【数１９】
2016217930.tif 000021
【０１０９】
本願発明者が鋭意研究したところによると、磁場Ｂは数式２０に示される多変数多項式（３変数で、変数に関する２次式の多項式）にて精度良く近似されることが判明した。数式２０において、ａij（ｉは１から３の整数、ｊは１から７の整数）は係数であり、ｘ、ｙ、ｚは磁場の近似値Ｂの空間座標であり、Ｂiは磁場の近似値Ｂの第ｉ成分である。なお、数式２０は、非線形多項式（３変数で、変数に関する１次項と２次項とを含む多項式）でもある。
【０１１０】
【数２０】
2016217930.tif 000022
【０１１１】
数式２０において、ｉ＝１が磁場ＢのＸ成分Ｂxを表し、ｉ＝２が磁場ＢのＹ成分Ｂyを表し、ｉ＝３が磁場ＢのＺ成分Ｂzを表す。具体的には、磁場Ｂの各成分は数式２１で表される。
【０１１２】
【数２１】
2016217930.tif 000023
【０１１３】
数式２０に示される多変数多項式の右辺の第１項（ａi1）は、計測対象空間全体での平行磁場（原点からのずれ、オフセット磁場）を表す。数式２０に示される多変数多項式の右辺の第２項（ａi2ｘ）から第４項（ａi4ｙ）は、線形磁場（磁場の傾き、傾斜磁場）を表す。数式２０に示される多変数多項式の右辺の第５項（ａi5ｘｙ）から第７項（ａi7ｚｘ）は交互磁場（磁場のねじれ、回転磁場）を表す。
【０１１４】
磁場Ｂには、ビオ・サバールの法則にて、電流要素ベクトルと位置ベクトルとの外積項に比例する成分が存在する。そこで、本願発明者は、磁場Ｂの近似式にねじれ成分を表すｘｙ項とｙｚ項とｚｘ項とを導入した。このような理由に基づき、数式２０に示される多変数多項式により、計測対象空間の任意の位置における磁場Ｂが高精度に近似される訳である。
【０１１５】
数式２０に示される多変数多項式では、ＸＹＺの３成分のそれぞれに対して７つの未知数ａij（ｊは１から７の整数）が存在するので、合計３×７＝２１個の未知数が存在することになる。これら２１個の未知数は計測対象空間に固有の値となり、計測対象空間に対応するこれらの未知数を特定すれば、計測対象空間での磁場Ｂを高精度に近似できることになる。続いて、これらの未知数ａijを特定する方法を説明する。
【０１１６】
まず、ノイズ磁気センサー項ベクトルＲ（以下では、磁気センサー項ベクトルＲとも称する）を数式２２と定義する。
【０１１７】
【数２２】
2016217930.tif 000024
【０１１８】
数式２２の磁気センサー項ベクトルＲを用いると、数式２１の磁場ベクトルＢは数式２３のように表される。
【０１１９】
【数２３】
2016217930.tif 000025
【０１２０】
数式２３の最後の等号は、未知数行列ａ（以下では、第１行列ａとも称する）の定義である。同様にして、数式２０で表される行列要素は数式２４にて表される。
【０１２１】
【数２４】
2016217930.tif 000026
【０１２２】
ノイズ磁気センサー３０の数をαで表す。上述したように、未知数ａijが２１個存在するので、本実施形態ではαは７以上の整数であり、一例としては、α＝８である。各ノイズ磁気センサー３０が磁場ＢのＸＹＺの３成分を計測するので、αが７以上であれば少なくとも２１個の未知数ａijを特定することができる。
【０１２３】
ｋ番目のノイズ磁気センサー３０の位置を、数式２５に示すように、ノイズ磁気センサー位置ベクトルｒk（以下では、磁気センサー位置ｒkとも称する）にて表す。ここで、ｋは１からαの整数であり、本実施形態の例では、ノイズ磁気センサー３０の数がα＝８個であるので、ｋは１から８の整数である。
【０１２４】
【数２５】
2016217930.tif 000027
【０１２５】
また、ｋ番目のノイズ磁気センサー３０の位置（磁気センサー位置ｒk）における磁場Ｂを、数式２６に示すように、ｋ番目の検出磁場ベクトルＢk（以下では、検出磁場Ｂkとも称する）にて表す。
【０１２６】
【数２６】
2016217930.tif 000028
【０１２７】
数式２６に数式２３を適応させると、数式２７が得られる。
【０１２８】
【数２７】
2016217930.tif 000029
【０１２９】
数式２７において、最後から２番目の等号は、磁気センサー位置ｒkにおける磁気センサー項ベクトルＲkの定義である。このとき、数式２４で表され、ｋ番目のノイズ磁気センサー３０が検出する検出磁場Ｂkのｉ行目の成分を表す行列要素Ｂikは、数式２８にて表される。
【０１３０】
【数２８】
2016217930.tif 000030
【０１３１】
次に、数式２７または数式２８を用いて、α個の検出磁場ベクトルＢk全体で形成する検出磁場行列Ｍ（第２行列Ｍとも称する）と、α個の磁気センサー項ベクトルＲk全体で形成する磁気センサー項行列Ｐ（第３行列Ｐとも称する）とを、それぞれ数式２９と数式３０とで表す。
【０１３２】
【数２９】
2016217930.tif 000031
【０１３３】
【数３０】
2016217930.tif 000032
【０１３４】
数式２９に示されているように、検出磁場行列Ｍは３行α列の行列であり、ｉ行ｋ列の行列要素ＭikはＭik＝Ｂikと表現される。すなわち、検出磁場行列Ｍは、３行１列の検出磁場ベクトルＢkをｋ＝１からαまでのα列並べた行列である。したがって、例えば、Ｂikはｋ番目のノイズ磁気センサー３０が検出する磁場のｉ行目の成分である。上述したように、ｉ＝１がＸ成分であり、ｉ＝２がＹ成分であり、ｉ＝３がＺ成分である。
【０１３５】
また、数式３０に示されているように、磁気センサー項行列Ｐは７行α列の行列であり、ｇ行ｋ列の行列要素ＰgkはＰgk＝Ｒgkと表現される。すなわち、磁気センサー項行列Ｐは７行１列の磁気センサー項ベクトルＲkをｋ＝１からαまで、α列並べた行列である。
【０１３６】
したがって、Ｒgkは、ｋ番目のノイズ磁気センサー３０（磁気センサー位置ｒk）における磁気センサー項ベクトルＲkのｇ行目の成分である。例えば、ｇ＝２のときＲgkはｘkであり、ｇ＝７のときＲgkはｚkｘkである。これらの検出磁場行列Ｍ（第２行列Ｍ）と磁気センサー項行列Ｐ（第３行列Ｐ）とは、未知数行列ａ（第１行列ａ）を用いて数式３１の関係で表現される。
【０１３７】
【数３１】
2016217930.tif 000033
【０１３８】
数式３１を行列要素にて表現すると数式３２となる。
【０１３９】
【数３２】
2016217930.tif 000034
【０１４０】
上述したように、ｉは１から３の整数であり、ノイズ磁気センサー３０を表すｋは１からαの整数である。したがって、例えばｋ番目のノイズ磁気センサー３０が検出する磁場ＢのＹ成分（ｉ＝２）であるＢ2kは、数式３３にて表現される。
【０１４１】
【数３３】
2016217930.tif 000035
【０１４２】
数式３１で求めるべきは、未知数行列ａ（第１行列ａ）である。ノイズ磁気センサー３０の数αが７であれば、磁気センサー項行列Ｐ（第３行列Ｐ）は７行７列の正方行列となり、その逆行列が存在する。この場合は、未知数行列ａ（第１行列ａ）は、数式３４に示すように、磁気センサー項行列Ｐ（第３行列Ｐ）の逆行列Ｐ-1を検出磁場行列Ｍ（第２行列Ｍ）に右から掛けて求められる。
【０１４３】
【数３４】
2016217930.tif 000036
【０１４４】
一方、ノイズ磁気センサー３０の数αが８以上であれば、磁気センサー項行列Ｐ（第３行列Ｐ）は正方行列にはならないので、その逆行列は存在しない。この場合は、未知数行列ａ（第１行列ａ）は、数式３５に示すように、磁気センサー項行列Ｐ（第３行列Ｐ）の疑似逆行列（一般化逆行列ともいう）Ｐ+を検出磁場行列Ｍ（第２行列Ｍ）に右から掛けて求められる。
【０１４５】
【数３５】
2016217930.tif 000037
【０１４６】
数式３５において、磁気センサー項行列Ｐ（第３行列Ｐ）の疑似逆行列Ｐ+は、数式３６で求められる。
【０１４７】
【数３６】
2016217930.tif 000038
【０１４８】
数式３６に示されるように、疑似逆行列Ｐ+は、磁気センサー項行列Ｐの転置行列ＰTと磁気センサー項行列Ｐとの積の逆行列に、磁気センサー項行列Ｐの転置行列ＰTを掛けたものである。なお、磁気センサー項行列Ｐの転置行列ＰTは、磁気センサー項行列Ｐの行列要素を行と列とに関して入れ替えたもので、α行７列の行列であり、数式３７にて表される。
【０１４９】
【数３７】
2016217930.tif 000039
【０１５０】
このように、疑似逆行列Ｐ+を用いると、最小二乗法の原理が働き、未知数行列ａ（第１行列ａ）に対して誤差を最小とする最適解が定まる。
【０１５１】
以上により、数式２０に示される多変数多項式における２１個の未知数が特定され、計測対象空間の任意の位置における磁場Ｂを高精度に近似できる。この結果、心磁センサー１０の位置における外部磁場の近似値（近似値ベクトルＡ）を高精度で算出することができる。
【０１５２】
具体的には、数式３４または３５にて求めた未知数行列ａを数式２３に適応する。その際に、磁気センサー項ベクトルＲには心磁センサー１０の座標を適応させる。その結果、数式２３の磁場ベクトルＢは、心磁センサー１０の位置における近似値ベクトルＡを表すことになる。そして、算出された外部磁場の近似値（近似値ベクトルＡ）を心磁センサー１０の計測値から除去することにより、計測対象とする心磁場を低ノイズで計測することができる。すなわち、心磁計測のような微弱信号を高Ｓ／Ｎ比（Signal to Noise Ratio）で取り出すことが可能になる。
【０１５３】
このように、本実施形態に係る磁気計測システム１によれば、計測対象とする磁場（例えば、心磁場）に対して大きな外部磁場（例えば、磁気ノイズ）が存在する場合でも、計測対象とする位置における外部磁場を演算により高精度に近似して計測値から除去できる。したがって、心磁場等の微弱な磁場を精度良く計測することができる。
【０１５４】
また、図２（ｂ）に示すように、心磁センサー１０がマトリックス状に配置された複数のセンサー素子１１で構成されていると、各センサー素子１１においてガスセル１２に入射するレーザー光１８ａの光軸がばらついてしまう場合がある。本実施形態に係る磁気計測システム１によれば、各センサー素子１１においてガスセル１２に入射するレーザー光１８ａの光軸にばらつきがあっても、各センサー素子１１に対して個別に外部磁場を高精度に近似して除去できるので、微弱な磁場の計測精度を高く維持することができる。
【０１５５】
さらに、本実施形態に係る磁気計測システム１によれば、心磁場計測処理を行列演算を用いて実行するので、磁気計測装置２の制御部をゲートアレイ等の簡易な素子で構成できる。これにより、磁気計測システム１を、より低コストで実現することができる。
【０１５６】
次に、ノイズ磁気センサー３０の配置の実施例を挙げて、第１の実施形態に係る外部磁場の近似値計算における磁気センサー項行列Ｐ（第３行列Ｐ）の計算方法を説明する。
【０１５７】
（実施例１−１）
図４および図５は、実施例１−１に係るノイズ磁気センサーの配置を説明する図である。詳しくは、図４（ａ）は斜視図であり、図４（ｂ）は図４（ａ）の＋Ｘ方向側から見た平面図である。図５（ａ）は図４（ａ）の＋Ｙ方向側から見た平面図であり、図５（ｂ）は図４（ａ）の＋Ｚ方向側から見た平面図である。
【０１５８】
図４および図５において、８個のノイズ磁気センサー３０を個別に識別するため、第１番目から第８番目までそれぞれノイズ磁気センサー３１，３２，３３，３４，３５，３６，３７，３８と称する。また、これらのノイズ磁気センサー３１，３２，３３，３４，３５，３６，３７，３８を総称する場合にノイズ磁気センサー３０という。
【０１５９】
磁気計測システム１において、ノイズ磁気センサー３１，３３，３６，３８の４個のセンサーは、例えば本体部６ａ（図１参照）等に取り付けられ、計測対象空間の上部におけるＸ−Ｙ平面に平行な平面の４隅に配置されている。ノイズ磁気センサー３２，３４，３５，３７の４個のセンサーは、例えばＹ方向テーブル４ｃ（図１参照）等に取り付けられ、計測対象空間の下部におけるＸ−Ｙ平面に平行であり、かつ、ノイズ磁気センサー３１，３３，３６，３８が配置される平面と平面視で重なる平面の４隅に配置されている。
【０１６０】
したがって、図４（ａ）に示すように、８個のノイズ磁気センサー３０（３１，３２，３３，３４，３５，３６，３７，３８）は、直方体３０ａの各頂点に配置されている。なお、直方体３０ａは、Ｘ−Ｙ平面に平行な２つの平面と、Ｙ−Ｚ平面に平行な２つの平面と、Ｘ−Ｚ平面に平行な２つの平面とで構成される六面体である。
【０１６１】
また、８個のノイズ磁気センサー３０は、直方体３０ａの中心３０ｃと心磁センサー１０の中心１０ｃとがほぼ一致するように配置される。本実施例では、直方体３０ａの中心３０ｃと心磁センサー１０の中心１０ｃとが一致し、図４（ａ），（ｂ）におけるＸＹＺ直交座標系の原点に直方体３０ａの中心３０ｃおよび心磁センサー１０の中心１０ｃが配置されるものとする。
【０１６２】
このように、２ｎ個（ｎは３以上の整数）のノイズ磁気センサー３０をｎ個の対とし、各対に関して、各対をなすノイズ磁気センサー３０の重心位置に心磁センサー１０を配置する。こうすると、各対からの計測値が心磁センサー１０の位置での近似値ベクトルＡに対して同じ重要度で貢献することになる。さらに、各対の長さを全て揃えると、２ｎ個（ｎは３以上の整数）のノイズ磁気センサー３０からの計測値の重要度が全て同じになるので、より高精度に心磁センサー１０の位置での近似値ベクトルＡを計測することができる。
【０１６３】
直方体３０ａのＸ軸に沿った辺の長さを２Ｌ1とし、直方体３０ａのＹ軸に沿った辺の長さを２Ｌ2とし、直方体３０ａのＺ軸に沿った辺の長さを２Ｌ3とする。図４（ａ），（ｂ）に示すように、第１のノイズ磁気センサー３１と第４のノイズ磁気センサー３４と第５のノイズ磁気センサー３５と第８のノイズ磁気センサー３８とは、Ｙ−Ｚ平面に平行なＸ＝Ｌ1の平面に配置される。また、第２のノイズ磁気センサー３２と第３のノイズ磁気センサー３３と第６のノイズ磁気センサー３６と第７のノイズ磁気センサー３７とは、Ｙ−Ｚ平面に平行なＸ＝−Ｌ1の平面に配置される。
【０１６４】
図４（ａ）および図５（ａ）に示すように、第１のノイズ磁気センサー３１と第３のノイズ磁気センサー３３と第５のノイズ磁気センサー３５と第７のノイズ磁気センサー３７とは、Ｘ−Ｚ平面に平行なＹ＝Ｌ2の平面に配置される。また、第２のノイズ磁気センサー３２と第４のノイズ磁気センサー３４と第６のノイズ磁気センサー３６と第８のノイズ磁気センサー３８とは、Ｘ−Ｚ平面に平行なＹ＝−Ｌ2の平面に配置される。
【０１６５】
図４（ａ）および図５（ｂ）に示すように、第１のノイズ磁気センサー３１と第３のノイズ磁気センサー３３と第６のノイズ磁気センサー３６と第８のノイズ磁気センサー３８とは、Ｘ−Ｙ平面に平行なＺ＝Ｌ3の平面に配置される。また、第２のノイズ磁気センサー３２と第４のノイズ磁気センサー３４と第５のノイズ磁気センサー３５と第７のノイズ磁気センサー３７とは、Ｘ−Ｙ平面に平行なＺ＝−Ｌ3の平面に配置される。
【０１６６】
８個のノイズ磁気センサー３０を配置する際に、最も好ましい直方体３０ａは、Ｌ1／２1/2＝Ｌ2＝Ｌ3＝Ｌとなる直方体３０ａである。すなわち、一辺の長さが２Ｌの正方形（Ｙ−Ｚ断面）を底面とし、高さ（Ｘ軸に沿った辺の長さ）が２×２1/2×Ｌの四角柱である。このような直方体３０ａでは、底面の法線はＸ軸と平行となり、二つの側面の法線はＹ軸に平行となり、残りの２つの側面の法線はＺ軸に平行となる。このような直方体３０ａの各頂点に配置されたノイズ磁気センサー３０の位置ベクトルｒk（磁気センサー位置ｒk）は、数式３８で表される。
【０１６７】
【数３８】
2016217930.tif 000040
【０１６８】
８個のノイズ磁気センサー３０をこのように配置すると、磁気センサー項行列Ｐ（第３行列Ｐ）は数式３９で示されるので、計算が比較的簡便となる。また、この場合、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２とを結ぶ線分と、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４とを結ぶ線分とが直交するので、原点付近（心磁センサー１０が位置する付近）の第２磁場の近似値を高精度に計算することができる。同様に、第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６とを結ぶ線分と、第７のノイズ磁気センサー３７と第８のノイズ磁気センサー３８とを結ぶ線分とが直交する。
【０１６９】
【数３９】
2016217930.tif 000041
【０１７０】
さらに、ＸＹＺ座標系の１単位をＬとすると、磁気センサー項行列Ｐ（第３行列Ｐ）は数式４０に示すように１と−１とで表現されるので、計算がより簡便となる。
【０１７１】
【数４０】
2016217930.tif 000042
【０１７２】
（実施例１−２）
実施例１−２では、８個のノイズ磁気センサー３０と心磁センサー１０との位置関係は実施例１−１と同じであるため図示を省略するが、最も好ましい直方体３０ａは、Ｌ1＝Ｌ2＝Ｌ3＝Ｌとなる直方体３０ａである。すなわち、直方体３０ａは一辺の長さが２Ｌの立方体である。このような一辺の長さが２Ｌの立方体の各頂点に配置されたノイズ磁気センサー３０の位置ベクトルｒk（磁気センサー位置ｒk）は、数式４１で表される。
【０１７３】
【数４１】
2016217930.tif 000043
【０１７４】
８個のノイズ磁気センサー３０をこのように配置すると、磁気センサー項行列Ｐ（第３行列Ｐ）は、数式４２に示すように１とＬと−Ｌとで表現されるので、計算が簡便となる。
【０１７５】
【数４２】
2016217930.tif 000044
【０１７６】
さらに、ＸＹＺ座標系の１単位をＬとすると、磁気センサー項行列Ｐ（第３行列Ｐ）は数式４３に示すように１と−１とで表現されるので、計算がより簡便となる。
【０１７７】
【数４３】
2016217930.tif 000045
【０１７８】
（実施例１−３）
図６および図７は、実施例１−３に係るノイズ磁気センサーの配置を説明する図である。詳しくは、図６（ａ）は斜視図であり、図６（ｂ）は図６（ａ）の＋Ｘ方向側から見た平面図である。図７（ａ）は図６（ａ）の＋Ｙ方向側から見た平面図であり、図７（ｂ）は図６（ａ）の＋Ｚ方向側から見た平面図である。
【０１７９】
実施例１−３では、図６（ａ）および図７（ｂ）に示すように、直方体３０ａは、一辺の長さが２1/2×Ｌの正方形（Ｘ−Ｙ断面）を底面とし、高さ（Ｚ軸に沿った辺の長さ）が２×Ｌの四角柱である。さらに、８個のノイズ磁気センサー３０（直方体３０ａ）と心磁センサー１０との位置関係が実施例１−１および実施例１−２と異なる。より具体的には、図７（ｂ）に示すように、＋Ｚ方向側から見た平面視で、直方体３０ａは心磁センサー１０に対してＺ軸を中心として４５°回転した位置関係となっている。８個のノイズ磁気センサー３０は、直方体３０ａの各頂点に配置され、直方体３０ａの中心３０ｃと心磁センサー１０の中心１０ｃとはほぼ一致する。
【０１８０】
実施例１−３では、図６（ａ）および図６（ｂ）に示すように、８個のノイズ磁気センサー３０のうち、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４と第７のノイズ磁気センサー３７と第８のノイズ磁気センサー３８との４個が、Ｙ−Ｚ平面に平行なＸ＝０の平面に配置されている。そして、これら４個のノイズ磁気センサー３０は、Ｘ＝０の平面内で正方形（第１正方形とも称する）の各頂点に配置されている。図６（ｂ）に示す平面視で、第１正方形の重心は原点にほぼ一致するように配置されており、心磁センサー１０はこの第１正方形の重心と原点とを含むように配置されている。
【０１８１】
また、図６（ａ）および図７（ａ）に示すように、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２と第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６との４個が、Ｘ−Ｚ平面に平行なＹ＝０の平面に配置されている。そして、これら４個のノイズ磁気センサー３０は、Ｙ＝０の平面内で正方形（第２正方形とも称する）の各頂点に配置されている。図７（ａ）に示す平面視で、第２正方形の重心は原点にほぼ一致するように配置されており、心磁センサー１０はこの第２正方形の重心と原点とを含むように配置されている。
【０１８２】
このように配置されたノイズ磁気センサー３０の位置ベクトルｒk（磁気センサー位置ｒk）は、数式４４で表される。
【０１８３】
【数４４】
2016217930.tif 000046
【０１８４】
８個のノイズ磁気センサー３０をこのように配置すると、磁気センサー項行列Ｐ（第３行列Ｐ）は、数式４５に示すように１とＬと−Ｌとで表現されるうえ、ゼロの行列要素が増えるので、計算がより簡便となる。
【０１８５】
【数４５】
2016217930.tif 000047
【０１８６】
さらに、ＸＹＺ座標系の１単位をＬとすると、磁気センサー項行列Ｐ（第３行列Ｐ）は数式４６に示すように１と−１とで表現されるので、計算がより一層簡便となる。
【０１８７】
【数４６】
2016217930.tif 000048
【０１８８】
以上述べたように、第１の実施形態に係る磁気計測システム１によれば、心磁場等の微弱な磁場に対して大きな外部磁場が存在する場合でも、計測対象とする位置における外部磁場を演算により高精度に近似して計測値から除去できるので、計測対象とする微弱な磁場を精度良く計測することができる。
【０１８９】
（第２の実施形態）
次に、第２の実施形態に係る磁気計測システムにおける外部磁場の近似値の計算方法を説明する。第２の実施形態に係る磁気計測システムは、第１の実施形態に対して、外部磁場の近似値計算における未知数行列等の表現方法が異なる点以外は、システムの構成を含め第１の実施形態とおおむね同様である。
【０１９０】
＜外部磁場の近似値の計算方法＞
第１の実施形態では、外部磁場の近似値計算において、未知数ａijから数式２３にて表現される未知数行列ａが定義され、この未知数行列ａを数式３４や数式３５を用いて解いた。これに対して、第２の実施形態では、未知数ａijから数式４７にて表現される未知数ベクトルｂ（第１ベクトルｂとも称する）を定義し、この未知数ベクトルｂを解く点が異なっている。
【０１９１】
【数４７】
2016217930.tif 000049
【０１９２】
以下に、第２の実施形態に係る磁気計測システムにおいて、未知数ベクトルｂを解く方法を説明する。数式４７に示されるように、未知数ベクトルｂ（第１ベクトルｂ）は、２１個の未知数ａijを１列に並べた、２１行１列の列ベクトルである。次に、α個のノイズ磁気センサー３０（磁気センサー位置ｒk）における検出磁場ベクトルＢk全体で形成する検出磁場ベクトルＮ（第２ベクトルＮとも称する）を数式４８にて表す。
【０１９３】
【数４８】
2016217930.tif 000050
【０１９４】
数式４８に示される様に、検出磁場ベクトルＮ（第２ベクトルＮ）は、３×α個の検出磁場行列要素Ｂikを１列に並べた３α行１列の列ベクトルである。次に、α個の磁場センサー項ベクトルＲk全体で形成する磁気センサー項行列Ｑ（第４行列Ｑ）を数式４９にて定義する。
【０１９５】
【数４９】
2016217930.tif 000051
【０１９６】
数式４９の中式で、行ベクトルＲkTは、磁気センサー項ベクトルＲkの転置行列であり、１行７列の行ベクトルである。また、数式４９の中式で、ゼロベクトル０は、行列要素が全てゼロの１行７列の行ベクトルである。したがって、第４行列Ｑは、３α行２１列の行列である。数式４９にて定義される第４行列Ｑを用いると、検出磁場ベクトルＮと未知数ベクトルｂとは数式５０にて表現される。
【０１９７】
【数５０】
2016217930.tif 000052
【０１９８】
ここで、第１の実施形態と同様に、数式５０を行列要素にて表現すると数式３２のようになる。したがって、数式５０は、第１の実施形態と同じ方程式系を表していることが判る。
【０１９９】
数式５０で求めるべきは、未知数ベクトルｂである。ノイズ磁気センサー３０の数αが７であれば、磁気センサー項行列Ｑ（第４行列Ｑ）は２１行２１列の正方行列となり、逆行列が存在する。この場合は、未知数ベクトルｂ（第１ベクトルｂ）は、数式５１に示すように、第４行列Ｑの逆行列Ｑ-1を検出磁場ベクトルＮ（第２ベクトルＮ）に左から掛けて求められる。
【０２００】
【数５１】
2016217930.tif 000053
【０２０１】
一方、ノイズ磁気センサー３０の数αが８以上であれば、磁気センサー項行列Ｐは正方にはならないので、逆行列は存在しない。この場合は、未知数ベクトルｂ（第１ベクトルｂ）は、数式５２に示すように、第４行列Ｑの疑似逆行列（一般化逆行列とも云う）Ｑ+を検出磁場ベクトルＮ（第２ベクトルＮ）に左から掛けて求められる。
【０２０２】
【数５２】
2016217930.tif 000054
【０２０３】
数式５２において、磁気センサー項行列Ｑ（第４行列Ｑ）の疑似逆行列Ｑ+は、数式５３で求められる。
【０２０４】
【数５３】
2016217930.tif 000055
【０２０５】
数式５３に示されるように、疑似逆行列Ｑ+は、第４行列Ｑの転置行列ＱTと第４行列Ｑとの積の逆行列に、第４行列Ｑの転置行列ＱTを掛けたものである。なお、第４行列Ｑの転置行列ＱTは、第４行列Ｑの行列要素を行と列とに関して入れ替えたもので、２１行３α列の行列となり、数式５４にて表される。
【０２０６】
【数５４】
2016217930.tif 000056
【０２０７】
未知数ベクトルｂを求める際に疑似逆行列Ｑ+を用いると、最小二乗法の原理が働き、誤差を最小とする最適解が定まる。このように、第２の実施形態では、未知数ベクトルｂの全体に対して最小二乗法が適用されるので、第１の実施形態で求められる最適解よりもさらに適した解を得ることができ、計測対象空間での磁場をより正確に近似することができる。
【０２０８】
（第３の実施形態）
次に、第３の実施形態に係る磁気計測システムにおける外部磁場の近似値の計算方法を説明する。第３の実施形態に係る磁気計測システムは、第１の実施形態に対してシステムの構成は同様であり、第２の実施形態に対して外部磁場の近似値計算における未知数ベクトルｂ等の表現方法が異なる点以外はおおむね同様である。
【０２０９】
＜外部磁場の近似値の計算方法＞
第２の実施形態では、外部磁場の近似値計算において、２１個の未知数ａijから数式４７にて表現される未知数ベクトルｂ（第１ベクトルｂ）が定義され、この未知数ベクトルｂを数式５１や数式５２を用いて解いた。これに対して、第３の実施形態では、磁場Ｂにマクスウェル方程式の第二式を適応させて、未知数ａijを解く点が異なっている。
【０２１０】
マクスウェル方程式の第二式は、磁場に関するガウスの法則を表し、磁場の発散がゼロであることを表している。マクスウェル方程式の第二式は、数式５５にて表現される。この数式５５に示される関係を、数式２０や数式２１、数式２３等に適応させる。
【０２１１】
【数５５】
2016217930.tif 000057
【０２１２】
数式５５に数式２０を代入すると数式５６が得られる。
【０２１３】
【数５６】
2016217930.tif 000058
【０２１４】
数式５６の右辺において、最初の括弧はｘに関する偏微分を表し、２番目の括弧はｙに関する偏微分を表し、３番目の括弧はｚに関する偏微分を表している。磁場に関するガウスの法則では、数式５６が恒等的にゼロでなければならない。したがって、数式５６の右辺における定数項も、ｘに対する比例項も、ｙに対する比例項も、ｚに対する比例項も、全てゼロでなければならない。このことから、数式５７が求められる。
【０２１５】
【数５７】
2016217930.tif 000059
【０２１６】
数式５７は４つの恒等式を含むので、未知数ａijは２１個から２１−４＝１７個に減少する。具体的には、２１個の未知数ａijを解く際に数式５８を適応する。
【０２１７】
【数５８】
2016217930.tif 000060
【０２１８】
数式５８の各左辺に現れる４つの未知数は解く必要がないので、本実施形態では、数式５９にて表現される未知数ベクトルｃ（第３ベクトルｃとも称する）を定義し、これを解く。
【０２１９】
【数５９】
2016217930.tif 000061
【０２２０】
数式５９に示されるように、未知数ベクトルｃ（第３ベクトルｃ）は、２１個の未知数ａijからａ26とａ34とａ36とａ37との４つを除いた１７個の未知数を１列に並べた、１７行１列の列ベクトルである。この第３ベクトルｃに対応する磁気センサー項行列Ｓ（第５行列Ｓ）を数式６０にて定義する。
【０２２１】
【数６０】
2016217930.tif 000062
【０２２２】
数式６０に示されるように、磁気センサー項行列Ｓ（第５行列Ｓ）は、３α行１７列の行列である。３α行の３行毎に、一つの磁気センサー項ベクトルＲkが対応する。具体的には、ｋ番目の磁気センサー位置ｒkにおける磁気センサー項ベクトルＲkに相当する成分は、第５行列Ｓの３ｋ−２行目と３ｋ−１行目と３ｋ行目とに現れる。第５行列Ｓの３ｋ−２行目は、ｋ番目のノイズ磁気センサー３０が検出する磁場の１行目の成分Ｂ1kを求めるのに用いられる。
【０２２３】
同様に、第５行列Ｓの３ｋ−１行目は、ｋ番目のノイズ磁気センサー３０が検出する磁場の２行目の成分Ｂ2kを求めるのに用いられ、第５行列Ｓの３ｋ行目は、ｋ番目のノイズ磁気センサー３０が検出する磁場の３行目の成分Ｂ3kを求めるのに用いられる。第５行列Ｓの３ｋ−２行目と、第５行列Ｓの３ｋ−１行目と、第５行列Ｓの３ｋ行目と、を３行１７列の部分行列Ｔkにて表すと、部分行列Ｔkは数式６１となる。
【０２２４】
【数６１】
2016217930.tif 000063
【０２２５】
第５行列Ｓは、部分行列Ｔkを用いると、ｋ＝１の部分行列Ｔ1からｋ＝αの部分行列Ｔαまでのα個をα行１列に並べた行列となり、数式６２にて表現される。
【０２２６】
【数６２】
2016217930.tif 000064
【０２２７】
数式５９にて定義される未知数ベクトルｃ（第３ベクトルｃ）と数式６０にて定義される磁気センサー項行列Ｓ（第５行列Ｓ）とを用いると、検出磁場ベクトルＮは数式６３にて表現される。
【０２２８】
【数６３】
2016217930.tif 000065
【０２２９】
ここで、第２の実施形態と同様に、マクスウェル方程式の第二式を表す数式５８の関係を考慮して、数式６３を行列要素にて表現すると数式３２となる。したがって、数式６３は、第１の実施形態や第２の実施形態と同じ方程式系を表していることが判る。
【０２３０】
数式６３で求めるべきは、未知数ベクトルｃである。ノイズ磁気センサーの数αが６以上であれば、未知数の個数１７に対して、検出磁場ベクトルＮ（第２ベクトルＮ）は１８行以上となるので、最小二乗法を適応して、未知数を特定することができる。この場合、磁気センサー項行列Ｓ（第５行列Ｓ）は正方行列にはならないので、逆行列は存在しない。この場合は、未知数ベクトルｃ（第３ベクトルｃ）は、数式６４に示すように、磁気センサー項行列Ｓ（第５行列Ｓ）の疑似逆行列（一般化逆行列とも云う）Ｓ+を検出磁場ベクトルＮ（第２ベクトルＮ）に左から掛けて求められる。
【０２３１】
【数６４】
2016217930.tif 000066
【０２３２】
数式６４において、第５行列Ｓの疑似逆行列Ｓ+は数式６５で求められる。
【０２３３】
【数６５】
2016217930.tif 000067
【０２３４】
数式６５に示されるように、疑似逆行列はＳ+は、第５行列Ｓの転置行列ＳTと第５行列Ｓとの積の逆行列に、第５行列Ｓの転置行列ＳTを掛けたものである。なお、第５行列Ｓの転置行列ＳTは、第５行列Ｓの行列要素を行と列とに関して入れ替えたものである。
【０２３５】
未知数ベクトルｃを求める際に疑似逆行列Ｓ+を用いると、最小二乗法の原理が働き、誤差を最小とする最適解が定まる。このように、第３の実施形態では、マクスウェル方程式の第二式を考慮するので、第２の実施形態よりも少ない数のノイズ磁気センサー３０にて計測対象空間での磁場を特定することができる。また、第２の実施形態と同じ数のノイズ磁気センサー３０を用いる場合は、未知数が４個少ない分だけ、第２の実施形態で求められる最適解よりもさらに適した解を得ることができ、計測対象空間での磁場をより正確に近似することができる。
【０２３６】
なお、本実施形態では磁場に関するガウスの法則の結果（数式５８）を第２の実施形態に適応した例を示したが、磁場に関するガウスの法則の結果（数式５８）を第１の実施形態に適応してもよい。
【０２３７】
次に、ノイズ磁気センサー３０の配置の実施例を挙げて、第３の実施形態に係る外部磁場の近似値計算における磁気センサー項行列Ｓ（第５行列Ｓ）の計算方法を説明する。
【０２３８】
（実施例３−１）
図８および図９は、実施例３−１に係るノイズ磁気センサーの配置を説明する図である。詳しくは、図８（ａ）は斜視図であり、図８（ｂ）は図８（ａ）の＋Ｘ方向側から見た平面図である。図９（ａ）は図８（ａ）の＋Ｙ方向側から見た平面図であり、図９（ｂ）は図８（ａ）の＋Ｚ方向側から見た平面図である。
【０２３９】
実施例３−１では、図８（ａ）に示すように、６個のノイズ磁気センサー３０を２個ずつ、Ｘ軸上とＹ軸上とＺ軸上とに、原点を中心として対象に配置する。したがって、実施例３−１では、６個のノイズ磁気センサー３０が、一辺の長さが２1/2×Ｌの正八面体３０ｂの各頂点に配置される。正八面体３０ｂの中心３０ｃと心磁センサー１０の中心１０ｃとはほぼ一致する。
【０２４０】
実施例３−１では、図８（ａ）および図８（ｂ）に示すように、６個のノイズ磁気センサー３０のうち、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４と第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６との４個が、Ｙ−Ｚ平面に平行なＸ＝０の平面に配置されている。
【０２４１】
また、図８（ａ）および図９（ａ）に示すように、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２と第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６との４個が、Ｘ−Ｚ平面に平行なＹ＝０の平面に配置されている。さらに、図８（ａ）および図９（ｂ）に示すように、第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２と第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４との４個が、Ｘ−Ｙ平面に平行なＺ＝０の平面に配置されている。
【０２４２】
第１のノイズ磁気センサー３１と第２のノイズ磁気センサー３２とを結ぶ線分と、第３のノイズ磁気センサー３３と第４のノイズ磁気センサー３４とを結ぶ線分と、第５のノイズ磁気センサー３５と第６のノイズ磁気センサー３６とを結ぶ線分と、のうちの少なくとも２本の線分は直交する。そして、この互いに直交する２本の線分に平行な平面に対して、残りの１本の線分が交差する配置である。心磁センサー１０は、互いに直交する２本の線分に平行な平面に残りの１本の線分が交差する交差部を含む位置に配置されている。
【０２４３】
このように配置されたノイズ磁気センサー３０の位置ベクトルｒk（磁気センサー位置ｒk）は、数式６６で表される。
【０２４４】
【数６６】
2016217930.tif 000068
【０２４５】
６個のノイズ磁気センサー３０をこのように配置すると、数式６７に示すように、磁気センサー項行列Ｓ（第５行列Ｓ）は１とＬと−Ｌとで表現されるうえ、ゼロの行列要素が増えるので、計算がより簡便となる。
【０２４６】
【数６７】
2016217930.tif 000069
【０２４７】
さらに、ＸＹＺ座標系の１単位をＬとすると、数式６８に示すように、磁気センサー項行列Ｓ（第５行列Ｓ）は１と−１とで表現されるので、計算がより一層簡便となる。
【０２４８】
【数６８】
2016217930.tif 000070
【０２４９】
上述した実施形態は、あくまでも本発明の一態様を示すものであり、本発明の範囲内で任意に変形および応用が可能である。変形例としては、例えば、以下のようなものが考えられる。
【０２５０】
（変形例）
上記実施形態では、磁場Ｂを数式２０に示される３変数で変数に関する２次式の多項式で近似したが、磁場Ｂの近似式は数式２０に限定されるものではない。例えば、外部磁場が空間的に高次数の勾配磁場を含む場合は、数式２０にさらに高次の項を追加してもよい。なお、この場合は、未知数となる係数ａijの個数以上の数のノイズ磁気センサー３０が必要となるので、係数ａijの個数が多くなると配置するノイズ磁気センサー３０の個数も多くなり、未知数行列ａ（第１行列ａ）もより大きなものとなる。
【符号の説明】
【０２５１】
１…磁気計測システム、２…磁気計測装置（処理装置）、１０…心磁センサー（第１磁気センサー）、３０…ノイズ磁気センサー（第２磁気センサー）、３１…第１のノイズ磁気センサー（第１の第２磁気センサー）、３２…第２のノイズ磁気センサー（第２の第２磁気センサー）、３３…第３のノイズ磁気センサー（第３の第２磁気センサー）、３４…第４のノイズ磁気センサー（第４の第２磁気センサー）。
【図１】
2016217930.tif 000088
【図２】
2016217930.tif 000089
【図３】
2016217930.tif 000090
【図４】
2016217930.tif 000091
【図５】
2016217930.tif 000092
【図６】
2016217930.tif 000093
【図７】
2016217930.tif 000094
【図８】
2016217930.tif 000095
【図９】
2016217930.tif 000096