JP2023550434A - 改良型音響源測位法 - Google Patents

改良型音響源測位法 Download PDF

Info

Publication number
JP2023550434A
JP2023550434A JP2023530282A JP2023530282A JP2023550434A JP 2023550434 A JP2023550434 A JP 2023550434A JP 2023530282 A JP2023530282 A JP 2023530282A JP 2023530282 A JP2023530282 A JP 2023530282A JP 2023550434 A JP2023550434 A JP 2023550434A
Authority
JP
Japan
Prior art keywords
vector
velocity vector
time
equation
sound
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
JP2023530282A
Other languages
English (en)
Inventor
ジェローム・ダニエル
スルダン・キティク
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Orange SA
Original Assignee
France Telecom SA
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by France Telecom SA filed Critical France Telecom SA
Publication of JP2023550434A publication Critical patent/JP2023550434A/ja
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S5/00Position-fixing by co-ordinating two or more direction or position line determinations; Position-fixing by co-ordinating two or more distance determinations
    • G01S5/18Position-fixing by co-ordinating two or more direction or position line determinations; Position-fixing by co-ordinating two or more distance determinations using ultrasonic, sonic, or infrasonic waves
    • G01S5/22Position of source determined by co-ordinating a plurality of position lines defined by path-difference measurements
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/802Systems for determining direction or deviation from predetermined direction
    • G01S3/803Systems for determining direction or deviation from predetermined direction using amplitude comparison of signals derived from receiving transducers or transducer systems having differently-oriented directivity characteristics
    • G01S3/8034Systems for determining direction or deviation from predetermined direction using amplitude comparison of signals derived from receiving transducers or transducer systems having differently-oriented directivity characteristics wherein the signals are derived simultaneously
    • G01S3/8036Systems for determining direction or deviation from predetermined direction using amplitude comparison of signals derived from receiving transducers or transducer systems having differently-oriented directivity characteristics wherein the signals are derived simultaneously derived directly from separate directional systems
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S11/00Systems for determining distance or velocity not using reflection or reradiation
    • G01S11/14Systems for determining distance or velocity not using reflection or reradiation using ultrasonic, sonic, or infrasonic waves

Abstract

本発明は、少なくとも1つの壁がある空間内の少なくとも1つの音源を測位するために、たとえば、アンビソニック型の少なくとも1つのマイクロフォンにより取得された音響信号の処理に関する。時間周波数変換が、取得された信号に適用され、取得された信号に基づき、実部および虚部のある一般複素速度ベクトルV(f)が、周波数領域にて表され、このベクトルは、全方向成分W(f)以外の成分を分母に有する。特に、このベクトルは、*発信源とマイクロフォンとの直線経路である、第1のベクトルU0として表される、第1の音響経路と、*壁での反射音によって生じ、第2のベクトルU1として表される、少なくとも第2の音響経路との合成を特徴づけ、第2の経路は、マイクロフォンにて、直線経路を基準として、第1の遅延TAU1がある。遅延TAU1、第1のベクトルU0、および第2のベクトルU1の関数として、直線経路の方向(DoA)と、発信源からマイクロフォンまでの距離d0と、発信源から壁までの距離z0とのうちの少なくとも1つのパラメータが決定される。

Description

本発明は、音響源測位法に関し、具体的には、小型マイクロフォンシステム(たとえば、以下に用いる「アンビオフォニック(ambiophonic)」または「アンビソニック(ambisonic)」表現にて、音を受信可能なマイクロフォン)により、音響方向、すなわち「DoA」(到来方向)を推定するための測位法に関する。
適用例としては、たとえば、ビーム形成がある。ビーム形成は、音源の空間分離を含み、それは具体的には、音声認識(たとえば、音声対話による仮想アシスタント)を改良するためのものである。また、このような処理は、立体音響コーディング(主要信号を個別に符号化するための音声場面の事前分析)として利用可能である。すなわち、没入感のある音声コンテンツを、場合によっては視聴覚的に(芸術、電子音楽、映画、または他の目的に)、空間編集可能とする。また、これにより、遠隔会議での人の発話をトラッキングすることや、音声事象(関連映像付きまたは関連映像なし)を検出することも、可能となる。
アンビソニック(またはこれと等価な)符号化に関連する最新技術では、ほとんどの手法が、周波数解析によって生じる空間要素(通例、短時間フーリエ変換、すなわち「STFT」での処理により生じる時間周波数表現、またはフィルタバンクにより生じる狭帯域時間信号の表現)に基づくものである。
1次アンビソニック信号は、補遺にある式1にしたがって、ベクトル形式で収集される。便宜上、式1の符号化規約が、ここに提示されているが、他の規約への変換が適切に実施可能であるため、それに限定されるわけではない。このように、単位ベクトルU1により記述される方向(したがって、発信源の方向DoA)から到来して発信信号s1(t)を担う単一平面波と同等な場が、式2(補遺)により記述され得る。
実際に、信号は、周波数領域にてフレームごとに分析されて、式3(補遺)が得られ、単一波の場合には、式4の形式、および拡張としてN波については、式5の形式が得られる。
ある種の方法は、式6および式7に表すように、速度ベクトルV(f)または強度ベクトルI(f)の分析に基づく(前者は、全方向基準成分のパワーで正規化された後者の代替)。
複素周波数サンプルを利用する方法は、位置推定を、本質的には、このようなベクトルの実部に含まれる情報に基礎づけている(能動強度(active intensity)、および、位相場の傾きに直接関連した波動伝播の特徴に関連)。
虚部(エネルギー勾配と関連づけられた反応部(reactive part))については、定常音響現象の特徴と考えられる。
単一平面波の場合、速度ベクトルはV=U1からなることが、よくわかるであろう。
既知の方法(「DirAC」として知られている)は、サブバンドへとフィルタ処理された時間サンプルであって、実数で強度ベクトルも同様の時間サンプルか、あるいは、到来方向(または、より正確には、その反対方向)を指定するのに利用する強度ベクトルの実部だけである複素周波数サンプルで動作する。さらに、ベクトルのノルムと音場のエネルギーとの比に関連した、いわゆる「拡散性」係数の計算により、該当周波数で得られる情報が、その代わりに、方向成分の特徴となるか(この場合、ベクトルの方向が位置を決定する)、あるいは、「雰囲気(ambience)」(拡散反響音、および/または未区分の第2音源の混合)の特徴となるかどうか、決定することができるようになる。
以下に「VVM」として示す他の方法は、速度ベクトル、および実部の角度方向の統計値に基づき、実部および虚部と、それらのノルムとの比に関するある係数により重み付けられている。球場図法(spherical cartography)(2Dヒストグラム、たとえば、正距円筒図(equirectangular))は、全周波数サンプルについて、および一定数の時間フレームについて値を収集することにより、確立されている。したがって、推定値は、実質的には、最大確率に基づき、ある待ち時間(latency)にしたがう。
「共分散」法と称する他の種類の方法は、第1の拡張として提示されることもあり、空間成分の共分散行列(パワースペクトル密度行列、すなわち「PSD」と称することもある)を、周波数サブバンドにて算出することを含む。ここでまた、虚部は、完全に無視されることもある。なお、この行列の第1行(または第1列)は、空間成分がアンビソニック型である場合には、強度ベクトルと等価である。これらの手法の多くは、「部分空間」法、および、特に、多くの周波数サブバンド上で稼働するときや、高い空間解像度を利用するときに、負担が大きくなることもあるアルゴリズムを含む。
これらの「ベクトルに基づく」または「行列に基づく」方法は、一方で、測位可能な音響源成分、すなわち経路と関連づけられ、他方で、環境成分と関連づけられた、「方向」成分を、得ようとするものである。
このような方法の観測された限界のなかで、これらの方法は、単一の同時音響源の場合であっても、反射音を伴う直接音(音響源の方向を示す)の干渉により妨害される。ある程度の室内効果があると、たとえば、適正な推定値を得ることができないことが多く、かつ/または、推定値は偏りすぎることが多い。音響位置決定および捕捉装置を含む物体(たとえば、アンビソニックマイクロフォン)が、たとえばテーブル上、または壁の近くに配置されると(および/または、これが音響源の場合)、このような反射面は、総合的な角度バイアスを誘発しやすい。
確かに、位置特定には、一般に、同一の音響源と関連づけられた、反射音を伴う直接音の総合的な干渉により、バイアスがかかってしまう。速度ベクトルに基づく場合、主に考慮されるのは、速度ベクトルの実部であり、虚部は、通例、無視される(あるいは、少なくとも、十分に使われることはない)。音響的反射音は、迷惑なものと考えられ、推定問題に組み込まれることはない。したがって、誘発された特定の干渉構造を考慮することなく、無視されてモデル化されない成分が残る。
このように、上述の種類の適用例では、音響測位は、一般に、角度項のみについて推定される。さらに、単一の捕捉点に基づく距離評価を提案した効果的な手法はなさそうである(同時的、またはより一般には「小型」マイクロフォンシステムについて、特徴的と考えられる。すなわち、音響源からの距離に比べて小さなサイズ、通例、アンビソニックマイクロフォンであれば約10センチの容量内に格納)。
しかしながら、ある種の適用状況では、発信源からの距離についての追加情報が、方向に加えて(したがって、XYZでの立体位置特定)、必要である。これらは、たとえば、
- 立体的に取り込まれた実環境での仮想ナビゲーション(発信源の角度および強度の適切な修正は、当該物体とマイクロフォンとの相対XYZ移動に依存するため)、
- 発話者を特定するための発信源の測位(特に、接続されたスピーカまたは他の装置について)、
- 住居または産業環境での監視、警報装置、
その他。
文献FR1911723に提示された、特に有益な手法は、音の速度ベクトルを用いて、特に音の到来方向、その遅延(それによる、発信源までの距離)、および反響音の遅延を取得することにより、障壁の位置を特定する。このような実装例により、直接波と少なくとも1つの間接波(反射による)との干渉をモデル化するとともに、速度ベクトル全体(その虚部および実部)につき、モデルの表現を利用することが、可能となる。
しかしながら、この技術は、既に用いられているものの、さらに改良されるべきものでもある。
本発明は、この状況を改良するものである。
FR1911723
Jerome Daniel(2000)
少なくとも1つのマイクロフォンにより取得された音響信号を処理する方法が、提案され、
少なくとも1つの壁がある空間内の少なくとも1つの音源を測位するために、
- 時間周波数変換が、取得された信号に適用され、
- 取得された信号に基づき、一般速度ベクトルV'(f)が、周波数領域にて表され、全方向成分W(f)以外の基準成分D(f)が分母に現れる、速度ベクトルV(f)の式から推定され、その式は、実部および虚部のある複素式であり、一般速度V'(f)は、
*発信源とマイクロフォンとの直線経路である、第1のベクトルU0として表される、第1の音響経路と、
*壁での反射音によって生じ、第2のベクトルU1として表される、少なくとも第2の音響経路との合成を特徴づけ、
第2の経路は、マイクロフォンにて、直線経路を基準として、第1の遅延TAU1があり、
- 遅延TAU1、第1のベクトルU0、および第2のベクトルU1の関数として、
*直線経路の方向(DoA)と、
*発信源からマイクロフォンまでの距離d0と、
*発信源から壁までの距離z0とのうちの少なくとも1つのパラメータが決定される、方法。
上述の一般速度ベクトルV'(f)は、全方向性についての分母の成分の関数として一般に表現される、速度ベクトルV(f)から、このように構成される。一般速度ベクトルV'(f)で、上述の文献FR1911723における意味の「従来の」速度ベクトルV(f)を置き換えて、分母における全方向成分以外の「基準」成分が得られる。基準成分は、音の到来方向に向かって、実際に、より「高分離」となり得る。さらに、図6Aから図6Dおよび図7を参照して以下に示す例示的実施形態では、基準成分を算出可能とする音の到来方向が、たとえば、正確なDoAに向かって徐々に収束する反復処理手順における第1の反復処理中に、従来の速度ベクトルV(f)を用いることにより、第1の近似値として得られる。
そして、上述のパラメータDoA、d0、z0の決定は、速度ベクトルV(f)に代えて、より関連性の高い基準成分を伴う、このような一般速度ベクトルV'(f)を用いることにより、具体的に、より正確に、かつ/または、より精確になることが、本発明の意図の範囲内で、観測された。具体的には、本方法は、より安定しており、特に、音響的に強い反射音が、マイクロフォン近傍または能動音源近傍に配置された障壁に由来する状況にて、より安定している。
既述のパラメータDoA、d0、z0は、上述のように引用されたものである。なお、典型的な実施形態では、一般速度ベクトルの式により、特に、上述の遅延TAU1を決定することができる。
一実施形態では、本方法は、上述のように、複数回の反復処理であって、以前の回の反復処理で取得された直線経路の方向(DoA)の近似値に基づいて決定された基準成分D(f)を、分母に含んだ、一般速度ベクトルV'(f)が、少なくとも一部に用いられた、複数回の反復処理を含む。ほぼ全ての状況下で、これらの反復処理は、より精確な方向DoAに収束する。
そして、このような方法は、一般速度ベクトルV'(f)の代わりに、「従来の」速度ベクトルV(f)が用いられる、第1の反復処理を含む。文献FR1911723に記述されているように、速度ベクトルV(f)は、周波数領域にて表現され、分母には、全方向成分W(f)が現れている。そして、第1の反復処理が終わると、少なくとも、直線経路の方向(DoA)の第1の近似値を決定することができる。
このように、少なくとも、第1の反復処理に続く第2の反復処理に、一般速度ベクトルV'(f)が、用いられ、分母の全方向成分W(f)が基準成分D(f)に置き換えられた、速度ベクトルV(f)の式から推定され、後者は、全方向成分W(f)よりも、空間的に高分離である。
たとえば、一実施形態では、基準成分D(f)は、直線経路の方向(DoA)の上述の第1の近似値に対応する方向にて、より高分離である。
反復処理は、所定の基準にしたがって収束に達するまで、繰り返され得る。特に、これは、少なくとも、マイクロフォンと発信源との間の音伝播環境における障害物(すなわち、上述の「障壁」)での第1の反射音を、適正な確実性で特定するために、因果関係基準であり得る。
具体的な一実施形態では、各反復処理は、
- 直線経路(DoA)に沿った音の到来に関連したピークに加えて、各々が少なくとも1つの壁での反射音に関連した一連のピークを、時間領域にて取得するために、周波数から時間への逆変換も、一般速度ベクトルV'(f)の式に対して適用され、
- 時間座標が、直線経路のもの未満であり、振幅が、選択閾値を超えた信号が、一連のピークに現れた場合には、新規の反復処理が実行され(場合により、適応可能)、
因果関係基準は、信号の振幅が、閾値未満である場合に満たされる。このような一連のピークを取得することは、通例、補遺にある式B4=39bに示す形式につながり得るものであり、図2を参照してさらに後述するが、ここで、一般速度ベクトルV'(f)に適用されることは当然である。
上述の方法の反復処理は、たとえば、
- 上述の信号の振幅が、選択閾値未満である、第1の場合、
- 反復処理の繰り返しが、この信号の振幅の大幅な減少につながらない、第2の場合に、
終了され得る。
一例示的実施形態では、第2の場合に続いて、以下の各ステップを実施し、取得された信号は、サンプルの連続フレームの形式として得られ、各ステップは、
- 各フレームに対し、そのフレーム内に音の開始があることについてのスコアが、推定されるステップ(たとえば、補遺の式53などの式)と、
- 閾値より高スコアのフレームが、該フレーム内で取得された音響信号を処理するために、選択されるステップとである。
確かに、DoAの解へ向けての収束が、直接的に第1の反射音が発生する障壁に近接していることにより、容易でない場合、音の開始に対するこれらの障壁での直接的反射音を得ることが好ましい(音の放出が始まった場合)。
「従来の」および「一般」速度ベクトルのそれぞれの式を考慮して、取得信号がアンビソニックマイクロフォンにより得られる実施形態では、「従来の」速度ベクトルV(f)は、周波数領域にて、以下の種類の形式の1次アンビソニック成分により表され得る。
V(f)=1/W(f)[X(f),Y(f),Z(f)]T,
W(f)は、全方向成分であり、
周波数領域にて、1次アンビソニック成分により表される、一般速度ベクトルV'(f)は、以下の種類の形式で提示される。
V(f)=1/D(f)[X(f),Y(f),Z(f)]T,
D(f)は、上述の全方向成分以外の基準成分である。
ここで考慮する次数は、速度ベクトルの成分を3次元座標系で表すことが可能な1次であるが、他の実装例も可能であり、特に、高次のアンビソニック次数にて可能である。
一実施形態では、第1のベクトルU0と同等の直線経路の方向の推定値は、周波数領域にて表された一般速度ベクトルV'(f)の実部の周波数の集合における平均値から、決定され得る(一般速度ベクトルV'(f)に対して、ここに当然に適用される式24の形式にしたがう)。
このように、周波数領域での速度ベクトルの式により、ベクトルU0の推定値を提供可能である。
ただし、さらに改良された実施形態では、
- 周波数から時間への逆変換が、一般速度ベクトルに対して、時間領域V'(t)内にて表すために適用され、
- 直線経路の時間後に、少なくとも1つの最小値が、一般速度ベクトルV'(t)maxの式にて、時間の関数として検索され、
- そこから、最大値V'(t)maxとなる時間に対応する第1の遅延TAU1が、推論される。
さらに、この実施形態では、
- 第2のベクトルU1は、ベクトルV1を以下のように定義するために、時間指標t=0、TAU1、および2xTAU1にて記録された、正規化速度ベクトルV'の値の関数として推定され、
V1=V'(TAU1)-((V'(TAU1)V'(2.TAU1))/||V'(TAU1)||2)V'(0)
そして、ベクトルU1は、U1=V1/||V1||。
そして、
- 第1のベクトルU0および第2のベクトルU1のそれぞれの角度PHI0およびPHI1は、壁について、以下のように決定可能であり、
PHI0=arcsin(U0.nR)およびPHI1=arcsin(U1.nR)、ここで、nRは、単位ベクトルであり、壁に垂直であり、
- 発信源とマイクロフォンとの距離d0は、第1の遅延TAU1の関数として、以下の種類の関係式
d0 = (TAU1xC)/((cosPHI0/cosPHI1)-1)、ここでCは音速、により決定可能である。
さらに、発信源から壁までの距離z0は、以下の種類の関係式
z0=d0(sinPHI0-sinPHI1)/2
により決定可能である。
このように、発信源の測位に関する全パラメータを決定可能であり(たとえば、図1におけるもの)、ここでは、単一の壁がある場合であるが、このモデルは、いくつかの壁がある場合へと一般化可能である。
このように、空間は、複数の壁を備える実施形態では、
- 周波数から時間への逆変換が、一般速度ベクトルに対して、時間領域V'(t)内にて、一連のピークの形式で表すために適用され(補遺の式39bに対する第1の手法に対応する形式)、
- 一連のピークにて、複数の壁のうちの1つの壁での反射音に関連したピークが、特定され、各々の特定されたピークは、直線経路を基準としたときの、対応する壁nでの反射音によって生じる音響経路の第1の遅延TAUnの関数である時間座標を有し、
- 壁nでの反射音よって生じる音響経路を表す第1のベクトルU0および各第2のベクトルUnの各第1の遅延TAUnの関数として、少なくとも1つのパラメータが、
*直線経路の方向(DoA)と、
*発信源からマイクロフォンまでの距離d0と、
*発信源から壁nまでの少なくとも1つの距離znとのうちから決定される。
「従来の」速度ベクトルに適用される例であるものの、一般速度ベクトルにも適応可能である、図5Bに示すように、逆変換(周波数から時間へ)後、速度ベクトルの式(従来および一般)は、一連のピークを示し、一連のピークは、説明のために図2にも示され、ここで、最大値は、直線経路と、壁での少なくとも1つの反射音により生じる経路との間の上述の遅延(TAU1、2TAU1など、TAU2、2TAU2など)の複数の値、およびこれらの遅延の組み合わせ(TAU1+TAU2、2TAU1+TAU2、TAU1+2TAU2など)に達する。
そして、これらのピークは、特に、壁nでの少なくとも1つの反射音に関するピークを特定するために利用可能であり、これにより、この壁nと関連づけられた遅延TAUnについて、複数の時間座標(x1、x2、x3など)がある。
様々な遅延が組み合わされることにより、単一の遅延(TAU1、TAU2、TAU3など)、および対応する壁の存在の特定が難しくなり得るため、最小の正の時間座標を有するピークの第1の部分が、この部分にて、壁での反射音と各々関連づけられたピークを特定するために、予め選択され得る(第1の遅延後に現れ得る、様々な遅延の組み合わせTAU1+TAU2、2TAU1+TAU2、TAU1+2TAU2など、とならないように)。しかしながら、このような実装例では、上述の因果関係基準が満たされていることが、想定されている(あるいは、「第2の」ピークも、負の乗数を伴う遅延の組み合わせで、「正の」遅延との組み合わせが、最終的に小さな正の時間座標となるものとしても、取得可能である)。
このように、壁nでの反射音に関連したピークが、この壁nと関連づけられた遅延TUnの倍数である時間座標を有することもある、理想的な場合には、最小の正の時間座標を有するピークの第1部分を、壁での単一反射音と関連づけられた当該部分の各ピークを特定するために、予め選択可能である。
マイクロフォンにより取得された信号が、一連のサンプルの形式となる実施形態では、時間的に減少する指数関数的変動を伴う重み付けウィンドウを、これらのサンプルに適用することが、より一般に可能である(図5Aを参照して、さらに後述)。
さらに、このウィンドウは、音の開始における最初の部分(あるいは、音の開始の前)に配置可能である。これにより、複数の反射音による困難性が、避けられる。
このような重み付けウィンドウの適用例により、時間領域での速度ベクトルの式を利用することによる、パラメータU0、d0などの、よりバイアスの少ない第1の推定値を得ることが可能であり、特に、「従来の」速度ベクトルが関係する場合、たとえば、本方法の第1の反復処理の状況にて、可能となる。確かに、反射音の累積的な強さが、直接音よりも大きくなる、ある種の状況では、上述のパラメータの推定値に、バイアスがかかり得る。これらの状況は、速度ベクトルの時間表現にて、ピークが、負の時間座標で観測される場合(図5Bの上部の曲線)、検出可能である。上述の種類の重み付けウィンドウの適用例により、これらのピークが、図5Bの曲線の底部に示すように、正の座標に戻ることが可能となり、よりバイアスの少ない推定値が得られる。
しかしながら、この実装例は、「従来の」速度ベクトルの代わりに一般速度ベクトルを用いることにより、この種の状況にて、パラメータU0、d0などの、ほぼバイアスのない推定値が既に可能であるなら、任意的となる。しかしながら、このような処理は、たとえば、「従来の」速度ベクトルを用いる本方法の第1の反復処理にて、あるいは、反復処理が非収束的な上述の第2の場合にも、実施され得る。
一実施形態では、さらに、重み付けq(f)が、周波数領域での速度ベクトル(一般または従来)に対し、以下の種類の式にしたがって、反復して、各々が周波数帯域fと関連づけられて、適用されてもよい。
q(f)=exp(-|Im(V(f)).m|/(||Im(V(f))||)
ここで、Im(V(f))は、速度ベクトル(従来または一般で、ここでは、単に「V(f)」と記す)の虚部であり、mは、ベクトルU0により定義される平面に垂直な単位ベクトルであり、壁の法線である(通例、以下にさらに詳述する図1のZ軸)。
このような実施形態により、上述のパラメータを決定するのに最も有効な周波数帯域が、選択される。
また、本発明は、上述の方法を実施するように構成された処理回路を備えた、音響信号処理装置に関する。
説明のために、図4は、このような処理回路を模式的に示し、この回路は、
- マイクロフォン(たとえば、アンビソニック的に信号を生成するためのいくつかの圧電板を備えてもよい)によって取得された信号SIGを受信するための入力インターフェースINと、
- ワークメモリMEMと協働して信号を処理し、特に、各々の値が出力インターフェースOUTによりもたらされ得る、所望のパラメータd0、U0などを導出するために、一般速度ベクトルの式を展開するプロセッサPROCとを、備えている。
このような装置は、立体環境にて音源を測位するモジュールの形態をとり得るものであり、このモジュールは、マイクロフォン(音響アンテナまたは他の種類のもの)と接続している。逆に、拡張現実で、仮想空間(1つまたは複数の壁を含む)における発信源の所与の位置に基づいて、音を生成するエンジンであってもよい。
また、本発明は、処理回路のプロセッサにより実行された場合に、上述の方法を実施する命令を含む、コンピュータプログラムに関する。
たとえば、図3A、図3Bおよび図7は、このようなプログラムのアルゴリズムとして考えられるフローチャート例を示す。
他の側面では、このようなプログラムが格納された非一時的コンピュータ可読媒体が、提供されている。
それによる詳細な説明では、一般速度ベクトルおよび「従来の」一般ベクトルが、「速度ベクトル」という同じ用語が区別されることなく、同じV表記(V(f)、V(t))で、特に、補遺に示す式において、指定されている。具体的に、一般速度ベクトルに関する場合、速度は、その用語で明示的に指定され、V'(V'(f)、V'(t))と表記される。図5Bまでの説明の最初の部分では、文献FR1911723に用いられている形式の基礎にあるとともに、ここで、補遺にある式に含まれている原理が、確認される。
より一般的には、他の特徴、詳細、および利点が、以下の詳細な説明を読み、以下の添付の図面を分析することにより、明らかになるであろう。
図1に、説明のために、一実施形態による音源測位に伴う様々なパラメータを示す。 図2に、説明のために、周波数から時間への逆変換(「IDFT」)後の速度ベクトルの時間表現として提示する様々な一連のピークを示す。 図3Aに、関連するパラメータU0、d0などを決定するアルゴリズム処理における最初の数ステップを示し、 図3Bに、図3Aの処理ステップの続きを示す。 図4に、一実施形態による本発明の意図の範囲での装置を模式的に示す。 図5Aに、一実施形態による、経時的に指数関数的に減少する取得信号のサンプルのための重み付けウィンドウを示す。 図5Bでは、以下の、速度ベクトルのIDFT後の時間表現を比較している。- 重み付けウィンドウによるサンプルの事前処理なし(上の曲線)、- ウィンドウ処理あり(下の曲線)。 図6Aに、図7を参照して後述する方法の反復処理が繰り返されることによる、一般速度ベクトルV'(t)の時間表現にある関連するピークの概要を示す。 図6Bに、図7を参照して後述する方法の反復処理が繰り返されることによる、一般速度ベクトルV'(t)の時間表現にある関連するピークの概要を示す。 図6Cに、図7を参照して後述する方法の反復処理が繰り返されることによる、一般速度ベクトルV'(t)の時間表現にある関連するピークの概要を示す。 図6Dに、高度に模式的に、説明のために、本方法を数回反復したときに、一般速度ベクトルV'(f)の式の分母に現れる基準成分D(f)の形式を示し、 図7に、例としてここに示す実施形態による、本発明の意図の範囲内での反復的方法のステップを、模式的に示す。
速度ベクトルは、周知の方法により、算出可能である。しかしながら、ある種の具体的パラメータ設定が、得られる最終結果を改良するために、推奨されてもよい。
通常、アンビソニック信号の周波数スペクトルB(f)は、最初に、短時間フーリエ変換(すなわちSTFT)により、一連の時間フレームb(t)について、一般的には重なって(たとえば、重複/追加で)、得られることが典型的である。ここで、アンビソニック成分の次数は、4つの成分について、m=1であってもよい(なお、一般性を失うことなく、計算は、高次にも適応可能)。
そして、各時間フレームについて、次に、速度ベクトルが、全周波数サンプルについて、方向成分X(f)、Y(f)およびZ(f)の
- 従来の速度ベクトルでは、全方向成分W(f)(補遺の式6)に対する比として、または
- 一般速度ベクトルでは、基準成分D(f)に対する比として、補遺の式6にて、W(f)がD(f)に置き換えられて、算出される。
以下に説明する加重和による時間的平滑化または統合を導入した実施形態を、考えることができる。
このような比(X(f)/W(f)、Y(f)/W(f)、Z(f)/W(f);X(f)/D(f)、Y(f)/D(f)、Z(f)/D(f))において、音響信号のスペクトル成分が、有益な周波数(たとえば、広い周波数帯域で)の実質的な量を実際に強化する場合、元の信号の特徴は、音響チャネルの特徴を強調するように、実質的に除去される。
上述の適用例では、安定した音響環境(反射壁および物体、場合によっては回折性など)での信号s(t)を発する(少なくともいくつかの連続したフレームにわたって、位置において、および放射において)安定した特徴の音響源の状況を考えることができるので、実際の「部屋」の中でなくとも、通常「室内効果」と称するものの原因となる。これらの信号は、アンビソニックマイクロフォンにより、捕捉される。アンビソニック信号b(t)は、直接および間接経路に沿った信号s(t)の様々なバージョンの空間符号化の組み合わせである、いわゆる「音響チャネル効果」に由来する。これは、結果的に、各チャネル(すなわち次元)が、補遺の式8として表されるアンビソニック成分と関連づけられた空間インパルス応答h(t)による信号の畳み込みとなる。
このインパルス応答は、SRIR、すなわち、「空間室内インパルス応答(Spatial Room Impulse Response)」称するものであり、一般に、以下の一連の時間的ピークとして表される。
- 直接音に対応する時刻t=TAU0(伝播時)で測位された第1のピーク、
- 第1の反射音に対応する時刻t=TAU1での第2のピーク、
などである。
このように、これらのピークにて、これらの波面の到来方向を、第1の近似値での式9-1で与えられるベクトルunの式で、読取可能である必要がある。実際には、空間インパルス応答は、未知のデータであるが、ここで、アンビソニック信号b(t)に基づいて算出された速度ベクトルを通じて、どのように、その特徴のいくつかに、間接的に至ることになるのか、説明する。
これを強調するため、まず、選択された観測時間にわたる、インパルス応答h(t)、送信信号s(t)、およびアンビソニック信号b(t)(式9-2)の関係を、説明する。正確には、この式は、測定ノイズが無く、当該時間にわたって、直接的または間接的に信号が補足される他の音響源も無いことを、想定している。したがって、発信源からの全ての直接および間接信号は、この時間にわたって補足される。
この時間の全体にわたり、フーリエ変換を実施して、結果として得られた速度ベクトルには、空間インパルス応答の一意的な特徴があることが、示されている。このいわゆるLT変換(STFTよりも「長期(longer term)」であるため)は、式10にしたがって、b(t)、s(t)およびh(t)を、B(f)、S(f)およびH(f)へと変換する。この時間的な補助は、いくつかの連続した信号フレームにわたって延びる時間ウィンドウに、対応し得る。
そして、式11にて算出された速度ベクトルの式は、周波数領域にて算出された畳み込みの式から、導かれる。当該時間にわたる各周波数fについて、非ゼロ・エネルギー(実際には、検出可能)であるとした場合、この式11は、音響チャネルで特徴的となり(換言すると、室内効果)、送信信号のものではなくなる。
実際には、上述したように、よくある手法は、観測信号が、式9の畳み込み積に完全にかつそれだけに由来することを、原則として確認することができない時間ウィンドウに対して、各短時間フーリエ変換が適用される、フレームごとの時間周波数解析を実行するためのものである。これは、厳密には、速度ベクトルが、音響チャネル(式11の右辺におけるものとして)を特徴づけるだけの形式では記述され得ないことを、意味している。しかしながら、ここで、この記述の状況に沿って、できる限り近似がとられるとともに(以下に導入する式20)、以下に示す短時間解析の利点が利用される。
後続のステップでは、発信源から発せられて、マイクロフォンにより受信された信号の直線経路を特徴づけている、一連のエネルギーピークを取得して、1つまたは複数の壁での第1の反射音を、これらの反射音が特定可能である限りにおいて取得する。そして、空間インパルス応答の開始を特徴づけるものに集中することができる。すなわち、最初に、直接音の方向が推論される第1の時間ピーク、および、場合によっては、第1の反射音を特徴づける後続の時間ピークに集中することができる。
このようにするために、直接音と、少なくとも1つの反射音との干渉の効果を、音源の位置を特定するための関連パラメータを推定するように、複素速度ベクトルの式で調べる。
式12に示すインパルス応答の開始について、N個の反射(n=1,…,N)を組み合わせた直線経路(n=0)の簡易なモデルを、導入する。この式において、gn、TAUnおよびunは、それぞれ、減衰、遅延、および、マイクロフォンシステムに到達するまでの(n回反射)インデックスnの波の到来方向である。以下に、簡潔にするために、しかし、一般性に制限を課すことなく、式13の項をn=0に設定することになる、直接音に対しての遅延および減衰を考慮する。
対応する周波数の式は、直接音については、gamma0=1の特別な場合として、式14で与えられる。変数gammanは、0より大きい任意のnについて、周波数fの関数であることは当然である。
アンビソニック場の周波数表現は、後半を無視するならば、式16で与えられることになる。
そして、短期速度ベクトル(short-term velocity vector)は、式17で表され、あるいは、分母の成分が(約)ゼロとなるときに、無限の値を(ほぼ)回避するように、非ゼロイプシロン項で正規化されたものとして、式18で表される。式17または式18では、成分W(従来の速度ベクトルに特徴的)は、一般速度ベクトルで表すために、基準成分Dで置き換えられ得る。実際には、一般的な場合、分母では、DによりWが置き換えられ、従来の速度ベクトルVの式は、D=Wである特別な場合に対応している。しかしながら、ここでの説明を簡潔にするために、D=Wの特別な場合に関連した表記は、従来の速度ベクトルについての、補遺の最初の式に提示されているが、DによりWを置き換えることに留意することで、簡単に一般速度ベクトルへと置き換えることが可能である。
短時間解析により、経時的に、かつ、元の信号の動的発展にしたがい、空間インパルス応答内の波面のサブミックスに特徴的な、周波数フィンガープリント(以下、「FDVV」と記す)を、観測することが可能となる。所与の観測について、特徴的なサブミックス(smxは「サブミックス(submix)である」)は、時間および周波数領域で、式19によりモデル化される。
以下に説明する手法では、周波数フィンガープリントFDVVを、式20に表す近似(一般に、空間状況以外では、厳密に等しいわけではない)による、サブミックスHsmxの暗黙的モデル(implicit model)に由来するものとして、特徴づけるようにしている。実際に従来の速度ベクトルにつき、式20での表現が、ここに与えられ、H_Wを、以下の行列DΘによる行列ベクトルHで置き換えることにより、一般速度ベクトルに適用可能である。
特に、信号開始時には、暗黙的モデルhsmx(t)は、空間インパルス応答の開始hearly(t)と、少なくとも、相対波面方向および遅延について、似ているようにみえる。暗黙的相対ゲインパラメータgnは、時間ウィンドウおよび信号の動的特徴に影響され、インパルス応答のものと、必ずしも一致するようにはならない。ここで、直接波(DoAを提供)、および1つまたはいくつかの初期反射音に、主に注目して、観測が特徴的である状況について、主に検討される。
特に、説明のために、従来の速度ベクトルを周波数領域で推定しつつ、単一の反射音のみを考慮する処理の一例を、以下に提示し、一般速度ベクトルの場合については、以後説明する。ここで(実質的には直接音と第1の反射音との間の)単一の干渉の場合について検討し、具体的な空間周波数構造を強調し、速度ベクトルの実部だけでなく虚部をも考慮することにより、所望のパラメータをどのように決定するのかを示す。確かに、アンビソニック場は、式21にしたがって記述され、これにより、式22にしたがって、速度ベクトルを推論する。この式から、実部および虚部は、周波数が式23に示す当該音響スペクトルを進む場合、立体空間中の平行な部分を進行する(それぞれ、アフィンおよび線形)ことになる。アフィン部分(実部)は、直接波および間接波をそれぞれ示す単位ベクトルU0およびU1を含む直線上にあり、この2つの部分は、これら2つのベクトルの中間面に垂直である(このように、ベクトルの虚部は、線形部分上にあるため、そのままである)。さらに、統計計算により、波同士の位相差が一様に分布すること(それにより、周波数の代表的な掃引)を想定し、速度ベクトルの実部の平均は、式24で表すようなベクトルU0に等しく、最大確率は、式25に示すような波の振幅で加重されたU0およびU1の平均である。したがって、最大確率に基づくDoAの検出は、規則的な角度バイアスにより妨げられ、直接音とその方向の中間的方向を与える。式23は、周波数を2つの波の間の遅延TAU1の逆数に周期的に等しくして、この空間の掃引がなされたことを示す。したがって、このような空間周波数構造が観測可能な場合、方向U0およびU1、ならびに観測からの遅延TAU1を、抽出することができる。さらに、これらのパラメータを時間領域で推定する他の方法を、以下に示す(図2との関連で説明する)。
マイクロフォンの基準系との関連で、反射面の姿勢について予想することで、U0、U1、TAU1の推定から、マイクロフォンに対する発信源の絶対距離d、および場合によっては両者の高度の情報を、推論可能である。実際に、図1に示すように、発信源S0からマイクロフォンMへの距離をd0とし、d1を、反射面Rについての鏡像S1とすると、面Rは、ベクトルU0およびU1がなす平面に垂直となる。これら3点(M、S0、S1)は、面Rと直交する同一平面にある。反射面の姿勢(すなわち傾斜)を規定するために決定すべきパラメータは、規定されたままである。床または天井による反射の場合(U1が、床または天井を指しているため、このように検出)、アンビソニックマイクロフォンの座標系のXY平面に水平かつ平行という前提を、利用することができる。そして、距離d0およびd1は、関係式26により関連づけられ、これにより、さらに、マイクロフォンMから軸(S0、S1)までの距離が、直接的に与えられ、PHI0およびPHI1は、ベクトルU0およびU1のそれぞれの仰角である。
また、直接音に対する反射音の遅延TAU1の推定値が得られ、これにより、距離と距離との他の関係式27が、利用できるようになる。これらの差が、係数cを音速として、音響経路の遅延を示すためである。
d1をd0の関数として表すことにより、この最後の量だけが未知となり、式28にしたがって推定可能である。発信源から反射面までの距離、すなわち、式29にしたがって、その高さ、あるいは地面を基準とした高度z0を取得するとともに、式30にてマイクロフォンについても取得する。
様々なパラメータU0、U1、PHI0、PHI1、d1、d0などを、床での反射音の例で、図1に示す。天井での反射音について、同様のパラメータが推論され得ることは当然である。同様に、類似のパラメータが、マイクロフォンの座標系を基準とした姿勢が既知の他の反射面Rでの反射音について、推論可能であり、姿勢は、法線nR(面Rに直交する単位ベクトル)により特徴づけられている。反射面Rに対して、角度PHI0およびPHI1を、一般にPHI0=arcsin(U0.nR)およびPHI1=arcsin(U1.nR)として改めて定義すれば十分である。このように、各反射音の場合と関連づけられたベクトルU1により、拡張現実またはロボット工学の用途に、音響検知による位置推定のため、これらの障害物のそれぞれの位置を決定することができる。
反射面の姿勢nRが最初は未知である場合、別々の回での観測により、同一の反射面による反射音が検出される少なくとも2つの発信源位置と関連づけられた波面パラメータの推定値があれば、その姿勢を完全に推定可能である。このように、第1の組のパラメータ(U0、U1、TAU1)、および少なくとも第2の組(U0'、U1'、TAU1')が得られる。U0およびU1が、面Rと直交する面を定義すると、それらのベクトル積が、この面Rの軸を定義し、同じことが、U0'およびU'1から引いたベクトル積にも適用される。
これらのベクトル積(同一線上にはない)は、面Rの姿勢を定義する。
しかしながら、干渉波が2つのみ(直接音および反射波)であるモデルの一つの制約は、障壁での様々な第1の反射音を区別することが困難となり得ることである。さらに、追加の反射音が導入されると、速度ベクトルの空間周波数的挙動は、すぐにより複雑となる。確かに、実部および虚部の経路が、非自明に、いくつかの軸に沿って、
- 直接波の平行面で、および2つの反射音について、
- または、一般に、空間全体で、
組み合わされる。
いくつかの反射面が考慮される場合、これらの複合空間周波数分布により、モデルのパラメータを決定することが、非常に複雑となる。
この問題に対する一解決策は、ある振幅(一時的、立ち上がり信号)の開始時に、より単純な音響ミックスの出現を特定する機会を得るために、時間的により高分離な(すなわち、より短い時間ウィンドウの)時間周波数解析を行うことである。すなわち、当該フレームにあるミックスでの直接音に干渉する反射音の数を減らすことである。しかしながら、ある種の状況では、一連の反射音と関連づけられた遅延が、直接音への干渉における第1の反射音の効果を隔離するのには、近すぎることもある。
したがって、複数の干渉の影響を簡単に分離して特徴づけることを可能にする、以下の処理が提案される。第1のステップは、速度ベクトルのフィンガープリントを、時間領域へと(すなわち「時間領域速度ベクトル(Time-Domain Velocity Vector)」である「TDVV」)、式31に示す逆フーリエ変換により変換することを含んでいる。これには、ある軸と関連づけられ速度ベクトルが複雑な遍歴として現れる周波数の循環性の効果をより疎らな、したがってより簡単に解析可能なデータに集約する効果がある。確かに、このような変換により、ピーク列が、規則的な時間間隔で現れるようになり、最も顕著なピークが、簡単に検出および抽出可能となる(たとえば、図5B参照)。
注目すべき特性として、構成により(逆フーリエ変換により)、t=0でのベクトルが、周波数領域にて速度ベクトルの平均(正の周波数の半スペクトルのみを考慮する場合には、実部の平均)と等しくなることがある。このような観測は、主DoA U0を推定することに関係している。
2つの干渉波(直接音および反射音)についての速度ベクトルの周波数モデルから始めて、式32のテイラー展開により、分母を、有用となるように、再度定式化することができる。式32により与えられるxおよびgamma1についての条件では、(従来の)速度ベクトル式33に至り、反射音の振幅が直接音よりも小さいという条件(g1<g0=1、これは、一般に、音の開始時)では、この式の逆フーリエ変換が収束し、式34として表すように定式化され、ここで、第1のピークは、t=0で特定され、U0(直接音の方向)が得られ、そして、反射音の直接音との干渉のピーク列の特徴が得られる。
これらのピークは、遅延TAU1の複数の時刻t=kTAU1(非ゼロ整数k>0)に位置し、ノルムにおける振幅が指数関数的に減少する(ゲインg1にしたがう)。従来の速度ベクトルを用いることにより、これらは、U0-U1の差と同一直線上にある方向と、全て関連づけられ、このため、これら2つのベクトルの中間面に直交し、方向が交互になる(符号)。速度ベクトルを時間領域へと変換することの利点は、所望のパラメータが、疎らで、ほぼ直接的に表現されていることである(図2)。
このように、主DoA U0に加えて、以下を決定することができるようになる。
- いくつかの区別された壁についてのこともある、遅延TAU1、
- 次に、単位ベクトルnへと正規化された、U0-U1と同一直線上にあるベクトルであり、たとえば式41で、以下のために、利用可能(従来の速度ベクトルに与えられる)、
- U1を、中間面に対してU0に対称なものとして推測、
- 任意であるが、減衰パラメータg1(時間周波数解析パラメータにより、特に、解析ウィンドウの形態により、そして、観測された音響事象を基準とした、その時間配置により、修正可能である。したがって、このパラメータの推定は、ここで検討する適用状況では、有用性がより低い)。
以下の時間ピークを複数観測することにより、これらが、同じく一連のものとして、実質的に一貫しているかどうか(複数の遅延TAU1、複数の遅延TAU2など)、それにより、同じ干渉の特徴を調べることが可能となり、そうでなければ、たとえば、複数の反射音があることを特定することが、必要となる。
以下に、「好適な条件」の場合を検討し、式35のガンマのNについての和が1未満であれば、反射音がN個の場合、テイラー展開を適用して、式35にしたがい、従来の速度ベクトルを得る。テイラー級数は、最初の式の分母を、変換するものであり、式36の多項式の法則を用いて、書き換え可能であり、従来の速度ベクトルVモデルの式を、いくつかの級数の和として、式37のSC項で表す「交差級数(cross series)」で、認識することができる。一般速度ベクトルV'について、補遺の最後方の式B2にある少しだけ異なる式を発見した。この式は、式35bとしても表され、それは、式35に対応するものの、一般速度ベクトルについて与えられているからである。さらに、補遺の最後の方で、一般速度ベクトルV'に特有の式が与えられており、式番号の後に「b」を示すことにより、従来の速度ベクトルに特有なものとして上述した式との対応がとられている(式xxb)。
従来の速度ベクトルについての式38、および任意の周波数f(一般速度ベクトルについての式B3=38b)に対する条件下で、遅延SARCが組み合わされた、以下の時間級数の式39(一般速度ベクトルについての式B4=39b)が、逆フーリエ変換によって推論される。U0(直接音の方向)を与えるt=0での第1のピークが特定されて、各反射音につき、直接音に対する反射音の干渉に特徴的なピーク列が、特定される。たとえば、図2では、これらのピークは、連続した正の時間座標TAU、2TAU、3TAUなどに位置する。これらは、壁での反射音と直線経路との間の遅延TAUの倍数である。
そして、いくつかの壁でのいくつかの反射音と、直接音との干渉の一連の特徴が(より大きな時間座標で)現れ、ここで、遅延は、様々な遅延の(正の倍数での)別の組み合わせである。
実際に、図2に、このような一連のものを、直接音に干渉する2つの反射音という単純な場合で示す。各マーカー(円、交差、ダイアモンド)は、ベクトルU0、U1、U2(直接音、第1の反射音、および第2の反射音のそれぞれの特徴)による、時間座標の関数としての時間フィンガープリントTDVVに対する寄与を縦座標にて示す。このように、直接音の受信は、円で示される、時刻がゼロで振幅が1の第1のピークにより特徴づけられることがわかる。直線経路に対する第1の反射音の干渉(遅延TAU1)により、TAU1、2xTAU1、3xTAU1などでの第1のピーク列が生じ、ここで、一端に交差(cross)が印され、他端に円が印される(上と下)。直線経路に対する第2の反射音の干渉(遅延TAU2)により、TAU2、2xTAU2、3xTAU2などでの第2のピーク列が生じ、ここで、一端にダイアモンドが印され、他端に円が印される。そして、「交差列」、すなわち、反射音と反射音の干渉(第1の遅延:TAU1+TAU2、そして、2TAU1+TAU2、さらにTAU1+2TAU2など)の要素となる。交差列は、一般の場合には式が得られるものの、書けば長くなるものであり、簡潔にするためにここでは明示的に書くことはしないが、これは、特に、ここに提示する処理での関連パラメータの推定に利用するのに必要ではないためである。
以下に、パラメータを連続的に推定することによる、時間フィンガープリントの解析について説明する。
算出した時系列に基づくモデルパラメータの推定は、既述の単一反射音の場合と同様になされる。まず、最も一般的な状況(後に扱う特殊な場合を除く)を検討すると、これは、遅延が「重複」していない好適な場合に対応し、上述の列は、時間的な一致がない。すなわち、特定可能なピークが、そのうちの1つにのみに属す。したがって、時間ピークを、t=0からの遅延が増大するとして示すことにより、遅延TAUnewを伴って検出された新たなピークは、既に特定した列のものであるか、または、新たな列の開始を規定するものである。確かに、既に特定された反射音の特徴がある遅延の組を考慮すると、正であるか一部ゼロである整数kが、式40にしたがって、TAUnewをもたらすのであれば、第1の場合が検出される。それ以外では、第2の場合となり、特定した反射音の組は、新たな遅延TAUN+1を導入することにより、増加する。これは、単一反射音の場合について説明したのと同様に、推定可能な方向と関連づけられている。
実際には、多くの時間ピークの説明を試みる必要はないであろう。観測した第1のピークに限定するのは、特に、振幅(すなわち、絶対値の大きさ)が、以下のものより大きいので、最も容易に検出可能であるためである。このように、遅延が公倍数であるものの、高い(または低くない)ランクKi;Kjである状況が、振幅の関数として、上述の処理により分析可能である。
暗黙的ゲインgn(n>0)の絶対値の和が1未満である限り(従来の速度ベクトルについて式38)、逆フーリエ変換(式31)により、時間が正で展開する単一方向時間フィンガープリントが得られる。
一方、暗黙的ゲインgn(n>0)の絶対値の和が1より大きい場合、逆フーリエ変換により、「双方向」時間フィンガープリントTDVVが得られ、列は、一般に正の時間および負の時間の双方へ向かって展開する(説明のための図5Bの上の曲線)。たとえば、直接波が、1つまたは複数の障壁での反射音に由来する波の振幅の和未満の振幅となる場合、このように1つまたは複数の反射ゲインが1より大きくなる状況が、生じ得る。この「好適でない場合」、時刻ゼロでの主要なピークは、もはやベクトルu0に厳密には対応しないものの、後者と、反射音の方向を示すベクトルのある程度大きな割合との組み合わせに、対応している。このことが、測位バイアスにつながる(「推定DoA」において)。他の兆候は、主要なピークが、一般に1とは異なるノルムを有し、さらに多くの場合、1未満となることである。
本発明は、特に、このような状況において、安定した方法を提案する。通常の全方向成分Wの代わりに分母に現れるD成分でのDoAに向けて、空間的分離を与えることにより、速度ベクトルの式を調整することが、提案される。
基準成分Dに直接的手法を提供することにより、インデックスnの各反射音と関連づけられた相対的減衰が、係数BETAnによる分母に対して割り当てられ、同時に、全体の係数として、NU0が算出され(補遺の式B1)、これが、一般速度ベクトルの場合として既に導入された、N回反射音のモデルとして、式B2=35bで与えられる一般速度ベクトルの式につながる。なお、式の右辺に提示するように、テイラー級数展開の条件は、式B3=38bにより与えられる。追加の減衰係数BETAnにより、この条件は、より多くの状況で、より容易に満たされることがわかる。この条件、および結果としてのモデルに則した指標は、時系列全体の因果的性質であることが、確認される。全周波数について満たされるこの条件下で、時間領域での一般速度ベクトルモデルは、式B4=39bの形式で特定され、これは、(式39での従来の速度ベクトルの場合と同様に)
- t=0での第1のピーク、これにより、その方向がDoAを与え、ベクトルU0は、式B5を正規化することにより得られ、
- 反射音と同じ数の列で、各々が、反射音と直接音との干渉と関連づけられ、これについて、一定の間隔TAUnで観察されるベクトルの値が、式B6に現れ、
- SARCとして示す遅延と組み合わされた列で、後続の推定手順では用いられないものを、示している。
補遺の最後の方の式B6から始めて、連続するベクトル列間の特定の関係、特に、最も顕著な2つの第1のベクトルV'(TAUn)とV'(2.TAUn)との間の関係が得られる。これにより、式B7に、係数(-Gn/BETAn)を示し、これは「RHO」として示され、その式B8により、同じ列の第1の2つのベクトルV'(TAUn)およびV'(2.TAUn)のスカラー積として、推定値が提供され、このスカラー積は、第1のもののノルムの2乗で除算される。RHO係数を、式B6へと復帰させることにより、得られた式B9が再編成されて、式B10が得られる。この式の右辺は、ベクトルUn(具体的には、第1の反射音と、それに関連する列とに注目した場合、ベクトルU1)を示し、係数NU0/Gnの影響を受け、これは正であり(位相反転を伴う反射音のような、おそらくは稀な状況は除く)、このため、左辺V'(TAUn)-RHO.V'(0)を正規化することにより、取得される。
さらに、全体の係数NU0は、影響する他の係数を、基準指向性よりも、統合する可能性が高く、たとえば、全体としての振幅の低減となり、これは、信号源の周波数帯域幅を制限することにより、かつ/または、ノイズを部分的にマスキングすることにより、起こりうる(後者の効果は、一般に、モデル化するには複雑であるが)。なお、興味深いことに、究極的には、ベクトルU1(または、より一般的にはUn)の方向は、同様に推定可能であり、この全体としての振幅低減NU0の理由となる。
また、この推定態様は、従来の速度ベクトルにも適用される(その場合、BETAn=1を考慮することが、単に必要)。
特に、一般速度ベクトルを用いてDoAなどのパラメータを決定する実施形態の実際的な例について、以下に説明する。
図6Aから図6Dおよび図7を説明のために参照して記述する本実施形態では、遅延の第1の推定値(図7のステップS71)が得られ、従来の速度ベクトル
V(f)=1/W(f)[X(f),Y(f),Z(f)]T
が、ここで、たとえば、1次アンビソニック成分について、「通常」算出される。
ステップS721では、上述の計算は、従来の速度ベクトルV(f)の周波数の式に基づいて、従来の速度ベクトルV(t)の時間の式の推定時まで、実行される。
ステップS731では、従来の速度ベクトルV(t)の時間の式の解析は、ピークの時系列として実行される。この解析では、特に、効率的に(低干渉で)特定可能かどうか判別し、単一方向時系列が、補遺の式39および式40に述べるように(真に「因果的反応」として)正の時間でのみ展開する。
従来の速度ベクトルV(t)の式の時間領域での解析にて、負の時間でのピーク構造(たとえば、エネルギーまたは振幅が、選択された閾値THRよりも高い)が、典型的には、図6Aの負の横座標にあるピークなどとして発見された。そして、これは、式39では適正に特定できず、そのため、ピークV(t=0)にて与えられるDoAの推定に、バイアスがかかることを意味している。さらに、これは、式39に至るテイラー級数の収束条件が、満たされていないことの兆候となり、そのため、V(f)の分母の計算で、直接音に混ざった間接波の量が、相対的に大きすぎることの兆候となる。ここに提案する改良では、この間接波の割合は、空間フィルタにより低減する。このことは、速度ベクトルV(f)の推定における指向性(分母で作用しているもの)を、改良することが必要であることを意味している。
そして、空間フィルタは、既に取得した速度ベクトルからのV(t==0)として、ステップS751にて推定された方向(DoA)に、ビーム形成するために取得された、アンビソニックデータに対して適用される(ステップS71)。確かに、おそらくは誤りを含むが、この第1の推定値により、近似値を提供することができる。この近似値は、明らかに粗いものの、ビームの方向を、直接音の発信源に向けるとともに、より離れた角度部分から到来する反射音を低減するのには、十分である。
修正速度ベクトルV'(f)、そして、時間領域ではV'(t)(ステップS781)が、これらのフィルタ処理されたアンビソニックデータに基づいて、算出される。
テストS732では、修正速度ベクトルV'(t)の時間式に、0未満の時間座標のピークが残っているかどうか、判別する。先行する図6Aに比べて改良されたものが観測可能であっても、図6Bに例として示すように、負の時間座標にある信号構造(たとえば、エネルギー(図7で「NRJ」として示す)として推定され、たとえば、負の時間での信号の積分により与えられる)が、閾値THRよりも大きく残っているかどうか、判別可能である。
この場合、本方法は、既に取得したDoA(n)の粗い推定値を得ることにより、繰り返され(S752)、指向性が直接音の方向を示す基準成分D(f)(ステップS762の方法のn回目の反復処理について、D(n)として示す)を決定し、D(f)は、先行する反復処理での推定値D(n-1)よりも高分離であり、一般速度ベクトルV'(f)の推定における後者のD(n-1)を置き換えて(S772)、ステップS782でV'(t)を得る。このように、基準成分の指向性が、直接音の推定方向を、先行する反復処理での基準成分よりも高分離に、「捕捉」する。本実施形態では、1を超えるアンビソニック次数とする必要がなく、指向性の姿勢および形状の双方を調整して、直接音をより良く補足することにより、たとえばある種の反射音を、さらに、無くしていくことができる。
このように、本方法は、図6Cに示すように、負の時間でのピークの振幅またはエネルギーが、選択された閾値THRよりも小さくなるまで、繰り返され得る。
このように、直接音に向かって向上している分離性が、反復処理が進行するにつれて、1次の式での速度ベクトル(従来、そして一般)の分母にある成分へと、継続的に伝えられる。図6Dでは、このように、球(薄い灰色)の形状における全方向W成分から、本例では、濃い灰色のスーパーカーディオイド(supercardioid)の形状である、より高分離な成分D(1)へと、次に、非常に濃い灰色のより狭いスーパーカーディオイドの形状のD(2)へと移行する。
図7に示す方法に戻り、さらに詳細に説明すると、最初のステップS71は、周波数領域にて従来の速度ベクトルV(f)を、全方向成分W(f)を分母に伴って構成して始まる。ステップS721では、時間領域での式V(t)が推定される。そして、テストS731では、信号構造が特定され、ピークがある従来の速度ベクトルV(t)の時間式を示し、負の時間座標(t<0)でこの信号構造のエネルギーNRJが、固定閾値THR(矢印NOK)より低くとどまり、現在の音響状況により、バイアスのないDoAが、従来の速度ベクトルから直接的に、既に抽出可能である。この場合、パラメータDoA、U0、U1などは、上述のように、ステップS741にて決定可能である。そうでなければ(テストS731から出る矢印OK)、従来の速度ベクトルによるDoAの直接推定にはバイアスがかかっており、以下に述べるように、一般速度ベクトルを決定するために速度ベクトルを調整する、少なくとも第1回の反復処理(n=1)を、実行する必要がある。
このDoAの測定値から、バイアスがかかっていても(ステップS751で取得)、ステップS761にて、周波数領域にて基準成分D(1)が推定され、ステップS771にて、今度は「一般」速度ベクトルV'(f)の式における全方向成分W(f)を置き換える。ステップS781では、一般速度V'(t)の時間式は、テストS732にて、顕著なエネルギー(閾値THRを超える)が、負の時間座標でのこの式V'(t)の信号構造に残るかどうか決定するために、推定される。そうでない場合(テストS732から出る矢印NOK)、処理は、パラメータDoAなどを与えることにより、この第1の反復処理で停止して、ステップS742へ移行する。そうでなければ、ステップS791での処理の反復インデックスnを更新することにより、処理が繰り返される(ここで、S79xで示すステップは、処理の反復処理を基準とし、インデックスnを増加させるか(ステップS793)、あるいは、処理S792~S794の終了を決定する)。
上述のように、先行する反復処理(ステップS752)にて推定された粗いDoAに基づき、新たな基準成分D(n)が、ステップS762にて推定されて、ステップS772にて、一般速度ベクトルV'(f)の分母における以前の基準成分D(n-1)を置き換える。周波数領域での一般速度ベクトルV'(f)のこの新たな式から、時間領域でのこの式V'(t)が、ステップS782にて決定される。その信号構造の比較(閾値THRを基準としたエネルギー)が、テストS733にて繰り返されて、これから推定可能な新たなDoAに、バイアスがかかっているかどうか、判別する。そうでない場合(テストS733から出る矢印NOK)、パラメータ、特にDoAなどが、ステップS743にて、図6A~図6Cの図示例にあるように、3回反復した後に、取得可能である。そうでなければ(テストS733から出る矢印OK)、粗く、おそらくはバイアスがかかっていても、最後に推定されたDoAを用いて、処理は、ステップS752から再度反復される必要がある。
以前に推定されたDoAから、基準成分D(f)を算出する可能な例を、以下に説明する。補遺にある式にみられる形式では、成分DΘ(f)が、以下のアンビソニック成分の行列をとることから(またはベクトルDΘによって重み付けられた和から)、通例、抽出される。
DΘ(f)=DΘ.B(f)、ここで、
- B(f)は、周波数領域でのアンビソニック場を記述する信号のベクトル、たとえば、
など、ここで、
であり、信号S(f)を乗せて、
により記述される方向から到来する平面波の場合であり、直接波と間接波との混合について、
であり、ここで、
が、球面調和関数係数
のベクトルであり、
- DΘは、「ステアリングベクトル(steering vector)」型のベクトルであってもよく、特定の方向Θ(単位ベクトル
によっても指定可能)に一般的に向いてビーム形成され、
「ステアリングゲイン(steering gain)」を規定することにより、
となり、「ステアリングゲイン」は、
となるように規定されている。
1次によっては、球関数
は、軸対称であり、とり得る自由度は、姿勢に加えて、(もしあるならば)正および負のローブ間の割合にのみ影響する。2次から始めて、ベクトルDΘの係数にしたがい、主ローブ(おそらくは、目標方向)は、必ずしも対称である必要はなく、2次ローブは、さらに様々な形状および姿勢をとり得る。したがって、DΘの表記に関わらず、ビーム形成は、目標の主方向Θだけによってパラメータ化されているわけではない。
指向性が軸対称に形成された特定の場合、これは、DΘ=Y(Θ).Diag(gbeamshape)
の形式となる。
そして、ゲイン
は、スカラー積
の多項式Pbeamshape(・)として表される(以下の、ルジャンドル多項式の変形)。
対角係数gbeamshapeは、以下を考慮に入れる。
- 一方では、アンビソニック符号変換の選択、それによる球面調和関数
の計算、
- 他方では、これらの対角係数を選択することによって選択されたビーム形状を、調整することにより、ビームの形状を扱う基準(たとえば、主ローブおよび2次ローブの精密さおよび割合)。
これらの側面について、Jerome Daniel(2000)の論文、特に、182~186ページおよび図3.14を参照することが有用である。ここで、空間復号のために提案されているツールは、ここでD基準成分として示す、アンビソニック信号からの指向性構成に対して、直接的に適用可能である。
なお、このような形式にしたがい、gbeamshapeの関数として上述した成分D(f)に対して分離性を付与する、ゲイン係数を、規定することができる。
なお、再び図7を参照すると、音響状況によって、良好な「因果」視点での時間領域にて、一般速度ベクトルV'(t)についての形式が発見されないため、良好なDoA推定値が得られないようになっている。再び図7を参照すると、たとえば、時間領域にて一般速度ベクトルV'(t)を示す信号の形状を、これ以上改善できない場合、処理の反復処理から離脱するように、終了基準が付加される。このように、前回の反復処理の最後に、テストS733にて、この信号のエネルギーが、負の時間座標(t<0)にて閾値THRよりも大きく、かつ、処理の反復が、一般速度ベクトルV'(t)の推定値を改善しない(すなわち、これ以上改善しない)場合、すなわち、以下の両者にて示し得る場合、
- 負の時間にて、閾値THRを超えている信号エネルギー、および
- ある反復処理n-1から次回の反復処理nへ進んでも低減しない信号エネルギー(テストS792から出る矢印NOK)、
処理の反復は、ステップS794にて停止し得る。
そうでなければ(テストS792から出る矢印OK)、処理は、ステップS793にて、次の反復処理のために、実行可能であり、反復カウンタnを増加させることにより開始する。
このような「最悪の」音響状況にて、バイアスなしのDoAが得られる収束に至らずに、処理がステップS794にて停止せざるを得ない場合が発生することを、最小化するために、上述の文献FR1911723での教示内容を採用することができ、これは、たとえば、バイアスなしの判定(たとえば、信号フレームの開始)の機会を増やすのに最良のフレームを取り出すことを可能とする解決策のためである。
確かに、文献FR1911723に記述されているように、この問題の相対的重要性により、ベクトルU0が、どの程度、DoAの適正な推定値(バイアスが弱い)を提供し得るのかを評価することができ、このように、推定のため、そして、あるフレームについて得られた推定値を優先的に保持することができるように、信頼係数を提供する。推定バイアスのリスクが、過剰であると判明した場合、図3Aおよび図3Bを参照して後述するように、最もこの問題に直面していないフレームが、選択され得る。
そして、以下に説明する実施形態は、特に、たとえば、図7を参照して上述した処理の第1の反復の間、従来の速度ベクトルの推定に適用され得る。既に文献FR1911723で説明されたように、このような処理が従来の速度ベクトルに適用されることは、以下で予め用いられる。
このため、時間サブフレームの周波数解析により、所与の空間についての第1のピークを観測することへと、進むことができる。信号の開始(エネルギーの立ち上がり、過渡など)があるフレームは、最も早い波面のみを含む音響の混合を観測することができ、それは、直接音、および1つまたは複数の反射音である(式38にしたがって、上述の「ガンマ和」が1未満にとどまるように)。
信号の開始を含むフレームについて、時間ウィンドウを、たとえば、非対称かつ全体的に減少する形状とすることにより、周波数解析用に(場合によっては、動的に)調整することができることで、ウィンドウの「凹凸」が、信号の立ち上がり端(開始、過渡)に対してより重みを与え、それにより、直接音に対して与えて次第に重みが軽くなる(たとえば、ほぼ指数関数的にであるが、これは必須ではない)。このように、初期の波面に比べて、以降の波面の振幅は、人工的に低減し、テイラー級数の収束条件に達し、単一方向の時間展開が保障される。
時間ウィンドウ処理の指数関数的減少型の例を、以下に示す。結果としての時間フィンガープリントを、波の到来方向の推定の際に実質的なバイアスなく、好適な場合へと導くために、この例は、解析された信号に適用される。便宜上、ゼロ時点として指定された時刻t0から有効な、この動作を実施する。この時刻は、好ましくは、無音から始まる信号の開始時に対応している。式42にあるように、ALPHA>0であり、s(t)およびh(t)を含む畳み込み形式を再び積分することにより、式43の形式を得る。
そして、式44は、この選択が正当となる指数関数を適切に含み、式45に与えられる形式が得られ、これは、式46を確立することを意味する。
したがって、直接音に加えられた反射音の組によるインパルス応答をモデル化すると、式47が得られる。
このように、ガンマの和が1以上である場合(「双方向級数」となることもある)、こうして「適応した」ゲインの和が(式48)が1未満となるように、減衰係数ALPHAを決定することが、常に可能である。
そして、時間フィンガープリントが、実際に単一方向となることが観測され、これは、ピークにより、減少指数関数ウィンドウ適用後、正の時間でのみ、明らかとなる(図5Bの下部)。実際には、観測された信号のエネルギーは、非常に迅速に指数関数的に減少することが、さらに観測され、信号の切り捨ての推定値に対する数的影響は、相対的に短い切り捨て時間を超えると、完全に無視できるようになる。換言すると、より短い時間で、全励起信号とその反響音の双方にわたる長期の解析の利点が得られる。確かに、観測された「TDVV」は、信号の動態によるエラーなく、干渉モデルに準拠している。したがって、このようなウィンドウによる重み付けには、利用可能な時間フィンガープリントを理想的に求めることができる、二重の特性がある。
実際には、予め反射音の振幅を知ることなく、ALPHAでの減衰を、決定することが必要であり、好ましくは、時間フィンガープリントの単一方向性を確実にするのに十分に低い値と、間接波を検出して推定する機会の減少を避けるために、低すぎはしない値との間の妥協点が探求される。たとえば、この値は、測定した現象を物理的に代表する時間tEXP(通例、5ms)による減衰係数aEXPについて、ALPHA=-(log aEXP)/tEXPのように決定され得る。
反復処理(たとえば、二分法による)が、減衰値を調整するために、実装され得る。閾値の減衰値を超えると、取得した時間フィンガープリントが、双方向で検出された場合、したがって、おそらくはバイアスのかかったベクトルU0で、より強い減衰で解析が繰り返され、そうでなければ、少なくとも推定値U0が適用され、以下のピークが、良好に識別可能でなければ(減衰による低減のため)、2つの先行するものの中間的な減衰で、解析が繰り返されるなどであり、必要であれば、ベクトルU1が推定できるまで繰り返される。
それでもなお、ウィンドウを指数関数的に減少させる手法は、特に、干渉が顕著に増幅されるウィンドウ処理の開始時に、干渉に敏感であり得る。ウィンドウ処理の開始時、起動が直前であれば、ノイズ以外の干渉は、単に発信源自体の反響音であることもある。そして、このような干渉を低減するために、処理のノイズ除去を導入可能である。
一般に、様々な形状および/または大きさの時間ウィンドウが、提供され得るものであり、「好適なフィンガープリント」を得る機会を最大化するために、ウィンドウ間で重なることもある。
当初のDFTの大きさは、この解析ウィンドウよりも一般的に大きいものとして選択される。
これは、サンプルの連続したブロック(すなわち「フレーム」)の形式で、所与のサンプリング周波数でサンプリングされる、デジタル音響信号処理の状況におけるものであることは、当然である。
また、開始や過渡などの検出に対して、事前処理を、そして、時間周波数のノイズ除去を、任意的に提供可能であり、これは、たとえば、他の環境からの要素の導入を回避し、かつ/または、場の発信源を干渉フィンガープリントへと拡散させるように、マスク(二値である場合もある時間周波数フィルタ)を規定することによるものである。マスク(逆変換の結果)のインパルス応答を算出して、ピークの解析に対するマスクの影響を検証する必要がある。あるいは、おそらくは同様に干渉する混合に対応する周波数フィンガープリントの加重平均を、次に算出するように、格納されるフレームのフィンガープリントの周波数加重へと統合可能である(通例、信号開始時に、該当する発信源が移動していないことを確認することによるものであり、これは、遅延の推定により推察可能)。
そして、このように、ピークの抽出および観測に進む、たとえば、ノルム|V(t)|:最大ピーク、そして、次に、TAU1(一般に)などを与えることによる。
そして、時間フィンガープリントの分析に進む。分析は、({tau_n}およびV(sum(k_n.tau_n))にしたがって):
- 短すぎる時間的補助(too-short temporal support)へのFFTの選択による、時間的再ループ(ある種の循環型「エイリアシング」)があるかどうか、
- 単方向増加級数(progressive unidirectional series)、もしくは逆に、双方向級数(bidirectional series)があるか、
または、目立った減衰のない特定の級数の場合(ゲイン和(gn)が1に近いままの場合)、もしくは、逆進級数(retrograde series)の場合(少なくとも暗黙的ゲインg_n>1の場合)、
を検出することによってなされる。
そして、
- 「良好なフレーム」または「良好なフィンガープリント」スコア(信頼性のある推定値、単方向性であるため、おそらくは、DoAバイアスのないものを、可能とする)を割り当てて保存し、
- 推定値(Un)を得て、
- 必要であれば、適切な時間ウィンドウを選択することにより、上流側解析を調整する。
時間フィンガープリントの解析については、上述したが、周波数解析は、以下のように、より簡潔に実行可能である。
時刻ゼロでのピークは、構成により、完全なスペクトルにわたり速度ベクトルの平均(ハミルトニアンの対称性により、実部がなくなる)に、あるいは、正の周波数のみを考慮する場合には、その実部にさえ等しいことが、数学的に容易に理解されるであろう。直接音にのみ関心がある場合、FDVVの逆変換を算出してDoAの推定値を得ることに意味がないことが、推定され得る。しかしながら、時間に基づくTDVVの推定により、このDoAに信頼性があるかどうか(正で進行する時間につれて発展する基準)が、検出可能となる。
この好適な場合は、混合が依然として非常に複雑というわけでなければ、元の信号の開始の際に、よりもっともらしく観測される。これは、一般に、これらの時期に、測定値を得るのに十分である。
さらに、実際には、VVの周波数および時間フィンガープリントは、干渉波の混合の理想的モデルにて、常に識別可能というわけではない。元の信号は、送信電力の欠如により、主要な各時点にて、周波数の顕著な範囲を十分にまたは常に励起するわけではなく、捕捉された音場の他の成分の同時並行性を考慮に入れる場合もある(不十分なSNRまたはSIR)。これは、多少拡散した背景音(他の音源)、マイクロフォンノイズにつながり得る。
そして、以下のうちの少なくとも1つ、あるいはいくつかの組み合わせにしたがって、処理を適用可能である。
- 改良されたアルゴリズムにしたがう信号開始検出を伴う、時間周波数サンプルの選択、
- いくつかのフレームにわたる(たとえば、おそらくは動的に、注目するフレームの|W(f)|2および忘却因子により重み付けされたV(f)の平均)速度ベクトルの平滑化、
- 信号開始フレームの選択について、|W(f)|2により重み付けされたV(f)の平均を算出して(抽出された遅延が同じ場合)、周波数フィンガープリントを補い、時間フィンガープリントを強化。
算出の際の効率化のため、情報においてより一貫したものとして検出されたフレームについてのみ、TDVV、または上流ではFDVVの算出を実行することも推奨され得るものであり、たとえば、この情報は、簡潔な処理により検出可能である状況で、信号開始フレームであり、この場合、信号の立ち上がり端で、解析ウィンドウを位置づけるのに依然として有利である。
非整数遅延の良好な推定値(時系列における分数遅延およびその倍数)について、サンプル間補間によるピーク推定、および/または局所周波数解析(時間的に制約のある近傍にてピークを分離することによる)が、考慮され得るものであり、遅延は、位相応答に基づいて精密に調整される。
時間ピークの事前選択は、級数の特徴のある遅延の現在の推定値に基づいて、なされ得る。
このように、実施されたステップは、図3Aおよび図3Bに示すような、実行可能な例示的実施形態では、要約されている。ステップS1では、アンビソニック信号の(時間から周波数への)フーリエ変換が、算出され、それは、一連の「フレーム」の形態(一連のサンプルのブロック)でなされる。各変換フレームk(ステップS2)について、動的マスクが、信号対ノイズ比が閾値を下回るいくつかの周波数帯域に対して適用可能である(いくつかの周波数帯域は、実際に、高雑音であり得るものであり、たとえば、マイクロフォン固有であるかまたは他のものでのノイズにより、この周波数帯にて取得された信号の利用が損なわれる)。特に、周波数帯域によるノイズを探索することが、ステップS3にて「全方向」成分Wについて優先的に実行され、ノイズにより変更された周波数帯域(たとえば、SNR<0dBというような閾値を超過)が、ステップS4にてマスキングされる(すなわち、ゼロに設定)。
そして、ステップS5では、速度ベクトルV(f)が、たとえば式6(または式11、式18もしくは式20の形式で)により、周波数領域にて算出される。
ここで、一例示的実施形態では、後述のように算出する重みq(f)を適用して、周波数帯域fに対して、何かしらの重要性を付与する。このような実装例により、速度ベクトルV(f)を、展開が顕著な周波数帯域にて表すことが可能となる。このようにするために、最適な重みが、U0およびV(f)の関数として、反復して算出される。このように、図3Aのアルゴリズム処理に戻り、ステップS6では、様々な重みq(f)が、1として設定される。ステップS7では、帯域ごとにV(f)に対して適用された重みq(f)を、適用して、Vbar(f)=q(f)V(f)とする。ステップS8では、U0は、各フレームkについて判別されて、
U0(k)=E(Re(Vbar(f)))となり、ここで、E(x)は、たとえば、xの期待値であり、したがって推定された速度ベクトルVbar(f)の実部の全周波数にわたる平均と、同等である。
この第1の推定値U0(k)が粗いのは、当然である。それは、先行する判別でのU0(k)について、重みが計算されることにより、ベクトルV(f)の虚部に基づき、式49を用いて反復して調整され、ここで、ベクトルmは、ベクトルU0により規定される平面に垂直な単位ベクトルであり、壁の法線である(たとえば、図1の方向z)。また、ベクトルmは、ステップS9でのU0の関数として反復して推定され、重みが、ステップS10での式49により算出される。得られた重みは、ステップS7にて適用され、U0の推定は、テストS11の出力にて収束に至るまで、調整される。この段階で、U0(k)は、様々なフレームについて推定済である。
ここから、U1は、上述の式41のような関係により、推論可能である。ここで説明した変形例では、U1は、式50~式52により決定される。これらの式は、ステップS7で得られたベクトルVbar(f)に対して、逆変換IDFT(周波数から時間へ)を、ステップS12にて前もって適用して、速度ベクトルの時間に基づく表現V(t)を取得する。このような実施形態により、図2を参照して理解されるように、異なる反射面のそれぞれの場合に、様々な遅延TAU1、TAU2などを特定可能である。第1の遅延TAU1が、直線経路の受信の瞬間に続く時刻での第1のピークV(t)であるため、特定される。このように式51では、tmax(k)は、フレームkについて算出されたV(t)kの絶対値を最大化する瞬間である。
テストS13は、各フレームについて、V(t=0)の絶対値が、実際に、t>0でV(t)よりも大きいことを確認する。この条件を満たさないフレームは、ステップS14にて廃棄される。様々な遅延TAU1そしてTAU2が、次に、ステップS15で決定される(遅延TAU1に対応するものを、式51にて比較されるV(t)kの絶対値から除去することによる)など。遅延TAUmは、時刻tおよびtmax(k)が、まずは、サンプルインデックス(時刻ゼロが、インデックスゼロの基準として用いられる)について表現されることを考慮して、各反復処理mにて得られた成分tmaxが、式52にしたがってサンプリング周波数fsにより除算されることにより、与えられる。そして、ベクトルU1、U2なども、式50により、算出可能である。
また、特に、ステップS16にて、式28で与えられるd0など、他のパラメータも決定可能である(次に、テストS17にて、d0min=0およびd0max=5mなどの従来の室内データの一貫性を検証することによる。そうでない場合、フレームは、エラーを含み、ステップS14にて廃棄され得る)。
ステップS18は、第1の反射音での音の開始を代表する「良好な」フレームを選択することからなる。このようなフレームを選択するための基準D(k)は、例として、式53により、説明され得る。ここで、C(f)i(k)は、フレームkの第1の変換(時間から周波数へ)に由来する時間周波数サンプル(t,f)において、アンビソニックチャネルiにて検出される大きさ(振幅の絶対値)を指定する。イプシロンは、ゼロでない正の値を指定して、信号がなくて分母がゼロになることを避ける。Fは、用いられる周波数サブバンドの総数を指定する。
このように、全フレームD(k)の基準のうち、式53から算出される基準D(k)が、ステップS21で得られる最大値Dmaxの90%よりも小さくはないフレームのみを、ステップS22にて保持可能である。
このように、ステップS18では、全フレームについて、D(k)の値が算出されて、ステップS19では、この処理は、U0(k)、d0(k)、D(k)を、様々なフレームに対して提供する。ステップS20では、D(k)の値は、ステップS21にて最高値を特定し、ステップS22にてD(k)値が0.9Dmax未満のフレームを廃棄するために、収集される。
最後に、ステップS23にて、ここで保持されているベクトルU0は、好ましくは、保持されている様々なフレームのベクトルU0の(平均値よりもむしろ)中央値である。また、保持されている距離d0は、保持されている様々なフレームの距離d0の中央値である。
本発明は、例として上述した実施形態に限定されないことは当然である。このことは、他の変形例にも適用される。
1次のアンビソニック(FOA)信号の処理に対する適用例について上述した。次数は、空間解像度を豊かにするために高くしてもよい。
確かに、1次のアンビソニック表現について説明したが、より高次とすることも可能である。このような場合、速度ベクトルの計算量は、より高次の方向成分の成分W(f)に対する比によって増加し、ベクトルUnは、次数が高くなると暗黙的に増加する。次元の増加(3を超えて)、およびこれに伴う空間解像度により、ベクトルU0、U1、…、Un相互の区別が良好となり、時間フィンガープリントについて、(U0-Un)に比例するピークV(k*TAUn)の検出がベクトルU0およびUnが角度的に近い場合であっても容易になり、この場合というのは、グレージング反射音のときに生じる(たとえば、発信源が地面から遠い、かつ/または近い場合)。したがって、これにより、所望のパラメータU0、U1、d0などをより正確に推定することが可能になる。また、3つの第1次の成分(X,Y,Z)のみが、ここで分子に保持されるということが、基準成分を分母に構築するのに使えるより高次の成分を利用可能であることから、独立するようになっていることが明記される。全ての場合(分母に関わらず)にて、上述の処理(およびFR1911723に示された処理)を、より高次の成分を分子に加えることにより向上させることが、考慮され得るものであり、このように、速度ベクトルの次元を増やし、特に、時間領域にてピークをより区別できるようになる。
より一般的には、速度ベクトルは、周波数領域にて「同時的」な空間的音響表現の成分間の比により、置き換え可能であり、上述の空間表現の特徴のある座標系で機能し得る。
複数の発信源の例の場合を軽減するため、より一般的に、ニューラルネットワークを含む人工知能の方法により、TDVVの計算が利用可能である。考えられるある種の学習戦略(たとえば、モデルからのフィンガープリントについて、またはウィンドウ処理SRIR、および必ずしも元の信号からである必要はない)により、ネットワークは、所与の室内状況を検出および推定するために、一連のフレームを利用することを学習可能である。
さらに、DoAの第1の粗い推定値を決定するために、従来の速度ベクトルV(f)を最初に推定する可能性、そして、より正確なDoAを得るために、この粗い推定値に基づいて、一般速度ベクトルV'(f)の推定値を調整する可能性について、図6A~図6Dおよび図7を参照して上述した。これは、単に、ある可能な実装例であることは当然である。たとえば、ある変形例では、空間を、直接、いくつかの部分に分割することが可能であり、このように、これらの部分の各々における分母D(f)内の成分に対して、指向性を与え、DoA算出にて、反復的方法(ステップS761から開始する図7の型であり、D(1)が当該角度部分の関数として簡潔に算出されてもよい)により収束を試み、これらの反復的方法は、特に、各部分に対して並行して実行される。これは、一つの可能な例示的実装例である。他には、「最良」(またはいくつかの最良の)角度部分が、上述の因果モデルの有効性の基準により選択され得るものであり、その方向またはそれらの方向に最適化された推定値が選択され、形式の変形例として分母に含まれている。確かに、より一般的に、第1のステップにて、かつ/または後続のステップの間に、様々に形成された指向性(分母にて)とそれぞれ関連づけられた複数の一般速度ベクトルを、評価することを、さらに予測することができる。
より一般的には、信頼性のある音の到来方向(すなわちDoA)を与えるビーム形成を探索することが、最適化の一般的な問題として探求され、この範囲にて、様々な戦略が要求され得る。したがって、以下のことを特定可能である。
- 因果モデルの有効性を表現/予測する最小化基準(最小化される関数)。ある程度単純で、そのため改良可能な方法にて、上記の説明では、これが、負の時間座標での信号の相対エネルギーの探索に注目している。
- 最適化されるべきパラメータは、ビーム形成パラメータであり、たとえば、式A1およびA4に含まれる行列Dの係数、またはシータ(Θ)方向として定式化され、軸対称を選ぶ場合、すなわち、1次に制限される場合には、式EQ.A5のビーム形状パラメータgbeamshapeとして定式化される。
- 最初のパラメータの1つ(または複数)の組は、通例、「全方向の」指向性(gbeamshape=[1 0…])、または、以前の使用の際に格納された好適な指向性、あるいはその他、空間を表す方向の組を指す複数の指向性、さらに、指向性の様々な形状をも表現。
- 通例、最新の推定DoAにおける音響ビームを再配置することが、常に十分に安定的な選択ではないため、テストされたパラメータを(反復処理の際)調整する原理。改良の欠如によってアルゴリズムを停止するよりは、格納された状況の1つ(おそらくは、最小化基準の視点から最良)から再度開始して、他の軸(たとえば、指向性パラメータの形式)または他の軸の組み合わせに沿ってパラメータを調整するために、必要とされている。
より一般的には、たとえば、確率勾配またはバッチ最適化を含む通常の手法が検討されるが、何度もなされる反復処理には、かなりのコストがかかることがある。
なお、共通の最適化タスクとは異なり、最終目的のパラメータ(通例、Unベクトル)は、直接的には最適化されたものではないが、それから結果として得られるものである。なお、因果モデルへの準拠に帰結するという意味で、潜在的には、全てが「最適な」ビーム形成パラメータの多くの異なる組がある。これにより、同じパラメータUnの組を、おそらくは正確に、推論することが可能である。
補遺
式13 g0=1;τ0=0
式15 γ0=1
式27 d1-d01c
式44 e-αt=e-α(t-τ).e-ατ
式A1 DΘ(f)=DΘ.B(f)
式A5 DΘ=Y(Θ).Diag(gbeamshape)
MEM ワークメモリ
PROC プロセッサ
OUT 出力インターフェース
IN 入力インターフェース

Claims (18)

  1. 少なくとも1つのマイクロフォンにより取得された音響信号を処理する方法であって、
    少なくとも1つの壁がある空間内の少なくとも1つの音源を測位するために、
    - 時間周波数変換が、取得された信号に適用され、
    - 取得された信号に基づき、一般速度ベクトルV'(f)が、周波数領域にて表され、全方向成分W(f)以外の基準成分D(f)が分母に現れる、速度ベクトルV(f)の式から推定され、その式は、実部および虚部のある複素式であり、前記一般速度ベクトルV'(f)は、
    *前記発信源と前記マイクロフォンとの直線経路である、第1のベクトルU0として表される、第1の音響経路と、
    *前記壁での反射音によって生じ、第2のベクトルU1として表される、少なくとも第2の音響経路との合成を特徴づけ、
    前記第2の経路は、前記マイクロフォンにて、前記直線経路を基準として、第1の遅延TAU1があり、
    - 前記遅延TAU1、前記第1のベクトルU0、および前記第2のベクトルU1の関数として、
    *前記直線経路の方向(DoA)と、
    *前記発信源から前記マイクロフォンまでの距離d0と、
    *前記発信源から前記壁までの距離z0とのうちの少なくとも1つのパラメータが決定される、方法。
  2. 複数回の反復処理であって、以前の回の反復処理で取得された前記直線経路の方向(DoA)の近似値に基づいて決定された基準成分D(f)を、分母に含んだ、前記一般速度ベクトルV'(f)が、少なくとも一部に用いられた、複数回の反復処理を含む、請求項1に記載の方法。
  3. 前記速度ベクトルV(f)が、前記一般速度ベクトルV'(f)の代わりに用いられた、第1の反復処理を含み、前記第1の反復処理の最後にて、少なくとも前記直線経路の方向(DoA)の第1の近似値を決定するために、前記速度ベクトルV(f)は、周波数領域で表され、分母に前記全方向成分W(f)を有し、
    少なくとも、前記第1の反復処理に続く第2の反復処理に、前記一般速度ベクトルV'(f)が、用いられ、分母の前記全方向成分W(f)が前記基準成分D(f)に置き換わった、前記速度ベクトルV(f)の式から推定され、前記基準成分D(f)は、前記全方向成分W(f)よりも、空間的に高分離である、請求項2に記載の方法。
  4. 前記基準成分D(f)は、前記直線経路の方向(DoA)の前記第1の近似値に対応した方向で、前記全方向成分W(f)よりも高分離である、請求項3に記載の方法。
  5. 前記反復処理は、因果関係基準にしたがって、収束に至るまで繰り返される、請求項2から4のいずれか一項に記載の方法。
  6. 各反復処理にて、
    - 前記直線経路(DoA)に沿った前記音の到来に関連したピークに加えて、各々が少なくとも1つの壁での反射音に関連した一連のピークを、時間領域にて取得するために、周波数から時間への逆変換も、前記一般速度ベクトルV'(f)の前記式に対して適用され、
    - 時間座標が、前記直線経路のもの未満であり、振幅が、選択閾値を超えた信号が、前記一連のピークに現れた場合には、新規の反復処理が実行され、
    前記因果関係基準は、前記信号の振幅が、前記閾値未満である場合に満たされる、請求項5に記載の方法。
  7. 前記反復処理は、
    - 前記信号の振幅が、選択閾値未満である、第1の場合、
    - 前記反復処理の繰り返しが、前記信号の振幅の大幅な減少につながらない、第2の場合に、
    終了する、請求項5または6に記載の方法。
  8. 前記第2の場合に続いて、以下の各ステップを実施し、前記取得された信号は、サンプルの連続フレームの形式として得られ、各ステップは、
    - 各フレームに対し、そのフレーム内に音の開始があること(式53)について、スコアが推定されるステップと(S18)、
    - 閾値より高スコアのフレームが、該フレーム内で取得された音響信号を処理するために、選択されるステップ(S22)とである、請求項7に記載の方法。
  9. 前記取得された信号は、アンビソニックマイクロフォンにより受信され、前記速度ベクトルV(f)は、周波数領域にて、以下の種類の形式の1次アンビソニック成分により表され、
    V(f)=1/W(f)[X(f),Y(f),Z(f)]T,
    W(f)は、前記全方向成分であり、
    前記一般速度ベクトルV'(f)は、周波数領域にて、以下の種類の形式の1次アンビソニック成分により表され、
    V(f)=1/D(f)[X(f),Y(f),Z(f)]T,
    D(f)は、前記全方向成分以外の前記基準成分である、請求項1から8のいずれか一項に記載の方法。
  10. 前記第1のベクトルU0と同等の前記直線経路の方向の推定値は、周波数領域にて表された前記一般速度ベクトルV'(f)の実部の周波数の集合における平均値から、決定される(式24)、請求項1から9のいずれか一項に記載の方法。
  11. - 周波数から時間への逆変換が、前記一般速度ベクトルに対して、時間領域V'(t)内にて表すために適用され、
    - 前記直線経路の時間後に、少なくとも1つの最大値が、前記一般速度ベクトルV'(t)maxの式にて、時間の関数として検索され、
    - そこから、最大値V'(t)maxとなる時間に対応する第1の遅延TAU1が、推論される、請求項1から10のいずれか一項に記載の方法。
  12. - 前記第2のベクトルU1は、ベクトルV1を以下のように定義するために、時間指標t=0、TAU1、および2xTAU1にて記録された、正規化速度ベクトルV'の値の関数として推定され、
    V1=V'(TAU1)-((V'(TAU1).V'(2.TAU1))/||V'(TAU1)||2)V'(0),
    ベクトルU1は、U1=V1/||V1||
    により与えられる、請求項11に記載の方法。
  13. - 第1のベクトルU0および第2のベクトルU1のそれぞれの角度PHI0およびPHI1は、前記壁について、以下のように決定され、
    PHI0=arcsin(U0.nR)およびPHI1=arcsin(U1.nR)、ここで、nRは、単位ベクトルであり、前記壁に垂直であり、
    - 前記発信源と前記マイクロフォンとの距離d0は、第1の遅延TAU1の関数として、以下の種類の関係式
    d0=(TAU1xC)/((cosPHI0/cosPHI1)-1)、ここでCは音速、により決定される、請求項12に記載の方法。
  14. 前記発信源から前記壁までの距離z0は、以下の種類の関係式
    z0=d0(sinPHI0-sinPHI1)/2
    により決定される、請求項13に記載の方法。
  15. 前記空間は、複数の壁を備え、
    - 周波数から時間への逆変換が、前記一般速度ベクトルに対して、時間領域V'(t)内にて、一連のピークの形式で表すために適用され(式39b、図2)、
    - 一連のピークにて、前記の複数の壁のうちの1つの壁での反射音に関連したピークが、特定され、各々の特定されたピークは、前記直線経路を基準としたときの、対応する壁nでの反射音によって生じる音響経路の第1の遅延TAUnの関数である時間座標を有し、
    - 壁nでの反射音よって生じる音響経路を表す前記第1のベクトルU0および各第2のベクトルUnの各第1の遅延TAUnの関数として、少なくとも1つのパラメータが、
    *前記直線経路の方向(DoA)と、
    *前記発信源から前記マイクロフォンまでの距離d0と、
    *前記発信源から前記壁nまでの少なくとも1つの距離znとのうちから決定される、請求項1から14のいずれか一項に記載の方法。
  16. 壁nでの反射音に関連した前記ピークは、当該壁nと関連づけられた遅延TUnの倍数である時間座標を有し、最小の正の時間座標を有するピークの第1部分が、壁での単一反射音と関連づけられた当該部分の前記各ピークを特定するために、予め選択される、請求項15に記載の方法。
  17. 請求項1から16のいずれか一項に記載の方法を実施する処理回路を備えた、音響信号処理装置。
  18. 処理回路のプロセッサにより実行された場合に、請求項1から16のいずれか一項に記載の方法を実施する命令を含む、コンピュータプログラム。
JP2023530282A 2020-11-19 2021-10-15 改良型音響源測位法 Pending JP2023550434A (ja)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
FR2011874 2020-11-19
FR2011874A FR3116348A1 (fr) 2020-11-19 2020-11-19 Localisation perfectionnée d’une source acoustique
PCT/FR2021/051801 WO2022106765A1 (fr) 2020-11-19 2021-10-15 Localisation perfectionnée d'une source acoustique

Publications (1)

Publication Number Publication Date
JP2023550434A true JP2023550434A (ja) 2023-12-01

Family

ID=75108412

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2023530282A Pending JP2023550434A (ja) 2020-11-19 2021-10-15 改良型音響源測位法

Country Status (7)

Country Link
US (1) US20240012093A1 (ja)
EP (1) EP4248231A1 (ja)
JP (1) JP2023550434A (ja)
KR (1) KR20230109670A (ja)
CN (1) CN116472471A (ja)
FR (1) FR3116348A1 (ja)
WO (1) WO2022106765A1 (ja)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR3102325A1 (fr) * 2019-10-18 2021-04-23 Orange Localisation perfectionnée d’une source acoustique

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR2998438A1 (fr) * 2012-11-16 2014-05-23 France Telecom Acquisition de donnees sonores spatialisees
FR3067511A1 (fr) * 2017-06-09 2018-12-14 Orange Traitement de donnees sonores pour une separation de sources sonores dans un signal multicanal
FR3081641A1 (fr) * 2018-06-13 2019-11-29 Orange Localisation de sources sonores dans un environnement acoustique donne.

Also Published As

Publication number Publication date
WO2022106765A1 (fr) 2022-05-27
FR3116348A1 (fr) 2022-05-20
EP4248231A1 (fr) 2023-09-27
US20240012093A1 (en) 2024-01-11
WO2022106765A8 (fr) 2023-05-04
CN116472471A (zh) 2023-07-21
KR20230109670A (ko) 2023-07-20

Similar Documents

Publication Publication Date Title
RU2589469C2 (ru) Устройство и способ для позиционирования микрофона, основываясь на пространственной плотности мощности
TWI530201B (zh) 經由自抵達方向估值提取幾何資訊之聲音擷取技術
RU2555188C2 (ru) Устройство, система (варианты), способ получения информации о направлении и компьютерный программный продукт
Chazan et al. Multi-microphone speaker separation based on deep DOA estimation
JP6467736B2 (ja) 音源位置推定装置、音源位置推定方法および音源位置推定プログラム
JPWO2006059806A1 (ja) 音声認識装置
JP2007033445A (ja) 信号源の軌跡をモデル化する方法及びシステム
US20150088497A1 (en) Speech processing apparatus, speech processing method, and speech processing program
Kumatani et al. Beamforming with a maximum negentropy criterion
US10393571B2 (en) Estimation of reverberant energy component from active audio source
US20230026881A1 (en) Improved Localization of an Acoustic Source
Bologni et al. Acoustic reflectors localization from stereo recordings using neural networks
JP2023550434A (ja) 改良型音響源測位法
KR20210137146A (ko) 큐의 클러스터링을 사용한 음성 증강
Pertilä Acoustic source localization in a room environment and at moderate distances
Dehghan Firoozabadi et al. A novel nested circular microphone array and subband processing-based system for counting and DOA estimation of multiple simultaneous speakers
Krueger et al. Bayesian Feature Enhancement for ASR of Noisy Reverberant Real-World Data.
JP2005258215A (ja) 信号処理方法及び信号処理装置
Luo et al. Fast Random Approximation of Multi-channel Room Impulse Response
KR102624195B1 (ko) 음성의 명시적 공간 필터링을 위한 지도 학습 방법 및 시스템
EP4171064A1 (en) Spatial dependent feature extraction in neural network based audio processing
US11835625B2 (en) Acoustic-environment mismatch and proximity detection with a novel set of acoustic relative features and adaptive filtering
Nguyen et al. Location Estimation of Receivers in an Audio Room using Deep Learning with a Convolution Neural Network.
Olgun et al. Data-driven threshold selection for direct path dominance test
CN117437930A (zh) 用于多通道语音信号的处理方法、装置、设备和存储介质