TOP > 国内特許検索 > 雑音抑圧装置および雑音抑圧方法 > 明細書

明細書 :雑音抑圧装置および雑音抑圧方法

発行国 日本国特許庁(JP)
公報種別 特許公報(B2)
特許番号 特許第5115952号 (P5115952)
公開番号 特開2008-236270 (P2008-236270A)
登録日 平成24年10月26日(2012.10.26)
発行日 平成25年1月9日(2013.1.9)
公開日 平成20年10月2日(2008.10.2)
発明の名称または考案の名称 雑音抑圧装置および雑音抑圧方法
国際特許分類 H03H  17/04        (2006.01)
G10L  21/0224      (2013.01)
G10L  21/0208      (2013.01)
H03H  21/00        (2006.01)
FI H03H 17/04 633Z
G10L 21/02 101A
G10L 21/02 102B
H03H 21/00
請求項の数または発明の数 9
全頁数 30
出願番号 特願2007-071688 (P2007-071688)
出願日 平成19年3月19日(2007.3.19)
新規性喪失の例外の表示 特許法第30条第1項適用 信号処理研究会(SIP)(社団法人電子情報通信学会、2006年9月26日~9月27日開催) 応用音響研究会(EA)(社団法人電子情報通信学会、2006年11月23日~11月24日開催) 応用音響研究会(EA)(社団法人電子情報通信学会、2006年11月23日~11月24日開催) 2006 インターナショナル シンポジウム オン インテリジェント シグナル プロセッシング アンド コミュニケーション システムズ(ISPACS2006)(トリプルイー 2006年12月12日発行)
審査請求日 平成22年3月16日(2010.3.16)
特許権者または実用新案権者 【識別番号】803000115
【氏名又は名称】学校法人東京理科大学
発明者または考案者 【氏名】田邉 造
【氏名】古川 利博
個別代理人の代理人 【識別番号】100105050、【弁理士】、【氏名又は名称】鷲田 公一
審査官 【審査官】佐藤 聡史
参考文献・文献 国際公開第2005/015737(WO,A1)
調査した分野 H03H 15/00 -21/00
G10L 21/02
H04B 1/76 - 3/44
H04B 7/005- 7/015
特許請求の範囲 【請求項1】
所望の情報に雑音が混在した観測情報のみから前記所望情報を推定する雑音抑圧装置であって、
前記観測情報を取得する取得手段と、
有色雑音を駆動源として用いるカルマンフィルタを用いて、取得された前記観測情報から前記雑音を除去して前記所望情報を抽出する抽出手段と、を有し、
前記カルマンフィルタは、
状態空間モデルの状態方程式において自己回帰モデルの係数を使用しないように構成されている、
雑音抑制装置。
【請求項2】
前記抽出手段は、
時刻nのみの観測情報に対して、時刻nまでの情報により前記所望情報を含む時刻n+1のシステムの状態量を推定した場合の推定誤差の第1相関値行列を算出する第1の相関演算部と、
時刻nのみの観測情報に対して、前記第1の相関演算部によって算出された前記推定誤差の第1相関値行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の最適推定値ベクトルと、時刻nまでの情報による時刻n+1での前記状態量の最適推定値ベクトルと、前記観測情報を含む観測量の推定誤差ベクトルと、の関係を規定するための重み係数行列を算出する重み係数算出部と、
時刻nのみの観測情報に対して、時刻nまでの情報による時刻n+1での前記状態量の第1最適推定値ベクトルを算出する第1の最適推定値算出部と、
時刻nのみの観測情報に対して、前記重み係数算出部によって算出された前記重み係数行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の第2最適推定値ベクトルを算出する第2の最適推定値算出部と、
時刻nのみの観測情報に対して、時刻n+1までの情報により当該時刻の前記状態量を推定した場合の推定誤差の第2相関値行列を算出する第2の相関演算部と、
を有する請求項1記載の雑音抑圧装置。
【請求項3】
前記第1の相関演算部は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第2の相関演算部によって算出された前記推定誤差の第2相関値行列を用いて、前記推定誤差の第1相関値行列の算出を行い、
前記重み係数算出部は、
前記第1の相関演算部によって算出された前記推定誤差の第1相関値行列、与えられた観測遷移行列、および与えられた雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第1の最適推定値算出部は、
前記状態遷移行列、および、与えられたまたは前回前記第2の最適推定値算出部によって算出された前記状態量の第2最適推定値ベクトルを用いて、前記状態量の第1最適推定値ベクトルの算出を行い、
前記第2の最適推定値算出部は、
前記第1の最適推定値算出部によって算出された前記状態量の第1最適推定値ベクトル、前記重み係数算出部によって算出された前記重み係数行列、前記観測遷移行列、および時刻n+1のみにおける観測量を用いて、前記状態量の第2最適推定値ベクトルの算出を行い、
前記第2の相関演算部は、
前記重み係数算出部によって算出された前記重み係数行列、前記観測遷移行列、および前記第1の相関演算部によって算出された前記推定誤差の第1相関値行列を用いて、前記推定誤差の第2相関値行列の算出を行う、
請求項2記載の雑音抑圧装置。
【請求項4】
所望の情報に雑音が混在した観測情報のみから前記所望情報を推定する雑音抑圧方法であって、
前記観測情報を取得する取得ステップと、
有色雑音を駆動源として用いるカルマンフィルタを用いて、取得した前記観測情報から前記雑音を除去して前記所望情報を抽出する抽出ステップと、を有し、
前記カルマンフィルタは、
状態空間モデルの状態方程式において自己回帰モデルの係数を使用しないように構成されている、
雑音抑圧方法。
【請求項5】
前記抽出ステップは、
時刻nのみの観測情報に対して、時刻nまでの情報により前記所望情報を含む時刻n+1のシステムの状態量を推定した場合の推定誤差の第1相関値行列を算出する第1の相関演算工程と、
時刻nのみの観測情報に対して、前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の最適推定値ベクトルと、時刻nまでの情報による時刻n+1での前記状態量の最適推定値ベクトルと、前記観測情報を含む観測量の推定誤差ベクトルと、の関係を規定するための重み係数行列を算出する重み係数算出工程と、
時刻nのみの観測情報に対して、時刻nまでの情報による時刻n+1での前記状態量の第1最適推定値ベクトルを算出する第1の最適推定値算出工程と、
時刻nのみの観測情報に対して、前記重み係数算出工程で算出した前記重み係数行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の第2最適推定値ベクトルを算出する第2の最適推定値算出工程と、
時刻nのみの観測情報に対して、時刻n+1までの情報により当該時刻の前記状態量を推定した場合の推定誤差の第2相関値行列を算出する第2の相関演算工程と、
を有する請求項4記載の雑音抑圧方法。
【請求項6】
前記第1の相関演算工程は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第2の相関演算工程で算出した前記推定誤差の第2相関値行列を用いて、前記推定誤差の第1相関値行列の算出を行い、
前記重み係数算出工程は、
前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列、与えられた観測遷移行列、および与えられた雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第1の最適推定値算出工程は、
前記状態遷移行列、および、与えられたまたは前回前記第2の最適推定値算出工程で算出した前記状態量の第2最適推定値ベクトルを用いて、前記状態量の第1最適推定値ベクトルの算出を行い、
前記第2の最適推定値算出工程は、
前記第1の最適推定値算出工程で算出した前記状態量の第1最適推定値ベクトル、前記重み係数算出工程で算出した前記重み係数行列、前記観測遷移行列、および時刻n+1のみにおける観測量を用いて、前記状態量の第2最適推定値ベクトルの算出を行い、
前記第2の相関演算工程は、
前記重み係数算出工程で算出した前記重み係数行列、前記観測遷移行列、および前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列を用いて、前記推定誤差の第2相関値行列の算出を行う、
請求項5記載の雑音抑圧方法。
【請求項7】
所望の情報に雑音が混在した観測情報のみから前記所望情報を推定するための雑音抑圧プログラムであって、
コンピュータに、
前記観測情報を取得する取得ステップと、
有色雑音を駆動源として用いるカルマンフィルタを用いて、取得した前記観測情報から前記雑音を除去して前記所望情報を抽出する抽出ステップと、前記カルマンフィルタは、状態空間モデルの状態方程式において自己回帰モデルの係数を使用しないように構成されている、
を実行させるための雑音抑制プログラム。
【請求項8】
前記抽出ステップは、
時刻nのみの観測情報に対して、時刻nまでの情報により前記所望情報を含む時刻n+1のシステムの状態量を推定した場合の推定誤差の第1相関値行列を算出する第1の相関演算工程と、
時刻nのみの観測情報に対して、前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の最適推定値ベクトルと、時刻nまでの情報による時刻n+1での前記状態量の最適推定値ベクトルと、前記観測情報を含む観測量の推定誤差ベクトルと、の関係を規定するための重み係数行列を算出する重み係数算出工程と、
時刻nのみの観測情報に対して、時刻nまでの情報による時刻n+1での前記状態量の第1最適推定値ベクトルを算出する第1の最適推定値算出工程と、
時刻nのみの観測情報に対して、前記重み係数算出工程で算出した前記重み係数行列を用いて、時刻n+1までの情報による当該時刻での前記状態量の第2最適推定値ベクトルを算出する第2の最適推定値算出工程と、
時刻nのみの観測情報に対して、時刻n+1までの情報により当該時刻の前記状態量を推定した場合の推定誤差の第2相関値行列を算出する第2の相関演算工程と、
を有する請求項7記載の雑音抑圧プログラム。
【請求項9】
前記第1の相関演算工程は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第2の相関演算工程で算出した前記推定誤差の第2相関値行列を用いて、前記推定誤差の第1相関値行列の算出を行い、
前記重み係数算出工程は、
前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列、与えられた観測遷移行列、および与えられた雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第1の最適推定値算出工程は、
前記状態遷移行列、および、与えられたまたは前回前記第2の最適推定値算出工程で算出した前記状態量の第2最適推定値ベクトルを用いて、前記状態量の第1最適推定値ベクトルの算出を行い、
前記第2の最適推定値算出工程は、
前記第1の最適推定値算出工程で算出した前記状態量の第1最適推定値ベクトル、前記重み係数算出工程で算出した前記重み係数行列、前記観測遷移行列、および時刻n+1のみにおける観測量を用いて、前記状態量の第2最適推定値ベクトルの算出を行い、
前記第2の相関演算工程は、
前記重み係数算出工程で算出した前記重み係数行列、前記観測遷移行列、および前記第1の相関演算工程で算出した前記推定誤差の第1相関値行列を用いて、前記推定誤差の第2相関値行列の算出を行う、
請求項8記載の雑音抑圧プログラム。
発明の詳細な説明 【技術分野】
【0001】
本発明は、カルマンフィルタに基づく雑音抑圧装置および雑音抑圧方法に関する。
【背景技術】
【0002】
所望の情報に不必要な情報が混在した受信情報(付加雑音などにより破損した情報)から不必要な情報を取り除き、所望情報のみを抽出することは、音声や無線通信、画像、姿勢制御などの分野における重要な技術であり、近年盛んに研究開発が行われている。
【0003】
例えば、音声分野における公知の雑音抑圧方法としては、単一のマイクロホンを用いた方法や、複数のマイクロホンから構成されるマイクロホンアレイを用いた方法が提案されている。
【0004】
しかしながら、マイクロホンアレイを用いた方法では、雑音信号の数が増大すると、マイクロホンの数も比例して増加することが避けらず、コストが増大する。そのため、現在は、単一のマイクロホンを用いた雑音抑圧方法の開発が主流となっている。
【0005】
単一のマイクロホンしか用いない従来の雑音抑圧方法のアルゴリズムとしては、以下のようなものが知られている。
【0006】
非特許文献1記載のANC(適応ノイズキャンセラ)アルゴリズムは、音声信号の周期性を利用してノイズ信号を低減する。
【0007】
非特許文献2には、線形予測に基づいた雑音抑圧アルゴリズムが記載されている。このアルゴリズムは、ピッチ推定や、雑音パワースペクトラム、雑音の平均方向に関する事前知識を必要としない。
【0008】
また、上記アルゴリズムとは別に、カルマンフィルタに基づいた雑音抑圧アルゴリズムが、非特許文献3に提案されている。このアルゴリズムは、音声信号を自己回帰(AR:AutoRegressive)過程でモデル化する。さらに、このアルゴリズムは、自己回帰モデルの係数(以下「AR係数」という)を推定し、推定したAR係数を用いたカルマンフィルタに基づいて雑音抑圧を実行する。このアルゴリズムは、AR係数を推定する必要があるため、上記他のアルゴリズムとは本質的に異なっている。
【0009】
カルマンフィルタに基づいたアルゴリズムの多くは、通常、2段階で動作する。すなわち、このようなアルゴリズムは、最初にAR係数を推定し、次に推定したAR係数を用いたカルマンフィルタに基づいて雑音抑圧を行う。

【非特許文献1】J.R. Deller, J.G. Proakis, J.H.L. Hansen, "Discrete-Time Processing of Speech Signals," Macmillan Press, 1993
【非特許文献2】A. Kawamura, K. Fujii, Y. Itoh and Y. Fukui, “A Noise Reduction Method Based on Linear Prediction Analysis,” IEICE Trans. Fundamentals, vol.J85-A, no.4, pp.415-423, May 2002
【非特許文献3】W. Kim and H. Ko, "Noise Variance Estimation for Kalman Filtering of Noise Speech," IEICE Trans. Inf. & syst., vol.E84-D, no.1, pp.155-160, Jan 2001
【発明の開示】
【発明が解決しようとする課題】
【0010】
しかしながら、非特許文献1に記載された公知のアルゴリズムは、音声信号のピッチ周期の正確な推定を必要とする。そのため、このアルゴリズムは、その雑音抑圧能力が付加雑音によって劣化してしまうという問題点を有している。
【0011】
この点、非特許文献2記載のアルゴリズムは、音声信号のピッチ周期の正確な推定を必要とせずに、雑音抑圧が可能である。さらに、このアルゴリズムは、その原理が単純であり、演算量を少なくすることができるといった長所を有している。しかし、このアルゴリズムは、その雑音抑圧能力が入力音声信号の周期性や線形性などの特性に依存している。言い換えると、このアルゴリズムは、雑音抑圧能力の中に音声信号に依存するパラメータが存在しているため、その実用には一定の限界がある。
【0012】
非特許文献3記載のアルゴリズムは、強力な雑音抑圧能力を有し、特に高い音質を得たい音響分野への応用に適した手法である。
【0013】
しかしながら、一方で、このアルゴリズムは、AR係数を必要とするため、AR係数の推定精度に雑音抑圧性能(つまり、当該カルマンフィルタアルゴリズムの性能)が大きく依存してしまうという問題点を有している。すなわち、AR係数が正確に推定されない場合、雑音を除去し切れないのみならず、場合によっては雑音に加えて音声信号まで除去してしまう可能性がある。これらは、雑音除去された音声信号の音質の劣化を引き起こす要因となりうる。
【0014】
この点、一般には、AR係数の正確な推定は困難である。AR係数の正確な推定は、例えば、雑音除去であれば、クリアな音声信号に依存しているからである。このことは、音声信号が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。また、仮に何らかの手法でリアルタイムにAR係数を正確に推定することが可能となったとしても、処理が増加するため演算量の問題は避けられない。
【0015】
本発明は、かかる点に鑑みてなされたものであり、AR係数の推定を必要としないシンプルで雑音抑圧能力が高い、カルマンフィルタに基づく雑音抑圧装置および雑音抑圧方法を提供することを目的とする。
【課題を解決するための手段】
【0016】
本発明の雑音抑圧装置は、所望の情報に不必要な情報が混在した情報を取得する取得手段と、カルマンフィルタのみを用いて、取得された情報から前記不必要な情報を除去して前記所望情報を抽出する抽出手段と、を有し、前記カルマンフィルタは、状態空間モデルの観測方程式において自己回帰モデルの係数を使用しないように構成されている、構成を採る。
【0017】
本発明の雑音抑圧方法は、所望の情報に不必要な情報が混在した情報を取得する取得ステップと、カルマンフィルタのみを用いて、取得した情報から前記不必要な情報を除去して前記所望情報を抽出する抽出ステップと、を有し、前記カルマンフィルタは、状態空間モデルの観測方程式において自己回帰モデルの係数を使用しないように構成されている、ようにした。
【発明の効果】
【0018】
本発明によれば、カルマンフィルタを用いつつ、AR係数の推定を必要とすることなく、シンプルな構成で、雑音抑圧能力を向上することができる。
【発明を実施するための最良の形態】
【0019】
以下、本発明の実施の形態について、図面を参照して詳細に説明する。
【0020】
(実施の形態1)
図1は、本発明の実施の形態1に係る雑音抑圧装置の構成を示すブロック図である。
【0021】
図1に示す雑音抑圧装置100は、入力部110、サンプリング部120、A/D変換部130、バッファ140、雑音抑圧処理部150、および出力部160を有する。
【0022】
入力部110は、観測信号を入力する。観測信号は、信号源からのクリアな信号と、付加雑音信号とが合わさった信号である。入力部110は、例えば、入力したアナログの観測信号を入力処理して、サンプリング部120に出力する。入力処理は、例えば、帯域制限処理や自動利得制御処理などである。
【0023】
サンプリング部120は、所定のサンプリング周波数(例えば、16kHz)で、入力されたアナログの観測信号をサンプリング処理し、A/D変換部130に出力する。A/D変換部130は、サンプリングされた観測信号の振幅値を所定の分解能(例えば、8bit)でA/D変換処理し、バッファ140に送る。バッファ140は、所定のサンプリング数Nの信号フレーム(ブロック)を雑音抑圧処理部150に出力する。
【0024】
雑音抑圧処理部140は、本発明の特徴的な構成要素であり、AR係数を用いないカルマンフィルタを内蔵している。すなわち、カルマンフィルタに基づく従来の雑音抑圧方法(非特許文献3参照)(以下「従来手法」という)では、線形予測を用いてAR係数を推定した後、その結果を用いてカルマンフィルタを実行することで雑音抑圧を実現しているのに対し、本発明の雑音抑圧方法(以下「本手法」という)では、カルマンフィルタのみを用いて雑音抑圧を実現している。そのため、本手法では、信号源からのクリアな信号のみを用いて状態方程式を構成し、そのクリアな信号と付加雑音信号を用いて観測方程式を構成している。以下、雑音抑圧処理部140で行われるカルマンフィルタに基づく雑音抑圧処理動作について詳細に説明する。
【0025】
雑音抑圧処理部140に入力される観測信号r(n)は、信号源からのクリアな信号(例えば、音声信号など)d(n)以外に付加雑音信号v(n)を含んでおり、次の式(1)を満たす。
【数1】
JP0005115952B2_000002t.gif
ここで、nとは、装置の時刻nである。時刻nは、サンプリング部120で生成された離散的な時間系列において、処理開始時刻を時刻0と仮定したときに、そこからn番目の時刻のことを意味する。
【0026】
従来のAR過程に基づくモデル化方法(従来手法)では、信号d(n)は、AR係数を用いて次の式(2)でモデル化されると仮定している。
【数2】
JP0005115952B2_000003t.gif
ここで、α(n)は時刻nでのAR係数、KはAR係数の次数、e(n)は信号d(n)が、式(2)に示すK次のAR過程でモデル化されるとした場合の予測誤差(モデルリング誤差)である。
【0027】
公知のように、従来手法では、付加雑音信号v(n)は、零平均であり、白色雑音であることが前提条件である。言い換えると、従来手法では、信号d(n)と付加雑音信号v(n)は無相関であり、つまり、次の式(3)を満たす。
【数3】
JP0005115952B2_000004t.gif
ここで、E[・]はアンサンブル平均を示す。
【0028】
AR係数は、AR係数推定アルゴリズムにより推定される。従来手法において最も重要な点は、カルマンフィルタを用いた高性能の雑音抑圧を達成するために、AR係数の正確な推定を必要とすることである。このことからも、カルマンフィルタの雑音抑圧能力がAR係数の推定精度に大きく依存しているため雑音抑圧能力が大きく劣化することは容易に想像可能である。
【0029】
本手法では、AR係数を用いずにカルマンフィルタの状態空間モデル、つまり音源信号からのクリアな信号のみを用いて状態方程式、およびそのクリアな信号と付加雑音信号とを用いて観測方程式を構成している。
【0030】
以下、本手法での状態空間モデルの構成法について説明する。表記を容易にするために、まずK×1次の信号ベクトルx(n)を次の式(4)で定義する。添え字“p”は本発明により考案された表現であることを示す。
【数4】
JP0005115952B2_000005t.gif

【0031】
次に、式(2)と同様にして信号d(n)のモデルを構成する。ここでは、AR係数を用いないとしても、従来手法と同様に、未知のK×K次の状態遷移行列Φ(n+1)と、K×1次の駆動雑音ベクトルδ(n+1)を導入する。これにより、本実施の形態の状態空間モデル(信号ベクトルにより記述される)の状態方程式の形として、次の式(5)が定まる。
【数5】
JP0005115952B2_000006t.gif

【0032】
式(5)が等式としてAR係数を用いずに成立するためには、K×K次の状態遷移行列Φ(n+1)、K×1次の駆動雑音ベクトルδ(n+1)は、次の式(6)および式(7)を満たすことが求められる。
【数6】
JP0005115952B2_000007t.gif
【数7】
JP0005115952B2_000008t.gif

【0033】
したがって、式(5)は次の式(8)に書き直せる。
【数8】
JP0005115952B2_000009t.gif

【0034】
式(6)で表されるK×K次の状態遷移行列Φ(n+1)の行列要素は0と1のみであり、特に第1行がすべて0である。このことは、本手法の状態空間モデルの特徴の一つである。また、上記から明らかなように、式(5)、式(6)、式(7)の導出において、駆動雑音ベクトルδ(n+1)に対して前提条件を与えていないことに注意すべきである。すなわち、駆動雑音ベクトルδ(n+1)は有色であってよい。この理由については後述する。
【0035】
本手法のカルマンフィルタアルゴリズムでは、特異値分解を避けるためにK×1次の観測ベクトルy(n+1)を次の式(9)で定義する。
【数9】
JP0005115952B2_000010t.gif

【0036】
よって、式(1)および式(9)から、本手法の状態空間モデルの観測方程式として、次の式(10)が導出される。
【数10】
JP0005115952B2_000011t.gif
ここで、MはK×K次の状態遷移行列、εはK×1次の付加雑音ベクトルであり、それぞれ次の式(11)および式(12)で定義される。
【数11】
JP0005115952B2_000012t.gif
【数12】
JP0005115952B2_000013t.gif

【0037】
以上より、本手法のカルマンフィルタアルゴリズムで用いる状態空間モデルの状態方程式および観測方程式として、式(5)および式(10)がそれぞれ求められた。
【0038】
図2は、本実施の形態に係る雑音抑圧装置の状態空間モデルを説明するブロック線図である。
【0039】
図2において、500は時刻nにおける信号ベクトルx(n)、501は時刻n+1における信号ベクトルx(n+1)、502は時刻nにおける観測ベクトルy(n)、503は時刻nにおける付加雑音ベクトルε(n)、504は時刻n+1における駆動雑音ベクトルδ(n+1)、505は状態遷移行列Φ、506は状態遷移行列Mである。
【0040】
図2から理解できるように、本手法のカルマンフィルタアルゴリズムは、駆動雑音が有色であるにもかかわらず、従来のカルマンフィルタアルゴリズムと同様の手順で実行できる(理由は後述する)。本手法のカルマンフィルタアルゴリズムに基づいた信号推定は、後述するフローチャートに従って行われる。
【0041】
図3は、本実施の形態に係る雑音抑圧装置のカルマンフィルタアルゴリズムに基づく信号推定手順を説明するフローチャートである。
【0042】
まず、推定信号ベクトルの初期値x(0|0)、推定誤差の共分散行列の初期値P(0|0)およびスカラー量である付加雑音ベクトルの共分散Rε(n)[i,j]の値を、次の式(13)に示すように設定する(ST100)。
【数13】
JP0005115952B2_000014t.gif

【0043】
ここで、Iは単位行列である。また、σは付加雑音信号v(n)の雑音分散であり、既知と仮定している。もし付加雑音信号v(n)が白色雑音であり零平均であれば、σは、次の式(14)で与えられる。
【数14】
JP0005115952B2_000015t.gif
ここで、Nは所定のサンプル数である。
【0044】
次に、カルマンフィルタアルゴリズムを実行する(ST101)。このカルマンフィルタアルゴリズムの処理手順は、後で詳細に説明する。
【0045】
次に、時刻nが所定のサンプル数Nに達したか否かを判定し(ST102)、所定のサンプル数Nに達していない場合は(ST102:NO)、ステップST101に戻って、カルマンフィルタアルゴリズムを繰り返し、所定のサンプル数Nに達した場合は(ST102:YES)、上記一連の処理を終了する。
【0046】
本手法では、AR係数を用いずに状態空間モデルを設定している。そのため、従来手法で必要であったAR係数を推定するステップを削減することができ、1段階処理で信号推定が可能となる。このことは本発明の大きな特徴の一つである。
【0047】
図4は、図3のカルマンフィルタアルゴリズムの処理内容を示すフローチャートである。
【0048】
ここで、開始時の時刻はn+1であり、以下のフローでは、時刻n+1での観測信号r(n+1)および観測ベクトルy(n+1)と、時刻nでの推定誤差の共分散行列P(n|n)とから、時刻n+1での推定信号ベクトルxpe(n+1|n)を推定し、同時に推定誤差の共分散行列P(n+1|n+1)を更新する。
【0049】
まず、次の式(15)を用いて、駆動雑音ベクトルの共分散Rδ(n)[i,j]の値を計算する(ST200)。
【数15】
JP0005115952B2_000016t.gif

【0050】
次に、推定誤差の共分散行列P(n+1|n)を計算し、更新する(ST201)。ここで、推定誤差の共分散行列P(n+1|n)の更新式は、次の式(16)である。ここでは、以下、添え字“pe”および“e”は、推定値であることを示す。
【数16】
JP0005115952B2_000017t.gif

【0051】
この式(16)を展開すると、次の式(17)となる。
【数17】
JP0005115952B2_000018t.gif

【0052】
ここで、式(16)の第3項および第4項に注目する。第3項および第4項は、推定信号ベクトルxpe(n|n)と駆動雑音ベクトルδ(n|n)のアンサンブル平均を含んでいる、したがって、推定信号ベクトルxpe(n|n)と駆動雑音ベクトルδ(n|n)とが無相関でない、言い換えれば、駆動雑音ベクトルδ(n|n)が有色のとき、第3項および第4項は0にならない。すなわち、本手法のカルマンフィルタアルゴリズムが、駆動雑音信号が有色信号の場合に、従来と同様なカルマンフィルタアルゴリズム手順で適用可能か否かは、第3項および第4項が、式(15)で表される推定誤差の共分散行列P(n+1|n)の更新にどの程度影響するかによる。
【0053】
そこで、式(17)の第3項に着目する。第3項は、次の式(18)のように、
【数18】
JP0005115952B2_000019t.gif
と書き表せる。
【0054】
ここで、次の式(19)のように、
【数19】
JP0005115952B2_000020t.gif
とすれば、次の式(20)となる。
【数20】
JP0005115952B2_000021t.gif

【0055】
式(5)および次の式(21)
【数21】
JP0005115952B2_000022t.gif
を考慮すれば、式(18)は、次の式(22)のように書き直すことができる。
【数22】
JP0005115952B2_000023t.gif
ただし、次の式(23)である。
【数23】
JP0005115952B2_000024t.gif
ここで、d(n)は時刻nでの推定信号である。
【0056】
式(23)の第1項は真値のアンサンブル平均、第2項は推定値のアンサンブル平均を示す。したがって、その差は小さく、e(l+1)は無視できる値である。よって、推定誤差の共分散行列P(n+1|n)の更新式(15)は、次の式(24)のように書き直される。このことから、駆動雑音が有色であっても、本手法の式(5)および式(9)から構成される状態空間モデル(状態方程式と観測方程式)であれば、従来手法と同様なカルマンフィルタアルゴリズム手順が適用可能となる。
【数24】
JP0005115952B2_000025t.gif

【0057】
次に、カルマンゲインK(n+1)を計算する(ST202)。ここで、カルマンゲインK(n+1)は、次の式(25)で与えられる。
【数25】
JP0005115952B2_000026t.gif

【0058】
そして、これらの値から、まず推定信号ベクトルxpe(n+1|n)を計算する(ST203)。ここで、推定信号ベクトルxpe(n+1|n)は、次の式(26)で与えられる。
【数26】
JP0005115952B2_000027t.gif

【0059】
そして、推定信号ベクトルxpe(n+1|n+1)を計算する(ST204)。ここで、推定信号ベクトルxpe(n+1|n+1)は、次の式(27)で与えられる。
【数27】
JP0005115952B2_000028t.gif

【0060】
最後に、推定誤差の共分散行列P(n+1|n+1)を計算、更新し(ST205)、フローを終了する。ここで、推定誤差の共分散行列P(n+1|n+1)は、次の式(28)で与えられる。
【数28】
JP0005115952B2_000029t.gif

【0061】
雑音抑圧処理部140は、上述した手順、アルゴリズムによって推定された推定信号をを出力部150に出力する。出力部150は、例えば、スピーカやディスプレイ、通信手段、記憶装置などで構成されている。
【0062】
AR係数の推定を必要とする従来のカルマンフィルタの問題点は、AR係数の推定精度にカルマンフィルタアルゴリズムの能力が依存していることである。これに対して、本手法は、カルマンフィルタアルゴリズムのみで実行しているため、AR係数の推定精度に依存しない。
【0063】
なお、本手法のカルマンアルゴリズムは、以下に説明する更に改良された方法(以下「改良手法」という)によっても実行できる。以下、この改良手法について説明する。
【0064】
改良手法のカルマンアルゴリズムの状態方程式は、上記式(5)と同一である。ただし、観測方程式を、次の式(29)のように書き直す。改良手法の観測方程式(29)は、本手法の観測方程式(10)の一部の要素を取り出してきたものである。
【数29】
JP0005115952B2_000030t.gif
ここで、y’(n+1)、m、ε’(n+1)は、それぞれ改良手法における観測信号、雑音信号、状態遷移ベクトルであり、次の式(30)を満たす。
【数30】
JP0005115952B2_000031t.gif

【0065】
改良手法の大きな特徴の一つは、本手法の観測方程式(10)と改良手法の観測方程式(30)の比較からわかるように、改良前の方法(本手法)における、観測ベクトルy(n+1)、雑音ベクトルε(n+1)、および状態遷移行列mが、それぞれスカラーとベクトルに変更されていることである。したがって、改良手法はより少ない演算量で実行できるであろうことは明らかである。また、上記から明らかなように、改良手法においても、本手法と同様に、駆動雑音ベクトルδ(n+1)に対して前提条件を与えていないことに注意すべきである。すなわち、改良手法においても駆動雑音ベクトルδ(n+1)は有色であってよい。
【0066】
すなわち、式(17)で与えられた推定誤差の共分散行列P(n+1|n)を再度、行列形式で書き下すと、次の式(31)となる。
【数31】
JP0005115952B2_000032t.gif
ここで、また以下の式で、行列式中、斜線網掛け部分で示される部分は駆動雑音が有色であるときに影響を受ける部分である。
【0067】
改良手法を用いたカルマンフィルタアルゴリズムは、図4に示す本手法のカルマンフィルタアルゴリズムとは、カルマンゲインK(n+1)の計算(ST202)および推定信号ベクトルxpe(n+1|n+1)の計算(ST204)における計算式が異なる。
【0068】
すなわち、カルマンゲインK(n+1)の計算式は、次の式(32)、つまり、
【数32】
JP0005115952B2_000033t.gif
となる。ただし、雑音信号ε’(n+1)の共分散をσとしたときに、次の式(33)が成り立つ。
【数33】
JP0005115952B2_000034t.gif

【0069】
また、推定信号ベクトルxpe(n+1|n+1)の計算式は、次の式(34)となる。
【数34】
JP0005115952B2_000035t.gif
となる
【0070】
ここで、式(34)のxpe(n+1|n+1)の第1要素に着目すると、斜線網掛けではない(つまり、有色雑音の影響を受けない)ことから、駆動雑音が有色であってもクリアな信号の推定信号を得られる。したがって、改良手法も、本手法と同様に、従来のカルマンフィルタアルゴリズムで実行可能である。既に述べたように、改良手法では、本手法に比べて、演算量を低減できる。しかしながら、本手法が観測ベクトルy(n+1)を用いているのに対して、改良手法方法では、スカラー量である観測信号y‘(n+1)を用いている。そのサイズの違いから、過去の観測量をより積極的に用いることができる本手法のほうが、より雑音抑圧能力が高いことは明らかである。
【0071】
以上をまとめると、本発明における上記二つの手法(本手法および改良手法)は、カルマンフィルタに必要な状態空間モデルの観測方程式を変更することによって、演算量を大幅に減少させることが可能である。より具体的には、本発明における上記二つの手法では、AR係数の推定を必要としないため、その演算量は従来手法に比べて少ない。また、改良手法では、演算に用いる変数のサイズが本手法に比べて小さいため、その演算量は本手法に比べて少ない。すなわち、演算量に関して、
改良手法<本手法<従来手法
である。
【0072】
また、本発明における上記二つの手法を、例えば、半導体集積回路や半導体回路などのハードウエアとして実施する場合、また、パーソナルコンピュータなどで実行可能なソフトウエアとして実施する場合のいずれにおいても、その構成は従来手法よりも単純化される。したがって、本発明における二つの手法を用いれば、回路規模やプログラム量が大幅に低減できるであろうことは明らかである。
【0073】
また、本発明における上記二つの手法は、カルマンフィルタに必要な状態空間モデルの観測方程式を変更することによって、雑音抑圧能力を大幅に向上させることが可能である。より具体的には、本発明における上記二つの手法では、AR係数の推定を必要としないため、雑音抑圧能力がAR係数の推定精度に依存する従来手法に比べて、雑音抑圧能力が高い。また、改良手法では、演算に用いる変数のサイズが本手法に比べて小さいため、過去の観測量を積極的に用いていない。したがって、その雑音抑圧能力は本手法に比べて低い。すなわち、雑音抑圧能力に関して、
従来手法<改良手法<本手法
である。
【0074】
例えば、音響分野において、音質を若干落としても(ただし、従来手法よりは高い)、演算速度を速くしたい場合に、改良手法は有効である。
【0075】
AR係数を用いない本発明における上記二つの手法は、状態方程式の駆動雑音が有色信号となっており、従来のカルマンフィルタの仮定、つまり、駆動雑音が白色信号であるという仮定に反している。しかしながら、本発明の状態空間モデルにおける状態方程式と観測方程式の構造的な性質により、駆動雑音が白色信号である従来のカルマンフィルタアルゴリズと全く同じアルゴリズムで実行可能であり、その理由は上記の通りである。また、本発明の改良手法においても、同様に実行可能であることは上記より明らかである。
【0076】
本発明における上記二つの方法は、AR係数を使わないようにカルマンフィルタを構成する。言い換えると、カルマンフィルタのみで実行可能な雑音抑圧方法であるため、従来の雑音抑圧能力の問題を解決していることは明らかである。
【0077】
(実施の形態2)
実施の形態2は、実施の形態1に示す雑音抑圧装置を音声に適用した場合である。
【0078】
図5は、本発明の実施の形態2に係る雑音抑圧装置の構成を示すブロック図である。
【0079】
図5に示す雑音抑圧装置200は、本実施の形態の雑音抑圧処理を実行できるパーソナルコンピュータ210、マイクロホン220、サンプリング部120、およびA/D変換部130を有する。
【0080】
パーソナルコンピュータ210は、操作装置211、ディスプレイ212、バスインタフェース213、記録装置214、主記憶メモリ215、および中央演算装置216を有する。操作装置211は、典型的にはキーボートやマウスなどであるが、音声認識装置などを用いてもよい。使用者は、操作装置211を用い、ディスプレイ212で確認をしながらコンピュータを操作できる。
【0081】
パーソナルコンピュータ210において本実施の形態の雑音抑圧処理を実行させるプログラムソフトウエアは、記録装置214に格納されていてもよいし、バスインタフェース213を介して外部からダウンロードされてきてもよい。記録装置214は、典型的にはハードディスク装置であるが、CD-ROM装置やDVD装置、フラッシュメモリなどの可搬性のあるものであってもよい。また、それらの組み合わせであってもよい。
【0082】
サンプリング部120およびA/D変換部130は、パーソナルコンピュータ210内部に格納された内蔵カード(ボード)であってもよいし、バスインタフェース213を経由して接続された外部設置型機器であってもよい。
【0083】
マイクロホン220からの観測音声信号は、サンプリング部120に入力される。サンプリング部120は、所定のサンプリング周波数(例えば、16kHz)で、入力されたアナログの観測音声信号をサンプリング処理し、A/D変換部130に出力する。A/D変換部130は、サンプリングされた観測音声信号の振幅値を所定の分解能(例えば、8bit)でA/D変換処理し、一時格納する。A/D変換部130は、所定のサンプリング数Nの音声フレーム単位で、ディジタル化した観測音声信号をパーソナルコンピュータ210のバスインタフェース213に出力する。
【0084】
パーソナルコンピュータ210はバスインタフェース213に出力された観測音声信号を一時、主記憶メモリ215に格納し、その後、所定の音声フレーム(サンプリング数)単位で、雑音抑圧処理を施した上で、再度主記憶メモリ215に格納する。雑音抑圧処理は、主記憶メモリ215や記録装置214に格納されたソフトウエアをバスインタフェース213経由で中央演算装置216に呼び出し、実行させることで行われる。
【0085】
パーソナルコンピュータ210は、使用者の操作により、処理を実行したり、中断、終了させたりする。また使用者の操作により、処理した推定音声信号を外部に出力してもよい。
【0086】
次に、本実施の形態における音声の雑音抑圧処理について、図面を参照しつつ説明する。
【0087】
本実施の形態における音声の雑音抑圧の性能評価の目的で音声波形の数値シミュレーションを行った。
【0088】
図6は、本実施の形態における音声波形シミュレーションの第1の例の結果を示す図であり、図7は、本実施の形態における音声波形シミュレーションの第2の例の結果を示す図である。また、図8は、本実施の形態における音声波形シミュレーションの第3の例の結果を示す図であり、図9は、本実施の形態における音声波形シミュレーションの第4の例の結果を示す図である。
【0089】
日本人成人男性および女性のオリジナル音声信号(クリアな音声信号)は、無響室において、16kHzのサンプリングレートでサンプリング、ディジタル化した。本例では、二つの音声信号サンプルを検討した。二つの音声信号サンプルは、(A-1)図6(A)に示す無声区間のないクリアな音声信号、(A-2)図7(A)に示す無声区間をもつクリアな音声信号である。それぞれの音声信号サンプルは、以後、音声(A-1)および音声(A-2)と参照する。
【0090】
雑音信号は、人工的な付加雑音信号であり、本例では、二つの雑音信号サンプルを検討した。二つの雑音信号サンプルは、(B-1)図6(B)および図7(B)に示す付加白色ガウス雑音、(B-2)図8(B)および図9(B)に示す付加有色バブル雑音である。それぞれの雑音信号サンプルは、以後、雑音(B-1)および雑音(B-2)と参照される。
【0091】
雑音信号の雑音分散σは既知として次の式(35)で表されるものとする。すなわち、
【数35】
JP0005115952B2_000036t.gif
である。ここで、Lは全音声信号サンプル数である。
【0092】
また、信号雑音比SNRinを式(36)で定義する。
【数36】
JP0005115952B2_000037t.gif

【0093】
従来手法と本発明の方法(本手法)とによる雑音抑圧の結果を、音声(A-1)および音声(A-2)と雑音(B-1)の組み合わせ、ならびに音声(A-1)および音声(A-2)と雑音(B-2)の組み合わせからなる観測音声信号の信号波形から比較した。
【0094】
図6(C)と図7(C)は、それぞれ音声(A-1)と雑音(B-1)、音声(A-2)と雑音(B-1)からなる観測音声信号波形である。
【0095】
図6(D)と図7(D)は、それぞれ音声(A-1)と雑音(B-1)、音声(A-2)と雑音(B-1)の合成波形に対して従来手法で雑音抑圧を行った後の推定音声信号波形である。
【0096】
一方、図6(E)と図7(E)は、それぞれ音声(A-1)と雑音(B-1)、音声(A-2)と雑音(B-1)の合成波形に対して従来手法で雑音抑圧を行った後の推定音声信号波形である。
【0097】
図6(D)と図7(D)に対して図6(A)と図7(A)に示すオリジナルの音声信号とを参照すると、図6(D)と図7(D)に対して図6(A)と図7(A)に示すクリアな音声信号とを比較すると、従来手法による雑音抑圧では、雑音抑圧後に推定音声信号の振幅が小さくなっており、クリアな音声信号が抑圧されていることがわかる。加えて、従来手法による雑音抑圧では、サンプリング数の増加とともに、雑音抑圧後の推定音声信号の波形が、図6(A)と図7(A)に示すクリアな音声信号の波形から変形していく。
【0098】
さらに、従来手法の雑音抑圧では、無声区間を有する音声(A-2)に対して、推定音声信号が抑圧されるだけでなく、無声区間においてオリジナルの雑音信号と異なる雑音が観察されている。これは、従来手法では、無声区間では信号d(n)がゼロであるにもかかわらず、式(2)でAR係数を求めようとするため、AR係数の値は発散し、不安定な状態を与えるからと推測される。
【0099】
また、このことから、雑音信号が有色の場合、従来手法の適用は困難であろうことは容易に推測される。
【0100】
従来手法とは対照的に、図6(E)および図7(E)に示す本発明の雑音抑圧では、雑音抑圧後の推定音声信号の波形が、図6(A)と図7(A)に示すクリアな音声信号の波形と非常に似通っている。
【0101】
次に、図8(D)と図9(D)に対して図8(A)と図9(A)に示すクリアな音声信号とを比較すると、従来手法による雑音抑圧では、雑音(B-2)を含む観測音声信号に対して、非常に劣った結果を与えていることがわかる。これは、従来手法では、有色雑音である雑音(B-2)を含んだ観測音声信号に対してAR係数を正確に推定することが困難であるためである。
【0102】
一方、本発明の雑音抑圧方法では、雑音(B-2)の場合も雑音(B-1)の場合と同程度の雑音抑圧が達成されている。
【0103】
すなわち、上述のように、本発明の雑音抑制方法は、白色、有色雑音、無声区間の有無に関わらず有効である。これは、本発明の雑音抑制方法の大きな特徴の一つである。
【0104】
本発明の雑音抑制方法が有色の付加雑音にも適用できることは、既述のカルマンアルゴリズムにおける推定誤差の共分散行列P(n+1|n)の更新式(23)の導出過程からも明らかであるが、上述の本実施の形態の数値シミュレーションは、このことを支持している。
【0105】
次に、雑音抑制性能を数値的に比較するため、雑音抑圧性能SNRoutを式(37)で定義し、その数値シミュレーションを行った。
【数37】
JP0005115952B2_000038t.gif
ここで、dei(n)は推定音声信号を表す。
【0106】
表1は、本実施の形態における雑音抑圧性能の数値シミュレーションの一例の結果を示す表であり、音声(A-1)と雑音(B-1)の組み合わせにおける、幾つかのSNRin、K(従来手法においては式(2)におけるAR係数の次数、本発明では状態遷移行列Φの次数)の値における、従来手法と本発明の方法の雑音抑圧性能SNRoutを比較して示している。
【0107】
【表1】
JP0005115952B2_000039t.gif

【0108】
表2は、本実施の形態における雑音抑圧性能の数値シミュレーションの他の例の結果を示す表であり、音声(A-1)と雑音(B-1)の組み合わせにおける、幾つかのSNRin、Kの値における、従来手法と本発明の方法の雑音抑圧性能SNRoutを比較して示している。
【0109】
【表2】
JP0005115952B2_000040t.gif

【0110】
表1および表2を参照すると、本発明の方法は、すべてのSNRin、Kの値において、従来の方法に比べて雑音抑圧能力を改善していることがわかる。
【0111】
特に、表2に示す有色雑音の場合には、従来手法は非常に劣った結果を与えているのにも関わらず、本発明の方法は表1に示す白色雑音の場合と同程度の結果を示している。すなわち、本発明の雑音抑圧方法は、白色雑音、有色雑音両者に効果的で雑音の性質に堅牢な雑音抑圧方法であるといえる。
【0112】
また、表1および表2に見られるように、本発明の方法では、Kの値に対して雑音抑圧性能SNRoutは安定であり、Kの値の増加に伴い増加する傾向にある。対照的に従来の方法では、表1に見られるように、Kの値に対して雑音抑圧性能SNRoutは不安定である。これは、従来手法では、最適なKの値、つまりAR係数の次数を決定することが困難であることを意味している。
【0113】
AR係数推定を必要とする従来手法において最も問題になることは、AR係数の次数の正確な推定は一般に困難であることである。なぜなら、AR係数の次数の正確な推定は,例えば雑音除去であれば、クリアな音声信号に依存しているからである。
【0114】
このことは、クリアな音声信号が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。AR係数の次数が正確でない場合には、カルマンフィルタアルゴリズの性能が劣化することは容易に想像可能である。また、何らかの手法でリアルタイムに推定することが可能となったとしても、処理が増加することより演算量などの問題を避けることは不可能である。
【0115】
次に、本実施の形態に係る推定音声信号の音声品質を評価するためにリスニングテストによる主観的評価を行った。
【0116】
音声品質評価に用いたオリジナル音声信号とオリジナル雑音信号は第2の実施の形態と同一であり、その説明は省略する。雑音信号は、異なるSNRin(0、5、および10[dB])で音声信号に加えた。
【0117】
音声品質評価は、ACR(絶対範疇評価)に基づいた5段階MOS(平均オピニオン値)を用いたリスニングテストにより行った。50人の聴取者が雑音抑圧により得られた推定音声信号のうち幾つかを評価した。各々の聴取者は、ポイント1からポイント5を決定する。ポイント5が最良である。
【0118】
図10は、本実施の形態における雑音抑圧後の音声品質の主観的評価結果の一つの例を示す図であり、音声(A-1)と雑音(B-1)の組み合わせにおける、従来手法と本発明の方法とリスニングテストの結果を比較して示している。
【0119】
また、図11は、本実施の形態における雑音抑圧後の音声品質の主観的評価結果の他の例を示す図であり、音声(A-1)と雑音(B-2)の組み合わせにおける、従来手法と本発明の方法とリスニングテストの結果を比較して示している。
【0120】
図10および図11から、本発明の方法で推定した音声信号のスコアは、すべてのSNRin値において従来手法のスコアより高いことがわかる。特にその差は、音声(A-1)と雑音(B-2)の組み合わせに対して大きい。
【0121】
以上より、本発明の雑音抑圧方法は、音声信号の音声品質を犠牲にすることのない、白色雑音、有色雑音に効果的な優れた雑音抑圧方法であるといえる。
【0122】
(実施の形態3)
本発明に係る雑音抑圧装置は、雑音抑圧以外の用途にも応用可能である。以下、その一つの応用例について説明する。
【0123】
図12は、本発明の雑音抑圧方法、つまりAR係数を必要としないカルマンフィルタが適用されたマルチキャリア受信装置の構成を示すブロック図である。
【0124】
図12に示すマルチキャリア受信装置300は、主に、検波部301、GI(ガードインターバル)除去部302、チャネル等化部303、チャネル推定部304、FFT部305、復調部306、および復号部307を有する。
【0125】
本実施の形態において、チャネル推定部304は、本発明の状態空間モデルを用いたカルマンフィルタに基づいたチャネル推定を実行できるように構成されている。より具体的には、チャネル推定部304は、サンプリング部120、A/D変換部130、バッファ140、雑音抑圧処理部150、出力部160を有する構成をとる。
【0126】
検波部301は伝送路上で周波数選択性フェージング等の影響を受けた受信信号を、中間周波数で直交検波し、GI除去部402に出力する。GI除去部302は、受信信号のガードインターバルを除去し、シンボル単位に連なった信号をチャネル等化部303と、チャネル推定部304に出力する。
【0127】
チャネル推定部304に入力された信号は、入力部110に出力される。入力部は信号に所定の入力信号処理を施し、サンプリング部120に出力する。サンプリング部120は、所定のサンプリング周波数(例えば16kHz)で、入力されたアナログの受信信号をサンプリング処理し、A/D変換部130に出力する。A/D変換部130は、サンプリングされた受信信号の振幅値を所定の分解能(例えば8bit)でA/D変換処理し、バッファ140に送る。バッファ140は所定のサンプリング数Nの受信信号フレーム(ブロック)を雑音抑圧処理部150に出力する。
【0128】
雑音抑圧処理部150は、カルマンフィルタに基づいたチャネル推定をサブキャリア全体、または各サブキャリアついて行い、推定された(サブ)チャネルゲインを出力部160に出力する。出力部160は、入力された推定チャネルゲインをチャネル等化部302に出力する。雑音処理部150での処理の詳細については後述する。
【0129】
チャネル等化部303は、入力された推定チャネルゲインを用いて、GI除去部302から入力した信号の同期検波を行い、結果をFFT部305に出力する。FFT部305は、フーリエ変換処理を行い、受信信号を各サブキャリア信号成分に分離して、復調部306に出力する。復調部306は信号の復調処理を行い、結果を復号部307に出力する。復号部307は、信号の復号を行い、デジタルデータを出力する。
【0130】
本実施の形態では、まず、フェージングチャンネルを既述の状態空間モデルで表現する。マルチキャリア通信において、第kサブキャリアのみに着目すると、その受信信号y(t)は、次の式(38)
【数38】
JP0005115952B2_000041t.gif
となる。ここでnは時刻、Nは総キャリア数、S(t)、H(t)は第kキャリアの送信信号、チャネルゲインである。
【0131】
表記を容易にするために、N×1次チャネルゲインベクトルhp(n)を次の式(39)で定義する。
【数39】
JP0005115952B2_000042t.gif
ここで、Kは後述の状態遷移行列の次数である。
【0132】
本実施の形態では、本発明の状態空間モデルを用いてカルマンフィルタに基づき、チャネルゲインの推定値を適応的に求める。状態空間モデルの状態ベクトルx(n+1)を次の式(40)で定義する。
【数40】
JP0005115952B2_000043t.gif

【0133】
伝送路特性の時間変動を表す状態方程式は、次の式(41)で記述される。
【数41】
JP0005115952B2_000044t.gif
また、観測方程式は、次の式(42)で記述される。
【数42】
JP0005115952B2_000045t.gif
ここで、δ(n+1)は駆動雑音ベクトル、ε(n+1)は雑音信号ベクトル、C、Dは状態遷移行列、Gは定数行列である。
【0134】
駆動雑音ベクトルδ(n+1)、雑音信号ベクトルε(n+1)、状態遷移行列C、D、定数行列Gは以下の式(43)から式(48)で定義される。
【数43】
JP0005115952B2_000046t.gif
【数44】
JP0005115952B2_000047t.gif
【数45】
JP0005115952B2_000048t.gif
【数46】
JP0005115952B2_000049t.gif
【数47】
JP0005115952B2_000050t.gif
【数48】
JP0005115952B2_000051t.gif
ここで、IN×N、0N×NはN×Nの単位行列、零行列である。
【0135】
上記の定式化は、従来Vector Kermanと呼ばれているチャネルゲイン推定法に、本発明の状態空間モデルを適用したものである。フェージングチャネル間の相関だけでなく、サブチャネル間の相関も利用しているため良い推定精度を有する。しかしながら、例えば、式(38)の状態ベクトルx(n+1)のサイズはK×Nであり、計算量が増大するといった問題がある。
【0136】
この問題を改良したものが、従来Per Subcarrier Kalmanと呼ばれている手法である。この方法は、上記方法をサブキャリア単位に分割して処理を行うものである。本発明の状態空間モデルを適用した場合、上記、式(40)から式(48)をサブキャリアk単位に分割することで所望の式が得られる。
【0137】
サブキャリア単位への分割は、状態遷移行列Cを例にとると、式(42)において、N×Nの単位行列、零行列であるIN×N、0N×Nを、スカラー量である1,0に置き換えることに相当する。したがって、行列のサイズは1/Nになり、計算量を低減することができる。無線伝送の場合には、音響分野とは異なり、推定される信号の品質によりも処理速度が要求される。そのため、本方法はより実用的である。
【0138】
次に、本実施の形態におけるチャネル推定精度の比較結果の一例について、図面を参照しつつ説明する。
【0139】
チャネル推定精度の比較検討は数値シミュレーションにより行った。数値シミュレーションの条件は次のように設定した。送信フレームは64のシンボルからなり、そのうち4シンボルがパイロットシンボル、60シンボルがデータシンボルである。総送信フレーム数は200、総送信データシンボル数は12000である.評価量はNMSE(Normarized Mean Square Error)を用いた。
【0140】
図13は、本実施の形態におけるチャネル推定精度数値シミュレーションの一つの例の結果を示す図であり、fT=0.045、(NT)=0.08μsにおける、信号雑音比SNRに対するNMSE特性を示している。ここでfは最大ドップラー周波数、TはOFDMシンボル周期、Tはサンプリング間隔、τmaxは最大遅延スプレッドである。図13から、本発明の方法の推定精度が従来手法に比べ向上していることがわか。これは、本発明の方法では、AR係数の推定誤差によるチャネルの推定精度の劣化が軽減できたためと考えられる。また、fが180Hzより大きい場合、本発明の方法の性能向上がより大きなものとなることは明らかである。
【0141】
図14は、本実施の形態におけるチャネル推定精度数値シミュレーションの他の例の結果を示す図であり、SNR=20dB、τmax/(NT)=0.08μsにおける、最大ドップラー周波数fに対するNMSE特性を示している。この結果より、本発明の方法は、NMSEが最大ドップラー周波数fに依存しないため、従来手法より良好な性能を有していることがわかる。すなわち、チャネル間干渉の影響が深刻となるフェージング変動の激しい環境においても、本発明の方法の有効性が確認できる。
【0142】
本発明は、上記各実施の形態に限定されるものではない。
【0143】
本発明の雑音抑圧装置は、ノイズが含まれた音声信号(得られた情報)からクリアな音声信号(必要な情報)を取り出すことが可能である。その一つの実施の形態として、カーナビゲーション装置などに必要不可欠な音声認識装置の前処理雑音除去装置が考えられる。
【0144】
また、画像分野においては、何らかの原因でぼけてしまったぼけ画像(得られた情報)からぼけのとれたクリアな画像(必要な情報)を取り出すことが可能であり、画像処理装置として活用可能である。
【0145】
さらに、従来、AR過程によるモデル化とカルマンフィルタアルゴリズムとを組み合わせを用いた通信・信号処理全般にわたり、本発明が適応可能であることはいうまでもない。
【0146】
また、医療分野では、従来、妊婦の胎児の状況を検査するには、個人が購入することができない高価な装置と高い専門知識とが必要があったが、本発明によれば、妊婦の体から母胎の心音や胎児の心音やその他の音(得られた情報)から不必要な音(情報)を取り除き胎児の心音(必要な情報)のみを取り出すことが可能になり、通院せずとも自宅で胎児の健康状態を、その心音から容易に確認することが可能となる。
【0147】
上記各実施の形態の説明に用いた各機能要素は、例えば、集積回路として実現される。これらは、個別に1チップ化されてもよいし、一部または全てを含むように1チップ化されてもよい。また、集積回路製造後にプログラムすることが可能なFPGA(Field Programmable Gate Array)や、回路を再構成可能なリコンフィギュラブル・プロセッサを利用してもよい。
【0148】
さらに、上記各実施の形態はハードウエアに限定されるものではなく、ソフトウエアによってもよい。その逆も真である。また、それらの組み合わせであってもよい。
【産業上の利用可能性】
【0149】
本発明に係る雑音抑圧装置および雑音抑圧方法は、カルマンフィルタを用いつつ、AR係数の推定を必要とすることなく、シンプルな構成で、雑音抑圧能力を向上することができる雑音抑圧装置および雑音抑圧方法として有用である。
【図面の簡単な説明】
【0150】
【図1】本発明の実施の形態1に係る雑音抑圧装置の構成を示すブロック図
【図2】本実施の形態に係る雑音抑圧装置の状態空間モデルを説明するブロック線図
【図3】本実施の形態に係る雑音抑圧装置のカルマンフィルタアルゴリズムに基づく信号推定手順を説明するフローチャート
【図4】図3のカルマンフィルタアルゴリズムの処理内容を示すフローチャート
【図5】本発明の実施の形態2に係る雑音抑圧装置の構成を示すブロック図
【図6】本実施の形態における音声波形数値シミュレーションの第1の例の結果を示す図
【図7】本実施の形態における音声波形数値シミュレーションの第2の例の結果を示す図
【図8】本実施の形態における音声波形数値シミュレーションの第3の例の結果を示す図
【図9】本実施の形態における音声波形数値シミュレーションの第4の例の結果を示す図
【図10】本実施の形態における雑音抑圧後の音声品質の主観的評価結果の一例を示す図
【図11】本実施の形態における雑音抑圧後の音声品質の主観的評価結果の他の例を示す図
【図12】本発明の雑音抑圧方法が適用されたマルチキャリア受信装置の構成を示すブロック図
【図13】本実施の形態におけるチャネル推定精度数値シミュレーションの一つの例の結果を示す図
【図14】本実施の形態におけるチャネル推定精度数値シミュレーションの他の例の結果を示す図
【符号の説明】
【0151】
100、200 雑音抑圧装置
110 入力部
120 サンプリング部
130 A/D変換部
140 バッファ
150 雑音抑圧処理部
160 出力部
210 パーソナルコンピュータ
211 操作装置
212 ディスプレイ
213 バスインタフェース
214 記録装置
215 主記憶メモリ
216 中央演算装置
220 マイクロホン
300 マルチキャリア受信装置
301 検波部
302 GI除去部
303 チャネル等化部
304 チャネル推定部
305 FFT部
306 復調部
307 復号部
500、501 信号ベクトル
502 観測ベクトル
503 付加雑音信号ベクトル
504 駆動雑音ベクトル
505 状態遷移行列
506 状態遷移行列
図面
【図1】
0
【図2】
1
【図3】
2
【図4】
3
【図5】
4
【図6】
5
【図7】
6
【図8】
7
【図9】
8
【図10】
9
【図11】
10
【図12】
11
【図13】
12
【図14】
13