(57)【要約】      （修正有）【課題】音源位置を特定する音源探索方法等に関し、従来の遅延和法と比べ分解能を向上させる。【解決手段】Ｚ＝ＺＳのｘ−ｙ平面上に存在する音源から発せられた音をＺ＝ＺＨのｘ−ｙ平面上に配列された複数のマイクロホンでピックアップすることにより生成された複数チャンネルの時系列データを取得し、取得された多チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより多チャンネルの所望の角周波数ωのスペクトルを生成し、生成した多チャンネルの角周波数ωのスペクトルに基づいて行列要素を算出し、算出された行列要素を連立方程式に代入してを算出する。【選択図】図１
2008286710.tif 000002
【特許請求の範囲】
【請求項１】
【数１】
2008286710.tif 000291
のｘ−ｙ平面上に存在する音源から発せられた音を
【数２】
2008286710.tif 000292
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された、複数チャンネルの時系列データを取得するデータ取得ステップと、
前記データ取得ステップで取得された複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【数３】
2008286710.tif 000293
のスペクトルを生成するスペクトル算出ステップと、
式
【数４】
2008286710.tif 000294
【数５】
2008286710.tif 000295
【数６】
2008286710.tif 000296
に、前記スペクトル生成ステップで算出された複数チャンネルの角周波数
【数７】
2008286710.tif 000297
のスペクトルを代入することにより、行列要素
【数８】
2008286710.tif 000298
【数９】
2008286710.tif 000299
【数１０】
2008286710.tif 000300
を算出する行列要素算出ステップと、
前記行列要素算出ステップで算出された行列要素
【数１１】
2008286710.tif 000301
【数１２】
2008286710.tif 000302
【数１３】
2008286710.tif 000303
を連立方程式
【数１４】
2008286710.tif 000304
に代入して該連立方程式（５．１）を解くことにより、
【数１５】
2008286710.tif 000305
を算出する音源探索ステップとを有することを特徴とする音源探索方法。
【請求項２】
観測点と計算点の離散化を同じ密度で行なうことにより、前記音源探索ステップが、前記連立方程式（５．１）に代えて、
【数１６】
2008286710.tif 000306
【数１７】
2008286710.tif 000307
【数１８】
2008286710.tif 000308
を連立方程式
【数１９】
2008286710.tif 000309
に代入して該連立方程式（５．２）を解くことにより、
【数２０】
2008286710.tif 000310
を算出するステップであることを特徴とする請求項１記載の音源探索方法。
【請求項３】
前記音源探索ステップが、
【数２１】
2008286710.tif 000311
を算出するステップであることを特徴とする請求項１又は２記載の音源探索方法。
【請求項４】
前記行列要素算出ステップが、行列要素
【数２２】
2008286710.tif 000312
を算出するにあたり、前記式（４．１２）に代えて、該式（４．１２）中の
【数２３】
2008286710.tif 000313
として、
【数２４】
2008286710.tif 000314
を採用した式
【数２５】
2008286710.tif 000315
に、前記スペクトル算出ステップで算出された複数チャンネルの角周波数
【数２６】
2008286710.tif 000316
のスペクトルを代入することにより、行列要素
【数２７】
2008286710.tif 000317
を算出するステップであることを特徴とする請求項１から３のうちいずれか１項記載の音源探索方法。
【請求項５】
前記行列要素算出ステップが、行列要素
【数２８】
2008286710.tif 000318
を算出するにあたり、前記式（４．１３）に代えて、該式（４．１３）中の
【数２９】
2008286710.tif 000319
として、
【数３０】
2008286710.tif 000320
を採用した式
【数３１】
2008286710.tif 000321
に、前記スペクトル算出ステップで算出された複数チャンネルの角周波数
【数３２】
2008286710.tif 000322
のスペクトルを代入することにより、行列要素
【数３３】
2008286710.tif 000323
を算出するステップであることを特徴とする請求項１から４のうちいずれか１項記載の音源探索方法。
【請求項６】
前記行列要素算出ステップが、行列要素
【数３４】
2008286710.tif 000324
を算出するにあたり、前記式（４．１４ａ）に代えて、該式（４．１４ａ）中の
【数３５】
2008286710.tif 000325
として
【数３６】
2008286710.tif 000326
を採用した式
【数３７】
2008286710.tif 000327
に、前記スペクトル算出ステップで算出された複数チャンネルの角周波数
【数３８】
2008286710.tif 000328
のスペクトルを代入することにより、行列要素
【数３９】
2008286710.tif 000329
を算出するステップであることを特徴とする請求項１から５のうちいずれか１項記載の音源探索方法。
【請求項７】
前記データ取得ステップが、
【数４０】
2008286710.tif 000330
のｘ−ｙ平面上に存在する音源から発せられた音を
【数４１】
2008286710.tif 000331
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するステップであることを特徴とする請求項１から５のうちのいずれか１項記載の音源探索方法。
【請求項８】
前記データ取得ステップが、
【数４２】
2008286710.tif 000332
のｘ−ｙ平面上に存在する音源から発せられた音を
【数４３】
2008286710.tif 000333
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するステップであることを特徴とする請求項１から５のうちのいずれか１項記載の音源探索方法。
【請求項９】
【数４４】
2008286710.tif 000334
のｘ−ｙ平面上に存在する音源から発せられた音を
【数４５】
2008286710.tif 000335
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された、複数チャンネルの時系列データを取得するデータ取得部と、
前記データ取得部で取得した複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【数４６】
2008286710.tif 000336
のスペクトルを算出するスペクトル算出部と、
式
【数４７】
2008286710.tif 000337
【数４８】
2008286710.tif 000338
【数４９】
2008286710.tif 000339
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数５０】
2008286710.tif 000340
のスペクトルを代入することにより、行列要素
【数５１】
2008286710.tif 000341
【数５２】
2008286710.tif 000342
【数５３】
2008286710.tif 000343
を算出する行列要素算出部と、
前記行列要素算出部で算出された行列要素
【数５４】
2008286710.tif 000344
【数５５】
2008286710.tif 000345
【数５６】
2008286710.tif 000346
を連立方程式
【数５７】
2008286710.tif 000347
に代入して該連立方程式（５．１）を解くことにより、
【数５８】
2008286710.tif 000348
を算出する音源探索部とを備えたことを特徴とする音源探索装置。
【請求項１０】
当該音源探索装置が、観測点と計算点の離散化を同じ密度で行なうものであって、前記音源探索部が、前記連立方程式（５．１）に代えて、
【数５９】
2008286710.tif 000349
【数６０】
2008286710.tif 000350
【数６１】
2008286710.tif 000351
を連立方程式
【数６２】
2008286710.tif 000352
に代入して該連立方程式（５．２）を解くことにより
【数６３】
2008286710.tif 000353
を算出するものであることを特徴とする請求項９記載の音源探索装置。
【請求項１１】
前記音源探索部が、
【数６４】
2008286710.tif 000354
を算出するものであることを特徴とする請求項９又は１０記載の音源探索装置。
【請求項１２】
前記行列要素算出部が、行列要素
【数６５】
2008286710.tif 000355
を算出するにあたり、前記式（４．１２）に代えて、該式（４．１２）中の
【数６６】
2008286710.tif 000356
として、
【数６７】
2008286710.tif 000357
を採用した式
【数６８】
2008286710.tif 000358
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数６９】
2008286710.tif 000359
のスペクトルを代入することにより、行列要素
【数７０】
2008286710.tif 000360
を算出するものであることを特徴とする請求項９から１１のうちいずれか１項記載の音源探索装置。
【請求項１３】
前記行列要素算出部が、行列要素
【数７１】
2008286710.tif 000361
を算出するにあたり、前記式（４．１３）に代えて、該式（４．１３）中の
【数７２】
2008286710.tif 000362
として、
【数７３】
2008286710.tif 000363
を採用した式
【数７４】
2008286710.tif 000364
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数７５】
2008286710.tif 000365
のスペクトルを代入することにより、行列要素
【数７６】
2008286710.tif 000366
を算出するものであることを特徴とする請求項９から１２のうちいずれか１項記載の音源探索装置。
【請求項１４】
前記行列要素算出部が、行列要素
【数７７】
2008286710.tif 000367
を算出するにあたり、前記式（４．１４ａ）に代えて、該式（４．１４ａ）中の
【数７８】
2008286710.tif 000368
として
【数７９】
2008286710.tif 000369
を採用した式
【数８０】
2008286710.tif 000370
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数８１】
2008286710.tif 000371
のスペクトルを代入することにより、行列要素
【数８２】
2008286710.tif 000372
を算出するものであることを特徴とする請求項９から１３のうちいずれか１項記載の音源探索装置。
【請求項１５】
前記データ取得部が、
【数８３】
2008286710.tif 000373
のｘ−ｙ平面上に存在する音源から発せられた音を
【数８４】
2008286710.tif 000374
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであることを特徴とする請求項９から１４のうちのいずれか１項記載の音源探索装置。
【請求項１６】
前記データ取得部が、
【数８５】
2008286710.tif 000375
のｘ−ｙ平面上に存在する音源から発せられた音を
【数８６】
2008286710.tif 000376
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであることを特徴とする請求項９から１４のうちのいずれか１項記載の音源探索装置。
【請求項１７】
プログラムを実行する演算処理装置内で実行され、該演算処理装置を、
【数８７】
2008286710.tif 000377
のｘ−ｙ平面上に存在する音源から発せられた音を
【数８８】
2008286710.tif 000378
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された、複数チャンネルの時系列データを取得するデータ取得部と、
前記データ取得部で取得した複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【数８９】
2008286710.tif 000379
のスペクトルを算出するスペクトル算出部と、
式
【数９０】
2008286710.tif 000380
【数９１】
2008286710.tif 000381
【数９２】
2008286710.tif 000382
に、前記スペクトル算出ステップで算出された複数チャンネルの角周波数
【数９３】
2008286710.tif 000383
のスペクトルを代入することにより、行列要素
【数９４】
2008286710.tif 000384
【数９５】
2008286710.tif 000385
【数９６】
2008286710.tif 000386
を算出する行列要素算出部と、
前記行列要素算出部で算出された行列要素
【数９７】
2008286710.tif 000387
【数９８】
2008286710.tif 000388
【数９９】
2008286710.tif 000389
を連立方程式
【数１００】
2008286710.tif 000390
に代入して該連立方程式（５．１）を解くことにより、
【数１０１】
2008286710.tif 000391
を算出する音源探索部とを備えた音源特定装置として動作させることを特徴とする音源探索プログラム。
【請求項１８】
当該音源探索プログラムが、観測点と計算点の離散化を同じ密度で行なわさせるものであって、前記音源探索部が、前記連立方程式（５．１）に代えて、
【数１０２】
2008286710.tif 000392
【数１０３】
2008286710.tif 000393
【数１０４】
2008286710.tif 000394
を連立方程式
【数１０５】
2008286710.tif 000395
に代入して該連立方程式（５．２）を解くことにより
【数１０６】
2008286710.tif 000396
を算出するものであることを特徴とする請求項１７記載の音源探索プログラム。
【請求項１９】
前記音源探索部が、
【数１０７】
2008286710.tif 000397
を算出するものであることを特徴とする請求項１７又は１８記載の音源探索プログラム。
【請求項２０】
前記行列要素算出部が、行列要素
【数１０８】
2008286710.tif 000398
を算出するにあたり、前記式（４．１２）に代えて、該式（４．１２）中の
【数１０９】
2008286710.tif 000399
として、
【数１１０】
2008286710.tif 000400
を採用した式
【数１１１】
2008286710.tif 000401
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数１１２】
2008286710.tif 000402
のスペクトルを代入することにより、行列要素
【数１１３】
2008286710.tif 000403
を算出するものであることを特徴とする請求項１７から１９のうちいずれか１項記載の音源探索プログラム。
【請求項２１】
前記行列要素算出部が、行列要素
【数１１４】
2008286710.tif 000404
を算出するにあたり、前記式（４．１３）に代えて、該式（４．１３）中の
【数１１５】
2008286710.tif 000405
として、
【数１１６】
2008286710.tif 000406
を採用した式
【数１１７】
2008286710.tif 000407
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数１１８】
2008286710.tif 000408
のスペクトルを代入することにより、行列要素
【数１１９】
2008286710.tif 000409
を算出するものであることを特徴とする請求項１７から２０のうちいずれか１項記載の音源探索プログラム。
【請求項２２】
前記行列要素算出部が、行列要素
【数１２０】
2008286710.tif 000410
を算出するにあたり、前記式（４．１４ａ）に代えて、該式（４．１４ａ）中の
【数１２１】
2008286710.tif 000411
として
【数１２２】
2008286710.tif 000412
を採用した式
【数１２３】
2008286710.tif 000413
に、前記スペクトル算出部で算出された複数チャンネルの角周波数
【数１２４】
2008286710.tif 000414
のスペクトルを代入することにより、行列要素
【数１２５】
2008286710.tif 000415
を算出するものであることを特徴とする請求項１７から２１のうちいずれか１項記載の音源探索プログラム。
【請求項２３】
前記データ取得部が、
【数１２６】
2008286710.tif 000416
のｘ−ｙ平面上に存在する音源から発せられた音を
【数１２７】
2008286710.tif 000417
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであることを特徴とする請求項１７から２２のうちのいずれか１項記載の音源探索プログラム。
【請求項２４】
前記データ取得部が、
【数１２８】
2008286710.tif 000418
のｘ−ｙ平面上に存在する音源から発せられた音を
【数１２９】
2008286710.tif 000419
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであることを特徴とする請求項１７から２２のうちのいずれか１項記載の音源探索プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音源位置を特定する音源探索方法および音源探索装置、並びに、コンピュータ等の演算処理装置内で実行されてその演算処理装置を音源探索装置として動作させる音源探索プログラムに関する。
【背景技術】
【０００２】
音源位置を特定する手法としては、古くから、遅延和法に代表されるマイクロホンアレイを利用した方法が知られている。この方法は、音源から離れた所で観測したデータを基に処理を行う便利さがある反面、分解能は波長のオーダー以下には望めないという欠点がある。この欠点を克服するために、音源のごく近く、すなわち、エバネッセント波が存在する近接場でデータを取得し、逆伝播により音源面での粒子速度分布を推定することにより放射音場特性を算出するという手法が提案されている（非特許文献１，２）。この手法は、分解能に波長の制限を受けないということでは、画期的な手法であるが、他方、粒子速度の算出点は、受音アレイ点上に限定されるとか、音波の波長が受音アレイの大きさに一致する周波数より上の周波数域でないと、アレイのゴーストが強く出て使用不可能になるという欠点（アレイ長による低周波限界）も合わせ持った手法である。
【０００３】
ここで、後述する式の導出過程や実施形態で参照するための文献として非特許文献３〜５、特許文献１を挙げておく。
【非特許文献１】Ｅ．Ｇ．Ｗｉｌｌｉａｍｓ，“Ｆｏｕｒｉｅｒ  Ａｃｏｕｓｔｉｃｓ；Ｓｏｕｎｄ  Ｒａｄｉａｔｉｏｎ  ａｎｄ  Ｎｅａｒｆｉｅｌｄ  Ａｃｏｕｓｔｉｃａｌ  Ｈｏｌｏｇｒａｐｈｙ，”（Ａｃａｄｅｍｉｃ  Ｐｒｅｓｓ，  Ｓａｎ  Ｄｉｅｇｏ，  １９９９  ）
【非特許文献２】Ｊ．Ｈａｌｄ，“Ｎｏｎ−ｓｔａｔｉｏｎａｒｙ  ＳＴＳＦ，”Ｂｒｕｅｌ  ＆  Ｋｊａｅｒ  Ｔｅｃｈｎｉｃａｌ  Ｒｅｖｉｅｗ，Ｎｏ．１，ｐｐ．１−３６（２０００）
【非特許文献３】Ａ．メシア、“量子力学１−３”、（東京図書、東京、１９８１−１９８２）（特に、１）
【非特許文献４】ディラック、“量子力学”、（岩波書店、東京、１９６８）
【非特許文献５】Ｌ．Ｓ．Ｒｏｄｂｅｒｇ  ａｎｄ  Ｒ．Ｍ．Ｔｈａｌｅｒ，“Ｉｎｔｒｏｄｕｃｔｉｏｎ  ｔｏ  ｔｈｅ  Ｑｕａｎｔｕｍ  Ｔｈｅｏｒｙ  ｏｆ  Ｓｃａｔｔｅｒｉｎｇ  ，”（Ａｃａｄｅｍｉｃ  Ｐｒｅｓｓ，Ｎｅｗ  Ｙｏｒｋ，１９６７）
【特許文献１】特開２００２−２７８９４８号公報
【発明の開示】
【発明が解決しようとする課題】
【０００４】
本発明は、上記事情に鑑み、エバネッセント波が存在する近接場よりも離れた場所でデータを取得し、かつ、従来の遅延和法と比べ分解能を向上させた音源探索方法および音源探索装置、並びに、コンピュータ等の演算処理装置内で実行されて、その演算処理装置を音源を高い分解能で特定する音源探索装置として動作させる音源探索プログラムを提供することを目的とする。
【課題を解決するための手段】
【０００５】
上記目的を達成する本発明の音源探索方法は、
【０００６】
【数１】
2008286710.tif 000003
【０００７】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０００８】
【数２】
2008286710.tif 000004
【０００９】
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された複数チャンネルの時系列データを取得するデータ取得ステップと、
データ取得ステップで取得された複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【００１０】
【数３】
2008286710.tif 000005
【００１１】
のスペクトルを生成するスペクトル算出ステップと、
式
【００１２】
【数４】
2008286710.tif 000006
【００１３】
【数５】
2008286710.tif 000007
【００１４】
【数６】
2008286710.tif 000008
【００１５】
に、スペクトル生成ステップで生成された複数チャンネルの角周波数
【００１６】
【数７】
2008286710.tif 000009
【００１７】
のスペクトルを代入することにより、行列要素
【００１８】
【数８】
2008286710.tif 000010
【００１９】
【数９】
2008286710.tif 000011
【００２０】
【数１０】
2008286710.tif 000012
【００２１】
を算出する行列要素算出ステップと、
行列要素算出ステップで算出された行列要素
【００２２】
【数１１】
2008286710.tif 000013
【００２３】
【数１２】
2008286710.tif 000014
【００２４】
【数１３】
2008286710.tif 000015
【００２５】
を連立方程式
【００２６】
【数１４】
2008286710.tif 000016
【００２７】
に代入して該連立方程式（５．１）を解くことにより、
【００２８】
【数１５】
2008286710.tif 000017
【００２９】
を算出する音源探索ステップとを有することを特徴とする。
【００３０】
ここで、観測点と計算点の離散化を同じ密度で行なうことにより、上記音源探索ステップが、連立方程式（５．１）に代えて、
【００３１】
【数１６】
2008286710.tif 000018
【００３２】
【数１７】
2008286710.tif 000019
【００３３】
【数１８】
2008286710.tif 000020
【００３４】
を連立方程式
【００３５】
【数１９】
2008286710.tif 000021
【００３６】
に代入してその連立方程式（５．２）を解くことにより
【００３７】
【数２０】
2008286710.tif 000022
【００３８】
を算出するステップであってもよい。
【００３９】
また、上記音源探索ステップが、
【００４０】
【数２１】
2008286710.tif 000023
【００４１】
を算出するステップであることが好ましい。
【００４２】
さらに、上記行列要素算出ステップが、行列要素
【００４３】
【数２２】
2008286710.tif 000024
【００４４】
を算出するにあたり、式（４．１２）に代えて、その式（４．１２）中の
【００４５】
【数２３】
2008286710.tif 000025
【００４６】
として、
【００４７】
【数２４】
2008286710.tif 000026
【００４８】
を採用した式
【００４９】
【数２５】
2008286710.tif 000027
【００５０】
に、スペクトル算出ステップで算出された複数チャンネルの角周波数
【００５１】
【数２６】
2008286710.tif 000028
【００５２】
のスペクトルを代入することにより、行列要素
【００５３】
【数２７】
2008286710.tif 000029
【００５４】
を算出するステップであることが好ましく、
また、上記行列要素算出ステップが、行列要素
【００５５】
【数２８】
2008286710.tif 000030
【００５６】
を算出するにあたり、式（４．１３）に代えて、その式（４．１３）中の
【００５７】
【数２９】
2008286710.tif 000031
【００５８】
として、
【００５９】
【数３０】
2008286710.tif 000032
【００６０】
を採用した式
【００６１】
【数３１】
2008286710.tif 000033
【００６２】
に、スペクトル算出ステップで算出された複数チャンネルの角周波数
【００６３】
【数３２】
2008286710.tif 000034
【００６４】
のスペクトルを代入することにより、行列要素
【００６５】
【数３３】
2008286710.tif 000035
【００６６】
を算出するステップであることが好ましく、
さらに、上記行列要素算出ステップが、行列要素
【００６７】
【数３４】
2008286710.tif 000036
【００６８】
を算出するにあたり、式（４．１４ａ）に代えて、その式（４．１４ａ）中の
【００６９】
【数３５】
2008286710.tif 000037
【００７０】
として
【００７１】
【数３６】
2008286710.tif 000038
【００７２】
を採用した式
【００７３】
【数３７】
2008286710.tif 000039
【００７４】
に、スペクトル算出ステップで算出された複数チャンネルの角周波数
【００７５】
【数３８】
2008286710.tif 000040
【００７６】
のスペクトルを代入することにより、行列要素
【００７７】
【数３９】
2008286710.tif 000041
【００７８】
を算出するステップであることが好ましい。
【００７９】
ここで、本発明の音源探索方法において、上記データ取得ステップは、
【００８０】
【数４０】
2008286710.tif 000042
【００８１】
のｘ−ｙ平面上に存在する音源から発せられた音を
【００８２】
【数４１】
2008286710.tif 000043
【００８３】
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するステップであってもよく、あるいは、
上記データ取得ステップは、
【００８４】
【数４２】
2008286710.tif 000044
【００８５】
のｘ−ｙ平面上に存在する音源から発せられた音を
【００８６】
【数４３】
2008286710.tif 000045
【００８７】
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するステップであってもよい。
【００８８】
また、上記目的を達成する本発明の音源探索装置は、
【００８９】
【数４４】
2008286710.tif 000046
【００９０】
のｘ−ｙ平面上に存在する音源から発せられた音を
【００９１】
【数４５】
2008286710.tif 000047
【００９２】
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された複数チャンネルの時系列データを取得するデータ取得部と、
データ取得部で取得した複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【００９３】
【数４６】
2008286710.tif 000048
【００９４】
のスペクトルを算出するスペクトル算出部と、
式
【００９５】
【数４７】
2008286710.tif 000049
【００９６】
【数４８】
2008286710.tif 000050
【００９７】
【数４９】
2008286710.tif 000051
【００９８】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【００９９】
【数５０】
2008286710.tif 000052
【０１００】
のスペクトルを代入することにより、行列要素
【０１０１】
【数５１】
2008286710.tif 000053
【０１０２】
【数５２】
2008286710.tif 000054
【０１０３】
【数５３】
2008286710.tif 000055
【０１０４】
を算出する行列要素算出部と、
行列要素算出部で算出された行列要素
【０１０５】
【数５４】
2008286710.tif 000056
【０１０６】
【数５５】
2008286710.tif 000057
【０１０７】
【数５６】
2008286710.tif 000058
【０１０８】
を連立方程式
【０１０９】
【数５７】
2008286710.tif 000059
【０１１０】
に代入して該連立方程式（５．１）を解くことにより、
【０１１１】
【数５８】
2008286710.tif 000060
【０１１２】
を算出する音源探索部とを備えたことを特徴とする。
【０１１３】
ここで、本発明の音源探索装置が、観測点と計算点の離散化を同じ密度で行なうものであって、上記音源探索部が、連立方程式（５．１）に代えて、
【０１１４】
【数５９】
2008286710.tif 000061
【０１１５】
【数６０】
2008286710.tif 000062
【０１１６】
【数６１】
2008286710.tif 000063
【０１１７】
を連立方程式
【０１１８】
【数６２】
2008286710.tif 000064
【０１１９】
に代入してその連立方程式（５．２）を解くことにより
【０１２０】
【数６３】
2008286710.tif 000065
【０１２１】
を算出するものであってもよい。
【０１２２】
また、上記音源探索部が、
【０１２３】
【数６４】
2008286710.tif 000066
【０１２４】
を算出するものであることが好ましい。
【０１２５】
さらに、上記行列要素算出部が、行列要素
【０１２６】
【数６５】
2008286710.tif 000067
【０１２７】
を算出するにあたり、記式（４．１２）に代えて、その式（４．１２）中の
【０１２８】
【数６６】
2008286710.tif 000068
【０１２９】
として、
【０１３０】
【数６７】
2008286710.tif 000069
【０１３１】
を採用した式
【０１３２】
【数６８】
2008286710.tif 000070
【０１３３】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０１３４】
【数６９】
2008286710.tif 000071
【０１３５】
のスペクトルを代入することにより、行列要素
【０１３６】
【数７０】
2008286710.tif 000072
【０１３７】
を算出するものであることが好ましく、
また、上記行列要素算出部が、行列要素
【０１３８】
【数７１】
2008286710.tif 000073
【０１３９】
を算出するにあたり、式（４．１３）に代えて、その式（４．１３）中の
【０１４０】
【数７２】
2008286710.tif 000074
【０１４１】
として、
【０１４２】
【数７３】
2008286710.tif 000075
【０１４３】
を採用した式
【０１４４】
【数７４】
2008286710.tif 000076
【０１４５】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０１４６】
【数７５】
2008286710.tif 000077
【０１４７】
のスペクトルを代入することにより、行列要素
【０１４８】
【数７６】
2008286710.tif 000078
【０１４９】
を算出するものであることが好ましく、
さらに、上記行列要素算出部が、行列要素
【０１５０】
【数７７】
2008286710.tif 000079
【０１５１】
を算出するにあたり、前記式（４．１４ａ）に代えて、その式（４．１４ａ）中の
【０１５２】
【数７８】
2008286710.tif 000080
【０１５３】
として
【０１５４】
【数７９】
2008286710.tif 000081
【０１５５】
を採用した式
【０１５６】
【数８０】
2008286710.tif 000082
【０１５７】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０１５８】
【数８１】
2008286710.tif 000083
【０１５９】
のスペクトルを代入することにより、行列要素
【０１６０】
【数８２】
2008286710.tif 000084
【０１６１】
を算出するものであることが好ましい。
【０１６２】
ここで、本発明の音源探索装置において、上記データ取得部は、
【０１６３】
【数８３】
2008286710.tif 000085
【０１６４】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０１６５】
【数８４】
2008286710.tif 000086
【０１６６】
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであってもよく、あるいは、
上記データ取得部は、
【０１６７】
【数８５】
2008286710.tif 000087
【０１６８】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０１６９】
【数８６】
2008286710.tif 000088
【０１７０】
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであってもよい。
【０１７１】
さらに、上記目的を達成する本発明の音源探索プログラムは、プログラムを実行する演算処理装置内で実行され、その演算処理装置を、
【０１７２】
【数８７】
2008286710.tif 000089
【０１７３】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０１７４】
【数８８】
2008286710.tif 000090
【０１７５】
のｘ−ｙ平面上に分散する複数の地点でピックアップすることにより生成された複数チャンネルの時系列データを取得するデータ取得部と、
データ取得部で取得した複数チャンネルの時系列データを各チャンネルごとに離散フーリエ変換することにより、複数チャンネルの所望の角周波数
【０１７６】
【数８９】
2008286710.tif 000091
【０１７７】
のスペクトルを算出するスペクトル算出部と、
式
【０１７８】
【数９０】
2008286710.tif 000092
【０１７９】
【数９１】
2008286710.tif 000093
【０１８０】
【数９２】
2008286710.tif 000094
【０１８１】
に、スペクトル生成ステップで生成された複数チャンネルの角周波数
【０１８２】
【数９３】
2008286710.tif 000095
【０１８３】
のスペクトルを代入することにより、行列要素
【０１８４】
【数９４】
2008286710.tif 000096
【０１８５】
【数９５】
2008286710.tif 000097
【０１８６】
【数９６】
2008286710.tif 000098
【０１８７】
を算出する行列要素算出部と、
行列要素算出部で算出された行列要素
【０１８８】
【数９７】
2008286710.tif 000099
【０１８９】
【数９８】
2008286710.tif 000100
【０１９０】
【数９９】
2008286710.tif 000101
【０１９１】
を連立方程式
【０１９２】
【数１００】
2008286710.tif 000102
【０１９３】
に代入して該連立方程式（５．１）を解くことにより、
【０１９４】
【数１０１】
2008286710.tif 000103
【０１９５】
を算出する音源探索部とを備えた音源特定装置として動作させることを特徴とする。
【０１９６】
ここで、本発明の音源探索プログラムが、観測点と計算点の離散化を同じ密度で行なわさせるものであって、上記音源探索部が、前記連立方程式（５．１）に代えて、
【０１９７】
【数１０２】
2008286710.tif 000104
【０１９８】
【数１０３】
2008286710.tif 000105
【０１９９】
【数１０４】
2008286710.tif 000106
【０２００】
を連立方程式
【０２０１】
【数１０５】
2008286710.tif 000107
【０２０２】
に代入してその連立方程式（５．２）を解くことにより
【０２０３】
【数１０６】
2008286710.tif 000108
【０２０４】
を算出するものであってもよい。
【０２０５】
また、上記音源探索部が、
【０２０６】
【数１０７】
2008286710.tif 000109
【０２０７】
を算出するものであることが好ましい。
【０２０８】
さらに、上記行列要素算出部が、行列要素
【０２０９】
【数１０８】
2008286710.tif 000110
【０２１０】
を算出するにあたり、式（４．１２）に代えて、その式（４．１２）中の
【０２１１】
【数１０９】
2008286710.tif 000111
【０２１２】
として、
【０２１３】
【数１１０】
2008286710.tif 000112
【０２１４】
を採用した式
【０２１５】
【数１１１】
2008286710.tif 000113
【０２１６】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０２１７】
【数１１２】
2008286710.tif 000114
【０２１８】
のスペクトルを代入することにより、行列要素
【０２１９】
【数１１３】
2008286710.tif 000115
【０２２０】
を算出するものであることが好ましく、
また、上記行列要素算出部が、行列要素
【０２２１】
【数１１４】
2008286710.tif 000116
【０２２２】
を算出するにあたり、式（４．１３）に代えて、その式（４．１３）中の
【０２２３】
【数１１５】
2008286710.tif 000117
【０２２４】
として、
【０２２５】
【数１１６】
2008286710.tif 000118
【０２２６】
を採用した式
【０２２７】
【数１１７】
2008286710.tif 000119
【０２２８】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０２２９】
【数１１８】
2008286710.tif 000120
【０２３０】
のスペクトルを代入することにより、行列要素
【０２３１】
【数１１９】
2008286710.tif 000121
【０２３２】
を算出するものであることが好ましく、
さらに、上記行列要素算出部が、行列要素
【０２３３】
【数１２０】
2008286710.tif 000122
【０２３４】
を算出するにあたり、前記式（４．１４ａ）に代えて、その式（４．１４ａ）中の
【０２３５】
【数１２１】
2008286710.tif 000123
【０２３６】
として
【０２３７】
【数１２２】
2008286710.tif 000124
【０２３８】
を採用した式
【０２３９】
【数１２３】
2008286710.tif 000125
【０２４０】
に、スペクトル算出部で算出された複数チャンネルの角周波数
【０２４１】
【数１２４】
2008286710.tif 000126
【０２４２】
のスペクトルを代入することにより、行列要素
【０２４３】
【数１２５】
2008286710.tif 000127
【０２４４】
を算出するものであることが好ましい。
【０２４５】
さらに、本発明の音源探索プログラムにおいても、上記データ取得部は、
【０２４６】
【数１２６】
2008286710.tif 000128
【０２４７】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０２４８】
【数１２７】
2008286710.tif 000129
【０２４９】
のｘ−ｙ平面上の前記複数の地点に置かれた複数のマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであってもよく、あるいは、
上記データ取得部は、
【０２５０】
【数１２８】
2008286710.tif 000130
【０２５１】
のｘ−ｙ平面上に存在する音源から発せられた音を
【０２５２】
【数１２９】
2008286710.tif 000131
【０２５３】
のｘ−ｙ平面上の前記複数の地点を経由して移動するマイクロホンでピックアップすることにより各地点に対応した複数チャンネルの時系列データを取得するものであってもよい。
【発明の効果】
【０２５４】
以上の本発明によれば、音源位置を高分解能で特定することができる。
【発明を実施するための最良の形態】
【０２５５】
１．はじめに
以下では、本発明の実施形態の説明に先立って、近接場音響ホログラフィーを逆問題として捉えた時の解法についてまとめておく。考え方は、一般的な作用素論によるものである。ここでは、重み付き二乗誤差を考えて定式化しておく。
【０２５６】
以下に、記号の一覧を示しておく。大体は、以下のように使い分けているが、文字にも限りがあるし、また、時に、同じ文字を別の意味で用いている所もあるかも知れない。が、前後関係から、混同は避けられるであろう。
【０２５７】
【数１３０】
2008286710.tif 000132
【０２５８】
：角周波数
【０２５９】
【数１３１】
2008286710.tif 000133
【０２６０】
：周波数
【０２６１】
【数１３２】
2008286710.tif 000134
【０２６２】
：微小振幅音波の伝播速度
【０２６３】
【数１３３】
2008286710.tif 000135
【０２６４】
：平衡状態での流体の密度
【０２６５】
【数１３４】
2008286710.tif 000136
【０２６６】
：自由空間中の微小振幅音波に対する波数
【０２６７】
【数１３５】
2008286710.tif 000137
【０２６８】
：流体密度の微小変動分、音圧
【０２６９】
【数１３６】
2008286710.tif 000138
【０２７０】
：粒子速度ベクトル
【０２７１】
【数１３７】
2008286710.tif 000139
【０２７２】
：ラグランジュ乗数
【０２７３】
【数１３８】
2008286710.tif 000140
【０２７４】
：２次元座標
【０２７５】
【数１３９】
2008286710.tif 000141
【０２７６】
を略記したもの
【０２７７】
【数１４０】
2008286710.tif 000142
【０２７８】
：２次元波数
【０２７９】
【数１４１】
2008286710.tif 000143
【０２８０】
を略記したもの
【０２８１】
【数１４２】
2008286710.tif 000144
【０２８２】
：２次元波数の絶対値、すなわち
【０２８３】
【数１４３】
2008286710.tif 000145
【０２８４】
【数１４４】
2008286710.tif 000146
【０２８５】
：作用素
【０２８６】
【数１４５】
2008286710.tif 000147
【０２８７】
のエルミート共役作用素
２．音波伝播の記述
今、
【０２８８】
【数１４６】
2008286710.tif 000148
【０２８９】
の面上に音源があるとする。
【０２９０】
【数１４７】
2008286710.tif 000149
【０２９１】
での音場は、
【０２９２】
【数１４８】
2008286710.tif 000150
【０２９３】
での速度ポテンシャル
【０２９４】
【数１４９】
2008286710.tif 000151
【０２９５】
を用いて、
【０２９６】
【数１５０】
2008286710.tif 000152
【０２９７】
【数１５１】
2008286710.tif 000153
【０２９８】
のように表される。ただし、
【０２９９】
【数１５２】
2008286710.tif 000154
【０３００】
は音圧、
【０３０１】
【数１５３】
2008286710.tif 000155
【０３０２】
は粒子速度、
【０３０３】
【数１５４】
2008286710.tif 000156
【０３０４】
は静止流体の密度、
【０３０５】
【数１５５】
2008286710.tif 000157
【０３０６】
は微小振幅音波に対する伝播速度、
【０３０７】
【数１５６】
2008286710.tif 000158
【０３０８】
は角周波数である。
【０３０９】
【数１５７】
2008286710.tif 000159
【０３１０】
での速度ポテンシャルは、
【０３１１】
【数１５８】
2008286710.tif 000160
【０３１２】
のように、
【０３１３】
【数１５９】
2008286710.tif 000161
【０３１４】
の面上の音圧、あるいは、粒子速度の
【０３１５】
【数１６０】
2008286710.tif 000162
【０３１６】
成分を用いてスペクトル表示される。ただし、
【０３１７】
【数１６１】
2008286710.tif 000163
【０３１８】
【数１６２】
2008286710.tif 000164
【０３１９】
、また、
【０３２０】
【数１６３】
2008286710.tif 000165
【０３２１】
である。従って、
【０３２２】
【数１６４】
2008286710.tif 000166
【０３２３】
での音場は、
【０３２４】
【数１６５】
2008286710.tif 000167
【０３２５】
の面上の音圧、あるいは、粒子速度の
【０３２６】
【数１６６】
2008286710.tif 000168
【０３２７】
成分が、何らかの方法で判れば、後は、計算で出せるということになる。しかし、
【０３２８】
【数１６７】
2008286710.tif 000169
【０３２９】
の面上の音圧、あるいは、粒子速度の
【０３３０】
【数１６８】
2008286710.tif 000170
【０３３１】
成分を直接測ることは難しいことが多く、通常は、
【０３３２】
【数１６９】
2008286710.tif 000171
【０３３３】
を満たす
【０３３４】
【数１７０】
2008286710.tif 000172
【０３３５】
の面上で音圧を観測し、このデータから、
【０３３６】
【数１７１】
2008286710.tif 000173
【０３３７】
の面上の音圧、あるいは、粒子速度の
【０３３８】
【数１７２】
2008286710.tif 000174
【０３３９】
成分を推定するという方法を採る。
【０３４０】
特別の場合として、（半）無限空間中の音場を考える限りにおいて、ある位置の音場は、他の平面上の音場を用いて、一意的に表される。今考えている、
【０３４１】
【数１７３】
2008286710.tif 000175
【０３４２】
での音圧は、
【０３４３】
【数１７４】
2008286710.tif 000176
【０３４４】
の平面上の音圧、あるいは、粒子速度の
【０３４５】
【数１７５】
2008286710.tif 000177
【０３４６】
成分を用いて、
【０３４７】
【数１７６】
2008286710.tif 000178
【０３４８】
のように表される。（２．６）は、スペクトルで議論するならば、
【０３４９】
【数１７７】
2008286710.tif 000179
【０３５０】
を意味している。すなわち、
【０３５１】
【数１７８】
2008286710.tif 000180
【０３５２】
から
【０３５３】
【数１７９】
2008286710.tif 000181
【０３５４】
への伝播の様子は、音圧・音圧の組み合わせなら、
【０３５５】
【数１８０】
2008286710.tif 000182
【０３５６】
で記述され、粒子速度・音圧の組み合わせならば、
【０３５７】
【数１８１】
2008286710.tif 000183
【０３５８】
で表されることが判る。すなわち、（２．７）は、
【０３５９】
【数１８２】
2008286710.tif 000184
【０３６０】
と書ける。（２．１０）の関係は、ブラ・ケットを用いた抽象表示（非特許文献３〜５参照）を採用するならば、
【０３６１】
【数１８３】
2008286710.tif 000185
【０３６２】
と書ける。（２．１０）は（２．１１）をスペクトル展開した時の成分を示したものである。（２．１１）は、まとめて、
【０３６３】
【数１８４】
2008286710.tif 000186
【０３６４】
と書ける。
【０３６５】
ここで、
【０３６６】
【数１８５】
2008286710.tif 000187
【０３６７】
は、
【０３６８】
【数１８６】
2008286710.tif 000188
【０３６９】
と
【０３７０】
【数１８７】
2008286710.tif 000189
【０３７１】
を代表させた表記である。
【０３７２】
（２．１２）の解は、
【０３７３】
【数１８８】
2008286710.tif 000190
【０３７４】
の一般化逆作用素
【０３７５】
【数１８９】
2008286710.tif 000191
【０３７６】
を用いて、
【０３７７】
【数１９０】
2008286710.tif 000192
【０３７８】
と表される。この一般化逆作用素
【０３７９】
【数１９１】
2008286710.tif 000193
【０３８０】
の表現として、次節以降では、条件付最小二乗法に基づく表現を与えることにする。
３．最小二乗解（一般論）
【０３８１】
【数１９２】
2008286710.tif 000194
【０３８２】
で音圧を観測するとして、観測値との誤差を考えるならば、（２．１２）から、
【０３８３】
【数１９３】
2008286710.tif 000195
【０３８４】
を得、重み付き二乗誤差は、
【０３８５】
【数１９４】
2008286710.tif 000196
【０３８６】
となる。ただし、
【０３８７】
【数１９５】
2008286710.tif 000197
【０３８８】
は重みを与える正定値作用素である。ここで記号
【０３８９】
【数１９６】
2008286710.tif 000198
【０３９０】
はエルミート共役を意味する。最小二乗問題とは、（３．２）の誤差を最小にするように、
【０３９１】
【数１９７】
2008286710.tif 000199
【０３９２】
での物理量
【０３９３】
【数１９８】
2008286710.tif 000200
【０３９４】
を推定するという問題である。以上は、何の制限もない場合であったが、ｉｌｌ−ｃｏｎｄｉｔｉｏｎに近い問題では、解の安定性を保証するために、最小ノルム解とか、滑らかな解とかいった条件付きでの解も考えられる。そのために、
【０３９５】
【数１９９】
2008286710.tif 000201
【０３９６】
なる作用素
【０３９７】
【数２００】
2008286710.tif 000202
【０３９８】
を考え、
【０３９９】
【数２０１】
2008286710.tif 000203
【０４００】
なる量を最小ならしめる解を考えることが逆問題ではよく行われる。ここで、作用素
【０４０１】
【数２０２】
2008286710.tif 000204
【０４０２】
は、非負定値作用素になっていることに注意されたい。逆問題では、（３．２）の値が一定になるという条件の下で、（３．４）の値を最小にするような解を考えるのである。この問題は、ラグランジアン、
【０４０３】
【数２０３】
2008286710.tif 000205
【０４０４】
あるいは、
【０４０５】
【数２０４】
2008286710.tif 000206
【０４０６】
に対する極値問題として定式化され、最小二乗解は、（３．６）を
【０４０７】
【数２０５】
2008286710.tif 000207
【０４０８】
で微分してゼロと置くことにより定式化される
【０４０９】
【数２０６】
2008286710.tif 000208
【０４１０】
の解、
【０４１１】
【数２０７】
2008286710.tif 000209
【０４１２】
として求められる。ここで出てきたラグランジュ乗数は、（３．７）の解（３．８）を（３．２）へ入れた値が一定値になるという条件、
【０４１３】
【数２０８】
2008286710.tif 000210
【０４１４】
、から定まることになるが、実際には、経験的に定数を用いることが多い。
４．積分方程式の具体的な導出
（３．７）を具体的に解くにあたり、問題をはっきりさせておく。今の問題は、
【０４１５】
【数２０９】
2008286710.tif 000211
【０４１６】
で音圧、
【０４１７】
【数２１０】
2008286710.tif 000212
【０４１８】
を観測値として得たとして、
【０４１９】
【数２１１】
2008286710.tif 000213
【０４２０】
での物理量、
【０４２１】
【数２１２】
2008286710.tif 000214
【０４２２】
を求める問題とする。そのため、ブラ・ケットに関する公式（非特許文献３〜５参照）、
【０４２３】
【数２１３】
2008286710.tif 000215
【０４２４】
【数２１４】
2008286710.tif 000216
【０４２５】
【数２１５】
2008286710.tif 000217
【０４２６】
【数２１６】
2008286710.tif 000218
【０４２７】
【数２１７】
2008286710.tif 000219
【０４２８】
【数２１８】
2008286710.tif 000220
【０４２９】
及び、
【０４３０】
【数２１９】
2008286710.tif 000221
【０４３１】
に注意する。これを用いて（３．７）を変形すると、まず、次式を得る。
【０４３２】
【数２２０】
2008286710.tif 000222
【０４３３】
（４．１０）は、更に、
【０４３４】
【数２２１】
2008286710.tif 000223
【０４３５】
と書き直しうる。（４．１１）は、未知関数
【０４３６】
【数２２２】
2008286710.tif 000224
【０４３７】
に関する第１種Ｆｒｅｄｈｏｌｍ型積分方程式である。（４．１１）に於ける積分核などは、
【０４３８】
【数２２３】
2008286710.tif 000225
【０４３９】
【数２２４】
2008286710.tif 000226
【０４４０】
【数２２５】
2008286710.tif 000227
【０４４１】
【数２２６】
2008286710.tif 000228
【０４４２】
で与えられるが、この具体的な形は、未知関数
【０４４３】
【数２２７】
2008286710.tif 000229
【０４４４】
として何を選ぶか、また、
【０４４５】
【数２２８】
2008286710.tif 000230
【０４４６】
として何を採用するかによって変わってくる。これについては後述する。ただ、重要な事は、（４．７）、（４．８）、
【０４４７】
【数２２９】
2008286710.tif 000231
【０４４８】
【数２３０】
2008286710.tif 000232
【０４４９】
を用いると、
【０４５０】
【数２３１】
2008286710.tif 000233
【０４５１】
であるから、
【０４５２】
【数２３２】
2008286710.tif 000234
【０４５３】
であるならば、
【０４５４】
【数２３３】
2008286710.tif 000235
【０４５５】
となり、（４．１１）は対称核の積分方程式になることである。この場合には、Ｈｉｌｂｅｒｔ−Ｓｃｈｍｉｄｔの定理、を利用できることになる。また、（３．４）の下で示したように、
【０４５６】
【数２３４】
2008286710.tif 000236
【０４５７】
であるから、
【０４５８】
【数２３５】
2008286710.tif 000237
【０４５９】
に選んでおけば、（４．１１）の積分核は正定値になること、従って、Ｍｅｒｃｅｒの定理を利用できることなども附記しておく。
５．積分方程式の解法（数値解法）
（４．１１）の積分方程式を数値的に解くのに最も単純な方法は、積分を離散近似で評価し、（４．１１）を、
【０４６０】
【数２３６】
2008286710.tif 000238
【０４６１】
とし、連立方程式の解として、近似解を考える方法である。この近似では微小面積比などの形で比例因子を必要とするので、得られた解に補正係数を掛けるなどの処置をする必要がある。もし、観測点と計算点の離散化を同じ密度で行うならば、（５．１）は、
【０４６２】
【数２３７】
2008286710.tif 000239
【０４６３】
と出来、補正係数は不要となる。ただし、この場合には、
【０４６４】
【数２３８】
2008286710.tif 000240
【０４６５】
となっている。
【０４６６】
ここでは、この数値解法のみを考えることとし、積分方程式の解法には触れないことにする。
６．積分核などの具体形
ここでは、第４節で定義した関数、
【０４６７】
【数２３９】
2008286710.tif 000241
【０４６８】
【数２４０】
2008286710.tif 000242
【０４６９】
【数２４１】
2008286710.tif 000243
【０４７０】
の具体形を与えておく。ここでは、考える未知関数
【０４７１】
【数２４２】
2008286710.tif 000244
【０４７２】
としては、
【０４７３】
【数２４３】
2008286710.tif 000245
【０４７４】
での粒子速度の
【０４７５】
【数２４４】
2008286710.tif 000246
【０４７６】
成分
【０４７７】
【数２４５】
2008286710.tif 000247
【０４７８】
を考えることにする。この時、
【０４７９】
【数２４６】
2008286710.tif 000248
【０４８０】
としては、（２．９）の
【０４８１】
【数２４７】
2008286710.tif 000249
【０４８２】
を用いるべきとなる。（２．４）、（２．５）、
【０４８３】
【数２４８】
2008286710.tif 000250
【０４８４】
【数２４９】
2008286710.tif 000251
【０４８５】
を思い出すと、
【０４８６】
【数２５０】
2008286710.tif 000252
【０４８７】
は
【０４８８】
【数２５１】
2008286710.tif 000253
【０４８９】
の関数であることが判る。（２．４）、（２．９）、（４．１２）の形を考えると、今の数値積分に適した関数
【０４９０】
【数２５２】
2008286710.tif 000254
【０４９１】
としては、
【０４９２】
【数２５３】
2008286710.tif 000255
【０４９３】
を採用すればよいことが判る。これを採用することにより積分が収束することになる。ただし、定数
【０４９４】
【数２５４】
2008286710.tif 000256
【０４９５】
は積分が無次元となるように便宜上付けた。この時、
【０４９６】
【数２５５】
2008286710.tif 000257
【０４９７】
【数２５６】
2008286710.tif 000258
【０４９８】
を得る。残る
【０４９９】
【数２５７】
2008286710.tif 000259
【０５００】
の計算であるが、作用素
【０５０１】
【数２５８】
2008286710.tif 000260
【０５０２】
の形には任意性があったことを考慮に入れると、
【０５０３】
【数２５９】
2008286710.tif 000261
【０５０４】
そのものを、その正値性さえ考慮に入れれば、任意に選べることが判る。ここでは、無限積分の収束性も考慮に入れて、
【０５０５】
【数２６０】
2008286710.tif 000262
【０５０６】
と選ぶことにする。この時、
【０５０７】
【数２６１】
2008286710.tif 000263
【０５０８】
を得る。ただし、
【０５０９】
【数２６２】
2008286710.tif 000264
【０５１０】
である。
【０５１１】
実際の逆問題解法では、以上の式を数値積分して係数行列を求め、連立方程式を解くことになる。
【０５１２】
以上で本発明の理論的な基礎となる、近接場音響ホログラフィーを逆問題として捉えたときの解法の説明を終了し、以下、本発明の実施形態について説明する。
【０５１３】
図１は、音源とマイクロホンアレイのレイアウトを示す図である。
【０５１４】
ここでは、
【０５１５】
【数２６３】
2008286710.tif 000265
【０５１６】
のｘ−ｙ平面上の互いに近接した位置に２つの点音源１１，１２が並んでいる。また、
【０５１７】
【数２６４】
2008286710.tif 000266
【０５１８】
のｘ−ｙ平面上の、所定のピッチで二次元的に分散する多数の地点にそれぞれ配置された多数のマイクロホン２１，２２，２３，…からなるマイクロホンアレイ２０が置かれている。
【０５１９】
ここでは、音源１１，１２からの音をこのマイクロホンアレイ２０でピックアップして得た信号が解析されて、音源１１，１２の位置が特定される。
【０５２０】
尚、ここでの例は、マイクロホンアレイ２０を備えているが、マイクロホンアレイ２０を備える代わりに、１本もしくは複数本のマイクロホンを上記の多数の地点を経由する経路で移動させながら、音源１１，１２からの音をピックアップしてもよい。
【０５２１】
図２は、本発明の音源探索装置の一実施形態として動作する音源探索用コンピュータの外観斜視図である。本発明の一実施形態としての音源探索装置は、この音源探索用コンピュータ１００のハードウェアとその内部で実行されるソフトウェアとの組合せにより実現されている。
【０５２２】
この音源探索用コンピュータ１００は、ＣＰＵ、ＲＡＭ、磁気ディスク、通信用ボード等を内蔵した本体１０１、本体１０１からの指示により表示画面１０２ａ上に画面表示を行なうディスプレイ１０２、この音源探索用コンピュータ１００内に、オペレータの指示や文字情報を入力するためのキーボード１０３、表示画面１０２ａ上の任意の位置を指定することによりその位置に表示されているアイコン等に応じた指示を入力するマウス１０４を備えている。
【０５２３】
本体１０１には、ＣＤ−ＲＯＭ１０５（図２参照）が取り出し自在に装填されて装填されたＣＤ−ＲＯＭ１０５をドライブするＣＤ−ＲＯＭドライブも内蔵されている。
【０５２４】
ここでは、ＣＤ−ＲＯＭ１０５に、本発明の一実施形態としての音源探索プログラムが記憶されており、そのＣＤ−ＲＯＭ１０５が本体１０１内に装填され、ＣＤ−ＲＯＭドライブによりそのＣＤ−ＲＯＭ１０５に記憶された音源探索プログラムが音源探索用コンピュータ１００の磁気ディスク内にインストールされる。音源探索用コンピュータ１００の磁気ディスク内にインストールされた音源探索プログラムが起動されると、この音源探索用コンピュータ１００は、本発明の音源探索装置の一実施形態として動作する。
【０５２５】
本実施形態では、この音源探索用コンピュータ１００には、信号線１０５を介して、図１に示すマイクロホンアレイ２０を構成する多数のマイクロホン２１，２２，２３，…が接続されている。
【０５２６】
図３は、図２に示す音源探索用コンピュータ１００のハードウェア構成図である。
【０５２７】
このハードウェア構成図には、中央演算処理装置（ＣＰＵ）１１１、ＲＡＭ１１２、磁気ディスクコントローラ１１３、ＣＤ−ＲＯＭドライブ１１５、マウスコントローラ１１６、キーボードコントローラ１１７、ディスプレイコントローラ１１８、およびＡ／Ｄ変換ボード１２０が示されており、それらはバス１１０で相互に接続されている。
【０５２８】
Ａ／Ｄ変換ボード１２０上には、フィルタ１２１およびＡ／Ｄ変換器１２２が搭載されており、このＡ／Ｄ変換ボード１２０は、図１に示す音源探索用コンピュータ１００の本体１０１内に装填されている。
【０５２９】
このＡ／Ｄ変換ボード１２０上のＡ／Ｄ変換器１２２は、マイクロホンアレイ２０を構成する各マイクロホン２１，２２，２３，…で音声をピックアップして得られたアナログの受音信号を所定の標本化時間間隔（サンプリング間隔）で標本化（サンプリング）してディジタルの時系列データを生成するものであり、またこのＡ／Ｄ変換ボード１２０に搭載されたフィルタ１２１は、Ａ／Ｄ変換器１２２により得られるディジタルの時系列データにエイリアジングノイズ（折り返し雑音）が混入するのを防止するための、低域通過型のアンチエイリアジングフィルタである。
【０５３０】
ＣＤ−ＲＯＭドライブ１１５は、図２を参照して説明したように、ＣＤ−ＲＯＭ１０５が装填され、装填されたＣＤ−ＲＯＭ１０５をアクセスするものである。
【０５３１】
また、図３には、磁気ディスクコントローラ１１３によりアクセスされる磁気ディスク１１４、マウスコントローラ１１６により制御されるマウス１０４、キーボードコントローラ１１７により制御されるキーボード１０３、およびディスプレイコントローラ１１８により制御されるディスプレイ１０２も示されている。
【０５３２】
図４は、本発明の一実施形態としての音源探索プログラムの構成を示す図である。
【０５３３】
ここでは、この音源探索プログラム２００は、ＣＤ−ＲＯＭ１０５に記憶されており、この音源探索プログラム２００は、データ取得部２１０、スペクトル算出部２２０、行列要素算出部２３０、および音源探索部２４０から構成されている。
【０５３４】
この音源探索プログラム２００を構成する各部の作用については後述する。
【０５３５】
図５は、本発明の一実施形態としての音源探索方法の構成を示すフローチャートである。
【０５３６】
この音源探索方法は、データ取得ステップ（ステップａ）、スペクトル算出ステップ（ステップｂ）、行列要素算出ステップ（ステップｃ）、および音源探索ステップ（ステップｄ）とから構成されている。
【０５３７】
図５に示す音源探索方法は、本実施形態では、図４に示す音源探索プログラム２００が図２、図３に示す音源探索用コンピュータ１００にインストールされて実行されることにより、その音源探索プログラム２００により実施される。
【０５３８】
この音源探索方法の各ステップの作用については、以下に説明する図６の音源探索装置の説明の際に、図４の音源探索プログラム２００の各部の作用も合わせて一緒に説明する。
【０５３９】
図６は、図２、図３に示す音源探索用コンピュータ１００に図４に示す音源探索プログラム２００がインストールされて実行されることにより実現された、本発明の音源探索装置の一実施形態の機能ブロック図である。
【０５４０】
この図６に示す音源探索装置３００には、アンチエイリアジングフィルタ３０１、アナログ・デジタル変換部３０２、データ取得部３１０、スペクトル算出部３２０、行列要素算出部３３０、および音源探索部３４０から構成されている。
【０５４１】
アンチエイリアジングフィルタ３０１およびアナログ・デジタル変換部３０２は、図３に示すＡ／Ｄ変換ボード１２０上のフィルタ１２１およびＡ／Ｄ変換器１２２に相当する。それ以外のデータ取得部３１０、スペクトル算出部３２０、行列要素算出部３３０、および音源探索部３４０は、音源探索用コンピュータ１００のＣＰＵ１１１（図１、図２参照）で実行される音源探索プログラム２００によって音源探索用コンピュータ１００内に実現される機能である。図４に示す音源探索プログラム２００および図５に示す音源探索方法との関係では、図６の音源探索装置３００のデータ取得部３１０は、図４の音源探索プログラム２００のデータ取得部２１０および図５の音源探索方法のデータ取得ステップ（ステップａ）に対応し、図６の音源探索装置３００のスペクトル算出部３２０は、図４の音源探索プログラム２００のスペクトル算出部２２０および図５の音源探索方法のスペクトル算出ステップ（ステップｂ）に対応し、図７の音源探索装置３００の行列要素算出部３３０は、図４の音源探索プログラム２００の行列要素算出部２３０および図５の音源探索方法の行列要素算出ステップ（ステップｃ）に対応し、図６の音源探索装置３００の音源探索部３４０は、図４の音源探索プログラム２００の音源探索部２４０および図５の音源探索方法の音源探索ステップ（ステップｄ）に対応する。ただし、図６の音源探索装置３００のデータ取得部３１０、スペクトル算出部３２０、行列要素算出部３３０、および音源探索部３４０は、それぞれがハードウェアとソフトウェアとの複合からなるものであるのに対し、図４の音源探索プログラム２００の、それらに対応する要素は、アプリケーションプログラムのみからなるものである。また図５の音源探索方法は、本実施形態では、図４の音源探索プログラムが、図２、図３の音源探索用コンピュータ１００にインストールされてその音源探索用コンピュータ１００内で実行されることにより、その音源探索用コンピュータ１００内で実施される。
【０５４２】
以下、この図６の音源探索装置３００の各部の作用を説明することで、図４の音源探索プログラム２００および図５の音源探索方法の説明を兼ねることとする。
【０５４３】
図６の音源探索装置のデータ取得部３１０は、図３に示すＡ／Ｄ変換ボード１２０を制御して、そのＡ／Ｄ変換ボード１２０におけるフィルタリングおよびＡ／Ｄ変換により得られたデジタルの時系列データを、図４の音源探索プログラム２００内に取り込む役割りを成すものであり、ハードウェア上は、図４の音源探索プログラム２００を構成するデータ取得部２１０を実行するＣＰＵ１１１等がこれに相当する。
【０５４４】
このデータ取得部３１０では、上記のようにして、図１に示すマイクロホンアレイ２０を構成する多数のマイクロホン２１，２２，２３，…で生成された多数チャンネルの受音信号に対応する多数チャンネルの時系列データが取得される。
【０５４５】
次に、スペクトル算出部３２０では、ＦＦＴ（Ｆａｓｔ  Ｆｏｕｒｉｅｒ  Ｔｒａｎｓｆｏｒｍ）などが利用されて、各チャンネルの時系列データが離散フーリエ変換される。このとき、所望の角周波数ωのスペクトルが十分な高精度で得られるように工夫することが望ましい。所望の角周波数ωのスペクトルを十分な高精度で得るためには、アナログ・デジタル変換部３０２での標本化時間間隔（サンプリング時間間隔）を十分細かくとるなどの工夫が考えられるが、前掲の特許文献１で提案されているアルゴリズムを利用することが好ましい。ただし、本実施形態は、所望の角周波数ωのスペクトルを高精度に求めることの手法の如何を問うものではなく、この点についての更なる詳細説明は省略する。
【０５４６】
次の行列要素算出部３３０では、前述の式（４．１２）、（４．１３）、（４．１４ａ）に、スペクトル算出部３２０で算出された多数チャンネルの角周波数ωのスペクトルを代入することにより、行列要素
【０５４７】
【数２６５】
2008286710.tif 000267
【０５４８】
【数２６６】
2008286710.tif 000268
【０５４９】
【数２６７】
2008286710.tif 000269
【０５５０】
が算出される。
【０５５１】
ここで、本実施形態では、行列要素
【０５５２】
【数２６８】
2008286710.tif 000270
【０５５３】
を算出するにあたっては、式（４．１２）中の
【０５５４】
【数２６９】
2008286710.tif 000271
【０５５５】
として、前述の式（６．１）が採用され、その式（６．１）を用いて式（４．１２）を変形した式（６．２）に、スペクトル算出部３２０で算出された多数チャンネルの角周波数ωのスペクトルが代入されて、行列要素
【０５５６】
【数２７０】
2008286710.tif 000272
【０５５７】
が算出される。
【０５５８】
また、これと同様に、本実施形態では、行列要素
【０５５９】
【数２７１】
2008286710.tif 000273
【０５６０】
を算出するにあたっては、式（４．１３）中の
【０５６１】
【数２７２】
2008286710.tif 000274
【０５６２】
として前述の式（６．４）が採用され、その式（６．４）を用いて式（４．１３）を変形した式（６．５）に、スペクトル算出部３２０で算出された多数チャンネルの角周波数ωのスペクトルが代入されて、行列要素
【０５６３】
【数２７３】
2008286710.tif 000275
【０５６４】
が算出される。
【０５６５】
さらに、これと同様に、本実施形態では、行列要素
【０５６６】
【数２７４】
2008286710.tif 000276
【０５６７】
を算出するにあたっては、行列要素
【０５６８】
【数２７５】
2008286710.tif 000277
【０５６９】
の算出と同様、式（４．１４ａ）中の
【０５７０】
【数２７６】
2008286710.tif 000278
【０５７１】
として前述の式（６．１）が採用され、その式（６．１）を用いて式（４．１４ａ）を変形した式（６．３）に、スペクトル算出部３２０で算出された多チャンネルの角周波数ωのスペクトルが代入されて、行列要素
【０５７２】
【数２７７】
2008286710.tif 000279
【０５７３】
が算出される。
【０５７４】
さらに、図６に示す音源探索装置３００の音源探索部３４０では、行列要素算出部３３０で算出された行列要素
【０５７５】
【数２７８】
2008286710.tif 000280
【０５７６】
【数２７９】
2008286710.tif 000281
【０５７７】
【数２８０】
2008286710.tif 000282
【０５７８】
が連立方程式（５．１）に代入されその連立方程式（５．１）を解くことにより
【０５７９】
【数２８１】
2008286710.tif 000283
【０５８０】
が算出される。
【０５８１】
本実施形態では、観測点の離散化（図１に示すマイクロホンアレイ２０に並ぶマイクロホン２１，２２，２３，…の間隔）と計算点の離散化が同じ密度で行なわれるため、この場合、式（５．１）は式（５．２）のように変形される。本実施形態では、行列要素算出部３３０で算出された行列要素
【０５８２】
【数２８２】
2008286710.tif 000284
【０５８３】
【数２８３】
2008286710.tif 000285
【０５８４】
【数２８４】
2008286710.tif 000286
【０５８５】
が連立方程式（５．２）に代入されその連立方程式（５．２）が解かれ
【０５８６】
【数２８５】
2008286710.tif 000287
【０５８７】
が算出される。
【０５８８】
さらに具体的には、本実施形態では、
【０５８９】
【数２８６】
2008286710.tif 000288
【０５９０】
として
【０５９１】
【数２８７】
2008286710.tif 000289
【０５９２】
が採用されて、
【０５９３】
【数２８８】
2008286710.tif 000290
【０５９４】
でのｘ−ｙ平面上の、角周波数ωのｚ方向の粒子速度分布が算出される。
【０５９５】
図７は、従来から広く行なわれている近接場音響ホログラフィーの解析手法であるフーリエ変換法による計算結果を示す図であり、本発明に対する比較例に相当する。
【０５９６】
また、図８は、本実施形態により求められた粒子速度分布を示す図であり、本発明の実施例に相当する。
【０５９７】
これら図７と図８を比較すると、明らかに分かるように、図８に示す実施例の方が音源位置を高精度に特定している。
【図面の簡単な説明】
【０５９８】
【図１】音源とマイクロホンアレイのレイアウトを示す図である。
【図２】本発明の音源探索装置の一実施形態として動作する音源探索用コンピュータの外観斜視図である。
【図３】図２に示す音源探索用コンピュータ１００のハードウェア構成図である。
【図４】本発明の一実施形態としてのピーク位置探索プログラムの構成を示す図である。
【図５】本発明の一実施形態としての音源探索方法の構成を示すフローチャートである。
【図６】本発明の音源探索装置の一実施形態の機能ブロック図である。
【図７】従来から広く行なわれている近接場音響ホログラフィーの解析手法であるフーリエ変換法による計算結果を示す図である。
【図８】本実施形態により求められた粒子速度分布を示す図である。
【符号の説明】
【０５９９】
１１，１２    音源センサ
２０    マイクロホンアレイ
２１，２２，２３，…    マイクロホン
１００    音源探索用コンピュータ
１０１    本体
１０２    ディスプレイ
１０３    キーボード
１０４    マウス
１０５    ＣＤ−ＲＯＭ
１１０    バス
１１１    ＣＰＵ
１２０    Ａ／Ｄ変換ボード
１２１    フィルタ
１２２    Ａ／Ｄ変換器
２００    音源探索プログラム
２１０    データ取得部
２２０    スペクトル算出部
２３０    行列要素算出部
２４０    音源探索部
３００    音源探索装置
３０１    アンチエイリアジングフィルタ
３０２    アナログ・デジタル変換部
３１０    データ取得部
３２０    スペクトル算出部
３３０    行列要素算出部
３４０    音源探索部
【図１】
2008286710.tif 000420
【図２】
2008286710.tif 000421
【図３】
2008286710.tif 000422
【図４】
2008286710.tif 000423
【図５】
2008286710.tif 000424
【図６】
2008286710.tif 000425
【図７】
2008286710.tif 000426
【図８】
2008286710.tif 000427