(57)【要約】      （修正有）【課題】簡素な構成で更新量の制御が可能であるとともに、系の変動への追従性の高い適応フィルタ装置を提供する。【解決手段】適応フィルタ装置は、総ステップ数２ｍ２の前半分は誤差信号ｅ（ｎ）と等しく、残りのステップは零とした信号ｅ´（ｎ）を算出し、同様にｘ（ｎ）の前半分はｘ（ｎ）と等しく、ｍ１〜（２ｍ２−１）ステップは零とした信号ｘ´（ｎ）を算出し、信号ｅ´（ｎ）と、信号ｘ´（ｎ）の共役複素数である信号ｘ＊´（ｎ）との相互相関関数を算出し、相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき誤差を処理して置換相互相関関数を算出し、置換相互相関関数を離散フーリエ変換して、置換クロススペクトルを算出し、置換入力信号のパワースペクトルを算出し、置換クロススペクトル及びパワースペクトルに基づきフィルタ係数の更新値を算出する。【選択図】図７
2007235354.tif 000002
【特許請求の範囲】
【請求項１】
伝達関数が未知である系へ入力された系入力信号と等しい入力信号が入力され、前記入力信号をフィルタ処理することにより出力信号を出力し、かつ、前記フィルタ処理におけるフィルタ係数を前記系からの出力される系出力信号と前記出力信号との差である誤差信号及び前記入力信号に基づき、前記系の伝達関数と等しくなるようにフィルタ係数を更新する適応フィルタ装置であって、
総ステップ数２ｍ２であり、０〜（ｍ２−１）ステップは誤差信号のＭｉ＋１〜Ｍｉ＋ｍ２ステップと等しく、ｍ２〜（２ｍ２−１）ステップは零とした置換誤差信号を出力する手段と、
総ステップ数２ｍ２であり、０〜（ｍ１−１）（ｍ１≦ｍ２）ステップは入力信号のＭｉ＋１〜Ｍｉ＋ｍ１ステップと等しく、ｍ１〜（２ｍ２−１）ステップは零とした置換入力信号を出力する手段と、
前記置換誤差信号と、前記置換入力信号の共役複素数である共役置換入力信号との相互相関関数を算出する手段と、
前記算出した相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき誤差を処理して置換相互相関関数を算出する手段と、
前記置換相互相関関数を離散フーリエ変換して、置換クロススペクトルを算出する手段と、
前記置換入力信号のパワースペクトルを算出する手段と、
前記置換クロススペクトル及び前記パワースペクトルに基づき前記フィルタ係数の更新値を算出する手段と、備えることを特徴とする適応フィルタ装置。
【請求項２】
前記相互相関関数を算出する手段は、
前記置換誤差信号と、前記置換入力信号とのクロススペクトルを逆離散フーリエ変換することにより相互相関関数を算出することを特徴とする請求項１記載の適応フィルタ装置。
【請求項３】
前記誤差を処理する手段は、
前記相互相関関数のｍ１〜（２ｍ２−１）ステップの値を零に置換する手段を備えることを特徴とする請求項１又は２記載の適応フィルタ装置。
【請求項４】
前記誤差を処理する手段は、
前記相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき、閾値を設定する手段と、
前記相互相関関数の０〜（２ｍ２−１）ステップの値が前記算出した閾値以下となる場合には、零に置換して、置換相互相関関数を算出する手段とを更に備えることを特徴とする請求項３記載の適応フィルタ装置。
【請求項５】
請求項１から４何れかに記載の適応フィルタ装置であって、
前記算出した置換クロススペクトルを加算平均する手段と、
前記算出したパワースペクトルを加算平均する手段と、を更に備え、
前記フィルタ係数の更新値を算出する手段は、加算平均した置換クロススペクトル及びパワースペクトルに基づき、フィルタ係数の更新値を算出することを特徴とする適応フィルタ装置。
【請求項６】
請求項１から５何れかに記載の適応フィルタ装置を用いたことを特徴とするエコーキャンセラ。
【請求項７】
伝達関数が未知である系へ入力された系入力信号と等しい入力信号が入力され、前記入力信号をフィルタ処理することにより出力信号を出力し、かつ、前記フィルタ処理におけるフィルタ係数を前記系からの出力される系出力信号と前記出力信号との差である誤差信号及び前記入力信号に基づき、前記系の伝達関数と等しくなるようにフィルタ係数を更新する適応フィルタ装置のフィルタ係数を更新する方法であって、
総ステップ数２ｍ２であり、０〜（ｍ２−１）ステップは誤差信号のＭｉ＋１〜Ｍｉ＋ｍ２ステップと等しく、ｍ２〜（２ｍ２−１）ステップは零とした置換誤差信号を出力するステップと、
総ステップ数２ｍ２であり、０〜（ｍ１−１）（ｍ１≦ｍ２）ステップは入力信号のＭｉ＋１〜Ｍｉ＋ｍ１ステップと等しく、ｍ１〜（２ｍ２−１）ステップは零とした置換入力信号を出力するステップと、
置換入力信号の共役を取り、共役置換入力信号を出力するステップと、
前記置換誤差信号と、前記共役置換入力信号の相互相関関数を算出するステップと、
前記算出した相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき誤差を処理して置換相互相関関数を算出するステップと、
前記置換相互相関関数を離散フーリエ変換して、置換クロススペクトルを算出するステップと、
前記置換入力信号のパワースペクトルを算出するステップと、
前記置換クロススペクトル及び前記パワースペクトルに基づき前記フィルタ係数の更新値を算出するステップと、備えることを特徴とするフィルタ係数の更新値を算出する方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明はクロススペクトル法を用いてフィルタ係数の更新量を算定する適応フィルタ装置、フィルタ係数の更新方法に関する。
【背景技術】
【０００２】
従来より、携帯電話や遠隔会議システムなどには、出力された音声が入力機器に拾われエコーやハウリングが起こるのを防止するため、エコーキャンセラが用いられている。図１６は、携帯電話などに用いられている従来のエコーキャンセラの構成の一例を示す図である。同図に示すように、エコーキャンセラにはＬＭＳアルゴリズムを用いた適応フィルタを備えている。マイクにより収音された音響信号には、携帯電話の利用者の発話音声の信号に加えて、スピーカからマイクに到達した入力信号が含まれる。エコーキャンセラは、適応フィルタのフィルタ係数をスピーカからマイクまでの伝達関数と等しくなるように更新し、このフィルタ係数を用いて入力信号を処理した信号を、マイクに収音された音響信号から減算する。これにより、エコーキャンセラより出力された信号（以下、誤差信号という）は、スピーカより出力され、マイクに収音された音声を除去することができる。
【０００３】
このような適応フィルタのフィルタ更新アルゴリズムであるＬＭＳアルゴリズムでは、非特許文献１に記載されているように、下記の式（１）のように、入力信号と誤差信号との畳込積分の解である推定量に対して、ステップサイズパラメータμを積算して更新量を算出することにより、システムの安定性を確保している。ただし、ｈｐ（ｔ）：適応前のフィルタ係数のベクトル、ｈｐ（ｔ+１）：適応後のフィルタ係数のベクトル、μ：ステップサイズパラメータ、ｅ（ｔ）：ｎ番目の繰り返し時の誤差信号、ｕ（ｔ）：ｔ番目の繰り返し時のタップ入力ベクトルとする。
【数１】
2007235354.tif 000003
【非特許文献１】Ｓ．ヘイキン著、武部幹  訳、“適応フィルタ入門”現代工学社、１９８７年９月１０日初版発行、ｐ．１０８−１０９
【非特許文献２】Manabu Fukushima, Hiroto Inoue, Ken'itiro Kamura, Hirofumi Yanagawa, Ken'iti Kido, "A method for the determination of noise factor in estimated transfer function - cross spectral technique by use of 1-0 and 1-000 windows -", Proc. of The 18th International Congress on Acoustics, pp.166-169, vol.25, no.2, 2004
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、ＬＭＳアルゴリズムを用いてフィルタ係数の更新する方法では、フィルタ係数の更新量を制御するために、推定精度によらずステップサイズパラメータμを一定として使用しているため、入力信号ｘ（ｎ）として白色信号が入力され、精度良く推定量（式（１）におけるｅ（ｎ）ｕ（ｎ））を算出できた場合でも、ステップサイズパラメータμを乗算してしまうため、系の変化への追従性が低下する一因となっている。
また、逆に、入力信号ｘ（ｎ）として周波数特性に偏りのある有色信号が入力され、推定量の算出精度が悪い場合でも、推定量にステップサイズパラメータμを乗算した分だけ更新してしまうため、フィルタ係数の精度が低くなってしまう。
【０００５】
ステップサイズパラメータμを誤差信号ｅ（ｎ）に基づき、動的に変動させることも可能であるが、誤差信号ｅ（ｎ）はスピーカからマイクまでの系の変動のみではなく、発話音声ｖ（ｔ）の影響を受けるため、発話音声が発生したのか、系が変動したのかを検知しなければならない。しかし、発話音声の発生を検知するためには複雑な構成が必要となるため、携帯電話等の小型の装置には適していない。
【０００６】
本発明は、上記の問題に鑑みなされたものであり、簡素な構成で更新量の制御が可能であるとともに、系の変動への追従性の高い適応フィルタ装置を提供することを目的とする。
【課題を解決するための手段】
【０００７】
本発明の適応フィルタ装置は、伝達関数が未知である系へ入力された系入力信号と等しい入力信号が入力され、前記入力信号をフィルタ処理することにより出力信号を出力し、かつ、前記フィルタ処理におけるフィルタ係数を前記系からの出力される系出力信号と前記出力信号との差である誤差信号及び前記入力信号に基づき、前記系の伝達関数と等しくなるようにフィルタ係数を更新する適応フィルタ装置であって、総ステップ数２ｍ２であり、０〜（ｍ２−１）ステップは誤差信号のＭｉ＋１〜Ｍｉ＋ｍ２ステップと等しく、ｍ２〜（２ｍ２−１）ステップは零とした置換誤差信号を出力する手段と、総ステップ数２ｍ２であり、０〜（ｍ１−１）（ｍ１≦ｍ２）ステップは入力信号のＭｉ＋１〜Ｍｉ＋ｍ１ステップと等しく、ｍ１〜（２ｍ２−１）ステップは零とした置換入力信号を出力する手段と、前記置換誤差信号と、前記置換入力信号の共役複素数である共役置換入力信号との相互相関関数を算出する手段と、前記算出した相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき誤差を処理して置換相互相関関数を算出する手段と、前記置換相互相関関数を離散フーリエ変換して、置換クロススペクトルを算出する手段と、前記置換入力信号のパワースペクトルを算出する手段と、前記置換クロススペクトル及び前記パワースペクトルに基づき前記フィルタ係数の更新値を算出する手段と、を備えることを特徴とする。
【０００８】
ここで、前記相互相関関数を算出する手段は、前記置換誤差信号と、前記置換入力信号とのクロススペクトルを逆離散フーリエ変換することにより相互相関関数を算出してもよい。
【０００９】
前記誤差を処理する手段は、前記相互相関関数のｍ１〜（２ｍ２−１）ステップの値を零に置換する手段を備えてもよい。また、前記誤差を処理する手段は、前記相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき、閾値を設定する手段と、前記相互相関関数の０〜（２ｍ２−１）ステップの値が前記算出した閾値以下となる場合には、零に置換して、置換相互相関関数を算出する手段とを更に備えてもよい。
【００１０】
また、上記の適応フィルタは、前記算出した置換クロススペクトルを加算平均する手段と、前記算出したパワースペクトルを加算平均する手段と、を更に備え、前記フィルタ係数の更新値を算出する手段は、加算平均した置換クロススペクトル及びパワースペクトルに基づき、フィルタ係数の更新値を算出してもよい。
【００１１】
本発明の適応フィルタ装置によれば、更新量の推定精度に基づき更新量の制御を行うためステップパラメータμが不要となり、安定性を確保しつつ、系の変動に対する追従性を向上することができる。また、簡素な演算装置からなるため、携帯電話等の小型の装置への適用も可能である。
【００１２】
また、本発明は、上記の適応フィルタを用いたエコーキャンセラを含むものとする。また、本発明は、伝達関数が未知である系へ入力された系入力信号と等しい入力信号が入力され、前記入力信号をフィルタ処理することにより出力信号を出力し、かつ、前記フィルタ処理におけるフィルタ係数を前記系からの出力される系出力信号と前記出力信号との差である誤差信号及び前記入力信号に基づき、前記系の伝達関数と等しくなるようにフィルタ係数を更新する適応フィルタ装置のフィルタ係数を更新する方法であって、総ステップ数２ｍ２であり、０〜（ｍ２−１）ステップは誤差信号のＭｉ＋１〜Ｍｉ＋ｍ２ステップと等しく、ｍ２〜（２ｍ２−１）ステップは零とした置換誤差信号を出力するステップと、総ステップ数２ｍ２であり、０〜（ｍ１−１）（ｍ１≦ｍ２）ステップは入力信号のＭｉ＋１〜Ｍｉ＋ｍ１ステップと等しく、ｍ１〜（２ｍ２−１）ステップは零とした置換入力信号を出力するステップと、置換入力信号の共役を取り、共役置換入力信号を出力するステップと、前記置換誤差信号と、前記共役置換入力信号の相互相関関数を算出するステップと、前記算出した相互相関関数のｍ２〜（２ｍ２−１）ステップの値に基づき誤差を処理して置換相互相関関数を算出するステップと、前記置換相互相関関数を離散フーリエ変換して、置換クロススペクトルを算出するステップと、前記置換入力信号のパワースペクトルを算出するステップと、前記置換クロススペクトル及び前記パワースペクトルに基づき前記フィルタ係数の更新値を算出するステップと、備えることを特徴とするフィルタ係数の更新値を算出する方法も含むものとする。
【発明の効果】
【００１３】
本発明によれば、適応フィルタとしての安定性を確保しつつ、未知の系の変動に対するフィルタ係数更新の追従性を向上することができる。
【発明を実施するための最良の形態】
【００１４】
以下、本発明の適応フィルタの一実施形態について、図面に基づき詳細に説明する。
図１は、本実施形態の適応フィルタ装置１１を用いた装置の一例であるエコーキャンセラ１０の構成を示す図である。エコーキャンセラ１０は、例えば、携帯電話などの、外部より入力された入力信号ｘ（ｔ）を音声出力するスピーカ１２と、発話者の発話音声ｖ（ｔ）を収集するマイク１３とを備えた装置において、スピーカ１２から音声出力された入力信号ｘ（ｔ）がマイク１３に再び収音されることにより生じるエコーやハウリングを防止するために設けられている。
【００１５】
ここで、スピーカ１２からマイク１３までの系を考えると、スピーカ１２に伝送される入力信号をｘ（ｔ）、スピーカ１２からマイク１３までの伝達関数をｈ（ｔ）、発話音声などの音声信号をｖ（ｔ）とすると、マイク１３の収音する音響信号ｙ（ｔ）は次式（２）で表される。
【数２】
2007235354.tif 000004
【００１６】
また、適応フィルタ装置１１のフィルタ係数をｈｐ（ｔ）とすると、適応フィルタ装置１１より出力される出力信号ｄ（ｔ）は次式（３）で表される。
【数３】
2007235354.tif 000005
【００１７】
ここで、適応フィルタ装置１１のフィルタ係数ｈｐ（ｔ）と、伝達関数ｈ（ｔ）との差をΔｈ（ｔ）とすると、エコーキャンセラ１０から出力される誤差信号ｅ（ｔ）は次式（４）で表される。
【数４】
2007235354.tif 000006
【００１８】
ここで、適応フィルタ装置１１のフィルタ係数ｈｐ（ｔ）が系の伝達関数ｈ（ｔ）と等しくなるように更新することにより、Δｈ（ｔ）が零に近づく。これにより、伝送されていく信号ｅ（ｔ）が音声信号ｖ（ｔ）と略等しくなる。
【００１９】
例えば、携帯電話のエコーキャンセラとして用いられる場合などのように、系の伝達関数ｈ（ｔ）は変動しやすい場合にも、適応フィルタ装置１１が系の入力信号ｘ（ｔ）、及び系の出力信号ｙ（ｔ）と出力信号ｄ（ｔ）との差である誤差信号ｅ（ｔ）に基づき、フィルタ係数ｈｐ（ｔ）を系の伝達関数ｈ（ｔ）と等しくなるように逐次更新することにより、エコーのない音声信号を送信することができる。
【００２０】
このように、フィルタ係数を更新する方法の一つに、ＬＭＳアルゴリズムがある。ＬＭＳアルゴリズムでは、次式（５）によりフィルタ係数ｈｐ（ｔ）を更新する。
【数５】
2007235354.tif 000007
式（５）の右辺第２項はフィルタ係数の更新量を示すが、このようにフィルタ係数の更新量は、推定量（ｙ（ｔ）ｕ（ｔ））にステップサイズパラメータμを積算することにより、算出している。このようにステップサイズパラメータμを積算するのは、推定量の算出精度は入力信号の周波数特性の影響を大きく受けるので、適応フィルタの安定性を確保するためである。しかし、ステップサイズパラメータμを一律に設定するために、推定量を高精度に算出できた場合にも、更新量はステップサイズパラメータμを積算してしまうためにその値が小さくなり、系の変化に対する追従性が低下してしまう。
【００２１】
信号の周波数特性の影響を受けずに安定な更新を行うことができる手法の一つに、信号をブロック化し、周波数領域で更新量を算出するブロックＬＭＳアルゴリズムがある。ブロックＬＭＳアルゴリズムでは、入力信号ｘ（ｔ）及び誤差信号ｅ（ｔ）を周波数領域に変換するためにＤＦＴを行い、周波数領域で更新量を算出したのち、ＩＤＦＴを行い、時間領域における更新量を算出する。しかし、ＤＦＴは有限長の信号を繰り返しにより無限長として信号を扱ってフーリエ変換を行うため、再び、更新量を周波数領域から時間領域に戻す際に時間窓外にはみ出した信号が、逆側に表れるという現象が生じ、計算誤差の生じる原因となっていた。
【００２２】
ここで、発明者らは、非特許文献２において、線形時不変な系の伝達関数をクロススペクトル法を用いて、より高精度に算出する方法を提案している。この方法によれば、伝達関数を推定する過程でＤＦＴの巡回性に起因して生じる誤差（以下、推定誤差という）の大きさに基づき、推定精度を求め、これに基づき推定誤差を処理することにより、より高精度に伝達関数を算出することができる。
そこで、発明者らは、適応フィルタのフィルタ係数の更新値をクロススペクトル法を用いて算出し、この際、上記の文献において発明者らが提案する推定誤差の処理方法を用いて推定誤差を見積もり、これに基づき誤差の処理を行うことにより、入力信号が周波数特性に偏りある場合でも、ステップサイズパラメータμを用いることなく、適応フィルタの安定性を確保することができると考えた。
【００２３】
そこで、まず、上記の文献において、発明者らの提案する線形時不変を仮定できる系の伝達関数をクロススペクトル法により算出する際に、ＤＦＴの巡回性に起因して生じる誤差を処理する方法について説明する。一例として、図２に示すような、入力信号ｘ（ｎ）が入力し、出力信号ｙ（ｎ）が出力される系の伝達関数ｈ（ｎ）を求める場合について説明する。同図に示す系の伝達関数ｈ（ｎ）はクロススペクトル法によれば、次式（５）で表される。なお、式中のＸ（ｋ）は系の入力信号ｘ（ｎ）をＤＦＴした信号を、Ｙ（ｋ）は系の出力信号ｙ（ｎ）をＤＦＴした信号を、Ｈ（ｋ）は系の推定伝達関数ｈ（ｎ）をＤＦＴした信号を、Ｇｘｙ（ｎ）はｘ（ｎ）とｙ（ｎ）のクロススペクトルを、Ｇｘ（ｎ）はｘ（ｎ）のパワースペクトルを示す。
【数６】
2007235354.tif 000008
なお、クロススペクトルＧｘｙ（ｎ）、及びパワースペクトルＧｘ（ｎ）は以下の式（７）及び（８）で算出される。
【数７】
2007235354.tif 000009
【００２４】
また、クロススペクトルＧｘｙ（ｎ）と、ｘ（ｎ）と、ｙ（ｎ）との相互相関関数Ｒｘｙ（ｎ）との間には、次式（９）が成立する。
【数８】
2007235354.tif 000010
【００２５】
同様に、パワースペクトルＧｘ（ｎ）と、ｘ（ｎ）の自己相関関数Ｒｘ（ｎ）との間には、以下の関係式が成立する。
【数９】
2007235354.tif 000011
【００２６】
ここで、例えば、総ステップ数が２Ｎである入力信号ｘ（ｎ）及び出力信号ｙ（ｎ）の、相互相関関数Ｒｘｙ（ｎ）は、次式（１１）で表される。
【数１０】
2007235354.tif 000012
ここで、ＤＦＴは、有限長の信号を無限に繰り返されているものと仮定することで、フーリエ変換できるようにしているため、周波数伝達関数Ｈ（ｋ）を周波数領域から実数領域に逆フーリエ変換（以下、ＩＤＦＴという）すると、時間窓より前及び後にはみ出す部分は、信号の逆側に表れてしまう。
このため、式（７）により、周波数領域でクロススペクトルを算出しているが、実数領域について考えると、上記の式（１１）は、次式（１２）のように計算していることになる。
【数１１】
2007235354.tif 000013
【００２７】
このため、クロススペクトル法により伝達関数を算出する際に、周波数領域で周波数伝達関数Ｈ（ｋ）を周波数領域から実数領域に逆フーリエ変換（以下、ＩＤＦＴという）すると、時間窓より前及び後にはみ出す部分は、信号の逆側に誤差として表れてしまう。
【００２８】
そこで、発明者らの提案する方法は、非特許文献２に記載されているように、クロススペクトル法の算出式の式（６）における入力信号ｘ（ｎ）及び出力信号ｙ（ｎ）に代えて、以下に説明する置換入力信号及び置換出力信号を用いて、ＤＦＴの巡回性に起因する誤差を分離し、これに基づき推定量の算出精度を見積もるものである。
【００２９】
置換出力信号ｙ´（ｎ）は総ステップ数が２ｍ２であって、０〜（ｍ２−１）ステップの値は出力信号の（Ｍ＋１）〜（Ｍ＋ｍ２）ステップの値と等しく、ｍ２〜（２ｍ２−１）ステップの値は、零である信号である。このような置換出力信号ｙ´（ｎ）は、まず、図３に示すように、出力信号ｙ（ｎ）のＭｉ＋１ステップよりステップ数ｍ２の信号を切り出し、次に、図４に示すように、この信号の後ろにステップ数ｍ２であり全ての値が零である信号を追加することにより得られる。
【００３０】
また、置換入力信号ｘ´（ｎ）は総ステップ数が２ｍ２であって、０〜（ｍ１−１）ステップの値は入力信号の（Ｍ＋１）〜（Ｍ＋ｍ１）ステップと等しく、ｍ１〜（２ｍ２−１）ステップの値は零である信号である。このような置換入力信号ｘ´（ｎ）は、まず、図３に示すように、入力信号ｘ（ｎ）のＭｉ＋１ステップよりｍ１（ｍ１≦ｍ２）ステップの信号を切り出し、次に、図４に示すように、この信号の後ろに（２ｍ２−ｍ１）ステップの値が零である信号を追加することにより得られる。
【００３１】
これらの置換入力信号ｘ´（ｎ）及び置換出力信号ｙ´（ｎ）を、式（５）における入力信号ｘ（ｎ）及び出力信号ｙ（ｎ）に代えて用いることにより次式（１３）が導かれる。なお、Ｘ´（ｋ）及びＹ´（ｋ）は、置換入力信号ｘ´（ｎ）及び置換出力信号ｙ´（ｎ）をＤＦＴした信号を表す。
【数１２】
2007235354.tif 000014
【００３２】
置換入力信号ｘ´（ｎ）及び置換出力信号ｙ´（ｎ）を用いた場合について、相互相関関数について考えると、次式（１４）のように表される。
【数１３】
2007235354.tif 000015
【００３３】
式（１４）において、ｍ２≦ｎ≦２ｍ２−１でｘ（ｎ）＝０であるため、右辺第３項、第４項は零となり、また、式（１４）における右辺第２項はｍ＞ｍ２と、ｍ＜ｍ２とで分割され次式（１５）のように表すことができる。
【数１４】
2007235354.tif 000016
【００３４】
ここで、ｙ（ｎ）はｍ２≦ｎ≦２ｍ２−１で零のため式（１５）で右辺第２項は零との積和となり零である。また，式（１５）第３項はｘ（ｎ）のインパルス応答を１つも含まないｙ（ｎ）との積和であるため、全て無相関な信号系列の積和、すなわちＤＦＴの巡回性に起因する誤差である。観測点数と同数でＤＦＴして推定したインパルス応答は、式（１５）右辺第１項と第３項の和である。したがって式（１５）右辺第３項は、ＤＦＴの巡回性に起因する誤差の総量の半分に相当する。このことから、式（１５）右辺第３項より誤差の推定を行う。この項は，推定値のｍ２〜２ｍ２−１の範囲に表われることから、この範囲の係数を用いて推定値の誤差の大きさを見積もることができる。
一例として、図５にｍ２＝５の場合における相互相関関数Ｒｘｙ（ｍ）を算出する過程を示す図を示す。同図に円で囲んだ部分がＤＦＴの巡回性に起因する誤差であり、式（１５）のように、ｍ２〜２ｍ２−１の範囲に分離されていることが確認できる。
【００３５】
さらに、発明者らは、“拡張型クロススペクトル法を用いた音場変化の短時間追従に関する検討”, 石丸善隆, 村上尚, 長尾優次, 加村健一郎, 福島学, 柳川博文, 日本音響学会春季発表会, pp.507-508,2004  に記載されているように、上記のようにして算定した推定値の精度に基づき、例えば、図６に示すように、算出した推定値のｍ２〜２ｍ２−１の範囲の最大値を閾値として設定し、推定値の０〜２ｍ２−１の範囲の値が、閾値以下となる場合には、零に置換するといった方法により、更新量を制御する試みを行っている。
【００３６】
本実施形態の適応フィルタ装置は、以上説明したクロススペクトル法の推定誤差を処理する技術を適応フィルタにおけるフィルタ係数の更新量の算出に応用したものであり、フィルタ係数の推定誤差を算出し、上記の方法と同様に、推定誤差の計算精度を見積り、この計算精度に基づき閾値を設定し、算出したフィルタ係数の更新量において、前記設定した閾値よりも小さい成分を零に置換することにより、ステップサイズパラメータμを用いなくとも安定なフィルタ係数の更新を可能としている。
【００３７】
次に、上記のクロススペクトル法における誤差の処理方法を応用した適応フィルタのフィルタ係数の更新量の算出方法について詳細に説明する。
ここでは、図１で説明したような、入力信号ｘ（ｔ）が入力し、フィルタ係数ｈｐ（ｎ）により処理することにより出力信号ｙ（ｔ）が出力し、かつ、フィルタ係数ｈｐを、入力信号ｘ（ｔ）及び系の出力信号ｄ（ｔ）と出力信号ｙ（ｔ）の差である誤差信号ｅ（ｔ）に基づき、更新する適応フィルタを例にとって説明する。
まず、クロススペクトル法を用いた適応フィルタにおけるフィルタ係数の更新量の算出方法について定式化する。
【００３８】
マイクの収音した信号をＡ／Ｄ変換した離散化信号ｙ（ｎ）は次式（１６）であらわされる。
【数１５】
2007235354.tif 000017
【００３９】
また、誤差信号ｅ（ｎ）は次式（１７）であらわされる。
【数１６】
2007235354.tif 000018
【００４０】
ここで、更新量Δｈ（ｎ）を次式（１８）で表すと、式（１７）は式（１９）のように表される。
【数１７】
2007235354.tif 000019
【００４１】
この式（１９）を、ＤＦＴにより周波数領域に変換すると、以下の式が導かれる。
【数１８】
2007235354.tif 000020
【００４２】
ここで、両辺にＸ（ｋ）の共役複素数であるＸ＊（ｋ）を乗算すると、次式（２５）が導かれる。
【数１９】
2007235354.tif 000021
【００４３】
さらに、入力信号を連続からブロックへと分割しているため、夫々のブロックの平均を考えると次式（２６）が導かれる。
【数２０】
2007235354.tif 000022
【００４４】
ここで、入力信号ｘ（ｎ）と、騒音信号ｖ（ｎ）とは、無相関であり、上記の式（２６）の右辺第二項は、加算平均することにより零に近づくため、次式（２７）が導かれる。
【数２１】
2007235354.tif 000023
【００４５】
よって、上記の式（２６）は次式（２８）のようになる。
【数２２】
2007235354.tif 000024
【００４６】
よって、ΔＨ（ｋ）は、次式（２９）で算出することができる。
【数２３】
2007235354.tif 000025
【００４７】
次式（２９）のように、ＩＤＦＴにより時間領域に変換することにより、フィルタ係数の更新量Δｈｐ（ｋ）は次式（３０）により算出される。
【数２４】
2007235354.tif 000026
【００４８】
式（３０）に示すように、更新量Δｈｐ（ｋ）は入力信号ｘ（ｎ）及び誤差信号ｅ（ｎ）のクロススペクトル及び入力信号ｘ（ｎ）のパワースペクトルに基づき算出できる。このため、上記説明したクロススペクトル法における推定精度の見積が可能になり、これに基づき更新量を制御することにより、ステップパラメータμを用いる必要がなくなる。よって、フィルタ係数の更新量は次式（３１）により算出することができる。
【数２５】
2007235354.tif 000027
【００４９】
上記のように定式化したフィルタ係数の更新量の算出式を用いた本実施形態の適応フィルタ装置におけるフィルタ更新量を算出する流れを以下説明する。
【００５０】
図７は、本実施形態の適応フィルタ装置における処理の流れを示す図である。適応フィルタ装置には、図１に示した通り、入力信号ｘ（ｔ）と、系の出力信号ｙ（ｔ）と適応フィルタの出力信号ｄ（ｔ）の差である誤差信号ｅ（ｔ）が入力される。
【００５１】
入力された誤差信号ｅ（ｔ）は、ステップ１１０において、Ａ／Ｄ変換され、離散信号のｅ（ｎ）に変換される。
次に、ステップ１１２において、図８に示すように、離散信号ｅ（ｎ）の、Ｍｉ＋１〜Ｍｉ＋ｍ２（ｉ＝１、…、Ｌ）ステップを切り出し、ｅｉ（ｎ）を出力する。
次に、ステップ１１４において、図９に示すように、切り出された信号ｅｉ（ｎ）の後部に長さがｍ２ステップであり、値が零である信号を接続し、ｅ´ｉ（ｎ）を出力する。
これにより、ｅ´ｉ（ｎ）は、総ステップ数２ｍ２であり、０〜ｍ２−１ステップは誤差信号ｅ（ｎ）のＭｉ＋１〜Ｍｉ＋ｍ２ステップと等しく、ｍ２〜２ｍ２−１ステップは零である信号（請求項における置換誤差信号に相当）となる。
次に、ステップ１１６において、ｅ´ｉ（ｎ）を次式（３２）のように、高速フーリエ変換（以下、ＦＦＴという）して、Ｅ´ｉ（ｋ）を算出する。
【数２６】
2007235354.tif 000028
【００５２】
上記のステップ１１０〜１１６と並行して、入力信号ｘ（ｔ）に対して以下の操作を行う。
入力された入力信号ｘ（ｔ）は、ステップ１００において、Ａ／Ｄ変換され、離散信号のｘ（ｎ）に変換される。
次に、ステップ１０２において、ｘ（ｎ）の、Ｍｉ＋１〜Ｍｉ＋ｍ１（ｉ＝１、…、Ｌ）ステップを切り出し、ｘｉ（ｎ）を出力する。なお、ここで信号を切り出す長さｍ１ステップは、ステップ１１２において誤差信号を切り出す長さｍ２ステップ以下（ｍ１≦ｍ２）であるものとする。
次に、ステップ１０４において、切り出された信号ｘｉ（ｎ）の後部に長さが（２ｍ２−ｍ１）ステップであり、値が零である信号を接続し、ｘ´ｉ（ｎ）を出力する。
これにより、ｘ´ｉ（ｎ）は総ステップ数２ｍ２であり、０〜ｍ１−１（ｍ１≦ｍ２）ステップは入力信号のＭｉ＋１〜Ｍｉ＋ｍ１ステップと等しく、ｍ１〜２ｍ２−１ステップは零である信号（すなわち、請求項における置換入力信号に相当）となる。
次に、ステップ１０６において、ｘ´ｉ（ｎ）を次式（３３）のように、ＦＦＴして、Ｘ´ｉ（ｋ）を算出する。
【数２７】
2007235354.tif 000029
次に、ステップ１０８において、Ｘ´ｉ（ｋ）の共役複素数を算出し、Ｘ´＊ｉ（ｋ）とする。
【００５３】
次に、ステップ１１８において、ステップ１０６及びステップ１０８で算出したＸ´ｉ（ｋ）とＸ´＊ｉ（ｋ）を積算する。
さらに、ステップ１２０において、Ｘ´ｉ（ｋ）Ｘ´＊ｉ（ｋ）をｉ＝１、…、Ｌについて加算平均する。これにより、次式（３４）により置換パワースペクトルＧｘ´が算出される。
【数２８】
2007235354.tif 000030
【００５４】
上記のステップ１１８〜１２０と並行して、以下の工程を行う。
ステップ１２２において、前記ステップ１０８で算出したＸ´＊ｉ（ｋ）と、ステップ１１６で算出したＥ´＊ｉ（ｋ）とを積算する。
次に、ステップ１２４において、前記算出したＸ´＊ｉ（ｋ）Ｅ´＊ｉ（ｋ）を逆高速フーリエ変換（以下、ＩＦＦＴという）する。次式（３５）に示すように、これにより得られる信号をｚｉ（ｎ）とする。
【数２９】
2007235354.tif 000031
次に、ステップ１２６において、前記算出した信号ｚｉ（ｎ）の推定誤差処理を行う。前記クロススペクトル法において推定誤差処理を行った場合と同様、図６に示すように、まず、信号ｚｉ（ｎ）のステップｍ２〜２ｍ２−１における絶対値の最大値を閾値として設定する。そして、ステップ０〜２ｍ２−１について、信号ｚｉ（ｎ）の絶対値が設定した閾値以下となる場合にはその値を零に置換する。なお、閾値はステップｍ２〜２ｍ２−１における絶対値の最大値であるので、ステップｍ２〜２ｍ２−１のｚｉ（ｎ）の値は必ず零に置換される。このように推定誤差処理を行った信号をｚ´ｉ（ｎ）（請求項における置換相互相関関数に相当）とする。
次に、ステップ１２８において、次式（３６）のように、前記算出した信号ｚ´ｉ（ｎ）をＦＦＴする。これにより得られた信号をＺ´ｉ（ｋ）とする。
【数３０】
2007235354.tif 000032
次に、ステップ１３０において、Ｚ´ｉ（ｋ）をｉ＝１、…、Ｌについて加算平均する。これにより、次式（３７）のように置換クロススペクトルＧｘ´ｅ´が算出される。
【数３１】
2007235354.tif 000033
【００５５】
次に、ステップ１３２において、次式（３８）に示すように、ステップ１３０で算出された置換クロススペクトルＧｘ´ｅ´をステップ１２０で算出された置換パワースペクトルＧｘ´で除算することにより、周波数領域における更新量ΔＨ（ｋ）が算出される。
【数３２】
2007235354.tif 000034
ステップ１３４において、次式（３９）のように上記算出した周波数領域における更新量ΔＨ（ｋ）をＩＦＦＴすることにより、フィルタ係数の更新量Δｈ（ｎ）が算出される。
【数３３】
2007235354.tif 000035
ステップ１３６において、次式（４０）のように、前記算出したフィルタ係数の更新量Δｈ（ｎ）を前ステップにおけるフィルタ係数ｈ（ｎ）に加えることにより、次のステップにおけるフィルタ係数Δｈ（ｎ+ＭＬ）が算出される。
【数３４】
2007235354.tif 000036
【００５６】
適応フィルタ装置は、以上説明したように更新したフィルタ係数を用いて、入力信号ｘ（ｎ）を処理することにより出力信号ｄ（ｎ）は算出し、算出したｄ（ｎ）をＤ／Ａ変換して、出力する。
【００５７】
以上説明したように、本実施形態の適応フィルタ装置によれば、クロススペクトルを算出する際に、離散フーリエ変換の巡回性により生じる推定誤差を分離し、この分離した推定誤差に基づき、更新量を制御するため、ステップサイズパラメータが不要となるため、安定性を確保しつつ、系の変動に対する追従性を向上することができる。また、簡素な演算装置からなるため、携帯電話等の小型の装置への適用も可能である。
なお、以上説明した実施形態では、ステップ１２６において、ＤＦＴの巡回性に起因する誤差の大きさに基づき閾値を設定して、信号ｚｉ（ｎ）の０〜２ｍ２−１ステップについて、値が閾値以下となる場合に零に置換することにより誤差の処理を行う構成としたが、これに限らず、クロススペクトル法と同様に、信号ｚｉ（ｎ）ｍ２〜２ｍ２−１ステップを零に置換するのみでも、ＤＦＴの巡回性に起因する誤差を低減できる。
【００５８】
以下、図１に示すエコーキャンセラに、ＬＭＳアルゴリズムにより更新を行う適応フィルタを用いた場合と、上記説明した適応フィルタを用いた場合について計算機を用いたシミュレーションによる検討を行ったので説明する。
【００５９】
本検討では、系の伝達関数ｈ（ｔ）として、防音室にて計測したインパルス応答（板振動スピーカから本体マイクまでの系）を用い、入力信号として白色雑音を用いて検討を行った。また、シミュレーションにおけるサンプリング周波数は８ｋＨｚとした。本発明の適応フィルタ及び比較対象のＬＭＳアルゴリズムの適応フィルタのフィルタ係数は、ともに１２８ステップとした。ＬＭＳアルゴリズムにおけるステップサイズパラメータμは０．０１とした。
【００６０】
図１０（Ａ）は、ＬＭＳアルゴリズムを用いた適応フィルタによるエコーキャンセル効果を示すグラフであり、同図（Ｂ）は、信号の切り出し長さを１００ｍｓとした本発明の適応フィルタによるエコーキャンセル効果を示す図である。（Ａ）及び（Ｂ）を比較するとわかるように、本発明の適応フィルタはＬＭＳアルゴリズムを用いた適応フィルタに比べ、短時間でエコー成分が低減している。これにより、系の変化の追従性が向上していることが確認された。
【００６１】
また、図（Ａ）は、ＬＭＳアルゴリズムを用いた適応フィルタによるエコーキャンセル効果を示すグラフであり、同図（Ｂ）は、信号の切り出し長さを１５０ｍｓとした本発明の適応フィルタによるエコーキャンセル効果を示す図である。信号の切り出し長さを１５０ｍｓとした場合も同様に、系の変化に対する追従性が向上していることがわかる。
【００６２】
図１１〜１４は、夫々、信号の切り出し長さを１００ｍｓとした本発明の適応フィルタの１６０ｍｓ、２２５ｍｓ、２８８ｍｓ、３５２ｍｓにおける入力信号ｘ（ｔ）と、誤差信号ｅ（ｔ）との周波数領域における伝達関数を示すグラフである。各図において、実線は入力信号ｘ（ｔ）と、誤差信号ｅ（ｔ）との周波数伝達関数を、破線は未知の系の伝達関数を示す。すなわち、破線と実線との差の大きさが、各周波数におけるエコーキャンセル効果を示す。
【００６３】
図１１〜１４に示すように、時間の経過とともに、全周波帯域で略等しくエコーキャンセル効果が得られていることがわかる。このように、本発明の適応フィルタは、周波数による偏りなく、全周波数に対して有効であることが確認された。
【００６４】
また、図１５は、入力信号が、音声（女性アナウンサ音声）である場合のエコーキャンセル効果を示すグラフであり、同図（Ａ）は入力信号を、同図（Ｂ）は、ＬＭＳアルゴリズムを用いた適応フィルタの誤差信号ｅ（ｔ）を、同図（Ｃ）は、信号の切り出し長さを１００ｍｓとした本発明の適応フィルタの誤差信号ｅ（ｔ）を、同図（Ｄ）は、信号の切り出し長さを１５０ｍｓとした本発明の適応フィルタの誤差信号ｅ（ｔ）を示す。
【００６５】
同図に示すように、ＬＭＳアルゴリズムでは、入力信号の有色性が推定結果に影響を及ぼすため、エコーキャンセル効果が低下してしまい、誤差信号に入力信号のエコーが大きく見られる。しかし、同図（Ｃ）、（Ｄ）に示すように、本発明の適応フィルタによれば、入力信号が有色性を持つ場合にも、推定結果にほとんど影響を及ぼさず、良好なエコーキャンセル効果が得られている。このことから、本発明の適応フィルタは、入力信号の周波数特性によらず、良好な推定を行うことができることが確認された。
【図面の簡単な説明】
【００６６】
【図１】適応フィルタ装置の用いられたエコーキャンセラの構成を示す図である。
【図２】クロススペクトル法による伝達関数の算定の対象となる未知の系を示す図である。
【図３】クロススペクトル法における置換入力信号及び置換出力信号を取得する方法を示す図（その１）である。
【図４】クロススペクトル法における置換入力信号及び置換出力信号を取得する方法を示す図（その２）である。
【図５】相互相関関数に表れるＤＦＴの巡回性に起因した誤差を説明するための図である。
【図６】分離した誤差に基づき、更新量を求める方法を説明するための図である。
【図７】本実施形態の適応フィルタ装置の処理の流れを示す図である。
【図８】誤差信号より置換誤差信号を取得する方法を説明するための図（その１）である。
【図９】誤差信号より置換誤差信号を取得する方法を説明するための図（その２）である。
【図１０】ＬＭＳアルゴリズムを用いた適応フィルタと、本実施形態の適応フィルタのエコーキャンセラ効果を比較したグラフである。
【図１１】本実施形態の適応フィルタの１６０ｍｓにおける誤差信号を周波数軸で表したグラフである。
【図１２】本実施形態の適応フィルタの２２４ｍｓにおける誤差信号を周波数軸で表したグラフである。
【図１３】本実施形態の適応フィルタの２８８ｍｓにおける誤差信号を周波数軸で表したグラフである。
【図１４】本実施形態の適応フィルタの３５２ｍｓにおける誤差信号を周波数軸で表したグラフである。
【図１５】入力信号が、音声（女性アナウンサ音声）である場合の、時間領域におけるエコーキャンセル効果を示すグラフである。
【図１６】従来のＬＭＳアルゴリズムを用いた適応フィルタを用いたエコーキャンセラの構成を示す図である。
【符号の説明】
【００６７】
１０  エコーキャンセラ
１１  適応フィルタ装置
１２  スピーカ
１３  マイク
２０  系
【図１】
2007235354.tif 000037
【図２】
2007235354.tif 000038
【図３】
2007235354.tif 000039
【図４】
2007235354.tif 000040
【図５】
2007235354.tif 000041
【図６】
2007235354.tif 000042
【図７】
2007235354.tif 000043
【図８】
2007235354.tif 000044
【図９】
2007235354.tif 000045
【図１０】
2007235354.tif 000046
【図１１】
2007235354.tif 000047
【図１２】
2007235354.tif 000048
【図１３】
2007235354.tif 000049
【図１４】
2007235354.tif 000050
【図１５】
2007235354.tif 000051
【図１６】
2007235354.tif 000052