(57)【要約】【課題】微弱な磁場であっても十分な精度で計測することが可能な磁場計測装置を提供すること。【解決手段】計測対象の磁場と環境磁場とを含む磁場を計測する少なくとも１つの第１磁気センサーと、前記環境磁場を計測する複数の第２磁気センサーと、前記第２磁気センサーの計測値に基づいて前記環境磁場の分布を推定し、前記第１磁気センサーの計測値と推定した前記環境磁場の分布とに基づいて、前記計測対象の磁場を計算する磁場計算部と、を含む、磁場計測装置。【選択図】図４
2017191040.tif 000002
【特許請求の範囲】
【請求項１】
計測対象の磁場と環境磁場とを含む磁場を計測する少なくとも１つの第１磁気センサーと、
前記環境磁場を計測する複数の第２磁気センサーと、
前記第２磁気センサーの計測値に基づいて前記環境磁場の分布を推定し、前記第１磁気センサーの計測値と推定した前記環境磁場の分布とに基づいて、前記計測対象の磁場を計算する磁場計算部と、を含む、磁場計測装置。
【請求項２】
前記磁場計算部は、
前記第１磁気センサーの検出ベクトルと前記第１磁気センサーの位置情報と推定した前記環境磁場の分布とに基づいて、前記第１磁気センサーによる前記環境磁場の計測値を推定し、前記第１磁気センサーの計測値と推定した前記環境磁場の計測値とに基づいて、前記計測対象の磁場を計算する、請求項１に記載の磁場計測装置。
【請求項３】
前記磁場計算部は、
前記第１磁気センサーの計測値と前記第１磁気センサーの利得とに基づいて前記第１磁気センサーの位置における磁場の近似値を計算し、推定した前記環境磁場の計測値と前記第１磁気センサーの利得とに基づいて前記第１磁気センサーの位置における前記環境磁場の近似値を計算し、前記第１磁気センサーの位置における磁場の近似値と前記第１磁気センサーの位置における前記環境磁場の近似値との差分により前記計測対象の磁場を計算する、請求項２に記載の磁場計測装置。
【請求項４】
前記磁場計算部は、
前記第２磁気センサーの検出ベクトルと前記第２磁気センサーの位置情報と前記第２磁気センサーの計測値とに基づいて、前記環境磁場の分布を推定する、請求項１乃至３のいずれか１項に記載の磁場計測装置。
【請求項５】
前記第１磁気センサーの計測値と前記第２磁気センサーの計測値とに基づいて前記環境磁場の分布を推定し、推定した前記環境磁場の分布に基づいて前記第１磁気センサーの検出ベクトル及び前記第２磁気センサーの検出ベクトルを算出する校正部を含む、請求項２乃至４のいずれか１項に記載の磁場計測装置。
【請求項６】
前記磁場計算部は、
前記環境磁場を前記第１磁気センサーの位置を変数とする多項式で近似し、前記第２磁気センサーの計測値に基づいて、前記多項式の係数を計算することにより、前記環境磁場の分布を推定する、請求項１乃至５のいずれか１項に記載の磁場計測装置。
【請求項７】
前記磁場計算部は、
前記環境磁場の発散がゼロであるものとして、前記多項式の係数を計算する、請求項６に記載の磁場計測装置。
【請求項８】
前記磁場計算部は、
前記環境磁場の回転がゼロであるものとして、前記多項式の係数を計算する、請求項６又は７に記載の磁場計測装置。
【請求項９】
前記第２磁気センサーの検出ベクトルのうちの少なくとも２つは互いに直交する、請求項１乃至８のいずれか１項に記載の磁場計測装置。
【請求項１０】
前記第１磁気センサー及び前記第２磁気センサーは、
アルカリ金属原子が収容され、直線偏光が入射するセルと、
前記セルを出射した光を第１軸方向の光と第２軸方向の光に分離する偏光分離器と、
前記第１軸方向の光を検出する第１光検出器と、
前記第２軸方向の光を検出する第２光検出器と、を有する、請求項１乃至９のいずれか１項に記載の磁場計測装置。
【請求項１１】
前記第２磁気センサーが有する前記セルは、互いに同一平面上に配置されている、請求項１０に記載の磁場計測装置。
【請求項１２】
計測対象の磁場と環境磁場とを含む磁場を計測する少なくとも１つの第１磁気センサーの計測値を取得する工程と、
前記環境磁場を計測する複数の第２磁気センサーの計測値を取得する工程と、
前記第２磁気センサーの計測値に基づいて前記環境磁場の分布を推定する工程と、
前記第１磁気センサーの計測値と推定した前記環境磁場の分布とに基づいて、前記計測対象の磁場を計算する工程と、を含む、磁場計測方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、磁場計測装置及び磁場計測方法に関する。
【背景技術】
【０００２】
地磁気に比べて微弱な心臓の磁場（心磁場）や脳の磁場（脳磁場）等の生体磁場を計測するための磁場計測装置が知られている。磁場計測装置は、非侵襲であるため、被検体（生体）に負荷をかけずに臓器の状態を計測することができる。この種の磁場計測装置では、地磁気等の磁気ノイズの影響を低減させるために、被検体（生体）を磁気シールド装置の内部に収容して生体磁場を計測するが、シールドしきれない磁気ノイズにより計測精度が低下するおそれがある。
【０００３】
特許文献１には、統計的な手法を用いて、磁気測定データに重畳している磁気ノイズを低減することにより、磁気素子のオフセット成分を算出し、磁気測定データを校正する手法が開示されている。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特開２０１４−２５７１５号公報
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、一般に、心磁場や脳磁場等の生体磁場は、環境磁場（磁気ノイズ）と比較して非常に弱いため、このような微弱な磁場を計測する計測装置に特許文献１の手法を適用しても、十分な計測精度を確保することは難しい。
【０００６】
本発明は、以上のような問題点に鑑みてなされたものであり、本発明のいくつかの態様によれば、微弱な磁場であっても十分な精度で計測することが可能な磁場計測装置及び磁場計測方法を提供することができる。
【課題を解決するための手段】
【０００７】
本発明は前述の課題の少なくとも一部を解決するためになされたものであり、以下の態様または適用例として実現することが可能である。
【０００８】
［適用例１］
本適用例に係る磁場計測装置は、計測対象の磁場と環境磁場とを含む磁場を計測する少なくとも１つの第１磁気センサーと、前記環境磁場を計測する複数の第２磁気センサーと、前記第２磁気センサーの計測値に基づいて前記環境磁場の分布を推定し、前記第１磁気センサーの計測値と推定した前記環境磁場の分布とに基づいて、前記計測対象の磁場を計算する磁場計算部と、を含む。
【０００９】
本適用例に係る磁場計測装置によれば、磁場計算部は、第１磁気センサーの計測値に基づいて第１磁気センサーの位置における磁場（計測対象の磁場と環境磁場とを含む磁場）を計算することができ、複数の第２磁気センサーの計測値に基づいて推定した環境磁場の分布から第１磁気センサーの位置における環境磁場を計算することができるので、例えば、これらの磁場の差分やその近似値を計測対象の磁場として計算することができる。このように、本適用例に係る磁場計測装置によれば、計測対象の磁場が微弱な磁場であっても、相対的に大きい環境磁場を精度良く推定することにより、十分な精度で計測することが
できる。
【００１０】
［適用例２］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記第１磁気センサーの検出ベクトルと前記第１磁気センサーの位置情報と推定した前記環境磁場の分布とに基づいて、前記第１磁気センサーによる前記環境磁場の計測値を推定し、前記第１磁気センサーの計測値と推定した前記環境磁場の計測値とに基づいて、前記計測対象の磁場を計算してもよい。
【００１１】
本適用例に係る磁場計測装置によれば、磁場計算部は、推定した環境磁場の分布とともに、第１磁気センサーの検出ベクトル（検出軸の向きと利得の情報）や位置情報も用いることにより、第１磁気センサーによる環境磁場の計測値を精度良く推定することができるので、例えば、第１磁気センサーの計測値と推定した環境磁場の計測値との差分やその近似値を計測対象の磁場として精度良く計算することができる。
【００１２】
［適用例３］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記第１磁気センサーの計測値と前記第１磁気センサーの利得とに基づいて前記第１磁気センサーの位置における磁場の近似値を計算し、推定した前記環境磁場の計測値と前記第１磁気センサーの利得とに基づいて前記第１磁気センサーの位置における前記環境磁場の近似値を計算し、前記第１磁気センサーの位置における磁場の近似値と前記第１磁気センサーの位置における前記環境磁場の近似値との差分により前記計測対象の磁場を計算してもよい。
【００１３】
本適用例に係る磁場計測装置によれば、磁場計算部は、例えば、第１磁気センサーの計測値を第１磁気センサーの利得で除算することにより第１磁気センサーの位置における検出軸方向の磁場（計測対象の磁場と環境磁場とを含む磁場）を算出することができ、推定した環境磁場の計測値を第１磁気センサーの利得で除算することにより第１磁気センサーの位置における検出軸方向の環境磁場を算出することができる。そして、第１磁気センサーの検出軸の方向と計測方向とのずれが小さい場合、第１磁気センサーの位置における検出軸方向の磁場の計算値を第１磁気センサーの位置における磁場の近似値とし、第１磁気センサーの位置における検出軸方向の環境磁場の計算値を第１磁気センサーの位置における環境磁場の近似値とすることができる。従って、本適用例に係る磁場計測装置によれば、すべての第１磁気センサーの検出軸の方向が揃っているために計測対象の磁場を正しく計算できないような場合であっても、磁場計算部は、第１磁気センサーの位置における磁場（計測対象の磁場と環境磁場とを含む磁場）の近似値と環境磁場の近似値との差分として、計測対象の磁場を近似計算することができる。
【００１４】
［適用例４］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記第２磁気センサーの検出ベクトルと前記第２磁気センサーの位置情報と前記第２磁気センサーの計測値とに基づいて、前記環境磁場の分布を推定してもよい。
【００１５】
本適用例に係る磁場計測装置によれば、磁場計算部は、第２磁気センサーの計測値とともに、第２磁気センサーの検出ベクトル（検出軸の向きと利得の情報）や位置情報も用いることにより、環境磁場の分布を精度良く推定することができる。
【００１６】
［適用例５］
上記適用例に係る磁場計測装置は、前記第１磁気センサーの計測値と前記第２磁気センサーの計測値とに基づいて前記環境磁場の分布を推定し、推定した前記環境磁場の分布に基づいて前記第１磁気センサーの検出ベクトル及び前記第２磁気センサーの検出ベクトル
を算出する校正部を含んでもよい。
【００１７】
本適用例に係る磁場計測装置によれば、校正部は、例えば、第１磁気センサーによって計測対象の磁場が計測されない状態において、第１磁気センサーの計測値と第２磁気センサーの計測値に基づいて、第１磁気センサーと第２磁気センサーを含む空間における環境磁場の分布を推定することにより、第１磁気センサーの検出ベクトル（検出軸の向きと利得の情報）や第２磁気センサーの検出ベクトル（検出軸の向きと利得の情報）を精度良く算出することができる。従って、本適用例に係る磁場計測装置によれば、磁場計算部は、精度良く算出された第１磁気センサーの検出ベクトル（検出軸の向きと利得の情報）及び第２磁気センサーの検出ベクトル（検出軸の向きと利得の情報）を用いて、計測対象の磁場を精度良く計算することができる。
【００１８】
［適用例６］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記環境磁場を前記第１磁気センサーの位置を変数とする多項式で近似し、前記第２磁気センサーの計測値に基づいて、前記多項式の係数を計算することにより、前記環境磁場の分布を推定してもよい。
【００１９】
本適用例に係る磁場計測装置によれば、磁場計算部は、第１磁気センサーの位置を変数とする多項式を用いて第１磁気センサーの位置における環境磁場を精度良く近似し、第２磁気センサーの計測値に基づいて計算した当該多項式の係数に対応づけて第１磁気センサーの位置での環境磁場の分布を精度良く推定することができる。
【００２０】
［適用例７］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記環境磁場の発散がゼロであるものとして、前記多項式の係数を計算してもよい。
【００２１】
本適用例に係る磁場計測装置によれば、環境磁場の発散がゼロであるとの条件により、多項式の係数の数を減らすことができるので、磁場計算部の計算量が軽減され、あるいは、計測対象の磁場の計算精度が向上する。
【００２２】
［適用例８］
上記適用例に係る磁場計測装置において、前記磁場計算部は、前記環境磁場の回転がゼロであるものとして、前記多項式の係数を計算してもよい。
【００２３】
本適用例に係る磁場計測装置によれば、環境磁場の回転がゼロであるとの条件により、多項式の係数の数を減らすことができるので、磁場計算部の計算量が軽減され、あるいは、計測対象の磁場の計算精度が向上する。
【００２４】
［適用例９］
上記適用例に係る磁場計測装置において、前記第２磁気センサーの検出ベクトルのうちの少なくとも２つは互いに直交してもよい。
【００２５】
本適用例に係る磁場計測装置によれば、第２磁気センサーの検出ベクトル（検出軸）がすべて同一方向を向いている場合と比較して、第２磁気センサーの周辺の平面的な環境磁場の分布あるいは空間的な環境磁場の分布を精度良く推定することができるので、計測対象の磁場の計算精度が向上する。
【００２６】
［適用例１０］
上記適用例に係る磁場計測装置において、前記第１磁気センサー及び前記第２磁気センサーは、アルカリ金属原子が収容され、直線偏光が入射するセルと、前記セルを出射した
光を第１軸方向の光と第２軸方向の光に分離する偏光分離器と、前記第１軸方向の光を検出する第１光検出器と、前記第２軸方向の光を検出する第２光検出器と、を有してもよい。
【００２７】
本適用例によれば、第１磁気センサー及び第２磁気センサーとして光ポンピング式磁気センサーを用いて、微弱な磁場であっても十分な精度で計測することができる磁場計測装置を実現することができる。
【００２８】
［適用例１１］
上記適用例に係る磁場計測装置において、前記第２磁気センサーが有する前記セルは、互いに同一平面上に配置されていてもよい。
【００２９】
本適用例によれば、複数の第２磁気センサーのそれぞれのセルを１つの容器（保温機構）に収容して保温することができるとともに、各セルへの光の分岐機構を簡易化することができるので、磁場計測装置の製造コストを削減することができる。
【００３０】
［適用例１２］
上記適用例に係る磁場計測方法は、計測対象の磁場と環境磁場とを含む磁場を計測する少なくとも１つの第１磁気センサーの計測値を取得する工程と、前記環境磁場を計測する複数の第２磁気センサーの計測値を取得する工程と、前記第２磁気センサーの計測値に基づいて前記環境磁場の分布を推定する工程と、前記第１磁気センサーの計測値と推定した前記環境磁場の分布とに基づいて、前記計測対象の磁場を計算する工程と、を含む。
【００３１】
本適用例に係る磁場計測方法によれば、第１磁気センサーの計測値に基づいて第１磁気センサーの位置における磁場（計測対象の磁場と環境磁場とを含む磁場）を計算することができ、複数の第２磁気センサーの計測値に基づいて推定した環境磁場の分布から第１磁気センサーの位置における環境磁場を計算することができるので、例えば、これらの磁場の差分やその近似値を計測対象の磁場として計算することができる。このように、本適用例に係る磁場計測方法によれば、計測対象の磁場が微弱な磁場であっても、相対的に大きい環境磁場を精度良く推定することにより、十分な精度で計測することができる。
【図面の簡単な説明】
【００３２】
【図１】本実施形態に係る磁場計測装置の構成例を示す概略側面図。
【図２】磁気センサーユニットの模式側面図。
【図３】磁気センサーユニットの模式平面図。
【図４】処理装置の構成例を示す図。
【図５】本実施形態に係る磁場計測装置の校正方法についての説明図。
【図６】処理装置の校正部が校正処理を行う手順の一例を示すフローチャート図。
【図７】図７は、図６の工程Ｓ３〜工程Ｓ７の処理に対応するブロック線図。
【図８】検出ベクトル行列の更新処理の手順の一例を示すフローチャート図。
【図９】本実施形態に係る磁場計測方法についての説明図。
【図１０】処理装置の磁場計算部が磁場計算処理を行う手順の一例を示すフローチャート図。
【発明を実施するための形態】
【００３３】
以下、本発明の好適な実施形態について図面を用いて詳細に説明する。なお、以下に説明する実施の形態は、特許請求の範囲に記載された本発明の内容を不当に限定するものではない。また以下で説明される構成の全てが本発明の必須構成要件であるとは限らない。
【００３４】
１．第１実施形態
１−１．磁場計測装置の構成
図１は、本実施形態に係る磁場計測装置の構成例を示す概略側面図である。図１に示すように、本実施形態の磁場計測装置１は、計測対象物としての被検体（生体）９の心臓から発せられる心磁場や被検体（生体）９の脳から発せられる脳磁場等を計測する装置である。図１に示すように、磁場計測装置１は、図示しない少なくとも１つの第１磁気センサー１１（図２、図３及び図５参照）を備えた磁気センサーユニット１０と、複数の第２磁気センサー３０と、図示しない処理装置２（図５参照）と、土台３と、テーブル４と、磁気シールド装置６とを備えている。
【００３５】
磁気センサーユニット１０が備える第１磁気センサー１１は、計測対象となる心磁場や脳磁場等の微弱な磁場（計測対象の磁場）と外部磁場（磁気ノイズ）等の環境磁場とを含む磁場を計測するセンサーであり、心磁計や脳磁計等として使用される。第２磁気センサー３０は、外部磁場（磁気ノイズ）等の環境磁場を計測するセンサーである。第１磁気センサー１１および第２磁気センサー３０としては、光ポンピング式磁気センサー、ＳＱＵＩＤ式磁気センサー、フラックスゲート磁気センサー、ＭＩセンサー、ホール素子等を用いることができる。
【００３６】
磁場計測装置１の高さ方向（図１における上下方向）をＺ方向とする。Ｚ方向は鉛直方向である。土台３、テーブル４の上面が延在する方向をＸ方向およびＹ方向とする。Ｘ方向およびＹ方向は水平方向であり、Ｘ方向とＹ方向とは直交する方向である。横たわった状態の被検体９の身長方向（図１における左右方向）をＸ方向とする。
【００３７】
土台３は磁気シールド装置６（本体部６ａ）の内側の底面上に配置され、本体部６ａの外側にまで、Ｘ方向（被検体９の移動可能方向）に沿って延在している。テーブル４は、Ｘ方向テーブル４ａと、Ｚ方向テーブル４ｂと、Ｙ方向テーブル４ｃとを有している。土台３上には、Ｘ方向直動機構３ａによりＸ方向に沿って移動するＸ方向テーブル４ａが設置されている。Ｘ方向テーブル４ａの上には、図示しない昇降装置によりＺ方向に沿って昇降するＺ方向テーブル４ｂが設置されている。Ｚ方向テーブル４ｂの上には、図示しないＹ方向直動機構によりレール上をＹ方向に沿って移動するＹ方向テーブル４ｃが設置されている。
【００３８】
磁気シールド装置６は、開口部６ｃを有する角筒状の本体部６ａを備えている。本体部６ａの内部は空洞となっており、Ｙ方向およびＺ方向を通る面（Ｙ−Ｚ断面でＸ方向に直交した平面）の断面形状は概ね四角形になっている。心磁場を計測する際は、本体部６ａの内部に被検体９がテーブル４上に横たわった状態で収容される。本体部６ａはＸ方向に延在しており、これ自体でパッシブ磁気シールドとして機能する。
【００３９】
磁気センサーユニット１０および第２磁気センサー３０は、磁気シールド装置６の本体部６ａの内部に配置されている。磁気シールド装置６は、地磁気等の外部磁場が、磁気センサーユニット１０が配置された空間へ流入する事態を抑制している。すなわち、磁気シールド装置６により、磁気センサーユニット１０が配置された空間は外部磁場に比べて著しく低磁場とされ、外部磁場の磁気センサーユニット１０への影響が抑制されている。
【００４０】
本体部６ａの開口部６ｃから＋Ｘ方向に土台３が突出している。磁気シールド装置６の大きさは、例えば、Ｘ方向の長さが約２００ｃｍ程度であり、開口部６ｃの一辺が９０ｃｍ程度である。そして、開口部６ｃから、磁気シールド装置６内に、テーブル４に横たわった被検体９がテーブル４と共に土台３上をＸ方向に沿って移動して出入することができる。
【００４１】
図示しない処理装置２は、磁気センサーユニット１０が備える第１磁気センサー１１か
らの電気信号と、第２磁気センサー３０からの電気信号を受け取って、心磁場や脳磁場等の磁場を計測する装置である。処理装置２が発生する電気信号により磁場や残留磁場が発生して磁気センサーユニット１０に検出されるとノイズとなる。そのため、処理装置２は、発生される磁場や残留する磁場が磁気センサーユニット１０に到達し難くなるように、磁気シールド装置６の開口部６ｃから離れた場所に設置されているのが好ましい。
【００４２】
磁気シールド装置６の本体部６ａは、比透磁率が例えば数千以上の強磁性体、または、高伝導率の導体によって形成される。強磁性体にはパーマロイ、フェライト、または鉄、クロムもしくはコバルト系のアモルファス等を用いることができる。高伝導率の導体には、例えば、アルミニウム等で、渦電流効果によって磁場低減効果を有するものを用いることができる。なお、強磁性体と高伝導率の導体とを交互に積層して本体部６ａを形成することも可能である。
【００４３】
本体部６ａおよび土台３の＋Ｘ方向側および−Ｘ方向側の端には補正コイル（ヘルムホルツコイル）６ｂが設置されている。補正コイル６ｂの形状は枠状であり、本体部６ａを囲むように配置されている。補正コイル６ｂは、本体部６ａの内部空間へ流入する流入磁場を補正するためのコイルである。流入磁場は、外部磁場が開口部６ｃを通過して内部空間に入り込む磁場を指す。流入磁場は開口部６ｃに対してＸ方向で最も強くなる。補正コイル６ｂは、処理装置２から供給される電流により流入磁場をキャンセルするように磁界を発生させる。
【００４４】
磁気センサーユニット１０は、本体部６ａの天井に支持部材７を介して固定されている。磁気センサーユニット１０は、Ｚ方向における磁場の強度成分を計測する。すなわち、磁気センサーユニット１０が備える第１磁気センサー１１のそれぞれの検出軸は、Ｚ方向を向いている。被検体９の心磁場を計測する際は、被検体９における計測位置である胸部９ａが磁気センサーユニット１０と対向する位置になるようにＸ方向テーブル４ａおよびＹ方向テーブル４ｃを移動させ、胸部９ａが磁気センサーユニット１０に接近するようにＺ方向テーブル４ｂを上昇させる。
【００４５】
第２磁気センサー３０は、磁気センサーユニット１０の周囲に複数個配置されている。第２磁気センサー３０のそれぞれは、Ｘ方向、Ｙ方向又はＺ方向における磁場の成分を計測する。すなわち、第２磁気センサー３０のそれぞれの検出軸は、Ｘ方向、Ｙ方向又はＺ方向を向いている。
【００４６】
１−２．磁気センサーユニットの構成
図２及び図３は、本実施形態に係る磁気センサーユニット１０の構造を示す模式図である。詳しくは、図２は磁気センサーユニット１０の模式側面図であり、図３は磁気センサーユニット１０の模式平面図である。
【００４７】
図３に示すように、磁気センサーユニット１０には、レーザー光源１８からレーザー光１８ａが供給される。レーザー光源１８から発せられたレーザー光１８ａは光ファイバー１９を通って磁気センサーユニット１０に供給される。磁気センサーユニット１０と光ファイバー１９とは、光コネクター２０を介して接続されている。
【００４８】
レーザー光源１８は、セシウムの吸収線に応じた波長のレーザー光１８ａを出力する。レーザー光１８ａの波長は特に限定されないが、本実施形態では、例えば、Ｄ１線に相当する８９４ｎｍの波長に設定している。レーザー光源１８はチューナブルレーザーであり、レーザー光源１８から出力されるレーザー光１８ａは一定の光量を有する連続光である。
【００４９】
光コネクター２０を介して供給されたレーザー光１８ａは、−Ｙ方向に進行して偏光板２１に入射する。偏光板２１を通過したレーザー光１８ａは、直線偏光になっている。そして、レーザー光１８ａは、第１ハーフミラー２２、第２ハーフミラー２３、第３ハーフミラー２４、第１反射ミラー２５に順次入射する。
【００５０】
第１ハーフミラー２２、第２ハーフミラー２３および第３ハーフミラー２４は、レーザー光１８ａの一部を反射して＋Ｘ方向に進行させ、一部のレーザー光１８ａを通過させて−Ｙ方向に進行させる。第１反射ミラー２５は、入射したレーザー光１８ａを全て＋Ｘ方向に反射する。第１ハーフミラー２２、第２ハーフミラー２３、第３ハーフミラー２４、第１反射ミラー２５により、レーザー光１８ａは４つの光路に分割される。各光路のレーザー光１８ａの光強度が同じ光強度になるように、各ミラーの反射率が設定されている。
【００５１】
次に、図２に示すように、レーザー光１８ａは第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９に順次入射する。第４ハーフミラー２６、第５ハーフミラー２７および第６ハーフミラー２８は、レーザー光１８ａの一部を反射して＋Ｚ方向に進行させ、一部のレーザー光１８ａを通過させて＋Ｘ方向に進行させる。第２反射ミラー２９は、入射したレーザー光１８ａを全て＋Ｚ方向に反射する。
【００５２】
第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９により、１つの光路のレーザー光１８ａは４つの光路に分割される。各光路のレーザー光１８ａの光強度が同じ光強度になるように、各ミラーの反射率が設定されている。したがって、レーザー光１８ａは１６個の光路に分離される。そして、各光路のレーザー光１８ａの光強度が同じ強度になるように、各ミラーの反射率が設定されている。
【００５３】
第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９の＋Ｚ方向側には、レーザー光１８ａの各光路に、４行４列の１６個のガスセル１２が設置されている。そして、第４ハーフミラー２６、第５ハーフミラー２７、第６ハーフミラー２８、第２反射ミラー２９にて反射したレーザー光１８ａは、ガスセル１２を通過する。ガスセル１２は、内部に空隙を有する箱であり、この空隙にはアルカリ金属のガスが封入されている。アルカリ金属は特に限定されず、カリウム、ルビジウムまたはセシウムを用いることができる。本実施形態では、例えばアルカリ金属にセシウムを用いている。
【００５４】
各ガスセル１２の＋Ｚ方向側には、偏光分離器１３が設置されている。偏光分離器１３は、入射したレーザー光１８ａを、互いに直交する２つの偏光成分のレーザー光１８ａに分離する素子である。偏光分離器１３には、例えば、ウォラストンプリズムまたは偏光ビームスプリッターを用いることができる。
【００５５】
偏光分離器１３の＋Ｚ方向側には第１光検出器１４が設置され、偏光分離器１３の＋Ｘ方向側には第２光検出器１５が設置されている。偏光分離器１３を通過したレーザー光１８ａは第１光検出器１４に入射し、偏光分離器１３にて反射したレーザー光１８ａは第２光検出器１５に入射する。第１光検出器１４および第２光検出器１５は、入射したレーザー光１８ａの光量に応じた電流を処理装置２に出力する。
【００５６】
第１光検出器１４および第２光検出器１５が磁場を発生すると測定に影響を与える可能性があるので、第１光検出器１４および第２光検出器１５は非磁性の材料で構成されることが望ましい。磁気センサーユニット１０は、Ｘ方向の両面およびＹ方向の両面に設置されたヒーター１６を有している。ヒーター１６は磁界を発生しない構造であることが好ましく、例えば、流路中に蒸気や熱風を通過させて加熱する方式のヒーターを用いることができる。ヒーターの代わりに、高周波電圧によりガスセル１２を誘電加熱してもよい。
【００５７】
磁気センサーユニット１０は、被検体９（図１参照）の＋Ｚ側に配置される。被検体９が発する磁気ベクトルは、−Ｚ方向側から磁気センサーユニット１０に入る。磁気ベクトルは、第４ハーフミラー２６〜第２反射ミラー２９を通過し、ガスセル１２を通過した後、偏光分離器１３を通過して磁気センサーユニット１０から出る。
【００５８】
ガスセル１２内のセシウムは、加熱されてガス状態になっている。そして、直線偏光になったレーザー光１８ａをセシウムガスに照射することにより、セシウム原子が励起され磁気モーメントの向きが揃えられる。この状態でガスセル１２を磁気ベクトルが通過するとき、セシウム原子の磁気モーメントが磁気ベクトルの磁場により歳差運動する。この歳差運動をラーモア歳差運動と称する。
【００５９】
ラーモア歳差運動の大きさは、磁気ベクトルの磁場の強さと正の相関を有している。ラーモア歳差運動は、レーザー光１８ａの偏向面を回転させる。ラーモア歳差運動の大きさとレーザー光１８ａの偏向面の回転角の変化量とは、正の相関を有する。したがって、磁場の強さとレーザー光１８ａの偏向面の回転角の変化量とは、正の相関を有している。
【００６０】
偏光分離器１３は、レーザー光１８ａを直交する２成分の直線偏光に分離する。そして、第１光検出器１４および第２光検出器１５は、直交する２成分の直線偏光の強さを検出する。これにより、第１光検出器１４および第２光検出器１５は、レーザー光１８ａの偏向面の回転角を検出することができる。そして、処理装置２は、レーザー光１８ａの偏向面の回転角の変化から、磁場を計算することができる。
【００６１】
ガスセル１２、偏光分離器１３、第１光検出器１４、および第２光検出器１５により第１磁気センサー１１が構成される。この第１磁気センサー１１は、光ポンピング式磁気センサーや光ポンピング原子磁気センサーと称されるセンサーである。第１磁気センサー１１の感度は、Ｚ方向において高く、Ｚ方向と直交する方向において低くなっている。図３に示すように、例えば、磁気センサーユニット１０には、第１磁気センサー１１が４行４列の１６個配置されている。磁気センサーユニット１０における第１磁気センサー１１の個数および配置は特に限定されない。第１磁気センサー１１は、３行以下でもよく５行以上でもよい。同様に第１磁気センサー１１は、３列以下でもよく５列以上でもよい。第１磁気センサー１１の個数が多い程空間分解能を高くすることができる。
【００６２】
磁気センサーユニット１０が配置される計測対象空間は、磁気シールド装置６（図１参照）により外部磁場の流入が抑制されているが、外部磁場の流入を皆無とすることは困難である。換言すれば、磁気センサーユニット１０には、計測対象の磁場（心磁場）と環境磁場（磁気ノイズ）とが印加される。そのため、第１磁気センサー１１で計測して得られる計測値には、計測対象の磁場（心磁場）に基づく信号成分と環境磁場（磁気ノイズ）に基づくノイズ成分とが含まれる。したがって、計測対象の磁場（心磁場）を正確に取得するためには、第１磁気センサー１１で得られる計測値からノイズ成分を高精度に除去する必要がある。
【００６３】
１−３．第２磁気センサーの構成
第２磁気センサー３０は、磁気センサーユニット１０が配置される計測対象空間における環境磁場（磁気ノイズ）を計測するためのものである。第２磁気センサー３０で得られる計測値から計測対象空間における環境磁場（磁気ノイズ）を特定することにより、第１磁気センサー１１で得られる計測値から環境磁場（磁気ノイズ）成分を除去することができる。第２磁気センサー３０は、環境磁場（磁気ノイズ）を検知し、計測対象の磁場（心磁場）を検知しないものとする。
【００６４】
複数の第２磁気センサー３０の検出軸はすべてＺ方向を向いていてもよいが、複数の第２磁気センサー３０の検出軸（後述する検出ベクトルｋβ１〜ｋβＮ）のうちの少なくとも２つは互いに直交していることが好ましい。例えば、少なくとも１つの第２磁気センサー３０の検出軸がＺ方向を向いており、他の第２磁気センサー３０の検出軸はＸ方向又はＹ方向を向いていてもよい。これにより、第２磁気センサー３０の検出軸がすべてＺ方向を向いている場合と比較して、第２磁気センサー３０の周辺の平面的な環境磁場の分布あるいは空間的な環境磁場の分布を精度良く推定することができるので、後述する本実施形態に係る磁場計測方法による磁場の計算精度が向上する。
【００６５】
第２磁気センサー３０として用いるセンサーの種類は限定されないが、例えば、上述の第１磁気センサー１１と同様の光ポンピング式磁気センサーを用いることができる。すなわち、第１磁気センサー１１と同様、第２磁気センサー３０は、アルカリ金属原子が収容され、直線偏光が入射するセル（ガスセル１２）と、当該セルを出射した光を第１軸方向の光と第２軸方向の光に分離する偏光分離器１３と、第１軸方向の光を検出する第１光検出器１４と、第２軸方向の光を検出する第２光検出器１５と、を有してもよい。
【００６６】
そして、第２磁気センサー３０が有するセル（図２及び図３のガスセル１２に相当する）は、互いに同一平面上に配置されているのが好ましい。このようにすれば、複数の第２磁気センサー３０のそれぞれのセルを１つの容器（保温機構）に収容して保温することができるとともに、各セルへのレーザー光の分岐機構を簡易化することができるので、磁場計測装置１の製造コストを削減することができる。
【００６７】
１−４．処理装置の構成
図４は、処理装置２の構成例を示す図である。図４に示すように、処理装置２は、演算部１００、記憶部１１０、操作部１２０及び表示部１３０を含んで構成されている。
【００６８】
操作部１２０は、演算部１００が行う処理に必要な情報（磁場の計測開始指示や計測条件等の各種指示等）を入力するためのものであり、例えば、ボタンスイッチやレバースイッチ、ダイヤルスイッチ等の各種スイッチ、タッチパネル、キーボード、マウス等であってもよい。
【００６９】
表示部１３０は、演算部１００の処理結果を文字、グラフ、表、アニメーション、その他の画像として表示するものであり、例えば、ＬＣＤ（Liquid Crystal Display）やＥＬディスプレイ（Electroluminescence display）等であってもよい。
【００７０】
なお、１つのタッチパネル型ディスプレイで操作部１２０と表示部１３０の機能を実現するようにしてもよい。
【００７１】
記憶部１１０は、演算部１００が各種の処理を行うためのプログラムやデータ等を記憶するためのものであり、例えば、ＲＯＭ（Read Only Memory）やフラッシュＲＯＭ、ＲＡＭ（Random Access Memory）等の各種ＩＣメモリーやハードディスクやメモリーカードなどの記録媒体等により構成される。
【００７２】
特に、本実施形態では、記憶部１１０には、演算部１００によって読み出され、磁場計測装置１の校正処理を実行するための校正プログラム１１１と、計測対象の磁場を計算する処理（磁場計算処理）を実行するための磁場計算プログラム１１２とが記憶されている校正プログラム１１１及び磁場計算プログラム１１２はあらかじめ記憶部１１０に記憶されていてもよいし、演算部１００がネットワークを介してサーバーから校正プログラム１１１や磁場計算プログラム１１２を受信して記憶部１１０に記憶させてもよい。
【００７３】
また、記憶部１１０は、演算部１００の作業領域として用いられ、演算部１００が各種プログラムに従って実行した演算結果等を一時的に記憶する。さらに、記憶部１１０は、演算部１００の処理により生成されたデータのうち、長期的な保存が必要なデータを記憶してもよい。
【００７４】
演算部１００は、例えば、ＣＰＵ（Central Processing Unit）等のマイクロプロセッサーで実現され、上述した校正処理や磁場計算処理等を行う。
【００７５】
本実施形態では、演算部１００は、校正プログラム１１１を実行することにより校正部１０１として機能する。すなわち、校正プログラム１１１は処理装置２（コンピューター）を校正部１０１として機能させる（あるいは、処理装置２に校正処理を実行させる）ためのプログラムである。校正部１０１は、第１磁気センサー１１の計測値と第２磁気センサー３０の計測値を取得して磁場計測装置１の校正処理を行う。この校正処理の詳細については後述する。
【００７６】
また、本実施形態では、演算部１００は、磁場計算プログラム１１２を実行することにより磁場計算部１０２として機能する。すなわち、磁場計算プログラム１１２は処理装置２（コンピューター）を磁場計算部１０２として機能させる（あるいは、処理装置２に校正処理を実行させる）ためのプログラムである。磁場計算部１０２は、第１磁気センサー１１の計測値と第２磁気センサー３０の計測値を取得して磁場計算処理を行う。この磁場計算処理の詳細については後述する。
【００７７】
１−５．磁場計測装置の校正処理
本実施形態に係る磁場計測装置の校正方法について詳細に説明した後、処理装置２の校正部１０１が当該校正方法に対応する校正処理を行う手順について説明する。
【００７８】
本実施形態に係る校正方法は、第１磁気センサー１１や第２磁気センサー３０に限らず任意の磁気センサーを備えた磁場計測装置に適用可能であり、以下では、より一般的な概念に拡張して説明するために、第１磁気センサー１１と第２磁気センサー３０を区別せず、単に「磁気センサー」と称することにする。
【００７９】
図５に示すように、磁気センサーの個数Ｗ、各磁気センサーｉ（ｉ＝１〜Ｗ）の検出ベクトルｋｉや位置ベクトルｒｉは任意の値をとるものとする。検出ベクトルｋｉは、各磁気センサーｉの検出軸方向の単位ベクトルと各磁気センサーｉの利得との積を表すベクトルであり、位置ベクトルｒｉは原点Ｏから各磁気センサーｉの位置までをベクトル表現したものである。磁場計測装置１の校正処理とは、Ｗ個の磁気センサーの検出ベクトルｋ１〜ｋＷを求める処理である。
【００８０】
図５に示すように、Ｗ個の磁気センサーには、磁場ｂが印加される。磁場ｂは一様磁場だけではなく、高次の勾配磁場も含むものとする。校正用の磁場ｂは、人工的に造れた磁場でもよいし、地磁気のような自然の磁場でもよい。
【００８１】
任意の点（ｘ，ｙ，ｚ）での時間ｔにおける計算上の磁場の各成分ｂｘ（ｔ），ｂｙ（ｔ），ｂｚ（ｔ）は、計測する磁場分布の次数に合わせることが理想であるが、本実施形態では次式（１）の二次非線形多項式で表されるものとする。
【００８２】
【数１】
2017191040.tif 000003
【００８３】
ここで、磁気センサーｉの位置での時間ｔにおける磁場計算値を（ｂｉｘ（ｔ），ｂｉｙ（ｔ），ｂｉｚ（ｔ））とすると、Ｗ個の磁気センサーの位置での時間ｔにおける磁場計算値ベクトルｂ（ｔ）は、次式（２）で表される。
【００８４】
【数２】
2017191040.tif 000004
【００８５】
時間ｔ＝１〜Ｔでの磁場計算値ベクトルｂ（１）〜ｂ（Ｔ）を統合して、磁場計算値行列Ｂは次式（３）で表される。なお、式（３）において、ｔｒはベクトルの転置を表す。
【００８６】
【数３】
2017191040.tif 000005
【００８７】
次に、次式（４）のように、多項式（１）の係数の組を３０次元列ベクトルａ（ｔ）で表す。なお、式（４）において、ｔｒはベクトルの転置を表す。
【００８８】
【数４】
2017191040.tif 000006
【００８９】
ベクトルａ（ｔ）は時系列で変化することを想定するため、時間ｔ＝１〜Ｔにおけるベクトルａ（１）〜ａ（Ｔ）を統合した多項式係数行列Ａを次式（５）のように定義する。なお、式（５）において、ｔｒはベクトルの転置を表す。
【００９０】
【数５】
2017191040.tif 000007
【００９１】
定義された二次非線形多項式（１）に基づいて、３Ｗ×３０の位置情報行列Ｐを次式（６）のように定義する。
【００９２】
【数６】
2017191040.tif 000008
【００９３】
ここで、次の関係式（７）が成り立つ。すなわち、位置情報行列Ｐは、Ｗ個の磁気センサーの位置ベクトルｒ１〜ｒＷを磁場計算値行列Ｂに変換するための行列である。
【００９４】
【数７】
2017191040.tif 000009
【００９５】
また、Ｗ個の磁気センサーの利得（感度に相当）ｇ１〜ｇＷを要素に含む利得行列Ｇを次式（８）のように定義する。利得行列Ｇは、Ｗ×Ｗの正方行列である。
【００９６】
【数８】
2017191040.tif 000010
【００９７】
また、磁気センサーｉの検出軸方位をＸＹＺ直交座標系上の単位ベクトル（ｓｉｘ，ｓｉｙ，ｓｉｚ）で表し、Ｗ個の磁気センサーの検出軸方位の単位ベクトルを統合したＷ×３Ｗの検出軸行列Ｓを次式（９）のように定義する。なお、ｓｉｘ２＋ｓｉｙ２＋ｓｉｚ２＝１である。
【００９８】
【数９】
2017191040.tif 000011
【００９９】
そして、次式（１０）に示すように、利得行列Ｇと検出軸行列Ｓとの積を検出ベクトル行列Ｋとする。ここで、磁気センサーｉの検出ベクトルｋｉ＝（ｇｉｓｉｘ，ｇｉｓｉｙ，ｇｉｓｉｚ）であり、検出ベクトル行列Ｋは、Ｗ個の磁気センサーの検出ベクトルｋ１〜ｋＷを要素に含んでいる。
【０１００】
【数１０】
2017191040.tif 000012
【０１０１】
また、時間tにおける磁気センサーｉの計算上の観測値（推定値）を観測値計算値ｌｉ（ｔ）’で表し、時間tにおけるＷ個の磁気センサーの観測値計算値ｌ１（ｔ）’〜ｌＷ（ｔ）’を統合した観測値計算値ベクトルｌ（ｔ）’を次式（１１）のように定義する。なお、式（１１）において、ｔｒはベクトルの転置を表す。
【０１０２】
【数１１】
2017191040.tif 000013
【０１０３】
このとき、時間ｔ＝１〜Ｔにおける観測値計算値ベクトルｌ（１）’〜ｌ（Ｔ）’を統合した観測値計算値行列Ｌ’、検出ベクトル行列Ｋ、位置情報行列Ｐ及び多項式係数行列Ａの間には次式（１２）の関係式が成り立つ。式（１２）において、Ｄ＝ＫＰである。
【０１０４】
【数１２】
2017191040.tif 000014
【０１０５】
同様に、時間tにおける磁気センサーｉの実際の観測値（実測値）をセンサー観測値ｌｉ（ｔ）で表し、時間tにおけるＷ個の磁気センサーのセンサー観測値ｌ１（ｔ）〜ｌＷ（ｔ）を統合したセンサー観測値ベクトルｌ（ｔ）を次式（１３）のように定義する。なお、式（１３）において、ｔｒはベクトルの転置を表す。
【０１０６】
【数１３】
2017191040.tif 000015
【０１０７】
そして、次式（１４）のように、センサー観測値ベクトルｌ（ｔ）と観測値計算値ベクトルｌ（ｔ）’との差を観測値誤差ベクトルｖ（ｔ）とする。
【０１０８】
【数１４】
2017191040.tif 000016
【０１０９】
また、時間ｔ＝１〜Ｔにおける観測値誤差ベクトルｖ（１）〜ｖ（Ｔ）を統合した行列を観測値誤差行列Ｖとすると、次式（１５）のように、観測値誤差行列Ｖは、時間ｔ＝１〜Ｔにおけるセンサー観測値ベクトルｌ（１）〜ｌ（Ｔ）を統合したセンサー観測値行列Ｌと観測値計算値行列Ｌ’との差で表される。
【０１１０】
【数１５】
2017191040.tif 000017
【０１１１】
そして、次式（１６）のように定義される観測値誤差行列Ｖのノルム||Ｖ||を最小にする検出ベクトル行列Ｋと多項式係数行列Ａを求める最適化問題を解けば、Ｗ個の磁気センサーの検出ベクトルｋ１〜ｋＷが得られる。ただし、この最適化問題の解が収束するのに多大な時間を要する場合もあるため、実際には、観測値誤差行列Ｖのノルム||Ｖ||が許容値εよりも小さくなるときの検出ベクトルｋ１〜ｋＷを求めるのが現実的である。
【０１１２】
【数１６】
2017191040.tif 000018
【０１１３】
図６は、処理装置２の校正部１０１（図４参照）が上述した磁場計測装置１の校正方法に対応する校正処理を行う手順の一例を示すフローチャート図である。なお、図６の校正処理は、被検体（生体）９がテーブル４上に横たわっていない状態（被検体（生体）９からの心磁場や脳磁場の影響がない状態）で行われる。
【０１１４】
図６の例では、まず、校正部１０１は、センサー観測値行列Ｌを取得する（工程Ｓ１）。具体的には、校正部１０１は、ｔ＝１〜Ｔにおける各第１磁気センサー１１の計測値及
び各第２磁気センサー３０の計測値を取得し、これらの計測値をセンサー観測値ベクトルｌ（１）〜ｌ（Ｔ）として統合したセンサー観測値行列Ｌを取得する。
【０１１５】
次に、校正部１０１は、検出ベクトル行列Ｋを初期値Ｋ０に設定する（工程Ｓ２）。後述する工程Ｓ３〜Ｓ８の処理により、検出ベクトル行列Ｋを真値に近い値に収束させるためには、初期値Ｋ０は、真値との差が小さい値が好ましく、例えば、設計値（各第１磁気センサー１１の配置及び各第２磁気センサーの配置から推定される値）であってもよい。なお、初期値Ｋ０の各要素の値は、あらかじめ記憶部１１０に記憶されている。
【０１１６】
次に、校正部１０１は、多項式係数行列Ａを導出する（工程Ｓ３）。具体的には、校正部１０１は、工程Ｓ１で取得したセンサー観測値行列Ｌ、検出ベクトル行列Ｋ（工程Ｓ２で設定した初期値Ｋ０）及び位置情報行列Ｐから次式（１７）によって多項式係数行列Ａを導出する。なお、位置情報行列Ｐの各要素の値は、あらかじめ記憶部１１０に記憶されている。
【０１１７】
【数１７】
2017191040.tif 000019
【０１１８】
式（１７）において、（ＫＰ）＋はＫＰの擬似逆行列であり、Ｄ＋はＤ（＝ＫＰ）の擬似逆行列である。擬似逆行列Ｄ＋（＝（ＫＰ）＋）は次式（１８）で定義される。なお、式（１８）において、Ｔは行列の転置を表す。
【０１１９】
【数１８】
2017191040.tif 000020
【０１２０】
次に、校正部１０１は、磁場計算値行列Ｂを導出する（工程Ｓ４）。具体的には、校正部１０１は、工程Ｓ３で導出した多項式係数行列Ａ及び位置情報行列Ｐから式（７）によって磁場計算値行列Ｂを導出する。
【０１２１】
次に、校正部１０１は、検出ベクトル行列Ｋを更新する（工程Ｓ５）。具体的には、校正部１０１は、工程Ｓ１で取得したセンサー観測値行列Ｌと工程Ｓ４で導出した磁場計算値行列Ｂの擬似逆行列Ｂ＋との行列積に対応する行列を導出して検出ベクトル行列Ｋを更新する。ただし、実際には、センサー観測値行列Ｌと擬似逆行列Ｂ＋との行列積を計算しても検出ベクトル行列Ｋが正しく導出されない場合もあるため、本実施形態では、センサー観測値行列Ｌの要素であるセンサー観測値ベクトルｌ１〜ｌＷと磁場計算値行列Ｂの要素である磁場計算値行列ｂ１〜ｂＷから、検出ベクトル行列Ｋの要素である検出ベクトルｋ１〜ｋＷを導出することにより、検出ベクトル行列Ｋを更新する。この検出ベクトル行列Ｋの更新処理の詳細については後述する。
【０１２２】
次に、校正部１０１は、観測値計算値行列Ｌ’を導出する（工程Ｓ６）。具体的には、校正部１０１は、工程Ｓ４で導出した磁場計算値行列Ｂ及び工程Ｓ５で更新した検出ベクトル行列Ｋから次式（１９）によって観測値計算値行列Ｌ’を導出する。
【０１２３】
【数１９】
2017191040.tif 000021
【０１２４】
次に、校正部１０１は、観測値誤差行列Ｖを導出する（工程Ｓ７）。具体的には、校正部１０１は、工程Ｓ１で取得したセンサー観測値行列Ｌ及び工程Ｓ３で導出した観測値計算値行列Ｌ’から式（１５）によって観測値誤差行列Ｖを導出する。
【０１２５】
次に、校正部１０１は、観測値誤差行列Ｖのノルム||Ｖ||が許容値εよりも小さいか否かを判定する（工程Ｓ８）。具体的には、校正部１０１は、工程Ｓ３で導出した観測値誤差行列Ｖから式（１６）によってノルム||Ｖ||を算出し、許容値εと比較する。
【０１２６】
そして、校正部１０１は、ノルム||Ｖ||が許容値ε以上である場合は（工程Ｓ８のＮ）、工程Ｓ３以降の処理を再び行い、ノルム||Ｖ||が許容値εよりも小さくなれば（工程Ｓ８のＹ）、校正処理を終了する。なお、図７は、図６の工程Ｓ３〜工程Ｓ７の処理に対応するブロック線図である。
【０１２７】
図８は、検出ベクトル行列Ｋの更新処理（図６の工程Ｓ５の処理）の手順の一例を示すフローチャート図である。
【０１２８】
図８の例では、まず、校正部１０１は、変数ｉを１に初期化し（工程Ｓ５１）する。
【０１２９】
次に、校正部１０１は、センサー観測値ベクトルｌｉと磁場計算値行列ｂｉから、検出ベクトルｋｉを算出する（工程Ｓ５２）。具体的には、校正部１０１は、センサー観測値ベクトルｌｉと磁場計算値行列ｂｉから、次式（２０）によって検出ベクトルｋｉを算出する。なお、式（２０）において、Ｔはベクトルの転置を表す。
【０１３０】
【数２０】
2017191040.tif 000022
【０１３１】
ここで、センサー観測値ベクトルｌｉは、次式（２１）のように、磁気センサーｉによる時間ｔ＝１〜Ｔにおけるセンサー観測値ｌｉ（１）〜ｌｉ（Ｔ）を統合したベクトルとして定義され、式（２２）のように、センサー観測値行列Ｌはセンサー観測値ベクトルｌ１〜ｌＷを統合した行列である。
【０１３２】
【数２１】
2017191040.tif 000023
【０１３３】
【数２２】
2017191040.tif 000024
【０１３４】
また、磁場計算値行列ｂｉは、次式（２３）のように、磁気センサーｉの位置での時間ｔ＝１〜Ｔにおける磁場計算値を統合したベクトルとして定義され、式（２４）のように
、磁場計算値行列Ｂは、磁場計算値行列ｂ１〜ｂＷを統合した行列である。
【０１３５】
【数２３】
2017191040.tif 000025
【０１３６】
【数２４】
2017191040.tif 000026
【０１３７】
次に、校正部１０１は、変数ｉを１だけ増加する（工程Ｓ５３）。
【０１３８】
次に、校正部１０１は、変数ｉがＷよりも大きいか否かを判定する（工程Ｓ５４）。すなわち、校正部１０１は、検出ベクトルｋ１〜ｋＷの算出をすべて終了したか否かを判定する。
【０１３９】
そして、校正部１０１は、変数ｉがＷ以下である場合は（工程Ｓ５４のＮ）、工程Ｓ５２以降の処理を再び行い、変数ｉがよりも小さくなれば（工程Ｓ５４のＹ）、検出ベクトルｋ１〜ｋＷから検出ベクトル行列Ｋを導出し、検出ベクトル行列Ｋの更新処理を終了する。上述したように、検出ベクトルｋｉ＝（ｇｉｓｉｘ，ｇｉｓｉｙ，ｇｉｓｉｚ）であり、校正部１０１は、式（１０）のように、検出ベクトルｋ１〜ｋＷを統合することで、検出ベクトル行列Ｋを導出する（工程Ｓ５５）。
【０１４０】
１−６．磁場計測処理
本実施形態に係る磁場計測方法について詳細に説明した後、処理装置２の磁場計算部１０２が当該磁場計測方法に対応する磁場計算処理を行う手順について説明する。
【０１４１】
図９に示すように、ＸＹＺ直交座標系（絶対座標系とする）で表される空間中に、計測対象の磁場ｂｓ（心磁場）と環境磁場ｂｎ（磁気ノイズ）とを同時に計測する磁気センサ
ー（第１磁気センサー１１に相当する）がＭ個アレイ状に並べて構成されるα群と、環境磁場ｂｎ（磁気ノイズ）を計測する磁気センサー（第２磁気センサー３０に相当する）がＮ個分散あるいはアレイ状に並べて構成されるβ群が存在するものとする。なお、β群の磁気センサーは複数個（Ｎ≧２）であるが、α群の磁気センサーは１つ（Ｍ＝１）でもよいし、複数個（Ｍ≧２）でもよい。
【０１４２】
各磁気センサーｉは一軸成分出力のベクトル型の磁気センサーであり、それぞれ固有の検出ベクトルｋｉを持つ。検出ベクトルｋｉの長さは系全体における利得を表し、その方位は射影軸を表す。各磁気センサーｉによる計測値であるセンサー観測値ベクトルｌｉは、磁気センサーｉに印加される磁場と検出ベクトルｋｉの内積で表せるものとする。
【０１４３】
また磁気センサーに印加される環境磁場ｂｎは時事刻々と変化し、一様磁場だけではなく高次の勾配磁場も含む。環境磁場ｂｎは、人工的に造れた磁場でもよいし、地磁気のような自然の磁場でもよい。また、磁気センサーの周囲の磁性体により磁場分布が歪むことも有りうる。
【０１４４】
任意の位置（ｘｉ，ｙｉ，ｚｉ）、任意の時点における計算上の環境磁場の各成分（ｂｎｉｘ，ｂｎｉｙ，ｂｎｉｚ）は、計測する磁場分布の次数に合わせることが理想であり、本実施形態では次式（２５）の二次非線形多項式で表されるものとする。
【０１４５】
【数２５】
2017191040.tif 000027
【０１４６】
b群の磁気センサーｊの位置（ｘｊ，ｙｊ，ｚｊ）における環境磁場ｂｎの磁場計算値を（ｂｎｊｘ，ｂｎｊｙ，ｂｎｊｚ）とすると、b群のＮ個の磁気センサーの位置における磁場計算値ベクトルＢｎβは、次式（２６）で表される。
【０１４７】
【数２６】
2017191040.tif 000028
【０１４８】
次に、次式（２７）のように、多項式（２６）の係数の組を３０次元列ベクトル（多項式係数ベクトル）ａで表す。なお、式（２７）において、ｔｒはベクトルの転置を表す。
【０１４９】
【数２７】
2017191040.tif 000029
【０１５０】
定義された二次非線形多項式（２５）に基づいて、３Ｎ×３０の位置情報行列Ｐβを次式（２８）のように定義する。
【０１５１】
【数２８】
2017191040.tif 000030
【０１５２】
ここで、次の関係式（２９）が成り立つ。すなわち、位置情報行列Ｐβは、β群のＮ個の磁気センサーの位置ベクトルｒ１＝（ｘ１，ｙ１，ｚ１）〜ｒＮ＝（ｘＮ，ｙＮ，ｚＮ）を磁場計算値ベクトルＢｎβに変換するための行列である。
【０１５３】
【数２９】
2017191040.tif 000031
【０１５４】
また、β群のＮ個の磁気センサーの利得（感度に相当）ｇβ１〜ｇβＮを要素に含む利得行列Ｇβを次式（３０）のように定義する。利得行列Ｇβは、Ｎ×Ｎの正方行列である。
【０１５５】
【数３０】
2017191040.tif 000032
【０１５６】
また、β群の磁気センサーｊの検出軸方位をＸＹＺ直交座標系上の単位ベクトルｓβｊ＝（ｓβｊｘ，ｓβｊｙ，ｓβｊｚ）で表し、β群のＮ個の磁気センサーの検出軸方位の単位ベクトルｓβ１〜ｓβＮを統合したＮ×３Ｎの検出軸行列Ｓβを次式（３１）のように定義する。なお、ｓβｊｘ２＋ｓβｊｙ２＋ｓβｊｚ２＝１である。
【０１５７】
【数３１】
2017191040.tif 000033
【０１５８】
そして、次式（３２）に示すように、利得行列Ｇβと検出軸行列Ｓβとの積を検出ベクトル行列Ｋβとする。ここで、β群の磁気センサーｊの検出ベクトルｋβｊ＝（ｇβｊｓβｊｘ，ｇβｊｓβｊｙ，ｇβｊｓβｊｚ）であり、検出ベクトル行列Ｋβは、β群のＮ個の磁気センサーの検出ベクトルｋβ１〜ｋβＮを要素に含んでいる。検出ベクトル行列Ｋβは、互いに直交する検出ベクトルを含むため、環境磁場のＸＹＺ方向の３軸成分の検出が可能となる。
【０１５９】
【数３２】
2017191040.tif 000034
【０１６０】
また、β群の磁気センサーｊの計算上の観測値（推定値）を観測値計算値ｌβｊ’で表し、β群のＮ個の磁気センサーの観測値計算値ｌβ１’〜ｌβＮ’を統合した観測値計算値ベクトルｌβ’を次式（３３）のように定義する。なお、式（３３）において、ｔｒはベクトルの転置を表す。
【０１６１】
【数３３】
2017191040.tif 000035
【０１６２】
このとき、観測値計算値ベクトルｌβ’、利得行列Ｇβ、検出軸行列Ｓβ、位置情報行列Ｐβ、検出ベクトル行列Ｋβ、多項式係数ベクトルａ及び磁場計算値ベクトルＢｎβの間には次式（３４）の関係式が成り立つ。式（３４）において、Ｄβ＝ＫβＰβである。
【０１６３】
【数３４】
2017191040.tif 000036
【０１６４】
同様に、β群の磁気センサーｊの実際の観測値（実測値）をセンサー観測値ｌβｊで表し、β群のＮ個の磁気センサーのセンサー観測値ｌβ１〜ｌβＮを統合したセンサー観測値ベクトルｌβを次式（３５）のように定義する。なお、式（３５）において、ｔｒはベクトルの転置を表す。
【０１６５】
【数３５】
2017191040.tif 000037
【０１６６】
そして、次式（３６）のように、センサー観測値ベクトルｌβと観測値計算値ベクトルｌβ’との差を観測値誤差ベクトルｖβとする。
【０１６７】
【数３６】
2017191040.tif 000038
【０１６８】
そして、次式（３７）のように定義される観測値誤差ベクトルｖβのノルム||ｖβ||を最小にする多項式係数ベクトルａを求める最適化問題を解く。なお、式（３７）において、ｖβｊ＝ｌβｊ−ｌβｊ’である。
【０１６９】
【数３７】
2017191040.tif 000039
【０１７０】
この最適化問題を最小二乗法を用いて解く場合、センサー観測値ベクトルｌβから導出される多項式係数ベクトルａの一般的な解は次式（３８）のようになる。
【０１７１】
【数３８】
2017191040.tif 000040
【０１７２】
また、多項式係数ベクトルａは、Ｄβの疑似逆行列Ｄβ＋を用いて次式（３９）のようにも表される。
【０１７３】
【数３９】
2017191040.tif 000041
【０１７４】
次に、式（３９）により定められた多項式係数ベクトルａを用いて、α群のＭ個の磁気センサーの位置における磁場計算値ベクトルＢｎαを行列計算により求める。
【０１７５】
α群の磁気センサーｊの位置（ｘｊ，ｙｊ，ｚｊ）における環境磁場ｂｎの磁場計算値を（ｂｎｊｘ，ｂｎｊｙ，ｂｎｊｚ）とすると、α群のＭ個の磁気センサーの位置における磁場計算値ベクトルＢｎαは、次式（４０）で表される。
【０１７６】
【数４０】
2017191040.tif 000042
【０１７７】
定義された二次非線形多項式（２５）に基づいて、３Ｍ×３０の位置情報行列Ｐαを次式（４１）のように定義する。
【０１７８】
【数４１】
2017191040.tif 000043
【０１７９】
ここで、次の関係式（４２）が成り立つ。すなわち、位置情報行列Ｐαは、α群のＭ個の磁気センサーの位置ベクトルｒ１＝（ｘ１，ｙ１，ｚ１）〜ｒＭ＝（ｘＭ，ｙＭ，ｚＭ）を磁場計算値ベクトルＢｎαに変換するための行列である。
【０１８０】
【数４２】
2017191040.tif 000044
【０１８１】
また、α群のＭ個の磁気センサーの利得（感度に相当）ｇα１〜ｇαＮを要素に含む利得行列Ｇαを次式（４３）のように定義する。利得行列Ｇαは、Ｍ×Ｍの正方行列である。
【０１８２】
【数４３】
2017191040.tif 000045
【０１８３】
また、α群の磁気センサーｊの検出軸方位をＸＹＺ直交座標系上の単位ベクトルｓαｊ＝（ｓαｊｘ，ｓαｊｙ，ｓαｊｚ）で表し、α群のＭ個の磁気センサーの検出軸方位の単位ベクトルｓα１〜ｓαＭを統合したＭ×３Ｍの検出軸行列Ｓαを次式（４４）のように定義する。なお、ｓαｊｘ２＋ｓαｊｙ２＋ｓαｊｚ２＝１である。
【０１８４】
【数４４】
2017191040.tif 000046
【０１８５】
そして、次式（４５）に示すように、利得行列Ｇαと検出軸行列Ｓαとの積を検出ベクトル行列Ｋαとする。ここで、α群の磁気センサーｊの検出ベクトルｋαｊ＝（ｇαｊｓαｊｘ，ｇαｊｓαｊｙ，ｇαｊｓαｊｚ）であり、検出ベクトル行列Ｋαは、α群のＭ個の磁気センサーの検出ベクトルｋα１〜ｋαＭを要素に含んでいる。
【０１８６】
【数４５】
2017191040.tif 000047
【０１８７】
また、α群の磁気センサーｊの計算上の観測値（計測対象の磁場と環境磁場を含む磁場の推定値）を観測値計算値ｌαｊ’で表し、α群のＭ個の磁気センサーの観測値計算値ｌα１’〜ｌαＭ’を統合した観測値計算値ベクトルｌα’を次式（４６）のように定義する。なお、式（４６）において、ｔｒはベクトルの転置を表す。
【０１８８】
【数４６】
2017191040.tif 000048
【０１８９】
同様に、α群の磁気センサーｊの実際の観測値（実測値）をセンサー観測値ｌαｊで表し、α群のＭ個の磁気センサーのセンサー観測値ｌα１〜ｌαＭを統合したセンサー観測値ベクトルｌαを次式（４７）のように定義する。なお、式（４７）において、ｔｒはベクトルの転置を表す。
【０１９０】
【数４７】
2017191040.tif 000049
【０１９１】
このとき、観測値計算値ベクトルｌα’、利得行列Ｇα、検出軸行列Ｓα、位置情報行列Ｐα、検出ベクトル行列Ｋα及び多項式係数ベクトルａの間には次式（４８）の関係式が成り立つ。
【０１９２】
【数４８】
2017191040.tif 000050
【０１９３】
ここで、α群のＭ個の磁気センサーは、第１磁気センサー１１に相当し、例えば、被検体（生体）９の生体胸壁に垂直な磁場成分のみを計測すればよいものとする。従って、図９に示すように、α群の各磁気センサーｊの検出ベクトルｋαｊはＺ軸方向に揃っており、Ｚ軸に対する検出ベクトルｋαｊのずれθは小さいものとする。すなわち、単位ベクトルｓα１〜ｓαＭはすべて（０，０，１）に近似可能と考える。ただし、利得ｇα１〜ｇαＭは互いに異なるものとする。
【０１９４】
このとき、式（４８）で示される観測値計算値ベクトルｌα’を利得ｇα１〜ｇαＭの分だけ補正することで、すなわち、次式（４９）のように、式（４８）の右辺にＧα−１を左から掛けることで、α群のＭ個の磁気センサーの位置における環境磁場（磁気ノイズ）の近似値ｌｎ’が得られる。
【０１９５】
【数４９】
2017191040.tif 000051
【０１９６】
また、次式（５０）により、α群のＭ個の磁気センサーの位置における全磁場（計測対象の磁場（心磁場）＋環境磁場（磁気ノイズ））の近似値ｌｓｎ’が得られる。
【０１９７】
【数５０】
2017191040.tif 000052
【０１９８】
そして、次式（５１）のように、α群のＭ個の磁気センサーの位置における全磁場（計測対象の磁場（心磁場）＋環境磁場（磁気ノイズ））の近似値ｌｓｎ’と環境磁場（磁気ノイズ）の近似値ｌｎ’との差を取ることで、計測対象の磁場（心磁場）の近似値ｌｓ’が得られる。
【０１９９】
【数５１】
2017191040.tif 000053
【０２００】
図１０は、処理装置２の磁場計算部１０２（図４参照）が上述した磁場計測方法に対応する磁場計算処理を行う手順の一例を示すフローチャート図である。なお、図１０の磁場計算処理は、被検体（生体）９がテーブル４上に横たわっている状態（被検体（生体）９からの心磁場や脳磁場の計測が可能な状態）で行われる。また、図１０の磁場計算処理に先立ち、図６に手順の一例を示した校正処理が行われているものとする。すなわち、校正処理において、センサー観測値ベクトルｌα，ｌβが統合されたセンサー観測値行列Ｌに基づいて、検出ベクトル行列Ｋα，Ｋβが統合された検出ベクトル行列Ｋが更新（算出）され、検出ベクトル行列Ｋα，Ｋβの各要素の値が記憶部１１０に記憶されているものとする。
【０２０１】
図１０の例では、まず、磁場計算部１０２は、センサー観測値ベクトルｌα，ｌβを取得する（工程Ｓ１０１）。具体的には、磁場計算部１０２は、各第１磁気センサー１１の計測値を取得し、これらの計測値を統合したセンサー観測値ベクトルｌαを取得する。また、磁場計算部１０２は、各第２磁気センサー３０の計測値を取得し、これらの計測値を統合したセンサー観測値ベクトルｌβを取得する。
【０２０２】
次に、磁場計算部１０２は、センサー観測値ベクトルｌβから多項式係数ベクトルａを算出する（工程Ｓ１０２）。具体的には、磁場計算部１０２は、工程Ｓ１で取得したセンサー観測値ベクトルｌβ、校正処理で得られた検出ベクトル行列Ｋβ及び位置情報行列Ｐβから式（３８）によって多項式係数ベクトルａを算出する。なお、式（３８）において、Ｄβ＝ＫβＰβである。また、位置情報行列Ｐβの各要素の値は、あらかじめ記憶部１１０に記憶されており、検出ベクトル行列Ｋβの各要素の値は、校正処理において記憶部１１０に記憶される。
【０２０３】
次に、磁場計算部１０２は、センサー観測値ベクトルｌαから第１磁気センサー１１の位置での全磁場（計測対象の磁場＋環境磁場)の近似値ｌｓｎ’を算出する（工程Ｓ１０
３）。具体的には、各第１磁気センサー１１の検出ベクトルｋαｊはすべてＺ軸方向に揃っており、Ｚ軸に対する検出ベクトルｋαｊのずれθは小さいものとして、単位ベクトルｓαｊ＝（ｓαｊｘ，ｓαｊｙ，ｓαｊｚ）≒（０，０，１）に近似する。そして、磁場計算部１０２は、（ｓαｊｘ，ｓαｊｙ，ｓαｊｚ）≒（０，０，１）と検出ベクトル行列Ｋαから式（４３）及び式（４５）によって利得行列Ｇαを算出し、その逆行列Ｇα−１を算出する。さらに、磁場計算部１０２は、工程Ｓ１で取得したセンサー観測値ベクトルｌα及び利得行列Ｇαの逆行列Ｇα−１から式（５０）によって第１磁気センサー１１の位置での全磁場（計測対象の磁場＋環境磁場)の近似値ｌｓｎ’を算出する。なお、検出ベクトル行列Ｋαの各要素の値は、校正処理において記憶部１１０に記憶される。
【０２０４】
次に、磁場計算部１０２は、多項式係数ベクトルａから第１磁気センサー１１の位置での環境磁場の近似値ｌｎ’を算出する（工程Ｓ１０４）。具体的には、磁場計算部１０２は、利得行列Ｇαの逆行列Ｇα−１、検出ベクトル行列Ｋα、位置情報行列Ｐα及び工程Ｓ１０２で算出した多項式係数ベクトルａから式（４９）によって第１磁気センサー１１の位置での環境磁場の近似値ｌｎ’を算出する。なお、位置情報行列Ｐαの各要素の値は、あらかじめ記憶部１１０に記憶されている。
【０２０５】
最後に、磁場計算部１０２は、工程Ｓ１０３で算出した第１磁気センサー１１の位置での全磁場（計測対象の磁場＋県環境磁場）の近似値ｌｓｎ’及び工程Ｓ１０４で算出した第１磁気センサー１１の位置での環境磁場の近似値ｌｎ’から式（５１）によって計測対象の磁場の近似値ｌｓ’を算出し（工程Ｓ１０５）、磁場計測処理を終了する。
【０２０６】
１−７．作用効果
以上に説明したように、本実施形態では、磁場計算部１０２は、Ｎ個の第２磁気センサー３０の計測値（センサー観測値ベクトルｌβ）に基づいて環境磁場の分布（多項式係数ベクトルａ）を推定し、Ｍ個の第１磁気センサー１１の計測値（センサー観測値ベクトルｌα）と推定した環境磁場の分布（多項式係数ベクトルａ）とに基づいて、計測対象の磁場（計測対象の磁場の近似値ｌｓ’）を計算する。
【０２０７】
本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の計測値に基づいてＭ個の第１磁気センサー１１の位置における磁場（計測対象の磁場と環境磁場とを含む磁場）を計算することができ、Ｎ個の第２磁気センサー３０の計測値に基づいて推定した環境磁場の分布からＭ個の第１磁気センサー１１の位置における環境磁場を計算することができるので、例えば、これらの磁場の差分やその近似値を計測対象の磁場として計算することができる。このように、本実施形態に係る磁場計測装置１によれば、計測対象の磁場が微弱な磁場であっても、相対的に大きい環境磁場を精度良く推定することにより、十分な精度で計測することができる。
【０２０８】
また、本実施形態では、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の検出ベクトル（検出ベクトル行列Ｋα）とＭ個の第１磁気センサー１１の位置情報（位置情報行列Ｐα）と推定した環境磁場の分布（多項式係数ベクトルａ）とに基づいて、第１磁気センサー１１による環境磁場の計測値を推定し（式（４８））、Ｍ個の第１磁気センサー１１の計測値（センサー観測値ベクトルｌα）と推定した環境磁場の計測値（観測値計算値ベクトルｌα’）とに基づいて、計測対象の磁場（計測対象の磁場の近似値ｌｓ’）を計算する。
【０２０９】
このように、本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、推定した環境磁場の分布とともに、Ｍ個の第１磁気センサー１１の検出ベクトル（検出軸の向きと利得の情報）や位置情報も用いることにより、Ｍ個の第１磁気センサー１１による環境磁場の計測値を精度良く推定することができるので、計測対象の磁場を精度良く計算する
ことができる。
【０２１０】
また、本実施形態では、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の計測値（センサー観測値ベクトルｌα）とＭ個の第１磁気センサー１１の利得（利得行列Ｇα）とに基づいてＭ個の第１磁気センサー１１の位置における磁場の近似値（ｌｓｎ’）を計算し（式（５０））、推定した環境磁場の計測値（観測値計算値ベクトルｌα’）とＭ個の第１磁気センサー１１の利得（利得行列Ｇα）とに基づいてＭ個の第１磁気センサー１１の位置における環境磁場の近似値（ｌｎ’）を計算し（式（４９））、Ｍ個の第１磁気センサー１１の位置における磁場の近似値（ｌｓｎ’）とＭ個の第１磁気センサー１１の位置における環境磁場の近似値（ｌｎ’）との差分により計測対象の磁場（計測対象の磁場の近似値ｌｓ’）を計算する（式（５１））。
【０２１１】
このように、本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の計測値をＭ個の第１磁気センサーの利得で除算することによりＭ個の第１磁気センサーの位置における検出軸方向の磁場（計測対象の磁場と環境磁場とを含む磁場）を算出することができ、推定した環境磁場の計測値をＭ個の第１磁気センサー１１の利得で除算することによりＭ個の第１磁気センサー１１の位置における検出軸方向の環境磁場を算出することができる。そして、本実施形態では、Ｍ個の第１磁気センサー１１の検出軸の方向と計測方向（Ｚ軸方向）とのずれが小さいので、Ｍ個の第１磁気センサー１１の位置における検出軸方向の磁場の計算値をＭ個の第１磁気センサー１１の位置における磁場の近似値とし、Ｍ個の第１磁気センサー１１の位置における検出軸方向の環境磁場の計算値をＭ個の第１磁気センサー１１の位置における環境磁場の近似値とすることができる。従って、本実施形態に係る磁場計測装置１によれば、すべての第１磁気センサー１１の検出軸の方向が揃っているために計測対象の磁場を正しく計算できなくても、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の位置における磁場（計測対象の磁場と環境磁場とを含む磁場）の近似値と環境磁場の近似値との差分として、計測対象の磁場を近似計算することができる。
【０２１２】
また、本実施形態では、磁場計算部１０２は、Ｎ個の第２磁気センサー３０の検出ベクトル（検出ベクトル行列Ｋβ）とＮ個の第２磁気センサー３０の位置情報（位置情報行列Ｐβ）とＮ個の第２磁気センサーの計測値（センサー観測値ベクトルｌβ）とに基づいて、環境磁場の分布（多項式係数ベクトルａ）を推定する（式（３８））。
【０２１３】
このように、本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、Ｎ個の第２磁気センサーの計測値とともに、Ｎ個の第２磁気センサーの検出ベクトル（検出軸の向きと利得の情報）や位置情報も用いることにより、環境磁場の分布を精度良く推定することができる。
【０２１４】
また、本実施形態では、磁場計算部１０２は、環境磁場をＭ個の第１磁気センサー１１の位置を変数とする多項式で近似し（式（４２））、Ｎ個の第２磁気センサー３０の計測値（センサー観測値ベクトルｌβ）に基づいて、多項式の係数を計算することにより、環境磁場の分布（多項式係数ベクトルａ）を推定する。
【０２１５】
このように、本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、Ｍ個の第１磁気センサー１１の位置を変数とする多項式を用いてＭ個の第１磁気センサー１１の位置における環境磁場を精度良く近似し、Ｎ個の第２磁気センサー３０の計測値に基づいて計算した当該多項式の係数に対応づけてＭ個の第１磁気センサー１１の位置での環境磁場の分布を精度良く推定することができる。
【０２１６】
また、本実施形態では、校正部１０１は、Ｍ個の第１磁気センサー１１の計測値とＮ個
の第２磁気センサー３０の計測値とに基づいて（センサー観測値ベクトルｌα，ｌβが統合されたセンサー観測値行列Ｌに基づいて）環境磁場の分布（多項式係数行列Ａ）を推定し、推定した環境磁場の分布（多項式係数行列Ａ）に基づいてＭ個の第１磁気センサー１１の検出ベクトル及びＮ個の第２磁気センサーの検出ベクトル（検出ベクトル行列Ｋαと検出ベクトル行列Ｋβを統合した検出ベクトル行列Ｋ）を算出する。
【０２１７】
このように、本実施形態に係る磁場計測装置１によれば、校正部１０１は、Ｍ個の第１磁気センサー１１によって計測対象の磁場が計測されない状態において、Ｍ個の第１磁気センサー１１の計測値とＮ個の第２磁気センサー３０の計測値に基づいて、Ｍ個の第１磁気センサー１１とＮ個の第２磁気センサー３０を含む空間における環境磁場の分布を推定することにより、Ｍ個の第１磁気センサー１１の検出ベクトル（検出軸の向きと利得の情報）やＮ個の第２磁気センサー３０の検出ベクトル（検出軸の向きと利得の情報）を精度良く算出することができる。従って、本実施形態に係る磁場計測装置１によれば、磁場計算部１０２は、精度良く算出されたＭ個の第１磁気センサー１１の検出ベクトル（検出軸の向きと利得の情報）及びＮ個の第２磁気センサー３０の検出ベクトル（検出軸の向きと利得の情報）を用いて、計測対象の磁場を精度良く計算することができる。
【０２１８】
２．第２実施形態
第１実施形態に係る磁場計測装置１（磁場計測方法）では、環境磁場の分布を表す非線形多項式（２５）は、環境磁場が本来持っている規則性を一切考慮せずに設定されている。これに対して、第２実施形態に係る磁場計測装置１（磁場計測方法）では、環境磁場の発散が零という法則を反映させる点で第１実施形態と異なり、その他は第１実施形態と同様である。すなわち、第２実施形態に係る磁場計測装置１（磁場計測方法）では、次式（５２）が成り立つものとする。そして、磁場計算部１０２は、環境磁場の発散がゼロであるものとして、多項式（２５）の係数（多項式係数ベクトルａ）を計算する。
【０２１９】
【数５２】
2017191040.tif 000054
【０２２０】
そして、式（５２）に式（２５）を代入し、各係数の関係を求めると次式（５３）が成り立つ。
【０２２１】
【数５３】
2017191040.tif 000055
【０２２２】
式（５３）が恒等式であることを利用して、次式（５４）の関係式が得られる。
【０２２３】
【数５４】
2017191040.tif 000056
【０２２４】
式（５４）に示す関係式が４個得られることから、第１実施形態では３０個あった係数ａｘ１〜ａｘ１０，ａｙ１〜ａｙ１０，ａｚ１〜ａｚ１０が２６個に減少する。そして、磁場計算部１０２は、図１０の手順に従い、計測対象の磁場（近似値ｌｓ’）を計算する。従って、第２実施形態に係る磁場計測装置１（磁場計測方法）によれば、磁場計算部１０２の計算量が軽減され、あるいは、計測対象の磁場の計算精度が向上する。
【０２２５】
３．第３実施形態
第１実施形態に係る磁場計測装置１（磁場計測方法）では、環境磁場の分布を表す非線形多項式（２５）は、環境磁場が本来持っている規則性を一切考慮せずに設定されている。これに対して、第３実施形態に係る磁場計測装置１（磁場計測方法）では、環境磁場の回転が零という法則を反映させる点で第１実施形態と異なり、その他は第１実施形態と同様である。すなわち、第３実施形態に係る磁場計測装置１（磁場計測方法）では、次式（５５）が成り立つものとする。そして、磁場計算部１０２は、環境磁場の回転がゼロであるものとして、多項式（２５）の係数（多項式係数ベクトルａ）を計算する。なお、環境磁場の回転が零になるには、計測する空間内で伝導電流および変位電流が零となる条件が必要となるが、この条件は満たされるものとする。
【０２２６】
【数５５】
2017191040.tif 000057
【０２２７】
そして、式（５５）に式（２５）を代入し、各係数の関係を求めると次式（５６）が成り立つ。
【０２２８】
【数５６】
2017191040.tif 000058
【０２２９】
式（５６）が恒等式であることを利用して、次式（５７）の関係式が得られる。
【０２３０】
【数５７】
2017191040.tif 000059
【０２３１】
式（５７）に示す関係式が１２個得られるが、このうち、ａｙ７＝ａｘ６は、ａｚ５＝ａｙ７とａｘ６＝ａｚ５から求められるため、実際には１１個の関係式が得られることから、第１実施形態では３０個あった係数ａｘ１〜ａｘ１０，ａｙ１〜ａｙ１０，ａｚ１〜ａｚ１０が１９個に減少する。そして、磁場計算部１０２は、図１０の手順に従い、計測対象の磁場（近似値ｌｓ’）を計算する。従って、第３実施形態に係る磁場計測装置１（磁場計測方法）によれば、磁場計算部１０２の計算量が軽減され、あるいは、計測対象の磁場の計算精度が向上する。
【０２３２】
なお、本実施形態において、磁場計算部１０２は、第２実施形態と同様、さらに、環境磁場の発散がゼロであるものとして、多項式（２５）の係数（多項式係数ベクトルａ）を計算してもよい。これにより、係数がさらに４個減って、１５個にまで減少するので、磁場計算部１０２の計算量がさらに軽減され、あるいは、計測対象の磁場の計算精度がさらに向上する。
【０２３３】
４．変形例
本発明は本実施形態に限定されず、本発明の要旨の範囲内で種々の変形実施が可能である。
【０２３４】
例えば、上記の各実施形態では、磁場計測装置１の校正部１０１が校正処理を行っているが、磁場計測装置１とは異なる校正装置が磁場計測装置１の校正処理を行ってもよい。すなわち、磁場計測装置１の処理装置２は校正部１０１を有していなくてもよい。この場合、当該校正装置は、校正処理によって得られた検出ベクトル行列Ｋα，Ｋβを処理装置２の記憶部１１０に書き込み、磁場計測装置１（磁場計算部１０２）は、記憶部１１０に書き込まれた検出ベクトル行列Ｋα，Ｋβを用いて磁場計算処理を行えばよい。
【０２３５】
また、例えば、上記の各実施形態では、磁場計測装置１が被検体９（生体）の心磁場や脳磁場を計測しているが、磁場計測装置１は、心磁場や脳磁場以外の生体磁場を計測してもよいし、生体磁場以外の磁場（微弱な磁場）を計測するものであってもよい。
【０２３６】
また、上記の各実施形態では、α群の磁気センサーの検出ベクトルはＺ軸方向に揃っているものとして、計測対象の磁場（心磁場）の近似値ｌｓ’の計算式（５１）を導いたが、α群は、お互いに直交する２軸または３軸方向に検出軸を持つ複数の磁気センサーを含むものとして、以下のようにして、計測対象の磁場（心磁場）を計算してもよい。
【０２３７】
まず、式（４７）に示したセンサー観測値ベクトルｌαから、α群のＭ個の磁気センサーの周辺の磁場分布を表す多項式（例えば、二次非線形多項式）の多項式係数ベクトルａαを、式（３８）と同様の次式（５８）によって求める。式（５８）において、Ｄα＝ＫαＰαである。なお、心磁分布は環境磁場の分布と比較してより高次の項を含むため、α群のＭ個の磁気センサーの周辺の磁場分布を表す多項式の項は式（２５）と同数とは限らない。
【０２３８】
【数５８】
2017191040.tif 000060
【０２３９】
次に、式（５８）を次式（５９）に代入して、α群のＭ個の磁気センサーの周辺磁場Ｂｓｎを算出する。
【０２４０】
【数５９】
2017191040.tif 000061
【０２４１】
この周辺磁場Ｂｓｎは、環境磁場Ｂｎα（式（４２）で求められる）と心磁場Ｂｓを含むので、次式（６０）により、心磁場Ｂｓを算出することができる。
【０２４２】
【数６０】
2017191040.tif 000062
【０２４３】
上述した実施形態および変形例は一例であって、これらに限定されるわけではない。例えば、各実施形態および各変形例を適宜組み合わせることも可能である。
【０２４４】
本発明は、実施の形態で説明した構成と実質的に同一の構成（例えば、機能、方法及び結果が同一の構成、あるいは目的及び効果が同一の構成）を含む。また、本発明は、実施の形態で説明した構成の本質的でない部分を置き換えた構成を含む。また、本発明は、実施の形態で説明した構成と同一の作用効果を奏する構成又は同一の目的を達成することができる構成を含む。また、本発明は、実施の形態で説明した構成に公知技術を付加した構成を含む。
【符号の説明】
【０２４５】
１…磁場計測装置、２…処理装置、３…土台、３ａ…Ｘ方向直動機構、４…テーブル、４ａ…Ｘ方向テーブル、４ｂ…Ｚ方向テーブル、４ｃ…Ｙ方向テーブル、６…磁気シールド装置、６ａ…本体部、６ｂ…補正コイル（ヘルムホルツコイル）、６ｃ…開口部、７…支持部材、９…被検体（生体）、９ａ…胸部、１０…磁気センサーユニット、１１…第１磁気センサー、１２…ガスセル、１３…偏光分離器、１４…第１光検出器、１５…第２光検出器、１６…ヒーター、１８…レーザー光源、１８ａ…レーザー光、１９…光ファイバー
、２０…光コネクター、２１…偏光板、２２…第１ハーフミラー、２３…第２ハーフミラー、２４…第３ハーフミラー、２５…第１反射ミラー、２６…第４ハーフミラー、２７…第５ハーフミラー、２８…第６ハーフミラー、２９…第２反射ミラー、３０…第２磁気センサー、１００…演算部、１０１…校正部、１０２…磁場計算部、１１０…記憶部、１１１…校正プログラム、１１２…磁場計算プログラム、１２０…操作部、１３０…表示部
【図１】
2017191040.tif 000063
【図２】
2017191040.tif 000064
【図３】
2017191040.tif 000065
【図４】
2017191040.tif 000066
【図５】
2017191040.tif 000067
【図６】
2017191040.tif 000068
【図７】
2017191040.tif 000069
【図８】
2017191040.tif 000070
【図９】
2017191040.tif 000071
【図１０】
2017191040.tif 000072