(57)【要約】      （修正有）【課題】従来技術よりも、音場推定が有効な空間領域が大きい音場推定装置、その方法及びプログラムを提供する。【解決手段】音場推定装置２００は、球体上の所定の位置にそれぞれマイクロホンを備えるM個の球面マイクロホンアレー１、２の周波数領域の収音信号u(ω,m,j)を用いて、音場を構成する平面波の強度からなるベクトルを推定する平面波分解部２１３と、平面波の強度からなるベクトルの推定値a(ω)と、仮想的なマイクロホンの位置r-とを用いて、仮想マイクロホンの位置r-での周波数領域の収音信号u^(ω,r-)を推定する補間推定部２１６とを含む。【選択図】図２
2017130899.tif 000002
【特許請求の範囲】
【請求項１】
m=1,2,…,Mを球面マイクロホンアレーのインデックス、jm=1,2,…,Jm、rmを極座標の動径、θj_m及びφj_mを極座標の偏角、ωを時間周波数のインデックスとし、d-mを中心とする半径rmの球体上の偏角θj_m及びφj_mのJm個の位置r-(m,jm)=d-m+[rmsinθj_mcosφj_m  rmsinθj_msinφj_m  rmcosθj_m]Tにそれぞれマイクロホンを備えるM個の球面マイクロホンアレーmの周波数領域の収音信号u(ω,m,jm)を用いて、音場を構成する平面波の強度からなるベクトルを推定する平面波分解部と、
前記平面波の強度からなるベクトルの推定値a(ω)と、仮想的なマイクロホンの位置r-pとを用いて、前記仮想マイクロホンの位置r-pでの周波数領域の収音信号u^(ω,r-p)を推定する補間推定部とを含む、
音場推定装置。
【請求項２】
請求項１の音場推定装置であって、
αを所定の実数、p=1,2,…,Pとし、rを極座標の動径、θp及びφpを極座標の偏角、収音信号u^(ω,r,θp,φp)の球面波スペクトルun,m(ω,r)の絶対値|un,m(ω,r)|と、収音信号u^(ω,αr,θp,φp)の球面波スペクトルun,m(ω,αr)の絶対値|un,m(ω,αr)|との大小関係に基づき、剛球型マイクロホンアレーの収音信号から算出される球面波スペクトルvn,m(ω,r)を推定するアレー型変換部を含む、
音場推定装置。
【請求項３】
請求項２の音場推定装置であって、仮想的な頭部の位置を剛球型マイクロホンアレーの位置として、
前記球面波スペクトルvn,m(ω,r)から得られるP個の時間領域の収音信号y(t,r-p)と仮想的な頭部の方向に基づき、仮想的な頭部の位置と方向におけるバイノーラル信号を生成するバイノーラル信号生成部を含む、
音場推定装置。
【請求項４】
請求項１から請求項３の何れかの音場推定装置であって、
Pを1以上の整数、p=1,2,…,P、a(ω)=[a1(ω),a2(ω),…,al(ω),…,aL'(ω)]T、iを虚数、kを波数とし、音場を構成するl番目の平面波の極座標の偏角をθl及びφlとし、前記補間推定部は、
【数２８】
2017130899.tif 000030
により、[dx dy dz]を中心とする半径rの球体上の偏角θp及びφpの仮想的なマイクロホンの周波数領域の収音信号u^(ω,r-p)を推定する、
音場推定装置
【請求項５】
請求項１から請求項４の何れかの音場推定装置であって、
前記M個の球面マイクロホンアレーmは、開球型のマイクロホンアレーであり、
||a(ω)||1をa(ω)のL1ノルムとし、正則化のための定数をλとし、前記平面波分解部は、
【数２９】
2017130899.tif 000031
をコスト関数とし、コスト関数Ｊを最小にする推定値a-(ω)を求める、
音場推定装置。
【請求項６】
請求項１から請求項４の何れかの音場推定装置であって、
前記M個の球面マイクロホンアレーは、剛球型のマイクロホンアレーであり、
||a(ω)||1をa(ω)のL1ノルムとし、正則化のための定数をλとし、Ynm'をオーダーn、次数m'の球調和関数とし、kを波数、bn(krm)をモード強度関数とし、前記平面波分解部は、
【数３０】
2017130899.tif 000032
をコスト関数とし、コスト関数Ｊを最小にする推定値a-(ω)を求める、
音場推定装置。
【請求項７】
m=1,2,…,Mを球面マイクロホンアレーのインデックス、jm=1,2,…,Jm、rmを極座標の動径、θj_m及びφj_mを極座標の偏角、ωを時間周波数のインデックスとし、平面波分解部が、d-mを中心とする半径rmの球体上の偏角θj_m及びφj_mのJm個の位置r-(m,jm)=d-m+[rmsinθj_mcosφj_m  rmsinθj_msinφj_m  rmcosθj_m]Tにそれぞれマイクロホンを備えるM個の球面マイクロホンアレーmの周波数領域の収音信号u(ω,m,jm)を用いて、音場を構成する平面波の強度からなるベクトルを推定する平面波分解ステップと、
補間推定部が、前記平面波の強度からなるベクトルの推定値a(ω)と、仮想的なマイクロホンの位置r-pとを用いて、前記仮想マイクロホンの位置r-pでの周波数領域の収音信号u^(ω,r-p)を推定する補間推定ステップとを含む、
音場推定方法。
【請求項８】
請求項１から請求項６の何れかの音場推定装置としてコンピュータを機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、ある位置に配置されたマイクロホンの収音信号を用いて、他の位置にマイクロホンが配置された場合に得られる収音信号を推定する技術に関する。
【背景技術】
【０００２】
近年、オーディオ再生技術は２チャネルステレオから５．１チャネル再生に拡大し、さらに２２．２チャネル再生や波面合成法の研究開発が進められ、再生そのものの臨場感を大きく向上させることと、臨場感の高い再生エリアをなるべく拡大することが図られている。
【０００３】
このような多チャネルオーディオ再生方法を評価検証するには、再生された音場を計測することが重要になる。例えば波面合成法では、実際に収録された音場と再現された音場とを比較し、その相違を把握する必要がある。その理由は、収録音場を再生信号へ変換する信号処理、収録した信号のエンコードとデコード、再生装置が設置された部屋の音響特性などの諸要因が音場の再現精度に影響するためであり、再現精度の高い手法を確立することが重要だからである。
【０００４】
（従来法１）
音場を計測する方法として、対象とする測定エリアの一部に局所的にマイクロホンを集中配置し、その測定結果から周辺エリアの音場を推定することが考えられる。一例として球面マイクロホンアレーの検討が進められている。球面マイクロホンアレーとは、数十以上のマイクロホン素子を半径raの球面上に配置して構成するマイクロホンアレーであり、raは数cmから十数cmの範囲にある。
【０００５】
図１に、従来技術における球面マイクロホンアレー１をもちいる音場推定処理の信号フローを示す。球面上に配置されたJ個のマイクロホンでそれぞれ収音された時間領域の信号y(t,ra,Ωj)は、短時間フーリエ変換部１１１により周波数領域の信号u(i,ω,ra,Ωj)に変換される。ただし、tは時刻、iはフレーム、Jは2以上の整数、ωは時間周波数、j=1,2,…,Jである。なお、以下の処理では、フレーム単位で処理を行うが、表記を簡略化するために、iを省略する。Ωjはj番目のマイクロホン素子の球面上の位置であり、elevation角θjとazimuth角φjのペアで指定される。Ωj=(θj,φj)である。
【０００６】
球面波スペクトル変換部１１２は、周波数ωごとに次式により球面波スペクトルun,m(ω,ra)を求める。
【０００７】
【数１】
2017130899.tif 000003
【０００８】
ただしαjは、式(1)の積和が次式で表される球調和関数の直交条件を満たすように適切に設定された重みである。
【０００９】
【数２】
2017130899.tif 000004
【００１０】
なお、Ynm(θj,φj)はオーダーn、次数mの球調和関数であり、*は複素共役を意味する。n=0,1,…,N、m=-n,-n+1,…,nである。δnn'はn=n'のときに１、n≠n'のときに0となる値であり、δmm'はm=m'のときに１、m≠m'のときに0となる値である。オーダー数Nまでの球面波スペクトルを得るには、(N+1)2個以上のマイクロホン素子が必要になる。
【００１１】
なお、これ以降は測定対象範囲の外側にある音源によって生成された音場を測定すること、すなわち内部問題を扱う。別の言い方をすると、球面マイクロホンアレーの成す球体の外側にある音源によって生成された音場を測定する。
【００１２】
球面マイクロホンアレーの中心を原点として、極座標系(r,Ω)=(r,θ,φ)で音場を考える。
【００１３】
外挿推定部１１６では、周波数ωで次式により、半径raの位置から半径rへと音場を外挿し、極座標系(r,Ω)=(r,θ,φ)における収音信号u(ω,r,Ω)を求める。言い換えると、球面マイクロホンアレー上に配置されたマイクロホンの収音信号を用いて球面マイクロホンアレーの成す球体の外側にある音場を推定する。
【００１４】
【数３】
2017130899.tif 000005
【００１５】
ただしkは波数k=ω/c(cは音速)であり、bn( )はモード強度関数である。
【００１６】
非特許文献１では、マイクロホン素子を中空で球面上に配置する開球型(open sphere)の球面マイクロホンアレーの場合が示されている。この場合、モード強度関数は次式で表される。
【００１７】
【数４】
2017130899.tif 000006
【００１８】
である。ただしiは虚数であり、jn( )はn次の球ベッセル関数である。球面マイクロホンアレーが剛球表面にマイク素子を配置して構成されている場合には、非特許文献２に基づき、モード強度関数は次式で表される。
【００１９】
【数５】
2017130899.tif 000007
【００２０】
ただしhn( )はn次の第1種ハンケル関数である。なお、「A'」はAの微分を意味する。
【００２１】
短時間逆フーリエ変換部１１８は、空間的に外挿した収音信号を周波数領域の信号u(ω,r,Ω)から時間領域の信号y(t,r,Ω)に戻し、出力する。
【００２２】
なお、式（３）は、球面波スペクトルにbn(kr)/bn(kra)を適用し、Ynm(θ,φ)で積和をとっている。このYnm(θ,φ)の積和は、逆球面波スペクトル変換に対応する。そのため、空間的に外挿した収音信号u(ω,r,Ω)は周波数領域の信号となっている。
【００２３】
開球型のマイクロホンアレーによる計測では、特異点の影響を避けられず、jn(kr)=0になるkおよびrで測定不能になる。具体的には、音場が存在してもｊn(kr)=0が満たされるとき、その出力が0になる。しかし剛球型のマイクロホンアレーには特異点がなく、測定不能にならない。そのため球面型マイクロホンアレーとしては、剛球型のマイクロホンアレーを使うことが主流である。
【先行技術文献】
【非特許文献】
【００２４】
【非特許文献１】T. Abhayapala and D. Ward, "Theory and design of high order  sound field microphones using spherical microphone array", in Acoustics, Speech, and Signal Processing (ICASSP), IEEE International Conference on, 2002, pp. II-1949.
【非特許文献２】Meyer, Jens; Elko, Gary, "A highly scalable spherical microphone array based on an orthonormal decomposition of the soundfield", Acoustics, Speech, and Signal Processing (ICASSP), 2002 IEEE International Conference on, 2002, II-1781 - II-1784.
【発明の概要】
【発明が解決しようとする課題】
【００２５】
従来技術では、半径raの位置から半径rへと音場を外挿する際に、ベッセル関数jn(kr)もしくはベッセル関数jn(kr)とハンケル関数hn(kr)とを使用している。参考文献１によれば、2関数とも大域的な傾向としてkrが増大すると、1/krのペースで減少する。
（参考文献１）Ｅ．Ｇ．ウィリアムズ、「フーリエ音響学」、シュプリンガー・フェアラーク、２００５、p.234-236．
例えばrがraの10倍になると、外挿の推定値は約1/10と急激に小さくなってしまう。そのため、外挿が有効な空間領域は球面マイクロホンアレー表面の周囲に限定されてしまう。また同じ理由により、周波数ωが高くなってもk=ω/cが大きくなり、外挿の推定値が急激に小さくなる。つまり周波数を高くすると外挿が有効な空間領域が急激に狭まってしまう。
【００２６】
本発明は、従来技術よりも、音場推定が有効な空間領域が大きい音場推定装置、その方法及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【００２７】
上記の課題を解決するために、本発明の一態様によれば、音場推定装置は、m=1,2,…,Mを球面マイクロホンアレーのインデックスとし、jm=1,2,…,Jm、rmを極座標の動径、θj_m及びφj_mを極座標の偏角、ωを時間周波数のインデックス、d-mを中心とする半径rmの球体上の偏角θj_m及びφj_mのJm個の位置r-(m,jm)=d-m+[rmsinθj_mcosφj_m  rmsinθj_msinφj_m  rmcosθj_m]Tにそれぞれマイクロホンを備えるM個の球面マイクロホンアレーmの周波数領域の収音信号u(ω,m,jm)を用いて、音場を構成する平面波の強度からなるベクトルを推定する平面波分解部と、平面波の強度からなるベクトルの推定値a(ω)と、仮想的なマイクロホンの位置r-pとを用いて、仮想マイクロホンの位置r-pでの周波数領域の収音信号u^(ω,r-p)を推定する補間推定部とを含む。
【００２８】
上記の課題を解決するために、本発明の他の態様によれば、音場推定方法は、m=1,2,…,Mを球面マイクロホンアレーのインデックス、jm=1,2,…,Jm、rmを極座標の動径、θj_m及びφj_mを極座標の偏角、ωを時間周波数のインデックスとし、平面波分解部が、d-mを中心とする半径rmの球体上の偏角θj_m及びφj_mのJm個の位置r-(m,jm)=d-m+[rmsinθj_mcosφj_m  rmsinθj_msinφj_m  rmcosθj_m]Tにそれぞれマイクロホンを備えるM個の球面マイクロホンアレーmの周波数領域の収音信号u(ω,m,jm)を用いて、音場を構成する平面波の強度からなるベクトルを推定する平面波分解ステップと、補間推定部が、平面波の強度からなるベクトルの推定値a(ω)と、仮想的なマイクロホンの位置r-pとを用いて、仮想マイクロホンの位置r-pでの周波数領域の収音信号u^(ω,r-p)を推定する補間推定ステップとを含む。
【発明の効果】
【００２９】
本発明によれば、従来技術よりも、音場推定が有効な空間領域が大きいという効果を奏する。
【図面の簡単な説明】
【００３０】
【図１】従来技術に係る音場推定装置の機能ブロック図。
【図２】第一実施形態に係る音場推定装置の機能ブロック図。
【図３】第一実施形態に係る音場推定装置の処理フローの例を示す図。
【図４】第一実施形態、その変形例１における仮想マイクロホンの位置の概要を示す図。
【図５】第二実施形態に係る音場推定装置の機能ブロック図。
【図６】第二実施形態に係る音場推定装置の処理フローの例を示す図。
【図７】第三実施形態に係る音場推定装置の機能ブロック図。
【図８】第三実施形態に係る音場推定装置の処理フローの例を示す図。
【発明を実施するための形態】
【００３１】
以下、本発明の実施形態について、説明する。なお、以下の説明に用いる図面では、同じ機能を持つ構成部や同じ処理を行うステップには同一の符号を記し、重複説明を省略する。以下の説明において、テキスト中で使用する記号「^」「-」等は、本来直前の文字の真上に記載されるべきものであるが、テキスト記法の制限により、当該文字の直後に記載する。式中においてはこれらの記号は本来の位置に記述している。また、ベクトルや行列の各要素単位で行われる処理は、特に断りが無い限り、そのベクトルやその行列の全ての要素に対して適用されるものとする。
【００３２】
＜第一実施形態のポイント＞
本実施形態では、(1)単一の球面マイクロホンアレーでなく、複数のマイクロホンアレーをもちい、(2)球面波スペクトルの代わりに、周波数領域の収音信号から音場を構成する平面波の集まりを直接求めて、この平面波をもちいて音場を推定する。(1)と(2)により、音場推定が有効な空間範囲を大きく広げることが可能になる。以下、その方法を説明する。
【００３３】
＜第一実施形態に係る音場推定装置２００＞
図２は第一実施形態に係る音場推定装置２００の機能ブロック図を、図３はその処理フローを示す。
【００３４】
音場推定装置２００は、短時間フーリエ変換部２１１、平面波分解部２１３、補間推定部２１６及び短時間逆フーリエ変換部２１８を含む。
【００３５】
音場推定装置２００は、球面マイクロホンアレー１から時間領域の収音信号y(t,1,j)(ただし、j=1,2,…,J)を受け取り、球面マイクロホンアレー２から時間領域の収音信号y(t,2,j)(ただし、j=1,2,…,J)を受け取り、仮想マイクロホンの位置情報r-を受け取り、仮想マイクロホンの位置r-における時間領域の収音信号y(t,r-)を推定し、出力する。
【００３６】
本実施形態では、球面マイクロホンアレー１、２が含む球面上のマイクロホン数、球面マイクロホンアレー１、２におけるマイクロホンの配置、球面マイクロホンアレー１、２の半径が同一であるとする。球面マイクロホンアレー１、２は、開球型の球面マイクロホンアレーであり、半径raの球面上にJ個のマイクロホンが配置され、球面上のマイクロホン配置はelevation角とazimuth角のペア(θj,φj)で指定される。球面マイクロホンアレー１の中心位置をd-1=[x1,y1,z1]とし、球面マイクロホンアレー２の中心位置をd-2=[x2,y2,z2]とする。球面マイクロホンアレーｍ（ｍ＝１，２）上の第ｊマイクロホンの３次元位置は、
【００３７】
【数６】
2017130899.tif 000008
【００３８】
で与えられる。このマイクロホンによる収音信号を、時間領域で収音信号y(t,m,j)(ただし、m=1,2、j=1,2,…,J)とする。
【００３９】
＜短時間フーリエ変換部２１１＞
短時間フーリエ変換部２１１は、時間領域の収音信号y(t,m,j)(ただし、m=1,2、j=1,2,…,J)を受け取り、短時間フーリエ変換により、時間領域の収音信号y(t,m,j)を周波数領域の収音信号u(i,ω,m,j)(ただし、iはフレーム番号、ω=1,2,…,F、j=1,2,…,J)に変換し（Ｓ２１１）、出力する。なお、以降の処理はフレームi毎に行うが、記載を簡略化するため、フレーム番号iを省略する。また、時間領域の信号を周波数領域の信号に変換する方法であれば、短時間フーリエ変換以外の方法を用いてもよい。
【００４０】
＜平面波分解部２１３＞
平面波分解部２１３は、周波数領域の収音信号u(ω,m,j)(ただし、ω=1,2,…,F、m=1,2、j=1,2,…,J)を受け取り、これらの値を用いて、音場を構成する平面波の強度からなるベクトルを推定し（Ｓ２１３）、推定値a-(ω)(ただし、ω=1,2,…,F)を出力する。例えば、平面波分解部２１３は、音場を構成する平面波の集まりを求めるために、次のコスト関数Ｊを最小にする解ベクトル(推定値)a-(ω)を求める。
【００４１】
【数７】
2017130899.tif 000009
【００４２】
なお、λは正則化のための定数であり、λを大きくするほどu(ω,m,j)にのる雑音に対し、推定をロバストにすることができる。また、式(7)中のD(ω)は辞書行列とよばれ、
【００４３】
【数８】
2017130899.tif 000010
【００４４】
であり、そのl番目の列ベクトルD(ω,l)は、elevation角とazimuth角のペア(θl,φl)で指定された方向から振幅１の平面波が入射したとき、原点での位相が0の状態での各球面マイクロホンアレー１，２での観測値からなるベクトルである。例えば、L'個の平面波が全方位から万遍なく取得できるようにL'個の入射角(θl,φl)を設定する。例えば、正多面体の頂点の方向から平面波が入射するようにL'個の入射角(θl,φl)を設定する。
【００４５】
辞書行列D(ω)のl番目の列ベクトルD(ω,l)は、入射角(θl,φl)の平面波が到来する際の球面マイクロホンアレーｍの第ｊマイクロホンの観測値
【００４６】
【数９】
2017130899.tif 000011
をもちいて、
【００４７】
【数１０】
2017130899.tif 000012
で表される。
【００４８】
また、解ベクトルa-(ω)のl番目の値は、l番目の平面波の振幅に対応する。a(ω)=[a1(ω),a2(ω),…,al(ω),…,aL'(ω)]Tとする。
【００４９】
上記のＬ１ノルムの正則化項をもちいる凸最適化は、解ベクトルa-(ω)として０を多く含むスパースなベクトルを導きだす。そのため、参考文献２で示されているように、あらかじめ想定する平面波の個数L'がマイクロホン数を大きく上回るような冗長な場合でも、平面波をうまく抽出することが可能である。
(参考文献２)A. Wabnitz, N. Epain, A. van.Shaik, C. Jin, "Reconstruction of spatial sound field using compressed sensing", in Acoustics, Speech, and Signal Processing (ICASSP), IEEE International Conference on, 2011.
【００５０】
＜補間推定部２１６＞
補間推定部２１６は、推定値a(ω)と仮想マイクロホンの位置情報r-=(rx,ry,rz)とを受け取り、仮想マイクロホンの位置r-での周波数領域の収音信号u^(ω,r-)(ただし、ω=1,2,…,F)を次式により推定し（Ｓ２１６）(言い換えると、解ベクトルa-(ω)をパラメータとする平面波モデルから仮想マイクロホンの出力値u^(ω,r-)を推定し)、出力する。
【００５１】
【数１１】
2017130899.tif 000013
【００５２】
ただし、●は内積を意味し、k-lはl番目の平面波の入射方向に対応する波数ベクトルであり、次式により表される。
【００５３】
【数１２】
2017130899.tif 000014
【００５４】
Tは転置を表す。なお、仮想マイクロホンの位置情報r-は、例えば、音場推定装置２００の利用者により入力される。
【００５５】
＜短時間逆フーリエ変換部２１８＞
短時間逆フーリエ変換部２１８は、周波数領域の収音信号u^(ω,r-)(ただし、ω=1,2,…,F)を受け取り、逆短時間フーリエ変換により、収音信号u^(ω,r-)を時間領域の収音信号y(t,r-)に変換し（Ｓ２１８）、出力する。なお、時間領域の信号を周波数領域の信号に変換する方法として、短時間フーリエ変換部２１１における変換方法に対応する方法を用いればよい。
【００５６】
＜効果＞
以上の構成により、従来技術よりも、音場推定が有効な空間領域が大きい音場推定装置を実現することができる。
【００５７】
＜変形例１＞
第一実施形態では、１つの仮想的なマイクロホンを想定し、その位置で収音される信号を推定した。しかし、当然、複数の位置に仮想的なマイクロホンを想定してもよい。また、仮想的なマイクロホンを同一の球面上に配置することで、半径rの開球型の仮想的なマイクロホンアレーを構成することができる。例えば、仮想的なマイクロホンアレーの中心が、原点から見てD=[dx dy dz]の位置にあり、P個の仮想的なマイクロホンを備える仮想的なマイクロホンアレーを構成した場合、音場推定装置２００は、P個の仮想マイクロホンの位置情報(r-p)(ただし、p=1,2,…,P)を受け取り、収音信号y(t,r-p)(ただし、p=1,2,…,P)を出力する。
【００５８】
仮想的なマイクロホンアレーの球面上のp番目のマイクロホン位置をr-p(ただし、p=1,2,…,P)とするとき、補間推定部２１６は、周波数領域の収音信号u^(ω,r-p)を次式により推定する。
【００５９】
【数１３】
2017130899.tif 000015
【００６０】
図４は、第一実施形態、その変形例１における仮想マイクロホンの位置の概要を示している。仮想的なマイクロホンアレーの中心[dx dy dz]が[rxry rz]、半径r=0、仮想的なマイクロホンアレーが備える仮想的なマイクロホンの個数Pが1のとき、第一実施形態となるため、第一実施形態は変形例１の一例とも言える（ただし、本変形例では半径rの球面上のマイクロホンが設けられる。一方、第一実施形態では半径r=0、つまり、球面上ではなく、点上にマイクロホンが設けられる）。
【００６１】
＜変形例２＞
第一実施形態では、開球型マイクロホンアレー２つを音場に設置する場合を説明した。本変形例では剛球型マイクロホンアレー２つを音場に設置する場合について説明する。
【００６２】
剛球型マイクロホンアレーは、半径ra、マイクロホン数Jとし、球面上のマイクロホン配置はelevation角とazimuth角のペア(θj,φj)で指定されるとする。elevation角とazimuth角のペア(θl、φl)で指定された方向から振幅１の平面波が入射するとき、音場は入射波と散乱波からなる。
【００６３】
剛球型マイクロホンアレーの中心と座標系の原点が一致するとき、第jマイクロホンで観測される信号は
【００６４】
【数１４】
2017130899.tif 000016
【００６５】
になる。第m(m=1,2)の剛球型マイクロホンアレーの中心が原点からdm離れている場合、その第jマイクロホンで観測される信号は、位相差を考慮して
【００６６】
【数１５】
2017130899.tif 000017
になる。
【００６７】
そこで式(9)の代わりに式(16)をもちいて、式(10)の辞書行列D(ω)のl番目の列ベクトルD(ω,l)を生成すれば、あとは同様に最適化問題を解くことで、剛球マイクロホンアレーの出力信号から平面波を抽出することができる。
【００６８】
なお式(15)は無限個の項をふくむため、実際には有限のnをもちい、数値計算によりvrigid(ω,l,m,j)を求める。ra=4cmのとき、n=10程度にとればよい。
【００６９】
＜その他の変形例＞
本実施形態では、マイクロホンアレー１、２におけるマイクロホンの配置、球面マイクロホンアレー１、２の半径が同一であるとしているが、異なってもよい。また、マイクロホンアレーの個数は２つである必要はなく、複数個であればよい。例えば、Mを2以上の整数の何れかとし、M個の球面マイクロホンアレーｍ（ｍ＝１，２，…，Ｍ）上の第ｊmマイクロホンの３次元位置は、球面マイクロホンアレーｍの中心をd-mとし、半径をrmとし、
【００７０】
【数１６】
2017130899.tif 000018
【００７１】
で与えられる。このマイクロホンによる収音信号を、時間領域で収音信号y(t,m,jm)(jm=1,2,…,Jmとし、Jmは球面マイクロホンアレーｍが備えるマイクロホンの個数)とする。また、式(7),(8)は次式に置き換える。
【００７２】
【数１７】
2017130899.tif 000019
【００７３】
辞書行列D(ω)のl番目の列ベクトルD(ω,l)は、入射角(θl,φl)の平面波が到来する際の球面マイクロホンアレーｍの第ｊmマイクロホンの観測値
【００７４】
【数１８】
2017130899.tif 000020
をもちいて、
【００７５】
【数１９】
2017130899.tif 000021
で表される。
【００７６】
また、第一実施形態の変形例２の式(15),(16)は次式に置き換える。
【００７７】
【数２０】
2017130899.tif 000022
【００７８】
＜第二実施形態＞
第一実施形態の変形例１と異なる部分を中心に説明する。
【００７９】
第一実施形態の変形例１では、仮想的に開球型のマイクロホンアレーを想定し、その収音信号を推定した。第二実施形態では、第一実施形態の変形例１の構成をベースとして、開球型のマイクロホンアレーの代わりに、仮想的に剛球型マイクロホンアレーを想定し、その収音信号を推定する。
【００８０】
図５は第二実施形態に係る音場推定装置３００の機能ブロック図を、図６はその処理フローを示す。
【００８１】
音場推定装置３００は、短時間フーリエ変換部２１１、平面波分解部２１３、補間推定部２１６及び短時間逆フーリエ変換部２１８を含み、さらに、アレー型変換部３１７を含む。
【００８２】
最初に、仮想的な球面マイクロホンアレーとして、参考文献３の二重式開球型の球面マイクロホンアレー(dual open sphere microphone array)による収音を想定する。このマイクロホンアレーではマイクロホン素子が半径rの球面もしくは半径αrの球面上に配置されており、α=1.2が推奨されている。
（参考文献３）I. Balmages, B. Rafaely, "Open-Sphere Designs for Spherical Microphone Arrays", IEEE Transactions on Audio, Speech, and Language Processing, vol. 15, no. 2, pp 727-732, 2007.
例えば、Q=P×2とし、Q個の仮想的なマイクロホン素子のうち、P個の仮想的なマイクロホン素子の位置を変形例１と同様とする。つまり、仮想的なマイクロホンアレーの中心が、原点から見てD=[dx dy dz]の位置にあり、仮想的なマイクロホンアレーの球面上のp番目の仮想的なマイクロホンの位置r-pは、
【００８３】
【数２１】
2017130899.tif 000023
【００８４】
である。Q個の仮想的なマイクロホン素子のうち、残りP個の仮想的なマイクロホンを中心が[dx dydz]、半径αrの球面上に配置し、q番目の仮想的なマイクロホンの位置を
【００８５】
【数２２】
2017130899.tif 000024
【００８６】
とする。また、Ωq=(θq,φq)=Ωp=(θp,φp)とする。つまり、q番目（ｑ＝P+p）のマイクロホンとp番目のマイクロホンとは、仮想的なマイクロホンアレーの中心からみて同じ方向にあり、q番目のマイクロホンは半径rの球面上にあり、q番目のマイクロホンは半径αrの球面上にある。
【００８７】
補間推定部２１６は、推定値Aと仮想的なマイクロホンアレーの中心Dと、仮想マイクロホンのP個の位置情報(r-p)(ただし、p=1,2,…,P)とP個の位置情報(r-q)(ただし、q=P+1,P+2,…,Q)とを受け取り、仮想マイクロホンの位置(r-p)及び(r-q)での周波数領域の収音信号u^(ω,r-p)(ただし、p=1,2,…,P)、u^(ω,r-q)(ただし、q=P+1,P+2,…,Q)を推定し（Ｓ２１６）、出力する。なお、P個の位置情報(r-q)(ただし、q=P+1,P+2,…,Q)に代えて、αのみを受け取る構成とし、P個の位置情報(r-p)とαとから、P個の位置情報(r-q)を計算して求めてもよい。
【００８８】
＜アレー型変換部３１７＞
アレー型変換部３１７は、周波数領域の収音信号u^(ω,r-p)(ただし、p=1,2,…,P)、u^(ω,r-q)(ただし、q=P+1,P+2,…,Q)を受け取り、次式により、球面波スペクトルun,m(ω,r)およびun,m(ω,αr)に変換する。
【００８９】
【数２３】
2017130899.tif 000025
【００９０】
開球型の球面マイクロホンアレーでは、特異点の影響によりjn(kr)=0になるkおよびrで測定が不可能になる。しかしun,m(ω,r)とun,m(ω,αr)のうち、絶対値の大きい方を選択することで、二重式開球型の球面マイクロホンアレーは特異点の影響を回避することができる。
【００９１】
そこで、アレー型変換部３１７は、
【００９２】
【数２４】
2017130899.tif 000026
【００９３】
を用い、
|un,m(ω,r)|>|un,m(ω,αr)|のとき
【００９４】
【数２５】
2017130899.tif 000027
【００９５】
とし、|un,m(ω,r)|≦|un,m(ω,αr)|のとき
【００９６】
【数２６】
2017130899.tif 000028
【００９７】
として、球面波スペクトルvn,m(ω,r)を求める。
【００９８】
アレー型変換部３１７は、最後に逆球面波スペクトル変換
【００９９】
【数２７】
2017130899.tif 000029
【０１００】
を適用する。これにより、最初に仮想的に設置した二重式開球型の球面マイクロホンアレーの位置に、半径ｒの剛球型マイクロホンアレーを設置した場合の収音信号を周波数領域で得ることができる。アレー型変換部３１７は周波数領域の信号v(ω,r-p)(ただし、p=1,2,…,P)を短時間逆フーリエ変換部２１８に出力する。
【０１０１】
＜効果＞
このような構成とすることで、第一実施形態の変形例１と同様の効果を得ることができる。さらに、剛球型のマイクロホンアレーを設置した場合の収音信号を仮想的に得ることができる。なお、本実施形態と第一実施形態の変形例２とを組合せてもよい。
【０１０２】
＜第三実施形態＞
剛球型マイクロホンアレーのバーチャルリアリティへの適用が、参考文献４で示されている。
（参考文献４）R. Duraiswami1, D. N. Zotkin1, Z. Li, E. Grassi, N. A. Gumerov, L. S. Davis, "High Order Spatial Audio Capture and Binaural Head-Tracked Playback over Headphones with HRTF Cues", Proceedings 119th convention of AES, 2005.
この参考文献４では、固定された剛球型マイクロホンアレーの収音信号および仮想的な頭部の方向を入力とし、指定方向に頭を向けたときに右耳および左耳に聞こえる信号（バイノーラル信号）を出力する方法が示されている。球面マイクロホンアレーが全方向に収音しているために、マイクロホン素子およびマイクロホンアレーを動かすことなく、指定された任意の方向に対応したバイノーラル信号を生成可能である。すなわち、受聴者の頭部回転をリアルタイムに計測して入力すると、その回転運動に追随してバイノーラル信号を生成して、受聴者に提示できる。
【０１０３】
第二実施形態では、仮想的に設置した剛球型マイクロホンアレーの収音信号を求める方法を示した。この収音信号に対して、図７のようにこのバイノーラル信号生成法を組み合わせる構成が、本実施形態の構成である。
【０１０４】
第二実施形態と異なる部分を中心に説明する。
【０１０５】
図７は第三実施形態に係る音場推定装置４００の機能ブロック図を、図８はその処理フローを示す。
【０１０６】
音場推定装置４００は、短時間フーリエ変換部２１１、平面波分解部２１３、補間推定部２１６、アレー型変換部３１７及び短時間逆フーリエ変換部２１８を含み、さらに、バイノーラル信号生成部４１９を含む。
【０１０７】
＜バイノーラル信号生成部４１９＞
バイノーラル信号生成部４１９は、仮想的な頭部の方向（姿勢）と時間領域の収音信号y(t,r-p)(ただし、p=1,2,…,P、剛球型球面マイクロホンアレーの収音信号に相当)とを受け取り、例えば参考文献４に記載の方法により、これらの信号から仮想的な頭部の位置と方向におけるバイノーラル信号y(t,R),y(t,L)を生成し（Ｓ４１９）、音場推定装置４００の出力値として出力する。なお、仮想的な頭部の位置は、仮想的なマイクロホンアレーの中心D=[dx dy dz]に相当し、時間領域の収音信号y(t,r-p)は、仮想的な頭部の位置における剛球型球面マイクロホンアレーの収音信号に相当する。そのため、バイノーラル信号生成部４１９では、仮想的な頭部の方向（姿勢）と時間領域の収音信号y(t,r-p)とから仮想的な頭部の位置と方向におけるバイノーラル信号y(t,R),y(t,L)を生成することができる。
【０１０８】
参考文献４の手法は頭での回転運動にしか追随できず、頭の並進運動には対応できない。しかし、本実施形態の構成では、剛球型球面マイクロホンアレーを仮想的に並進移動させることができる。そのために本実施形態は、頭部の回転運動および並進運動の両方に追随してバイノーラル信号を生成することを可能にする。
【０１０９】
＜その他の変形例＞
本発明は上記の実施形態及び変形例に限定されるものではない。例えば、上述の各種の処理は、記載に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。その他、本発明の趣旨を逸脱しない範囲で適宜変更が可能である。
【０１１０】
＜プログラム及び記録媒体＞
また、上記の実施形態及び変形例で説明した各装置における各種の処理機能をコンピュータによって実現してもよい。その場合、各装置が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、上記各装置における各種の処理機能がコンピュータ上で実現される。
【０１１１】
この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。
【０１１２】
また、このプログラムの流通は、例えば、そのプログラムを記録したＤＶＤ、ＣＤ−ＲＯＭ等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させてもよい。
【０１１３】
このようなプログラムを実行するコンピュータは、例えば、まず、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、一旦、自己の記憶部に格納する。そして、処理の実行時、このコンピュータは、自己の記憶部に格納されたプログラムを読み取り、読み取ったプログラムに従った処理を実行する。また、このプログラムの別の実施形態として、コンピュータが可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することとしてもよい。さらに、このコンピュータにサーバコンピュータからプログラムが転送されるたびに、逐次、受け取ったプログラムに従った処理を実行することとしてもよい。また、サーバコンピュータから、このコンピュータへのプログラムの転送は行わず、その実行指示と結果取得のみによって処理機能を実現する、いわゆるＡＳＰ（Application Service Provider）型のサービスによって、上述の処理を実行する構成としてもよい。なお、プログラムには、電子計算機による処理の用に供する情報であってプログラムに準ずるもの（コンピュータに対する直接の指令ではないがコンピュータの処理を規定する性質を有するデータ等）を含むものとする。
【０１１４】
また、コンピュータ上で所定のプログラムを実行させることにより、各装置を構成することとしたが、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。
【図１】
2017130899.tif 000033
【図２】
2017130899.tif 000034
【図３】
2017130899.tif 000035
【図４】
2017130899.tif 000036
【図５】
2017130899.tif 000037
【図６】
2017130899.tif 000038
【図７】
2017130899.tif 000039
【図８】
2017130899.tif 000040