(57)【要約】  測位システム１は、複数の受信機２０−Ａ、２０−Ｂ等と、各受信機から信号波の受信強度を示すデータを受け取って行列Ｑ・ベクトルｐ＝ベクトルｓ（ベクトルｐ：所定の基準位置からの波源の位置の各座標成分を各成分として含むベクトル、ベクトルｓ：所定の基準位置から各受信機までの距離の二乗を並べたベクトル）の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列を算出し、算出した逆行列をベクトルｓに作用させてベクトルｐに含まれる波源の各座標成分を算出する演算回路３１とを備える。行列Ｑおよびベクトルｐは、非線形項を線形項に置き換えた変数を含む。演算回路３１は、所定の受信機の受信強度と他の受信機の各受信強度との比の値を、第２線形成分と乗算される行または列の成分として含む行列Ｑの逆行列を求める。
2018062433.tif 000002
【特許請求の範囲】
【請求項１】
波源から出力された信号波をそれぞれ受信する４台または５台以上の複数の受信機であって、各々が、前記波源から出力された信号波を受信して前記信号波の受信強度を示すデータを出力する複数の受信機と、
前記複数の受信機の各々から前記信号波の受信強度を示すデータを受け取って、
行列Ｑ・ベクトルｐ＝ベクトルｓ
（ただし、ベクトルｐ：所定の基準位置からの前記波源の位置の各座標成分を各成分として含むベクトル、ベクトルｓ：前記所定の基準位置から前記各受信機までの距離の二乗を並べたベクトル）
の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列Ｑ-1を算出し、算出した前記逆行列Ｑ-1を前記ベクトルｓに作用させて、前記ベクトルｐに含まれる前記波源の位置の各座標成分を算出する演算回路を備え、
各受信機の位置の各座標成分と前記波源の位置の各座標成分との差の二乗の和が、前記各受信機と前記波源との距離の二乗に等しいという第１の関係式  を、
前記波源の位置の各座標成分を含む二次の項、一次の項、および、前記各受信機と前記波源との距離の二乗の項、の和が、所定の基準位置から前記各受信機までの距離の二乗に等しいという第２の関係式  に変形し、かつ、
前記ベクトルｐが、前記二次の項を第１線形成分に置換した成分、および、所定の受信機と前記波源との距離の二乗の項を第２線形成分に置換した成分をさらに含むよう定義したときにおいて、
前記演算回路は、前記所定の受信機の受信強度と他の受信機の各受信強度との比の値を、前記第２線形成分と乗算される行または列の成分として含む前記行列Ｑの逆行列Ｑ-1を求める、測位システム。
【請求項２】
前記複数の受信機が５台以上であり、前記行列Ｑが五次元の正則正方行列であるときにおいて、前記基準位置を（０，０，０）、前記波源の位置を（ｘ，ｙ，ｚ）、前記各受信機の位置を（ｘi，ｙi，ｚi）（ただしｉ＝ａ，ｂ，ｃ，ｄ，ｅ）とすると、
前記行列Ｑ、前記ベクトルｐ、および前記ベクトルｓは数１のように表現される、請求項１に記載の測位システム。
【数１】
2018062433.tif 000021
【請求項３】
前記ベクトルｐに含まれる前記ωは、前記基準位置から前記波源までの距離の二乗を表し、前記γは、前記所定の受信機から前記波源までの距離の二乗を表し、
前記行列Ｑに含まれる前記βb、βc、βdおよびβeは、前記所定の受信機の受信強度と他の各受信機の各受信強度との比である、請求項２に記載の測位システム。
【請求項４】
前記複数の受信機が４台以上であり、前記行列Ｑが四次元の正則正方行列であるときにおいて、前記基準位置を（０，０）、前記波源の位置を（ｘ，ｙ）、前記受信機の位置を（ｘi，ｙi）（ただしｉ＝ａ，ｂ，ｃ，ｄ）とすると、
前記行列Ｑ、前記ベクトルｐ、および前記ベクトルｓは数２のように表現される、請求項１に記載の測位システム。
【数２】
2018062433.tif 000022
【請求項５】
前記ベクトルｐに含まれる前記ωは、前記基準位置から前記波源までの距離の二乗を表し、前記γは、前記所定の受信機から前記波源までの距離の二乗を表し、
前記行列Ｑに含まれる前記βb、βcおよびβdは、前記所定の受信機の受信強度と他の各受信機の各受信強度との比である、請求項４に記載の測位システム。
【請求項６】
前記波源は、前記信号波として電磁波または音波を放射する、請求項１から５のいずれかに記載の測位システム。
【請求項７】
前記演算回路が算出した前記波源の座標を出力するインタフェースをさらに備えた、請求項１から６のいずれかに記載の測位システム。
【請求項８】
測位システムに用いられる測位装置であって、
前記測位システムは、波源から出力された信号波をそれぞれ受信する４台または５台以上の複数の受信機であって、各々が、前記波源から出力された信号波を受信して前記信号波の受信強度を示すデータを出力する複数の受信機を有し、
前記測位装置は、
前記複数の受信機の各々から前記信号波の受信強度を示すデータを受け取る入力端子と、
行列Ｑ・ベクトルｐ＝ベクトルｓ
（ただし、ベクトルｐ：所定の基準位置からの前記波源の位置の各座標成分を各成分として含むベクトル、ベクトルｓ：前記所定の基準位置から前記各受信機までの距離の二乗を並べたベクトル）
の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列Ｑ-1を算出し、算出した前記逆行列Ｑ-1を前記ベクトルｓに作用させて、前記ベクトルｐに含まれる前記波源の位置の各座標成分を算出する演算回路と
前記演算回路が算出した前記波源の位置の各座標成分のデータを出力する出力端子と
を備え、
各受信機の位置の各座標成分と前記波源の位置の各座標成分との差の二乗の和が、前記各受信機と前記波源との距離の二乗に等しいという第１の関係式  を、
前記波源の位置の各座標成分を含む二次の項、一次の項、および、前記各受信機と前記波源との距離の二乗の項、の和が、所定の基準位置から前記各受信機までの距離の二乗に等しいという第２の関係式  に変形し、かつ、
前記ベクトルｐが、前記二次の項を第１線形成分に置換した成分、および、所定の受信機と前記波源との距離の二乗の項を第２線形成分に置換した成分をさらに含むよう定義したときにおいて、
前記演算回路は、前記所定の受信機の受信強度と他の受信機の各受信強度との比の値を、前記第２線形成分と乗算される行または列の成分として含む前記行列Ｑの逆行列Ｑ-1を求める、測位装置。
【請求項９】
測位システムの測位装置に設けられたコンピュータによって実行されるコンピュータプログラムであって、
前記測位システムは、波源から出力された信号波をそれぞれ受信する４台または５台以上の複数の受信機であって、各々が、前記波源から出力された信号波を受信して前記信号波の受信強度を示すデータを出力する複数の受信機を有し、
前記コンピュータプログラムは、前記コンピュータに対し、
前記複数の受信機の各々から取得された前記信号波の受信強度を示すデータを読み出すステップと、
行列Ｑ・ベクトルｐ＝ベクトルｓ
（ただし、ベクトルｐ：所定の基準位置からの前記波源の位置の各座標成分を各成分として含むベクトル、ベクトルｓ：前記所定の基準位置から前記各受信機までの距離の二乗を並べたベクトル）
の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列Ｑ-1を算出するステップと、
算出した前記逆行列Ｑ-1を前記ベクトルｓに作用させて、前記ベクトルｐに含まれる前記波源の位置の各座標成分を算出するステップと、
前記演算回路が算出した前記波源の位置の各座標成分のデータを出力するステップと
を実行させ、
各受信機の位置の各座標成分と前記波源の位置の各座標成分との差の二乗の和が、前記各受信機と前記波源との距離の二乗に等しいという第１の関係式  を、
前記波源の位置の各座標成分を含む二次の項、一次の項、および、前記各受信機と前記波源との距離の二乗の項、の和が、所定の基準位置から前記各受信機までの距離の二乗に等しいという第２の関係式  に変形し、かつ、
前記ベクトルｐが、前記二次の項を第１線形成分に置換した成分、および、所定の受信機と前記波源との距離の二乗の項を第２線形成分に置換した成分をさらに含むよう定義したときにおいて、
前記逆行列Ｑ-1を算出するステップは、前記所定の受信機の受信強度と他の受信機の各受信強度との比の値を、前記第２線形成分と乗算される行または列の成分として含む前記行列Ｑの逆行列Ｑ-1を求める、コンピュータプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本開示は、測位システム、測位装置およびコンピュータプログラムに関する。
【背景技術】
【０００２】
従来、対象物を測位する測位システムの開発が行われてきた（たとえば特許文献１および２）。【０００３】
近年、Bluetooth（登録商標。以下同じ。）規格を利用した測位システムが普及しつつある。より具体的には、Bluetooth規格のうち、Bluetooth Low Energy（以下「ＢＬＥ」と略記する。）規格の信号を発生させる信号発生器を利用して、当該信号発生器の位置を推定するＢＬＥ測位システムが普及しつつある。【０００４】
ＢＬＥ測位システムの多くは３点測位アルゴリズムによる計測処理を採用している。３点測位アルゴリズムによる計測処理は、３つの受信機の各々の設置位置を示す座標情報、および、信号発生器と各受信機との各距離の情報を用いる。測位システムは、所定の演算を行い、３元連立非線形方程式の解として信号発生器の絶対位置座標を出力する。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特開２０１２−１９８０６６号公報
【特許文献２】特開平１１−１６０４０９号公報
【発明の概要】
【発明が解決しようとする課題】
【０００６】
３点測位アルゴリズムでは３元連立非線形方程式の解を算出することが必要となる。しかしながら３元連立非線形方程式を解くことには困難が伴う。そのため、通常は偏微分行列の逆行列演算とその繰り返し演算を利用するニュートン法が用いられる。【０００７】
ニュートン法では適性値への収束判定が必要であり、期待しない別解へ収束することがあり得る。また、ハードウェアやソフトウェアの負担が大きい等の問題がある。【０００８】
本願の、限定的ではない例示的なある実施形態は、ニュートン法とは異なる測位アルゴリズムを利用した測位システムを提供する。
【課題を解決するための手段】
【０００９】
本発明の例示的な実施形態によれば、測位システムは、波源から出力された信号波をそれぞれ受信する４台または５台以上の複数の受信機であって、各々が、波源から出力された信号波を受信して前記信号波の受信強度を示すデータを出力する複数の受信機と、前記複数の受信機の各々から前記信号波の受信強度を示すデータを受け取って、行列Ｑ・ベクトルｐ＝ベクトルｓ（ただし、ベクトルｐ：所定の基準位置からの前記波源の位置の各座標成分を各成分として含むベクトル、ベクトルｓ：前記所定の基準位置から前記各受信機までの距離の二乗を並べたベクトル）の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列Ｑ-1を算出し、算出した前記逆行列Ｑ-1を前記ベクトルｓに作用させて、前記ベクトルｐに含まれる前記波源の位置の各座標成分を算出する演算回路を備え、各受信機の位置の各座標成分と前記波源の位置の各座標成分との差の二乗の和が、前記各受信機と前記波源との距離の二乗に等しいという第１の関係式  を、前記波源の位置の各座標成分を含む二次の項、一次の項、および、前記各受信機と前記波源との距離の二乗の項、の和が、所定の基準位置から前記各受信機までの距離の二乗に等しいという第２の関係式  に変形し、かつ、前記ベクトルｐが、前記二次の項を第１線形成分に置換した成分、および、所定の受信機と前記波源との距離の二乗の項を第２線形成分に置換した成分をさらに含むよう定義したときにおいて、前記演算回路は、前記所定の受信機の受信強度と他の受信機の各受信強度との比の値を、前記第２線形成分と乗算される行または列の成分として含む前記行列Ｑの逆行列Ｑ-1を求める。
【発明の効果】
【００１０】
本発明の例示的な実施形態によれば、波源の位置を推定するためのアルゴリズムにおいて、演算に含まれる非線形項を線形化する。これにより、ニュートン法を利用する場合と比較して、多重解を回避でき、繰り返し演算を行う必要がなくなる。よってハードウェアへの負荷やソフトウェアの処理を軽減することができる。
【図面の簡単な説明】
【００１１】
【図１】図１は、本発明の例示的な実施形態にかかる測位システム１の構成を模式的に示す図である。
【図２】図２は、測位システム１が導入された環境の第一の例を示す図である。
【図３】図３は、測位システム１が導入された環境の第二の例を示す図である。
【図４】図４は、図２および図３における、棚２２および２３を含む４つの棚を示す図である。
【図５】図５は、棚２５に載置された商品２４ａに貼り付けられた信号発生器１０を示す図である。
【図６】図６は、棚２２および２３の間を走行する自走車両２４ｂに実装された信号発生器１０を示す図である。
【図７】図７は、棚２２および２３の間を歩行する買い物客のスマートフォン２４ｃに実装された信号発生器１０を示す図である。
【図８】図８は、信号発生器１０のハードウェアの構成を示す図である。
【図９】図９は、測位装置３０のハードウェアの構成を示す図である。
【図１０】図１０は、波源位置Ｐと各受信機との位置関係を示す図である。
【図１１】図１１は、ＣＰＵ３１の処理ブロックの構成を示す図である。
【図１２】図１２は、ＣＰＵ３１の処理手順を示すフローチャートである。
【図１３】図１３は、応用例にかかる処理を行うＣＰＵ３１の処理ブロックの構成を示す図である。
【図１４】図１４は、応用例にかかる処理を行うＣＰＵ３１の処理手順を示すフローチャートである。
【発明を実施するための形態】
【００１２】
上述した通り、３点測位アルゴリズムにニュートン法を利用した測位システムでは、適正な解に収束しないことがあるほか、さらに、別の問題も存在し得る。すなわち、信号発生器と各受信機との各距離の情報の取得方法である。
【００１３】
信号発生器と各受信機との各距離の情報は通常は既知ではない。直接、距離情報を得ることはできないため、ＢＬＥ測位システムでは、各受信機が受信した信号発生器からの信号の受信信号レベルを利用する。具体的には、事前に信号発生器と受信機とを１ｍ離れた距離に置き、その時の受信パワー値 Ｐ＿ｏ を予め取得しておく。その後、信号発生器から信号を受信すると、その時の受信信号レベルＰ＿ｉを利用して計測時の信号発生器と受信機との距離情報を得る。具体的には、計測時の距離情報（ｒ＿ｉ）2＝Ｐ＿ｏ／Ｐ＿ｉにより、求めることができる。
【００１４】
しかしながら、上述の処理には複数の問題が内包されている。まず、受信パワー値 Ｐ＿ｏ を予め取得しておくことが現実的には困難である場合がある。たとえば、施設内を移動する買い物客の位置を推定するために測位システムが用いられる場合を想定する。買い物客が所有する電子機器（たとえばスマートフォン）を信号発生器とすると、すべての電子機器を受信機の１ｍの位置に設置して受信パワー値を事前に求めることは現実的ではない。
【００１５】
利用される信号発生器の種類や数が予め決まっている施設（たとえば工場）内に設置された測位システムの場合には、上述した１ｍの位置の受信パワー値を事前に求めることが可能な場合はある。しかしながら、バッテリの残量等に起因して、事前計測時の信号発生器の送信出力と、測位時の信号発生器の送信出力とが異なる場合があり得る。受信機についても、事前計測時と測位時とで計測状況が異なる場合があり得る。
【００１６】
上述した問題は、全て測位誤差の要因になる。
【００１７】
本願発明者は、演算処理の負荷を軽減しつつ、測位誤差を抑え得る測位アルゴリズムについて検討を重ねた。その結果、新たな測位アルゴリズムを構築し、当該測位アルゴリズムを実装した測位システムを完成させるに至った。
【００１８】
本開示による測位システムでは、４台または５台以上の複数の受信機で、波源から出力された信号波をそれぞれ受信する。波源は、電磁波、音波等の信号波を出力する信号発生器である。
【００１９】
各受信機は、波源から出力された信号波を受信して信号波の受信強度を示すデータ（受信強度データ）を出力する。演算回路は、複数の受信機の各々から信号波の受信強度データを受け取って、所定の演算を行い、波源の位置の各座標成分を算出する。この演算は、所定の逆行列を求める演算を含むが、ニュートン法で利用される偏微分行列の逆行列演算とは異なる。
【００２０】
いま、２つのベクトルｐおよびｓを以下のように定義する。
【００２１】
ベクトルｐ：所定の基準位置からの波源の位置の各座標成分を各成分として含むベクトル
ベクトルｓ：所定の基準位置から各受信機までの距離の二乗を並べたベクトル
【００２２】
所定の基準位置は、測位システムが設置される空間における任意の位置であり、原点の位置として定義され得る。「位置」は、二次元空間内の位置であってもよいし、三次元空間内の位置であってもよい。
【００２３】
演算回路は、
行列Ｑ・ベクトルｐ＝ベクトルｓ
の関係式を満たす四次元または五次元の正則正方行列Ｑの逆行列Ｑ-1を算出し、算出した逆行列Ｑ-1をベクトルｓに作用させて、ベクトルｐに含まれる波源の位置の各座標成分を算出する。
【００２４】
なお、行列Ｑが正則正方行列である、という条件は、後に詳述するように、各受信機の位置を調整することにより容易に満足され得る。
【００２５】
以下、添付の図面を参照しながら、測位システム、測位装置およびコンピュータプログラムの構成例を説明する。ただし、必要以上に詳細な説明は省略する場合がある。たとえば、既によく知られた事項の詳細説明や実質的に同一の構成に対する重複説明を省略する場合がある。これは、以下の説明が不必要に冗長になるのを避け、当業者の理解を容易にするためである。なお、発明者は、当業者が本開示を十分に理解するために添付図面および以下の説明を提供するのであって、これらによって特許請求の範囲に記載の主題を限定することを意図するものではない。以下の説明においては、同一または類似する構成要素には、同一の参照符号を付している。
【００２６】
図１は、本発明の例示的な実施形態にかかる測位システム１の構成を模式的に示す。測位システム１は、信号発生器１０と、複数の受信機２０−Ａ、２０−Ｂ、・・・、２０−Ｎ（Ｎ：４または５以上の整数）と、測位装置３０とを備えている。
【００２７】
測位システム１の測位装置３０は、信号発生器１０から出力された信号波を利用して、信号発生器１０の位置を推定する。
【００２８】
信号発生器１０は、内蔵電池または外部から供給される電力を利用して信号を生成し、信号波として出力する電子回路である。信号発生器１０は、信号波として電磁波または音波を放射する。本明細書では、信号発生器１０を「波源」と呼ぶこともある。本明細書では、信号発生器１０は、ＢＬＥ規格に準拠する２．４ＧＨｚ帯域の電磁波を出力するとして説明する。
【００２９】
信号発生器１０は、測位対象が人である場合には、人が所有する電子機器に内蔵され得る。または、信号発生器１０は、測位対象が物である場合には当該物に貼り付けられてもよいし、当該物に内蔵されてもよい。
【００３０】
測位システム１は、測位システム１が設置された空間内で、信号発生器１０の位置を推定する。「空間」とは、本実施形態では主として三次元空間を想定している。図には、Ｘ軸、Ｙ軸およびＺ軸が示されている。ただし、「空間」は二次元空間であってもよい。たとえば、床面を走行する自走車両を測位する場合には、測位システム１は、床面である二次元空間上において、信号発生器１０の位置を推定してもよい。後述のように、測位システム１が二次元空間の位置を推定する場合には、受信機は少なくとも４台存在すればよい。測位システム１が三次元空間の位置を推定する場合には、受信機は少なくとも５台存在すればよい。
【００３１】
複数の受信機２０−Ａ、２０−Ｂ、・・・、２０−Ｎは、各々が、アンテナ装置（図示せず）を内蔵し、当該アンテナ装置を用いて波源から出力された信号波を受信する。アンテナ装置は、上述のＢＬＥ規格に準拠する電磁波を受信可能である。各受信機は、信号波の受信強度データを出力する。なお、アンテナ装置の構成は公知であるため、具体的な構成の説明は省略する。なお、信号発生器１０が音波を発生させる場合には、各受信機はマイクを有する装置であり得る。
【００３２】
測位装置３０は演算回路３１を有している。演算回路３１は、複数の受信機２０−Ａ、２０−Ｂ、・・・、２０−Ｎの各々から信号波の受信強度データを受け取って所定の演算を行い、波源の位置を示す各座標成分を算出する。演算の詳細は後に詳述する。
【００３３】
図２は、測位システム１が導入された環境の第一の例を示す。図２の例では、測位システム１は、物が載置される棚２２および２３等を有する工場、書店等に構築されている。複数の受信機の各々の位置が、「●」によって示されている。複数の受信機は、概ね分散して配置されている。なお、図２はＸＹ平面図であり、ｚ軸方向の関係は示されていない。
【００３４】
本実施形態では、複数の受信機は、予め定められた位置に設置されているとする。ある基準位置Ｏを原点とし、図示されるような３つの軸（Ｘ軸、Ｙ軸およびＺ軸）を定めたとき、各受信機の位置は、Ｘ軸、Ｙ軸およびＺ軸の各値によって表すことができる。以下では、受信機２０−Ａの座標をＡ（ｘa，ｙa，ｚa）などと表記し、ｘa，ｙaおよびｚaをそれぞれ「ｘ座標」「ｙ座標」「ｚ座標」と呼ぶことがある。図面上に例示された受信機２０−Ｊのｘ座標はｘJ、ｙ座標はｙJ、ｚ座標はｚJである。
【００３５】
なお、後述の理由から、本実施形態では演算に利用する複数の受信機の各ｘ座標が全て同じ値になることはないとする。換言すると、測位装置３０は、ｘ座標が少なくとも１つ異なる複数の受信機を利用して演算を行う。ｙ座標、および、ｚ座標についても同様である。そのため、たとえばｚ座標については、高さが少なくとも１つ異なる複数の受信機を利用することになる。
【００３６】
図３は、測位システム１が導入された環境の第二の例を示す。図面から理解されるように、図３では、複数の受信機が環境内の一部の区画Ｒにまとまって設けられている。なお図３もまたＸＹ平面図であり、ｚ軸方向の関係は示されていない。図３の構成は、上述した位置の条件を満たすように位置が調整されて一体化された複数の受信機のユニットを予め製造しておき、測位システム１の導入環境の天井部分に固定することによって容易に実現される。
【００３７】
本実施形態にかかる測位処理は、図２および図３のいずれの態様でも利用可能である。
【００３８】
次に、図４から図７を参照しながら、測位システム１における測位の対象を説明する。
【００３９】
信号発生器１０は、固定された位置に設置された物に設けられてもよいし、位置が変化し得る物、たとえば自走車両、持ち運び可能な電子機器に設けられてもよい。持ち運び可能な電子機器は、たとえば、携帯電話、スマートフォン、電子タグ装置である。
【００４０】
図４は、図２および図３における、棚２２および２３を含む４つの棚を示す。それぞれの棚は、物を載置するための複数の棚２５を有している。本実施形態では、複数の棚２５はＹ軸方向およびＺ軸方向に複数存在し得るとする。
【００４１】
図５は、棚２５に載置された商品２４ａに貼り付けられた信号発生器１０を示す。信号発生器１０から出力された信号波は、４または５台以上の受信機によって受信される。図５では、受信機２０−Ｓおよび２０−Ｔが、信号発生器１０から出力された信号波をそれぞれ受信している状態が示されている。測位システム１は、商品２４ａが載置された位置を推定し得る。
【００４２】
図６は、棚２２および２３の間を走行する自走車両２４ｂに実装された信号発生器１０を示す。本例でも、受信機２０−Ｓおよび２０−Ｔが、信号発生器１０から出力された信号波をそれぞれ受信している状態が示されている。測位システム１は、自走車両２４ｂの位置をリアルタイムで推定し得る。
【００４３】
図７は、棚２２および２３の間を歩行する買い物客のスマートフォン２４ｃに実装された信号発生器１０を示す。本例でも、受信機２０−Ｓおよび２０−Ｔが、信号発生器１０から出力された信号波をそれぞれ受信している状態が示されている。測位システム１は、買い物客の位置をリアルタイムで推定し得る。
【００４４】
図５〜図７の例では、信号発生器１０は１つだけ示されているが、複数存在してもよい。測位システム１は、並列的に処理を行うことにより、受信した信号波に含まれる、各信号発生器１０を一意に特定する識別情報から各信号発生器１０の位置を推定することができる。
【００４５】
次に、信号発生器１０および測位装置３０の構成を説明する。
【００４６】
図８は、信号発生器１０のハードウェアの構成を示す。信号発生器１０は、高周波信号を生成するためのＩＣ１１と、記憶装置１２と、アンテナ１４とを有する。記憶装置１２はたとえばフラッシュＲＯＭであり、信号発生器１０毎に一意の識別情報１３が格納されている。ＩＣ１１は、アンテナ１４を利用して識別情報を周期的に送信する。
【００４７】
図９は、測位装置３０のハードウェアの構成を示す。
【００４８】
測位装置３０は、ＣＰＵ３１と、メモリ３２と、通信回路３３とを有しており、これらは内部バスで接続されている。ＣＰＵ３１は、後述の処理により、個々の信号発生器１０の位置を推定し、推定した位置を示す位置情報を生成する演算回路である。メモリ３２は、たとえばＤＲＡＭであり、ＣＰＵ３１の処理に関連して利用されるワークメモリである。
【００４９】
通信回路３３は、たとえば、１または複数の通信コネクタを有する通信回路である。通信コネクタは、入力インタフェース３４ａ、および、波源の座標のデータを出力する出力インタフェース３４ｂを含む。入力インタフェース３４ａおよび出力インタフェース３４ｂは一体化され１つの通信コネクタとして実装されてもよい。
【００５０】
入力インタフェース３４ａは、いわゆる信号入力端子であり、受信機２０−Ａ〜２０−Ｎの各々から高周波電気信号を受信する。高周波電気信号は、各受信機が受信した電磁波（信号波）を変換して生成された信号である。
【００５１】
出力インタフェース３４ｂは、たとえばイーサネット（登録商標）規格の有線通信を行う通信端子であり、波源の座標のデータを出力する。波源の座標のデータに代えて、出力インタフェース３４ｂは、波源の座標を画像化した映像信号を出力してもよい。このとき、出力インタフェース３４ｂは、ＤＶＩ端子等の画像信号の出力端子であり得る。図９には、表示装置３５と接続された出力インタフェース３４ｂが示されている。
【００５２】
次に、本願発明者が完成させた、本実施形態による測位アルゴリズムを説明する。
【００５３】
三次元空間での測位のため、以下では５台の受信機を利用する。二次元空間での測位の場合には、４台の受信機があればよい。上述のとおり、受信機２０−Ａの位置Ａを（ｘa，ｙa，ｚa）などと表す。また、波源位置Ｐから受信機２０−Ａまでの距離をｒaなどと表す。測位すべき波源（信号発生器１０）の位置Ｐを（ｘ，ｙ，ｚ）と表す。
【００５４】
図１０は、波源位置Ｐと各受信機との位置関係を示す。図１０には、波源位置Ｐ、各受信機の位置Ａ〜Ｄが示されている。また、図１０には基準位置として原点Ｏが示されている。原点Ｏは任意に定め得る。各受信機の位置は、原点Ｏの位置（０，０，０）を基準として決定されている。
【００５５】
波源位置Ｐと受信機２０−Ａの位置Ａとの間には以下の関係が成り立つ。
【数１】
2018062433.tif 000003
【００５６】
同様に、波源位置Ｐと他の各受信機の位置との間にも以下の関係が成り立つ。
【数２】
2018062433.tif 000004
【数３】
2018062433.tif 000005
【数４】
2018062433.tif 000006
【数５】
2018062433.tif 000007
【００５７】
ここでパラメータβb〜βeを導入する。パラメータβbは、波源位置Ｐから受信機２０−Ａまでの距離ｒaに対する、波源位置Ｐから受信機２０−Ｂまでの距離ｒbの大きさの二乗である。よって以下のように表記される。
（数６）
βb＝（ｒb／ｒa）2
同様に、βc、βdおよびβeも、波源位置Ｐから受信機２０−Ａまでの距離ｒaと、波源位置Ｐから各受信機までの距離を用いて、まとめて以下のように表すことができる。
（数７）
βc＝（ｒc／ｒa）2
βd＝（ｒd／ｒa）2
βe＝（ｒe／ｒa）2
【００５８】
数６および７を利用して数２〜５を変形すると、以下の数８が得られる。
【数８】
2018062433.tif 000008
【００５９】
数１および数８の各式を展開して変形すると、以下の式が得られる。
【数９】
2018062433.tif 000009
【００６０】
いま、数９に含まれる未知数を、以下のように置き換える。
【数１０】
2018062433.tif 000010
【００６１】
上述の変数の置き換えは、二次の項（非線形項）を一次の項（線形項）に変更する操作である。つまり非線形方程式を線形化することに相当する。
【００６２】
すると、数９は行列を用いて以下のように表現することができる。
【数１１】
2018062433.tif 000011
【００６３】
数１１は行列Ｑ、列ベクトルｐおよびｓを用いて以下のように表現できる。
（数１２）
Ｑ・ｐ＝ｓ
【００６４】
ベクトルｐに含まれる要素のうち、ｘ、ｙおよびｚが、波源Ｐの位置を表す未知数である。行列Ｑの逆行列をＱ-1と表すと、ベクトルｐは以下の式により求めることができる。
（数１３）
ｐ＝Ｑ-1・ｓ
【００６５】
数１３は具体的には以下のように表現される。
【数１４】
2018062433.tif 000012
【００６６】
数１３および数１４が成立する条件として、行列Ｑの逆行列が存在すること、換言すれば、行列Ｑが正則正方行列であることが必要である。行列式が０でないことが条件となる。当該条件を満たすためには、（ｉ）行列Ｑのある行と他の行との間には整数倍の関係がないこと、および、（ｉｉ）行列Ｑのある列と他の列との間には整数倍の関係がないこと、が必要である。数１４に示されるように、行列Ｑの行に注目すると、上記関係（ｉ）は明らかに満たさない。上記関係を満たす場合には、５つの受信機の位置が少なくとも同じになる必要があるが、そのような配置はありえないからである。
【００６７】
そこで行列Ｑの列の関係を検討する。５つの受信機のすべてのｘ座標、ｙ座標またはｚ座標が一致する場合には、上記関係（ｉｉ）を満たす。そのため、５つの受信機を選択する際に、そのような条件を満たさない受信機の組を採用しなければよい。たとえば、ｚ座標については、高さが少なくとも１つ異なる５台の受信機を選択する。ｘ座標の関係およびｙ座標の関係も同様である。
【００６８】
次に、行列Ｑの逆行列を求める際の、パラメータβb〜βeの取り扱いを説明する。パラメータβb〜βeは、受信機２０−Ａ〜２０−Ｅの受信強度から求めることができる。
【００６９】
まず、フリスの伝達公式によれば、下記の式が成り立つことが知られている。
（数１５）
Ｐｒ＝（λ／４πｒ）2・Ｐｔ・Ｇｔ・Ｇｒ
【００７０】
ここでＰｒは受信電力、λは定数、ｒは送信機と受信機との距離、Ｐｔは信号発生器１０の送信電力、Ｇｔは信号発生器１０の送信アンテナ１４（図８）の利得、Ｇｒは受信アンテナの利得である。ここで、ＰｔおよびＧｔは、共通の波源であるから同じ値である。また、受信アンテナの利得Ｇｒは、すべての受信機２０−Ａ〜２０−Ｅで共通化しておくことにより、同じ値であると仮定することができる。
【００７１】
受信機２０−Ａ〜２０−Ｅの各々波源との関係は以下のように表現できる。なお、受信機２０−Ａ〜２０−Ｅの受信強度をＰａ〜Ｐｅと表す。
（数１６）
Ｐａ＝（λ／４πｒa）2・Ｐｔ・Ｇｔ・Ｇｒ
Ｐｂ＝（λ／４πｒb）2・Ｐｔ・Ｇｔ・Ｇｒ
Ｐｃ＝（λ／４πｒc）2・Ｐｔ・Ｇｔ・Ｇｒ
Ｐｄ＝（λ／４πｒd）2・Ｐｔ・Ｇｔ・Ｇｒ
Ｐｅ＝（λ／４πｒe）2・Ｐｔ・Ｇｔ・Ｇｒ
【００７２】
この数式と、パラメータβｂ〜βｅの定義式とを用いると、以下の数１７が得られる。
（数１７）
βb＝Ｐａ／Ｐｂ
βc＝Ｐａ／Ｐｃ
βd＝Ｐａ／Ｐｄ
βe＝Ｐａ／Ｐｅ
【００７３】
数１７から明らかなように、βb〜βeは、受信機２０−Ａ〜２０−Ｅの各受信強度を用いて具体的に求めることができる。
【００７４】
また、数１１の右辺は、基準位置である原点Ｏから、各受信機までの距離の二乗である。各受信機の位置は予め定められている。よって、原点Ｏからの距離の二乗もまた予め取得しておくことができる。数１１の右辺のベクトルは、および数１４に含まれる。
【００７５】
以上のとおり、数１４の右辺に示される、行列Ｑの逆行列Ｑ-1は演算によって求めることができる。また、原点位置から各受信機までの距離の二乗を並べたベクトルｓの各成分も予め求めることができる。これにより、数１４の左辺の成分ｘ、ｙ、ｚ、つまり波源の位置（ｘ，ｙ，ｚ）を推定することができる。
【００７６】
測位装置３０のＣＰＵ３１は、上述の原理に従った演算を実行する。以下、図１１および図１２を参照しながら、ＣＰＵ３１の動作を説明する。
【００７７】
図１１は、ＣＰＵ３１の処理ブロックの構成を示す。また図１２は、ＣＰＵ３１の処理手順を示すフローチャートである。
【００７８】
ＣＰＵ３１は、パラメータ演算器５１と、逆行列演算器５２と、ベクトル乗算器５３として機能する。図１１では、あたかも３つの構成要素が存在するように示されているが、これらは処理の単位を意味しているに過ぎない。図示された数の構成要素が含まれている必要はない。本実施形態では、ＣＰＵ３１は、図１２に示すフローチャートによる処理を行うコンピュータプログラムにしたがって動作する。ＣＰＵ３１は、コンピュータプログラムの命令に従い、時間によって、パラメータ演算器５１として動作し、逆行列演算器５２として動作し、ベクトル乗算器５３として動作する。処理ブロックから処理ブロックへの矢印は、データが次の演算に利用されることを意味している。
【００７９】
ただし、パラメータ演算器５１、逆行列演算器５２およびベクトル乗算器５３の少なくとも１つをハードウェアによって実現してもよい。たとえばＦＰＧＡ（Field-programmable gate array）を用いることにより、パラメータ演算器５１、逆行列演算器５２およびベクトル乗算器５３を１つの集積回路上に実装することが可能である。
【００８０】
以下では説明の便宜のため、パラメータ演算器５１、逆行列演算器５２およびベクトル乗算器５３が独立に設けられた構成要素であるとして説明する。
【００８１】
パラメータ演算器５１は、受信機２０−Ａ〜２０−Ｅの各々から受信強度の値（受信強度データ）を受け取り、数１７に示す演算を行うことにより、パラメータβb〜βeを算出する（ステップＳ１）。算出した各パラメータは逆行列演算器５２に送られる。
【００８２】
逆行列演算器５２は、パラメータ演算器５１からパラメータβb〜βeを受け取る。また逆行列演算器５２は、メモリ３２に格納されている各受信機の位置座標を示すデータを読み出す（ステップＳ２）。そして逆行列演算器５２は、パラメータβb〜βeおよび各受信機の位置座標を示すデータを用いて行列Ｑの逆行列を計算する（ステップＳ３）。逆行列の演算方法は公知であり、たとえば掃き出し法(sweep out method)を用いて計算することができる。なお、行列Ｑを入力すると逆行列Ｑ-1を出力するコンピュータプログラムのライブラリもまた公知であり容易に入手可能である。そのようなライブラリプログラムを用いてＣＰＵ３１を逆行列演算器５２として動作させてもよい。
【００８３】
ベクトル乗算器５３は、逆行列演算器５２によって計算された逆行列Ｑ-1およびメモリ３２に格納されている各受信機の位置座標を示すデータを受け取る。後者は、原点位置から各受信機までの距離の二乗を並べたベクトルｓを求めるために利用される。受け取ったデータに基づいて、ベクトル乗算器５３は、演算Ｑ-1・ｓを行う（ステップＳ４）。これにより、数１４左辺のベクトルｐを求めることができる。ベクトル乗算器５３は、求めたベクトルｐの成分であるｘ，ｙ，ｚを、波源の位置として出力する（ステップＳ５）。
【００８４】
以上の処理により、波源である信号発生器１０の位置を推定することができる。
【００８５】
上述の数１０に示す操作により、結果的に線形項が増えた。つまり、２つの非線形項を線形化したことにより、波源の位置を示すｘ、ｙ、ｚに加えて未知数が２つ増加した。５つの未知数を求めるためには、５つの連立方程式が必要となる。これが、５台の受信機からの信号を利用して、５行５列の行列Ｑが必要になった理由である。逆行列を求めるために、本願発明者は受信機の数を増やして行列のランクの減少を補った。
【００８６】
アルゴリズムを線形化することにより、ニュートン法を利用する場合と比較して、多重解を回避でき、繰り返し演算を行う必要もないため、ハードウェアへの負荷やソフトウェアの処理を軽減することができる。ＢＬＥ測位システムの導入は比較的手軽に行うことができる。そのような性格上、導入のコストを抑えることが求められる。線形計測処理は、入力に対して出力が一意的に決定するためテーブル化も可能でありシステムの低コスト化に貢献する。
【００８７】
「テーブル化」とは、予め測定されたβb、βc、βd、βeの組と、それらから推定されるｘ、ｙ、ｚの組とを対応付けたテーブルを用意することを言う。実際に受信強度データが得られ、βb、βc、βd、βeの組が決定されたとき、テーブルを参照して一致する組を検索する。一致する組がテーブルに存在した場合には、対応付けられたｘ、ｙ、ｚの組を読み出して出力する。行列演算等が不要になるため、非常に高速に結果を出力することができる。さらにＣＰＵ３１の負荷を大きく抑制できる。テーブルのエントリ数が増えるほど、実測したβb、βc、βd、βeの組に一致する組がテーブル上に存在する可能性が高くなる。
【００８８】
なお、テーブルのエントリ数を増加させることは極めて容易である。測位システム１が導入された環境で、波源となる信号発生器１０の位置を変えながら、都度、各受信機で信号波を受信してβb、βc、βd、βeの組を取得し、波源の位置ｘ、ｙ、ｚの組を推定しておけばよいからである。検索のヒット率が十分高いテーブルを設けておけば、比較的処理能力が低いＣＰＵを採用しても、十分高速に波源の位置を示す情報を出力することができる。
【００８９】
上述の説明から類推されるように、二次元空間の位置を推定する場合には、４台の受信機から受信強度データを所得できればよい。二次元空間の位置を推定するために必要な行列Ｑ、ベクトルｐおよびｓは以下のとおりである。
【数１８】
2018062433.tif 000013
【００９０】
なお、行列Ｑの逆行列が存在するための受信機の位置に関する条件は同じである。
【００９１】
信号発生器１０の位置を推定した後は、測位システム１は、たとえば図示されない誘導装置を利用して、自走車両２４ｂおよび／またはスマートフォン２４ｃを所定の位置まで誘導してもよい。誘導装置は移動先の位置を通知して、当該位置まで自走車両２４ｂまたはスマートフォン２４ｃを所有する人を誘導することができる。このとき、リアルタイムに、または間欠的に信号発生器１０の位置の推定処理を行い、適切に誘導を行うことができる。
【００９２】
なお、上述の説明では、二次元空間の位置を推定するためには４台の受信機からの受信強度データを利用すると説明したが、５台以上の受信機からの受信強度データを利用してもよい。また、三次元空間の位置を推定するためには５台の受信機からの受信強度データを利用すると説明したが、６台以上の受信機からの受信強度データを利用してもよい。
【００９３】
（応用例）
以下、応用例を説明する。本応用例では、たとえば三次元空間において６台以上の受信機から受信強度データが得られる例を想定して、数１１および数１２に示す行列演算を利用して波源の位置を推定する。ただし、受信機が５台で５個の受信強度データが得られる例は排除されない。
【００９４】
数１２に示される行列演算  Ｑ・ｐ＝ｓ  の行列Ｑ、ベクトルｐおよびｓは以下のとおりである。
【数１９】
2018062433.tif 000014
【００９５】
いま、Ｎ台の受信機の受信強度データから重複しない５個のデータを抽出して１つの組とし、それぞれの組について、図１２に示す処理を行うとする。５個のデータを１組として第ｉ組のデータから得られる行列Ｑおよびｓをそれぞれ、Ｑiおよびベクトルｓiと記述する。すると、Ｑi・ｐ＝ｓiが成り立つ。
【００９６】
さらに、第ｉ組のデータを得た時の計測状況などによって生じる誤差成分のベクトルをεiと表すと、以下の式が成り立つ。なお、誤差成分εiは、ベクトルｐおよびｓと同じ５行１列の縦ベクトルである。
（数２０）
Ｑi・ｐ＋εi＝ｓi
【００９７】
数２０を変形すると数２１が得られる。
（数２１）
εi＝ｓi−Ｑi・ｐ
【００９８】
誤差成分εiの二乗誤差はεTiεiと表すことができる。すべての組について二乗誤差が得られる。二乗誤差の総和が最も小さい値、すなわち０を取ると仮定して、その状況下でベクトルｐを求めることを考える。
【数２２】
2018062433.tif 000015
【００９９】
数２２に数２１を代入すると数２３が得られる。
【数２３】
2018062433.tif 000016
【０１００】
数２３の両辺をｐで微分すると数２４が得られ、さらに変形すると数２５が得られる。
【数２４】
2018062433.tif 000017
【数２５】
2018062433.tif 000018
【０１０１】
数２５から数２６が導出され、その結果、ベクトルｐを求めるための数２７が得られる。
【数２６】
2018062433.tif 000019
【数２７】
2018062433.tif 000020
【０１０２】
数２７からベクトルｐが得られ、ベクトルｐに含まれるｘ，ｙ，ｚの各値が、ベクトルｐの成分であるｘ，ｙ，ｚを、波源の位置として得ることができる。
【０１０３】
図１３は、数２７の演算を行うＣＰＵ３１の処理ブロックの構成を示す。また図１４は、ＣＰＵ３１の処理手順を示すフローチャートである。なお、図１３の処理のうち、ステップＳ２以降がＣＰＵ３１の処理に相当する。
【０１０４】
ＣＰＵ３１は、組合せパターン発生器６１と、ベクトル演算器６２−１〜６２−Ｋと、加算器６３−１および６３−２と、逆行列演算器６４と、ベクトル乗算器６５として機能する。図１１の例と同様、図１３の構成はＣＰＵ３１が行う処理の単位を意味している。本応用例でも、ＣＰＵ３１は、図１４に示すフローチャートによる処理を行うコンピュータプログラムにしたがって動作する。なお図１３のブロックの少なくとも１つをハードウェアによって実現してもよいことも図１１と同様である。
【０１０５】
測位装置３０は、Ｎ台の受信機の各々から出力された受信強度データを、通信回路３３（図９）を介して受信し、メモリ３２に格納する（ステップＳ１１）。
【０１０６】
組合せパターン発生器６１は、受け取ったＮ個のデータから５個１組のデータ群をＫ組抽出する（ステップＳ１２）。組数Ｋは組み合わせの数であり、NＣ5通りである。たとえばＮ＝９の場合には１２６通りのデータの組が得られる。組合せパターン発生器６１は、組み合わせごとにデータ群を読み出してベクトル演算器６２−１〜６２−Ｋの各々に入力する。
【０１０７】
ベクトル演算器６２−１〜６２−Ｋの各々は、メモリ３２から５個の受信強度データと、各受信強度データを出力した各受信機と原点Ｏとの距離の２乗の値を受け取る。これにより、数１９に示す行列Ｑおよびベクトルｓの各要素を得ることができる。
【０１０８】
ベクトル演算器６２は２種類の演算を行う（ステップＳ１３）。１つは行列演算ＱiTＱiであり、他の１つはベクトル演算Ｑiｓiである。行列演算ＱiTＱiの結果は加算器６３−１に送られる。ベクトル演算Ｑiｓiの結果は加算器６３−２に送られる。
【０１０９】
加算器６３−１は、ベクトル演算器６２−１〜６２−Ｋの各々から受け取った行列を加算して、数２７の右辺の括弧内の行列、つまり逆行列算出の対象となる行列、を求める（ステップＳ１４）。加算器６３−１は、求めた行列を逆行列演算器６４に送る。
【０１１０】
一方、加算器６３−２は、ベクトル演算器６２−１〜６２−Ｋの各々から受け取ったベクトルを加算して、和ベクトルを求める（ステップＳ１５）。得られた和ベクトルは、数２７の右辺の２つ目の括弧内のベクトルに対応する。加算器６３−１は、求めた和ベクトルをベクトル乗算器６５に送る。
【０１１１】
逆行列演算器６４は、受け取った行列の逆行列を求める（ステップＳ１６）。逆行列の演算方法は公知であり、たとえば掃き出し法(sweep out method)を用いて計算することができる。なお、ある行列を入力するとその逆行列を出力するコンピュータプログラムのライブラリもまた公知であり容易に入手可能である。そのようなライブラリプログラムを用いてＣＰＵ３１を逆行列演算器６２として動作させてもよい。
【０１１２】
ベクトル乗算器６３は、逆行列演算器６２によって計算された逆行列および加算器６３−２によって計算された和ベクトルを受け取り、その積を求める（ステップＳ１７）。この演算は、数２７に示すベクトル演算に対応する。これにより、数２７左辺のベクトルｐを求めることができる。ベクトル乗算器６３は、求めたベクトルｐの成分であるｘ，ｙ，ｚを、波源の位置として出力する（ステップＳ１８）。
【０１１３】
以上の処理により、波源である信号発生器１０の位置を推定することができる。
【０１１４】
上述した実施形態（応用例を含む）では、ベクトルｐおよびｓを列ベクトルとしたときの行列Ｑを説明した。しかしながら、ベクトルｐおよびｓを行ベクトルとして行列Ｑを構築してもよい。転置の関係にあれば数学的には同じであることは明らかである。
【産業上の利用可能性】
【０１１５】
本開示の誘導システムは、屋内または屋外を移動する移動体の位置の推定に用いられ得る。また測位した結果を利用する移動体の位置の制御に用いられ得る。
【符号の説明】
【０１１６】
１  測位システム
１０  信号発生器
２０−Ａ、２０−Ｂ、・・・、２０−Ｎ  受信機
３０  測位装置
３１  ＣＰＵ（演算回路）
３２  メモリ
３３  通信回路
３４ａ  入力インタフェース
３４ｂ  出力インタフェース
５１  パラメータ演算器
５２  逆行列演算器
５３  ベクトル乗算器
【図１】
2018062433.tif 000023
【図２】
2018062433.tif 000024
【図３】
2018062433.tif 000025
【図４】
2018062433.tif 000026
【図５】
2018062433.tif 000027
【図６】
2018062433.tif 000028
【図７】
2018062433.tif 000029
【図８】
2018062433.tif 000030
【図９】
2018062433.tif 000031
【図１０】
2018062433.tif 000032
【図１１】
2018062433.tif 000033
【図１２】
2018062433.tif 000034
【図１３】
2018062433.tif 000035
【図１４】
2018062433.tif 000036
【国際調査報告】
2018062433.tif 000037
2018062433.tif 000038