(57)【要約】【課題】音源の三次元方向の位置を精度良く推定することができる音響処理装置、音響処理方法、およびプログラムを提供することを目的とする。【解決手段】音響処理装置は、位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位部と、２個の収音部の組ごとに当該収音部のそれぞれから音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を推定音源位置として推定する音源位置推定部と、を備える。【選択図】図１
2020038123.tif 000002
【特許請求の範囲】
【請求項１】
位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位部と、
２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を前記推定音源位置として推定する音源位置推定部と、
を備える音響処理装置。
【請求項２】
前記音源位置推定部は、
前記中間点の分布を複数のクラスタに分類し、前記推定音源位置が前記音源に対応するクラスタに分類される確率である推定確率が高くなるように前記推定音源位置を更新する、請求項１に記載の音響処理装置。
【請求項３】
前記推定確率は、前記定位音源方向が定められるとき前記推定音源方向が得られる確率である第１確率と、前記中間点が定められるとき前記推定音源位置が得られる確率である第２確率と、前記中間点が分類されるクラスタの出現確率である第３確率と、をそれぞれ因子とする積である、請求項２に記載の音響処理装置。
【請求項４】
前記第１確率は、前記定位音源方向を基準とするフォン・ミーゼス分布に従い、前記第２確率は、前記中間点の位置を基準とする多次元ガウス関数に従い、
前記音源位置推定部は、
前記推定確率が高くなるように、前記フォン・ミーゼス分布の形状パラメータと、前記多次元ガウス関数の平均ならびに分散と、を更新する、請求項３に記載の音響処理装置。
【請求項５】
前記音源位置推定部は、
前記収音部の３個から定められる３個の前記中間点の重心を前記推定音源位置の初期値として定める、請求項１から請求項４のいずれか一項に記載の音響処理装置。
【請求項６】
前記複数チャネルの音響信号から音源ごとの音源別信号に分離する音源分離部と、
前記音源別信号のスペクトルを算出する周波数分析部と、
前記スペクトルを複数の第２クラスタに分類し、前記第２クラスタのそれぞれに分類される各スペクトルに係る音源が同一であるか否かを判定し、
同一と判定した音源の前記推定音源位置を、同一でないと判定した音源よりも優先して選択する音源特定部と、
を備える請求項１から請求項５のいずれか一項に記載の音響処理装置。
【請求項７】
前記音源特定部は、
前記第２クラスタのそれぞれに分類されるスペクトルに係る音源の前記推定音源位置の分散に基づいて当該第２クラスタの安定性を評価し、
前記安定性が高い第２クラスタほど当該第２クラスタにスペクトルが分類される音源の前記推定音源位置を優先して選択する、請求項６に記載の音響処理装置。
【請求項８】
音響処理装置における音響処理方法であって、
前記音響処理装置が、
位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位過程と、
２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を前記推定音源位置として推定する音源位置推定過程と、
を有する音響処理方法。
【請求項９】
音響処理装置のコンピュータに、
位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位手順と、
２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を音源位置推定手順と、
を実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音響処理装置、音響処理方法、およびプログラムに関する。
【背景技術】
【０００２】
環境理解において音環境の情報を取得することは重要である。従来から、音環境における種々の音源や雑音から特定の音源を検出するために、音源定位、音源分離、音源同定などの要素技術が提案されている。特定の音源は、例えば、鳥の鳴き声や人の発話など、装置の利用者である受聴者にとって有用な音である。音源定位とは、音源の方向や位置を推定することを意味する。推定された音源の方向や位置は、音源分離や音源同定の手掛かりとなる。また、音源分離とは、複数の音それぞれを分離することを意味する。音源同定とは、音源を同定することを意味する。
【０００３】
このような音源分離や音源定位においては、複数のマイクロホンを有するマイクロホンアレイを用いて音声信号の収音を行う。このようなマイクロホンアレイでは、マイクロホン毎に音の到来時間差が異なる。これにより、装置は音源方向を得ることができる。
【０００４】
一方、現在、どのような鳥の種類やどの方向から鳴いているかなどの識別等の処理を、自動的に行う研究が行われている。このような鳥の定位においては、鳥の大きさが小さく、鳥の動きが速いため、方向推定よりも位置推定が好ましい。
このため、非特許文献１には、方位角データを用いて二次元音源位置を得るための三角測量手法が開示されている。
【先行技術文献】
【非特許文献】
【０００５】
【非特許文献１】Daniel Gabriel, et.al., "Iterative Outlier Removal Method Using In-Cluster Variance Changes in Multi-Microphone Array Sound Source Localization", The Proceedings of The 79th National Convention of IPSJ, pp. 2-229-230, Nagoya University, Nagoya, Mar. 16-18,2017.
【発明の概要】
【発明が解決しようとする課題】
【０００６】
しかしながら、非特許文献１に記載の技術では、音源が同じ方向に位置し且つ高さが異なる位置にある音源同士を区別することができなかった。
【０００７】
本発明は、上記の問題点に鑑みてなされたものであって、音源の三次元方向の位置を精度良く推定することができる音響処理装置、音響処理方法、およびプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００８】
（１）上記目的を達成するため、本発明の一態様に係る音響処理装置（１）は、位置が異なるＭ（Ｍは、３以上の整数）個の収音部（２０）のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位部（１２０）と、２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を前記推定音源位置として推定する音源位置推定部（１４）と、を備える。
【０００９】
（２）また、本発明の一態様に係る音響処理装置において、前記音源位置推定部は、前記中間点の分布を複数のクラスタに分類し、前記推定音源位置が前記音源に対応するクラスタに分類される確率である推定確率が高くなるように前記推定音源位置を更新するようにしてもよい。
【００１０】
（３）また、本発明の一態様に係る音響処理装置において、前記推定確率は、前記定位音源方向が定められるとき前記推定音源方向が得られる確率である第１確率と、前記中間点が定められるとき前記推定音源位置が得られる確率である第２確率と、前記中間点が分類されるクラスタの出現確率である第３確率と、をそれぞれ因子とする積であるようにしてもよい。
【００１１】
（４）また、本発明の一態様に係る音響処理装置において、前記第１確率は、前記定位音源方向を基準とするフォン・ミーゼス分布に従い、前記第２確率は、前記中間点の位置を基準とする多次元ガウス関数に従い、前記音源位置推定部は、前記推定確率が高くなるように、前記フォン・ミーゼス分布の形状パラメータと、前記多次元ガウス関数の平均ならびに分散と、を更新するようにしてもよい。
【００１２】
（５）また、本発明の一態様に係る音響処理装置において、前記音源位置推定部は、前記収音部の３個から定められる３個の前記中間点の重心を前記推定音源位置の初期値として定めるようにしてもよい。
【００１３】
（６）また、本発明の一態様に係る音響処理装置において、前記複数チャネルの音響信号から音源ごとの音源別信号に分離する音源分離部（１２２）と、前記音源別信号のスペクトルを算出する周波数分析部（１２４）と、前記スペクトルを複数の第２クラスタに分類し、前記第２クラスタのそれぞれに分類される各スペクトルに係る音源が同一であるか否かを判定し、同一と判定した音源の前記推定音源位置を、同一でないと判定した音源よりも優先して選択する音源特定部（１６）と、を備えるようにしてもよい。
【００１４】
（７）また、本発明の一態様に係る音響処理装置において、前記音源特定部は、前記第２クラスタのそれぞれに分類されるスペクトルに係る音源の前記推定音源位置の分散に基づいて当該第２クラスタの安定性を評価し、前記安定性が高い第２クラスタほど当該第２クラスタにスペクトルが分類される音源の前記推定音源位置を優先して選択するようにしてもよい。
【００１５】
（８）上記目的を達成するため、本発明の一態様に係る音響処理方法は、音響処理装置における音響処理方法であって、前記音響処理装置が、位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位過程と、２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を前記推定音源位置として推定する音源位置推定過程と、を有する。
【００１６】
（９）上記目的を達成するため、本発明の一態様に係るプログラムは、音響処理装置のコンピュータに、位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位手順と、２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を音源位置推定手順と、を実行させる。
【発明の効果】
【００１７】
上述した（１）または（８）あるいは（９）によれば、それぞれ異なる収音部からの定位音源方向による半直線が交差しなくても半直線双方に対して垂直となる線分の中間点によって定まる点を推定音源位置として推定することができる。
【００１８】
上述した（２）によれば、それぞれ異なる収音部からの定位音源方向により定まる中間点が分類されるクラスタの範囲内に、対応する音源の推定音源位置が分類される可能性が高くなるように推定音源位置が調整される。クラスタの範囲内に音源が存在する可能性が高くなるので、調整される推定音源位置がより正確な音源位置として得られる。
【００１９】
一般に、定位音源方向、推定音源位置及び中間点は相互に依存するが、（３）の音源位置推定部は、第１確率、第２確率及び第３確率をそれぞれ独立な推定確率の因子として推定音源位置を定めることができる。そのため、（３）の構成によれば、推定音源位置の調整に係る計算負荷が低減する。
【００２０】
上述した（４）の構成によれば、第１確率の推定音源方向の関数、第２確率の推定音源位置の関数が、それぞれ形状パラメータ、平均ならびに分散といった少数のパラメータで表される。そのため、推定音源位置の調整に係る計算負荷がさらに低減する。
【００２１】
上述した（５）の構成によれば、推定音源位置の初期値を、音源が存在する可能性が高い３個の中間点をそれぞれ頂点とする三角形の領域内に設定することができる。そのため、調整による推定音源位置の変化が収束するまでの計算負荷が低減する。
【００２２】
上述した（６）の構成によれば、スペクトルに基づいて同一と判定されなかった音源の定位音源方向の中間点に基づいて推定された推定音源位置が棄却される可能性が高くなる。
そのため、互いに異なる音源の推定音源方向の中間点に基づいて推定音源位置が虚像として誤って選択される可能性を低くすることができる。
【００２３】
上述した（７）の構成によれば、推定音源位置が定常的な音源のスペクトルが分類される第２クラスタに対応する音源の推定音源位置が選択される可能性が高くなる。即ち、推定音源位置が選択される第２クラスタには、偶発的に互いに異なる音源の推定音源方向の中間点に基づいて推定される推定音源位置が含まれる可能性が低くなる。そのため、互いに異なる音源の推定音源方向の中間点に基づいて推定音源位置が虚像として誤って選択される可能性をさらに低くすることができる。
【図面の簡単な説明】
【００２４】
【図１】本実施形態に係る音響処理システムの構成を示すブロック図である。
【図２】マイクロホンアレイ（収音部２０）の外形例を示す図である。
【図３】二次元空間における２つのマイクロホンアレイそれぞれによって定位された音源方向を示す図である。
【図４】三次元空間における２つのマイクロホンアレイそれぞれによって定位された音源方向を示す図である。
【図５】本実施形態に係る３つのマイクロホンアレイから音源位置を推定する例を示す図である。
【図６】マイクロホンアレイの配置と推定される音源方向の一例を示す図である。
【図７】各マイクロホンアレイから推定される音源方向の組に基づく中間点（含む交点）の一例を示す図である。
【図８】本実施形態に係る音源位置の推定の初期設定処理の一例を示すフローチャートである。
【図９】４個のマイクロホンアレイＭＡ１〜ＭＡ４のうち、３個のマイクロホンアレイＭＡ１〜ＭＡ３をマイクロホンアレイｍ１〜ｍ３として選択し、それぞれ推定された定位音源方向ｄ’ｍ１、ｄ’ｍ２、ｄ’ｍ３の組み合わせから推定音源位置の初期値ｘｎを定める場合を示す図である。
【図１０】本実施形態に係る確率モデルの概念図である。
【図１１】本実施形態に係る音源方向探索の説明図である。
【図１２】本実施形態に係る音源位置更新処理の一例を示すフローチャートである。
【図１３】虚像の検出例を示す図である。
【図１４】本実施形態に係る周波数分析処理の一例を示すフローチャートである。
【図１５】本実施形態に係るスコア算出処理の一例を示すフローチャートである。
【図１６】本実施形態に係る音源選択処理の一例を示すフローチャートである。
【図１７】本実施形態に係る音響処理の一例を示すフローチャートである。
【図１８】処理対象のデータ区間の例を示す図である。
【図１９】マイクロホンアレイの配置例を示す図である。
【図２０】収音した第１の音響信号例を示す図である。
【図２１】マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第１の音響信号それぞれに対して解析した結果である。
【図２２】図２０の音響信号から推定した鳥の位置を示す例である。
【図２３】第２の音響信号から推定した鳥の位置を示す例である。
【図２４】従来技術による２次元位置推定結果と本実施形態の３次元位置推定結果を比較した結果を示す図である。
【図２５】収音した第２の音響信号の１チャネルの例を示す図である。
【図２６】マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第２の音響信号それぞれに対して方位角方向を解析した結果である。
【図２７】マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第２の音響信号それぞれに対して仰角方向を解析した結果である。
【発明を実施するための形態】
【００２５】
以下、本発明の実施の形態について図面を参照しながら説明する。なお、以下の説明に用いる図面では、各部材を認識可能な大きさとするため、各部材の縮尺を適宜変更している。
【００２６】
（二次元配置のマイクロホンによるマイクロホンアレイを用いた音源定位の問題点）
ここで、二次元配置のマイクロホンによるマイクロホンアレイを用いた音源定位の問題点を説明する。
二次元の位置推定を行うマイクロホンアレイでは、マイクロホン同士の仰角に差が無い水平方向に複数のマイクロホンが配置されている。
このようなマイクロホンアレイを用いた収音した場合は、音源の高さ方向の情報を取得できない。
【００２７】
（音響処理システム）
図１は、本実施形態に係る音響処理システムＳ１の構成を示すブロック図である。図１に示すように、音響処理システムＳ１は、音響処理装置１と、Ｍ個の収音部２０と、を含んで構成される。また、収音部２０は、収音部２０−１，２０−２，・・・２０−Ｍを備えている。
なお、以下の説明において、収音する音響信号の対象は野鳥の鳴き声であり、収音場所が野外であるとする。なお、本実施形態の高さ方向も含む位置推定を三次元位置推定という。
【００２８】
音響処理装置１は、Ｍ個の収音部２０のそれぞれから取得した複数チャネルの音響信号について音源定位を行い、各音源の音源方向である定位音源方向を推定する。音響処理装置１は、Ｍ個の収音部２０のうち２個の収音部２０の組ごとに、収音部の位置から各音源の推定音源方向への半直線双方に対して垂直となる線分（二直線間の距離を表す線分）の中間点（含む交点）を定める。推定音源方向は、それぞれの収音部２０から推定される音源の方向を意味する。推定される音源の位置を推定音源位置と呼ぶ。音響処理装置１は、定めた中間点（含む交点）の分布についてクラスタリングを行い複数のクラスタに分類する。音響処理装置１は、推定音源位置が、その音源に対応するクラスタに分類される確率である推定確率が高くなるように推定音源位置を更新する。音響処理装置１が行う処理については、後述する。
【００２９】
Ｍ個の収音部２０は、それぞれ異なる位置に配置される。個々の収音部２０は、それぞれ自部に到来した音を収音し、収音した音からＱ（Ｑは、２以上の整数）チャネルの音響信号を生成する。個々の収音部２０は、例えば、所定の領域内にそれぞれ異なる位置に配置されたＱ個のマイクロホン（電気音響変換素子）を含んで構成されるマイクロホンアレイである。個々の収音部２０について、各マイクロホンが配置される領域の形状は任意である。領域の形状は、四角形、円形、球形、楕円形、など、いずれであってもよい。個々の収音部２０は、取得したＱチャネルの音響信号を音響処理装置１に出力する。個々の収音部２０は、Ｑチャネルの音響信号を無線又は有線で送信するための入出力インタフェースを備えてもよい。個々の収音部２０は一定の空間を占めるが、特に断らない限り、収音部２０の位置とは、その空間を代表する一点（例えば、重心）の位置を意味する。また、収音部２０は、Ｑチャネルの音響信号を同じサンプリング周波数の信号を用いてアナログ信号からデジタル信号に変換し、デジタル信号に変換したＱチャネルの音響信号を音響処理装置１に出力する。
なお、収音部２０をマイクロホンアレイｍと呼ぶことがある。また、個々のマイクロホンアレイｍを、マイクロホンアレイｍｋ等と、インデックスｋ等を付して区別することがある。
【００３０】
（マイクロホンアレイ）
図２は、マイクロホンアレイ（収音部２０）の外形例を示す図である。
図２に示すマイクロホンアレイは、例えば外形が球状であり、１６個のマイクロホン２０１を球体の側面上に備えている。なお、マイクロホン２０１は、二次元方向の音源方向推定に用いられる従来の水平方向に加え、垂直方向にも配置されている。
【００３１】
（音響処理装置）
次に、音響処理装置１の構成例について説明する。
図１に示すように、音響処理装置１は、入力部１０、初期処理部１２、音源位置推定部１４、音源特定部１６、出力部１８、および記憶部１９を備えている。
初期処理部１２は、音源定位部１２０、音源分離部１２２、および周波数分析部１２４を備えている。
音源位置推定部１４は、初期値設定部１３０、および音源位置更新部１４２を備える。
音源特定部１６は、分散算出部１６０、スコア算出部１６２、および音源選択部１６４を備える。
【００３２】
入力部１０は、各マイクロホンアレイｍから入力されるＱチャネルの音響信号を初期処理部１２に出力する。入力部１０は、例えば、入出力インタフェースを含んで構成される。
【００３３】
音源定位部１２０は、入力部１０から入力され、各マイクロホンアレイｍｋから取得されたＱチャネルの音響信号に基づいて音源定位を行って、各音源の方向を予め定めた長さのフレーム（例えば、１００ｍｓ）ごとに推定する。音源定位部１２０は、音源定位において、例えば、ＭＵＳＩＣ（Ｍｕｌｔｉｐｌｅ  Ｓｉｇｎａｌ  Ｃｌａｓｓｉｆｉｃａｔｉｏｎ；多重信号分類）法を用いて方向ごとのパワーを示す空間スペクトルを算出する。
音源定位部１２０は、空間スペクトルに基づいて音源ごとの音源方向を定める。音源定位部１２０は、マイクロホンアレイｍごとに定めた各音源の音源方向を示す音源方向情報と、そのマイクロホンアレイｍが取得したＱチャネルの音響信号を対応付けて音源分離部１２２に出力する。ＭＵＳＩＣ法については、後述する。
【００３４】
音源分離部１２２には、音源定位部１２０からマイクロホンアレイｍごとの音源方向情報とＱチャネルの音響信号が入力される。音源分離部１２２は、各マイクロホンアレイｍについて、Ｑチャネルの音響信号を音源方向情報が示す定位音源方向に基づいて音源ごとの成分を示す音源別音響信号に分離する。音源分離部１２２は、音源別音響信号に分離する際、例えば、ＧＨＤＳＳ（Ｇｅｏｍｅｔｒｉｃ−ｃｏｎｓｔｒａｉｎｅｄ  Ｈｉｇｈ−ｏｒｄｅｒ  Ｄｅｃｏｒｒｅｌａｔｉｏｎ−ｂａｓｅｄ  Ｓｏｕｒｃｅ  Ｓｅｐａｒａｔｉｏｎ）法を用いる。音源分離部１２２は、各マイクロホンアレイｍについて、分離した音源ごとの音源別音響信号とその音源の定位音源方向を示す音源方向情報を対応付けて周波数分析部１２４と音源位置推定部１４に出力する。ＧＨＤＳＳ法については、後述する。
【００３５】
周波数分析部１２４には、各マイクロホンアレイｍについて音源ごとの音源別音響信号と音源方向情報が対応付けて入力される。周波数分析部１２４は、個々のマイクロホンアレイｍに係る音響信号から分離された各音源の音源別音響信号を所定の時間長（例えば、１２８点）のフレームごとに周波数分析を行ってスペクトル［Ｆｍ，１］、［Ｆｍ，２］〜［Ｆｍ，ｓｍ］を算出する。［…］は、ベクトル、行列など複数の値からなるセットを示す。ｓｍは、マイクロホンアレイｍが取得した音響信号から音源定位ひいては音源分離により推定された音源の音源数を示す。ここで、スペクトル［Ｆｍ，１］、［Ｆｍ，２］〜［Ｆｍ，ｓｍ］は、それぞれ行ベクトルである。周波数分析において、周波数分析部１２４は、例えば、各音源別音響信号に１２８点のハミング窓を作用して得られる信号に短時間フーリエ変換（ＳＴＦＴ：Ｓｈｏｒｔ  Ｔｅｒｍ  Ｆｏｕｒｉｅｒ  Ｔｒａｎｓｆｏｒｍ）を行う。周波数分析部１２４は、時間的に隣接するフレームを重複させ、分析対象の区間をなすフレームを逐次にシフトさせる。
周波数分析部１２４は、音源ごとのスペクトルを行間で統合して次式（１）に示すマイクロホンアレイｍごとのスペクトル行列［Ｆｍ］（ｍは、１からＭまでの間の整数）を構成する。周波数分析部１２４は、構成したスペクトル行列［Ｆ１］、［Ｆ２］、〜［ＦＭ］を、さらに行間で統合して次式（２）に示すスペクトル行列［Ｆ］を構成する。周波数分析部１２４は、構成したスペクトル行列［Ｆ］と、各音源の定位音源方向を示す音源方向情報とを対応付けて音源特定部１６に出力する。
【００３６】
【数１】
2020038123.tif 000003
【００３７】
【数２】
2020038123.tif 000004
【００３８】
初期値設定部１４０は、音源分離部１２２から入力されるマイクロホンアレイｍごとの音源方向情報に基づいて三角分割法（ｔｒｉａｎｇｕｌａｔｉｏｎ）を用いて音源の候補として推定される位置である推定音源位置の初期値を定める。三角分割法は、Ｍ個のうち３個のマイクロホンアレイの組から定められ、ある音源の候補に係る３個の中間点の重心を、その音源の推定音源位置の初期値として定める手法である。以下の説明では、音源の候補を音源候補と呼ぶ。中間点は、３個のマイクロホンアレイｍのうち２個のマイクロホンアレイｍの組ごとに、各マイクロホンアレイｍの位置を通り、そのマイクロホンアレイｍが取得された音響信号に基づいて推定された定位音源方向への半直線双方に対して垂直となる線分の中間点（半直線同士の交点も含む）である。初期値設定部１４０は、音源候補ごとの推定音源位置の初期値を示す初期推定音源位置情報を音源位置更新部１４２に出力する。初期値設定処理の例については、後述する。
【００３９】
音源位置更新部１４２は、各２個のマイクロホンアレイｍの組ごとに、それぞれのマイクロホンアレイｍから、そのマイクロホンアレイｍに基づく定位音源方向に係る音源候補の推定音源方向への半直線双方に対して垂直になる線分の中間点を定める。推定音源方向とは、推定音源位置への方向を意味する。音源位置更新部１４２は、定めた中間点の空間分布についてクラスタリングを行い複数のクラスタ（群）に分類する。音源位置更新部１４２は、音源候補ごとの推定音源位置がそれぞれの音源候補に対応するクラスタに分類される確率である推定確率が高くなるように、その推定音源位置を更新する。
音源位置更新部１４２は、音源候補ごとの推定音源位置の初期値として、初期値設定部１４０から入力される初期推定音源位置情報が示す推定音源位置の初期値を用いる。音源位置更新部１４２は、推定音源位置もしくは推定音源方向の更新量が所定の更新量の閾値未満となったとき、推定音源位置もしくは推定音源方向の変化が収束したと判定し、推定音源位置の更新を停止する。音源位置更新部１４２は、音源候補ごとの推定音源位置を示す推定音源位置情報を音源特定部１６に出力する。更新量が所定の更新量の閾値以上であるとき、音源位置更新部１４２は、音源候補ごとの推定音源位置を更新する処理を継続する。推定音源位置の更新処理の例については、後述する。
【００４０】
音源特定部１６は、分散算出部１６０と、スコア算出部１６２と、音源選択部１６４と、を含んで構成される。
分散算出部１６０には、周波数分析部１２４からスペクトル行列［Ｆ］と音源方向情報が入力され、音源位置推定部１４から推定音源位置情報が入力される。
分散算出部１６０は、次に説明する処理を所定の回数繰り返す。繰り返し回数Ｒは、予め分散算出部１６０に設定しておく。
【００４１】
分散算出部１６０は、スペクトル行列［Ｆ］が示す収音部２０ごとの各音源のスペクトルについてクラスタリングを行い、複数のクラスタ（群）に分類する。分散算出部１６０が実行するクラスタリングは、音源位置更新部１４２が実行するクラスタリングと独立である。分散算出部１６０は、クラスタリングの手法として、例えば、ｋ−平均法（ｋ−ｍｅａｎｓ  ｃｌｕｓｔｅｒｉｎｇ）を用いる。ｋ−平均法では、クラスタリングの対象とする複数のデータのそれぞれをランダムにｋ個のクラスタに割り当てる。分散算出部１６０は、各繰り返し回数ｒにおいてスペクトルごとの初期値として、割り当てられるクラスタを変更する。以下の説明では、分散算出部１６０が分類したクラスタを第２クラスタと呼ぶ。分散算出部１６０は、第２クラスタのそれぞれに属する複数のスペクトルの類似度を示す指標値を算出する。分散算出部１６０は、算出した指標値が所定の類似度を示す指標値よりも高いか否かにより、各スペクトルに係る音源候補が同一であるか否かを判定する。
【００４２】
音源候補が同一と判定した第２クラスタに対応する音源候補について、分散算出部１６０は、その推定音源位置情報が示すその音源候補の推定音源位置の分散を算出する。後述するように、この段階では、第２クラスタの個数よりも、音源位置更新部１４２が音源位置を更新する音源候補の個数の方が多くなる可能性があるためである。分散算出部１６０は、例えば、第２クラスタについて現在の繰り返し回数ｒにおいて算出した分散が、前回の繰り返し回数ｒ−１において算出した分散より大きいとき、スコアを０とする。分散算出部１６０は、その第２クラスタについて現在の繰り返し回数ｒにおいて算出した分散が、前回の繰り返し回数ｒ−１において算出した分散と等しいか、より小さいとき、スコアをεとする。εは、例えば、所定の正の実数である。分散の増加の頻度が多いほど、第２クラスタに分類される推定音源位置が繰り返し回数により異なる、つまり、第２クラスタとしての安定性が低くなる。言い換えれば、設定されるスコアは、第２クラスタの安定性を示す。音源選択部１６４において、スコアが高い第２クラスタほど対応する音源候補の推定音源位置が優先して選択される。
【００４３】
他方、音源候補が同一ではないと判定した第２クラスタについて、分散算出部１６０は、対応する音源候補がないと判定し、その推定音源位置の分散が有効ではないと判定し、スコアをδとする。δは、例えば、０より小さい負の実数である。これにより、音源選択部１６４において、音源候補が同一と判定した音源候補に係る推定音源位置が、同一と判定しなかった音源候補よりも優先して選択される。
分散算出部１６０は、第２クラスタごとの各繰り返し回数のスコアと推定音源位置を示すスコア算出情報をスコア算出部１６２に出力する。
【００４４】
スコア算出部１６２は、分散算出部１６０から入力されるスコア算出情報に基づいて第２クラスタに対応する音源候補ごとの最終スコアを算出する。ここで、スコア算出部１６２は、第２クラスタごとに有効な分散を定めた回数である有効を計数し、各回のスコアの合計値を算出する。スコアの合計値は、各回で分散が増加する回数である有効回数が多いほど大きくなる。即ち、第２クラスタの安定性が高いほど、スコアの合計値が大きくなる。なお、この段階では、１個の推定音源位置が複数の第２クラスタにまたがる場合がある。そこで、スコア算出部１６２は、推定音源位置ごとのスコアの合計値の総和を、計数した有効回数の総和で除算してその推定音源位置に対応する音源候補の最終スコアを算出する。スコア算出部１６２は、算出した音源候補の最終スコアと推定音源位置を示す最終スコア情報を音源選択部１６４に出力する。
【００４５】
音源選択部１６４は、スコア算出部１６２から入力される最終スコア情報が示す音源候補の最終スコアが、所定の最終スコアの閾値θ２以上となる音源候補を音源として選択する。音源選択部１６４は、最終スコアが、閾値θ２未満となる音源候補を棄却する。音源選択部１６４は、選択した音源について、音源ごとの推定音源位置を示す出力音源位置情報を出力部１８に出力する。
【００４６】
出力部１８は、音源選択部１６４から入力される出力音源位置情報を、音響処理装置１の外部に出力する。出力部１８は、例えば、入出力インタフェースを含んで構成される。出力部１８と入力部１０とは、共通のハードウェアで構成されてもよい。出力部１８は、出力音源位置情報を表示する表示部（例えば、ディスプレイ）を備えてもよい。音響処理装置１は、出力部１８とともに、又は出力部１８に代えて、出力音源位置情報を記憶する記憶媒体を含んで構成されてもよい。
【００４７】
記憶部１９は、予め所定の間隔（例えば、５°）で分布した方向ψごとの伝達関数を記憶する。記憶部１９は、各マイクロホンアレイが収音した音響信号を記憶しておいてもよい。記憶部１９は、後述するように所定期間の音響信号と導出されるデータを記憶する。
【００４８】
（ＭＵＳＩＣ法）
ここで、音源定位の一手法であるＭＵＳＩＣ法について説明する。
ＭＵＳＩＣ法は、以下に説明する空間スペクトルのパワーＰｅｘｔ（ψ）が極大であって、所定のレベルよりも高い方向ψを定位音源方向として定める手法である。音源定位部１２０が備える記憶部１９には、予め所定の間隔（例えば、５°）で分布した方向ψごとの伝達関数を記憶させておく。本実施形態では、次に説明する処理をマイクロホンアレイｍ毎に実行する。
【００４９】
音源定位部１２０は、音源から各チャネルｑ（ｑは、１以上Ｑ以下の整数）に対応するマイクロホンまでの伝達関数Ｄ［ｑ］（ω）を要素とする伝達関数ベクトル［Ｄ（ψ）］を方向ψごとに生成する。
音源定位部１２０は、各チャネルｑの音響信号ξｑを所定の要素数からなるフレームごとに周波数領域に変換することによって変換係数ξｑ（ω）を算出する。音源定位部１２０は、算出した変換係数を要素として含む入力ベクトル［ξ（ω）］から次式（３）に示す入力相関行列［Ｒξξ］を算出する。
【００５０】
【数３】
2020038123.tif 000005
【００５１】
式（３）において、Ｅ［…］は、…の期待値を示す。［…］は、…が行列又はベクトルであることを示す。［…］＊は、行列又はベクトルの共役転置（ｃｏｎｊｕｇａｔｅ  ｔｒａｎｓｐｏｓｅ）を示す。
音源定位部１２０は、入力相関行列［Ｒξξ］の固有値δｐ及び固有ベクトル［εｐ］を算出する。入力相関行列［Ｒξξ］、固有値δｐ、及び固有ベクトル［ξｐ］は、次式（４）に示す関係を有する。
【００５２】
【数４】
2020038123.tif 000006
【００５３】
式（４）において、ｐは、１以上Ｑ以下の整数である。インデックスｐの順序は、固有値δｐの降順である。
音源定位部１２０は、伝達関数ベクトル［Ｄ（ψ）］と算出した固有ベクトル［εｐ］に基づいて、次式（５）に示す周波数別空間スペクトルのパワーＰｓｐ（ψ）を算出する。
【００５４】
【数５】
2020038123.tif 000007
【００５５】
式（５）において、Ｄｍは、検出可能とする音源の最大個数（例えば、２）であって、Ｑよりも小さい予め定めた自然数である。
音源定位部１２０は、Ｓ／Ｎ比が予め定めた閾値（例えば、２０ｄＢ）よりも大きい周波数帯域における空間スペクトルＰｓｐ（ψ）の総和を全帯域の空間スペクトルのパワーＰｅｘｔ（ψ）として算出する。
【００５６】
なお、音源定位部１２０は、ＭＵＳＩＣ法に代えて、その他の手法を用いて定位音源方向を算出してもよい。例えば、重み付き遅延和ビームフォーミング（ＷＤＳ−ＢＦ：Ｗｅｉｇｈｔｅｄ  Ｄｅｌａｙ  ａｎｄ  Ｓｕｍ  Ｂｅａｍ  Ｆｏｒｍｉｎｇ）法が利用可能である。ＷＤＳ−ＢＦ法は、次式（６）に示すように各チャネルｑの全帯域の音響信号ξｑ（ｔ）の遅延和の二乗値を空間スペクトルのパワーＰｅｘｔ（ψ）として算出し、空間スペクトルのパワーＰｅｘｔ（ψ）が極大となる定位音源方向ψを探索する手法である。
【００５７】
【数６】
2020038123.tif 000008
【００５８】
式（６）において［Ｄ（ψ）］の各要素が示す伝達関数は、音源から各チャネルｑ（ｑは、１以上Ｑ以下の整数）に対応するマイクロホンまでの位相の遅延による寄与を示す。［ξ（ｔ）］は、時刻ｔの時点における各チャネルｑの音響信号ξｑ（ｔ）の信号値を要素とするベクトルである。
【００５９】
（ＧＨＤＳＳ法）
次に、音源分離の一手法であるＧＨＤＳＳ法について説明する。
ＧＨＤＳＳ法は、２つのコスト関数（ｃｏｓｔ  ｆｕｎｃｔｉｏｎ）として、分離尖鋭度（Ｓｅｐａｒａｔｉｏｎ  Ｓｈａｒｐｎｅｓｓ）ＪＳＳ（［Ｖ（ω）］）と幾何制約度（Ｇｅｏｍｅｔｒｉｃ  Ｃｏｎｓｔｒａｉｎｔ）ＪＧＣ（［Ｖ（ω）］）が、それぞれ減少するように分離行列［Ｖ（ω）］を適応的に算出する方法である。本実施形態では、各マイクロホンアレイｍが取得した音響信号のそれぞれから音源別音響信号を分離する。
【００６０】
分離行列［Ｖ（ω）］は、音源定位部１２０から入力されたＱチャネルの音響信号［ξ（ω）］に乗じることによって、検出される最大Ｄｍ個の音源それぞれの音源別音響信号（推定値ベクトル）［ｕ’（ω）］を算出するために用いられる行列である。ここで、［…］Ｔは、行列又はベクトルの転置を示す。
【００６１】
分離尖鋭度ＪＳＳ（［Ｖ（ω）］）、幾何制約度ＪＧＣ（［Ｖ（ω）］）は、それぞれ、次式（７）、（８）のように表される。
【００６２】
【数７】
2020038123.tif 000009
【００６３】
【数８】
2020038123.tif 000010
【００６４】
式（７）、（８）において、｜｜…｜｜２は、行列…のフロベニウスノルム（Ｆｒｏｂｅｎｉｕｓ  ｎｏｒｍ）である。フロベニウスノルムとは、行列を構成する各要素値の二乗和（スカラー値）である。φ（［ｕ’（ω）］）は、音源別音響信号［ｕ’（ω）］の非線形関数、例えば、双曲線正接関数（ｈｙｐｅｒｂｏｌｉｃ  ｔａｎｇｅｎｔ  ｆｕｎｃｔｉｏｎ）である。ｄｉａｇ［…］は、行列…の対角成分の総和を示す。従って、分離尖鋭度ＪＳＳ（［Ｖ（ω）］）は、音源別音響信号（推定値）のスペクトルのチャネル間非対角成分の大きさ、つまり、ある１つの音源が他の音源として誤って分離される度合いを表す指標値である。また、式（８）において、［Ｉ］は、単位行列を示す。従って、幾何制約度ＪＧＣ（［Ｖ（ω）］）とは、音源別音響信号（推定値）のスペクトルと音源別音響信号（音源）のスペクトルとの誤差の度合いを表す指標値である。
【００６５】
（マイクロホンアレイによる中間点（含む交点））
次に、マイクロホンアレイによる中間点（含む交点）について説明する。
まず、二次元空間におけるマイクロホンアレイによる交点について説明する。
図３は、二次元空間における２つのマイクロホンアレイそれぞれによって定位された音源方向を示す図である。図３において、横方向をｘ軸方向、縦方向をｙ軸方向とする。なお、ｘｙ平面は、例えば地面に水平な平面である。
この例では、マイクロホンアレイＭＡ１とＭＡ２の位置を通り、それぞれのマイクロホンアレイが取得した音響信号に基づいて推定された定位音源方向への半直線が定められる。これらの２本の半直線は、音源Ｓの位置において一点に交わる。
【００６６】
次に、三次元平面におけるマイクロホンアレイによる中間点（含む交点）について説明する。
図４は、三次元空間における２つのマイクロホンアレイそれぞれによって定位された音源方向を示す図である。図４において、横方向をｘ軸方向、縦方向をｙ軸方向、高さ方向をｚ軸方向とする。なお、ｘｙ平面は、例えば地面に水平な水平平面である。
【００６７】
本実施形態では、上述したように収音対象が野鳥の鳴き声であり、収音場所が野外であるため、反響音の影響が室内での収音と比べて少ない。このため、本実施形態では、推定音源位置を、互いに最も近い２つの方位線上の点として求める。
この例では、マイクロホンアレイＭＡ１とＭＡ２の位置を通り、それぞれのマイクロホンアレイが取得した音響信号に基づいて推定された定位音源方向への半直線が定められる。これらの２本の半直線は、三次元空間において、１点に交わらず互いに平行な場合がある。しかしながら、実際の環境で収音された音響信号を用いて解析する場合は、ある時間に単一の音源を示す全ての方向が互いに歪んだ線である仮定することができる。
このため、本実施形態では、推定音源位置を最短距離ｄの中点Ｐを交点とみなす。
【００６８】
マイクロホンアレイＭＡ１とＭＡ２の２つの方位角線上にある２つの点ＰとＱとを、次式（９）と次式（１０）のように表す。
【００６９】
【数９】
2020038123.tif 000011
【００７０】
【数１０】
2020038123.tif 000012
【００７１】
式（９）と式（１０）において、ｍｊとｍｋは２つのマイクロホンアレイの（ｘ，ｙ，ｚ）座標、ａｊとａｋは２つの方向線の方向余弦（単位ベクトル）、ｒｊとｒｋは対応するマイクロホンから点ＰとＱまでの距離である。
式（９）と式（１０）を用いて、点Ｐと点Ｑとのベクトルは、次式（１１）のように表される。
【００７２】
【数１１】
2020038123.tif 000013
【００７３】
ここで、点Ｐと点Ｑを通る２つの線の間において最も距離が近い線は、両方の線に直角の線である。
そして、ｍｊ、ｍｋ、ａｊ、ａｋの値は既知であるので、距離ｒｊとｒｋは、次式（１２）のように、ＰＱベクトルを用いて両方の線の方向ベクトルの内積を計算することで求めることができる。
【００７４】
【数１２】
2020038123.tif 000014
【００７５】
音源位置を推定するために、マイクロホンアレイから得られた各方向線のペアについて上記の計算を行った。
推定音源位置は、複数のマイクロホンアレイのペアの点Ｐと点Ｑを通る２つの線の間において最も距離が短い垂線の中間点の重心となる。このため、各ペアのＰとＱの中間点と、中間点の重心を次式（１３）によって計算する。
【００７６】
【数１３】
2020038123.tif 000015
【００７７】
式（１３）は、図５のように一度に３つのマイクロホンアレイから得られたデータに対して実行する。図５は、本実施形態に係る３つのマイクロホンアレイから音源位置を推定する例を示す図である。図５の座標系は図４と同じである。マイクロホンアレイＭＡ１とＭＡ２による両方の線に直角の線の中心点がＣ１である。マイクロホンアレイＭＡ２とＭＡ３による両方の線に直角の線の中心点がＣ２である。マイクロホンアレイＭＡ１とＭＡ３による両方の線に直角の線の中心点がＣ３である。また、定位音源方向ｄ’ｍ１とｄ’’ｍ１は、マイクロホンアレイＭＡ１による定位結果である。定位音源方向ｄ’ｍ２とｄ’’ｍ２は、マイクロホンアレイＭＡ２による定位結果である。定位音源方向ｄ’ｍ３とｄ’’ｍ３は、マイクロホンアレイＭＡ３による定位結果である。
【００７８】
なお、３つ以上のマイクロホンアレイを使用する場合でも、３つのマイクロホンアレイのすべての組み合わせに対して式（１３）を実施して例えばそれらの重心を求めることで、推定音源位置ｘを求める。
【００７９】
（初期値の設定）
次に、初期値設定部１４０（図１）が行う初期値の設定の例について説明する。各２個のマイクロホンアレイｍに基づいて定められる中間点（含む交点）は、理想的には各音源の音源位置と等しくなるはずである。図６を参照して、互いに異なる位置に設置されたマイクロホンアレイＭＡ１、ＭＡ２、ＭＡ３のそれぞれが取得した音響信号に基づいて音源Ｓの定位音源方向が推定される場合を例にする。図６は、マイクロホンアレイの配置と推定される音源方向の一例を示す図である。この例では、マイクロホンアレイＭＡ１、ＭＡ２、ＭＡ３の位置を通り、それぞれのマイクロホンアレイが取得した音響信号に基づいて推定された定位音源方向への半直線が定められる。これらの３本の半直線は、理想的には音源Ｓの位置において一点に交わる。
【００８０】
しかしながら、音源Ｓの定位音源方向には誤差が含まれる。現実的には、図７に示すように１つの音源に係る中間点（含む交点）Ｐ１、Ｐ２、Ｐ３の位置が互いに異なる。図７は、各マイクロホンアレイから推定される音源方向の組に基づく中間点の一例を示す図である。中間点Ｐ１は、マイクロホンアレイＭＡ１、ＭＡ２の位置を通り、それぞれのマイクロホンアレイＭＡ１、ＭＡ２が取得した音響信号から推定された音源Ｓの定位音源方向の半直線双方に対して垂直になる線分の中間点である。中間点Ｐ２は、マイクロホンアレイＭＡ２、ＭＡ３の位置を通り、それぞれのマイクロホンアレイＭＡ２、ＭＡ３が取得した音響信号から推定された音源Ｓの定位音源方向の半直線双方に対して垂直になる線分の中間点である。中間点Ｐ３は、マイクロホンアレイＭＡ１、ＭＡ３の位置を通り、それぞれのマイクロホンアレイＭＡ１、ＭＡ３が取得した音響信号から推定された音源Ｓの定位音源方向の半直線双方に対して垂直になる線分の中間点である。同一の音源Ｓについて、各マイクロホンアレイが取得した音響信号から推定される定位音源方向の誤差がランダムであれば、真の音源位置は、中間点Ｐ１、Ｐ２、Ｐ３のそれぞれを頂点とする三角形の内部の領域にあることが期待される。そこで、初期値設定部１４０は、中間点Ｐ１、Ｐ２、Ｐ３間の重心を、音源Ｓの候補である音源候補の推定音源位置の初期値ｘｎとして定める。
【００８１】
ただし、音源定位部１２０が各マイクロホンアレイｍから取得した音響信号から推定する音源方向の数は、１個には限らず、複数になることがある。そのため、中間点Ｐ１、Ｐ２、Ｐ３は、互いに同一の音源Ｓの方向に基づいて定められるとは限らない。そこで、初期値設定部１４０は、３個の中間点Ｐ１、Ｐ２、Ｐ３のうち、各２個の中間点間の距離Ｌ１２、Ｌ２３、Ｌ１３が、いずれも予め定めた距離の閾値θ１未満であるか、少なくとも中間点間の距離のいずれかが、その閾値θ１以上となる距離が存在するか否かを判定する。いずれも閾値θ１未満と判定するとき、初期値設定部１４０は、それらの中間点Ｐ１、Ｐ２、Ｐ３の重心を音源候補ｎの音源位置の初期値ｘｎとして採用する。初期値設定部１４０は、少なくとも中間点間の距離のいずれかが、その閾値θ１以上となる場合、中間点Ｐ１、Ｐ２、Ｐ３の重心を音源位置の初期値ｘｎとして定めずに、棄却する。
【００８２】
ここで、音源位置推定部１４には、Ｍ個のマイクロホンアレイＭＡ１，ＭＡ２，…，ＭＡＭのそれぞれの位置ｕＭＡ１，ｕＭＡ２，…，ｕＭＡＭを、予め設定させておく。個々のマイクロホンアレイｍの位置ｕＭＡ１，ｕＭＡ２，…，ｕＭＡＭを要素とする位置ベクトル［ｕ］は、次式（１４）で表わされる。
【００８３】
【数１４】
2020038123.tif 000016
【００８４】
式（１４）において、マイクロホンアレイｍの位置ｕＭＡｍ（ｍは、１からＭの間の整数）は、ｘ座標ｕＭＡｘｍ、ｙ座標ｕＭＡｙｍを要素値とする２次元の座標［ｕＭＡｘｍ，ｕＭＡｙｍ］である。
上述したように、音源定位部１２０は、各マイクロホンアレイＭＡｍが取得したＱチャネルの音響信号から、それぞれ最大Ｄｍ個の定位音源方向ｄ’ｍ（１），ｄ’ｍ（２），…，ｄ’ｍ（Ｄｍ）をフレームごとに定める。定位音源方向ｄ’ｍ（１），ｄ’ｍ（２），…，ｄ’ｍ（Ｄｍ）を要素とするベクトル［ｄ’］は、次式（１５）で表わされる。
【００８５】
【数１５】
2020038123.tif 000017
【００８６】
次に、本実施形態に係る音源位置の推定の初期設定処理の一例について説明する。
図８は、本実施形態に係る音源位置の推定の初期設定処理の一例を示すフローチャートである。
（ステップＳ１０１）初期値設定部１４０は、三角分割法においてＭ個のマイクロホンアレイから互いに異なる３個のマイクロホンアレイｍ１、ｍ２、ｍ３の組（ｔｒｉｐｌｅｔ）を選択する。処理後、ステップＳ１０２の処理に進む。
【００８７】
（ステップＳ１０２）初期値設定部１４０は、選択した３個の組のマイクロホンアレイｍ１、ｍ２、ｍ３のそれぞれについて、それぞれのマイクロホンアレイが取得した音響信号に基づいて推定された最大Ｄｍ個の音源から各１個の音源δ１、δ２、δ３の定位音源方向ｄ’ｍ１（δ１）、ｄ’ｍ２（δ２）、ｄ’ｍ３（δ３）を選択する。続けて、初期値設定部１４０は、３個のマイクロホンアレイのうち各２つのマイクロホンアレイの組（対；ｐａｉｒ）について、式（１３）を用いて中間点Ｐ１、Ｐ２、Ｐ３の重心を算出する。ここで、マイクロホンアレイｍ１とｍ２による両方の線に直角の線の中心点がＰ１、マイクロホンアレイｍ２とｍ３による両方の線に直角の線の中心点がＰ２、マイクロホンアレイｍ１とｍ３による両方の線に直角の線の中心点がＰ３である。次式（１６）に示すように、中間点Ｐ１は、マイクロホンアレイｍ１、ｍ２の位置と、定位音源方向ｄ’ｍ１（δ１）、ｄ’ｍ２（δ２）により定まる。中間点Ｐ２は、マイクロホンアレイｍ２、ｍ３の位置と、定位音源方向ｄ’ｍ２（δ２）、ｄ’ｍ３（δ３）により定まる。中間点ＰＰ３は、マイクロホンアレイｍ１、ｍ３の位置と、定位音源方向ｄ’ｍ１（δ１）、ｄ’ｍ３（δ３）により定まる。なお、中間点には、２つのマイクロホンによる定位音源方向への半直線同士の交点も含まれる。その後、ステップＳ１０３の処理に進む。
【００８８】
【数１６】
2020038123.tif 000018
【００８９】
（ステップＳ１０３）初期値設定部１４０は、マイクロホンアレイの組ごとに特定した音源ごとの中心点の重心を推定音源位置の初期値として定める。
算出した距離Ｌ１２、Ｌ２３、Ｌ１３がいずれも閾値θ１以下となる場合、初期値設定部１４０は、３個の中間点Ｐ１、Ｐ２、Ｐ３の組み合わせを、音源候補ｎに係る組み合わせとして選択する。その場合、初期値設定部１４０は、次式（１７）に示すように、中間点Ｐ１、Ｐ２、Ｐ３の重心を音源候補ｎの音源推定位置の初期値ｘｎとして定める。
他方、距離Ｌ１２、Ｌ２３、Ｌ１３の少なくともいずれか１つが閾値θ１より大きいとなる場合、初期値設定部１４０は、これらの中間点の組み合わせを棄却し、初期値ｘｎを定めない。式（１７）において、φは空集合を示す。その後、図６に示す処理を終了する。
【００９０】
【数１７】
2020038123.tif 000019
【００９１】
初期値設定部１４０は、マイクロホンアレイｍ１、ｍ２、ｍ３ごとに推定される定位音源方向の組み合わせｄ’ｍ１（δ１）、ｄ’ｍ２（δ２）、ｄ’ｍ３（δ３）ごとに、ステップＳ１０１〜Ｓ１０３の処理を実行する。これにより、音源候補として不適切な中間点の組み合わせが棄却され、音源候補ｎごとに音源推定位置の初期値ｘｎが定められる。なお、以下の説明では音源候補数を、Ｎで表す。
また、初期値設定部１４０は、Ｍ個のマイクロホンアレイのうち、３個のマイクロホンアレイの組ごとに、ステップＳ１０１〜Ｓ１０３の処理を実行してもよい。これにより、音源の候補ｎの検出漏れを少なくすることができる。
【００９２】
図９は、４個のマイクロホンアレイＭＡ１〜ＭＡ４のうち、３個のマイクロホンアレイＭＡ１〜ＭＡ３をマイクロホンアレイｍ１〜ｍ３として選択し、それぞれ推定された定位音源方向ｄ’ｍ１、ｄ’ｍ２、ｄ’ｍ３の組み合わせから推定音源位置の初期値ｘｎを定める場合を示す図である。
【００９３】
中間点Ｐ１の方向は、それぞれマイクロホンアレイｍ１、ｍ２の位置を基準とする定位音源方向ｄ’ｍ１、ｄ’ｍ２と同一の方向となる。中間点Ｐ２の方向は、それぞれマイクロホンアレイｍ２、ｍ３の位置を基準とする音源方向ｄ’ｍ２、ｄ’ｍ３と同一の方向となる。中間点Ｐ３の方向は、それぞれマイクロホンアレイｍ１、ｍ３の位置を基準とする定位音源方向ｄ’ｍ１、ｄ’ｍ３と同一の方向となる。定められた初期値ｘｎの方向は、それぞれマイクロホンアレイｍ１、ｍ２、ｍ３の位置を基準とする方向ｄ”ｍ１、ｄ”ｍ２、ｄ”ｍ３となる。よって、音源定位により推定される定位音源方向ｄ’ｍ１、ｄ’ｍ２、ｄ’ｍ３が、それぞれ推定音源方向ｄ”ｍ１、ｄ”ｍ２、ｄ”ｍ３に修正される。なお、定位音源方向ｄ’ｍ１（δ１）、ｄ’ｍ２（δ２）、ｄ’ｍ３（δ３）を要素とする方向ベクトル［ｄ”］は、次式（１８）で表される。なお、δ１、δ２、δ３は、それぞれ１からＤｍの間の整数である。
【００９４】
【数１８】
2020038123.tif 000020
【００９５】
（推定音源位置の更新処理）
次に、推定音源位置の更新処理について説明する。音源定位により推定される音源方向は誤差を含むため、音源方向間の中間点から推定される候補音源ごとの推定音源位置も誤差を含む。これらの誤差がランダムであれば、推定音源位置ならびに中間点は、各音源の真の音源位置の周囲に分布することが期待される。そこで、本実施形態に係る音源位置更新部１４２は、各２個のマイクロホンアレイ、推定音源方向間の中間点についてクラスタリングを行い、これらの中間点の分布を複数のクラスタに分類する。ここで、推定音源方向とは、推定音源位置の方向を意味する。クラスタリングの手法として、音源位置更新部１４２は、例えば、ｋ−平均法を用いる。音源位置更新部１４２は、音源候補ごとの推定音源位置がそれぞれの音源候補に対応するクラスタに分類される可能性の度合いである推定確率が高くなるように、その推定音源位置を更新する。
【００９６】
（確率モデル）
推定音源位置を算出する際、音源位置更新部１４２は、三角分割法に基づく確率モデルを用いる。この確率モデルでは、音源候補ごとの推定音源位置がそれぞれの音源候補に対応するクラスタに分類される推定確率が、第１確率と、第２確率と、第３確率と、をそれぞれ因子とする積で表されるように分解されるように近似できるものと仮定する。第１確率は、音源定位により定位音源方向が定められるとき、その音源に対応する音源候補の推定音源位置の方向である推定音源方向が得られる確率である。第２確率は、２つのマイクロホンアレイそれぞれの位置からその推定音源方向への半直線の中間点が定められるとき、その推定音源位置が得られる確率である。第３確率は、その中間点の分類されるクラスタへの出現確率である。
【００９７】
より具体的には、第１確率は、それぞれ定位音源方向ｄ’ｍｊ、ｄ’ｍｋを基準とするフォン・ミーゼス分布（ｖｏｎ−Ｍｉｓｅｓ  ｄｉｓｔｒｉｂｕｔｉｏｎ）に従うものと仮定する。つまり、第１確率は、音源定位により各マイクロホンアレイｍｊ、ｍｋが取得される音響信号から推定される定位音源方向ｄ’ｍｊ、ｄ’ｍｋに、確率分布がフォン・ミーゼス分布となる誤差が含まれるとの仮定に基づく。理想的には、図１０に示す例では、誤差がなければ、定位音源方向ｄ’ｍｊ、ｄ’ｍｋとして真の音源方向ｄｍｊ、ｄｍｋが得られる。図１０は、本実施形態に係る確率モデルの概念図である。
【００９８】
第２確率は、マイクロホンアレイｍｊ、ｍｋ、推定音源方向ｄｍｊ、ｄｍｋ間の中間点ｓｊ，ｋの位置を基準とする多次元ガウス関数に従うものと仮定する。つまり、第２確率は、各マイクロホンアレイｍｊ、ｍｋのそれぞれを通り、それぞれの方向が推定音源方向ｄｍｊ、ｄｍｋとなる半直線の中間点ｓｊ，ｋとなる推定音源位置に、確率分布が多次元ガウス分布となる誤差としてガウス雑音が含まれているとの仮定に基づく。理想的には、中間点ｓｊ，ｋの座標が多次元ガウス関数の平均値μｃｊ，ｋとなる。
従って、音源位置更新部１４２は、音源定位により得られた定位音源方向ｄ’ｍｊ、ｄ’ｍｋに基づいて、音源候補の推定音源方向を与える中間点ｓｊ，ｋの座標が、中間点ｓｊ，ｋの分布を近似する多次元ガウス関数の平均値μｃｊ，ｋに極力近づくように推定音源方向ｄｍｊ、ｄｍｋを推定する。
【００９９】
第３確率は、マイクロホンアレイｍｊ、ｍｋのそれぞれを通り、それぞれの方向が推定音源方向ｄｍｊ、ｄｍｋとなる半直線の中間点ｓｊ，ｋが分類されるクラスタｃｊ，ｋの出現確率を示す。つまり、第３確率は、その中間点ｓｊ，ｋに相当する推定音源位置のクラスタｃｊ，ｋへの出現確率を示す。
各クラスタと音源を対応付けるため、音源位置更新部１４２は、音源候補ごとの推定音源位置ｘｎの初期値について、初期クラスタリング（ｉｎｉｔｉａｌ  ｃｌｕｓｔｅｒｉｎｇ）を行ってクラスタの個数Ｃを定める。
【０１００】
初期クラスタリングでは、音源位置更新部１４２は、次式（１９）に示すように、音源候補ごとの推定音源位置ｘｎについて所定のユークリッド距離の閾値φをパラメータとしてそれぞれ用いて階層クラスタリング（ｈｉｅｒａｒｃｈｉｃａｌ  ｃｌｕｓｔｅｒｉｎｇ）を行って複数のクラスタに分類する。階層クラスタリングとは、１個の対象データだけを含む複数のクラスタを初期状態として生成し、それぞれ異なる対応データを含む２つのクラスタ間のユークリッド距離を算出し、算出したユークリッド距離が最も小さいクラスタ同士を逐次に併合して、新たなクラスタを形成する手法である。クラスタを併合する処理は、ユークリッド距離が閾値φに達するまで繰り返す。閾値φとして、例えば、音源位置の推定誤差よりも大きい値を予め設定しておけばよい。従って、閾値φより距離が小さい複数の音源候補同士が１つのクラスタに集約され、それぞれのクラスタが音源に対応付けられる。そして、クラスタリングにより得られるクラスタの数Ｃが音源数として推定される。
【０１０１】
【数１９】
2020038123.tif 000021
【０１０２】
式（１９）において、ｈｉｅｒａｒｃｈｙとは、階層クラスタリングを示す。ｃｎは、クラスタリングに得られる各クラスタのインデックスｃｎを示す。ｍａｘ（…）は、…の最大値を示す。
【０１０３】
次に、確率モデルの適用例について説明する。上述したように、各マイクロホンアレイｍｉについて、定位音源方向ｄ’ｍｉが定められるとき推定音源方向ｄｍｉが得られる第１確率（ｄ’ｍｉ，ｄｍｉ；βｍｉ）は、次式（２０）に示すフォン・ミーゼス分布に従うものと仮定する。
【０１０４】
【数２０】
2020038123.tif 000022
【０１０５】
フォン・ミーゼス分布は、最大値、最小値を、それぞれ１、０とする連続関数であり、定位音源方向ｄ’ｍｉと推定音源方向ｄｍｉが等しいときに最大値１をとり、定位音源方向ｄ’ｍｉと推定音源方向ｄｍｉのなす角が大きいほど関数値が小さくなる。式（２０）において、音源方向ｄ’ｍｉ、推定音源方向ｄｍｉは、それぞれ大きさが１に正規化された単位ベクトルで示されている。βｍｉは、関数値の広がりを示す形状パラメータを示す。形状パラメータβｍｉが大きいほど、第１の確率は正規分布に近似し、形状パラメータβｍｉが小さいほど、第２の確率は一様分布に近似する。Ｉ０（βｍｉ）は、第０次の第一種変形ベッセル関数を示す。フォン・ミーゼス分布は、音源方向のように角度に加わったノイズの分布をモデル化するうえで好適である。確率モデルでは、形状パラメータβｍｉをモデルパラメータの１つとする。
【０１０６】
音響処理システムＳ１全体として、定位音源方向［ｄ’］のもとで推定音源方向［ｄ］が得られる確率ｐ（［ｄ’］｜［ｄ］）は、次式（２１）に示すようにマイクロホンアレイｍｉ間での第１確率ｆ（ｄ’ｍｉ，ｄｍｉ；βｍｉ）の総乗と仮定する。
【０１０７】
【数２１】
2020038123.tif 000023
【０１０８】
ここで、定位音源方向［ｄ’］、推定音源方向［ｄ］は、それぞれ定位音源方向ｄ’ｍｉ、推定音源方向ｄｍｉを要素として含むベクトルである。
また、確率モデルでは、マイクロホンアレイｍｊ、ｍｋ、推定音源方向ｄｍｊ、ｄｍｋ間の中間点ｓｊ，ｋが得られるとき、その中間点ｓｊ，ｋが分類されるクラスタｃｊ，ｋに対応する推定音源位置が得られる第２確率ｐ（ｓｊ，ｋ｜ｃｊ，ｋ）が、次式（２２）に示す多変量ガウス分布Ｎ（ｓｊ，ｋ；μｃｊ，ｋ，Σｃｊ，ｋ）に従うことを仮定する。μｃｊ，ｋ、Σｃｊ，ｋは、それぞれ多変量ガウス分布の平均、分散を示す。この平均は、推定音源位置の分布の大きさや偏りを示す。中間点ｓｊ，ｋは、上述したように、マイクロホンアレイｍｊ、ｍｋそれぞれの位置ｕｊ、ｕｋと、推定音源方向ｄｍｊ、ｄｍｋとから定まる関数である。以下の説明では、中間点の位置を、ｇ（ｄｍｊ、ｄｍｋ）と示すことがある。確率モデルでは、平均μｃｊ，ｋ、分散Σｃｊ，ｋをモデルパラメータの一部とする。
【０１０９】
【数２２】
2020038123.tif 000024
【０１１０】
音響処理システムＳ１全体として、各２つのマイクロホンアレイ、推定音源方向［ｄ］間の中間点の分布が得られるとき、それぞれの候補音源に対応するクラスタ［ｃ］が得られる確率ｐ（［ｄ］｜［ｃ］）は、次式（２２）に示すように中間点間での第２確率ｐ（ｓｊ，ｋ｜ｃｊ，ｋ）の総乗に近似されるものと仮定する。［ｃ］は、クラスタｃｊ，ｋを要素として含むベクトルである。
【０１１１】
【数２３】
2020038123.tif 000025
【０１１２】
また、確率モデルでは、第３確率として、２つのマイクロホンアレイｍｊ、ｍｋ、推定音源方向ｄｍｊ、ｄｍｋ間の中間点ｓｊ，ｋが分類されるクラスタｃｊ，ｋの出現確率ｐ（ｃｊ，ｋ）をモデルパラメータの１つとする。このパラメータをπｃｊ，ｋと表すことがある。
【０１１３】
（音源位置の更新）
次に、上述した確率モデルを用いた音源位置の更新処理について説明する。
音源位置更新部１４２は、音源定位により定位音源方向［ｄ’］が得られるとき、音源候補ごとの推定音源位置［ｄ］がそれぞれの音源候補に対応するクラスタ［ｃ］に分類される推定確率ｐ（［ｃ］，［ｄ］，［ｄ’］）が高くなるように、推定音源位置［ｄ］を再帰的に更新する。音源位置更新部１４２は、各２つのマイクロホンアレイ、推定音源方向間の中間点の分布についてクラスタリングを行ってクラスタ［ｃ］に分類する。
推定音源位置［ｄ］を更新するため、音源位置更新部１４２は、ビタビ学習法（Ｖｉｔｅｒｂｉ  Ｔｒａｉｎｉｎｇ）を応用した手法を用いる。
【０１１４】
音源位置更新部１４２は、次式（２３）に示すように、モデルパラメータ［μ＊］，［Σ＊］，［β＊］を一定として、推定確率ｐ（［ｃ］，［ｄ］，［ｄ’］；［μ＊］，［Σ］＊，［β＊］）を最大化する推定音源位置［ｄ＊］、クラスタ［ｃ＊］を算出する処理と、次式（２４）に示すように、算出した推定音源位置［ｄ＊］、クラスタ［ｃ＊］を一定として、推定確率ｐ（［ｃ＊］，［ｄ＊］，［ｄ’］；［μ］，［Σ］，［β］）を最大化するモデルパラメータ［π＊］、［μ＊］、［Σ＊］、［β＊］を算出する処理と、を逐次に繰り返す。…＊は、最大化したパラメータ…を示す。ここで、最大化とは、巨視的に増加させること、もしくはそのための処理を意味し、その処理により一時的もしくは局所的に減少する場合もありうる。
【０１１５】
【数２４】
2020038123.tif 000026
【０１１６】
【数２５】
2020038123.tif 000027
【０１１７】
式（２３）の右辺は、式（２０）〜（２２）を代入して、次式（２５）に示すように変形される。
【０１１８】
【数２６】
2020038123.tif 000028
【０１１９】
式（２５）に示すように、推定確率ｐ（［ｃ］，［ｄ］，［ｄ’］）は、上述の第１確率と、第２確率と、第３確率と、をそれぞれ因子とする積で表される。但し、式（２５）において値がゼロ以下となる因子を、乗算対象としない。
式（２５）の右辺は、次式（２６）、（２７）に示すようにクラスタｃｊ，ｋの関数と音源方向［ｄ］の関数に分解される。従って、クラスタｃｊ，ｋと推定音源方向［ｄ］は、個々に更新可能となる。
【０１２０】
【数２７】
2020038123.tif 000029
【０１２１】
【数２８】
2020038123.tif 000030
【０１２２】
音源位置更新部１４２は、式（２６）の右辺の値をより大きくするように全ての中間点（含む交点）ｇ（ｄ＊ｍｊ，ｄ＊ｍｋ）をクラスタｃ＊ｊ，ｋを要素とするクラスタ［ｃ＊］に分類する。音源位置更新部１４２は、クラスタｃ＊ｊ，ｋを定める際、階層クラスタリングを行う。
階層クラスタリングは、各２つのクラスタ間の距離を算出し、最も距離が小さい２つのクラスタを併合して新たなクラスタを生成する処理を逐次に繰り返す手法である。このとき、音源位置更新部１４２は、２つのクラスタ間の距離として、一方のクラスタに分類される中間点ｇ（ｄ＊ｍｊ，ｄ＊ｍｋ）と他方のクラスタｃｊ’，ｋ’の中心である平均μｃｊ’，ｋ’との間の距離のうち最も小さい距離を用いる。
【０１２３】
一般に、推定音源方向［ｄ］は、他の変数との依存性が高いため解析的に最適値を算出することは困難である。そこで、式（２７）の右辺を次式（２８）に示すように近似的に推定音源方向ｄｍｉの関数に分解する。音源位置更新部１４２は、式（２４）の右辺第３〜５行に示す値をコスト関数としてより大きくするように個々の推定音源方向ｄｍｉを更新する。
【０１２４】
【数２９】
2020038123.tif 000031
【０１２５】
推定音源方向ｄｍｉを更新する際、音源位置更新部１４２は、次に説明する制約条件（ｃ１）、（ｃ２）のもとで、最急降下法（ｇｒａｄｉｅｎｔ  ｄｅｓｃｅｎｔ  ｍｅｔｈｏｄ）を用いて推定音源方向ｄ＊ｍｉを探索する。
（ｃ１）音源定位により推定された定位音源方向［ｄ’］のそれぞれが、それぞれ対応する真の音源方向［ｄ］に近似している。
（ｃ２）推定音源位置に相当する平均μｃｊ，ｋが、直前に更新された推定音源方向ｄ＊ｍｊ、ｄ＊ｍｋ、ｄ＊ｍｉに基づく３つの中間点Ｐｊ、Ｐｋ、Ｐｉを頂点とする三角形の領域内にある。但し、マイクロホンアレイｍｉは、マイクロホンアレイｍｊ、ｍｋとは別個のマイクロホンアレイである。
【０１２６】
例えば、推定音源方向ｄｍ３を更新する際、図１１に示すように、音源位置更新部１４２は、マイクロホンアレイｍ３から中間点Ｐ２の方向を起点ｄｍｉｎ（ｍ３）とし、マイクロホンアレイｍ３から中間点Ｐ１の方向を終点ｄｍａｘ（ｍ３）とする方向の範囲内で、上述のコスト関数が最も大きくなる推定音源方向ｄｍ３を推定音源方向ｄ＊ｍ３として定める。図１１は、本実施形態に係る音源方向探索の説明図である。他の音源方向ｄｍ１、ｄｍ２等を更新する際も、音源位置更新部１４２は、同様の制約条件を課してコスト関数が最も大きくなる推定音源方向ｄｍ１、ｄｍ２を探索する。即ち、音源位置更新部１４２は、マイクロホンアレイｍ１から中間点Ｐ３の方向を起点ｄｍｉｎ（ｍ１）とし、中間点Ｐ２の方向を終点ｄｍａｘ（ｍ１）とする方向の範囲内で、コスト関数が最も大きくなる推定音源方向ｄ＊ｍ１を探索する。音源位置更新部１４２は、マイクロホンアレイｍ２から中間点Ｐ１の方向を起点ｄｍｉｎ（ｍ２）とし、中間点Ｐ３の方向を終点ｄｍａｘ（ｍ２）とする方向の範囲内で、コスト関数が最も大きくなる推定音源方向ｄ＊ｍ２を探索する。従って、推定音源方向の探索領域が、直前に更新された推定音源方向ｄ＊ｍ１等に基づいて定めた探索領域内に制限されるので、計算量が低減することができる。また、コスト関数の非線形性による解の不安定性が回避される。
【０１２７】
なお、式（２４）の右辺は、式（２０）〜（２２）を代入して、次式（２９）に示すように変形される。音源位置更新部１４２は、式（２９）の右辺の値を大きくするように、モデルパラメータのセット［π＊］、［μ＊］、［Σ＊］、［β＊］を更新する。
【０１２８】
【数３０】
2020038123.tif 000032
【０１２９】
音源位置更新部１４２は、式（２９）の右辺の値をより大きくするため、次式（３０）に示す関係を用いて、定位音源方向［ｄ’］、更新された推定音源方向［ｄ＊］及び更新されたクラスタ［ｃ＊］に基づいて、各クラスタｃのモデルパラメータπ＊ｃ、μ＊ｃ、Σ＊ｃと各マイクロホンアレイｍのモデルパラメータβ＊ｍを算出することができる。
【０１３０】
【数３１】
2020038123.tif 000033
【０１３１】
式（３０）において、モデルパラメータπ＊ｃは、音源候補数Ｎに対する、推定音源位置がクラスタｃに属する音源候補数Ｎｃの割合、即ち、推定音源が分類されるクラスタｃへの出現確率を示す。モデルパラメータμ＊ｃは、クラスタｃに属する中間点ｓｊ，ｋ（＝ｇ（ｄ＊ｍｊ，ｄ＊ｍｋ））の座標の平均値、即ち、クラスタｃの中心を示す。モデルパラメータμ＊ｃは、クラスタｃに属する中間点ｓｊ，ｋの座標の分散を示す。モデルパラメータβ＊ｍは、マイクロホンアレイｉについての定位音源方向ｄ’ｍｉと推定音源方向ｄ＊ｍｉとの内積の平均値を示す。
【０１３２】
次に、本実施形態に係る音源位置更新処理の一例について説明する。
図１２は、本実施形態に係る音源位置更新処理の一例を示すフローチャートである。
【０１３３】
（ステップＳ１８２）音源位置更新部１４２は、更新処理に係る各種の初期値を設定する。音源位置更新部１４２は、初期値設定部１４０から入力された初期推定音源位置情報が示す音源候補ごとの推定音源位置の初期値を設定する。また、音源位置更新部１４２は、推定音源位置の初期値［ｄ］、クラスタの初期値［ｃ］、出現確率の初期値π＊ｃ、平均の初期値μ＊ｃ、分散の初期値Σ＊ｃ、形状パラメータの初期値β＊ｍ、を、それぞれ次式（３１）に示すように設定する。推定音源方向の初期値［ｄ］として、定位音源方向［ｄ’］が設定される。クラスタの初期値ｃｊ，ｋとして、音源推定位置の初期値ｘｎが属するクラスタｃｎが設定される。出現確率の初期値π＊ｃとして、クラスタ数Ｃの逆数が設定される。平均の初期値μ＊ｃとして、クラスタｃに属する音源推定位置の初期値ｘｎの平均値が設定される。分散の初期値Σ＊ｃとして、単位行列が設定される。形状パラメータの初期値β＊ｍとして、１が設定される。その後、ステップＳ１８４の処理に進む。
【０１３４】
【数３２】
2020038123.tif 000034
【０１３５】
（ステップＳ１８４）音源位置更新部１４２は、上述の制約条件のもとで、式（２４）の右辺に示すコスト関数が大きくなるように推定音源方向ｄ＊ｍｉを更新する。その後、ステップＳ１８６の処理に進む。
【０１３６】
（ステップＳ１８６）音源位置更新部１４２は、式（２６）に示す関係を用いて各クラスタｃの出現確率π＊ｃ、平均μ＊ｃ、分散Σ＊ｃと各マイクロホンアレイｍの形状パラメータβ＊ｍを算出する。その後、ステップＳ１８８の処理に進む。
【０１３７】
（ステップＳ１８８）音源位置更新部１４２は、更新した推定音源方向ｄ＊ｍｊ、ｄ＊ｍｋから中間点ｇ（ｄ＊ｍｊ，ｄ＊ｍｋ）を定める。音源位置更新部１４２は、式（２２）の右辺に示すコスト関数の値が大きくなるように、中間点（ｄ＊ｍｊ，ｄ＊ｍｋ）の分布についてクラスタリングを行って複数のクラスタｃｊ，ｋに分類する。その後、ステップＳ１９０の処理に進む。
【０１３８】
（ステップＳ１９０）音源位置更新部１４２は、音源方向ｄ＊ｍｉと推定音源位置ｘ＊ｎとする平均μｃｊ，ｋのいずれか又は両方の更新量を算出し、算出した更新量が所定の更新量よりも小さいか否かにより、収束したか否かを判定する。更新量は、例えば、更新前後の音源方向ｄ＊ｍｉの差分のマイクロホンアレイ間ｍｉ間の二乗和、平均μｃｊ，ｋの更新前後の差分のクラスタｃ間の二乗和の一方又はそれらの重み付き和のいずれであってもよい。収束したと判定する場合（ステップＳ１９０  ＹＥＳ）、ステップＳ１９２の処理に進む。収束していないと判定する場合（ステップＳ１９０  ＮＯ）、ステップＳ１８４の処理に戻る。
【０１３９】
（ステップＳ１９２）音源位置更新部１４２は、更新された推定音源位置ｘ＊ｎ最確（ｍｏｓｔ  ｐｒｏｂａｂｌｅ）音源位置として定める。音源位置更新部１４２は、音源候補ごとの推定音源位置を示す推定音源位置情報を音源特定部１６に出力する。音源位置更新部１４２は、更新された推定音源方向［ｄ＊］を最確音源方向として定め、音源候補ごとの推定音源方向を示す推定音源位置情報を音源特定部１６に出力してもよい。また、音源位置更新部１４２は、音源候補ごとの音源識別情報をさらに推定音源位置情報に含めて出力してもよい。音源識別情報には、各音源候補の推定音源位置の初期値に係る３個のマイクロホンアレイを示すインデックスの少なくともいずれか１つと、マイクロホンアレイごとの音源定位により推定された音源を示すインデックスの少なくともいずれか１つが含まれればよい。その後、図１２に示す処理を終了する。
【０１４０】
（音源特定部の処理）
次に、本実施形態に係る音源特定部１６の処理について説明する。音源位置更新部１４２は、３個のマイクロホンアレイのうち、各２つのマイクロホンアレイにより取得された音源方向の３つの中間点に基づいて、推定音源位置を定めていた。しかしながら、各マイクロホンアレイから取得されて音響信号により独立に音源方向が推定されうる。そのため、音源位置更新部１４２は、２個のマイクロホンアレイのそれぞれについて、互いに異なる音源の音源方向同士で中間点を定めてしまうことがある。その中間点は、音源が実在している位置とは異なる位置に生じるため、いわゆるゴースト（虚像）として検出されることがある。例えば、図１３に示す例では、マイクロホンアレイＭＡ１、ＭＡ２、ＭＡ３により、それぞれ音源Ｓ１、Ｓ２、Ｓ１の方向に音源方向が推定される。図１３は、虚像の検出例を示す図である。この場合、マイクロホンアレイＭＡ１、ＭＡ３による中間点Ｐ３は、いずれも音源Ｓ１の方向に基づいて定められるため、音源Ｓ１の位置に近似する。しかしながら、マイクロホンアレイＭＡ２、ＭＡ３による中間点Ｐ２は、それぞれ音源Ｓ２、Ｓ１の方向に基づいて定められるため、音源Ｓ１、Ｓ２のいずれの位置からも離れた位置となる。
【０１４１】
そこで、音源特定部１６は、マイクロホンアレイごとの各音源の音源別信号のスペクトルを複数の第２のクラスタに分類し、前記第２のクラスタのそれぞれに属する各スペクトルに係る音源が同一であるか否かを判定する。音源特定部１６は、同一と判定した音源の前記推定音源位置を、同一でないと判定した音源よりも優先して選択する。これにより、虚像の検出により音源位置が誤って推定されることが防止される。
【０１４２】
（周波数分析）
周波数分析部１２４は、音源ごとに分離された音源別音響信号について周波数分析を行う。図１４は、本実施形態に係る周波数分析処理の一例を示すフローチャートである。
【０１４３】
（ステップＳ２０２）周波数分析部１２４は、各マイクロホンアレイｍで取得された音響信号から分離された各音源の音源別音響信号をフレームごとに短時間フーリエ変換を行ってスペクトル［Ｆｍ，１］、［Ｆｍ，２］〜［Ｆｍ，ｓｍ］を算出する。その後、ステップＳ２０４の処理に進む。
【０１４４】
（ステップＳ２０４）周波数分析部１２４は、音源ごとに算出した周波数スペクトルをマイクロホンアレイｍごとに行間で統合して、スペクトル行列［Ｆｍ］を構成する。周波数分析部１２４は、マイクロホンアレイｍごとのスペクトル行列［Ｆｍ］を行間で統合してスペクトル行列［Ｆ］を構成する。周波数分析部１２４は、構成したスペクトル行列［Ｆ］と音源方向情報とを対応付けて音源特定部１６に出力する。その後、図１０に示す処理を終了する。
【０１４５】
（スコア算出）
音源特定部１６の分散算出部１６０とスコア算出部１６２は、次に例示するスコア算出処理を行う。
図１５は、本実施形態に係るスコア算出処理の一例を示すフローチャートである。
【０１４６】
（ステップＳ２２２）分散算出部１６０は、周波数分析部１２４から入力されるスペクトル行列［Ｆ］が示すマイクロホンアレイｍならびに音源の組ごとのスペクトルについてｋ−平均法を用いてクラスタリングを行い、複数の第２クラスタに分類する。クラスタ数Ｋは、予め分散算出部１６０に設定しておく。但し、分散算出部１６０は、スペクトルごとのクラスタの初期値を、繰り返し回数ｒごとに変更する。クラスタ数Ｋは、音源候補数Ｎと等しくしてもよい。分散算出部１６０は、スペクトルごとに分類される第２クラスタのインデックスｃｉ，ｘ＊ｎを要素として含むクラスタ行列［ｃ＊］を構成する。クラスタ行列［ｃ＊］の各列、各行は、それぞれマイクロホンアレイｉ、音源ｘ＊ｎに対応付けられる。マイクロホンアレイの数Ｍが３である場合、クラスタ行列［ｃ＊］は、次式（３２）に示すように、Ｎ行３列の行列となる。
【０１４７】
【数３３】
2020038123.tif 000035
【０１４８】
分散算出部１６０は、音源位置更新部１４２から入力される推定音源位置情報が示す音源候補ごとの音源識別情報に基づいて、各音源候補に対応する第２クラスタを特定する。
分散算出部１６０は、例えば、クラスタ行列において音源識別情報が示すマイクロホンアレイの列と音源列のうち、クラスタ行列に含まれるマイクロホンアレイの列と音源の行に配置された、インデックスが示す第２クラスタを特定することができる。
分散算出部１６０は、第２クラスタに対応する音源候補ごとの推定音源位置の分散Ｖｘ＊ｎを算出する。その後、ステップＳ２２４の処理に進む。
【０１４９】
（ステップＳ２２４）分散算出部１６０は、第２クラスタｃｘ＊ｎのそれぞれについて、分類された複数のスペクトルに係る音源が互いに同一の音源であるか否かを判定する。分散算出部１６０は、例えば、複数のスペクトルのうち、各２つのスペクトル間の類似度を示す指標が示す類似度が、いずれも所定の類似度よりも高いとき、同一の音源であると判定する。分散算出部１６０は、少なくとも１組のスペクトル間の類似度を示す指標が所定の類似度以下となるとき、同一の音源ではないと判定する。類似度の指標として、例えば、内積、ユークリッド距離、などを用いることができる。内積は、その値が大きいほど類似度が高いことを示す。ユークリッド距離は、その値が小さいほど類似度が低いことを示す。なお、分散算出部１６０は、複数のスペクトルの類似度の指標として、それらの分散を算出してもよい。分散算出部１６０は、分散が所定の分散の閾値よりも小さいとき、同一の音源であると判定し、分散がその閾値以上であるとき、同一の音源ではないと判定してもよい。同一の音源であると判定する場合（ステップＳ２２４；ＹＥＳ）、ステップＳ２２６の処理に進む。同一の音源ではないと判定する場合（ステップＳ２２４；ＮＯ）、ステップＳ２２８の処理に進む。
【０１５０】
（ステップＳ２２６）分散算出部１６０は、現在の繰り返し回数ｒにおいて第２クラスタｃｘ＊ｎについて算出した分散Ｖｘ＊ｎ（ｒ）が、前回の繰り返し回数ｒ−１に算出した分散Ｖｘ＊ｎ（ｒ−１）以下になったか否かを判定する。分散Ｖｘ＊ｎ（ｒ−１）以下になったと判定する場合（ステップＳ２２６；ＹＥＳ）、ステップＳ２３２の処理に進む。分散Ｖｘ＊ｎ（ｒ−１）より大きいと判定する場合（ステップＳ２２６；ＮＯ）、ステップＳ２３０の処理に進む。
【０１５１】
（ステップＳ２２８）分散算出部１６０は、現在の繰り返し回数ｒの第２クラスタｃｘ＊ｎの分散Ｖｘ＊ｎ（ｒ）をＮａＮと設定し、スコアｅｎ,ｒをδとする。ＮａＮは、分散が無効であることを示す記号（ｎｏｔ  ａ  ｎｕｍｂｅｒ）である。δは、０よりも小さい所定の実数である。その後、ステップＳ２３４の処理に進む。
【０１５２】
（ステップＳ２３０）分散算出部１６０は、現在の繰り返し回数ｒの第２クラスタｃｘ＊ｎのスコアｅｎ,ｒを０とする。その後、ステップＳ２３４の処理に進む。
【０１５３】
（ステップＳ２３２）分散算出部１６０は、現在の繰り返し回数ｒの第２クラスタｃｘ＊ｎのスコアｅｎ,ｒをεとする。その後、ステップＳ２３４の処理に進む。
【０１５４】
（ステップＳ２３４）分散算出部１６０は、現在の繰り返し回数ｒが所定の繰り返し回数Ｒに達したか否かを判定する。達していないと判定するとき（ステップＳ２３４；ＮＯ）、ステップＳ２３６の処理に進む。達したと判定するとき（ステップＳ２３４；ＹＥＳ）、分散算出部１６０は、第２クラスタごとの各回のスコアと推定音源位置を示すスコア算出情報をスコア算出部１６２に出力し、ステップＳ２３８の処理に進む。
【０１５５】
（ステップＳ２３６）分散算出部１６０は、現在の繰り返し回数ｒを、１増加させる。その後、ステップＳ２２２の処理に戻る。
【０１５６】
（ステップＳ２３８）スコア算出部１６２は、次式（３３）に示すように分散算出部１６０から入力されるスコア算出情報に基づいて、第２クラスタｃｘ＊ｎごとにスコアｅｎ,ｒの合計値ｅｎを算出する。スコア算出部１６２は、座標値ｘｎが相互に所定の範囲内にある推定音源位置ｘｉにそれぞれ対応する第２クラスタｉの合計値ｅｉの総和ｅ’ｎを算出する。これは、相互に座標値が等しいもしくは所定の範囲内にある推定音源位置に対応する第２クラスタを、１個の第２クラスタとして統合するためである。相互に座標値が等しいもしくは所定の範囲内にある推定音源位置に対応する第２クラスタが生じるのは、一般に各１個の音源からの発音期間の方が周波数分析に係るフレーム長よりも長いうえ、周波数特性が変動するためである。
【０１５７】
【数３４】
2020038123.tif 000036
【０１５８】
スコア算出部１６２は、次式（３４）に示すように分散算出部１６０から入力されるスコア算出情報に基づいて、第２クラスタｃｘ＊ｎごとに有効な分散が算出された回数を存在度数ａｎとして計数する。スコア算出部１６２は、有効な分散が算出されていないか否かを、分散Ｖｘ＊ｎ（ｒ）にＮａＮが設定されたか否かにより判定することができる。式（３０）の第１行の右辺のａｎ，ｒは、ＮａＮが設定された繰り返し回数ｒについて０、ＮａＮが設定されていない繰り返し回数ｒについて１となる。
スコア算出部１６２は、座標値ｘｎが相互に所定の範囲内にある推定音源位置ｘｉにそれぞれ対応する第２クラスタｉの存在度数ａｉの総和ａ’ｎを算出する。その後、ステップＳ２４０の処理に進む。
【０１５９】
【数３５】
2020038123.tif 000037
【０１６０】
（ステップＳ２４０）スコア算出部１６２は、次式（３５）に示すように、統合した第２クラスタｎのそれぞれについてスコアの総和ｅ’ｎを存在度数の総和ａ’ｎで除算して最終スコアｅ＊ｎを算出する。統合した第２クラスタｎは、個々の音源候補に対応する。スコア算出部１６２は、算出した音源候補ごとの最終スコアと推定音源位置を示す最終スコア情報を音源選択部１６４に出力する。その後、図１５に示す処理を終了する。
【０１６１】
【数３６】
2020038123.tif 000038
【０１６２】
上述の例では、ステップＳ２２８、Ｓ２３０、Ｓ２３２においてスコアｅｎ,ｒをそれぞれδ、０、εとする場合を例にしたが、これには限られない。ステップＳ２２８、Ｓ２３０、Ｓ２３２において定められるスコアｅｎ,ｒの値の大小関係は、その昇順であればよい。
【０１６３】
（音源選択）
音源選択部１６４は、次に例示する音源選択処理を行う。
図１６は、本実施形態に係る音源選択処理の一例を示すフローチャートである。
【０１６４】
（ステップＳ２４２）音源選択部１６４は、スコア算出部１６２から入力された最終スコア情報が示す音源候補の最終スコアｅ＊ｎが所定の最終スコアの閾値θ２以上であるか否かを判定する。閾値θ２以上と判定する場合（ステップＳ２４２；ＹＥＳ）、ステップＳ２４４の処理に進む。閾値θ２未満と判定する場合（ステップＳ２４２；ＮＯ）、ステップＳ２４６の処理に進む。
【０１６５】
（ステップＳ２４４）音源選択部１６４は、最終スコアｅ＊ｎが正常値（Ｉｎｌｉｅｒ）と判定し、その音源候補を音源として選択する。音源選択部１６４は、選択した音源に対応する推定音源位置を示す出力音源位置情報を音響処理装置１の外部に出力部１８を介して出力する。
【０１６６】
（ステップＳ２４６）音源選択部１６４は、最終スコアｅ＊ｎが異常値（Ｏｕｔｌｉｅｒ）と判定し、対応する音源候補を音源として選択せずに棄却する。その後、図１２に示す処理を終了する。
【０１６７】
（音響処理）
音響処理装置１は、全体として次に例示する音響処理を行う。
図１７は、本実施形態に係る音響処理の一例を示すフローチャートである。
【０１６８】
（ステップＳ１２）音源定位部１２０は、入力部１０から入力され、各マイクロホンアレイから取得された複数チャネルの音響信号に基づいて各音源の定位音源方向を予め定めた長さのフレームごとに推定する（音源定位）。音源定位部１２０は、音源定位において、例えば、ＭＵＳＩＣ法を用いる。その後、ステップＳ１４の処理に進む。
【０１６９】
（ステップＳ１４）音源分離部１２２は、音源ごとの定位音源方向に基づいて、各マイクロホンアレイから取得された音響信号を音源ごとの音源別音響信号に分離する。音源分離部１２２は、音源分離部において、例えば、ＧＨＤＳＳ法を用いる。その後、ステップＳ１６の処理に進む。
【０１７０】
（ステップＳ１６）初期値設定部１４０は、三角分割法により、３個のマイクロホンアレイのうち、各２個のマイクロホンアレイの組ごとに推定された定位音源方向に基づいて中間点（含む交点）を定める。初期値設定部１４０は、定めた中間点を音源候補の推定音源位置の初期値として定める。その後、ステップＳ１８の処理に進む。
【０１７１】
（ステップＳ１８）音源位置更新部１４２は、各２個のマイクロホンアレイの組ごとに推定音源方向に基づいて定められる中間点の分布を複数のクラスタに分類する。音源位置更新部１４２は、音源候補ごとの推定音源位置が、それぞれの音源候補に対応するクラスタに属する確率が高くなるように推定音源位置を更新する。ここで、音源位置更新部１４２は、上述の音源位置更新処理を行う。その後、ステップＳ２０の処理に進む。
【０１７２】
（ステップＳ２０）周波数分析部１２４は、各マイクロホンアレイについて音源ごとに分離した音源別音響信号について周波数分析を行い、スペクトルを算出する。その後、ステップＳ２２の処理に進む。
【０１７３】
（ステップＳ２２）分散算出部１６０は、算出したスペクトルを複数の第２クラスタに分類し、分類した第２クラスタに属するスペクトルに係る音源が相互に同一であるか否かを判定する。分散算出部１６０は、第２クラスタに属するスペクトルに係る音源候補ごとの推定音源位置の分散を算出する。スコア算出部１６２は、同一と判定された音源に係る第２クラスタを、同一でないと判定した音源に係る第２クラスタよりも大きくなるように第２クラスタごとの最終スコアを定める。スコア算出部１６２は、クラスタの安定性として、繰り返しごとの推定音源位置の分散の増加が稀な第２クラスタほど大きくなるように最終スコアを定める。ここで、分散算出部１６０とスコア算出部１６２は、上述のスコア算出処理を行う。その後、ステップＳ２４の処理に進む。
【０１７４】
（ステップＳ２４）音源選択部１６４は、最終スコアが、所定の最終スコアの閾値以上となる第２クラスタに対応する音源候補を音源として選択し、最終スコアの閾値未満となる第２クラスタに対応する音源候補を棄却する。音源選択部１６４は、選択した音源に係る推定音源位置を出力する。その後、図１７に示す処理を終了する。
【０１７５】
（フレームデータ解析）
音響処理システムＳ１は、記憶部１９を備え、図１７に示す音響処理を行う前に、各マイクロホンアレイが収音した音響信号を記憶しておいてもよい。記憶部１９は、音響処理装置１の一部として構成されてもよいし、音響処理装置１とは別個の外部機器に設置されてもよい。音響処理装置１は、記憶部から読み出した音響信号を用いて図１７に示す音響処理を行ってもよい（バッチ処理）。
【０１７６】
上述の図１７の音響処理のうち、音源位置更新処理（ステップＳ１８）、スコア算出処理（ステップＳ２２）は、複数のフレームの音響信号に基づく各種のデータを要するうえ、処理時間が長い。オンライン処理において、あるフレームについて図１７の処理を完了した後で、次のフレームの処理を開始すると、出力が間欠的となるため現実的ではない。
そこで、オンライン処理において、初期処理部１２によるステップＳ１２、Ｓ１４、Ｓ２０の処理が、音源位置推定部１４と音源特定部１６によるステップＳ１６、Ｓ１８、Ｓ２２、Ｓ２４の処理と並列に行われてもよい。但し、ステップＳ１２〜Ｓ１４、Ｓ２０の処理において、現時点ｔ０までの第１区間内の音響信号もしくは音響信号から導出された各種のデータを処理対象とする。ステップＳ１２、Ｓ１４、Ｓ２０の処理において、現時点ｔ０までの第１区間内の音響信号もしくは音響信号から導出された各種のデータを処理対象とする。ステップＳ１６、Ｓ１８、Ｓ２２、Ｓ２４の処理において、第１区間よりも過去の第２区間内の音響信号もしくは各種のデータを処理対象とする。
【０１７７】
図１８は、処理対象のデータ区間の例を示す図である。
図１８において、左右方向は時刻を示す。右上のｔ０は、現時点を示す。ｗｌは、個々のフレームｗ１、ｗ２、…のフレーム長を示す。音響処理装置１の入力部１０には、フレームごとに最新の音響信号が入力され、音響処理装置１の記憶部１９は、期間がｎｅ・ｗｌの音響信号と導出されるデータを記憶する。そして、記憶部は、フレームごとに最も過去の音響信号とデータを棄却する。ｎｅは、記憶される全データのフレーム数を示す。初期処理部１２は、全データのうち最新の第１区間内のデータを用いて、ステップＳ１２〜Ｓ１４、Ｓ２０の処理を行う。第１区間の長さが、初期処理長ｎｔ・ｗｌに相当する。ｎｔは、予め定めた初期処理長のフレーム数を示す。音源位置推定部１４と音源特定部１６は、全データのうち第１区間の終期よりも後の第２区間のデータを用いて、ステップＳ１６、Ｓ１８、Ｓ２２、Ｓ２４の処理を行う。第２区間の長さが、バッチ長ｎｂ・ｗｌに相当する。ｎｂは、予め定めたバッチ長のフレーム数を示す。第１区間、第２区間には、フレームごとに、それぞれ最新のフレームの音響信号、第ｎｔ＋１フレームの音響信号と導出されるデータが加入される。他方、第１区間、第２区間には、フレームごとに第ｎｔフレームの音響信号とその音響信号から導出されるデータと、第ｎｅフレームの音響信号と導出されるデータが棄却される。このように、初期処理部１２と、音源位置推定部１４ならびに音源特定部１６は、それぞれ第１区間内のデータと、第２区間内のデータとを使い分けることで、出力がフレーム間で継続するように図１３に示す音響処理がオンラインで実行可能となる。
【０１７８】
以上に説明したように、本実施形態に係る音響処理装置１は、位置が異なるＭ（Ｍは、３以上の整数）個の収音部のそれぞれから取得した複数チャネルの音響信号に基づいて音源の方向である定位音源方向を定める音源定位部と、２個の前記収音部の組ごとに当該収音部のそれぞれから前記音源の推定音源位置への方向である推定音源方向への半直線双方に対して垂直となる線分の中間点を前記推定音源位置として推定する音源位置推定部と、を備える。
この構成により、それぞれ異なる収音部からの定位音源方向による半直線が交差しなくても半直線双方に対して垂直となる線分の中間点（含む交点）によって定まる点を推定音源位置として推定することができる。
【０１７９】
また、本実施形態に係る音響処理装置１において、前記音源位置推定部は、前記中間点の分布を複数のクラスタに分類し、前記推定音源位置が前記音源に対応する。クラスタに分類される確率である推定確率が高くなるように前記推定音源位置を更新する
この構成により、それぞれ異なる収音部２０からの定位音源方向により定まる中間点が分類されるクラスタの範囲内に、対応する音源の推定音源位置が分類される可能性が高くなるように推定音源位置が調整される。クラスタの範囲内には音源が存在する可能性が高くなるため、調整される推定音源位置がより正確な音源位置として得られる。
【０１８０】
また、推定確率は、定位音源方向が定められるとき推定音源方向が得られる確率である第１確率と、中間点が定められるとき推定音源位置が得られる確率である第２確率と、中間点が分類されるクラスタの出現確率である第３確率と、をそれぞれ因子とする積である。
一般に、定位音源方向、推定音源位置及び中間点は相互に依存するが、音源位置推定部１４は、第１確率、第２確率及び第３確率をそれぞれ独立な推定確率の因子として推定音源位置を定めることができる。そのため、推定音源位置の調整に係る計算負荷が低減する。
【０１８１】
また、第１確率は、定位音源方向を基準とするフォン・ミーゼス分布に従い、第２確率は、中間点の位置を基準とする多次元ガウス関数に従う。音源位置推定部１４は、推定確率が高くなるように、フォン・ミーゼス分布の形状パラメータと、多次元ガウス関数の平均ならびに分散と、を更新する。
この構成により、第１確率の推定音源方向の関数、第２確率の推定音源位置の関数が、それぞれ形状パラメータ、平均ならびに分散といった少数のパラメータで表される。そのため、推定音源位置の調整に係る計算負荷がさらに低減する。
【０１８２】
また、音源位置推定部１４は、収音部２０の３個から定められる３個の中間点（含む交点）の重心を推定音源位置の初期値として定める。
この構成により、推定音源位置の初期値を、音源が存在する可能性が高い３個の中間点をそれぞれ頂点とする三角形の領域内に設定することができる。そのため、調整による推定音源位置の変化が収束するまでの計算負荷が低減する。
【０１８３】
また、音響処理装置１は、複数チャネルの音響信号から音源ごとの音源別信号に分離する音源分離部１２２と、音源別信号のスペクトルを算出する周波数分析部１２４を備える。音響処理装置１は、算出したスペクトルを複数の第２クラスタに分類し、第２クラスタのそれぞれに分類される各スペクトルに係る音源が同一であるか否かを判定し、同一と判定した音源の推定音源位置を、同一でないと判定した音源よりも優先して選択する音源特定部１６を備える。
この構成により、スペクトルに基づいて同一と判定されなかった音源の定位音源方向の中間点に基づいて推定された推定音源位置が棄却される可能性が高くなる。そのため、互いに異なる音源の推定音源方向の中間点に基づいて推定音源位置が虚像（ゴースト）として誤って選択される可能性を低くすることができる。
【０１８４】
音源特定部１６は、第２クラスタのそれぞれに分類されるスペクトルに係る音源の推定音源位置の分散に基づいて当該第２クラスタの安定性を評価し、安定性が高い第２クラスタほど当該第２クラスタにスペクトルが分類される音源の推定音源位置を優先して選択する。
この構成により、推定音源位置が定常的な音源のスペクトルが分類される第２クラスタに対応する音源の推定音源位置が選択される可能性が高くなる。即ち、推定音源位置が選択される第２クラスタには、偶発的に互いに異なる音源の推定音源方向の中間点に基づいて推定される推定音源位置が含まれる可能性が低くなる。そのため、互いに異なる音源の推定音源方向の中間点に基づいて推定音源位置が虚像として誤って選択される可能性をさらに低くすることができる。
【０１８５】
（測定結果）
次に、本実施形態の音響処理装置１を用いて測定した結果例を説明する。
まず、マイクアレイの配置例を説明する。
図１９は、マイクロホンアレイの配置例を示す図である。図１９において、符号ｍ１，ｍ２，ｍ３，ｍ４それぞれがマイクロホンアレイである。
また、符号ｇ１１が示す画像は、実空間におけるマイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４の配置を示す図である。また、符号ｇ１２は、ｘｙ平面におけるマイクロホンアレイ設置座標を示す図である。なお、図１９において、水平方向・北向きにｘ軸、垂直方向・西向きにｙ軸、鉛直方向上向きにｚ軸を取る。また、原点はマイクロホンアレイｍ４の位置とする。高さｚについて、ｘ，ｙにかかわらず常にＴ１の位置の地表面を基準（ｚ＝０）とする。また、北をｘｙ平面における方位の０度方向とする。また、各マイクロホンアレイの地上高さは１．３ｍである。また、使用したマイクロホンアレイそれぞれは、１６個のマイクロホンがデバイスの表面上の異なる垂直および水平位置に分散されて配置されている。この環境で数時間の鳥の歌を録音した。
【０１８６】
図２０は、収音した第１の音響信号例を示す図である。図２０において、横軸は時刻であり、縦軸は音響信号の大きさである。また、収音した音響信号は、野鳥の鳴き声である。録音時間は２５秒である。この収音した音響信号に対して、上述した三次元の位置推定と外れ値除去（ゴースト除去）処理を行った。
図２０に示すように、このデータは、４つの大きな鳴き声イベントから構成され、通過する車や鳥が遠くで鳴くなどのノイズなどが含まれている。
【０１８７】
図２１は、マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第１の音響信号それぞれに対して解析した結果である。符号ｇ２１はマイクロホンアレイｍ１が収音した音響信号に対して解析した結果である。符号ｇ２２はマイクロホンアレイｍ２が収音した音響信号に対して解析した結果である。符号ｇ２３はマイクロホンアレイｍ３が収音した音響信号に対して解析した結果である。符号ｇ２４はマイクロホンアレイｍ４が収音した音響信号に対して解析した結果である。符号ｇ２１〜ｇ２４において、横軸は時刻であり、縦軸は方位角（ｒａｄ）である。また、符号ｇ３１等の縦線は、各イベントのおおよその開始時間を示している。
【０１８８】
図２１に示したように、外れ値除去アルゴリズムが消去することができなかったいくつかの雑音イベントを伴うマイクロホンアレイｍ１、ｍ２およびｍ３によって事象を適切に推定することができた。
マイクロホンアレイｍ４では、いくつかの主要イベントが省略されたが、他の3つのマイクロホンアレイからの情報を使用してローカライズを実行することができる。
【０１８９】
図２２は、図２０の音響信号から推定した鳥の位置を示す例である。
図２２において、符号ｇ３１はＸＹ平面における推定位置であり、横軸がＸ（ｍ）軸、縦軸がＹ（ｍ）軸である。符号ｇ３２はＸＺ平面における推定位置であり、横軸がＸ（ｍ）軸、縦軸がＺ（ｍ）軸である。符号ｇ３３はＹＺ平面における推定位置であり、横軸がＹ（ｍ）軸、縦軸がＺ（ｍ）軸である。
また丸印は、正常値である。バツ印は、異常値である。三角印は、マイクロホンアレイの位置を表す。
【０１９０】
２番目のデータでは、最初のデータと同じ設定で２２秒ごとに異なる鳥の曲の録音を選択した。
図２３は、第２の音響信号（不図示）から推定した鳥の位置を示す例である。
図２３において、符号ｇ４１は、本実施形態の３次元位置推定によるＸＹ平面における推定位置であり、横軸がＸ（ｍ）軸、縦軸がＹ（ｍ）軸である。符号ｇ４２は、本実施形態の３次元位置推定によるＸＺ平面における推定位置であり、横軸がＸ（ｍ）軸、縦軸がＺ（ｍ）軸である。符号ｇ４３は、本実施形態の３次元位置推定によるＹＺ平面における推定位置であり、横軸がＹ（ｍ）軸、縦軸がＺ（ｍ）軸である。符号ｇ４４は、従来技術による２次元位置推定結果と本実施形態の３次元位置推定結果を比較したものであり、横軸がＸ（ｍ）軸、縦軸がＹ（ｍ）軸である。
また、符号ｇ４１〜ｇ４３において、丸印は正常値であり、バツ印は異常値であり、三角印はマイクロホンアレイの位置を表す。
符号ｇ４４において、丸印は従来技術による２次元位置推定結果であり、バツ印は本実施形態の３次元位置推定結果であり、三角印はマイクロホンアレイの位置を表す。
【０１９１】
図２３に示す例では、いくつかの鳥が鳴いているように見えますが、それぞれ同じＸＹ位置にあり、高さはそれぞれ異なっている。このため、従来技術による二次元位置推定では方位角方向のみを使用するため、結果は鳴いている１つの鳥の位置だけが推定される。一方、本実施形態の三次元位置推定を行った結果は、鳥の位置の高度の差も分析することができる。
【０１９２】
図２４は、従来技術による２次元位置推定結果と本実施形態の３次元位置推定結果を比較した結果を示す図である。
符号ｇ５１は、ＸＹ平面における従来技術による２次元位置推定結果であり、横軸がＸ（ｍ）、縦軸がＹ（ｍ）である。符号ｇ４１とｇ４２は、図２３と同じである。符号ｇ５１において、丸印は正常値であり、バツ印は異常値であり、三角印はマイクロホンアレイの位置を表す。
符号ｇ５１と符号ｇ４１を重ねた結果が図２３の符号ｇ４４の図である。符号ｇ５１と符号ｇ４１のように、上から見ると（ＸＹ平面）、鳥の推定位置は同じように見える。しかしながら、符号ｇ４２のように、違う高さで鳥が鳴いていることがわかる。
【０１９３】
図２５は、収音した第２の音響信号の１チャネルの例を示す図である。図２５において、横軸は時刻であり、縦軸は音響信号の大きさである。また、収音した音響信号は、野鳥の鳴き声である。録音時間は２２秒である。この収音した音響信号に対して、上述した三次元の位置推定と外れ値除去（ゴースト除去）処理を行った。
【０１９４】
図２６は、マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第２の音響信号それぞれに対して方位角方向を解析した結果である。符号ｇ６１はマイクロホンアレイｍ１が収音した音響信号に対して解析した結果である。符号ｇ６２はマイクロホンアレイｍ２が収音した音響信号に対して解析した結果である。符号ｇ６３はマイクロホンアレイｍ３が収音した音響信号に対して解析した結果である。符号ｇ６４はマイクロホンアレイｍ４が収音した音響信号に対して解析した結果である。符号ｇ６１〜ｇ６４において、横軸は時刻であり、縦軸は方位角（ｒａｄ）である。
【０１９５】
図２７は、マイクロホンアレイｍ１，ｍ２，ｍ３，ｍ４それぞれが収音した第２の音響信号それぞれに対して仰角方向を解析した結果である。符号ｇ７１はマイクロホンアレイｍ１が収音した音響信号に対して解析した結果である。符号ｇ７２はマイクロホンアレイｍ２が収音した音響信号に対して解析した結果である。符号ｇ７３はマイクロホンアレイｍ３が収音した音響信号に対して解析した結果である。符号ｇ７４はマイクロホンアレイｍ４が収音した音響信号に対して解析した結果である。符号ｇ７１〜ｇ７４において、横軸は時刻であり、縦軸は仰角（ｒａｄ）である。
【０１９６】
図２６と図１７に示す例では、鳥のイベントがほぼ一定であり、これらのイベントの正確な開始時間を決定することは困難であった。
また、ほとんどの場合、情報は記録されたデータの主な音響信号のイベントと同時に取得される。
【０１９７】
以上の測定結果に示したように、本実施形態によれば、三次元位置推定により、高さの異なる鳥を区別することができる。また、測定結果より、樹木の位置や高さから、推定された鳥の位置が妥当と思われる。
【０１９８】
以上、図面を参照してこの発明の実施形態について説明してきたが、具体的な構成は上述のものに限られることはなく、この発明の要旨を逸脱しない範囲内において様々な設計変更等をすることが可能である。
【０１９９】
例えば、分散算出部１６０は、図１５の処理のうちステップＳ２２２、Ｓ２２４の処理を行い、ステップＳ２２６〜Ｓ２４０の処理を行わなくてもよい。その場合には、スコア算出部１６２が省略されてもよい。その場合、音源選択部１６４は、第２クラスタに分類されるスペクトルに係る音源が互いに同一と判定された第２クラスタに対応する候補音源を音源として選択し、同一と判定されない第２クラスタに対応する候補音源を棄却してもよい。音源選択部１６４は、選択した音源に対応する推定音源位置を示す出力音源位置情報を音響処理装置１の外部に出力する。
また、音響処理装置１において、周波数分析部１２４と音源特定部１６が省略されてもよい。その場合、音源位置更新部１４２は、音源候補ごとの推定音源位置を示す推定音源位置情報を出力部１８に出力する。
【０２００】
音響処理装置１は、収音部２０−１〜２０−Ｍと一体化した単一の装置として構成されてもよい。
収音部２０の数Ｍは、３個に限られず４個以上であってもよい。また、収音部２０ごとに収音可能とする音響信号のチャネル数が異なってもよいし、それぞれの音響信号から推定可能な音源数が異なってもよい。
第１確率が従う確率分布は、フォン・ミーゼス分布に限られず、ロジスティック関数の導関数など、１次元空間内のある基準値に対する最大値を与える１次元の確率分布であればよい。
第２確率が従う確率分布は、多次元ガウス関数に限られず、多次元ロジスティック関数の一次導関数など、多次元空間内のある基準値に対する最大値を与える多次元の確率分布であればよい。
【０２０１】
なお、本発明における音響処理装置１の機能の全てまたは一部を実現するためのプログラムをコンピュータ読み取り可能な記録媒体に記録して、この記録媒体に記録されたプログラムをコンピュータシステムに読み込ませ、実行することにより音響処理装置１が行う処理の全てまたは一部を行ってもよい。なお、ここでいう「コンピュータシステム」とは、ＯＳや周辺機器等のハードウェアを含むものとする。また、「コンピュータシステム」は、ホームページ提供環境（あるいは表示環境）を備えたＷＷＷシステムも含むものとする。また、「コンピュータ読み取り可能な記録媒体」とは、フレキシブルディスク、光磁気ディスク、ＲＯＭ、ＣＤ−ＲＯＭ等の可搬媒体、コンピュータシステムに内蔵されるハードディスク等の記憶装置のことをいう。さらに「コンピュータ読み取り可能な記録媒体」とは、インターネット等のネットワークや電話回線等の通信回線を介してプログラムが送信された場合のサーバやクライアントとなるコンピュータシステム内部の揮発性メモリ（ＲＡＭ）のように、一定時間プログラムを保持しているものも含むものとする。
【０２０２】
また、上記プログラムは、このプログラムを記憶装置等に格納したコンピュータシステムから、伝送媒体を介して、あるいは、伝送媒体中の伝送波により他のコンピュータシステムに伝送されてもよい。ここで、プログラムを伝送する「伝送媒体」は、インターネット等のネットワーク（通信網）や電話回線等の通信回線（通信線）のように情報を伝送する機能を有する媒体のことをいう。また、上記プログラムは、前述した機能の一部を実現するためのものであってもよい。さらに、前述した機能をコンピュータシステムにすでに記録されているプログラムとの組み合わせで実現できるもの、いわゆる差分ファイル（差分プログラム）であってもよい。
【０２０３】
以上、本発明を実施するための形態について実施形態を用いて説明したが、本発明はこうした実施形態に何等限定されるものではなく、本発明の要旨を逸脱しない範囲内において種々の変形および置換を加えることができる。
【符号の説明】
【０２０４】
Ｓ１…音響処理システム、１…音響処理装置、２０，２０−１，２０−２，・・・２０−Ｍ…収音部、１０…入力部、１２…初期処理部、１４…音源位置推定部、１６…音源特定部、１８…出力部、１９…記憶部、１２０…音源定位部、１２２…音源分離部、１２４…周波数分析部、１３０…初期値設定部、１４２…音源位置更新部、１６０…分散算出部、１６２…スコア算出部、１６４…音源選択部
【図１】
2020038123.tif 000039
【図２】
2020038123.tif 000040
【図３】
2020038123.tif 000041
【図４】
2020038123.tif 000042
【図５】
2020038123.tif 000043
【図６】
2020038123.tif 000044
【図７】
2020038123.tif 000045
【図８】
2020038123.tif 000046
【図９】
2020038123.tif 000047
【図１０】
2020038123.tif 000048
【図１１】
2020038123.tif 000049
【図１２】
2020038123.tif 000050
【図１３】
2020038123.tif 000051
【図１４】
2020038123.tif 000052
【図１５】
2020038123.tif 000053
【図１６】
2020038123.tif 000054
【図１７】
2020038123.tif 000055
【図１８】
2020038123.tif 000056
【図１９】
2020038123.tif 000057
【図２０】
2020038123.tif 000058
【図２１】
2020038123.tif 000059
【図２２】
2020038123.tif 000060
【図２３】
2020038123.tif 000061
【図２４】
2020038123.tif 000062
【図２５】
2020038123.tif 000063
【図２６】
2020038123.tif 000064
【図２７】
2020038123.tif 000065