(57)【要約】【課題】  既知信号を用いずに伝搬路を推定すること。【解決手段】  本発明は、伝搬路のインパルス応答が有限長である条件のもとに伝搬路のインパルス応答を推定する方法を提供する。受信信号は、送信信号と伝播路の伝達関数の積で与えられる。したがって、受信信号に送信信号の逆の操作を行えば伝搬路の伝達関数が復元できる。まず、ＤＦＴ演算部１００は、伝搬路を介したベースバンド信号を送信側で供給するシンボル間隔の1/Mの間隔でサンプルする。さらに、このサンプル値系列を用いて受信信号に対してDFT演算を行う。つぎに、次数推定部２００は、伝搬路のインパルス応答の長さ、すなわち伝播路の次数を求める。最後に、インパルス応答推定部３００は、次数推定部２００で求めた伝播路の次数LとＤＦＴ演算部１００で求めたDFT演算の結果を使って、伝播路のインパルス応答を推定する。【選択図】  図１２
2007020126.tif 000002
【特許請求の範囲】
【請求項１】
伝搬路を介して信号を受信するステップと、
前記受信するステップで受信した信号を離散フーリエ変換(DFT)領域に変換した値と、前記受信した信号に送信信号の逆数に対応する未知数を乗じ、さらに、乗じて得た値を時間領域に変換することによって、インパルス応答を与える方程式を導出するステップと、
伝搬路を次数ＬのＦＩＲフィルタと仮定し、前記インパルス応答のうち、Ｌ番目以降の成分を０とすることによって、前記方程式の解を導出して前記送信信号を推定するステップと、
前記推定された送信信号と、前記受信した信号とを前記方程式に代入して、前記インパルス応答を導出し、伝搬路を推定するステップと、
を含むことを特徴とする伝搬路推定方法。
【請求項２】
前記方程式を導出するステップは、
前記受信した信号を高頻度サンプリングした信号の係数ベクトルを導出するステップと、
前記係数ベクトルに対し、Ｌ以上の正の整数であるＫを変換長としたDFTを施して、ベクトルＸ［ｋ］（ｋ＝０、１、・・・、Ｋ−１）を求めるステップと、
前記Ｘ［ｋ］と、前記受信した信号の逆数に対応する値を未知数Ｙ（ｋ）とし、それを乗じた系列に対し、逆離散フーリエ変換を施して、
【数１】
2007020126.tif 000043
で示される前記インパルス応答ｈ（ｎ）についての方程式を導出するステップと、
を含むことを特徴とする請求項１に記載の伝搬路推定方法。
【請求項３】
前記送信信号を推定するステップは、
前記方程式の左辺において、Ｌ番目以降のインパルス応答に対応する、
【数２】
2007020126.tif 000044
で定義されるＡＬと、Ｙ＝［Ｙ（０）  Ｙ（１）  ・・・  Ｙ（Ｋ−１）］Ｔで定義されるＹとを用いて得たＡＬＹをインパルス応答のL番目以降が０に対応させ、
【数３】
2007020126.tif 000045
で与えられる方程式の解を導出することによって、前記送信信号を推定するステップと、
を含むことを特徴とする請求項１に記載の伝搬路推定方法。
【請求項４】
前記送信信号を推定するステップは、
【数４】
2007020126.tif 000046
で定義されるｙと、
【数５】
2007020126.tif 000047
で定義されるＡＬと、前記ＡＬのｋ番目の列ベクトルをαｋで表したとき、
【数６】
2007020126.tif 000048
で定義されるＢＬと、
【数７】
2007020126.tif 000049
で定義されるＰＬとで構成され、
【数８】
2007020126.tif 000050
として定義した方程式を解くことによって、前記送信信号を推定することを特徴とする請求項１に記載の伝搬路推定方法。
【請求項５】
伝搬路を介して信号を受信する受信部と、
前記受信するステップで受信した信号をDFT領域に変換した値と、前記受信した信号に対応する送信信号の逆数と、を乗じ、さらに、乗じて得た値を時間領域に変換することによって、インパルス応答の方程式を導出部と、
伝搬路を次数ＬのＦＩＲフィルタと仮定し、前記インパルス応答のうち、Ｌ番目以降の成分を０とすることによって、前記方程式の解を導出して前記送信信号を推定する推定部と、
前記推定された送信信号と、前記受信した信号とを前記方程式に代入して、前記インパルス応答のうち前記ＦＩＲフィルタにかかる成分を導出し、伝搬路を推定する推定部と、
を備えることを特徴とする伝搬路推定装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、推定技術に関し、特に無線伝搬路の状態を推定する伝搬路推定方法および伝搬路推定装置に関する。
【背景技術】
【０００２】
一般に、無線通信においては、伝搬路において多重通信経路による線形歪が生じ、受信状態が劣化するといわれている。受信状態の劣化を回避するためには、伝搬路の状態を推定し、その状態に応じて、送信側、もしくは受信側において、通信の対象となる信号に対して、伝搬路で発生した歪を除去する処理を施す必要がある。その為には、伝送路の推定が必要である。従来、伝搬路の状態を推定する方法としては、信号に既知信号を含めて送信し、受信側で、既知信号がどの程度劣化したかを検査することによって、伝搬路の状態を推定していた（たとえば、特許文献１参照。）。
【特許文献１】特開２００３−３３２９４３号公報  （第１図）
【発明の開示】
【発明が解決しようとする課題】
【０００３】
しかしながら、送信信号に既知信号を含めることによって、送信すべき信号の量が圧迫され、伝送効率が低下するといった課題があった。
本発明はこうした状況に鑑みてなされたものであり、その目的は、既知信号を用いずに伝搬路を推定できる伝搬路推定方法および伝搬路推定装置を提供することにある。
【課題を解決するための手段】
【０００４】
上記課題を解決するために、本発明の伝搬路推定方法は、伝搬路を介して信号を受信するステップと、受信するステップで受信した信号をDFT領域に変換した値と、受信した信号に、送信信号の逆数に対応する値を未知数Ｙ（ｋ）とし、それを乗じ、さらに、乗じて得た値を時間領域に変換することによって、インパルス応答の方程式を導出するステップと、伝搬路を次数ＬのＦＩＲフィルタと仮定し、インパルス応答のうち、Ｌ番目以降の成分を０とすることによって、方程式の解を導出して送信信号を推定するステップと、推定された送信信号と、受信した信号から、インパルス応答のうちインパルス応答のうちＦＩＲフィルタにかかる成分を導出し、伝搬路を推定するステップと、を含む。
【０００５】
また、方程式を導出するステップは、受信した信号を高頻度サンプリングした信号の係数ベクトルを導出するステップと、係数ベクトルに対し、Ｌ以上の正の整数であるＫを変換長としたDFTを施して、ベクトルＸ［ｋ］（ｋ＝０、１、・・・、Ｋ−１）を求めるステップと、Ｘ［ｋ］と、それに送信信号の逆数に対応する値を未知数Ｙ（ｋ）を乗じた系列に対し、逆離散フーリエ変換（IDFT）を施して、
【数１】
2007020126.tif 000003
で示される前記インパルス応答ｈ（ｎ）についての方程式を導出するステップと、を含んでもよい。
【０００６】
また、送信信号を推定するステップは、方程式の左辺において、Ｌ番目以降のインパルス応答に対応する、
【数２】
2007020126.tif 000004
で定義されるＡＬと、Ｙ＝［Ｙ（０）  Ｙ（１）  ・・・  Ｙ（Ｋ−１）］Ｔで定義されるＹとを用いて得たＡＬＹをインパルス応答のL番目以降が０に対応させ、
【数３】
2007020126.tif 000005
で与えられる方程式の解を導出することによって、送信信号を推定するステップと、を含んでもよい。
【０００７】
また、送信信号を推定するステップは、
【数４】
2007020126.tif 000006
で定義されるｙと、
【数５】
2007020126.tif 000007
で定義されるＡＬと、ＡＬのｋ番目の列ベクトルをαｋで表したとき、
【数６】
2007020126.tif 000008
で定義されるＢＬと、
【数７】
2007020126.tif 000009
で定義されるＰＬとで構成され、
【数８】
2007020126.tif 000010
として定義された方程式を解くことによって、送信信号を推定してもよい。
【０００８】
この態様によると、Ｌ番目以降のインパルス応答は０になるという特性を利用することによって、既知信号を用いずに、所定の方程式の解にもとづいて伝搬路の伝達関数を求めることができる。また、伝搬路の伝達関数を求める際の計算量を低減することができる。
【０００９】
本発明のさらに別の態様は、伝搬路推定装置である。伝搬路推定装置は、伝搬路を介して信号を受信する受信部と、受信するステップで受信した信号をDFT領域に変換した値と、受信した信号に、送信信号の逆数に対応する値を未知数Ｙ（ｋ）とし、それを乗じ、さらに、乗じて得た値を時間領域に変換することによって、インパルス応答の方程式を得る導出部と、伝搬路を次数ＬのＦＩＲフィルタと仮定し、インパルス応答のうち、Ｌ番目以降の成分を０とすることによって、方程式の解を導出して送信信号を推定する推定部と、推定された送信信号と、受信した信号とを方程式に代入して、インパルス応答のうちＦＩＲフィルタにかかる成分を導出し、伝搬路を推定する推定部と、を備える。
【００１０】
この態様によると、Ｌ番目以降のインパルス応答は０になるという特性を利用することによって、既知信号を用いずに、所定の方程式の解にもとづいて伝搬路の伝達関数を求めることができる。また、伝搬路の伝達関数を求める際の計算量を低減することができる。
【００１１】
なお、以上の構成要素の任意の組合せ、本発明の表現を方法、装置、システム、記録媒体、コンピュータプログラムなどの間で変換したものもまた、本発明の態様として有効である。
【発明の効果】
【００１２】
本発明によれば、既知信号を用いずに伝搬路を推定できる伝搬路推定方法、伝搬路推定装置を提供することができる。
【発明を実施するための最良の形態】
【００１３】
（実施例）
概要
初めに本発明に係るブラインド伝搬路推定（Ｂｌｉｎｄ  Ｃｈａｎｎｅｌ  Ｅｓｔｉｍａｔｉｏｎ）の原理を概説する。ブラインド伝搬路推定は、伝搬路が有限インパルス応答（Ｆｉｎｉｔｅ  Ｉｍｐｕｌｓｅ  Ｒｅｓｐｏｎｓｅ。以下、「ＦＩＲ」と略称する。）であると仮定したときの直交周波数分割多重方式（Ｏｒｔｈｏｇｏｎａｌ  Ｆｒｅｑｕｅｎｃｙ  Ｄｉｖｉｓｉｏｎ  Ｍｕｌｔｉｐｌｅｘ。以下、「ＯＦＤＭ」と略称する。）を用いた送信について、伝搬路の状態を推定するものである。［１］では、本推定方法を適用する伝搬路モデルについて説明する。［２］では、既知の伝搬路の次数とインパルス応答零値法（Impulse Response Zeros Method。以下、「IRZ法」と略称する。）について説明する。［３］では伝搬路の次数を推定する方法を説明する。［４］においては、上述した方法の性能を比較したときのシミュレーション結果を示す。［５］においては、概要のまとめを示す。
【００１４】
［１］伝搬路モデル
ＯＦＤＭとして知られているマルチキャリア変調は、送信側において、逆離散フーリエ変換（IDFT）を用いてシンボルを多重し、受信側において、離散フーリエ変換（DFT）を用いてシンボルを復元する。ＯＦＤＭは、欧州デジタルオーディオ広帯域に関する規格として採用されており、また、さまざまな電化製品に用いられている。
【００１５】
ＯＦＤＭシステムのブラインド推定に使われている従来技術としてサブスペース法が最も重要である。この方法ではまず受信信号からその自己相関行列を数値平均により推定する。受信信号の自己相関行列を特異値分解し、信号部分と雑音部分とに分解する。雑音部分は伝搬路インパルス応答と直交しているという性質を利用して線形連立方程式つくり、それを解くことでインパルス応答を推定する。
【００１６】
他の推定方法として、伝搬路インパルス応答と送信シンボルの両者の関数としてコスト関数を定義し、反復法か最ゆう法のいずれかによって、このコスト関数を最小化することで問題を解決する方法がある。
この発明で想定している通信形態を説明する。一つのブロックで、Ｎ個のシンボルＳ（ｋ）、ｋ＝０〜Ｎ−１、の信号を送信するものとする。これらのシンボルの変換長がＮのＩＤＦＴを取り、その系列のz-変換をS（ｚ）とする。S（ｚ）を低域通過フィルタに通しパルス波形に整形する。それを変調し、整形されたパルスの系列として送信する。図１は、高頻度サンプリングを使ったＯＦＤＭのキャリア周波数帯通信モデルを示す図である。送信側のパルス整形フィルタは、割り当てられた帯域幅に送信信号を制限する役目を果たす。受信機は、受信信号を復調し、伝送レートよりＭ倍早いレートでベースバンド信号をサンプリングし、高頻度サンプリング（Ｆｒａｃｔｉｎｏａｌｌｙ  Ｓａｍｐｌｉｎｇ）された受信信号を得る。
【００１７】
すなわち、送信機ではＳ（ｚ）のサンプルがＴの間隔で低域通過フィルタに供給され、受信機では、復調した後に、Ｔ／Ｍ間隔で受信信号をサンプルしている。受信機に設けられた剰余演算を示すｍｏｄ（１−ｚ（―ＭＮ））の演算子については、後述する。
【００１８】
図２は、図１に対応するベースバンド送信モデルを示す図である。図２においては、伝搬路での歪は線形フィルタと付加雑音であると仮定している。図２において、Ｈ（ｚ）、Ｗ（ｚ）は、それぞれベースバンド伝搬路伝達関数と付加雑音を示す。伝達関数は、パルス整形低域通過フィルタと伝送路の線形フィルタの両者の影響を反映している。上向き矢印とＭを含むブロックは、受信機におけるサンプリングレートを上げるためのアップサンプリング処理を示す。
【００１９】
ここで、伝搬路をＦＩＲとし、信号ブロック間にガード間隔あるいはIDFT出力の最後のサンプルの複製を前段に配置する巡回変換（Ｃｙｃｌｉｃ  Ｐｒｅｆｉｘ。以下、「ＣＰ」と略称する。）を設けることで、ブロック間干渉が生じないようになっているものとする。受信信号を式（１）のベクトルz変換で示すことにする。
【数９】
2007020126.tif 000011
・・・式（１）
ここで、ベクトルx（ｎ）は、式（２）である。
【数１０】
2007020126.tif 000012
・・・式（２）
このベクトルx（ｎ）の定義式で使われている上付きのTは転置を示す。そうすると、受信信号は、式（３）で与えられる。
【数１１】
2007020126.tif 000013
・・・式（３）
【００２０】
この式に於いて、Ｈ（ｚ）およびＷ（ｚ）は、Ｘ（ｚ）と同じように、それぞれ伝送路の伝達関数と付加雑音に対するベクトルz変換である。また、ｍｏｄ（１−ｚ（−Ｎ））の演算処理はz変換に於いて、ｚ（−N）＝１を代入することと同じである。したがって、積Ｈ（ｚ）Ｓ（ｚ）の長さがＮを超えていた場合でも、X（ｚ）は、長さＮを超えることはない。
【００２１】
本発明は、伝搬路がＦＩＲである条件のもとにおけるＯＦＤＭシステムのためのブラインド伝搬路推定方法を提供する。以下、その基本的な考え方の概略を説明する。受信信号は、z変換領域では、付加的な雑音がないとすると、送信信号と伝送路の伝達関数の積で与えられる。したがって、受信信号に送信信号の逆のz変換を掛けて得られたz変換は伝達関数になる。この場合、伝搬路はFIRとしているので、復元した受信信号のインパルス応答のサンプルは、伝搬路のインパルス応答における０の領域において、０となる。この関係を利用して、送信信号の逆のz変換を未知とした、線形連立方程式を立てる。このようにして、ブラインド伝搬路推定の問題を線形連立方程式を解く問題に置き換えて解く。
ＯＦＤＭシステムに従来使われているサブスペース法と比較してみる。このサブスペース法は受信信号の自己相関行列を求め、その行列の特異値分解を求めてブラインド伝搬路推定の問題を解決している。この自己相関行列を数値平均で推定するのに相当数のブロック受信信号が必要であるが、本発明では一つのブロック受信信号だけでよい。また一般に、線形連立方程式を解く計算量は、行列の特異値分解の要する計算量よりも格段に少なくてよい。
【００２２】
［２］剰余多項式を用いたブラインド推定
まず、ボードレートでの伝送路インパルス応答長をＬとする。このインパルス応答長のことを伝搬路の次数とも呼ぶ。このとき伝搬路の伝達関数は式（４）のように書ける。
【数１２】
2007020126.tif 000014
・・・式（４）
ここで、ｈ（ｎ）＝[ｈ（ｎＭ）  ｈ（ｎＭ＋１）  ・・・  ｈ（ｎＭ＋Ｍ−１）]Ｔである。ここでは、伝搬路の次数は既知であると仮定する。なお、この伝搬路の次数を推定する方法については、［３］の所で説明する。
【００２３】
ある正の整数Ｋとして、Ｌ以上であり、Ｎの約数であるように選ぶことにする。すると、１−ｚ（―Ｋ）は１−ｚ（−Ｎ）の因数多項式になる。式（３）の両辺をｍｏｄ（１−ｚ（―Ｋ））によって剰余を計算すると、ＫがＬ以上なので、Ｈ（ｚ）ｍｏｄ（１−ｚ（―Ｋ））はＨ（ｚ）と等しい。これを式（５）で表す。
【数１３】
2007020126.tif 000015
・・・式（５）
ここで、Ｘ０（ｚ）、Ｗ０（ｚ）、及びＳ０（ｚ）は、それぞれ、Ｘ０（ｚ）≡Ｘ（ｚ）ｍｏｄ（１−ｚ（―Ｋ））、Ｗ０（ｚ）≡Ｗ（ｚ）ｍｏｄ（１−ｚ（―Ｋ））、及びＳ０（ｚ）≡Ｓ（ｚ）ｍｏｄ（１−ｚ（―Ｋ））である。Ｘ０（ｚ）の係数ベクトルは、受信信号Ｘ（ｚ）から式（６）によって計算できる。
【数１４】
2007020126.tif 000016
・・・式（６）
【００２４】
多項式（１−ｚ（―Ｋ））は、（１−ｅ（ｊ２πｋ/Ｋ）ｚ（−１））、ｋ＝０〜Ｋ−１の積に等しい。また、ｍｏｄ（１−ｅ（ｊ２πｋ/Ｋ）ｚ（−１））の演算を行うことは、zにｚ＝ｅ（ｊ２πｋ/Ｋ）の代入を行うことと等価になる。これらにより、式（５）は下記の式（７）と等価になる。
【数１５】
2007020126.tif 000017
・・・式（７）
ここで、Ｘ［ｋ］は式（８）で与えられる。
【数１６】
2007020126.tif 000018
・・・式（８）
この式（８）では、ｘ０（ｎ)、ｎ＝０、１、・・・、Ｋ−１は、Ｘ０（ｚ）の係数ベクトルである。また、Ｈ［ｋ］とＷ［ｋ］は、Ｘ［ｋ］の定義と同様の形式で定義されているものとする。送信信号は、Ｓ（ｋ）≡Ｓ（ｚ）ｍｏｄ（１−ｅ（―ｊ２πｋ/Ｎ）)によって与えられるため、ｚ＝ｅ（―ｊ２πｋ/Ｋ）をＳ（ｚ）に代入することによって得られる値は、Ｓ（ｋＮ／Ｋ）となる。ベクトルＸ［ｋ］ｋ＝０、１、・・・、Ｋ−１は、Ｘ0（ｚ）の係数に対して変換長ＫのＤＦＴを取ることで計算できる。
【００２５】
付加的な雑音がないと仮定すると、Ｈ［ｋ］はＸ［ｋ］と式（９）のように関連付けられる。
【数１７】
2007020126.tif 000019
・・・式（９）
ここで、Ｙ（ｋ）＝１／Ｓ（ｋＮ／Ｋ）である。一般的に、信号空間（signal constellation）は、ゼロもしくは無限大を含まないため、Ｙ（ｋ）は、ゼロでない有限値となる。送信シンボルＳ（ｋＮ／Ｋ）は未知なので、Ｙ（ｋ）、ｋ＝０〜Ｋ−１、を未知の値とする。このＹ（ｋ）を求めることが出来れば、式（８）により、伝搬路インパルス応答は、Ｘ［ｋ］Ｙ（ｋ）、ｋ＝０〜Ｋ−１、に対しＩＤＦＴ演算をすることによって得られる。すなわち、伝搬路インパルス応答は式（１０）によって求めることが出来る。
【数１８】
2007020126.tif 000020
・・・式（１０）
【００２６】
伝搬路の伝達関数は式（４）で与えられることより、ｎ＝Ｌ〜Ｋ−１の範囲において、ｈ（ｎ）＝０となる。これらのゼロを使うことにより、未知のＹ（ｋ）、ｋ＝０〜Ｋ−１に関する線形連立方程式を構成する。
【００２７】
その為に、式（１１）で与えられる（Ｋ−Ｌ）Ｍ×Ｋの行列を導入する。
【数１９】
2007020126.tif 000021
・・・式（１１）
すると、ｎ＝Ｌ〜Ｋ−１の範囲において、ｈ（ｎ）＝０となることから、式（１２）に示す線形連立方程式が得られる。
【数２０】
2007020126.tif 000022
・・・式（１２）
ここで、Ｙ＝［Ｙ（０）  Ｙ（１）  ・・・  Ｙ（Ｋ−１）］Ｔである。この式（１２）が、Ｙを求めるための基本方程式となる。
【００２８】
実際の伝播路では付加雑音がある。そのような場合には、式（１２）の連立方程式を解くよりも、ＡＬＹのノルムを最小化するベクトルＹを見つける方が自然である。通常、信号空間にはゼロや無限大の値を含まないように設計されている。Yの要素が取りえる値は信号空間の値の逆数値であり、したがって、ゼロや無限大にならない。そこで、Ｙ（０）＝１と置いて、ノルムを最小化するベクトルＹを求めることにする。ここで、Yから最初の要素Y(0)を除いた（Ｋ−１）次元のベクトルｙを式（１３）のように定義する。
【数２１】
2007020126.tif 000023
・・・式（１３）
さらに、行列ＢＬを式（１４）で定義する。
【数２２】
2007020126.tif 000024
・・・式（１４）
ここで、αｋは、ＡＬのｋ番目の列ベクトルである。次に、（Ｋ−１）×（Ｋ−１）の行列PLを式（１５）で定義する。
【数２３】
2007020126.tif 000025
・・・式（１５）
ここで、＊は、複素共役転置の演算子を表す。Ｙ（０）＝１の条件下においては、||ＡＬＹ||２＝||α０＋ＢＬｙ||２となるが、これを最小化するベクトルｙは、式（１６）の線形連立方程式の解として求められる。
【数２４】
2007020126.tif 000026
・・・式（１６）
この線形連立方程式は、（Ｋ−１）個の未知数と（Ｋ−１）個の式とから成り立っており、通常のサブルーチンによって、容易に解くことができる。
【００２９】
ＰＬの定義により、行列ＰＬの（ｎ―１、ｋ−１）番目（ｎおよびkは１以上でＫ−１以下。）の要素は、式（１７）で与えられる。
【数２５】
2007020126.tif 000027
・・・式（１７）
この式で、積Ｘ［ｎ］＊Ｘ［ｋ］の部分は、Ｌに依存しない。一方、Σによる加算は、受信信号とは独立であり、前もって計算しておくことが出来る。Ｑを、式（１８）で与えられる（Ｋ−１）×（Ｋ−１）の行列とすると、Ｑの（ｎ−１、ｋ−１）番目の要素は、積Ｘ［ｎ］＊Ｘ［ｋ］と等しくなる。したがって、PLの（ｎ−１、ｋ−１）番目の要素は、式（１９）で与えられる。
【数２６】
2007020126.tif 000028
・・・式（１８）
【数２７】
2007020126.tif 000029
・・・式（１９）
ここで、ｑ（ｎー１、ｋ−１）は、行列Ｑの（ｎ−１、ｋ−１）番目の要素である。
また、式（１６）の右辺のベクトルは、式（２０）によって与えられる。
【数２８】
2007020126.tif 000030
・・・式（２０）
【００３０】
式（１９）によって得たＰＬと式（２０）によって得た−ＢＬ＊α０を使って、ｙに関する線形連立方程式（１６）を解くことができ、||ＡＬＹ||２を最小化するベクトルＹを得ることができる。そのベクトルＹの解を使って、伝搬路インパルス応答は、式（１０）で推定することができる。この計算は、Ｘ［ｋ］Ｙ（ｋ）、ｋ＝０〜Ｋ−１に対する長さがＫのＩＤＦＴである。
【００３１】
上述のブラインド伝搬路推定法を、インパルス応答零値法（Impulse Response Zeros Method。以下、「IRZ法」と略称する。）と呼ぶことにする。これまでの議論をまとめると、IRZ法は以下のステップで実現できる。
１．式（６）と式（８）によって、Ｘ（ｚ）から、Ｘ［ｋ］、ｋ＝０〜Ｋ−１を計算する。
２．線形方程式（１６）を解くことによって、式（１２）の連立線形方程式の解Ｙを見つける。
３．伝搬路のインパルス応答ｈ（ｎ）の推定は、ｃを定数として、ｃｇ（ｎ）で与えられる。ここで、ｇ（ｎ）、ｎ＝０〜Ｋ−１は、Ｘ［ｋ］Ｙ（ｋ）、ｋ＝０〜Ｋ−１に対する変換長がＫのＩＤＦＴである。
【００３２】
［３］伝搬路の次数の推定について
通常、伝搬路の次数は既知でないため、それを推定する必要がある。ここでは、伝搬路の次数の推定とIRZ法が適用できる条件の問題について議論する。伝搬路の次数の推定においては、式（１１）で定義した行列でLのかわりにJを使ったＡJを用いる。Ｊは、固定の整数であり、インパルス応答長Ｌとは異なった値でもよい。まず、伝搬路において付加的な雑音がないものと仮定して、ＡＪのランクを評価する。式（１０）から、式（２１）が導かれる。
【数２９】
2007020126.tif 000031
・・・式（２１）
これらの方程式は、式（２２）を与える。
【数３０】
2007020126.tif 000032
・・・式（２２）
ここで、（ｎ＋ｍ）Ｋ≡ｎ＋ｍ  ｍｏｄ（Ｋ）であり、また、［ＤＦＴ］は、長さがＫのＤＦＴ行列である。つぎに、式（２３）で定義される（Ｋ−Ｊ）Ｍ×Ｋの行列CJを導入する。
【数３１】
2007020126.tif 000033
・・・式（２３）
【００３３】
式（１１）と式（２２）によって、ＣＪ［ＤＦＴ］の行列のｋ番目の列は、Ｙ（ｋ）αｋと等しくなる。ここで、αｋは、行列ＡＪのｋ番目の列である。［ＤＦＴ］は正則行列であり、Ｙ（ｋ）、ｋ＝０〜Ｋ−１は、ゼロでも無限大でもないので、ＣＪのランクは、ＡＪのランクと等しくなる。さらに、式（２４）で定義される、別の（Ｋ−Ｊ）Ｍ×Ｋの行列DJを導入する。
【数３２】
2007020126.tif 000034
・・・式（２４）
この行列は、ＣＪの列を入れ替えて得た行列であるため、この行列DJとＡＪのランクは等しくなる。
種々の伝搬路推定法で、その適用条件として伝搬路パラメータ行列(channel parameter matrix)の列が全て独立であることが要求されているが、その行列は式（２５）で定義されている。
【数３３】
2007020126.tif 000035
・・・式（２５）
【００３４】
他の伝搬路推定法で要求されているように、ここでも伝搬路パラメータ行列の列は全て独立と仮定する。この仮定のもとにおいて、ＪがＫ−Ｌ以下でＬ−１以上の領域においては、ＤＪの最初のＫ−Ｊ＋Ｌ−１列は独立しており、残りのＪ−Ｌ＋１列はゼロベクトルとなる。それゆえに、ＤＪのランクは、Ｋ−Ｊ＋Ｌ−１と等しくなる。また、ＤＪのランクは、ＡＪのランクと等しいため、伝搬路の次数Ｌは、ＡＪのランクから推定することができる。
【００３５】
付加的な雑音があった場合には、伝搬路の次数の推定にＡＪのランクを直接使うことは適していない。このような場合、式（２６）で定義されるＫ×Ｋの行列RJを用いる。
【数３４】
2007020126.tif 000036
・・・式（２６）
ＲＪの（ｎ、ｋ）番目の要素は、式（２７）で与えられる。
【数３５】
2007020126.tif 000037
・・・式（２７）
【００３６】
ＲＪの要素は、式（１９）と同様に、Ｋ×Ｋの行列［Ｘ[０]  Ｘ[１]  Ｘ[Ｋ−１]］＊［Ｘ[０]  Ｘ[１]  Ｘ[Ｋ−１]］から計算される。行列ＲＪは、エルミート行列（Ｈｅｒｍｉｔｉａｎ  Ｍａｔｒｉｘ）であるため、この行列の固有値は正となる。ＲＪの固有値をλＪ（ｎ）、ｎ＝０〜Ｋ−１、とする。ここで、λＪ（０）≦λＪ（１）≦・・≦λＪ（Ｋ−１）を満たすような順に並べておく。雑音電力がインパルス応答に比べて小さいとすると、λＪ（ｎ）≦δ、０≦ｎ≦Ｊ−Ｌで、残りのＫ−Ｊ＋Ｌ−１の固有値はδより大きい、となるような正定数δが存在する。この定数δは、雑音電力より計算することができる。ｎ０をλＪ（ｎ０）≦δを満たす最も大きな整数とすると、伝搬路の次数は、Ｌ＝Ｊ−ｎ０によって推定できる。
【００３７】
この固有値λＪ（ｎ）、ｎ＝０・・Ｋ−１、の系列は、ｎ０＝Ｊ−Ｌにおいて、急激に増加するが、この特徴を利用して伝搬路の次数を推定することも出来る。その為に、系列log(λＪ（ｎ）)に対して式（２８）のように定義される数値2次導関を利用する。数値2次導関数は、急激に増加する場所ｎ０においてピーク値を持つ。従って、そのピーク値をもつ場所ｎ０を求めれば、インパルス応答長は、Ｌ＝Ｊ−ｎ０の関係より推定できる。この方法は、雑音電力から定数δを計算する必要がない点で有利である。
【数３６】
2007020126.tif 000038
・・・式（２８）
【００３８】
Jは固定の整数であるが、Ｊ＝Ｋ／２と決めておくことにする。伝搬路の次数の推定アルゴリズムは下記のようにまとめられる。
１．Ｋとして、２Ｌ以上であり、かつ、Ｎの約数である整数を用い、ＲＫ／２の固有値を計算する。
２．固有値の系列から、その数値２次導関数ｄ２（log（λＪ（ｎ）））のピークｎ０を見つける。
３．Ｌ＝Ｋ／２―ｎ０により伝搬路の次数を推定する。
【００３９】
ここで、IRZ法が適用できる条件について述べる。
【００４０】
IRZ法の導出とそこで仮定された条件から、全てのシンボルＳ（ｋＮ／Ｋ）、ｋ＝０・・・Ｋ−１、がゼロでも無限大でもなく、また、ＡＬのランクがＫ−１である場合に、IRZ法が適用できることが分かる。ゼロもしくは無限大を含まないように信号空間を設計すれば、第１の条件は満足される。第２の条件に関しては、Ｓ（ｋＮ／Ｋ）、ｋ＝０・・・Ｋ−１、がゼロや無限大でないとき、ＡＬのランクがＤＬのランクと等しいことをすでに確認した。また、ＤＬのランクは、伝搬路パラメータ行列の列が全て独立であり、Ｋが２Ｌ以上のとき、Ｋ−１と等しくなる。この事より、以下のことが言える。
【００４１】
ブロックサイズがＮ、サンプルレートの増加率がＭ、伝搬路の次数がＬであるＯＦＤＭ通信システムに於いて、下記の条件を満たすとき、伝搬路インパルス応答はIRZ法により推定される。
１．IRZ法に用いられる整数ＫはＮの約数であり、かつＫが２Ｌ以上である。
２．全てのシンボルＳ（ｋＮ／Ｋ）、ｋ＝０・・・Ｋ−１、はゼロや無限大でない。
３．伝搬路パラメータ行列の列が全て独立である。
【００４２】
伝搬路パラメータ行列のサイズは、ＬＭ×（２Ｌ―１）であるため、上記条件３は、ＬＭ≧２Ｌ−１であることが要求される。すなわち、Ｍが２以上であることが要求される。
【００４３】
［４］シミュレーション
図３に、本発明の実施の形態にかかるシミュレーションで用いたインパルス応答を示す。この図に示すように、このインパルス応答を持つ伝搬路の次数Lは５である。信号空間としては、直角位相変調（Ｑｕａｄｒａｔｕｒｅ  Ｐｈａｓｅ  Ｓｈｉｆｔ  Ｋｅｙｉｎｇ。以下、「ＱＰＳＫ」と略称する。）を用いることにする。付加雑音としては、電力がσ２＝（１／２）Ｅ｛｜ｗ（ｎ）｜２｝で平均ゼロの白色ガウス雑音とする。DFTの定義によって、信号対雑音比（Ｓｉｇｎａｌ  Ｎｏｉｓｅ  Ｒａｔｉｏ。以下、「ＳＮＲ」と略称する。）は、ＳＮＲ＝１０ｌｏｇ（｜Ｓ（ｋ）｜２／（Ｎσ２））で与えられる。ＱＰＳＫにおいて、この比は、ＳＮＲ＝１０ｌｏｇ（１／（Ｎσ２））で与えられる。
【００４４】
推定法によって得られた解ベクトルをｇとすると、インパルス応答の推定はｃｇとなる。||ｃｇ−ｈ||２を最小化する定数ｃは、ｃ＝ｇ＊ｈ／||ｇ||２  で与えられる。この定数を使って、推定法の性能を、式（２９）で定義される正規化平均二乗誤差法（Ｎｏｒｍａｌｉｚｅｄ  Ｍｅａｎ  Ｓｑｕａｒｅ  Ｅｒｒｏｒ．以下、「ＮＭＳＥ」と略称する。）によって評価するものとする。
【数３７】
2007020126.tif 000039
・・・式（２９）
ここで、ｃｔｇｔはｔ番目の推定値であり、Nｔは推定値の個数である。
【００４５】
例１．シミュレーションによる性能の評価の前に、まず手計算で動作を確認できる単純な例に対してIRZ法を適用してみる。具体的には、付加雑音がないものとし、伝搬路伝達関数をＨ（ｚ）＝ｈ（０）＋ｈ（１）ｚ−１＝［１  １］Ｔ＋［１  −１］Ｔｚ−１とし、Ｎ＝８、Ｍ＝２、Ｌ＝２とした場合について考える。送信シンボルとしてはＳ（０）＝ｊ、Ｓ（ｋ）＝ｅｘｐ（−ｊπｋ／２）、ｋ＝１、２、・・・、７、を用いる。これらのシンボルの変換長８のＤＦＴは、ｓ（ｎ）＝（−１＋ｊ）／８、ｎ＝０、１、３、４、５、６、７とｓ（２）＝（７＋ｊ）／８となる。受信信号は、Ｘ（ｚ）≡Ｈ（ｚ）Ｓ（ｚ）ｍｏｄ（１−ｚ−８）で与えられる。すなわち、式（３０）で与えられる。
【数３８】
2007020126.tif 000040
・・・式（３０）
IRZ法で使う正整数Kとしては、２L以上とNの因数であるという条件を満たすように、Ｋ＝４を選ぶことにする。
【００４６】
まず、Ｘ（ｚ）から伝搬路次数を推定する。式（６）と式（８）を使うと、Ｘ［０］＝［２ｊ  ０］Ｔ、Ｘ［１］＝［−１＋ｊ  −１−ｊ］Ｔ、Ｘ［２］＝［０  ２］Ｔ、Ｘ［３］＝［−１−ｊ  −１＋ｊ］Ｔとなる。式（１１）により、行列ＡＫ／２は、式（３１）で得られる。
【数３９】
2007020126.tif 000041
・・・式（３１）
この行列のランクは３であり、ＲＫ／２＝Ａ＊Ｋ／２ＡＫ／２／（ＭＫ）の一つの固有値だけが０となる。それゆえに、λＫ／２（ｎ）＝０となる最大の整数ｎ０はｎ０＝０となり、伝搬路の次数は、Ｌ＝Ｋ／２−ｎ０＝２と正しく推定される。
【００４７】
次に、IRZ法を使って伝搬路インパルス応答を推定する。ＡＬ＝ＡＫ／２なので、IRZ法で必要なＡＬとしても式（３１）によって与えられるＡＫ／２を用いる。線形方程式（１６）に対応する式は、式（３２）のようになる。
【数４０】
2007020126.tif 000042
・・・式（３２）
【００４８】
この方程式の解はｙ＝［−ｊ  ｊ  −ｊ］Ｔとなるが、これはＹ＝［１  −ｊ  ｊ  −ｊ］Ｔに対応する。式（１０）を用いて、Ｘ［ｋ］Ｙ［ｋ］、ｋ＝０、１、２、３の変換長４のIDFTは、ｇ（０）＝［ｊ  ｊ］Ｔ、ｇ（１）＝［ｊ  −ｊ］Ｔ、ｇ（２）＝［０  ０］Ｔ、ｇ（３）＝［０  ０］Ｔとして得られる。伝搬路インパルス応答と得られた推定値には、ｈ（ｎ）＝−ｊｇ（ｎ）、ｎ＝０、１、２、３の関係が成り立つことが確認できる。
【００４９】
例２．（伝搬路次数の推定）
伝搬路次数の推定法には、RＫ／２の固有値λＫ／２（ｎ）、ｎ＝０、１、・・・、Ｋ−１、の系列の計算が必要である。N=256、Ｋ＝１６とした場合、SNRが２０、３０、４０、５０ｄＢのそれぞれについて、ｌｏｇ（λＫ／２（ｎ））、ｎ＝０、１、・・・、Ｋ−１、の系列は図４で示すようになる。図４は、本発明の実施の形態にかかるシミュレーションで用いた固有値の系列を示す図である。図から見て取れるように、SNRが大きいとき、ｌｏｇ（λＫ／２（３））からｌｏｇ（λＫ／２（４））に移るときに大きく増加している。したがって、伝搬路次数は、n0=3として、Ｌ＝Ｋ／２−n0＝５によって推定される。系列の数値2次導関数ｄ２（ｌｏｇ（λＫ／２（ｎ）））、ｎ＝１、２、・・・、Ｋ−２は、図５で与えられる。数値２次導関数は、SNRが大きいとき、ｎ＝３でピークとなり、伝搬路次数の推定に用いることができることが確認できる。
【００５０】
例３．（ブロック長Nを変化させたときのNMSE）
図６は、本発明の実施の形態にかかる、一つのブロックで送信するシンボル数N（ブロック長）を変化させたときのNMSEの変化の様子を示す図である。図６は、本発明の実施の形態にかかる第１の推定結果を示す図である。このシミュレーションでは、ＳＮＲ=40dＢと設定し、ｌｏｇ２Ｎを４から10まで変化させた。
【００５１】
この図から、IRZ法の性能は、ブロックサイズNとは無関係であることが確認できる。
【００５２】
例４．（ＳＮＲを変化させたときのNMSE）
ブロック長N=256として、ＳＮＲを変化させたときのＮＭＳＥのグラフを図７に示す。図７は、本発明の実施の形態にかかる第２の推定結果を示す図である。ブロック長N=５１２として、ＳＮＲを変化させたときのＮＭＳＥのグラフを図８に示す。図８は、本発明の実施の形態にかかる第３の推定結果を示す図である。
【００５３】
例５．（伝搬路次数の推定誤差による性能の劣化）
先の例においては、性能を、伝搬路推定が正しいとの仮定の下において評価した。この例は、伝搬路次数の推定誤差による性能の劣化の程度について評価する。IRZ法で使う伝搬路次数の推定Lを２から７に変化したときのＮＭＳＥの変化の様子を図９に示す。図９は、本発明の実施の形態にかかる第４の推定結果を示す図である。ＮＭＳＥは、Ｌ＝５のときのみ小さな値であり、伝搬路次数の正しい推定が重要であることが分かる。
【００５４】
例６．（受信シンボルの分布）
信号空間としてＱＰＳＫを使っているので、各送信シンボルS(n)は{1, j, -1, -j}の中のいずれかの値を取る。受信側では、IRZ法で推定された伝搬路のインパルス応答に対する等化器を受信信号に施して、送信シンボルを推定する。その推定されたシンボルを複素平面上にプロットしたときの様子を図１０及び図１１に示す。図１０はSNR=40dB（NMSE=0.0046）の時、図１１はSNR=30dB（NMSE=0.0040）の時の受信シンボルの分布を示す。これらの図は、本発明の実施の形態にかかる受信シンボルの分布を示す図である。
【００５５】
［５］概要のまとめ
本発明では、伝搬路がＦＩＲと仮定したＯＦＤＭ送信システムのためのブラインド伝搬路推定法を導出した。付加雑音がないと仮定すると、受信信号のz変換は送信信号のz変換と伝搬路の伝達関数の積として与えられる。したがって、送信信号のz変換の逆を受信信号のz変換に掛けて得られた式は伝搬路の伝達関数になる。伝搬路がＦＩＲ、すなわち、伝搬路のインパルス応答は有限長と仮定しているので、送信信号のz変換の逆を受信信号のz変換に掛けて得られた式で定まるインパルス応答に於いて、実際のインパルス応答長以降に対応する部分については０にならなければならない。本発明で提案するIRZ法は、上記の観測に基づき、送信信号のz変換の逆を未知とした線形方程式を立てそれを解くことで、ブラインド伝搬路推定問題を解決している。サブスペース法などの従来の方法では、受信信号から数値平均により受信信号の相互相関行列を近似することや、その行列の特異値分解の操作が必要である。IRZ法ではそれらの操作を必要とせず、計算効率がよい。
【００５６】
具体的構成
以上述べた概要に基づき、本実施の形態の具体的構成を次に述べる。図１２は、本発明の実施の形態にかかる伝搬路推定装置の構成例を示す図である。図１２に示す伝搬路推定装置は、ＤＦＴ演算部１００と、次数推定部２００と、インパルス応答推定部３００とを含む。ＤＦＴ演算部１００は、受信信号に対し、ダウンサンプル処理、ＤＦＴ演算などを行う。次数推定部２００は、ＤＦＴ演算部１００の結果を用いて、伝搬路の次数Ｌを推定する。インパルス応答推定部３００は、ＤＦＴ演算部１００の結果と、次数推定部２００で推定した伝搬路次数Ｌを用いてインパルス応答ｈ（ｎ）を推定する。なお、ＤＦＴ演算部１００とインパルス応答推定部３００の基本となる動作原理は、概要の［２］で述べたとおりである。また、次数推定部２００の基本となる動作原理は、概要の［３］で述べたとおりである。
【００５７】
図１３は、図１２のＤＦＴ演算部１００の構成例を示す図である。ＤＦＴ演算部１００は、ダウンサンプリング部１０で代表される第１ダウンサンプリング部１１、第２ダウンサンプリング部１２、第Ｍダウンサンプリング部１３と、ベクトル表現部２０と、変換部２２と、Ｋ点ＤＦＴ部２４とを含む。ＤＦＴ演算部１００は、受信側で、復調したベースバンド信号をT/Mの間隔でサンプルして得られた信号からX(z)、伝搬路の次数推定やインパルス応答の推定の際に必要な受信信号のDFTベクトルＸ［ｋ］ｋ＝０、１、・・・、Ｋ−１を求めるまでの手順をブロック図として示している。ダウンサンプリング部１０は、下向き矢印とＭを含むブロックは、サンプリングレートを下げるためのダウンサンプリング処理を示す。m番目のダウンサンプリングとその直前のzmによる掛け算により、その出力は図で示してあるようにx(nM+m)となる。ベクトル表現部２０は、式（１）と式（２）に従って、ダウンサンプリングの出力を長さがNのベクトルz変換表現X(z)に置き換える操作を示す。変換部２２は、式（６）に従って長さNの信号X(z)から長さKの受信信号X0(z)を求める操作を示す。最後に、Ｋ点ＤＦＴ部２４は、式（８）で定義される変換長KのDFT演算を行う操作を示す。
【００５８】
図１４は、図１２の次数推定部２００の構成例を示す図である。次数推定部２００は、行列計算部２６と、固有値計算部２８と、次数導出部３０とを含む。次数推定部２００は、図１３のＤＦＴ演算部１００で得られたＸ［ｋ］ｋ＝０、１、・・・、Ｋ−１を使って伝播路の次数を推定するまでの手順をブロック図として示している。行列計算部２６では、式（１１）に於いて、Lの代わりにK/2を使って行列AK/2を求める操作を示す。固有値計算部２８は、式（２６）に従いK×Kの行列RK/2を求め、その行列の固有値を求め、それらを小さい順に配列する操作を示す。次数導出部３０では固有値系列から式（２８）で定義される数値２次導関数のピーク値を与える場所n0を求め、その値からL=K/2-n0により次数の推定を行う。
【００５９】
図１５は、図１２のインパルス応答推定部３００の構成例を示す図である。インパルス応答推定部３００は、ベクトル計算部３２と、方程式求解部３４と、インパルス応答導出部３６とを含む。また、インパルス応答推定部３００は、次数推定部２００で推定した伝搬路の次数の推定LとＸ［ｋ］ｋ＝０、１、・・・、Ｋ−１を使って伝播路のインパルス応答を推定するまでの手順をブロック図として示している。ベクトル計算部３２は、式（１８）と式（１９）を使ってPLを求め、式（２０）を使ってBLa0を求める操作を示す。方程式求解部３４は、線形連立方程式の求解部では行列PL及びベクトル-BLa0で定まる線形連立方程式（１６）の解y=[Y(1) Y(2) … Y(K-1)]Tを求め、その解よりY=[1 Y(1) … Y(K-1)]Tを求める操作を示す。インパルス応答導出部３６ではＸ［ｋ］Ｙ（ｋ）、ｋ＝０〜Ｋ−１のＩＤＦＴを取り、それをｃｇ（ｎ）、n＝０〜Ｋ−１とし、ｃをある定数としてインパルス応答をｈ（ｎ）=ｃｇ（ｎ）、n＝０〜Ｋ−１により推定する操作を示す。
【００６０】
以上の構成、手順により、既知信号を用いずに、伝搬路の伝達関数を推定することができる。本発明の実施の形態により得られる主な性能、効果は、概要の［４］などで述べたとおりである。図１２のＤＦＴ演算部１００における「受信信号のベクトルDFT計算の手順」で必要な演算は、実質的にはDFTの演算の所のみで、それ以前は単にサンプル値の並べ替えと加算の処理である。図１２の次数推定部２００における「伝播路の次数推定の手順」では、K×K行列の固有値を求める処理が演算の付加の大部分を占める。従来の手法では数値平均法を使って受信信号の自己相関行列を推定し、その行列の特異値分解により得られた信号成分と雑音成分の分布より次数の推定を行っている。数値平均を取る際に、多数の受信信号が必要であるが、提案する手法では１つの受信信号のみで、伝播路の次数推定を行って、計算量のみならず推定に要する時間遅れの面でも有利である。図１２のインパルス応答推定部３００における「インパルス応答推定の手順」では、K-1次の線形連立方程式を解く処理が演算の付加の大部分を占める。従来の手法のこの部分に相当する処理では、行列の固有値を求めることが必要であるが、線形連立方程式を解くための演算量は固有値を求めるための演算量数よりも一般的に少ない。
【００６１】
以上、本発明を実施の形態をもとに説明した。この実施の形態は例示であり、それらの各構成要素や各処理プロセスの組合せにいろいろな変形例が可能なこと、またそうした変形例も本発明の範囲にあることは当業者に理解されるところである。
【図面の簡単な説明】
【００６２】
【図１】本発明の実施の形態にかかる高頻度サンプリングを使ったＯＦＤＭのキャリア周波数帯通信モデルを示す図である。
【図２】図１に対応するベースバンド送信モデルを示す図である。
【図３】本発明の実施の形態にかかるシミュレーションで用いたインパルス応答を示す図である。
【図４】本発明の実施の形態にかかるシミュレーションで用いた固有値の系列を示す図である。
【図５】本発明の実施の形態のシミュレーションにかかる固有値系列の数値２次導関数を示す図である。
【図６】本発明の実施の形態にかかる第１の推定結果を示す図である。
【図７】本発明の実施の形態にかかる第２の推定結果を示す図である。
【図８】本発明の実施の形態にかかる第３の推定結果を示す図である。
【図９】本発明の実施の形態にかかる第４の推定結果を示す図である。
【図１０】本発明の実施の形態にかかる第１の受信シンボルの分布を示す図である。
【図１１】本発明の実施の形態にかかる第２の受信シンボルの分布を示す図である。
【図１２】本発明の実施の形態にかかる伝搬路推定装置の構成例を示す図である。
【図１３】図１２のＤＦＴ演算部の構成例を示す図である。
【図１４】図１２の次数推定部の構成例を示す図である。
【図１５】図１２のインパルス応答推定部の構成例を示す図である。
【符号の説明】
【００６３】
１０  ダウンサンプリング部、  １１  第１ダウンサンプリング部、  １２  第２ダウンサンプリング部、  １３  第Ｍダウンサンプリング部、  ２０  ベクトル表現部、  ２２  変換部、  ２４  Ｋ点ＤＦＴ部、  ２６  行列計算部、  ２８  固有値計算部、  ３０  次数導出部、  ３２  ベクトル計算部、  ３４  方程式求解部、  ３６  インパルス応答導出部、  １００  ＤＦＴ演算部、  ２００  次数推定部、  ３００  インパルス応答推定部。
【図１】
2007020126.tif 000051
【図２】
2007020126.tif 000052
【図３】
2007020126.tif 000053
【図４】
2007020126.tif 000054
【図５】
2007020126.tif 000055
【図６】
2007020126.tif 000056
【図７】
2007020126.tif 000057
【図８】
2007020126.tif 000058
【図９】
2007020126.tif 000059
【図１０】
2007020126.tif 000060
【図１１】
2007020126.tif 000061
【図１２】
2007020126.tif 000062
【図１３】
2007020126.tif 000063
【図１４】
2007020126.tif 000064
【図１５】
2007020126.tif 000065