(57)【要約】      （修正有）【課題】衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して３次元の電離層電子密度分布をより正確に推定する。【解決手段】アレイアンテナが受信した電波の到来方位と仰角とを算出する第１算出部と、衛星信号の通過経路の総電子数を算出する総電子数算出部と、電子密度を推定する３次元空間を設定し当該３次元空間を複数の領域に分割する３次元電子密度推定空間設定部３６と、領域毎の電子密度を推定する電子密度推定部３８と、アレイアンテナが受信した電波のホップ数を算出するホップ数算出部３９と、ホップ数に基づいて算出したアレイアンテナから送信源までの距離と実際の距離との差が所定のしきい値以下である場合に電子密度推定部による推定結果を正しいと判定する判定部４１と、判定部４１の判定結果に応じて電波の伝搬経路における電子密度を修正する修正部４３とを備える。【選択図】図２
2011185665.tif 000002
【特許請求の範囲】
【請求項１】
複数の測位衛星から送信される衛星信号を受信する１以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定装置であって、
既知の送信源により送信される電波を受信するアレイアンテナと、
前記アレイアンテナが受信した電波の到来方位と仰角とを算出する第１算出部と、
前記１以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出部と、
電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する空間設定部と、
前記総電子数算出部により算出された総電子数に基づいて前記空間設定部により分割された領域毎の電子密度を推定する電子密度推定部と、
前記第１算出部により算出された電波の到来方位及び仰角と前記電子密度推定部により推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第２算出部と、
前記第２算出部により算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定部による推定結果を正しいと判定する判定部と、
前記判定部により前記電子密度推定部の推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正部と、
を備えることを特徴とする電子密度推定装置。
【請求項２】
前記電子密度推定部は、前記修正部により修正された電子密度に基づいて再度前記領域毎の電子密度を推定することを特徴とする請求項１記載の電子密度推定装置。
【請求項３】
前記修正部は、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する際に、前記伝搬経路としてレイトレーシング手法により特定されたものを用いることを特徴とする請求項１又は請求項２記載の電子密度推定装置。
【請求項４】
複数の測位衛星から送信される衛星信号を受信する１以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定方法であって、
アレイアンテナで既知の送信源により送信される電波を受信し、受信した電波の到来方位と仰角とを算出する第１算出ステップと、
前記１以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出ステップと、
電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する空間設定ステップと、
前記総電子数算出ステップにより算出された総電子数に基づいて前記空間設定ステップにより分割された領域毎の電子密度を推定する電子密度推定ステップと、
前記第１算出ステップにより算出された電波の到来方位及び仰角と前記電子密度推定ステップにより推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第２算出ステップと、
前記第２算出ステップにより算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定ステップによる推定結果を正しいと判定する判定ステップと、
前記判定ステップにより前記電子密度推定ステップの推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正ステップと、
を備えることを特徴とする電子密度推定方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、衛星から送信される複数の異なる周波数信号による通過経路の電離層電子密度分布を推定する電子密度推定装置及び電子密度推定方法に関する。
【背景技術】
【０００２】
衛星ナビゲーションシステム「Ｇａｌｉｌｅｏ」やＧＰＳ（Ｇｌｏｂａｌ  Ｐｏｓｉｔｉｏｎｉｎｇ  Ｓｙｓｔｅｍ）、ＱＺＳＳ（Ｑｕａｓｉ  Ｚｅｎｉｔｈ  Ｓａｔｅｌｌｉｔｅｓ  Ｓｙｓｔｅｍ；準天頂衛星システム）等の衛星から送信される信号は、地球上の観測点の測位を行うのに利用される。これらの衛星により放送される複数周波数信号の擬似距離差は、電波の伝搬経路上の電子数に関係していることが知られている。
【０００３】
衛星から地上に向けて送信される信号は、電離層を通過する際に伝搬遅延量を生じる。この電離層伝搬遅延量は、伝搬する信号の周波数に依存している。そのため、衛星からの複数の周波数信号を受信することにより、信号の通過経路における総電子数（ＴＥＣ：Ｔｏｔａｌ  Ｅｌｅｃｔｒｏｎ  Ｃｏｎｔｅｎｔ）を得ることができる。電波の伝搬経路上の総電子数を推定する場合には、複数周波数信号の送信・受信時刻による擬似距離（コード距離）差を使用する方法や複数周波数のキャリア位相擬似距離（フェーズ距離）の差を利用する方法が知られている。
【０００４】
非特許文献１や非特許文献２には、電離層における電子密度モデル関数であるＩＲＩ（Ｉｎｔｅｒｎａｔｉｏｎａｌ  Ｒｅｆｅｒｅｎｃｅ  Ｉｏｎｏｓｐｈｅｒｅ）モデルが記載されている。さらに、非特許文献４には、やはり電離層電子密度モデル関数の一つであるＧａｌｌａｇｈｅｒのモデルが記載されている。これらのモデル関数を用いて計算によりＴＥＣを求めることもできる。
【０００５】
非特許文献５や非特許文献６には、トモグラフィ手法による衛星ＴＥＣを用いた３次元電子密度推定が記載されており、ＭＡＲＴ（Ｍｕｌｔｉｐｌｉｃａｔｉｖｅ  Ａｌｇｅｂｒａｉｃ  Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ  Ｔｅｃｈｎｉｑｕｅ）と呼ばれる代数的再構成法が適用されている。一方、非特許文献７には、ＧＳＶＤ（Ｇｅｎｅｒａｌｉｚｅｄ  Ｓｉｎｇｕｌａｒ  Ｖａｌｕｅ  Ｄｅｃｏｍｐｏｓｉｔｉｏｎ）を用いたトモグラフィック再構成による電子密度推定が記載されている。
【０００６】
また、非特許文献８及び非特許文献９には、電離層を球核の薄い層とみなした電子密度計算について記載されており、ＧＰＳを使用した観測結果のみに基づいて３次元電子密度推定を試みている。さらに、非特許文献１０には、３次元方向の電子密度分布を高さ方向に計算した場合について記載されている。
【０００７】
ＨＦ帯を用いた電波通信において、受信した電波の到来方法及び仰角を測定し、この到来方法及び仰角と、ＩＲＩモデル等の電離層電子密度分布モデル等とを用いて、レイトレーシング手法等により電波の伝搬経路を計算することができる。このような手法については、非特許文献１１，１２，１３に記載されている。また、非特許文献１４には、電波の伝搬特性を予測する手法が記載されている。
【先行技術文献】
【非特許文献】
【０００８】
【非特許文献１】Ｄｉｅｔｅｒ  Ｂｉｌｉｔｚａ：“Ｉｎｔｅｒｎａｔｉｏｎａｌ  Ｒｅｆｅｒｅｎｃｅ  Ｉｏｎｏｓｐｈｅｒｅ  ２０００”，Ｒａｄｉｏ．Ｓｃｉｅｎｃｅ，Ｖｏｌ．３６，Ｎｕｍｂｅｒ２，ＰＰ２６１−２７５，Ｍａｒｃｈ／Ａｐｒｉｌ，２００１
【非特許文献２】Ｄｉｅｔｅｒ  Ｂｉｌｉｔｚａ，ｅｔ，ａｌ．，：“Ｉｎｔｅｒｎａｔｉｏｎａｌ  Ｒｅｆｅｒｅｎｃｅ  Ｉｏｎｏｓｐｈｅｒｅ  １９９０”，Ｎｏｖｅｍｂｅｒ，１９９０．
【非特許文献３】Ａ．Ｋｏｍｊａｔｈｙ：‘Ｇｌｏｂａｌ  Ｉｏｎｏｓｐｈｅｒｉｃ  Ｔｏｔａｌ  Ｅｌｅｃｔｒｏｎ  Ｃｏｎｔｅｎｔ  Ｍａｐｐｉｎｇ  Ｕｓｉｎｇ  ｔｈｅ  Ｇｌｏｂａｌ  Ｐｏｓｉｔｉｏｎｉｎｇ  Ｓｙｓｔｅｍ’  ，ＵＮＢ，Ｔｅｃｈｎｉｃａｌ  Ｒｅｐｏｒｔ  Ｎｏ．１８８，Ｓｅｐ．１９９７．
【非特許文献４】Ｇａｌｌａｇｈｅｒ，Ｄ．Ｌ．，Ｐ．Ｄ．Ｃｒａｖｅｎ，ａｎｄ  Ｒ．Ｈ．Ｃｏｍｆｏｒｔ，Ｇｌｏｂａｌ  ｃｏｒｅ  ｐｌａｓｍａ  ｍｏｄｅｌ，Ｊ．Ｇｅｏｐｈｙｓ．Ｒｅｓ．１０５，Ａ８，１８，８１９−１８，８３３，２０００．
【非特許文献５】Ｒａｙｍｕｎｄ，Ｔ．Ｄ．，Ａｕｓｔｉｎ，Ｊ．Ｒ．，Ｆｒａｎｋｅ，Ｓ．Ｊ．，Ｌｉｕ，Ｃ．Ｈ．，Ｋｌｏｂｕｃｈａｒ，Ｊ．Ａ．，ａｎｄ  Ｓｔａｌｋｅｒ，Ｊ．，：“Ａｐｐｌｉｃａｔｉｏｎ  ｏｆ  Ｃｏｍｕｔｅｒｉｚｅｄ  Ｔｏｍｏｇｒａｐｈｙ  ｔｏ  ｔｈｅ  Ｉｎｖｅｓｔｉｇａｔｉｏｎ  ｏｆ  Ｉｏｎｏｓｐｈｅｒｉｃ  Ｓｔｒｕｃｔｕｒｅｓ”，Ｒａｄｉｏ  Ｓｃｉ．，２５（５），７７１−７８９，１９９０．
【非特許文献６】Ｓｔｅｆａｎ  Ｓｃｈｌｕｅｒ，ｅｔ  ａｌ．：“Ｍｏｎｉｔｏｒｉｎｇ  ｔｈｅ  ３  Ｄｉｍｅｎｓｉｏｎａｌ  Ｉｏｎｏｓｐｈｅｒｉｃ  Ｅｌｅｃｔｒｏｎ  Ｄｉｓｔｒｉｂｕｔｉｏｎ  ｂａｓｅｄ  ｏｎ  ＧＰＳ  Ｍｅａｓｕｒｅｍｅｎｔｓ”
【非特許文献７】Ｋ．Ｂｈｕｙａｎ，ｅｔ  ａｌ：“Ｔｏｍｏｇｒａｐｈｉｃｒｅｃｏｎｓｔｒｕｃｔｉｏｎ  ｏｆ  ｔｈｅ  Ｉｏｎｏｓｐｈｅｒｅ  ｕｓｉｎｇ  ｇｅｎｅｒａｌｉｚｅｄ  ｓｉｎｇｕｌａｒ  ｖａｌｕｅ  ｄｅｃｏｍｐｏｓｉｔｉｏｎ”  ，Ｃｕｒｒｅｎｔ  Ｓｃｉ．，Ｖｏｌ１．８３，Ｎｏ．９，１０  Ｎｏｖ．，２００２．
【非特許文献８】Ｏｔｕｋａ，Ｙ．，ｅｔ．ａｌ．：“Ａ  ｎｅｗ  ｔｅｃｈｎｉｑｕｅ  ｆｏｒ  ｍａｐｐｉｎｇ  ｏｆ  ｔｏｔａｌ  ｅｌｅｃｔｒｏｎ  ｃｏｎｔｅｎｔ  ｕｓｉｎｇ  ＧＰＳ  ｎｅｔｗｏｒｋ  ｉｎ  Ｊａｐａｎ，”，Ｅａｒｔｈ  Ｐｌａｎｅｔｓ  Ｓｐａｃｅ，１１１−１２０，２００１．
【非特許文献９】Ｍａ，Ｇ．，ａｎｄ  Ｔ．Ｍａｒｕｙａｍａ（２００３），Ｄｅｒｉｖａｔｉｏｎ  ｏｆ  ＴＥＣ  ａｎｄ  ｅｓｔｉｍａｔｉｏｎ  ｏｆ  ｉｎｓｔｒｕｍｅｎｔａｌ  ｂｉａｓｅｓ  ｆｒｏｍ  ＧＥＯＮＥＴ  ｉｎ  Ｊａｐａｎ，Ａｎｎ．Ｇｅｏｐｈｙｓ．，２１，２０８３−２０９３．
【非特許文献１０】五十嵐来良，斉藤昭則，大塚雄一：“標準電離層モデル（ＩＲＩ）で導出したＴＥＣとＧＰＳ受信機網で観測したＴＥＣの比較解析”，第１回電離圏の利用と影響に関するシンポジウム講和集，ｐｐ２４−１〜８，２００３．
【非特許文献１１】前田憲一，後藤三男：“電波伝搬”，岩波全書，１９５３年２月
【非特許文献１２】Ｋ．Ｇ．Ｂｕｄｅｎ：“Ｔｈｅ  Ｐｒｏｐａｇａｔｉｏｎ  ｏｆ  ｒａｄｉｏ  ｗａｖｅｓ”．Ｃａｍｂｒｉｄｇｅ  Ｕｎｉｖｅｒｓｉｔｙ  Ｐｒｅｓｓ，１９８８．
【非特許文献１３】Ｉｗａｎｅ  Ｋｉｍｕｒａ：“Ｅｆｆｅｃｔｓ  ｏｆ  Ｉｏｎｓ  ｏｎ  Ｗｈｉｓｔｌｅｒ−Ｍｏｄｅ  Ｒａｙ  Ｔｒａｃｉｎｇ”，Ｒａｄｉｏ  Ｓｃｉｅｎｃｅ．Ｖｏｌ．Ｉ，Ｎｏ．３，２６９−２８３，Ｍａｒｃｈ  １９６６．
【非特許文献１４】“ＨＦ  ｐｒｏｐａｇａｔｉｏｎ  ｐｒｅｄｉｃｔｉｏｎ  ｍｅｔｈｏｄ”，ＩＴＵ−Ｒ  Ｐ．５３３−７．
【発明の概要】
【発明が解決しようとする課題】
【０００９】
しかしながら、地上で観測するＧＰＳ衛星電波の通過経路が限定されるので、非特許文献８及び非特許文献９に記載された手法のようにＧＰＳを使用した観測結果のみに基づいて３次元電子密度推定を行う場合には、十分な３次元空間の情報を得ることができないという問題点がある。
【００１０】
すなわち、従来の衛星観測結果を用いた技術は、ＧＰＳ衛星から受信機までの通過経路上の総電子数しか求められないので、観測できる３次元空間が限定されてしまい、３次元の電離層電子密度分布を十分推定しきれていない。その結果、電波が通過していない空間の３次元電子密度推定値は、誤差が多く含まれるものとなる可能性が大きい。このような問題が生じる根本的な要因は、観測値が少ない点にあると考えられる。
【００１１】
本発明は上述した従来技術の問題点を解決するもので、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して３次元の電離層電子密度分布をより正確に推定する電子密度推定装置及び電子密度推定方法を提供することを課題とする。
【課題を解決するための手段】
【００１２】
本発明に係る電子密度推定装置は、上記課題を解決するために、複数の測位衛星から送信される衛星信号を受信する１以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定装置であって、既知の送信源により送信される電波を受信するアレイアンテナと、前記アレイアンテナが受信した電波の到来方位と仰角とを算出する第１算出部と、前記１以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出部と、電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する空間設定部と、前記総電子数算出部により算出された総電子数に基づいて前記空間設定部により分割された領域毎の電子密度を推定する電子密度推定部と、前記第１算出部により算出された電波の到来方位及び仰角と前記電子密度推定部により推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第２算出部と、前記第２算出部により算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定部による推定結果を正しいと判定する判定部と、前記判定部により前記電子密度推定部の推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正部とを備えることを特徴とする。
【００１３】
本発明に係る電子密度推定方法は、上記課題を解決するために、複数の測位衛星から送信される衛星信号を受信する１以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定方法であって、アレイアンテナで既知の送信源により送信される電波を受信し、受信した電波の到来方位と仰角とを算出する第１算出ステップと、前記１以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出ステップと、電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する空間設定ステップと、前記総電子数算出ステップにより算出された総電子数に基づいて前記空間設定ステップにより分割された領域毎の電子密度を推定する電子密度推定ステップと、前記第１算出ステップにより算出された電波の到来方位及び仰角と前記電子密度推定ステップにより推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第２算出ステップと、前記第２算出ステップにより算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定ステップによる推定結果を正しいと判定する判定ステップと、前記判定ステップにより前記電子密度推定ステップの推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正ステップとを備えることを特徴とする。
【発明の効果】
【００１４】
本発明によれば、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して３次元の電離層電子密度分布をより正確に推定する電子密度推定装置及び電子密度推定方法を提供することができる。
【図面の簡単な説明】
【００１５】
【図１】本発明の実施例１の形態の電子密度推定装置及び測位システムの構成を示すブロック図である。
【図２】本発明の実施例１の形態の電子密度推定装置の詳細なブロック図である。
【図３】本発明の実施例１の形態の電子密度推定装置の動作を示すフローチャート図である。
【図４】本発明の実施例１の形態の電子密度推定装置内の３次元電子密度推定空間設定部による３次元空間設定を説明する図である。
【図５】本発明の実施例１の形態の電子密度推定装置内のホップ数算出部によるホップ数の算出方法を説明する図である。
【発明を実施するための形態】
【００１６】
以下、本発明の電子密度推定装置及び電子密度推定方法の実施の形態について図面を参照して詳細に説明する。
【実施例１】
【００１７】
以下、本発明の実施例について図面を参照しながら説明する。図１は、本発明の実施例１の電子密度推定装置３２及び測位システムの構成を示すブロック図である。この電子密度推定装置３２は、複数の測位衛星（測位衛星１０ａ、１０ｂ、１０ｃ）から送信される衛星信号を受信する１以上の衛星信号受信機２４を有し且つ衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用して電子密度を推定する。
【００１８】
まず、本実施の形態の構成を説明する。本実施例の電子密度推定装置３２及び測位システムは、図１に示すように、信号処理部１、測位衛星１０ａ、測位衛星１０ｂ、測位衛星１０ｃ、衛星信号処理系２０、電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０、及びインターネットデータ処理系４０で構成されている。ただし、本実施例において信号処理部１は、電子密度推定装置の構成の一部であるものとする。
【００１９】
測位衛星１０ａ、測位衛星１０ｂ、測位衛星１０ｃ、及び衛星信号処理系２０は、本発明の測位システムに対応する。なお、インターネットデータ処理系４０も含めて測位システムとすることも可能であるが、インターネットデータ処理系４０は、必ずしも必須のものではなく、付加的なものである。また、信号処理部１、衛星信号処理系２０、電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０、及びインターネットデータ処理系４０は、お互いに通信回線５０で接続されている。
【００２０】
測位衛星１０ａ、測位衛星１０ｂ、測位衛星１０ｃは、ＧＰＳ、Ｇａｌｉｌｅｏ、準天頂衛星（ＱＺＳＳ）等の航法衛星であり、複数の周波数の衛星信号を送信する。
【００２１】
衛星信号処理系２０は、アンテナ２２、衛星信号受信機２４、及び衛星信号処理装置２６から構成される。衛星信号受信機２４は、本発明の受信機に対応し、複数の測位衛星（測位衛星１０ａ、１０ｂ、１０ｃ）から送信される複数の周波数の衛星信号をアンテナ２２を介して受信する。また、衛星信号処理装置２６は、衛星信号に含まれる測位情報を用いて位置情報を得る。
【００２２】
信号処理部１は、複数のアンテナ１１−１，１１−２，…，１１−ｎからなるアレイアンテナ１１と、アンプ１２−１，１２−２，…，１２−ｎと、周波数変換部１３−１，１３−２，…，１３−ｎと、Ａ／Ｄ変換部１４−１，１４−２，…，１４−ｎと、演算部１５とを備える。
【００２３】
アレイアンテナ１１は、既知の送信源により送信される電波を受信する。例えば、アレイアンテナ１１は、既存の放送局から送信され、電離層を通過して到来する放送信号（電波）を受信して受信信号を出力する。アンプ１２−１〜１２−ｎは、それぞれアンテナ１１−１〜１１−ｎからの信号を増幅する周波数変換部１３−１〜１３−ｎは、アンプ１２−１〜１２−ｎにより増幅された信号をベースバンドの信号に変換する。Ａ／Ｄ変換部１４−１〜１４−ｎは、周波数変換部１３−１〜１３−ｎにより変換され出力された信号をＡ／Ｄ変換し、デジタル受信信号として出力する。
【００２４】
演算部１５は、本発明の第１算出部に対応し、Ａ／Ｄ変換部１４−１〜１４−ｎの各々により出力されたデジタル受信信号に基づいて、アレイアンテナ１１が受信した電波の到来方位と仰角とを算出する。例えば、演算部１５は、ＭＵＳＩＣ（Ｍｕｌｔｉｐｌｅ  Ｓｉｇｎａｌ  Ｃｌａｓｓｉｆｉｃａｔｉｏｎ）法や独立成分分析(ＩＣＡ：Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ)の手法などを用いて、既存の放送局からの電波の到来方位角および到来仰角を算出する。
【００２５】
電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０は、電離層電子密度分布モデルと演算部１５で算出された方位角及び仰角とを用いて、既存の放送局からの電波の伝搬経路をレイトレーシング手法により算出し、算出した伝搬経路に基づいて、電波の到来方位角および到来仰角を算出する。なお、レイトレーシング手法の詳細については、例えば非特許文献１２に記載されている。
【００２６】
あるいは、電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０は、電子密度推定装置３２を有し、受信した信号に基づき衛星位置、信号通過経路の電離層総電子数（ＴＥＣ）、周波数間バイアス等を算出する。受信したデータ・処理結果は、データサーバ系へＬＡＮ経由で伝送され保存される。
【００２７】
インターネットデータ処理系４０は、ルータ４２、衛星データ処理装置４４、及び外部インターネット網４６で構成されている。ルータ４２は、スイッチングハブでもよい。インターネットデータ処理系４０は、公開されている電離層関連の情報や国土地理院が公開しているＧＰＳ観測データ（ＧＥＯＮＥＴデータ）や国際的に観測結果を公開しているＩＧＳ（Ｉｎｔｅｒｎａｔｉｏｎａｌ  ＧＰＳ  Ｓｅｒｖｉｃｅ  ｆｏｒ  Ｇｅｏｄｙｎａｍｉｃｓ）データ等をインターネット経由で収集する装置である。インターネットデータ処理系４０により収集された結果は、電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０で処理される。ルータ４２は、セキュリティを考慮して設けられ、ファイアウォールとする。
【００２８】
図２は、本発明の実施例１の電子密度推定装置３２の詳細なブロック図である。図２に示すように、電子密度推定装置３２は、衛星ＴＥＣ算出部３３、ＴＥＣバイアス補正部３４、周波数間バイアス推定部３５、３次元電子密度推定空間設定部３６、電子密度モデル算出部３７、電子密度推定部３８、ホップ数算出部３９、判定部４１、及び修正部４３で構成されている。
【００２９】
衛星ＴＥＣ算出部３３、ＴＥＣバイアス補正部３４、及び周波数間バイアス推定部３５は、本発明の総電子数算出部に対応し、１以上の受信機の各々から複数の測位衛星の各々に対する擬似距離及びキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する。具体的には、衛星ＴＥＣ算出部３３は、衛星信号受信機２４から測位衛星１０ａ、１０ｂ、１０ｃの各々に対する擬似距離及びキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する。
【００３０】
周波数間バイアス推定部３５は、周波数間バイアスを算出し、ＴＥＣバイアス補正部３４に出力する。周波数間バイアス推定部３５による周波数間バイアスの推定方法は、どのようなものでもよく、従来手法を用いてよい。
【００３１】
ＴＥＣバイアス補正部３４は、周波数間バイアス推定部３５により推定された周波数間バイアスを用いて補正を行い、最終的な総電子数を算出する。
【００３２】
３次元電子密度推定空間設定部３６は、本発明の空間設定部に対応し、電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する。
【００３３】
電子密度モデル算出部３７は、上述したようなＩＲＩモデルやＧａｌｌａｇｈｅｒモデル等の電離層電子密度モデル関数に基づき電子密度モデルを算出する。
【００３４】
電子密度推定部３８は、総電子数算出部（ＴＥＣバイアス補正部３４）により算出された総電子数に基づいて３次元電子密度推定空間設定部３６により分割された領域毎の電子密度を推定する。本実施例において、電子密度推定部３８は、非特許文献５や非特許文献に記載されているようなＭＡＲＴやＧＳＶＤ等のトモグラフィの手法を用いて電子密度を推定するものとする。
【００３５】
ホップ数算出部３９は、本発明の第２算出部に対応し、通信回線５０を介して信号処理部１内の演算部１５から情報を得ることにより、演算部１５により算出された電波の到来方位及び仰角と電子密度推定部３８により推定された電子密度とに基づいて、アレイアンテナ１１が受信した電波のホップ数を算出する。
【００３６】
判定部４１は、ホップ数算出部３９により算出されたホップ数に基づいて算出したアレイアンテナ１１から送信源（電波を送信した放送局）までの距離と実際の距離との差が所定のしきい値以下である場合に、電子密度推定部３８による推定結果を正しいと判定する。
【００３７】
修正部４３は、判定部４１により電子密度推定部３８の推定結果が正しいと判断されない場合に、アレイアンテナ１１が受信した電波の伝搬経路における電子密度を修正する。すなわち、ホップ数算出部３９により算出されたホップ数に基づいて判定部４１が算出したアレイアンテナ１１から送信源までの距離と実際の距離との差の大きさが所定のしきい値よりも大きい場合に、修正部４３は、アレイアンテナ１１が受信した電波の伝搬経路における電子密度を修正する。なお、修正部４３は、アレイアンテナ１１が受信した電波の伝搬経路における電子密度を修正する際に、伝搬経路としてレイトレーシング手法により特定された伝搬経路を用いる。
【００３８】
また、電子密度推定部３８は、修正部４３により修正された電子密度に基づいて再度領域毎の電子密度を推定する。
【００３９】
次に、上述のように構成された本実施の形態の作用を説明する。図３は、本実施例の電子密度推定装置３２の動作を示すフローチャート図である。まず、衛星ＴＥＣ算出部３３は、１以上の受信機の各々から複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する（ステップＳ１）。具体的には、衛星ＴＥＣ算出部３３は、広範囲に設置した衛星受信機収集データからＧＰＳ等の衛星のＴＥＣを計算する。ここで、本実施例の衛星ＴＥＣ算出部３３は、衛星から放送される２周波の信号（ＧＰＳ衛星の場合Ｌ１周波数（１５７５．４２ＭＨｚ）とＬ２周波数（１２２７．６０ＭＨｚ）の信号）に基づいて信号伝搬経路のＴＥＣを算出するものとする。ＧＰＳ以外でも異なる２周波を放送している場合には同様な手法が使用できる。
【００４０】
ここでは、測位衛星１０ａと衛星信号受信機２４との間の擬似距離について考える。衛星信号通過経路のＴＥＣは、測位衛星１０ａから送信される衛星信号（ここではＬ１とＬ２で示す）に基づき求められる。まず、擬似距離（コード距離、シュードレンジ）とキャリア位相擬似距離（フェーズ距離）は、以下のように表すことができる。
【数１】
2011185665.tif 000003
【００４１】
ここで、ρは、遅延時間による距離（擬似距離）を示す。また、Φは、キャリア位相を示す。ｒは、真の距離を示す。ｃは光速である。また、δｔｕは、受信機時刻誤差を示し、δｔｓは、衛星時刻誤差を示す。なお、本実施例において、右下の添え字は基本的に地上（受信機等）に関連する項を示し、右上の添え字は基本的に上空（衛星等）に関連する項を示す。δｔｕ、Ｌ１ｏｒＬ２  ｂｉａｓは、受信機周波数依存ハードウェア依存バイアスを示し、δｔｓＬ１ｏｒＬ２  ｂｉａｓは、衛星周波数依存ハードウェア依存バイアスを示す。さらに、Ｉは、電離層伝搬遅延量を示し、Ｔは、対流圏伝搬遅延量を示す。ＮＬ１ｏｒＬ２，ａｍｂは、整数不確定値を示し、εは、観測誤差を示す。
【００４２】
衛星ＴＥＣ算出部３３は、２つの周波数の観測値の差をとることにより、ＴＥＣを求めることができる。
【数２】
2011185665.tif 000004
【００４３】
（３）式は、２つの周波数の擬似距離の差を示す。（４）式は、位相距離の差を示す。ここで、ｆは周波数である。また、λは波長を示す。Δｂｉａｓは、周波数間バイアスを示す。また、ＴＥＣｔｒｕｅは、真の総電子数を示す。
【００４４】
次に、ＴＥＣバイアス補正部３４は、周波数間バイアス推定部３５により事前に推定された周波数間バイアスΔｂｉａｓによる補正をＴＥＣに実施する（ステップＳ２）。まず、バイアス補正後の２つの周波数の擬似距離差及び位相距離差は、次式により表される。
【数３】
2011185665.tif 000005
【００４５】
ここで、（６）式には、観測値から分からない不確定値（λＬ１ΔＮＬ１，ａｍｂ−λＬ２ΔＮＬ２，ａｍｂ）が含まれている。そこで、ＴＥＣバイアス補正部３４は、（５）式と（６）式とを組み合わせて、以下に示す式により整数不確定値を削除する。この手法は、非特許文献３にも記載されている。
【数４】
2011185665.tif 000006
【００４６】
ここで、添え字ｋは、データの番号を示す。また、Ｍは、連続的に収集できたサンプルの合計を示す。さらに、ＴＥＣｋ，ｕは、衛星ｕの伝搬経路上の総電子数（個／ｍ２）を示す。
【００４７】
電子密度モデル算出部３７は、電離層電子密度モデル関数に基づき電子密度モデルを算出し、算出した電子密度モデルを３次元電子密度推定空間設定部３６に出力する。なお、電子密度モデル算出部３７による電離層モデルは、例えば非特許文献１，２，４等に記載されている。また、ステップＳ１，Ｓ２（特にＳ１）は、本発明の総電子数算出ステップに対応する。
【００４８】
次に、３次元電子密度推定空間設定部３６は、電子密度を推定する３次元空間を設定し、当該３次元空間を複数の領域に分割する。この３次元電子密度推定空間設定部３６の動作は、本発明の空間設定ステップに対応する。また、３次元電子密度推定空間設定部３６は、推定する３次元空間の電離層電子密度モデル値あるいは推定値を設定する。
【００４９】
具体的には、３次元電子密度推定空間設定部３６は、衛星から衛星受信機が配置されている地点まで３次元空間における電離層を、緯度・経度・高さ方向に分割し、電離層モデル値を使用して分割した領域内あるいは格子点の位置の電子密度を設定する。非特許文献１，２等に記載されたＩＲＩモデルは、高さ１０００ｋｍまでが有効範囲であるため、衛星が存在する空間までは設定できない。そこで、３次元電子密度推定空間設定部３６は、例えば非特許文献４に記載されたＧＣＰＭ（Ｇｌｏｂａｌ  Ｃｏｒｅ  Ｐｌａｓｍａ  Ｍｏｄｅｌ）を用いる手法や、ＩＲＩモデルの約４００ｋｍ〜１０００ｋｍを外挿する手法を採用することができる。
【００５０】
図４は、本実施例の電子密度推定装置３２内の３次元電子密度推定空間設定部３６による３次元空間設定を説明する図である。本実施例の３次元電子密度推定空間設定部３６は、図４に示すように、３次元空間をＸ１〜Ｘ１２の１２個の領域（ブロック）に分割しているが、実際には膨大な数のブロック数となる。図４において、測位衛星Ｓ１〜Ｓ４と衛星信号受信機Ｒｘ１〜Ｒｘ３との間における信号通過経路は、矢印により表されている。図４の衛星・ＴＥＣ・電子密度の関係は、次式により表される。
【数５】
2011185665.tif 000007
【００５１】
（９）式において、Ｘ１〜Ｘ１２は、空間を分割した領域（メッシュ）の電子密度を示し、現時点では不明である。また、ＴＥＣｉｊは、ｉ衛星・ｊ受信機間で決まる総電子数を示す。例えば、図４中の太い矢印ＴＥＣ２１は、測位衛星Ｓ２と衛星信号受信機Ｒｘ１との間の信号通過経路における総電子数を示す。ＴＥＣバイアス補正部３４により出力されたＴＥＣの値は、（９）式の右辺に代入される。
【００５２】
ここで、（９）式の行列式をＡＸ＝Ｅと表すことにする。行列Ａの中のゼロでない要素（１で表されている要素）は、ｉ衛星・ｊ受信機間の信号が通過したメッシュを示す。電子密度推定部３８は、ＡＸ＝Ｅの式における電子密度Ｘ（Ｘ１〜Ｘ１２）をトモグラフィ手法を用いて推定する（ステップＳ３）。
【００５３】
具体的には、電子密度推定部３８は、ＭＡＲＴあるいはＧＳＶＤ手法等により、衛星信号通過経路の電子密度を推定する。すなわち、電子密度推定部３８は、総電子数算出部（ＴＥＣバイアス補正部３４）により算出された総電子数に基づいて３次元電子密度推定空間設定部３６により分割された領域毎の電子密度を推定する。この電子密度推定部３８の動作は、本発明の電子密度推定ステップに対応する。
【００５４】
なお、本実施例において、電子密度推定部３８は、ＭＡＲＴあるいはＧＳＶＤ手法により電子密度を推定するが、本発明を実現するためには、必ずしもいずれかの手法を使う必要はなく、衛星信号受信機が得た観測値に基づいて暫定的に各領域の電子密度を推定できればよい。
【００５５】
ＭＡＲＴの手法は、上述したように非特許文献５，６等に記載されているが、ここで簡潔に説明する。
【数６】
2011185665.tif 000008
【００５６】
（１０）式は、ＭＡＲＴのアルゴリズムを示す式である。ａｉは、（９）式のＡ行列ｉ行ベクトルを表す。また、ａｉｊは、Ａ行列のｉｊ成分を表す。Ｘｊは、Ｘ列行列のｊ成分を表す。さらに、Ｅｉは、（９）式のＥ行列のｉ成分を表している。また、λは、収束の度合いを決めるパラメータである。
【００５７】
電子密度推定部３８は、（１０）式に示すＭＡＲＴのアルゴリズムを使用して、（１０）式の値が収束するまで繰り返し処理を実施し、３次元電子密度推定空間設定部３６により分割された領域毎の電子密度を推定する。その際に、電子密度推定部３８は、ＴＥＣバイアス補正部３４により出力されたＴＥＣの値をＥ行列の値に使用するとともに、電子密度モデル算出部３７により算出されたモデル値をＸの初期値として使用する。
【００５８】
なお、電子密度推定部３８は、信号が全く通過していない空間の電子密度を推定することはできない。
【００５９】
ＧＳＶＤ手法を採用する場合には、電子密度推定部３８は、以下に示す方式で電子密度を求める。なお、ＧＳＶＤ手法については、非特許文献７に記載がある。
【数７】
2011185665.tif 000009
【００６０】
（１１）式のＬは、Ｔｉｋｈｏｎｏｖ行列と呼ばれる。０次のＴｉｋｈｏｎｏｖ正規化では、Ｌを単位行列とする。αは、正規化パラメータで、α＞０とする。また、Ｌは、正値２次形式行列（固有値が全て正）である。
【００６１】
（１１）式を満足する電子密度Ｘは、次式で表すことができる。
【数８】
2011185665.tif 000010
【００６２】
なお、電子密度推定部３８は、ＡＸ＝Ｅの式における電子密度Ｘ（Ｘ１〜Ｘ１２）を推定する際に、本発明の発明者により発明され既に出願されている特願２００９−２９７３０９に記載された手法を用いて電子密度を推定してもよい。
【００６３】
一方、信号処理部１は、既知の送信源（放送局等）により送信された電波の受信状況を確認する（ステップＳ４）。なお、信号処理部１の代わりに電子密度推定装置３２内のホップ数算出部３９が行ってもよい。その場合には、ホップ数算出部３９は、演算部１５から受信した信号強度等の情報を受け取って当該処理を行う。本実施例において、この受信状況の確認処理を行うのは演算部１５であるものとする。また、ステップＳ４以降の処理を行うホップ数算出部３９、判定部４１、及び修正部４３は、その作用効果において本発明の特徴を示すものであり、従来手法と異なる点である。
【００６４】
アレイアンテナ１１は、既知の送信源により送信される電波を受信する。演算部１５は、Ａ／Ｄ変換部１４−１〜１４−ｎの各々により出力されたデジタル受信信号に基づいて、アレイアンテナ１１で受信した電波の到来方位と仰角とを算出する（本発明の第１算出ステップ）とともに、例えば、信号強度等の情報に基づいて次に示す（１３）式及び（１４）式を満たすか否かを判断することにより受信状況の確認を行う。
【数９】
2011185665.tif 000011
【００６５】
ここで、Ｐｒｃｖは、アレイアンテナ１１が受信した電波の信号強度である。受信した信号強度Ｐｒｃｖが予め設定したしきい値εｐ以上であるか否か判断される。
【００６６】
また、Ａｚｒｃｖは、受信した電波の到来方位である。さらに、Ａｚ０は、既知の送信源の位置とアレイアンテナ１１の設置位置とで決まる幾何学的方位であり、次式に示す球面三角法により求めたものである。
【数１０】
2011185665.tif 000012
【００６７】
ここで、θｓｎｄは、送信源余緯度である。また、θｒｃｖは、受信源余緯度である。λｓｎｄは、送信源経度である。λｒｃｖは、受信源経度である。
【００６８】
演算部１５は、（１４）式に示すように受信した信号の方位Ａｚｒｃｖが、幾何学的方位（Ａｚ０）±εＡＺの範囲内にあるか否かを判断する。なお、εＡＺは、予め演算部１５に設定された方位に対するしきい値である。
【００６９】
演算部１５は、受信した電波に基づく情報が（１３）式及び（１４）式を満たす場合には、既知の送信源から送信された電波を受信できたと判断する。ただし、既知の送信源の送信周波数が他の送信源で使用されていない場合や他の送信源の信号を受信していないことが明らかである場合等には、必ずしも方位の条件は必須ではなく、演算部１５は（１３）式のみを用いて受信状況の確認を行うことができる。
【００７０】
受信した電波に基づく情報が（１３）式及び（１４）式を満たさない場合には、演算部１５は、既知の送信源から送信された電波を受信できていないと判断し、既知送信源位置による処理は行わない（ステップＳ５）。
【００７１】
受信した電波に基づく情報が（１３）式及び（１４）式を満たす場合には、演算部１５は、既知の送信源から送信された電波を受信できていると判断し、その旨をホップ数算出部３９に伝達する。ホップ数算出部３９は、電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系３０により算出された電波の到来方位及び仰角と電子密度推定部３８により推定された電子密度とに基づいて、アレイアンテナ１１が受信した電波のホップ数を算出する（ステップＳ６）。このホップ数算出部３９の動作は、本発明の第２算出ステップに対応する。
【００７２】
図５は、本実施例の電子密度推定装置３２内のホップ数算出部３９によるホップ数の算出方法を説明する図である。図５に示すように、既知の送信源から送信された電波は、所定の電子密度を有する電離層と地上との間で反射し、ホップするようにして受信位置（アレイアンテナ１１）に到達する。
【００７３】
ホップ数算出部３９によるホップ数の算出方法には様々なものが考えられるので、いくつかの方法について例示する。まず、簡易伝搬距離推定について説明する。受信位置（すなわちアレイアンテナ１１の設置位置）で観測した信号の受信仰角をθとする。受信位置上空の電子密度で仰角θにより反射する際の電子密度は、非特許文献１１に記載の方式から以下に示す式で推定される。
【数１１】
2011185665.tif 000013
【００７４】
ここで、ｈは高さを示す。ε０は真空の誘電率を示す。ｍｅは電子の質量を示す。ｆは受信（送信）周波数を示す。ｑは電子の電荷を示す。ホップ数算出部３９は、電子密度推定部３８により推定された電子密度に基づいて、（１７）式で求められた密度ρに相当する電子密度を有する高さｈを特定する。なお、ホップ数算出部３９は、高さｈを特定する際に、非特許文献１〜３に関連した電子密度モデルを参照してもよいし、衛星を使用したトモグラフィ手法による電子密度モデルを参照してもよい。
【００７５】
次に、ホップ数算出部３９は、特定した高さｈをｈｒとして使用し、概略の受信位置からの１ホップ推定距離ｄ（図５参照）を求める。反射点高さｈｒと概略距離ｄとは、次に示す（１８）式あるいは（１９）式のような関係を有する。
【数１２】
2011185665.tif 000014
【００７６】
ここで、Ｒは受信位置における地球半径を示す。（１８）式は、電離層と地表が平面である場合の簡易的なモデルである。一方、（１９）式は、電離層が平面と考え、地表の湾曲を考慮したモデルである。
【００７７】
ホップ数算出部３９は、複数ホップを想定し、（１８）式あるいは（１９）式を用いて概略距離ｄを求め、ホップ数倍する。
【００７８】
Ｄ＝ｎｈｏｐ×ｄ      …（２０）
ここで、ｎｈｏｐはホップ数である。この時点でホップ数ｎｈｏｐは不明であるため、ホップ数算出部３９は、ｎｈｏｐに１，２，３，…の数値を仮に代入し、それぞれについてのＤの値を算出する。
【００７９】
次に、ホップ数算出部３９は、計算した距離Ｄと、球面三角法による（２１）式及び（２２）式を連立して解くことにより、距離Ｄの経度（λ１）及び緯度（φ１）を計算する。
【数１３】
2011185665.tif 000015
【００８０】
次に、ホップ数算出部３９は、以下に示す球面三角法の式を用いて、距離Ｄに相当する位置と既知の送信位置との間の距離を求める。
【数１４】
2011185665.tif 000016
【００８１】
ホップ数算出部３９は、ホップ数ｎｈｏｐの数値を例えば１，２，３，…と変えた場合のそれぞれについて、（２３）式に示す距離Ｘｎｈｏｐを算出し、算出したＸｎｈｏｐの中で最小の距離となるものを選択し、そのＸｎｈｏｐに相当するホップ数ｎｈｏｐを実際に伝搬したホップ数とする。なお、ホップ数算出部３９は、算出したホップ数の信頼性を高めるために、複数の観測結果を使用してホップ数毎の平均を求め、平均値を使用してＸｎｈｏｐの最小値を求めることもできる。この場合には、ホップ数算出部３９は、次式を用いて計算する。
【数１５】
2011185665.tif 000017
【００８２】
ホップ数算出部３９は、ホップ数ｎｈｏｐ毎に（２４）式の平均値を求める。ここで、（２４）式の左辺の添え字ｍは、平均を示す。ホップ数算出部３９は、この平均値の中で最小となる距離に相当するホップ数ｎｈｏｐを実際に伝搬してきたホップ数とする。
【００８３】
以上述べた簡易伝搬距離推定の手法の代わりに、ホップ数算出部３９は、レイトレーシングの手法を用いてホップ数を算出することもできる。この場合に、ホップ数算出部３９は、非特許文献１２，１３のレイトレーシング方式と電子密度モデル（非特許文献１〜３に関連した電子密度モデル、あるいは衛星を使用したトモグラフィ手法により求めたもの）を使用し、受信位置から複数ホップ数分レイトレーシングを実施する。その後、ホップ数算出部３９は、レイトレーシングの結果推定された推定送信源位置と（２１）〜（２４）式とを使用し、最小となる距離に相当するホップ数を実際に伝搬してきたホップ数とする。
【００８４】
また、ホップ数算出部３９は、簡易伝搬距離推定及びレイトレーシングの手法の代わりに、受信電力による推定でホップ数を算出することもできる。この場合に、ホップ数算出部３９は、非特許文献１４による信号強度の手法により推定する。
【数１６】
2011185665.tif 000018
【００８５】
ここで、Ｅは、非特許文献１４の信号強度式であり、ホップ数に関係しない項である。ｆは送信（受信）周波数を示す。ｐ´（ｎｈｏｐ）は仮想斜距離（ｖｉｒｔｕａｌ  ｓｌａｎｔ  ｒａｎｇｅ）を示す。Ｌｉ（ｎｈｏｐ）は吸収による減衰を示す。Ｐｔは既知送信局の送信電力を示す。Ｇｔは送信アンテナゲインを示す。ＬｍはＭＵＦ（Ｍａｘｉｍｕｍ  Ｕｓｅａｂｌｅ  Ｆｒｅｑｕｅｎｃｙ）を超える場合の減衰を示す。Ｌｈはオーロラあるいは他の減衰項を示す。Ｌｚはその他の減衰項を示す。
【００８６】
ホップ数算出部３９は、予め（２６）式の値を計算し、その後、ホップ数を変えてＥｒ式の値を（２５）式から計算する。このとき、Ｅｒは、サンプルの中央値を表している。ホップ数算出部３９は、何回かの受信強度の観測結果を演算部１５から受け取り、その中央値を代表値として（２５）式で計算したＥｒと比較する。その結果、ホップ数算出部３９は、最も近い信号強度に相当するホップ数を実際に伝搬してきたホップ数とする。
【００８７】
次に判定部４１（あるいはホップ数算出部３９又は修正部４３）は、レイトレーシングを実施し、既知の送信源位置からの電波を仮想的に追跡する（ステップＳ７）。図５の実線部分は、実際の電波の伝搬経路を示す。一方、１点鎖線は、計算した電子密度モデル（非特許文献１〜３に関連した電子密度モデル、あるいは衛星を使用したトモグラフィ手法により求められたもの）を使用したレイトレーシングによる伝搬経路を示す。
【００８８】
また、判定部４１は、ホップ数算出部３９により算出されたホップ数に基づいて算出したアレイアンテナ１１から送信源（電波を送信した放送局）までの距離と実際の距離との差の大きさが所定のしきい値ε以下であるか否かを判断し、所定のしきい値ε以下である場合に、電子密度推定部３８による推定結果を正しいと判定する（ステップＳ８）。この判定部４１の動作は、本発明の判定ステップに対応する。
【００８９】
すなわち、判定部４１は、ホップ数を特定した場合に（２３）式で算出したＸｎｈｏｐの値がＸｎｈｏｐ≦εを満たすか否かを判断し、満たす場合には電子密度推定部３８による推定結果を正しい（真値に近い）と判定して終了する。
【００９０】
Ｘｎｈｏｐ≦εを満たさない場合には、判定部４１は、修正部４３に通過経路の電子密度を修正させる（ステップＳ９）。すなわち、修正部４３は、判定部４１により電子密度推定部３８の推定結果が正しいと判断されない場合（Ｘｎｈｏｐ＞εの場合）に、アレイアンテナ１１が受信した電波の伝搬経路における電子密度を修正する。この修正部４３の動作は、本発明の修正ステップに対応する。
【００９１】
修正部４３は、アレイアンテナ１１が受信した電波の伝搬経路における電子密度を修正する際に、伝搬経路としてレイトレーシング手法により特定されたものを用いる。修正部４３は、レイトレーシングを自ら実施してもよいし、ホップ数算出部３９あるいは判定部４１からレイトレーシング結果を受け取ってもよい。
【００９２】
修正部４３は、レイトレーシングで計算した位置と受信した位置との距離をＤｒａｙとし、既知の受信位置と送信局位置との距離をＤ０として差分を計算する。
【００９３】
ΔＤ＝Ｄｒａｙ−Ｄ０      …（２７）
修正部４３は、計算した差分ΔＤに基づいて、レイトレーシングにより計算した伝搬経路の電子密度を修正する。ΔＤ＞０の場合には、レイトレーシングによる伝搬経路は、実際の伝搬経路よりも高い場所で反射した経路になっていると考えられる。逆に、ΔＤ＜０の場合には、レイトレーシングによる伝搬経路は、実際の伝搬経路よりも低い場所で反射した経路となっていると考えられる。
【００９４】
電子密度を修正する方法は、一定値割合で変化させる方法と、ΔＤの大きさに応じて変化させる方法とが考えられる。
【数１７】
2011185665.tif 000019
【００９５】
（２８）式は、一定割合で変化させる方式である。また、（２９）式は、ΔＤの大きさに依存して変化させる方式である。ここで、φは緯度を示す。λは経度を示す。ｈは高さを示す。また、α，βは、変化率に対する係数を表しており、修正部に予め設定されている値である。
【００９６】
修正部４３は、伝搬経路上の電子密度を（２８）式あるいは（２９）式を用いて修正し、修正した結果を電子密度推定部３８に出力する。電子密度推定部３８は、修正部４３により修正された電子密度に基づいて再度領域毎の電子密度を推定する。これにより、電子密度推定装置３２は、ステップＳ８においてＸｎｈｏｐ≦εを満たすまで、ステップＳ３〜Ｓ９の工程を繰り返す。
【００９７】
上述のとおり、本発明の実施例１の形態に係る電子密度推定装置及び電子密度推定方法によれば、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して３次元の電離層電子密度分布をより正確に推定することができる。
【００９８】
すなわち、本実施例の電子密度推定装置及び電子密度推定方法は、既知の送信源（放送局等）から放送される信号の伝搬経路を推定し、経路上の電子密度を求めるので、従来の衛星観測による３次元電子密度分布のデータを補強し修正することができる。特に、本発明は、既知の送信源として日本国内の放送局のみならず、海外の放送局等も利用することができるので、より容易且つ正確に３次元電子密度分布を推定することができる。
【００９９】
さらに、電子密度推定部３８は、修正部４３による修正結果を反映させて、再度電子密度の推定を行うので、より正確な３次元電子密度分布を求めることができる。
【０１００】
また、修正部４３は、電離層電子密度モデルを使用したレイトレーシング手法により特定された伝搬経路上の電子密度を修正するので、既知の送信源による電波を電子密度推定に十分活用することができ、衛星観測結果のみを用いた場合の観測値の少なさという問題点を解決することができる。
【０１０１】
さらに、本実施例の電子密度推定装置及び電子密度推定方法は、既知の送信源の位置と計算による送信源の位置との差が所定のしきい値ε以下になったか否かを判断基準とし、ステップＳ３からＳ９までの処理を繰り返すことにより３次元電子密度推定の近似精度を向上させることができる。
【産業上の利用可能性】
【０１０２】
本発明に係る電子密度推定装置及び電子密度推定方法は、衛星から送信される複数の異なる周波数信号による通過経路の電離層電子密度分布を推定する電子密度推定装置に利用可能である。
【符号の説明】
【０１０３】
１      信号処理部
１０ａ、１０ｂ、１０ｃ  測位衛星
１１    アレイアンテナ
１２    アンプ
１３    周波数変換部
１４    Ａ／Ｄ変換部
１５    演算部
２０    衛星信号処理系
２２    アンテナ
２４    衛星信号受信機
２６    衛星信号処理装置
３０    電離層ＴＥＣ・バイアス推定・３次元電子密度推定処理系
３２    電子密度推定装置
３３    衛星ＴＥＣ算出部
３４    ＴＥＣバイアス補正部
３５    周波数間バイアス推定部
３６    ３次元電子密度推定空間設定部
３７    電子密度モデル算出部
３８    電子密度推定部
３９    ホップ数算出部
４０    インターネットデータ処理系
４１    判定部
４２    ルータ
４３    修正部
４４    ＧＥＯＮＥＴ収集データ処理装置
４６    外部インターネット網
５０    通信回線
Ｓ１，Ｓ２，Ｓ３，Ｓ４  測位衛星
Ｒｘ１，Ｒｘ２，Ｒｘ３  衛星信号受信機
【図１】
2011185665.tif 000020
【図２】
2011185665.tif 000021
【図３】
2011185665.tif 000022
【図４】
2011185665.tif 000023
【図５】
2011185665.tif 000024