JP7614872B2

JP7614872B2 - Sound field reproduction device and program

Info

Publication number: JP7614872B2
Application number: JP2021019629A
Authority: JP
Inventors: 陽佐々木
Original assignee: Japan Broadcasting Corp
Current assignee: Japan Broadcasting Corp
Priority date: 2021-02-10
Filing date: 2021-02-10
Publication date: 2025-01-16
Anticipated expiration: 2041-02-10
Also published as: JP2022122414A

Description

本発明は、スピーカアレイを用いて、移動する音源によって形成される音場を再現する音場再現装置及びプログラムに関する。 The present invention relates to a sound field reproduction device and program that uses a speaker array to reproduce a sound field formed by a moving sound source.

従来、多数のスピーカを用いて、ある空間に任意の音場を形成する音場再現技術に関する研究が進められている。音場再現の主要な技術としては、波面合成法（ＷＦＳ：Wave Field Synthesis）、境界音場制御（ＢｏＳＣ：Boundary Surface Control）、高次アンビソニックス（ＨＯＡ：Higher Order Ambisonics）、スペクトル除算法（ＳＤＭ：Spectral Division Method）等が知られている。 Research has been ongoing into sound field reproduction technology that uses multiple speakers to create an arbitrary sound field in a certain space. Known sound field reproduction technologies include Wave Field Synthesis (WFS), Boundary Surface Control (BoSC), Higher Order Ambisonics (HOA), and Spectral Division Method (SDM).

図８は、波面合成法（ＷＦＳ）におけるスピーカアレイの配置例を示す図である。ＷＦＳは、レイリー積分に基づいた手法であり、直線状または平面状に配置されたスピーカアレイを用いて、所望音場の境界平面上の音圧または音圧勾配を再現することで、境界外部から到来する音源による波面を形成するものである。 Figure 8 shows an example of speaker array placement in wave field synthesis (WFS). WFS is a method based on the Rayleigh integral, and uses a speaker array placed in a line or plane to reproduce the sound pressure or sound pressure gradient on the boundary plane of the desired sound field, forming a wavefront due to a sound source coming from outside the boundary.

図９は、境界音場制御（ＢｏＳＣ）におけるスピーカアレイの配置例を示す図である。ＢｏＳＣは、キルヒホッフ・ヘルムホルツ積分方程式に基づいた手法である。ＢｏＳＣは、音場を再現したい領域の外部に、領域内部に向けて聴取者を囲むようにスピーカアレイを配置し、逆システムを用いてその領域の境界面上の音圧及び音圧勾配を制御することで、領域内部に所望音場を再現するものである。 Figure 9 shows an example of speaker array placement in boundary sound field control (BoSC). BoSC is a method based on the Kirchhoff-Helmholtz integral equation. BoSC places speaker arrays outside the area where you want to reproduce the sound field, surrounding the listener toward the inside of the area, and uses an inverse system to control the sound pressure and sound pressure gradient on the boundary surface of the area, thereby reproducing the desired sound field inside the area.

図１０は、高次アンビソニックス（ＨＯＡ）におけるスピーカアレイの配置例を示す図である。ＨＯＡは、球面調和関数を用いて所望音場を表現し、球面内部に向けて聴取者を囲むように配置した球面スピーカアレイを用いて、再現音場の球面調和係数を所望音場の球面調和係数とマッチングさせることで、音場を再現する手法である。 Figure 10 shows an example of speaker array placement in Higher Order Ambisonics (HOA). HOA is a method of reproducing a sound field by expressing a desired sound field using spherical harmonic functions, using a spherical speaker array arranged to surround the listener toward the inside of the sphere, and matching the spherical harmonic coefficients of the reproduced sound field with the spherical harmonic coefficients of the desired sound field.

図１１は、スペクトル除算法（ＳＤＭ）におけるスピーカアレイの配置例を示す図である。ＳＤＭは、角度スペクトルを用いて所望音場を表現し、直線状または平面状に配置されたスピーカアレイを用いて、再現音場の角度スペクトルを所望音場の角度スペクトルとマッチングさせることで、音場を再現する手法である。ＳＤＭの詳細については、例えば非特許文献１を参照されたい。 Figure 11 shows an example of speaker array placement in the Spectral Division Method (SDM). SDM is a technique for expressing a desired sound field using an angular spectrum, and reproducing the sound field by matching the angular spectrum of the reproduced sound field with the angular spectrum of the desired sound field using a speaker array arranged in a line or a plane. For details on SDM, see, for example, Non-Patent Document 1.

〔ＳＤＭにおける駆動信号の算出方法〕
以下、一般的なＳＤＭにおける駆動信号の算出方法について説明する。図１２は、一般的なＳＤＭを説明する図である。 [Method of calculating drive signal in SDM]
A method for calculating a drive signal in a typical SDM will be described below with reference to FIG.

ｘｙｚ空間上のｙ＝ｙ₀，ｚ＝ｚ₀に配置された無限直線音源１０１上の座標ｒ0＝（ｘ₀, ｙ₀，ｚ₀）の点において、角周波数ωの駆動信号をＤ（ｒ0，ω）とする。この場合の座標ｒ＝（ｘ，ｙ_ref，ｚ₀）の点における音圧Ｐ（ｒ，ω）は、以下の式にて表される。

A driving signal with an angular frequency ω is defined as D(r0,ω) at a point of coordinate r0=( _x0 , _y0 , _z0 ) on an infinite linear sound source 101 placed at y= _y0 , z= _z0 in the xyz space. In this case, the sound pressure P(r,ω) at the point of coordinate r=(x, _yref , _z0 ) is expressed by the following formula.

ここで、Ｇ（ｒ－ｒ0，ω）は、ベクトル（ｒ－ｒ0）で表現される方向及び距離の伝達関数であり、再現音場として自由音場を想定した場合、３次元の自由音場グリーン関数で表現され、以下の式にて表される。

ｋは波数である。 Here, G(r-r0, ω) is a transfer function of direction and distance expressed by the vector (r-r0). When a free sound field is assumed as the reproduced sound field, it is expressed by a three-dimensional free sound field Green's function and is given by the following formula.

k is the wave number.

音圧Ｐ（ｒ，ω）は、駆動信号Ｄ（ｒ0，ω）及び伝達関数Ｇ（ｒ－ｒ0，ω）の空間上の畳み込みによって得られるものと解釈でき、ｙ＝ｙ_refにおいてｘ軸の方向に沿ってフーリエ変換することにより、以下の角度スペクトル表現Ｐ^（ｋ_x，ｙ_ref，ｚ₀，ω）が得られる。

The sound pressure P(r,ω) can be interpreted as being obtained by spatially convoluting the drive signal D(r0,ω) and the transfer function G(r-r0,ω), and by Fourier transforming along the x-axis direction at y = y _ref , the following angular spectrum representation P^(k _x , y _ref , z ₀ ,ω) is obtained.

Ｐ^（ｋ_x，ｙ_ref，ｚ₀，ω）は、音圧Ｐ（ｒ，ω）をｙ＝ｙ_refで空間フーリエ変換することにより得られる角度スペクトルである。また、Ｄ^（ｋ_x，ｙ₀，ｚ₀，ω）は、駆動信号Ｄ（ｒ0，ω）をｙ＝ｙ₀で空間フーリエ変換することにより得られる角度スペクトルである。さらに、Ｇ^（ｋ_x，ｙ_ref－ｙ₀，ｚ₀，ω）は、伝達関数Ｇ（ｒ－ｒ0，ω）をｙ＝ｙ_refで空間フーリエ変換することにより得られる角度スペクトルである。ｋ_xはｘ軸方向の波数である。 P^( _kx , _yref , _z0 , ω) is the angular spectrum obtained by spatial Fourier transforming the sound pressure P(r, ω) with y = _yref . Furthermore, D^( _kx , _y0 , _z0 , ω) is the angular spectrum obtained by spatial Fourier transforming the drive signal D(r0, ω) with y = _y0 . Furthermore, G^( _kx , _yref - _y0 , _z0 , ω) is the angular spectrum obtained by spatial Fourier transforming the transfer function G(r - r0, ω) with y = _yref . _kx is the wave number in the x-axis direction.

ここで、波面を再現したい所望音場のｙ＝ｙ_refにおける音圧分布を空間フーリエ変換することにより得られる角度スペクトルをＰ^_d（ｋ_x，ｙ_ref，ｚ₀，ω）とすると、所望音場の波面を再現するための駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、以下の式にて表される。

Here, if the angular spectrum obtained by spatial Fourier transforming the sound pressure distribution at y = _yref of the desired sound field in which the wavefront is to be reproduced is denoted as P^ _d ( _kx , _yref , _z0 , ω), the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal for reproducing the wavefront of the desired sound field is expressed by the following equation.

座標ｒ_s＝（ｘ_s，ｙ_s，ｚ₀）に配置された点音源が周波数ωにおいてＳ（ω）で駆動された場合、所望音場のｙ＝ｙ_refでの角度スペクトルＰ^_d（ｋ_x，ｙ_ref，ｚ₀，ω）は、以下の式にて表される。

When a point sound source located at coordinates r _s = (x _s , y _s , z ₀ ) is driven by S(ω) at frequency ω, the angular spectrum P^ _d (k _x , y _ref , z ₀ , ω) of the desired sound field at y = y _ref is expressed by the following equation.

ここで、Ｈ₀ ⁽²⁾は第２種０次ハンケル関数であり、Ｋ₀は０次変形ベッセル関数である。また、０≦ｋ²－ｋ_x ²の場合の音波は、遠方まで伝搬する音波を示し、０＞ｋ²－ｋ_x ²の場合の音波はエバネッセント波と呼ばれ、振幅が急速に減衰する遠方まで伝搬しない音波を示す。 Here, H ₀ ⁽²⁾ is the 0th order Hankel function of the second kind, and K ₀ is the 0th order modified Bessel function. When 0≦k ² −k _x ² , the sound wave propagates far, and when 0>k ² −k _x ² , the sound wave is called an evanescent wave, and the amplitude attenuates rapidly and the sound wave does not propagate far.

同様に、角度スペクトルＧ^（ｋ_x，ｙ_ref－ｙ₀，ｚ₀，ω）は、以下の式にて表される。

ただし、０＞ｋ²－ｋ_x ²のときは、音源近傍のみに伝搬するエバネッセント波を表しており、エバネッセント波の再現は数値的に不安定になり易いため、一般には省かれることが多い。 Similarly, the angular spectrum G^(k _x , y _ref −y ₀ , z ₀ , ω) is expressed by the following equation.

However, when 0>k ² -k _x ² , it represents an evanescent wave that propagates only in the vicinity of the sound source, and since reproduction of the evanescent wave tends to be numerically unstable, it is generally omitted.

したがって、駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、以下の式にて表す。

以下、特に断りなく０＞ｋ²－ｋ_x ²の場合はＤ^（ｋ_x，ｙ₀，ｚ₀，ω）＝０とし、０≦ｋ²－ｋ_x ²の場合についてのみ説明するものとする。 Therefore, the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal is expressed by the following equation.

In the following, unless otherwise specified, when 0>k ² -k _x ^2, D^(k _x , y ₀ , z ₀ , ω)=0, and only the case when 0≦k ² -k _x ² will be explained.

図１３は、再現領域及び再現境界を説明する図である。直線状に配置されたスピーカアレイ１００が無限直線音源である場合、当該無限直線音源を、前記式（７）の駆動信号を用いて駆動することにより、ｙ≧ｙ_refの領域において音場を再現することができる。このため、ここではｙ＝ｙ_refを再現境界、ｙ≧ｙ_refの領域を再現領域と呼ぶこととする。 13 is a diagram for explaining the reproduction area and the reproduction boundary. When the speaker array 100 arranged in a line is an infinite linear sound source, the infinite linear sound source can be driven using the drive signal of the above formula (7) to reproduce a sound field in the area of y≧ _yref . For this reason, y= _yref is referred to as the reproduction boundary, and the area of y≧ _yref is referred to as the reproduction area.

前記式（７）の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、無限直線音源を用いて波面合成を行うことを想定しているが、実際には、直線状に配置した有限個のスピーカユニットからなるスピーカアレイ１００を用いて実現する。 The angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal in equation (7) is assumed to be produced by wave field synthesis using an infinite linear sound source, but in practice it is realized using a speaker array 100 consisting of a finite number of speaker units arranged in a line.

このため、前記式（７）の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）において、周波数、角度スペクトルを離散化すると共に、有限長で打ち切り、さらに離散逆フーリエ変換（ＩＤＦＴ：Inverse Discrete Fourier Transform）を行うことで、スピーカユニット毎の時間周波数領域の駆動信号を算出する必要がある。 For this reason, in the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal in equation (7) above, it is necessary to discretize the frequency and angular spectrum, truncate it to a finite length, and then perform an Inverse Discrete Fourier Transform (IDFT) to calculate the drive signal in the time-frequency domain for each speaker unit.

スピーカアレイ１００を構成するスピーカユニット数をＭ（自然数）、スピーカユニットの間隔をΔｘ、スピーカアレイ長をＬ（＝（Ｍ－１）Δｘ）とする。前記式（７）の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）に対し、離散逆フーリエ変換を行うことで、ｘに配置されたスピーカユニットの時間周波数領域の駆動信号Ｄ~（ｘ，ω）は、以下の式で表される。

ここで、ｍは、０～Ｍ－１（ｍは整数）であり、ｋ_x ^m＝２πｍ／（ΔｘＭ）は、ｘ軸方向の波数である。尚、スピーカアレイ１００は、ｘｙｚ空間において、予め設定されたｙ＝ｙ₀，ｚ＝ｚ₀のｘ軸上に配置される。 The number of speaker units constituting the speaker array 100 is M (a natural number), the spacing between the speaker units is Δx, and the length of the speaker array is L (=(M-1)Δx). By performing an inverse discrete Fourier transform on the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal in the time-frequency domain of the speaker unit arranged at x, the drive signal D~(x, ω) in the time-frequency domain of the speaker unit arranged at x is expressed by the following equation.

Here, m is an integer between 0 and M-1, and _kxm = ^2πm /(ΔxM) is the wave number in the x-axis direction. The speaker array 100 is disposed on the x-axis at a preset position of y = _y0 , z = _z0 in the xyz space.

そして、スピーカアレイ１００を構成するスピーカユニット毎に、時間領域の駆動信号が算出される。 Then, a time domain drive signal is calculated for each speaker unit that makes up the speaker array 100.

このようなＳＤＭ等の音場再現技術を用いることにより、ユーザが指定する座標に配置される点音源によって形成される音場を再現することができる。しかし、従来のスピーカ駆動信号算出方法によれば、静止した音源によって形成される音場を再現することは可能であるものの、音源が連続的に移動する状況は考慮されていない。 By using sound field reproduction technology such as SDM, it is possible to reproduce a sound field formed by a point sound source located at coordinates specified by the user. However, while conventional speaker drive signal calculation methods can reproduce a sound field formed by a stationary sound source, they do not take into account situations in which the sound source moves continuously.

計算機を用いて音声信号処理を行う場合、バッファリングされた信号に対して処理を行うことが一般的であり、静止した音源の座標をバッファ毎に移動させることで、断続的に移動する音源によって音場を再現することができる。しかし、この場合に再現される音場は、連続的に移動する音源によって形成される音場とは異なるため、所望音場とは音色が異なるという課題がある。 When processing audio signals using a computer, it is common to process buffered signals, and by moving the coordinates of a stationary sound source for each buffer, a sound field can be reproduced using an intermittently moving sound source. However, the reproduced sound field in this case differs from the sound field formed by a continuously moving sound source, and there is an issue that the timbre differs from the desired sound field.

この課題を解決するために、非特許文献２には、ＳＤＭをベースとし、移動する音源によって形成される音場の再現方法が提案されている。この音場再現方法では、時刻ｔ＝０のときに音源がｒ_s＝（ｘ_s，ｙ_s）に位置する。そして、時刻が－∞≦ｔ≦∞のときに、ｘ軸に対してφ方向に速度ｖで進行する音源により、ｙ＝ｙ_refで形成される音場を再現するための駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、以下の式にて表される。

ここで、ｋ^＝（ω－ω₀）／ｃ、ｋ^_t ²＝ｋ^²－（ω／ｃ）²である。ω₀は音源の角周波数であり、ｃは音速である。 To solve this problem, Non-Patent Document 2 proposes a method for reproducing a sound field formed by a moving sound source based on SDM. In this sound field reproduction method, the sound source is located at r _s = (x _s , y _s ) at time t = 0. Then, when the time is -∞≦t≦∞, the angular spectrum D^(k _x , y 0 , z 0 , ω) of the drive signal for reproducing a sound field formed at y = y _ref by a sound source traveling at a speed v in the _φ direction with respect to the _x axis is expressed by the following formula.

Here, k^ = (ω - ω ₀ )/c, and k^ _t ² = k^ ² - (ω/c) ^2. ω ₀ is the angular frequency of the sound source, and c is the speed of sound.

J.Ahrens, and S.Spors,“Sound field reproduction using planar and linear arrays of loudspeakers”, IEEE Trans. Audio, Speech, Lang. Process., vol.18, no.8, pp.2038-2050, Nov.2010.J. Ahrens, and S. Spors, “Sound field reproduction using planar and linear arrays of loudspeakers”, IEEE Trans. Audio, Speech, Lang. Process., vol.18, no.8, pp.2038-2050, Nov. 2010. Firtha Gergely and Fiala Peter,“Sound Field Synthesis of Uniformly Moving Virtual Monopoles”,JAES Volume 63 Issue 1/2 pp. 46-53, January 2015Firtha Gergely and Fiala Peter, “Sound Field Synthesis of Uniformly Moving Virtual Monopoles”, JAES Volume 63 Issue 1/2 pp. 46-53, January 2015

しかしながら、前述の移動する音源によって形成される音場の再現方法における前記式（９）の駆動信号の角度スペクトルでは、音源が無限遠（ｔ＝－∞）から無限遠（ｔ＝∞）まで等速直線運動が行われることが想定されている。 However, in the angular spectrum of the drive signal in equation (9) in the method for reproducing a sound field formed by a moving sound source described above, it is assumed that the sound source moves at a constant speed in a straight line from infinity (t = -∞) to infinity (t = ∞).

このため、前記式（９）の駆動信号の角度スペクトルにおいて、単純に空間方向に周波数を離散化すると共に、有限長で周波数の打ち切りを行うと、無限個の点音源が、異なる座標に時間エイリアスとして出現してしまう。 For this reason, if the frequency is simply discretized in the spatial direction and truncated to a finite length in the angular spectrum of the drive signal in equation (9), an infinite number of point sound sources will appear as time aliases at different coordinates.

本来は、１つの移動する音源によって形成される音場を再現するものであるが、同じ速度で移動する等間隔に整列した無限個の点音源によって形成される音場が再現されてしまう。 Originally, it was meant to reproduce the sound field formed by a single moving sound source, but it ended up reproducing the sound field formed by an infinite number of equally spaced point sound sources moving at the same speed.

図１４は、非特許文献２の音場再現方法により算出された所望音場の例を示す図である。図１４に示すように、前述の離散化及び打ち切り処理により、再現音場には、ｙ＝－２．５上をｘ＝－∞からｘ＝∞まで、同じ速度で移動する等間隔に整列した無限個の点音源が時間エイリアスとして出現してしまう。 Figure 14 shows an example of a desired sound field calculated by the sound field reproduction method of Non-Patent Document 2. As shown in Figure 14, due to the discretization and truncation process described above, an infinite number of equally spaced point sound sources moving at the same speed from x = -∞ to x = ∞ on y = -2.5 appear as time aliases in the reproduced sound field.

本来の所望音場は、ｙ＝－２．５上をｘ＝－∞からｘ＝∞まで等速直線運動で移動する１つの点音源によって形成されるべきである。 The desired sound field should be formed by a single point sound source moving at a constant speed in a straight line along y=-2.5 from x=-∞ to x=∞.

そこで、本発明は前記課題を解決するためになされたものであり、その目的は、ＳＤＭによるスピーカアレイの駆動信号を生成する際に、移動する点音源により形成される音場を、時間エイリアスを出現させることなく再現可能な音場再現装置及びプログラムを提供することにある。 The present invention has been made to solve the above problems, and its purpose is to provide a sound field reproduction device and program that can reproduce a sound field formed by a moving point sound source without causing time aliasing when generating a drive signal for a speaker array using SDM.

前記課題を解決するために、請求項１の音場再現装置は、複数のスピーカユニットからなるスピーカアレイを用いて、移動する点音源によって形成される音場を再現するための駆動信号を生成する音場再現装置において、前記スピーカアレイが、ｘｙｚ空間上のｚ＝ｚ₀のｙ＝ｙ₀上に配置されているものとし、時間サンプルをｎ、所定のバッファ区間０≦ｎ≦Ｎ－１における時間サンプルｎ毎の音源信号をｓ（ｎ）及び仮想音源座標を（ｘ_s（ｎ），ｙ_s（ｎ））、再現境界をｙ＝ｙ_ref、サンプリング周波数をｆ_s、スピーカユニット数をＭ（Ｍは自然数）、スピーカユニット間隔をΔｘ、前記音源信号ｓ（ｎ）及び前記仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））が格納されるバッファのバッファ長をＮ、角周波数をω、前記角周波数ωを音速ｃで除算することで得られる波数をｋ、ｘ軸方向の波数をｋ_ｘ、ｍ＝０～Ｍ－１（ｍは整数）、第２種０次ハンケル関数をＨ₀ ⁽²⁾、前記駆動信号の角度スペクトルをＤ^（ｋ_x，ｙ₀，ｚ₀，ω）として、予め設定された前記スピーカユニット数Ｍ及び前記スピーカユニット間隔Δｘ、並びに前記ｍを用いて、式：ｋ_ｘ＝２πｍ／(ΔｘＭ)により、前記波数ｋ_ｘを算出し、前記所定のバッファ区間０≦ｎ≦Ｎ－１における前記音源信号ｓ（ｎ）及び前記仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））、予め設定された前記再現境界ｙ_ref、前記サンプリング周波数ｆ_s、前記波数ｋ及び前記値ｙ₀，ｚ₀、並びに前記波数ｋ_ｘ及び前記角周波数ωを用いて、以下の式：

の演算を行うことで、前記駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を求め、前記駆動信号の角度スペクトルＤ^（ｋ _x ，ｙ ₀ ，ｚ ₀ ，ω）に対して空間周波数領域で逆フーリエ変換を行うことにより、前記スピーカアレイを構成する前記スピーカユニット毎の時間周波数領域の駆動信号を求め、前記スピーカアレイを構成する前記スピーカユニット毎の時間周波数領域の駆動信号に対して時間周波数領域で逆フーリエ変換を行うことにより、前記スピーカアレイを構成する前記スピーカユニット毎の時間領域の駆動信号を求める、ことを特徴とする。 In order to solve the above problem, the sound field reproducing device of claim 1 generates a drive signal for reproducing a sound field formed by a moving point sound source using a speaker array consisting of a plurality of speaker units, the speaker array being arranged on y= _y0 of z= _z0 in an xyz space, a time sample is n, a sound source signal for each time sample n in a predetermined buffer section 0≦n≦N-1 is s(n) and virtual sound source coordinates are ( _xs (n), _ys (n)), a reproduction boundary is y= _yref , a sampling frequency is _fs , the number of speaker units is M (M is a natural number), a speaker unit interval is Δx, a buffer length of a buffer in which the sound source signal s(n) and the virtual sound source coordinates ( _xs (n), _ys (n)) are stored is N, an angular frequency is ω, a wave number obtained by dividing the angular frequency ω by the sound speed c is k, and a wave number in the x-axis direction is kx _. , m=0 to M-1 (m is an integer), a zeroth-order Hankel function of the second kind is H ₀ ⁽²⁾ , and the angular spectrum of the drive signal is D^(k _x , y ₀ , z ₀ , ω) is used to calculate the wave number k _x by the equation: k _x = 2πm/(ΔxM) using the preset number of speaker units M and the speaker unit interval Δx, and m. Then, using the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) in the predetermined buffer section 0≦n≦N-1, the preset reproduction boundary y _ref , the sampling frequency f _s , the wave number k and the values y ₀ , z ₀ , and the wave number k _x and the angular frequency ω, the following equation:

the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal is obtained by performing an inverse Fourier transform in the spatial frequency domain on the angular spectrum D^(kx, y0, z0, ω) of the drive signal _to obtain _a drive _signal in the time frequency domain for each of the speaker units that constitute the speaker array, and the drive signal in the time domain for each of the speaker units that constitute the speaker array is obtained by performing an inverse Fourier transform in the time frequency domain on the drive signals in the time frequency domain for each of the speaker units that constitute the speaker array .

また、請求項２の音場再現装置は、複数のスピーカユニットからなるスピーカアレイを用いて、移動する点音源によって形成される音場を再現するための駆動信号を生成する音場再現装置において、前記スピーカアレイが、ｘｙｚ空間上におけるｚ＝ｚ₀のｙ＝ｙ₀上に配置されているものとし、時間サンプルをｎ、所定のバッファ区間０≦ｎ≦Ｎ－１における時間サンプルｎ毎の音源信号をｓ（ｎ）及び仮想音源座標を（ｘ_s（ｎ），ｙ_s（ｎ））、再現境界をｙ＝ｙ_ref、サンプリング周波数をｆ_s、スピーカユニット数をＭ（Ｍは自然数）、スピーカユニット間隔をΔｘ、前記音源信号ｓ（ｎ）及び前記仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））が格納されるバッファのバッファ長をＮ、角周波数をω、前記角周波数ωを音速ｃで除算することで得られる波数をｋ、ｘ軸方向の波数をｋ_ｘ、ｍ＝０～Ｍ－１（ｍは整数）、第２種０次ハンケル関数をＨ₀ ⁽²⁾、前記駆動信号の角度スペクトルをＤ^（ｋ_x，ｙ₀，ｚ₀，ω）として、予め設定された前記スピーカユニット数Ｍ及び前記スピーカユニット間隔Δｘ、並びに前記ｍを用いて、式：ｋ_ｘ＝２πｍ／（ΔｘＭ）により、前記波数ｋ_ｘを算出する波数算出部と、前記波数ｋの２乗値から、前記波数算出部により算出された前記波数ｋ_ｘの２乗値を減算し、当該減算結果の平方根を求め、前記再現境界ｙ_refから前記値ｙ₀を減算し、当該減算結果の絶対値を求め、前記平方根及び前記絶対値を乗算し、当該乗算結果の前記第２種０次ハンケル関数Ｈ₀ ⁽²⁾を演算することで、再現音場角度スペクトルを求める再現音場角度スペクトル算出部と、前記時間サンプルｎ毎に、前記音源信号ｓ（ｎ）及び前記仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））を入力し、前記所定のバッファ区間０≦ｎ≦Ｎ－１におけるＮ個の前記音源信号ｓ（ｎ）及び前記仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））をバッファに格納し、前記所定のバッファ区間０≦ｎ≦Ｎ－１について、前記波数ｋの２乗値から、前記波数算出部により算出された前記波数ｋ_ｘの２乗値を減算し、当該減算結果の平方根を求め、前記再現境界ｙ_refから前記仮想音源座標の値ｙ_s（ｎ）を減算し、当該減算結果の絶対値を求め、前記平方根及び前記絶対値を乗算し、当該乗算結果の前記第２種０次ハンケル関数Ｈ₀ ⁽²⁾を演算し、前記角周波数ω及び当該時間サンプルｎの乗算結果を前記サンプリング周波数ｆ_sで除算して除算結果を求め、前記波数ｋ_ｘ及び前記仮想音源座標の値ｘ_s（ｎ）を乗算し、当該乗算結果から前記除算結果を減算し、当該減算結果を偏角θとし、複素数平面上における前記偏角θの複素数ｅ^jθを算出し、前記第２種０次ハンケル関数Ｈ₀ ⁽²⁾の演算結果、前記音源信号ｓ（ｎ）及び前記複素数ｅ^jθを乗算することで、当該時間サンプルｎの角度スペクトルを求め、前記角度スペクトルをメモリに格納し、前記メモリから、前記所定のバッファ区間０≦ｎ≦Ｎ－１における前記時間サンプルｎのそれぞれの前記角度スペクトルを読み出し、前記時間サンプルｎのそれぞれの前記角度スペクトルを加算し、当該加算結果を前記サンプリング周波数ｆ_sで除算することで、所望音場角度スペクトルを求める所望音場角度スペクトル算出部と、前記所望音場角度スペクトル算出部により求めた前記所望音場角度スペクトルを、前記再現音場角度スペクトル算出部により求めた前記再現音場角度スペクトルで除算することで、前記駆動信号の角度スペクトルを求める除算部と、前記除算部により求めた前記駆動信号の角度スペクトルに対して空間周波数領域で離散逆フーリエ変換を行うことにより、前記スピーカアレイを構成する前記スピーカユニット毎の時間周波数領域の駆動信号を求める空間周波数領域離散逆フーリエ変換部と、前記空間周波数領域離散逆フーリエ変換部により求めた前記スピーカアレイを構成する前記スピーカユニット毎の時間周波数領域の駆動信号に対して時間周波数領域で離散逆フーリエ変換を行うことにより、前記スピーカアレイを構成する前記スピーカユニット毎の時間領域の駆動信号を求める時間周波数領域離散逆フーリエ変換部と、を備えたことを特徴とする。 A sound field reproduction device according to claim 2 generates a drive signal for reproducing a sound field formed by a moving point sound source using a speaker array consisting of a plurality of speaker units, the speaker array being arranged on y= _y0 of z= _z0 in an xyz space, a time sample is n, a sound source signal for each time sample n in a predetermined buffer section 0≦n≦N-1 is s(n) and virtual sound source coordinates are ( _xs (n), _ys (n)), a reproduction boundary is y= _yref , a sampling frequency is _fs , the number of speaker units is M (M is a natural number), a speaker unit interval is Δx, a buffer length of a buffer in which the sound source signal s(n) and the virtual sound source coordinates ( _xs (n), _ys (n)) are stored is N, an angular frequency is ω, a wave number obtained by dividing the angular frequency ω by the sound speed c is k, and a wave number in the x-axis direction is kx _. a wave number calculation unit that calculates the wave number kx according to the equation: _kx = _2πm / ^{(ΔxM) using the number M of speaker units and the speaker unit interval Δx that are set in advance, m = 0 to M-1 (m is an integer), a zeroth order Hankel function of the second kind as H 0 (2)} _, _and _an angular spectrum of the drive signal as D^(k _{x ,} y 0 , z 0 , ω), a wave number calculation unit that subtracts the squared value of the wave number kx calculated by the wave number calculation unit from the squared value of the wave number _k , calculates the square root of the subtraction result, subtracts the value y ₀ from the reproduction boundary y _ref , calculates the absolute value of the subtraction result, multiplies the square root and the absolute value, and calculates the zeroth order Hankel function of the second kind H ₀ of the multiplication result. a reproduced sound field angular spectrum calculation unit that calculates a reproduced sound field angular spectrum by calculating ⁽²⁾ ; and a reproduced sound field angular spectrum calculation unit that inputs the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) for each time sample n, stores N of the sound source signals s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) in the predetermined buffer section 0≦n≦N-1 in a buffer, subtracts the squared value of the wave number _kx calculated by the wave number calculation unit from the squared value of the wave number k for the predetermined buffer section 0≦n≦N-1, calculates the square root of the subtraction result, subtracts the value y _s (n) of the virtual sound source coordinates from the reproduction boundary y _ref , calculates the absolute value of the subtraction result, multiplies the square root and the absolute value, and calculates the second kind zeroth order Hankel function H ₀ of the multiplication result. ⁽²⁾ , divide the multiplication result of the angular frequency ω and the time sample n by the sampling frequency f _s to obtain a division result, multiply the wave number k _x and the value x _s (n) of the virtual sound source coordinates, subtract the division result from the multiplication result, set the subtraction result as an argument θ, calculate a complex number e ^jθ of the argument θ on a complex plane, multiply the calculation result of the second kind zeroth-order Hankel function H ₀ ⁽²⁾ by the sound source signal s(n) and the complex number e ^jθ to obtain an angular spectrum of the time sample n, store the angular spectrum in a memory, read out the angular spectrum of each of the time samples n in the predetermined buffer section 0≦n≦N−1 from the memory, add up the angular spectrum of each of the time samples n, and store the sum result in the sampling frequency f a division unit that calculates the angular spectrum of the drive signal by dividing the desired sound field angular spectrum calculated by the desired sound field angular spectrum calculation unit by the reproduced sound field angular spectrum calculated by the reproduced sound field angular spectrum calculation unit; _{a spatial} frequency domain discrete inverse Fourier transform unit that calculates a time frequency domain drive signal for each of the speaker units constituting the speaker array by performing a discrete inverse Fourier transform in the spatial frequency domain on the angular spectrum of the drive signal calculated by the division unit; and a time frequency domain discrete inverse Fourier transform unit that calculates a time domain drive signal for each of the speaker units constituting the speaker array by performing a discrete inverse Fourier transform in the time frequency domain on the time frequency domain drive signal for each of the speaker units constituting the speaker array calculated by the spatial frequency domain discrete inverse Fourier transform unit .

さらに、請求項３のプログラムは、コンピュータを、請求項１または２に記載の音場再現装置として機能させることを特徴とする。 Furthermore, the program of claim 3 is characterized in that it causes a computer to function as the sound field reproduction device of claim 1 or 2.

以上のように、本発明によれば、ＳＤＭによるスピーカアレイの駆動信号を生成する際に、移動する点音源により形成される音場を、時間エイリアスを出現させることなく再現することができる。 As described above, according to the present invention, when generating a driving signal for a speaker array using SDM, it is possible to reproduce a sound field formed by a moving point sound source without introducing time aliasing.

本発明の実施形態において想定する音場の構成例を説明する図である。FIG. 2 is a diagram illustrating an example of a configuration of a sound field assumed in an embodiment of the present invention. 本発明の実施形態による音場再現装置の構成例を示すブロック図である。1 is a block diagram showing an example of the configuration of a sound field reproduction device according to an embodiment of the present invention. 駆動信号算出部の構成例を示すブロック図である。4 is a block diagram showing an example of the configuration of a drive signal calculation unit; FIG. 駆動信号算出部の処理例を示すフローチャートである。10 is a flowchart showing an example of processing by a drive signal calculation unit. 再現音場角度スペクトル算出部の構成例を示すブロック図である。11 is a block diagram showing an example of the configuration of a reproduced sound field angular spectrum calculation unit. FIG. 所望音場角度スペクトル算出部の構成例を示すブロック図である。4 is a block diagram showing an example of the configuration of a desired sound field angular spectrum calculation unit. FIG. バッファに格納されたデータの構成例を示す図である。FIG. 13 is a diagram illustrating an example of the configuration of data stored in a buffer. 波面合成法（ＷＦＳ）におけるスピーカアレイの配置例を示す図である。FIG. 1 is a diagram showing an example of speaker array arrangement in wave field synthesis (WFS). 境界音場制御（ＢｏＳＣ）におけるスピーカアレイの配置例を示す図である。FIG. 1 is a diagram showing an example of speaker array arrangement in boundary sound field control (BoSC). 高次アンビソニックス（ＨＯＡ）におけるスピーカアレイの配置例を示す図である。FIG. 1 is a diagram showing an example of speaker array arrangement in Higher Order Ambisonics (HOA). スペクトル除算法（ＳＤＭ）におけるスピーカアレイの配置例を示す図である。FIG. 1 is a diagram showing an example of speaker array arrangement in the spectral division method (SDM). 一般的なＳＤＭを説明する図である。FIG. 1 is a diagram illustrating a typical SDM. 再現領域及び再現境界を説明する図である。FIG. 2 is a diagram illustrating a reproduction area and a reproduction boundary. 非特許文献２の音場再現方法により算出された所望音場の例を示す図である。FIG. 11 is a diagram showing an example of a desired sound field calculated by the sound field reproduction method of Non-Patent Document 2.

以下、本発明を実施するための形態について図面を用いて詳細に説明する。本発明は、ＳＤＭによるスピーカアレイの駆動信号を生成する際に、任意の軌道を描いて移動する点音源によって形成される音場を定式化すると共に、時間窓を考慮することを特徴とする。 The following describes in detail the embodiments of the present invention with reference to the drawings. The present invention is characterized by formulating a sound field formed by a point sound source that moves along an arbitrary trajectory and taking into account a time window when generating a driving signal for a speaker array using SDM.

これにより、滑らかに移動する音源によって形成される音場について、ドップラー効果を含めて再現することができるとともに、周波数領域の離散化によって生じる時間エイリアスの出現を防ぐことができる。 This makes it possible to reproduce the sound field formed by a smoothly moving sound source, including the Doppler effect, while preventing the appearance of time aliasing caused by discretization of the frequency domain.

〔スピーカアレイの駆動信号〕
まず、本発明の実施形態による音場再現装置が生成する、スピーカアレイの駆動信号について説明する。 [Speaker array drive signal]
First, a speaker array driving signal generated by the sound field reproduction device according to the embodiment of the present invention will be described.

音源信号の時間波形ｓ（ｔ）に対しフーリエ変換を行うことにより、音源信号の周波数応答Ｓ（ω₀）が得られるものとする。すなわち、以下の関係が成り立つものとする。

It is assumed that the frequency response S(ω ₀ ) of the sound source signal can be obtained by performing a Fourier transform on the time waveform s(t) of the sound source signal. That is, it is assumed that the following relationship holds:

この音源信号の時間波形ｓ（ｔ）に対し、以下で定義される窓掛けを行う。

Ｔは、時間窓長を示す。 The time waveform s(t) of this sound source signal is subjected to windowing defined as follows.

T denotes the time window length.

そうすると、以下の式が導出される。

Then, the following formula is derived:

このようにして得られた音源のモデルに対し、移動する音源に対する畳み込みを想定する。音源が移動することにより、音源と受音点との間のインパルス応答が時間と共に変化するため、この畳み込みは、時変インパルス応答と音源との畳み込みであると言える。 The model of the sound source obtained in this way is then convolved with a moving sound source. As the sound source moves, the impulse response between the sound source and the sound receiving point changes over time, so this convolution can be said to be a convolution of the time-varying impulse response with the sound source.

図１は、本発明の実施形態において想定する音場の構成例を説明する図である。ｘｙ平面（ｚ＝０とする）において、時刻ｔにおける音源の座標をｒ_s（ｔ）＝（ｘ_s（ｔ），ｙ_s（ｔ））とし、音源から受音点の座標ｒ＝（ｘ，ｙ）までのインパルス応答をｇ（ｒ－ｒ_s（ｔ），ｔ）とする。 1 is a diagram for explaining an example of the configuration of a sound field assumed in an embodiment of the present invention. In the xy plane (z=0), the coordinates of a sound source at time t are defined as r _s (t) = (x _s (t), y _s (t)), and the impulse response from the sound source to the coordinates r = (x, y) of a sound receiving point is defined as g(r - r _s (t), t).

受音点の座標ｒ＝（ｘ，ｙ）における応答は、以下の式のとおり、音源信号の時間波形ｓ（ｔ）と時変インパルス応答との畳み込みで表すことができる。

The response at the coordinates r=(x, y) of the sound receiving point can be expressed as a convolution of the time waveform s(t) of the sound source signal and the time-varying impulse response, as shown in the following equation.

そして、音源から受音点までの間のインパルス応答ｇ（ｒ－ｒ_s（τ），ｔ－τ）は、時間周波数軸にてフーリエ変換することにより、以下の式で表される。

The impulse response g(r−r _s (τ), t−τ) from the sound source to the sound receiving point is expressed by the following equation by Fourier transforming it on the time-frequency axis.

さらに、音源から受音点までの間の伝達関数Ｇ（ｒ－ｒ_s（τ），ω’）は、ｘ－ｘ_s（τ）でフーリエ変換することで、以下の式のように、音源から受音点までの間の角度スペクトルＧ^（ｋ_x，ｙ－ｙ_s（τ），ω’）を用いて表現することができる。

Furthermore, the transfer function G(r-r _s (τ), ω') from the sound source to the sound receiving point can be expressed by Fourier transforming it with x-x _s (τ) using the angular spectrum G^(k _x , y-y _s (τ), ω') from the sound source to the sound receiving point as shown in the following equation.

また、受音点の座標ｒ＝（ｘ，ｙ）の時刻ｔにおける音圧ｐ（ｒ，ｔ）は、時間方向にフーリエ変換することにより、以下の式に変形することができる。

Furthermore, the sound pressure p(r, t) at time t of the sound receiving point with coordinates r=(x, y) can be transformed into the following equation by Fourier transforming it in the time direction.

したがって、所望音場の角度スペクトルＰ^（ｋ_x）は、以下の式で表される。

Therefore, the angular spectrum P^(k _x ) of the desired sound field is expressed by the following equation:

ここで、音源から受音点までの間の角度スペクトルＧ^（ｋ_x，ｙ－ｙ_s（τ），ω）は、以下の式で表される。

Here, the angular spectrum G^(k _x , y−y _s (τ), ω) from the sound source to the sound receiving point is expressed by the following formula.

したがって、点音源の座標ｒ_s（ｔ）＝（ｘ_s（ｔ），ｙ_s（ｔ））が任意の起動で移動する場合の受音点の座標ｒ＝（ｘ，ｙ）における音圧Ｐ（ｒ，ω）は、以下の式にて表される。

Therefore, when the coordinates r _s (t) = (x _s (t), y _s (t)) of a point sound source move with an arbitrary velocity, the sound pressure P(r, ω) at the coordinates r = (x, y) of the sound receiving point is expressed by the following equation.

また、所望音場の角度スペクトルＰ^（ｋ_x）は、以下の式にて表される。

Moreover, the angular spectrum P^(k _x ) of the desired sound field is expressed by the following equation.

したがって、点音源の座標ｒ_s（ｔ）＝（ｘ_s（ｔ），ｙ_s（ｔ））が任意の起動で移動する場合において、この起動により形成される音場を再現する無限直線音源の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、以下の式にて表される。尚、無限直線音源は、ｘｙｚ空間においてｚ＝ｚ₀のｙ＝ｙ₀上に配置されるものとする。

Therefore, when the coordinates r _s (t) = (x _s (t), y _s (t)) of a point sound source moves with an arbitrary motion, the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal of an infinite linear sound source reproducing the sound field formed by this motion is expressed by the following formula: Note that the infinite linear sound source is assumed to be placed on y = y ₀ at z = z ₀ in the xyz space.

実際にコンピュータを用いて駆動信号を算出する場合、音源信号の時間波形ｓ（ｔ）は、通常時間方向に離散化されている。このため、前記式（２１）を時間方向に離散化すると、無限直線音源の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、以下の式で表される。

ここで、ｆ_sはサンプリング周波数であり、ｎ（０≦ｎ≦Ｎ－１）は時間サンプルであり、τ＝ｎ／ｆ_sである。また、Ｎはバッファ長であり、時間窓長はＴ＝（Ｎ－１）／ｆ_sである。 When actually calculating the excitation signal using a computer, the time waveform s(t) of the sound source signal is usually discretized in the time direction. Therefore, when the above formula (21) is discretized in the time direction, the angular spectrum D^( _kx , _y0 , _z0 , ω) of the excitation signal of the infinite linear sound source is expressed by the following formula.

where _fs is the sampling frequency, n (0≦n≦N−1) are the time samples, τ=n/ _fs , N is the buffer length, and the time window length is T=(N−1)/ _fs .

前記式（２２）において、ｙ₀，ｚ₀は、予め設定される。つまり、前記式（２２）の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）は、ｘ軸方向の波数ｋ_x及び所定数の角周波数ωを変数として算出される。 In the above formula (22), _y0 and _z0 are set in advance. That is, the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal in the above formula (22) is calculated using the wave number _kx in the x-axis direction and a predetermined number of angular frequencies ω as variables.

このように、本発明の実施形態による音場再現装置は、前記式（２２）の演算にて、時間周波数及び角度スペクトルを離散化すると共に、有限長での打ち切り（以下、「離散化及び打ち切り処理」という。）を行うことで、ｘ軸方向の波数ｋ_ｘ毎の（スピーカユニット数Ｍ分の）空間周波数領域における駆動信号を求める。 In this way, the sound field reproducing device according to the embodiment of the present invention discretizes the time frequency and the angular spectrum and truncates them to a finite length in the calculation of equation (22) (hereinafter referred to as "discretization and truncation processing"), thereby determining a drive signal in the spatial frequency domain (for the number M of speaker units) for each wave number _kx in the x-axis direction.

音場再現装置は、駆動信号の角度スペクトルに対し離散逆フーリエ変換を行うことで、スピーカアレイ１００を構成するスピーカユニット毎に、時間周波数領域の駆動信号を算出する。そして、音場再現装置は、時間周波数領域の駆動信号に対し離散逆フーリエ変換を行うことで、スピーカユニット毎に、時間領域のスピーカ駆動信号を算出する。これにより、移動する点音源によって形成される音場を再現するためのスピーカ駆動信号が生成される。 The sound field reproduction device calculates a time-domain drive signal for each speaker unit constituting the speaker array 100 by performing a discrete inverse Fourier transform on the angular spectrum of the drive signal. The sound field reproduction device then calculates a time-domain speaker drive signal for each speaker unit by performing a discrete inverse Fourier transform on the time-domain drive signal. This generates a speaker drive signal for reproducing a sound field formed by a moving point sound source.

〔音場再現装置〕
次に、本発明の実施形態による音場再現装置について説明する。図２は、本発明の実施形態による音場再現装置の構成例を示すブロック図である。この音場再現装置１は、駆動信号算出部１０、空間周波数領域離散逆フーリエ変換部１１及び時間周波数領域離散逆フーリエ変換部１２を備えている。 [Sound field reproduction device]
Next, a sound field reproduction device according to an embodiment of the present invention will be described. Fig. 2 is a block diagram showing an example of the configuration of a sound field reproduction device according to an embodiment of the present invention. This sound field reproduction device 1 includes a drive signal calculation unit 10, a spatial frequency domain discrete inverse Fourier transform unit 11, and a time frequency domain discrete inverse Fourier transform unit 12.

駆動信号算出部１０は、予め設定された再現境界ｙ_ref、サンプリング周波数ｆ_s、スピーカユニット数Ｍ、スピーカユニット間隔Δｘ及びバッファ長Ｎを入力する。また、駆動信号算出部１０は、時間サンプルｎ毎に、仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を入力する。仮想音源（点音源）は、任意の軌跡を描いて移動するものとし、その座標が時間サンプルｎ毎の仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））である。 The drive signal calculation unit 10 inputs a preset reproduction boundary y _ref , a sampling frequency f _s , the number of speaker units M, a speaker unit interval Δx, and a buffer length N. The drive signal calculation unit 10 also inputs virtual sound source coordinates (x _s (n), y _s (n)) and a sound source signal s(n) for each time sample n. The virtual sound source (point sound source) is assumed to move along an arbitrary trajectory, and its coordinates are the virtual sound source coordinates (x _s (n), y _s (n)) for each time sample n.

駆動信号算出部１０は、前記式（２２）の演算を行う際に、ｘ軸方向の波数ｋ_ｘ毎に（スピーカユニット数Ｍ分の）駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を求める。そして、駆動信号算出部１０は、ｘ軸方向の波数ｋ_ｘ毎の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を空間周波数領域離散逆フーリエ変換部１１に出力する。駆動信号算出部１０の詳細については後述する。 When performing the calculation of formula (22), the drive signal calculation unit 10 obtains the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal for each wave number _kx in the x-axis direction (for the number _M of speaker units). Then, the drive signal calculation unit 10 outputs the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal for each wave number kx in the x-axis direction to the spatial frequency domain discrete inverse Fourier transform unit 11. Details of the drive signal calculation unit 10 will be described later.

空間周波数領域離散逆フーリエ変換部１１は、予め設定されたスピーカユニット数Ｍ及びスピーカユニット間隔Δｘを入力すると共に、駆動信号算出部１０から、ｘ軸方向の波数ｋ_ｘ毎の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を入力する。 The spatial frequency domain discrete inverse Fourier transform unit 11 inputs the number M of speaker units and the speaker unit spacing Δx that are set in advance, and also inputs the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal for each wave number k _x in the x-axis direction from the drive signal calculation unit 10.

空間周波数領域離散逆フーリエ変換部１１は、駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）に対し離散逆フーリエ変換を行うことで、スピーカアレイ１００を構成するスピーカユニット毎の時間周波数領域の駆動信号Ｄ~（ｘ，ω）を求める。具体的には、空間周波数領域離散逆フーリエ変換部１１は、スピーカユニット数Ｍ及びスピーカユニット間隔Δｘからスピーカアレイ長Ｌ（＝ＭΔｘ）を算出し、前記式（８）の演算を行う。 The spatial frequency domain discrete inverse Fourier transform unit 11 performs a discrete inverse Fourier transform on the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal to obtain the time-frequency domain drive signal D~(x, ω) for each speaker unit constituting the speaker array 100. Specifically, the spatial frequency domain discrete inverse Fourier transform unit 11 calculates the speaker array length L (=MΔx) from the number of speaker units M and the speaker unit spacing Δx, and performs the calculation of the above formula (8).

空間周波数領域離散逆フーリエ変換部１１は、スピーカユニット毎の時間周波数領域の駆動信号Ｄ~（ｘ，ω）を時間周波数領域離散逆フーリエ変換部１２に出力する。 The spatial frequency domain discrete inverse Fourier transform unit 11 outputs the time frequency domain drive signal D~(x, ω) for each speaker unit to the time frequency domain discrete inverse Fourier transform unit 12.

時間周波数領域離散逆フーリエ変換部１２は、空間周波数領域離散逆フーリエ変換部１１からスピーカユニット毎の時間周波数領域の駆動信号Ｄ~（ｘ，ω）を入力する。そして、時間周波数領域離散逆フーリエ変換部１２は、時間周波数領域の駆動信号Ｄ~（ｘ，ω）に対し離散逆フーリエ変換を行うことで、スピーカアレイ１００を構成するスピーカユニット毎の時間領域のスピーカ駆動信号を求める。 The time-frequency domain discrete inverse Fourier transform unit 12 inputs the time-frequency domain drive signal D~(x,ω) for each speaker unit from the spatial frequency domain discrete inverse Fourier transform unit 11. The time-frequency domain discrete inverse Fourier transform unit 12 then performs a discrete inverse Fourier transform on the time-frequency domain drive signal D~(x,ω) to obtain a time-domain speaker drive signal for each speaker unit that constitutes the speaker array 100.

時間周波数領域離散逆フーリエ変換部１２は、スピーカユニット毎の時間領域のスピーカ駆動信号を、対応するスピーカアレイ１００を構成するスピーカユニットへ出力する。 The time-frequency domain discrete inverse Fourier transform unit 12 outputs the time domain speaker drive signal for each speaker unit to the corresponding speaker unit that constitutes the speaker array 100.

（駆動信号算出部１０）
次に、図２に示した駆動信号算出部１０について詳細に説明する。図３は、駆動信号算出部１０の構成例を示すブロック図であり、図４は、駆動信号算出部１０の処理例を示すフローチャートである。 (Drive signal calculation unit 10)
Next, a detailed description will be given of the drive signal calculation unit 10 shown in Fig. 2. Fig. 3 is a block diagram showing an example of the configuration of the drive signal calculation unit 10, and Fig. 4 is a flowchart showing an example of the processing performed by the drive signal calculation unit 10.

この駆動信号算出部１０は、ｘ軸方向波数算出部２０、波数算出部２４、所望音場角度スペクトル算出部２１、再現音場角度スペクトル算出部２２及び除算部２３を備えている。これらの構成部について、図４のフローチャートを用いて説明する。 This drive signal calculation unit 10 includes an x-axis direction wave number calculation unit 20, a wave number calculation unit 24, a desired sound field angle spectrum calculation unit 21, a reproduced sound field angle spectrum calculation unit 22, and a division unit 23. These components will be explained using the flowchart in Figure 4.

駆動信号算出部１０は、予め設定された再現境界ｙ_ref、サンプリング周波数ｆ_s、スピーカユニット数Ｍ、スピーカユニット間隔Δｘ及びバッファ長Ｎを入力する（ステップＳ４０１）。 The drive signal calculation unit 10 receives as input a preset reproduction boundary y _ref , a sampling frequency f _s , the number of speaker units M, a speaker unit interval Δx, and a buffer length N (step S401).

ｘ軸方向波数算出部２０は、予め設定されたスピーカユニット数Ｍ及びスピーカユニット間隔Δｘを入力する。そして、ｘ軸方向波数算出部２０は、スピーカユニット数Ｍ及びスピーカユニット間隔Δｘ、並びにインデックスｍを用いて、ｘ軸方向の波数ｋ_ｘ（＝２πｍ／（ΔｘＭ））を算出する（ステップＳ４０２）。そして、ｘ軸方向波数算出部２０は、ｘ軸方向の波数ｋ_ｘを、所望音場角度スペクトル算出部２１及び再現音場角度スペクトル算出部２２に出力する。波数算出部２４は、予め設定されたサンプリング周波数ｆs及びＤＦＴ点数Ｆを入力する。そして、波数算出部２４は、サンプリング周波数ｆs、ＤＦＴ点数Ｆ及び周波数インデックスｌ（０～Ｆ－１）を用いて、角周波数ω（＝２πｆsｌ／Ｆ）及び波数ｋ（＝ω／ｃ）を算出する（ステップＳ４０２）。そして、波数算出部２４は、波数ｋを所望音場角度スペクトル算出部２１及び再現音場角度スペクトル算出部２２に出力する。 The x-axis direction wave number calculation unit 20 inputs the number M of speaker units and the speaker unit interval Δx that are set in advance. Then, the x-axis direction wave number calculation unit 20 calculates the wave number k _x (=2πm/(ΔxM)) in the x-axis direction using the number M of speaker units and the speaker unit interval Δx, and the index m (step S402). Then, the x-axis direction wave number calculation unit 20 outputs the wave number k _x in the x-axis direction to the desired sound field angular spectrum calculation unit 21 and the reproduced sound field angular spectrum calculation unit 22. The wave number calculation unit 24 inputs the sampling frequency fs and the number of DFT points F that are set in advance. Then, the wave number calculation unit 24 calculates the angular frequency ω (=2πfsl/F) and the wave number k (=ω/c) using the sampling frequency fs, the number of DFT points F, and the frequency index l (0 to F-1) (step S402). Then, the wave number calculation unit 24 outputs the wave number k to the desired sound field angular spectrum calculation unit 21 and the reproduced sound field angular spectrum calculation unit 22 .

ここで、インデックスｍは、前記式（８）にて説明したとおり、ｍ＝０～Ｍ－１（ｍは整数）である。 Here, index m is 0 to M-1 (m is an integer) as explained in formula (8) above.

再現音場角度スペクトル算出部２２は、予め設定された再現境界ｙ_refを入力すると共に、波数算出部２４から波数ｋ、ｘ軸方向波数算出部２０からｘ軸方向の波数ｋ_ｘを入力する。 The reproduced sound field angular spectrum calculation unit 22 receives a preset reproduction boundary y _ref as an input, and also receives the wave number k from the wave number calculation unit 24 and the wave number _kx in the x-axis direction from the x-axis direction wave number calculation unit 20 .

再現音場角度スペクトル算出部２２は、再現境界ｙ_ref、波数ｋ、ｘ軸方向の波数ｋ_ｘ、予め設定された値ｙ₀、及び角周波数ωを用いて、前記式（２２）の分母に相当する再現音場角度スペクトルを算出する（ステップＳ４０３）。前述のとおり、スピーカアレイ１００は、ｘｙｚ空間において予め設定されたｚ＝ｚ₀のｙ＝ｙ₀上に配置されており、値ｙ₀は、スピーカアレイ１００が配置されているｘｙｚ空間のｙ値を示す。そして、再現音場角度スペクトル算出部２２は、再現音場角度スペクトルを除算部２３に出力する。 The reproduced sound field angular spectrum calculation unit 22 calculates the reproduced sound field angular spectrum corresponding to the denominator of the above formula (22) using the reproduction boundary y _ref , the wave number k, the wave number k _x in the x-axis direction, the preset value y ₀ , and the angular frequency ω (step S403). As described above, the speaker array 100 is placed on y=y ₀ of the preset z=z ₀ in the xyz space, and the value y ₀ indicates the y value of the xyz space in which the speaker array 100 is placed. Then, the reproduced sound field angular spectrum calculation unit 22 outputs the reproduced sound field angular spectrum to the division unit 23.

図５は、再現音場角度スペクトル算出部２２の構成例を示すブロック図である。この再現音場角度スペクトル算出部２２は、算出部４０及びメモリ４１を備えている。 Figure 5 is a block diagram showing an example of the configuration of the reproduced sound field angle spectrum calculation unit 22. This reproduced sound field angle spectrum calculation unit 22 includes a calculation unit 40 and a memory 41.

算出部４０は、予め設定された再現境界ｙ_refを入力すると共に、波数算出部２４から波数ｋ、ｘ軸方向波数算出部２０から（スピーカユニット数Ｍ分の）ｘ軸方向の波数ｋ_ｘを入力する。 The calculation unit 40 receives a preset reproduction boundary y _ref , as well as the wave number k from the wave number calculation unit 24 and the wave number _kx in the x-axis direction (for the number M of speaker units) from the x-axis direction wave number calculation unit 20 .

算出部４０は、以下の式にて、前記式（２２）の分母の演算を行うことで、再現音場角度スペクトルを算出し、再現音場角度スペクトルをメモリ４１に格納する。

The calculation unit 40 calculates the reproduced sound field angular spectrum by performing an operation on the denominator of the above-mentioned equation (22) using the following equation, and stores the reproduced sound field angular spectrum in the memory 41.

具体的には、算出部４０は、波数ｋの２乗値からｘ軸方向の波数ｋ_ｘの２乗値を減算し、減算結果の平方根を求め、再現境界ｙ_refから値ｙ₀を減算し、減算結果の絶対値を求め、前記平方根及び前記絶対値を乗算する。そして、算出部４０は、乗算結果の第２種０次ハンケル関数Ｈ₀ ⁽²⁾を演算することで、再現音場角度スペクトルを求める。この再現音場角度スペクトルは、前記式（１８）の音源から受音点までの間の角度スペクトルＧ^（ｋ_x，ｙ－ｙ_s（τ），ω）に対応するデータである。 Specifically, the calculation unit 40 subtracts the squared value of the wave number _kx in the x-axis direction from the squared value of the wave number k, calculates the square root of the subtraction result, subtracts the value _y0 from the reproduction boundary _yref , calculates the absolute value of the subtraction result, and multiplies the square root and the absolute value. The calculation unit 40 then calculates the second kind zeroth order Hankel function _H0 ⁽²⁾ of the multiplication result to obtain a reproduced sound field angular spectrum. This reproduced sound field angular spectrum is data corresponding to the angular spectrum G^( _kx , y- _ys (τ),ω) between the sound source and the sound receiving point in the above equation (18).

このように、再現音場角度スペクトルは、再現境界ｙ_ref等の予め設定されたパラメータのみを用いて算出される。このため、再現音場角度スペクトルは、当該駆動信号算出部１０による駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）の算出処理に先立って予め算出され、メモリ４１に格納される。 In this way, the reproduced sound field angular spectrum is calculated using only preset parameters such as the reproduction boundary _yref , etc. For this reason, the reproduced sound field angular spectrum is calculated in advance and stored in the memory 41 prior to the calculation process of the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal by the drive signal calculation unit 10.

図３及び図４に戻って、所望音場角度スペクトル算出部２１は、予め設定された再現境界ｙ_ref、サンプリング周波数ｆ_s及びバッファ長Ｎを入力すると共に、波数算出部２４から波数ｋ、x軸方向波数算出部２０からｘ軸方向の波数ｋ_ｘを入力する。 Returning to Figures 3 and 4, the desired sound field angular spectrum calculation unit 21 inputs the preset reproduction boundary _yref , sampling frequency _fs , and buffer length N, as well as the wave number k from the wave number calculation unit 24 and the wave number _kx in the x-axis direction from the x-axis direction wave number calculation unit 20.

所望音場角度スペクトル算出部２１は、バッファ長Ｎの時間サンプルｎ毎に、仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を入力する（ステップＳ４０４）。 The desired sound field angular spectrum calculation unit 21 receives the virtual sound source coordinates (x _s (n), y _s (n)) and the sound source signal s(n) for each time sample n of the buffer length N (step S404).

所望音場角度スペクトル算出部２１は、再現境界ｙ_ref、サンプリング周波数ｆ_s、ｘ軸方向の波数ｋ_ｘ及び波数ｋ、既に入力済みのバッファ長Ｎの時間サンプルｎ毎の仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）、並びに角周波数ω等を用いて、前記式（２２）の分子の演算を行うことで、所望音場角度スペクトルを算出する（ステップＳ４０５）。 The desired sound field angular spectrum calculation unit 21 calculates the desired sound field angular spectrum by performing an operation on the numerator of the above equation (22) using the reproduction boundary y _ref , the sampling frequency f _s , the wave number _kx and wave number k in the x-axis direction, the virtual sound source coordinates (x _s (n), y _s (n)) and sound source signal s(n) for each time sample n of the buffer length N that have already been input, and the angular frequency ω, etc. (step S405).

そして、所望音場角度スペクトル算出部２１は、所望音場角度スペクトルを除算部２３に出力する。 Then, the desired sound field angle spectrum calculation unit 21 outputs the desired sound field angle spectrum to the division unit 23.

図６は、所望音場角度スペクトル算出部２１の構成例を示すブロック図である。この所望音場角度スペクトル算出部２１は、入力部３０、バッファ３１、算出部３２及びメモリ３３を備えている。 Figure 6 is a block diagram showing an example of the configuration of the desired sound field angle spectrum calculation unit 21. This desired sound field angle spectrum calculation unit 21 includes an input unit 30, a buffer 31, a calculation unit 32, and a memory 33.

入力部３０は、時間サンプルｎ毎の仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を入力し、これらのデータをバッファ３１に格納する。この場合、バッファ３１には、常に最新のバッファ長Ｎの時間サンプルｎ毎の仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）が格納されるように、入力部３０により格納処理が行われる。 The input unit 30 inputs the virtual sound source coordinates ( _xs (n), _ys (n)) and sound source signal s(n) for each time sample n, and stores these data in the buffer 31. In this case, the input unit 30 performs a storage process so that the buffer 31 always stores the virtual sound source coordinates ( _xs (n), _ys (n)) and sound source signal s(n) for each time sample n of the latest buffer length N.

図７は、バッファ３１に格納されたデータの構成例を示す図である。バッファ３１には、バッファ区間０≦ｎ≦Ｎ－１のバッファ長Ｎの各区間において、バッファ区間０に仮想音源座標（ｘ_s（０），ｙ_s（０））及び音源信号ｓ（０）が格納される。同様に、バッファ３１には、バッファ区間Ｎ－１に仮想音源座標（ｘ_s（Ｎ－１），ｙ_s（Ｎ－１））及び音源信号ｓ（Ｎ－１）が格納される。 7 is a diagram showing an example of the structure of data stored in the buffer 31. In the buffer 31, in each section of a buffer length N, where 0≦n≦N-1, virtual sound source coordinates (x _s (0), y _s (0)) and sound source signal s(0) are stored in buffer section 0. Similarly, in the buffer 31, virtual sound source coordinates (x _s (N-1), y _s (N-1)) and sound source signal s(N-1) are stored in buffer section N-1.

入力部３０は、新たな仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を入力する毎に、バッファ３１に格納されたバッファ区間０～Ｎ－２のデータをバッファ区間１～Ｎ－１へシフトする。そして、入力部３０は、入力した新たな仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を仮想音源座標（ｘ_s（０），ｙ_s（０））及び音源信号ｓ（０）としてバッファ区間０に格納する。つまり、バッファ３１は、時間サンプルｎ毎に更新されるバッファである。 Every time new virtual sound source coordinates (x _s (n), y _s (n)) and sound source signal s(n) are input, the input unit 30 shifts the data in buffer section 0 to N-2 stored in the buffer 31 to buffer section 1 to N-1. Then, the input unit 30 stores the input new virtual sound source coordinates (x _s (n), y _s (n)) and sound source signal s(n) in buffer section 0 as virtual sound source coordinates (x _s (0), y _s (0)) and sound source signal s(0). In other words, the buffer 31 is a buffer that is updated every time sample n.

図６に戻って、算出部３２は、予め設定された再現境界ｙ_ref、サンプリング周波数ｆ_s及びバッファ長Ｎを入力すると共に、波数算出部２４から波数ｋ、ｘ軸方向波数算出部２０からｘ軸方向の波数ｋ_ｘを入力する。そして、算出部３２は、バッファ３１から、バッファ長Ｎ毎に、仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を読み出す。 6 , the calculation unit 32 receives a preset reproduction boundary y _ref , a sampling frequency f _s and a buffer length N, as well as a wave number k from the wave number calculation unit 24 and a wave number _kx in the x-axis direction from the x-axis direction wave number calculation unit 20. The calculation unit 32 then reads out virtual sound source coordinates (x _s (n), y _s (n)) and a sound source signal s(n) from the buffer 31 for each buffer length N.

算出部３２は、バッファ長Ｎ毎に、仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を用いた以下の式（前記式（２２）の分子）にて演算を行い、演算結果をメモリ３３に格納する。そして、算出部３２は、メモリ３３から演算結果を読み出して加算することで、所望音場角度スペクトルを算出する。そして、算出部３２は、所望音場角度スペクトルを除算部２３に出力する。

The calculation unit 32 performs a calculation for each buffer length N using the following equation (the numerator of the above equation (22)) using the virtual sound source coordinates ( _xs (n), _ys (n)) and sound source signal s(n), and stores the calculation results in the memory 33. The calculation unit 32 then reads out the calculation results from the memory 33 and adds them up to calculate the desired sound field angular spectrum. The calculation unit 32 then outputs the desired sound field angular spectrum to the division unit 23.

具体的には、算出部３２は、波数ｋの２乗値からｘ軸方向の波数ｋ_ｘの２乗値を減算し、減算結果の平方根を求め、再現境界ｙ_refから値ｙ_s（ｎ）を減算し、減算結果の絶対値を求め、前記平方根及び前記絶対値を乗算する。そして、算出部３２は、乗算結果の第２種０次ハンケル関数Ｈ₀ ⁽²⁾を演算する。 Specifically, the calculation unit 32 subtracts the squared value of the wave number _kx in the x-axis direction from the squared value of the wave number k, calculates the square root of the subtraction result, subtracts the value _ys (n) from the reproduction boundary _yref , calculates the absolute value of the subtraction result, and multiplies the square root and the absolute value. Then, the calculation unit 32 calculates the second kind zeroth order Hankel function _H0 ⁽²⁾ of the multiplication result.

算出部３２は、角周波数ω及び時間サンプルｎの乗算結果をサンプリング周波数ｆ_sで除算して除算結果を求め、ｘ軸方向の波数ｋ_ｘ及び値ｘ_s（ｎ）を乗算し、この乗算結果から前記除算結果を減算し、減算結果を偏角θ（＝ｋ_xｘ_s（ｎ）－ωｎ／ｆ_s）とする。そして、算出部３２は、複素数平面上における偏角θの複素数ｅ^jθを算出する。 The calculation unit 32 divides the multiplication result of the angular frequency ω and the time sample n by the sampling frequency f _s to obtain the division result, multiplies it by the wave number k _x in the x-axis direction and the value x _s (n), subtracts the division result from this multiplication result, and sets the subtraction result as the argument θ (= k _x x _s (n) - ωn/f _s ).The calculation unit 32 then calculates a complex number e ^jθ of the argument θ on the complex plane.

算出部３２は、第２種０次ハンケル関数Ｈ₀ ⁽²⁾の演算結果、音源信号ｓ（ｎ）及び複素数ｅ^jθを乗算することで、バッファ長Ｎの時間サンプルｎ毎に角度スペクトルを求め、これをメモリ３３に格納する。 The calculation unit 32 multiplies the calculation result of the second kind zeroth order Hankel function H ₀ ⁽²⁾ , the sound source signal s(n) and a complex number e ^jθ to obtain an angular spectrum for each time sample n of the buffer length N, and stores this in the memory 33 .

算出部３２は、バッファ長Ｎの全ての時間サンプルｎについての角度スペクトルの算出処理が完了すると、メモリ３３からバッファ長Ｎのバッファ区間０≦ｎ≦Ｎ－１における時間サンプルｎ毎の角度スペクトルを読み出す。そして、算出部３２は、これらを加算し、サンプリング周波数ｆsで除算することで、所望音場角度スペクトルを求める。この所望音場角度スペクトルは、前記式（２０）の所望音場の角度スペクトルＰ^（ｋ_x）に対応するデータである。 When the calculation process of the angular spectrum for all time samples n of the buffer length N is completed, the calculation unit 32 reads out the angular spectrum for each time sample n in the buffer section 0≦n≦N-1 of the buffer length N from the memory 33. The calculation unit 32 then adds these and divides by the sampling frequency fs to obtain the desired sound field angular spectrum. This desired sound field angular spectrum is data corresponding to the angular spectrum P^(k _x ) of the desired sound field in the above formula (20).

このように、所望音場角度スペクトルは、再現境界ｙ_ref等の予め設定されたパラメータ、並びにバッファ長Ｎの仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）を用いて、バッファ長Ｎ毎に算出される。 In this way, the desired sound field angular spectrum is calculated for each buffer length N using preset parameters such as the reproduction boundary y _ref , as well as the virtual sound source coordinates (x _s (n), y _s (n)) and sound source signal s(n) for the buffer length N.

図３及び図４に戻って、除算部２３は、所望音場角度スペクトル算出部２１から所望音場角度スペクトルを入力すると共に、再現音場角度スペクトル算出部２２から再現音場角度スペクトルを入力する。 Returning to Figures 3 and 4, the division unit 23 inputs the desired sound field angle spectrum from the desired sound field angle spectrum calculation unit 21 and also inputs the reproduced sound field angle spectrum from the reproduced sound field angle spectrum calculation unit 22.

具体的には、除算部２３は、所望音場角度スペクトル算出部２１から、バッファ長Ｎ毎に所望音場角度スペクトルを入力する。また、除算部２３は、当該除算部２３の除算による駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）の算出処理に先立って、再現音場角度スペクトル算出部２２のメモリ４１から再現音場角度スペクトルを読み出す。 Specifically, the division unit 23 inputs the desired sound field angular spectrum for each buffer length N from the desired sound field angular spectrum calculation unit 21. Furthermore, the division unit 23 reads out the reproduced sound field angular spectrum from the memory 41 of the reproduced sound field angular spectrum calculation unit 22 prior to the calculation process of the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal by division by the division unit 23.

除算部２３は、バッファ長Ｎ毎の所望音場角度スペクトルを再現音場角度スペクトルで除算し、駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を求め、これを空間周波数領域離散逆フーリエ変換部１１に出力する（ステップＳ４０６）。 The division unit 23 divides the desired sound field angular spectrum for each buffer length N by the reproduced sound field angular spectrum to obtain the angular spectrum D^( _kx , _y0 , _z0 , ω) of the drive signal, and outputs this to the spatial frequency domain discrete inverse Fourier transform unit 11 (step S406).

駆動信号算出部１０は、所定の終了条件を満たさない限り（ステップＳ４０７：Ｎ）、ステップＳ４０４へ移行する。そして、駆動信号算出部１０は、次のバッファ長Ｎについて駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を求めるステップＳ４０４～Ｓ４０６の処理を繰り返す。一方、駆動信号算出部１０は、所定の終了条件を満たす場合（ステップＳ４０７：Ｙ）、処理を終了する。 The drive signal calculation unit 10 proceeds to step S404 unless a predetermined termination condition is satisfied (step S407: N). The drive signal calculation unit 10 then repeats the processes of steps S404 to S406 to obtain the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal for the next buffer length N. On the other hand, if the predetermined termination condition is satisfied (step S407: Y), the drive signal calculation unit 10 terminates the process.

以上のように、本発明の実施形態の音場再現装置１によれば、駆動信号算出部１０は、予め設定されたスピーカユニット数Ｍ及びスピーカユニット間隔Δｘに基づいて、ｘ軸方向の波数ｋ_ｘ（＝２πｍ／(ΔｘＭ)）を算出し、予め設定されたサンプリング周波数ｆs及びＤＦＴ点数Ｆ等に基づいて、角周波数ω（＝２πｆsｌ／Ｆ）及び波数ｋ（＝ω／ｃ）を算出する。そして、駆動信号算出部１０は、ｘ軸方向の波数ｋ_ｘ（＝２πｍ／(ΔｘＭ)）、予め設定された再現境界ｙ_ref、サンプリング周波数ｆ_s、バッファ長Ｎ及び波数ｋ（＝ω／ｃ）、並びにバッファ長Ｎの時間サンプルｎ毎の仮想音源座標（ｘ_s（ｎ），ｙ_s（ｎ））及び音源信号ｓ（ｎ）等を用いて、前記式（２２）の演算を行うことで、ｘ軸方向の波数ｋ_ｘ毎の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）を求める。 As described above, according to the sound field reproduction device 1 of the embodiment of the present invention, the drive signal calculation unit 10 calculates the wave number _kx (=2πm/(ΔxM)) in the x-axis direction based on the preset number of speaker units M and speaker unit spacing Δx, and calculates the angular frequency ω (=2πfsl/F) and the wave number k (=ω/c) based on the preset sampling frequency fs and number of DFT points F, etc. Then, the drive signal calculation unit 10 calculates the angular spectrum D^( _{kx, y0, z0, ω) of the drive signal for each wave number kx in the x-axis direction by calculating equation (22) using the wave number kx (=2πm/(ΔxM)) in the x-axis direction, the preset reproduction boundary yref} _, _the _sampling _frequency fs, the buffer length N and the wave number k (=ω/c ₎ , and the virtual sound source coordinates ( _xs (n), ys( _n )) and sound source signal s(n) for each time sample _n of the buffer length N.

空間周波数領域離散逆フーリエ変換部１１は、ｘ軸方向の波数ｋ_ｘ毎の駆動信号の角度スペクトルＤ^（ｋ_x，ｙ₀，ｚ₀，ω）に対し空間周波数領域離散逆フーリエ変換を行うことで、スピーカユニット毎の時間周波数領域の駆動信号Ｄ~（ｘ，ω）を求める。 The spatial frequency domain discrete inverse Fourier transform unit 11 performs a spatial frequency domain discrete inverse Fourier transform on the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal for each wave number k _x in the x-axis direction, to obtain a drive signal D~(x, ω) in the time frequency domain for each speaker unit.

時間周波数領域離散逆フーリエ変換部１２は、スピーカユニット毎の時間周波数領域の駆動信号Ｄ~（ｘ，ω）に対し時間周波数領域離散逆フーリエ変換を行うことで、スピーカユニット毎の時間領域のスピーカ駆動信号を求める。そして、時間周波数領域離散逆フーリエ変換部１２は、スピーカ駆動信号を、スピーカアレイ１００を構成するスピーカユニットへ出力する。 The time-frequency domain discrete inverse Fourier transform unit 12 performs a time-frequency domain discrete inverse Fourier transform on the time-frequency domain drive signal D~(x, ω) for each speaker unit to obtain a time-domain speaker drive signal for each speaker unit. The time-frequency domain discrete inverse Fourier transform unit 12 then outputs the speaker drive signal to the speaker units that make up the speaker array 100.

このように、本発明の実施形態では、ＳＤＭによるスピーカアレイ１００の駆動信号を生成する際に、移動する点音源によって形成される音場を定式化し、時間サンプルｎ（＝０～Ｎ－１）毎の窓関数を考慮した前記式（２２）を用いるようにした。 In this way, in an embodiment of the present invention, when generating a drive signal for the speaker array 100 using SDM, the sound field formed by a moving point sound source is formulated, and the above equation (22) is used, which takes into account the window function for each time sample n (= 0 to N-1).

前述のとおり、従来の前記式（９）の駆動信号を用いて音場を再現した場合、空間方向に周波数を離散化したときに時間エイリアスが出現する。これに対し、本発明の実施形態の前記式（２２）を用いた場合には、窓関数を考慮しているため、時間エイリアスが生じる時刻で無音とすることができる。 As mentioned above, when a sound field is reproduced using the conventional driving signal of the above formula (9), time aliasing appears when the frequency is discretized in the spatial direction. In contrast, when the above formula (22) of the embodiment of the present invention is used, it is possible to silence the time when time aliasing occurs because the window function is taken into account.

本発明の実施形態では、ＳＤＭによるスピーカアレイの駆動信号を生成する際に、移動する点音源により形成される音場を、時間エイリアスを出現させることなく再現することができる。 In an embodiment of the present invention, when generating a driving signal for a speaker array using SDM, it is possible to reproduce a sound field formed by a moving point sound source without introducing time aliasing.

以上、実施形態を挙げて本発明を説明したが、本発明は前記実施形態に限定されるものではなく、その技術思想を逸脱しない範囲で種々変形可能である。 The present invention has been described above using embodiments, but the present invention is not limited to the above embodiments and can be modified in various ways without departing from the technical concept thereof.

尚、本発明の実施形態による音場再現装置１のハードウェア構成としては、通常のコンピュータを使用することができる。音場再現装置１は、ＣＰＵ、ＲＡＭ等の揮発性の記憶媒体、ＲＯＭ等の不揮発性の記憶媒体、及びインターフェース等を備えたコンピュータによって構成される。 In addition, a normal computer can be used as the hardware configuration of the sound field reproduction device 1 according to the embodiment of the present invention. The sound field reproduction device 1 is configured by a computer equipped with a CPU, a volatile storage medium such as RAM, a non-volatile storage medium such as ROM, an interface, etc.

音場再現装置１に備えた駆動信号算出部１０、空間周波数領域離散逆フーリエ変換部１１及び時間周波数領域離散逆フーリエ変換部１２の各機能は、これらの機能を記述したプログラムをＣＰＵに実行させることによりそれぞれ実現される。 The functions of the drive signal calculation unit 10, the spatial frequency domain discrete inverse Fourier transform unit 11, and the time frequency domain discrete inverse Fourier transform unit 12 provided in the sound field reproduction device 1 are each realized by causing a CPU to execute a program that describes these functions.

これらのプログラムは、前記記憶媒体に格納されており、ＣＰＵに読み出されて実行される。また、これらのプログラムは、磁気ディスク（フロッピー（登録商標）ディスク、ハードディスク等）、光ディスク（ＣＤ－ＲＯＭ、ＤＶＤ等）、半導体メモリ等の記憶媒体に格納して頒布することもでき、ネットワークを介して送受信することもできる。 These programs are stored in the storage medium and are read and executed by the CPU. In addition, these programs can be distributed by storing them on storage media such as magnetic disks (floppy disks, hard disks, etc.), optical disks (CD-ROMs, DVDs, etc.), and semiconductor memories, and can also be transmitted and received via a network.

１音場再現装置
１０駆動信号算出部
１１空間周波数領域離散逆フーリエ変換部
１２時間周波数領域離散逆フーリエ変換部
２０ｘ軸方向波数算出部
２１所望音場角度スペクトル算出部
２２再現音場角度スペクトル算出部
２３除算部
２４波数算出部
３０入力部
３１バッファ
３２，４０算出部
３３，４１メモリ
１００－１，１００－２，１００－３，１００－４，１００スピーカアレイ
１０１無限直線音源
Ｄ^（ｋ_x，ｙ₀，ｚ₀，ω）駆動信号の角度スペクトル
Ｐ^（ｋ_x）所望音場の角度スペクトル
Ｇ^（ｋ_x，ｙ－ｙ_s（τ），ω）音源から受音点までの間の角度スペクトル
Ｄ~（ｘ，ω）時間周波数領域の駆動信号
ｎ時間サンプル
ｓ（ｎ）音源信号
（ｘ_s（ｎ），ｙ_s（ｎ））仮想音源座標
ｙ_ref 再現境界
ｆ_s サンプリング周波数
Ｍスピーカユニット数
Δｘスピーカユニット間隔
Ｎバッファ長
ｋ_ｘｘ軸方向の波数
ω 角周波数
ｋ波数
Ｈ₀ ⁽²⁾ 第２種０次ハンケル関数
Ｋ₀ ０次変形ベッセル関数 1 Sound field reproduction device 10 Drive signal calculation unit 11 Spatial frequency domain discrete inverse Fourier transform unit 12 Time frequency domain discrete inverse Fourier transform unit 20 x-axis direction wave number calculation unit 21 Desired sound field angular spectrum calculation unit 22 Reproduced sound field angular spectrum calculation unit 23 Division unit 24 Wave number calculation unit 30 Input unit 31 Buffer 32, 40 Calculation units 33, 41 Memories 100-1, 100-2, 100-3, 100-4, 100 Speaker array 101 Infinite linear sound source D^(k _x , y ₀ , z ₀ , ω) Angular spectrum of drive signal P^(k _x ) Angular spectrum of desired sound field G^(k _x , y -y _s (τ), ω) Angular spectrum D~(x, ω) from sound source to sound receiving point Drive signal n in time frequency domain Time sample s(n) Sound source signal (x _s (n), y _s (n)) Virtual sound source coordinate y _ref Reproduction boundary f _s Sampling frequency M Number of speaker units Δx Speaker unit interval N Buffer length k _x Wave number in the x-axis direction ω Angular frequency k Wave number H ₀ ⁽²⁾ Second kind zeroth order Hankel function K ₀ Zeroth order modified Bessel function

Claims

A sound field reproduction device that generates a drive signal for reproducing a sound field formed by a moving point sound source using a speaker array consisting of a plurality of speaker units,
The speaker array is assumed to be arranged at y= _y0 at z= _z0 in the xyz space,
A time sample is n, a sound source signal for each time sample n in a predetermined buffer section 0≦n≦N−1 is s(n), virtual sound source coordinates are (x _s (n), y _s (n)), a reproduction boundary is y=y _ref , a sampling frequency is f _s , the number of speaker units is M (M is a natural number), a speaker unit interval is Δx, a buffer length of a buffer in which the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) are stored is N, an angular frequency is ω, a wave number obtained by dividing the angular frequency ω by the sound speed c is k, a wave number in the x-axis direction is k _x , m=0 to M−1 (m is an integer), a second kind zeroth order Hankel function is H ₀ ⁽²⁾ , and an angular spectrum of the drive signal is D^(k _x , y ₀ , z ₀ , ω),
Calculate the wave number _kx by the formula _kx = 2πm/(ΔxM) using the number M of speaker units and the speaker unit interval Δx that are set in advance, and m;
Using the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) in the predetermined buffer section 0≦n≦N−1, the preset reproduction boundary y _ref , the sampling frequency f _s , the wave number k and the values y ₀ and z ₀ , and the wave number _kx and the angular frequency ω, the following equation:

The angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal is calculated by the above calculation .
performing an inverse Fourier transform on the angular spectrum D^(k _x , y ₀ , z ₀ , ω) of the drive signal in the spatial frequency domain to obtain a drive signal in the time frequency domain for each of the speaker units constituting the speaker array;
A sound field reproduction device characterized in that a time-domain drive signal for each of the speaker units constituting the speaker array is obtained by performing an inverse Fourier transform in the time-frequency domain on the time-frequency domain drive signal for each of the speaker units constituting the speaker array .

A sound field reproduction device that generates a drive signal for reproducing a sound field formed by a moving point sound source using a speaker array consisting of a plurality of speaker units,
The speaker array is assumed to be arranged at y= _y0 at z= _z0 in the xyz space,
A time sample is n, a sound source signal for each time sample n in a predetermined buffer section 0≦n≦N−1 is s(n), virtual sound source coordinates are (x _s (n), y _s (n)), a reproduction boundary is y=y _ref , a sampling frequency is f _s , the number of speaker units is M (M is a natural number), a speaker unit interval is Δx, a buffer length of a buffer in which the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) are stored is N, an angular frequency is ω, a wave number obtained by dividing the angular frequency ω by the sound speed c is k, a wave number in the x-axis direction is k _x , m=0 to M−1 (m is an integer), a second kind zeroth order Hankel function is H ₀ ⁽²⁾ , and an angular spectrum of the drive signal is D^(k _x , y ₀ , z ₀ , ω),
a wave number calculation unit that calculates the wave number _kx by the formula _kx = 2πm/(ΔxM) using the number M of speaker units and the speaker unit interval Δx that are set in advance, and the m;
a reproduced sound field angular spectrum calculation unit that calculates a reproduced sound field angular spectrum by subtracting the squared value of the wave number _kx calculated by the wave number calculation unit from the squared value of the wave number k, calculating the square root of the subtraction result, subtracting the value _y0 from the reproduction boundary _yref , calculating the absolute value of the subtraction result, multiplying the square root and the absolute value, and calculating the second kind zeroth order Hankel function _H0 ⁽²⁾ of the multiplication result;
inputting the sound source signal s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) for each time sample n, and storing N sound source signals s(n) and the virtual sound source coordinates (x _s (n), y _s (n)) in the predetermined buffer section 0≦n≦N−1 in a buffer;
For the predetermined buffer interval 0≦n≦N−1,
subtracting the squared value of the wave number _kx calculated by the wave number calculation unit from the squared value of the wave number k, calculating the square root of the subtraction result, subtracting the value _ys (n) of the virtual sound source coordinates from the reproduction boundary _yref , calculating the absolute value of the subtraction result, multiplying the square root and the absolute value, and calculating the second kind zeroth-order Hankel function _H0 ⁽²⁾ of the multiplication result,
dividing a multiplication result of the angular frequency ω and the time sample n by the sampling frequency _fs to obtain a division result, multiplying the multiplication result by the wave number _kx and the value _xs (n) of the virtual sound source coordinates, subtracting the division result from the multiplication result, setting the subtraction result as a deflection angle θ, and calculating a complex number ^{ejθ of the deflection angle θ on a complex plane} ;
multiplying the calculation result of the second kind zeroth order Hankel function H ₀ ⁽²⁾ , the sound source signal s(n) and the complex number e ^jθ to obtain an angular spectrum of the time sample n, and storing the angular spectrum in a memory;
a desired sound field angular spectrum calculation unit that reads out from the memory the angular spectrum of each of the time samples n in the predetermined buffer section 0≦n≦N−1, adds up the angular spectrum of each of the time samples n, and divides the sum by the sampling frequency f _s to obtain a desired sound field angular spectrum;
a division unit that calculates an angular spectrum of the drive signal by dividing the desired sound field angular spectrum calculated by the desired sound field angular spectrum calculation unit by the reproduced sound field angular spectrum calculated by the reproduced sound field angular spectrum calculation unit;
a spatial frequency domain discrete inverse Fourier transform unit that performs a discrete inverse Fourier transform in a spatial frequency domain on the angular spectrum of the drive signal obtained by the division unit to obtain a drive signal in a time frequency domain for each of the speaker units that constitute the speaker array;
a time-frequency domain discrete inverse Fourier transform unit that performs a discrete inverse Fourier transform in a time-frequency domain on the time-frequency domain drive signals for each of the speaker units constituting the speaker array obtained by the spatial frequency domain discrete inverse Fourier transform unit, to obtain a time-domain drive signal for each of the speaker units constituting the speaker array;
A sound field reproduction device comprising:

A program for causing a computer to function as the sound field reproduction device according to claim 1 or 2.