説明

電子密度推定装置及び電子密度推定方法

【課題】衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して3次元の電離層電子密度分布をより正確に推定する。
【解決手段】アレイアンテナが受信した電波の到来方位と仰角とを算出する第1算出部と、衛星信号の通過経路の総電子数を算出する総電子数算出部と、電子密度を推定する3次元空間を設定し当該3次元空間を複数の領域に分割する3次元電子密度推定空間設定部36と、領域毎の電子密度を推定する電子密度推定部38と、アレイアンテナが受信した電波のホップ数を算出するホップ数算出部39と、ホップ数に基づいて算出したアレイアンテナから送信源までの距離と実際の距離との差が所定のしきい値以下である場合に電子密度推定部による推定結果を正しいと判定する判定部41と、判定部41の判定結果に応じて電波の伝搬経路における電子密度を修正する修正部43とを備える。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、衛星から送信される複数の異なる周波数信号による通過経路の電離層電子密度分布を推定する電子密度推定装置及び電子密度推定方法に関する。
【背景技術】
【0002】
衛星ナビゲーションシステム「Galileo」やGPS(Global Positioning System)、QZSS(Quasi Zenith Satellites System;準天頂衛星システム)等の衛星から送信される信号は、地球上の観測点の測位を行うのに利用される。これらの衛星により放送される複数周波数信号の擬似距離差は、電波の伝搬経路上の電子数に関係していることが知られている。
【0003】
衛星から地上に向けて送信される信号は、電離層を通過する際に伝搬遅延量を生じる。この電離層伝搬遅延量は、伝搬する信号の周波数に依存している。そのため、衛星からの複数の周波数信号を受信することにより、信号の通過経路における総電子数(TEC:Total Electron Content)を得ることができる。電波の伝搬経路上の総電子数を推定する場合には、複数周波数信号の送信・受信時刻による擬似距離(コード距離)差を使用する方法や複数周波数のキャリア位相擬似距離(フェーズ距離)の差を利用する方法が知られている。
【0004】
非特許文献1や非特許文献2には、電離層における電子密度モデル関数であるIRI(International Reference Ionosphere)モデルが記載されている。さらに、非特許文献4には、やはり電離層電子密度モデル関数の一つであるGallagherのモデルが記載されている。これらのモデル関数を用いて計算によりTECを求めることもできる。
【0005】
非特許文献5や非特許文献6には、トモグラフィ手法による衛星TECを用いた3次元電子密度推定が記載されており、MART(Multiplicative Algebraic Reconstruction Technique)と呼ばれる代数的再構成法が適用されている。一方、非特許文献7には、GSVD(Generalized Singular Value Decomposition)を用いたトモグラフィック再構成による電子密度推定が記載されている。
【0006】
また、非特許文献8及び非特許文献9には、電離層を球核の薄い層とみなした電子密度計算について記載されており、GPSを使用した観測結果のみに基づいて3次元電子密度推定を試みている。さらに、非特許文献10には、3次元方向の電子密度分布を高さ方向に計算した場合について記載されている。
【0007】
HF帯を用いた電波通信において、受信した電波の到来方法及び仰角を測定し、この到来方法及び仰角と、IRIモデル等の電離層電子密度分布モデル等とを用いて、レイトレーシング手法等により電波の伝搬経路を計算することができる。このような手法については、非特許文献11,12,13に記載されている。また、非特許文献14には、電波の伝搬特性を予測する手法が記載されている。
【先行技術文献】
【非特許文献】
【0008】
【非特許文献1】Dieter Bilitza:“International Reference Ionosphere 2000”,Radio.Science,Vol.36,Number2,PP261−275,March/April,2001
【非特許文献2】Dieter Bilitza,et,al.,:“International Reference Ionosphere 1990”,November,1990.
【非特許文献3】A.Komjathy:‘Global Ionospheric Total Electron Content Mapping Using the Global Positioning System’ ,UNB,Technical Report No.188,Sep.1997.
【非特許文献4】Gallagher,D.L.,P.D.Craven,and R.H.Comfort,Global core plasma model,J.Geophys.Res.105,A8,18,819−18,833,2000.
【非特許文献5】Raymund,T.D.,Austin,J.R.,Franke,S.J.,Liu,C.H.,Klobuchar,J.A.,and Stalker,J.,:“Application of Comuterized Tomography to the Investigation of Ionospheric Structures”,Radio Sci.,25(5),771−789,1990.
【非特許文献6】Stefan Schluer,et al.:“Monitoring the 3 Dimensional Ionospheric Electron Distribution based on GPS Measurements”
【非特許文献7】K.Bhuyan,et al:“Tomographicreconstruction of the Ionosphere using generalized singular value decomposition” ,Current Sci.,Vol1.83,No.9,10 Nov.,2002.
【非特許文献8】Otuka,Y.,et.al.:“A new technique for mapping of total electron content using GPS network in Japan,”,Earth Planets Space,111−120,2001.
【非特許文献9】Ma,G.,and T.Maruyama(2003),Derivation of TEC and estimation of instrumental biases from GEONET in Japan,Ann.Geophys.,21,2083−2093.
【非特許文献10】五十嵐来良,斉藤昭則,大塚雄一:“標準電離層モデル(IRI)で導出したTECとGPS受信機網で観測したTECの比較解析”,第1回電離圏の利用と影響に関するシンポジウム講和集,pp24−1〜8,2003.
【非特許文献11】前田憲一,後藤三男:“電波伝搬”,岩波全書,1953年2月
【非特許文献12】K.G.Buden:“The Propagation of radio waves”.Cambridge University Press,1988.
【非特許文献13】Iwane Kimura:“Effects of Ions on Whistler−Mode Ray Tracing”,Radio Science.Vol.I,No.3,269−283,March 1966.
【非特許文献14】“HF propagation prediction method”,ITU−R P.533−7.
【発明の概要】
【発明が解決しようとする課題】
【0009】
しかしながら、地上で観測するGPS衛星電波の通過経路が限定されるので、非特許文献8及び非特許文献9に記載された手法のようにGPSを使用した観測結果のみに基づいて3次元電子密度推定を行う場合には、十分な3次元空間の情報を得ることができないという問題点がある。
【0010】
すなわち、従来の衛星観測結果を用いた技術は、GPS衛星から受信機までの通過経路上の総電子数しか求められないので、観測できる3次元空間が限定されてしまい、3次元の電離層電子密度分布を十分推定しきれていない。その結果、電波が通過していない空間の3次元電子密度推定値は、誤差が多く含まれるものとなる可能性が大きい。このような問題が生じる根本的な要因は、観測値が少ない点にあると考えられる。
【0011】
本発明は上述した従来技術の問題点を解決するもので、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して3次元の電離層電子密度分布をより正確に推定する電子密度推定装置及び電子密度推定方法を提供することを課題とする。
【課題を解決するための手段】
【0012】
本発明に係る電子密度推定装置は、上記課題を解決するために、複数の測位衛星から送信される衛星信号を受信する1以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定装置であって、既知の送信源により送信される電波を受信するアレイアンテナと、前記アレイアンテナが受信した電波の到来方位と仰角とを算出する第1算出部と、前記1以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出部と、電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する空間設定部と、前記総電子数算出部により算出された総電子数に基づいて前記空間設定部により分割された領域毎の電子密度を推定する電子密度推定部と、前記第1算出部により算出された電波の到来方位及び仰角と前記電子密度推定部により推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第2算出部と、前記第2算出部により算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定部による推定結果を正しいと判定する判定部と、前記判定部により前記電子密度推定部の推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正部とを備えることを特徴とする。
【0013】
本発明に係る電子密度推定方法は、上記課題を解決するために、複数の測位衛星から送信される衛星信号を受信する1以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定方法であって、アレイアンテナで既知の送信源により送信される電波を受信し、受信した電波の到来方位と仰角とを算出する第1算出ステップと、前記1以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出ステップと、電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する空間設定ステップと、前記総電子数算出ステップにより算出された総電子数に基づいて前記空間設定ステップにより分割された領域毎の電子密度を推定する電子密度推定ステップと、前記第1算出ステップにより算出された電波の到来方位及び仰角と前記電子密度推定ステップにより推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第2算出ステップと、前記第2算出ステップにより算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定ステップによる推定結果を正しいと判定する判定ステップと、前記判定ステップにより前記電子密度推定ステップの推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正ステップとを備えることを特徴とする。
【発明の効果】
【0014】
本発明によれば、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して3次元の電離層電子密度分布をより正確に推定する電子密度推定装置及び電子密度推定方法を提供することができる。
【図面の簡単な説明】
【0015】
【図1】本発明の実施例1の形態の電子密度推定装置及び測位システムの構成を示すブロック図である。
【図2】本発明の実施例1の形態の電子密度推定装置の詳細なブロック図である。
【図3】本発明の実施例1の形態の電子密度推定装置の動作を示すフローチャート図である。
【図4】本発明の実施例1の形態の電子密度推定装置内の3次元電子密度推定空間設定部による3次元空間設定を説明する図である。
【図5】本発明の実施例1の形態の電子密度推定装置内のホップ数算出部によるホップ数の算出方法を説明する図である。
【発明を実施するための形態】
【0016】
以下、本発明の電子密度推定装置及び電子密度推定方法の実施の形態について図面を参照して詳細に説明する。
【実施例1】
【0017】
以下、本発明の実施例について図面を参照しながら説明する。図1は、本発明の実施例1の電子密度推定装置32及び測位システムの構成を示すブロック図である。この電子密度推定装置32は、複数の測位衛星(測位衛星10a、10b、10c)から送信される衛星信号を受信する1以上の衛星信号受信機24を有し且つ衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用して電子密度を推定する。
【0018】
まず、本実施の形態の構成を説明する。本実施例の電子密度推定装置32及び測位システムは、図1に示すように、信号処理部1、測位衛星10a、測位衛星10b、測位衛星10c、衛星信号処理系20、電離層TEC・バイアス推定・3次元電子密度推定処理系30、及びインターネットデータ処理系40で構成されている。ただし、本実施例において信号処理部1は、電子密度推定装置の構成の一部であるものとする。
【0019】
測位衛星10a、測位衛星10b、測位衛星10c、及び衛星信号処理系20は、本発明の測位システムに対応する。なお、インターネットデータ処理系40も含めて測位システムとすることも可能であるが、インターネットデータ処理系40は、必ずしも必須のものではなく、付加的なものである。また、信号処理部1、衛星信号処理系20、電離層TEC・バイアス推定・3次元電子密度推定処理系30、及びインターネットデータ処理系40は、お互いに通信回線50で接続されている。
【0020】
測位衛星10a、測位衛星10b、測位衛星10cは、GPS、Galileo、準天頂衛星(QZSS)等の航法衛星であり、複数の周波数の衛星信号を送信する。
【0021】
衛星信号処理系20は、アンテナ22、衛星信号受信機24、及び衛星信号処理装置26から構成される。衛星信号受信機24は、本発明の受信機に対応し、複数の測位衛星(測位衛星10a、10b、10c)から送信される複数の周波数の衛星信号をアンテナ22を介して受信する。また、衛星信号処理装置26は、衛星信号に含まれる測位情報を用いて位置情報を得る。
【0022】
信号処理部1は、複数のアンテナ11−1,11−2,…,11−nからなるアレイアンテナ11と、アンプ12−1,12−2,…,12−nと、周波数変換部13−1,13−2,…,13−nと、A/D変換部14−1,14−2,…,14−nと、演算部15とを備える。
【0023】
アレイアンテナ11は、既知の送信源により送信される電波を受信する。例えば、アレイアンテナ11は、既存の放送局から送信され、電離層を通過して到来する放送信号(電波)を受信して受信信号を出力する。アンプ12−1〜12−nは、それぞれアンテナ11−1〜11−nからの信号を増幅する周波数変換部13−1〜13−nは、アンプ12−1〜12−nにより増幅された信号をベースバンドの信号に変換する。A/D変換部14−1〜14−nは、周波数変換部13−1〜13−nにより変換され出力された信号をA/D変換し、デジタル受信信号として出力する。
【0024】
演算部15は、本発明の第1算出部に対応し、A/D変換部14−1〜14−nの各々により出力されたデジタル受信信号に基づいて、アレイアンテナ11が受信した電波の到来方位と仰角とを算出する。例えば、演算部15は、MUSIC(Multiple Signal Classification)法や独立成分分析(ICA:Independent Component Analysis)の手法などを用いて、既存の放送局からの電波の到来方位角および到来仰角を算出する。
【0025】
電離層TEC・バイアス推定・3次元電子密度推定処理系30は、電離層電子密度分布モデルと演算部15で算出された方位角及び仰角とを用いて、既存の放送局からの電波の伝搬経路をレイトレーシング手法により算出し、算出した伝搬経路に基づいて、電波の到来方位角および到来仰角を算出する。なお、レイトレーシング手法の詳細については、例えば非特許文献12に記載されている。
【0026】
あるいは、電離層TEC・バイアス推定・3次元電子密度推定処理系30は、電子密度推定装置32を有し、受信した信号に基づき衛星位置、信号通過経路の電離層総電子数(TEC)、周波数間バイアス等を算出する。受信したデータ・処理結果は、データサーバ系へLAN経由で伝送され保存される。
【0027】
インターネットデータ処理系40は、ルータ42、衛星データ処理装置44、及び外部インターネット網46で構成されている。ルータ42は、スイッチングハブでもよい。インターネットデータ処理系40は、公開されている電離層関連の情報や国土地理院が公開しているGPS観測データ(GEONETデータ)や国際的に観測結果を公開しているIGS(International GPS Service for Geodynamics)データ等をインターネット経由で収集する装置である。インターネットデータ処理系40により収集された結果は、電離層TEC・バイアス推定・3次元電子密度推定処理系30で処理される。ルータ42は、セキュリティを考慮して設けられ、ファイアウォールとする。
【0028】
図2は、本発明の実施例1の電子密度推定装置32の詳細なブロック図である。図2に示すように、電子密度推定装置32は、衛星TEC算出部33、TECバイアス補正部34、周波数間バイアス推定部35、3次元電子密度推定空間設定部36、電子密度モデル算出部37、電子密度推定部38、ホップ数算出部39、判定部41、及び修正部43で構成されている。
【0029】
衛星TEC算出部33、TECバイアス補正部34、及び周波数間バイアス推定部35は、本発明の総電子数算出部に対応し、1以上の受信機の各々から複数の測位衛星の各々に対する擬似距離及びキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する。具体的には、衛星TEC算出部33は、衛星信号受信機24から測位衛星10a、10b、10cの各々に対する擬似距離及びキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する。
【0030】
周波数間バイアス推定部35は、周波数間バイアスを算出し、TECバイアス補正部34に出力する。周波数間バイアス推定部35による周波数間バイアスの推定方法は、どのようなものでもよく、従来手法を用いてよい。
【0031】
TECバイアス補正部34は、周波数間バイアス推定部35により推定された周波数間バイアスを用いて補正を行い、最終的な総電子数を算出する。
【0032】
3次元電子密度推定空間設定部36は、本発明の空間設定部に対応し、電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する。
【0033】
電子密度モデル算出部37は、上述したようなIRIモデルやGallagherモデル等の電離層電子密度モデル関数に基づき電子密度モデルを算出する。
【0034】
電子密度推定部38は、総電子数算出部(TECバイアス補正部34)により算出された総電子数に基づいて3次元電子密度推定空間設定部36により分割された領域毎の電子密度を推定する。本実施例において、電子密度推定部38は、非特許文献5や非特許文献に記載されているようなMARTやGSVD等のトモグラフィの手法を用いて電子密度を推定するものとする。
【0035】
ホップ数算出部39は、本発明の第2算出部に対応し、通信回線50を介して信号処理部1内の演算部15から情報を得ることにより、演算部15により算出された電波の到来方位及び仰角と電子密度推定部38により推定された電子密度とに基づいて、アレイアンテナ11が受信した電波のホップ数を算出する。
【0036】
判定部41は、ホップ数算出部39により算出されたホップ数に基づいて算出したアレイアンテナ11から送信源(電波を送信した放送局)までの距離と実際の距離との差が所定のしきい値以下である場合に、電子密度推定部38による推定結果を正しいと判定する。
【0037】
修正部43は、判定部41により電子密度推定部38の推定結果が正しいと判断されない場合に、アレイアンテナ11が受信した電波の伝搬経路における電子密度を修正する。すなわち、ホップ数算出部39により算出されたホップ数に基づいて判定部41が算出したアレイアンテナ11から送信源までの距離と実際の距離との差の大きさが所定のしきい値よりも大きい場合に、修正部43は、アレイアンテナ11が受信した電波の伝搬経路における電子密度を修正する。なお、修正部43は、アレイアンテナ11が受信した電波の伝搬経路における電子密度を修正する際に、伝搬経路としてレイトレーシング手法により特定された伝搬経路を用いる。
【0038】
また、電子密度推定部38は、修正部43により修正された電子密度に基づいて再度領域毎の電子密度を推定する。
【0039】
次に、上述のように構成された本実施の形態の作用を説明する。図3は、本実施例の電子密度推定装置32の動作を示すフローチャート図である。まず、衛星TEC算出部33は、1以上の受信機の各々から複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき衛星信号の通過経路の総電子数を算出する(ステップS1)。具体的には、衛星TEC算出部33は、広範囲に設置した衛星受信機収集データからGPS等の衛星のTECを計算する。ここで、本実施例の衛星TEC算出部33は、衛星から放送される2周波の信号(GPS衛星の場合L1周波数(1575.42MHz)とL2周波数(1227.60MHz)の信号)に基づいて信号伝搬経路のTECを算出するものとする。GPS以外でも異なる2周波を放送している場合には同様な手法が使用できる。
【0040】
ここでは、測位衛星10aと衛星信号受信機24との間の擬似距離について考える。衛星信号通過経路のTECは、測位衛星10aから送信される衛星信号(ここではL1とL2で示す)に基づき求められる。まず、擬似距離(コード距離、シュードレンジ)とキャリア位相擬似距離(フェーズ距離)は、以下のように表すことができる。
【数1】

【0041】
ここで、ρは、遅延時間による距離(擬似距離)を示す。また、Φは、キャリア位相を示す。rは、真の距離を示す。cは光速である。また、δtは、受信機時刻誤差を示し、δtは、衛星時刻誤差を示す。なお、本実施例において、右下の添え字は基本的に地上(受信機等)に関連する項を示し、右上の添え字は基本的に上空(衛星等)に関連する項を示す。δtu、L1orL2 biasは、受信機周波数依存ハードウェア依存バイアスを示し、δtL1orL2 biasは、衛星周波数依存ハードウェア依存バイアスを示す。さらに、Iは、電離層伝搬遅延量を示し、Tは、対流圏伝搬遅延量を示す。NL1orL2,ambは、整数不確定値を示し、εは、観測誤差を示す。
【0042】
衛星TEC算出部33は、2つの周波数の観測値の差をとることにより、TECを求めることができる。
【数2】

【0043】
(3)式は、2つの周波数の擬似距離の差を示す。(4)式は、位相距離の差を示す。ここで、fは周波数である。また、λは波長を示す。Δbiasは、周波数間バイアスを示す。また、TECtrueは、真の総電子数を示す。
【0044】
次に、TECバイアス補正部34は、周波数間バイアス推定部35により事前に推定された周波数間バイアスΔbiasによる補正をTECに実施する(ステップS2)。まず、バイアス補正後の2つの周波数の擬似距離差及び位相距離差は、次式により表される。
【数3】

【0045】
ここで、(6)式には、観測値から分からない不確定値(λL1ΔNL1,amb−λL2ΔNL2,amb)が含まれている。そこで、TECバイアス補正部34は、(5)式と(6)式とを組み合わせて、以下に示す式により整数不確定値を削除する。この手法は、非特許文献3にも記載されている。
【数4】

【0046】
ここで、添え字kは、データの番号を示す。また、Mは、連続的に収集できたサンプルの合計を示す。さらに、TECk,uは、衛星uの伝搬経路上の総電子数(個/m)を示す。
【0047】
電子密度モデル算出部37は、電離層電子密度モデル関数に基づき電子密度モデルを算出し、算出した電子密度モデルを3次元電子密度推定空間設定部36に出力する。なお、電子密度モデル算出部37による電離層モデルは、例えば非特許文献1,2,4等に記載されている。また、ステップS1,S2(特にS1)は、本発明の総電子数算出ステップに対応する。
【0048】
次に、3次元電子密度推定空間設定部36は、電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する。この3次元電子密度推定空間設定部36の動作は、本発明の空間設定ステップに対応する。また、3次元電子密度推定空間設定部36は、推定する3次元空間の電離層電子密度モデル値あるいは推定値を設定する。
【0049】
具体的には、3次元電子密度推定空間設定部36は、衛星から衛星受信機が配置されている地点まで3次元空間における電離層を、緯度・経度・高さ方向に分割し、電離層モデル値を使用して分割した領域内あるいは格子点の位置の電子密度を設定する。非特許文献1,2等に記載されたIRIモデルは、高さ1000kmまでが有効範囲であるため、衛星が存在する空間までは設定できない。そこで、3次元電子密度推定空間設定部36は、例えば非特許文献4に記載されたGCPM(Global Core Plasma Model)を用いる手法や、IRIモデルの約400km〜1000kmを外挿する手法を採用することができる。
【0050】
図4は、本実施例の電子密度推定装置32内の3次元電子密度推定空間設定部36による3次元空間設定を説明する図である。本実施例の3次元電子密度推定空間設定部36は、図4に示すように、3次元空間をX1〜X12の12個の領域(ブロック)に分割しているが、実際には膨大な数のブロック数となる。図4において、測位衛星S1〜S4と衛星信号受信機Rx1〜Rx3との間における信号通過経路は、矢印により表されている。図4の衛星・TEC・電子密度の関係は、次式により表される。
【数5】

【0051】
(9)式において、X1〜X12は、空間を分割した領域(メッシュ)の電子密度を示し、現時点では不明である。また、TECは、i衛星・j受信機間で決まる総電子数を示す。例えば、図4中の太い矢印TECは、測位衛星S2と衛星信号受信機Rx1との間の信号通過経路における総電子数を示す。TECバイアス補正部34により出力されたTECの値は、(9)式の右辺に代入される。
【0052】
ここで、(9)式の行列式をAX=Eと表すことにする。行列Aの中のゼロでない要素(1で表されている要素)は、i衛星・j受信機間の信号が通過したメッシュを示す。電子密度推定部38は、AX=Eの式における電子密度X(X1〜X12)をトモグラフィ手法を用いて推定する(ステップS3)。
【0053】
具体的には、電子密度推定部38は、MARTあるいはGSVD手法等により、衛星信号通過経路の電子密度を推定する。すなわち、電子密度推定部38は、総電子数算出部(TECバイアス補正部34)により算出された総電子数に基づいて3次元電子密度推定空間設定部36により分割された領域毎の電子密度を推定する。この電子密度推定部38の動作は、本発明の電子密度推定ステップに対応する。
【0054】
なお、本実施例において、電子密度推定部38は、MARTあるいはGSVD手法により電子密度を推定するが、本発明を実現するためには、必ずしもいずれかの手法を使う必要はなく、衛星信号受信機が得た観測値に基づいて暫定的に各領域の電子密度を推定できればよい。
【0055】
MARTの手法は、上述したように非特許文献5,6等に記載されているが、ここで簡潔に説明する。
【数6】

【0056】
(10)式は、MARTのアルゴリズムを示す式である。aは、(9)式のA行列i行ベクトルを表す。また、aijは、A行列のij成分を表す。Xは、X列行列のj成分を表す。さらに、Eは、(9)式のE行列のi成分を表している。また、λは、収束の度合いを決めるパラメータである。
【0057】
電子密度推定部38は、(10)式に示すMARTのアルゴリズムを使用して、(10)式の値が収束するまで繰り返し処理を実施し、3次元電子密度推定空間設定部36により分割された領域毎の電子密度を推定する。その際に、電子密度推定部38は、TECバイアス補正部34により出力されたTECの値をE行列の値に使用するとともに、電子密度モデル算出部37により算出されたモデル値をXの初期値として使用する。
【0058】
なお、電子密度推定部38は、信号が全く通過していない空間の電子密度を推定することはできない。
【0059】
GSVD手法を採用する場合には、電子密度推定部38は、以下に示す方式で電子密度を求める。なお、GSVD手法については、非特許文献7に記載がある。
【数7】

【0060】
(11)式のLは、Tikhonov行列と呼ばれる。0次のTikhonov正規化では、Lを単位行列とする。αは、正規化パラメータで、α>0とする。また、Lは、正値2次形式行列(固有値が全て正)である。
【0061】
(11)式を満足する電子密度Xは、次式で表すことができる。
【数8】

【0062】
なお、電子密度推定部38は、AX=Eの式における電子密度X(X1〜X12)を推定する際に、本発明の発明者により発明され既に出願されている特願2009−297309に記載された手法を用いて電子密度を推定してもよい。
【0063】
一方、信号処理部1は、既知の送信源(放送局等)により送信された電波の受信状況を確認する(ステップS4)。なお、信号処理部1の代わりに電子密度推定装置32内のホップ数算出部39が行ってもよい。その場合には、ホップ数算出部39は、演算部15から受信した信号強度等の情報を受け取って当該処理を行う。本実施例において、この受信状況の確認処理を行うのは演算部15であるものとする。また、ステップS4以降の処理を行うホップ数算出部39、判定部41、及び修正部43は、その作用効果において本発明の特徴を示すものであり、従来手法と異なる点である。
【0064】
アレイアンテナ11は、既知の送信源により送信される電波を受信する。演算部15は、A/D変換部14−1〜14−nの各々により出力されたデジタル受信信号に基づいて、アレイアンテナ11で受信した電波の到来方位と仰角とを算出する(本発明の第1算出ステップ)とともに、例えば、信号強度等の情報に基づいて次に示す(13)式及び(14)式を満たすか否かを判断することにより受信状況の確認を行う。
【数9】

【0065】
ここで、Prcvは、アレイアンテナ11が受信した電波の信号強度である。受信した信号強度Prcvが予め設定したしきい値ε以上であるか否か判断される。
【0066】
また、Azrcvは、受信した電波の到来方位である。さらに、Azは、既知の送信源の位置とアレイアンテナ11の設置位置とで決まる幾何学的方位であり、次式に示す球面三角法により求めたものである。
【数10】

【0067】
ここで、θsndは、送信源余緯度である。また、θrcvは、受信源余緯度である。λsndは、送信源経度である。λrcvは、受信源経度である。
【0068】
演算部15は、(14)式に示すように受信した信号の方位Azrcvが、幾何学的方位(Az)±εAZの範囲内にあるか否かを判断する。なお、εAZは、予め演算部15に設定された方位に対するしきい値である。
【0069】
演算部15は、受信した電波に基づく情報が(13)式及び(14)式を満たす場合には、既知の送信源から送信された電波を受信できたと判断する。ただし、既知の送信源の送信周波数が他の送信源で使用されていない場合や他の送信源の信号を受信していないことが明らかである場合等には、必ずしも方位の条件は必須ではなく、演算部15は(13)式のみを用いて受信状況の確認を行うことができる。
【0070】
受信した電波に基づく情報が(13)式及び(14)式を満たさない場合には、演算部15は、既知の送信源から送信された電波を受信できていないと判断し、既知送信源位置による処理は行わない(ステップS5)。
【0071】
受信した電波に基づく情報が(13)式及び(14)式を満たす場合には、演算部15は、既知の送信源から送信された電波を受信できていると判断し、その旨をホップ数算出部39に伝達する。ホップ数算出部39は、電離層TEC・バイアス推定・3次元電子密度推定処理系30により算出された電波の到来方位及び仰角と電子密度推定部38により推定された電子密度とに基づいて、アレイアンテナ11が受信した電波のホップ数を算出する(ステップS6)。このホップ数算出部39の動作は、本発明の第2算出ステップに対応する。
【0072】
図5は、本実施例の電子密度推定装置32内のホップ数算出部39によるホップ数の算出方法を説明する図である。図5に示すように、既知の送信源から送信された電波は、所定の電子密度を有する電離層と地上との間で反射し、ホップするようにして受信位置(アレイアンテナ11)に到達する。
【0073】
ホップ数算出部39によるホップ数の算出方法には様々なものが考えられるので、いくつかの方法について例示する。まず、簡易伝搬距離推定について説明する。受信位置(すなわちアレイアンテナ11の設置位置)で観測した信号の受信仰角をθとする。受信位置上空の電子密度で仰角θにより反射する際の電子密度は、非特許文献11に記載の方式から以下に示す式で推定される。
【数11】

【0074】
ここで、hは高さを示す。εは真空の誘電率を示す。mは電子の質量を示す。fは受信(送信)周波数を示す。qは電子の電荷を示す。ホップ数算出部39は、電子密度推定部38により推定された電子密度に基づいて、(17)式で求められた密度ρに相当する電子密度を有する高さhを特定する。なお、ホップ数算出部39は、高さhを特定する際に、非特許文献1〜3に関連した電子密度モデルを参照してもよいし、衛星を使用したトモグラフィ手法による電子密度モデルを参照してもよい。
【0075】
次に、ホップ数算出部39は、特定した高さhをhrとして使用し、概略の受信位置からの1ホップ推定距離d(図5参照)を求める。反射点高さhrと概略距離dとは、次に示す(18)式あるいは(19)式のような関係を有する。
【数12】

【0076】
ここで、Rは受信位置における地球半径を示す。(18)式は、電離層と地表が平面である場合の簡易的なモデルである。一方、(19)式は、電離層が平面と考え、地表の湾曲を考慮したモデルである。
【0077】
ホップ数算出部39は、複数ホップを想定し、(18)式あるいは(19)式を用いて概略距離dを求め、ホップ数倍する。
【0078】
D=nhop×d …(20)
ここで、nhopはホップ数である。この時点でホップ数nhopは不明であるため、ホップ数算出部39は、nhopに1,2,3,…の数値を仮に代入し、それぞれについてのDの値を算出する。
【0079】
次に、ホップ数算出部39は、計算した距離Dと、球面三角法による(21)式及び(22)式を連立して解くことにより、距離Dの経度(λ1)及び緯度(φ1)を計算する。
【数13】

【0080】
次に、ホップ数算出部39は、以下に示す球面三角法の式を用いて、距離Dに相当する位置と既知の送信位置との間の距離を求める。
【数14】

【0081】
ホップ数算出部39は、ホップ数nhopの数値を例えば1,2,3,…と変えた場合のそれぞれについて、(23)式に示す距離Xnhopを算出し、算出したXnhopの中で最小の距離となるものを選択し、そのXnhopに相当するホップ数nhopを実際に伝搬したホップ数とする。なお、ホップ数算出部39は、算出したホップ数の信頼性を高めるために、複数の観測結果を使用してホップ数毎の平均を求め、平均値を使用してXnhopの最小値を求めることもできる。この場合には、ホップ数算出部39は、次式を用いて計算する。
【数15】

【0082】
ホップ数算出部39は、ホップ数nhop毎に(24)式の平均値を求める。ここで、(24)式の左辺の添え字mは、平均を示す。ホップ数算出部39は、この平均値の中で最小となる距離に相当するホップ数nhopを実際に伝搬してきたホップ数とする。
【0083】
以上述べた簡易伝搬距離推定の手法の代わりに、ホップ数算出部39は、レイトレーシングの手法を用いてホップ数を算出することもできる。この場合に、ホップ数算出部39は、非特許文献12,13のレイトレーシング方式と電子密度モデル(非特許文献1〜3に関連した電子密度モデル、あるいは衛星を使用したトモグラフィ手法により求めたもの)を使用し、受信位置から複数ホップ数分レイトレーシングを実施する。その後、ホップ数算出部39は、レイトレーシングの結果推定された推定送信源位置と(21)〜(24)式とを使用し、最小となる距離に相当するホップ数を実際に伝搬してきたホップ数とする。
【0084】
また、ホップ数算出部39は、簡易伝搬距離推定及びレイトレーシングの手法の代わりに、受信電力による推定でホップ数を算出することもできる。この場合に、ホップ数算出部39は、非特許文献14による信号強度の手法により推定する。
【数16】

【0085】
ここで、Eは、非特許文献14の信号強度式であり、ホップ数に関係しない項である。fは送信(受信)周波数を示す。p´(nhop)は仮想斜距離(virtual slant range)を示す。Li(nhop)は吸収による減衰を示す。Ptは既知送信局の送信電力を示す。Gtは送信アンテナゲインを示す。LmはMUF(Maximum Useable Frequency)を超える場合の減衰を示す。Lhはオーロラあるいは他の減衰項を示す。Lzはその他の減衰項を示す。
【0086】
ホップ数算出部39は、予め(26)式の値を計算し、その後、ホップ数を変えてEr式の値を(25)式から計算する。このとき、Erは、サンプルの中央値を表している。ホップ数算出部39は、何回かの受信強度の観測結果を演算部15から受け取り、その中央値を代表値として(25)式で計算したErと比較する。その結果、ホップ数算出部39は、最も近い信号強度に相当するホップ数を実際に伝搬してきたホップ数とする。
【0087】
次に判定部41(あるいはホップ数算出部39又は修正部43)は、レイトレーシングを実施し、既知の送信源位置からの電波を仮想的に追跡する(ステップS7)。図5の実線部分は、実際の電波の伝搬経路を示す。一方、1点鎖線は、計算した電子密度モデル(非特許文献1〜3に関連した電子密度モデル、あるいは衛星を使用したトモグラフィ手法により求められたもの)を使用したレイトレーシングによる伝搬経路を示す。
【0088】
また、判定部41は、ホップ数算出部39により算出されたホップ数に基づいて算出したアレイアンテナ11から送信源(電波を送信した放送局)までの距離と実際の距離との差の大きさが所定のしきい値ε以下であるか否かを判断し、所定のしきい値ε以下である場合に、電子密度推定部38による推定結果を正しいと判定する(ステップS8)。この判定部41の動作は、本発明の判定ステップに対応する。
【0089】
すなわち、判定部41は、ホップ数を特定した場合に(23)式で算出したXnhopの値がXnhop≦εを満たすか否かを判断し、満たす場合には電子密度推定部38による推定結果を正しい(真値に近い)と判定して終了する。
【0090】
nhop≦εを満たさない場合には、判定部41は、修正部43に通過経路の電子密度を修正させる(ステップS9)。すなわち、修正部43は、判定部41により電子密度推定部38の推定結果が正しいと判断されない場合(Xnhop>εの場合)に、アレイアンテナ11が受信した電波の伝搬経路における電子密度を修正する。この修正部43の動作は、本発明の修正ステップに対応する。
【0091】
修正部43は、アレイアンテナ11が受信した電波の伝搬経路における電子密度を修正する際に、伝搬経路としてレイトレーシング手法により特定されたものを用いる。修正部43は、レイトレーシングを自ら実施してもよいし、ホップ数算出部39あるいは判定部41からレイトレーシング結果を受け取ってもよい。
【0092】
修正部43は、レイトレーシングで計算した位置と受信した位置との距離をDrayとし、既知の受信位置と送信局位置との距離をDとして差分を計算する。
【0093】
ΔD=Dray−D …(27)
修正部43は、計算した差分ΔDに基づいて、レイトレーシングにより計算した伝搬経路の電子密度を修正する。ΔD>0の場合には、レイトレーシングによる伝搬経路は、実際の伝搬経路よりも高い場所で反射した経路になっていると考えられる。逆に、ΔD<0の場合には、レイトレーシングによる伝搬経路は、実際の伝搬経路よりも低い場所で反射した経路となっていると考えられる。
【0094】
電子密度を修正する方法は、一定値割合で変化させる方法と、ΔDの大きさに応じて変化させる方法とが考えられる。
【数17】

【0095】
(28)式は、一定割合で変化させる方式である。また、(29)式は、ΔDの大きさに依存して変化させる方式である。ここで、φは緯度を示す。λは経度を示す。hは高さを示す。また、α,βは、変化率に対する係数を表しており、修正部に予め設定されている値である。
【0096】
修正部43は、伝搬経路上の電子密度を(28)式あるいは(29)式を用いて修正し、修正した結果を電子密度推定部38に出力する。電子密度推定部38は、修正部43により修正された電子密度に基づいて再度領域毎の電子密度を推定する。これにより、電子密度推定装置32は、ステップS8においてXnhop≦εを満たすまで、ステップS3〜S9の工程を繰り返す。
【0097】
上述のとおり、本発明の実施例1の形態に係る電子密度推定装置及び電子密度推定方法によれば、衛星から送信される信号による通過経路の電離層電子密度分布を推定する際に、既知の放送局からの信号を利用して3次元の電離層電子密度分布をより正確に推定することができる。
【0098】
すなわち、本実施例の電子密度推定装置及び電子密度推定方法は、既知の送信源(放送局等)から放送される信号の伝搬経路を推定し、経路上の電子密度を求めるので、従来の衛星観測による3次元電子密度分布のデータを補強し修正することができる。特に、本発明は、既知の送信源として日本国内の放送局のみならず、海外の放送局等も利用することができるので、より容易且つ正確に3次元電子密度分布を推定することができる。
【0099】
さらに、電子密度推定部38は、修正部43による修正結果を反映させて、再度電子密度の推定を行うので、より正確な3次元電子密度分布を求めることができる。
【0100】
また、修正部43は、電離層電子密度モデルを使用したレイトレーシング手法により特定された伝搬経路上の電子密度を修正するので、既知の送信源による電波を電子密度推定に十分活用することができ、衛星観測結果のみを用いた場合の観測値の少なさという問題点を解決することができる。
【0101】
さらに、本実施例の電子密度推定装置及び電子密度推定方法は、既知の送信源の位置と計算による送信源の位置との差が所定のしきい値ε以下になったか否かを判断基準とし、ステップS3からS9までの処理を繰り返すことにより3次元電子密度推定の近似精度を向上させることができる。
【産業上の利用可能性】
【0102】
本発明に係る電子密度推定装置及び電子密度推定方法は、衛星から送信される複数の異なる周波数信号による通過経路の電離層電子密度分布を推定する電子密度推定装置に利用可能である。
【符号の説明】
【0103】
1 信号処理部
10a、10b、10c 測位衛星
11 アレイアンテナ
12 アンプ
13 周波数変換部
14 A/D変換部
15 演算部
20 衛星信号処理系
22 アンテナ
24 衛星信号受信機
26 衛星信号処理装置
30 電離層TEC・バイアス推定・3次元電子密度推定処理系
32 電子密度推定装置
33 衛星TEC算出部
34 TECバイアス補正部
35 周波数間バイアス推定部
36 3次元電子密度推定空間設定部
37 電子密度モデル算出部
38 電子密度推定部
39 ホップ数算出部
40 インターネットデータ処理系
41 判定部
42 ルータ
43 修正部
44 GEONET収集データ処理装置
46 外部インターネット網
50 通信回線
S1,S2,S3,S4 測位衛星
Rx1,Rx2,Rx3 衛星信号受信機

【特許請求の範囲】
【請求項1】
複数の測位衛星から送信される衛星信号を受信する1以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定装置であって、
既知の送信源により送信される電波を受信するアレイアンテナと、
前記アレイアンテナが受信した電波の到来方位と仰角とを算出する第1算出部と、
前記1以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出部と、
電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する空間設定部と、
前記総電子数算出部により算出された総電子数に基づいて前記空間設定部により分割された領域毎の電子密度を推定する電子密度推定部と、
前記第1算出部により算出された電波の到来方位及び仰角と前記電子密度推定部により推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第2算出部と、
前記第2算出部により算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定部による推定結果を正しいと判定する判定部と、
前記判定部により前記電子密度推定部の推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正部と、
を備えることを特徴とする電子密度推定装置。
【請求項2】
前記電子密度推定部は、前記修正部により修正された電子密度に基づいて再度前記領域毎の電子密度を推定することを特徴とする請求項1記載の電子密度推定装置。
【請求項3】
前記修正部は、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する際に、前記伝搬経路としてレイトレーシング手法により特定されたものを用いることを特徴とする請求項1又は請求項2記載の電子密度推定装置。
【請求項4】
複数の測位衛星から送信される衛星信号を受信する1以上の受信機を有し、前記衛星信号に含まれる測位情報を用いて位置情報を得る測位システムを利用する電子密度推定方法であって、
アレイアンテナで既知の送信源により送信される電波を受信し、受信した電波の到来方位と仰角とを算出する第1算出ステップと、
前記1以上の受信機の各々から前記複数の測位衛星の各々に対する擬似距離およびキャリア位相擬似距離に基づき前記衛星信号の通過経路の総電子数を算出する総電子数算出ステップと、
電子密度を推定する3次元空間を設定し、当該3次元空間を複数の領域に分割する空間設定ステップと、
前記総電子数算出ステップにより算出された総電子数に基づいて前記空間設定ステップにより分割された領域毎の電子密度を推定する電子密度推定ステップと、
前記第1算出ステップにより算出された電波の到来方位及び仰角と前記電子密度推定ステップにより推定された電子密度とに基づいて前記アレイアンテナが受信した電波のホップ数を算出する第2算出ステップと、
前記第2算出ステップにより算出されたホップ数に基づいて算出した前記アレイアンテナから前記送信源までの距離と実際の距離との差が所定のしきい値以下である場合に前記電子密度推定ステップによる推定結果を正しいと判定する判定ステップと、
前記判定ステップにより前記電子密度推定ステップの推定結果が正しいと判断されない場合に、前記アレイアンテナが受信した電波の伝搬経路における電子密度を修正する修正ステップと、
を備えることを特徴とする電子密度推定方法。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate


【公開番号】特開2011−185665(P2011−185665A)
【公開日】平成23年9月22日(2011.9.22)
【国際特許分類】
【出願番号】特願2010−49484(P2010−49484)
【出願日】平成22年3月5日(2010.3.5)
【出願人】(000003078)株式会社東芝 (54,554)
【Fターム(参考)】