"複数のスピーカを用いてエバネッセント波を再生する局所再生装置のフィルタ係数を決定するフィルタ係数決定装置であって、前記スピーカは円筒面上に配置されるものとし、前記エバネッセント波の再生方向の波数情報を算出する波数算出部と、前記スピーカの配置情報と前記波数情報とを用いて前記フィルタ係数を算出するフィルタ係数算出部と、を含む、フィルタ係数決定装置。請求項1記載のフィルタ係数決定装置であって、前記エバネッセント波の再生方向を円筒面の軸方向とし、kmaxを入力信号が含む波数の最大値とし、α1を1より大きい所定の値とし、前記波数算出部は、波数情報kz=α1kmaxを算出し、前記フィルタ係数算出部は円筒面の軸方向に配置された複数のスピーカの配置情報と前記波数情報kzとを用いて前記フィルタ係数を算出する、フィルタ係数決定装置。請求項1記載のフィルタ係数決定装置であって、前記エバネッセント波の再生方向を円筒面の軸方向とし、kを入力信号の波数とし、β1を1より大きい所定の値とし、前記波数算出部は、波数情報kzk=β1kを算出し、前記フィルタ係数算出部は円筒面の軸方向に配置された複数のスピーカの配置情報と前記波数情報kzkとを用いて前記フィルタ係数を算出する、フィルタ係数決定装置。請求項1記載のフィルタ係数決定装置であって、前記エバネッセント波の再生方向を円筒面の軸方向とし、kaを1より大きい所定の値とし、kを入力信号の波数とし、前記波数算出部は、波数情報kzk=±√ka2k2を算出し、前記フィルタ係数算出部は円筒面の軸方向に配置された複数のスピーカの配置情報と前記波数情報kzkとを用いて前記フィルタ係数を算出する、フィルタ係数決定装置。請求項1記載のフィルタ係数決定装置であって、前記エバネッセント波の再生方向を円筒面の周方向とし、α2を1より大きい所定の値とし、kmaxを入力信号が含む波数の最大値とし、aを前記円筒面の成す円の半径とし、前記波数算出部は、波数情報n=α2kmaxaを算出し、前記フィルタ係数算出部は円筒面の周方向に配置された複数のスピーカの配置情報と前記波数情報nとを用いて前記フィルタ係数を算出する、フィルタ係数決定装置。請求項1記載のフィルタ係数決定装置であって、前記エバネッセント波の再生方向を円筒面の周方向とし、β2を1より大きい所定の値とし、kを入力信号の波数とし、aを前記円筒面の成す円の半径とし、前記波数算出部は、波数情報nk=β2kaを算出し、前記フィルタ係数算出部は円筒面の周方向に配置された複数のスピーカの配置情報と前記波数情報nkとを用いて前記フィルタ係数を算出する、フィルタ係数決定装置。請求項1から6の何れかに記載のフィルタ係数決定装置であって、a,φm,zlを前記スピーカの配置情報とし、kを入力信号の波数とし、kzを軸方向の前記波数情報とし、nを周方向の前記波数情報とし、前記フィルタ係数算出部は、フィルタ係数をとして算出する、フィルタ係数決定装置。請求項1から請求項7の何れかに記載のフィルタ係数決定装置で決定されたフィルタ係数を複数のスピーカに設定し該スピーカを用いてエバネッセント波を再生する局所再生装置であって、前記フィルタ係数と前記入力信号とを乗じて前記スピーカを駆動するための駆動信号を生成するフィルタ乗算部と、をさらに含む、局所再生装置。複数のスピーカを用いてエバネッセント波を再生する局所再生方法のフィルタ係数を決定するフィルタ係数決定方法であって、前記スピーカは円筒面上に配置されるものとし、前記エバネッセント波の再生方向の波数情報を算出する波数算出ステップと、前記スピーカの配置情報と前記波数情報とを用いて前記フィルタ係数を算出するフィルタ係数算出ステップと、を含む、フィルタ係数決定方法。コンピュータを請求項1から請求項7の何れかに記載のフィルタ係数決定装置または請求項8に記載の局所再生装置として機能させるためのプログラム。しかしながら、従来技術は平面スピーカアレーを用いているため、その受聴領域91は平面スピーカアレーの前面の平板状に限定される。言い換えると、平板状の厚みはフィルタの設計方法により変更することができるが、その形状は平板状に限定される。本発明では、円筒面上に配置された複数のスピーカを用いて、その円筒面から遠ざかる方向円筒面の円の中心から拡がる方向なので以下放射方向ともいうに急峻に減衰するエバネッセント波を再生し、受聴領域を円筒状に形成する局所再生装置、及びその装置で用いるフィルタの係数を決定する技術を提供することを目的とする。上記の課題を解決するために、本発明の第一の態様によれば、円筒面上に配置される複数のスピーカを用いてエバネッセント波を再生する局所再生装置のフィルタ係数を決定する。エバネッセント波の再生方向の波数情報を算出し、スピーカの配置情報と波数情報とを用いてフィルタ係数を算出する。本発明に係る局所再生装置、及びその装置に用いるフィルタの係数を決定する技術であれば、従来とは異なる円筒状の受聴領域を形成することができる。以下、本発明の実施形態について説明する。なお、以下の説明に用いる図面では、同じ機能を持つ構成部や同じ処理を行うステップには同一の符号を記し、重複説明を省略する。<エバネッセント波の理論>以下、本実施形態において利用するエバネッセント波の理論について説明する。1円筒座標系における波動方程式の一般解音波はある法則に従って伝搬し、その音波の伝搬の秩序を表す方程式は波動方程式と呼ばれ、式1で表される。xは音波の進行方向における位置を、∂は偏微分を、tは時刻を、cは音速を、pは音圧を表す。図2の円筒座標系r,φ,z但し、rは半径を、φは角度を、zは高さを表すにおける波動方程式の一般解周波数領域は、式2で与えられる。但し、jは虚数単位、ωは周波数、kr、kzはそれぞれr、z方向なお、r方向は前述の放射方向である。以下、z方向を軸方向ともいうの波数、H1nkrrはn次の第1種ハンケル関数、H2nkrrはn次の第2種ハンケル関数、Ankz,ω、Bnkz,ωは任意の定数を表す。上記の任意の定数Ankz,ω、Bnkz,ωを決定するために、ある境界面r=aにおける境界条件を規定する必要がある。図3のように境界面94r=aが全ての波源93を包含する場合を考えると但し、図中、95は音圧の解析領域を示す、式2の一般解は式3のように書き換えられる。このとき、定数Bnkz,ωの値に関わらず音圧を表すことができる。境界面94r=a上の測定音圧は式3を満たすので、以下のように表すことができる。ここで円筒座標系における空間スペクトルPnr,kzと測定音圧Pr,φ,z,ωの関係を考えると、以下のように表すことができる。従って、境界条件境界面r=aが全ての波源を包含するとの条件を満たす場合には、式4の右辺と式Bにr=aを代入したときの右辺から以下のように表すことができる。これを式3に代入すると、以下のようになる。2軸方向におけるエバネッセント波発生条件軸方向zに沿った波長λzλz=2π/kzが音源の波長λλ=2π/k、但しkは波源の波数より短い場合、放射方向にエバネッセント波が発生する図4参照、なお、図中、96は受聴領域を示す。このときkz>kである。また波源の波数kとkr、kzの間には、k2rk2k2zという関係が成り立つ。kz>kの条件下では、となる。krが純虚数となる場合、式6中のハンケル関数H1nは変形ベッセル関数Knk’rrに書き換えられる。ハンケル関数と変形ベッセル関数の関係式は次式で表される。ここで変形ベッセル関数の変数k’rrが大きいと仮定すると、以下の近似が成り立つ。従って式6中のハンケル関数の比は以下のように書き換えられる。以上より、軸方向における発生条件を満たす場合のエバネッセント波の式は式12で与えられ、kz>kにおいて成り立つ。後述の波数算出部では、エバネッセント波の再生方向が軸方向と入力された場合に、上記の条件kz>kを満たすような波数kzを出力する。なお、エバネッセント波の再生方向とは、円筒面上に配置されたスピーカの内、再生されるスピーカの配置された方向を意味する。例えば、エバネッセント波の再生方向は軸方向とする場合には、図4の斜線部分97に位置するスピーカが再生される。3周方向におけるエバネッセント波発生条件ここではkz≦kとし、軸方向におけるエバネッセント波発生条件を満たさない場合を考える。周方向φに沿った波長λφλφ=2πa/nが音源の波長λλ=2π/kより短い場合も円筒面から遠ざかる方向つまり放射方向にエバネッセント波が発生する図5参照、なお、図中、98は受聴領域を示す。ここで2πaは円周の長さであり、nは周方向の波数である。このとき、2πa/n<λであるので、nが十分に大きい場合を考える。n→∞の場合のハンケル関数の漸近展開は、となる。ここでζ=krr=krであるkz=0と仮定。ζ/n<1のとき式の実部は無視することができ、このときハンケル関数はkrnに比例して減少する。従って式6中のハンケル関数の比は以下のように書き換えられる。以上より、周方向における発生条件を満たす場合のエバネッセント波の式は、式15で与えられ、kr<nにおいて成り立つ。後述の波数算出部180では、エバネッセント波の再生方向が周方向と入力された場合に、上記の条件kr<nを満たすような波数nを出力する。例えば、図5の斜線部分99に位置するスピーカが再生される。円筒座標系における波動方程式の一般解において放射方向に急峻に減衰するエバネッセント波の式に着目し、円筒面上に並べた複数のスピーカスピーカアレーによってエバネッセント波を再生する装置を構成した。エバネッセント波発生条件は上述の2種類が存在し、それぞれの条件に則ったフィルタ係数決定方法を以下の実施形態において示す。<第一実施形態に係る局所再生装置100>図6、7、8に局所再生装置100の構成を示す。図6は局所再生装置100の外観図であり、図7は局所再生装置100を横から見た概念図であり、図8は局所再生装置100を上から見た概念図である。局所再生装置100は、M×L個のスピーカ110lm但し、m=1,...,M、lはローマ字のエルでありl=1,...,Lとする。Mは周方向のスピーカ総数であり、Lは軸方向のスピーカ総数である。と、M×L個のA/D変換器160lmと、M×L個のフィルタ乗算部120lmと、M×L個のD/A変換器150lmと、M×L個のアンプ140lmと、フィルタ係数算出部170と、波数算出部180とを含む図7参照。なお、フィルタ乗算部120lmは、ディジタルフィルタによって構成される。スピーカ110lmは円筒面の周方向と軸方向に離散的に配置される図6参照。なお、スピーカ間の間隔は一定でなくてもよいが、各方向に対して一定であるほうが好ましい。図6のスピーカ110lmの位置をa,φm,zlとする。aは半径言い換えると、円筒面の軸からの距離を、φmは角度を、zlは高さを表し、aは全てのスピーカ110lmに共通であり、φm及びzlは各スピーカにより異なる値となる。製造者または利用者によって、これらのスピーカ110lmの位置以下配置情報ともいうa,φm,zlは図示しない記憶部に予め記憶されているものとする。局所再生装置100は信号入力端子30を介して、アナログ入力信号Skを受け取る図9参照、s1、但し、kは波源の波数を表す。A/D変換器160lmはアナログ入力信号Skを受け取り、ディジタル入力信号Skに変換し、フィルタ乗算部120lmに出力する。一方、波数算出部180は、再生に先立ち製造者や利用者等によって予め設定されたエバネッセント波の再生方向と波数決定パターンに従って、波数周方向の波数nまたは軸方向の波数kzの何れかの波数を算出しs2、フィルタ係数算出部170に出力する。フィルタ係数算出部170は、スピーカ110lmの位置aφm,zlと波数kとを用いて、フィルタ係数Haφm,zl,kを算出しs3、フィルタ乗算部120lmに設定する。フィルタ乗算部120lmは、ディジタル入力信号Skとフィルタ係数Haφm,zl,kとを乗じて駆動信号Da,φm,zl,kを生成しs4、D/A変換器150lmに出力する。D/A変換器150lmは、ディジタル駆動信号Daφm,zl,kをアナログ駆動信号Daφm,zl,kに変換し、アンプ140lmに出力する。アンプ140lmはアナログ駆動信号Daφm,zl,kを増幅し、スピーカ110lmに供給する。M×L個のスピーカ110lmは、入力信号Skに応じたアナログ駆動信号Daφm,zl,kを再生しs5、音響信号を出力する。以下にフィルタ乗算部120lm、フィルタ係数算出部170、波数算出部180の詳細を説明する。<フィルタ乗算部120lm>フィルタ乗算部120lmでは、フィルタ係数算出部170から出力されるフィルタ係数Ha,φm,zl,kとAD変換器160lmにより変換されたディジタル入力信号Skとを入力として、スピーカ110lmを駆動するための信号である駆動信号Da,φm,zl,kを生成し、出力する。フィルタ乗算部120lmには、再生に先立ち後述のフィルタ係数算出部170で計算しておいたフィルタ係数Ha,φm,zl,kを設定しておく。フィルタ乗算部120lmにディジタル入力信号Skが入力されると、フィルタ乗算部120lmは入力されたディジタル入力信号Skにフィルタ係数Ha,φm,zl,kを畳み込むことでエバネッセント波を発生するためのディジタル駆動信号Da,φm,zl,kを生成し、この駆動信号Da,φm,zl,kを対応するD/A変換器150lmに出力する。なお、駆動信号Da,φm,zl,kは以下の式であらわされる。フィルタ係数の算出方法についてはフィルタ係数算出部170で説明する。<フィルタ係数算出部170>フィルタ係数算出部170では、スピーカ110lmの配置情報a,φm,zlとエバネッセント波を再生するための波数以下波数情報ともいうkzまたはnの一方を入力とし、フィルタ乗算部120lmに与えるフィルタ係数Ha,φm,zl,kを算出し、フィルタ乗算部120lmへ出力する。フィルタ係数算出部170は軸方向の波数kz、周方向の波数nとスピーカ位置a,φm,zlを用いて、により、フィルタ係数Ha,φm,zl,kを算出する。式17は周波数領域におけるフィルタ係数Ha,φm,zl,kである。なお、入力信号の波数kと周波数fの関係は、音速をcとしてと表す。A/D変換器160lmにおける入力信号のサンプリング周波数をfs、周波数の帯域分割数をQとし、q=0,1,,Q-1とした場合に、各周波数帯域をと表す。式Cより、kは各周波数帯域fqにより異なる値となる。式17で用いる波数情報はkzまたはnは、後述する波数算出部180において波数kとエバネッセント波の再生方向、波数決定パターンを用いて算出されるので、フィルタ係数算出部170は各周波数帯域fqにおけるフィルタ係数を算出することができる。フィルタ係数算出部170は、周波数帯域fqにおけるフィルタ係数を逆フーリエ変換したものをフィルタ乗算部120lmにおけるフィルタ係数として設定する。つまり、Q個の周波数領域のフィルタ係数を算出し、それを逆フーリエ変換により時間領域のフィルタ係数に変換し、これをフィルタ乗算部120lmに設定する。なお、周波数の帯域分割数Qの値を大きくとればフィルタとしての精度は高くなるが計算コストが大きくなる。Qは任意に設定することができ、例えば数百程度の値とする。また2のべき乗数に設定すると演算が高速化できる。<波数算出部180>波数算出部180はエバネッセント波の再生方向軸方向または周方向と波数決定パターン波数決定パターンの種類については後述するを入力として、エバネッセント波の再生方向の波数kz、nを算出し、フィルタ係数算出部170へ出力する。以下では各入力条件における波数算出方法を説明する。1周方向にエバネッセント波を発生させる場合ここではエバネッセント波の再生方向が軸方向と入力された場合の波数算出方法を示す。軸方向にエバネッセント波を発生させるため、周方向の波数nはn=0と決定される。なお、これは周方向のスピーカが同相で駆動することを表すため、装置構成における周方向のスピーカ個数Mは1個でもよい。M=1の場合は円筒面状スピーカアレーではなく、軸方向の直線状スピーカアレーとして実現される。ここでは軸方向にエバネッセント波が発生する条件であるkz>kを満たすような波数kzを決定する。決定された波数を用いることで、式12で表されるエバネッセント波が再生される。これは放射方向に対して、指数関数的に音圧が減衰する特性をもつ。従って単一スピーカより急峻に減衰する特徴をもつが、波数の設定次第では可聴範囲が狭まりすぎてしまうという欠点もある。以下では3種類の波数kzの決定方法波数決定パターンを示す。i波数kzを音源の波数kの最大値kmaxを用いて決定する方法入力信号が含む波数kの最大値をkmaxとすると、kz=α1kmaxとして指定する。但し、α1は1より大きい任意の定数である。このようにして決定したkzは音源の波数kに依存しない。言い換えると、式C’に示す最大値kmax=kQ-1以外の波数k0,k1,...,kQ-2に依存しない。ここでのαの値は1~数十程度が目安であり、値を変化させると減衰しない周波数の位置が変わる傾向があるため、入力信号の周波数特性の応じた減衰特性の制御が可能である。なお、式C’及び式18より、kmax=2πfmax/cfmax=fQ-1である。ii波数kzを音源の波数kの定数倍として決定する方法β1を1より大きい任意の定数として、kzk=β1kとして指定する。ここでのβ1の値は1~数十程度が目安である。なお、式C’よりkは各周波数帯域fqに応じて変化する値なので、kzも各周波数帯域fq及びkに応じて変化する値となる。iii波数kzを減衰制御項から決定する方法エバネッセント波の式12において、Pr,φ,zのr方向の指数減衰項がekar但し、添え字kaはkaを表すとなるとすると、式20は、任意の減衰量と音源の波数から波数kzを求めることができることを示している。指定するkaの値については、あまりに大きいと装置近傍ですぐに減衰してしまう制御となるため、ka=1~数十として設計するのが目安である。なお、前述の通り、式C’よりkは各周波数帯域fqに応じて変化する値なので、iiの場合と同様に、kzも各周波数帯域fq及びkに応じて変化する値となる。2周方向にエバネッセント波を発生させる場合ここではエバネッセント波の再生方向が周方向と入力された場合の波数算出方法を示す。周方向に発生させるため、軸方向の波数kzはkz=0と決定される。なお、これは軸方向のスピーカが同相で駆動することを表すため、装置構成における軸方向のスピーカ個数Lは1個でもよい。但しL=1の場合は円筒面状スピーカアレーではなく、円状スピーカアレーとして実現される。周方向にエバネッセント波が発生する条件であるn>kaを満たすような波数nを決定する。決定された波数nを用いることで、式15で表されるエバネッセント波が再生される。これは放射方向に対して、べき乗関数に従って音圧が減衰する特性をもつ。従って軸方向に再生する場合と比較して減衰特性は緩やかになるが、波数を調整することで減衰特性を制御することが可能であり、可聴範囲を制御するのが容易であるといえる。以下では2種類の波数nの決定方法波数決定パターンを示す。i波数nを音源の波数の最大値を用いて決定する方法入力信号が含む波数kの最大値をkmaxとすると、n=α2kmaxaとして指定する。但し、α2は1より大きい任意の定数である。このようにして決定したnは音源の波数kに依存しない言い換えると最大値kmax以外の波数kに依存しない。ここでのα2の値は1~数十程度が目安であり、値を変化させると減衰しない周波数の位置が変わる傾向があるため、入力信号の周波数特性の応じた減衰特性の制御が可能である。なお、kmaxについては、1iで説明した値と同様である。ii波数nを音源の波数の定数倍として決定する方法β2を1より大きい任意の定数として、nk=β2kaとして指定する。ここでのβ2の値は1~数十程度が目安である。kについては1iiで説明した値と同様である。<効果>このような構成とすることで、円筒面スピーカアレーを用いてエバネッセント波を再生し、その受聴領域を従来とは異なる円筒状とし、減衰効果範囲を拡大できる、さらに指数減衰だけでなくべき乗関数に従うような減衰特性を有するエバネッセント波を再生することができる。<その他の変形例>第一実施形態では、局所再生装置100内部にフィルタ係数算出部170と波数算出部180を含む構成としたが、フィルタ係数算出部170と波数算出部180を含むフィルタ係数決定装置190を別装置として設けてもよい。フィルタ係数決定装置190は、局所再生装置100の内部に組込まれてもよいし、別装置として外部に設けてもよい。フィルタ係数決定装置190は、発生方向と波数決定パターンとスピーカの位置情報を用いて、フィルタ係数を決定し、局所再生装置100内のフィルタ乗算部120lmに出力する構成とする。各処理内容は第一実施形態において説明した通りである。本発明は上記の実施形態及び変形例に限定されるものではない。例えば、上述の各種の処理は、記載に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。その他、本発明の趣旨を逸脱しない範囲で適宜変更が可能である。<プログラム及び記録媒体>上述した局所再生装置またはフィルタ係数決定装置は、コンピュータにより機能させることもできる。この場合はコンピュータに、目的とする装置各種実施例で図に示した機能構成をもつ装置として機能させるためのプログラム、またはその処理手順各実施例で示したものの各過程をコンピュータに実行させるためのプログラムを、CD-ROM、磁気ディスク、半導体記憶装置などの記録媒体から、あるいは通信回線を介してそのコンピュータ内にダウンロードし、そのプログラムを実行させればよい。"