Sound field reproduction device and program

Abstract

To provide a sound field reproduction device capable of reproducing a sound field by using a speaker array that reproduces the wavefront taking into account the influence of the wall when installing a wall that has an opening between the sound source of the desired sound field and the desired boundary when generating a driving signal for speaker array by spectral division method (SDM), and a program.SOLUTION: A drive signal calculator 10 of a sound field reproduction device calculates the sound pressure angle spectrum at wall position yw, using virtual sound source coordinates (xs,ys) with a given buffer length, a sound source signal s(n), and a wall position yw, etc., calculates the sound pressure angle spectrum taking into account the influence of the wall by using y axis direction particle velocity transmittance distribution angle spectrum H(kx), and calculates the angle spectrum D^(kx, y0, z0, ω) of the driving signal by the SDM division formula.SELECTED DRAWING: Figure 3

Description

[Technical Field]

[0001] The present invention relates to a sound field reproduction device and a program for reproducing a sound field using a speaker array. [Background technology]

[0002] Research has been conducted on 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).

[0003] 7 is a diagram showing an example of speaker array arrangement in wave field synthesis (WFS). WFS is a method based on the Rayleigh integral, and uses a speaker array 100-1 arranged linearly or planarly to reproduce the sound pressure or sound pressure gradient on the boundary plane of a desired sound field, thereby forming a wavefront due to a sound source arriving from outside the boundary.

[0004] Figure 8 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 array 100-2 outside the area where a sound field is to be reproduced, facing the inside of the area so as to surround the listener, 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.

[0005] 9 is a diagram showing an example of speaker array arrangement in Higher Order Ambisonics (HOA). HOA is a technique for reproducing a sound field by expressing a desired sound field using spherical harmonic functions and using a spherical speaker array 100-3 arranged to surround a listener facing 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.

[0006] 10 is a diagram showing an example of speaker array arrangement 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 using a speaker array 100-4 arranged linearly or planarly to match the angular spectrum of the reproduced sound field with the angular spectrum of the desired sound field. For details of SDM, see, for example, Non-Patent Document 1.

[0007] [Calculation method of drive signal in SDM] A method for calculating a drive signal in a typical SDM will be described below. Fig. 11 is a diagram for explaining a typical SDM.

[0008] At the point of coordinate r0=(x0, y0, z0) on the infinite linear sound source 101 placed in the xyz space at y=y0, z=z0, the driving signal of angular frequency ω is assumed to be D(r0, ω). ref The sound pressure P(r,ω) at the point z(r, z0) is expressed by the following equation:

number

[0009] 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 as a three-dimensional free-field Green's function, as shown in the following equation.

number

[0010] The sound pressure P(r,ω) can be interpreted as being obtained by convolving the drive signal D(r0,ω) and the transfer function G(r-r0,ω) in space, and y=y ref By performing a Fourier transform along the x-axis direction at x ,y ref ,z0,ω) is obtained.

number

[0011] P^(k x ,y ref ,z0,ω) is the sound pressure P(r,ω) at y=y ref The angular spectrum is obtained by spatial Fourier transform. x ,y0,z0,ω) is the angular spectrum obtained by spatial Fourier transform of the drive signal D(r0,ω) at y=y0. Furthermore, G^(k x ,y ref -y0,z0,ω) is the transfer function G(r-r0,ω) as y=y ref The angular spectrum is obtained by spatial Fourier transform with k x is the wave number in the x-axis direction.

[0012] Here, the wavefront of the desired sound field to be reproduced is y=y des The angular spectrum obtained by spatial Fourier transform of the sound pressure distribution at the desired boundary (hereinafter referred to as the desired boundary) is P^ d (k x ,y des , z0, ω), the angular spectrum D^(k x , y0, z0, ω) are expressed by the following equation.

number

[0013] Coordinate r s =(x s ,y s , z0) is driven by S(ω) at angular frequency ω, the desired sound field y=y des Angular spectrum P^ atd (k x ,y des , z0,ω), i.e., the desired sound field angular spectrum, is expressed by the following equation:

number

[0014] where H0 (2) is the second kind Hankel function of degree 0, and K0 is the 0th order modified Bessel function. 2 -k x 2 The sound waves in the case of 0>k are sound waves that propagate far away. 2 -k x 2 In this case, the sound wave is called an evanescent wave, and the amplitude of the sound wave rapidly decays and does not propagate far.

[0015] Similarly, the angular spectrum G^(k x ,y ref -y0,z0,ω), i.e., the reproduced sound field angular spectrum, is expressed by the following equation:

number

[0016] Therefore, the angular spectrum of the drive signal D^(k x , y0, z0, ω) are expressed by the following formula.

number

[0017] 12 is a diagram illustrating 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 is driven using the drive signal of the above-mentioned equation (7), and thereby y≧y ref Therefore, in this case, y=y ref The reproduction boundary, y ≥ y ref This area is called the clipping point.

[0018] The angular spectrum D^(k x , y0, z0, ω) is assumed to be performed using wave field synthesis using an infinite linear sound source, but in reality it is realized using a speaker array 100 consisting of a finite number of speaker units arranged in a line.

[0019] Therefore, the angular spectrum D^(k x , y0, z0, ω), the frequency and angular spectrum must be discretized, truncated to a finite length, and then an inverse discrete Fourier transform (IDFT) must be performed to calculate the drive signal in the time-frequency domain for each speaker unit.

[0020] 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). x , y0, z0, ω), the driving signal D~(x, ω) in the time-frequency domain for the speaker unit placed at x is expressed by the following equation.

number

[0021] Then, a time domain drive signal is calculated for each speaker unit that constitutes the speaker array 100.

[0022] 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.

[0023] According to conventional speaker drive signal calculation methods, it is assumed that the sound source is installed in free space in the desired sound field in which the wavefront is to be reproduced, and the presence of obstructions such as walls is not taken into consideration. [Prior art documents] [Non-patent literature]

[0024] [Non-Patent Document 1] 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. Summary of the Invention [Problem to be solved by the invention]

[0025] However, when producing audio content using SDM, it is sometimes necessary to take into account occlusions caused by obstacles such as walls or changes in tone.

[0026] Therefore, the present invention has been made to solve the above-mentioned problems, and its purpose is to provide a sound field reproduction device and program that can reproduce a wavefront that takes into account the influence of a wall with an opening when a wall with an opening is placed between the sound source of the desired sound field and the desired boundary when generating a drive signal for a speaker array using SDM. [Means for solving the problem]

[0027] In order to solve the above problem, the sound field reproduction device of claim 1 is a sound field reproduction device that generates a drive signal for reproducing a sound field formed by a desired sound source using a speaker array consisting of a plurality of speaker units, and the speaker array is arranged on y=y0 of z=z0 in an xyz space, where k is a wave number and k x is the wave number in the x-axis direction, Δk x The wave number k x is the discretization interval, ω is the angular frequency, c is the speed of sound, (x s ,y s ) is the virtual sound source coordinates, which are the coordinates of the desired sound source, s(n) is the sound source signal of the desired sound source for each time sample n, y w is the wall position, H(k x ) at the wall position y w y-axis particle velocity transmittance distribution angular spectrum, y des is the desired boundary, y ref is the reproduction boundary, M is the number of speaker units, Δx is the speaker unit spacing, f s is the sampling frequency, F is the number of DFT points, H0 (2) is a second-order Hankel function, i is an imaginary unit, and the number M of speaker units and the interval Δx between speaker units, which are set in advance, and an index m=0 to M−1 (m is an integer), are used to calculate the formula: k x = 2πm / (ΔxM), the wave number k x (=k 1 x ,k 2 x ,···,k M x ) and calculate the discretization interval Δk when the speaker array is configured by the plurality of speaker units at equal intervals. xand a wave number calculation unit in the x-axis direction that calculates the predetermined sampling frequency f s and using the DFT point number F and frequency index l=0 to F-1 (l is an integer), the formula: ω=2πf s a wave number calculation unit that calculates the angular frequency ω from l / F, and calculates the wave number k from the equation: k=ω / c using the angular frequency ω and the sound speed c; and a wave number calculation unit that calculates the wave number k calculated by the x-axis wave number calculation unit. x , the wave number k calculated by the wave number calculation unit, and the preset reproduction boundary y ref and using the value y0, the formula: TIFF0007792263000009.tif13170, a desired sound field angle spectrum calculation unit that calculates a desired sound field angle spectrum, and a division unit that obtains the angular spectrum of the drive signal by dividing the desired sound field angle spectrum calculated by the desired sound field angle spectrum calculation unit by the reproduced sound field angle spectrum calculated by the reproduced sound field angle spectrum calculation unit, wherein the desired sound field angle spectrum calculation unit calculates the angular spectrum of the drive signal by dividing the virtual sound source coordinates (x s ,y s ) and the sound source signal s(n), the wave number k x At each wall position y w Sound pressure angular spectrum P^ d (k x ,y w ) is calculated, and the wave number k calculated by the wave number calculation unit and the wave number k calculated by the x-axis direction wave number calculation unit are calculated. x , and the sound pressure angular spectrum P^ d (k x ,y w ) to obtain the formula: TIFF0007792263000010.tif48170, calculate the vector P^, and x , the previously set y-axis particle velocity transmittance distribution angular spectrum H(k x ), the desired boundary y des and the wall position y w Using the formula: TIFF0007792263000011.tif23170(k n x and k m x where n, m = 1 to M (n, m are integers) and the angular spectrum Ĥ(k n x ,k m x ) and calculate the formula: TIFF0007792263000012.tif30170, the matrix H is calculated, and the discretization interval Δk calculated by the x-axis direction wave number calculation unit is calculated. x , and multiplying the vector P^ and the matrix H to obtain the desired sound field angular spectrum.

[0028] The program of claim 2 is a computer program for configuring a sound field reproduction device that generates a drive signal for reproducing a sound field formed by a desired sound source using a speaker array consisting of a plurality of speaker units, where the speaker array is assumed to be arranged on y=y0 of z=z0 in an xyz space, k is a wave number, k x is the wave number in the x-axis direction, Δk x The wave number k x is the discretization interval, ω is the angular frequency, c is the speed of sound, (x s ,y s ) is the virtual sound source coordinates, which are the coordinates of the desired sound source, s(n) is the sound source signal of the desired sound source for each time sample n, y w is the wall position, H(k x ) at the wall position y w y-axis particle velocity transmittance distribution angular spectrum, y des is the desired boundary, y ref is the reproduction boundary, M is the number of speaker units, Δx is the speaker unit spacing, f s is the sampling frequency, F is the number of DFT points, H0 (2) is a second-order Hankel function, i is an imaginary unit, and the number M of speaker units and the interval Δx between speaker units, which are set in advance, and an index m=0 to M−1 (m is an integer), are used to calculate the formula: k x = 2πm / (ΔxM), the wave number k x (=k1 x ,k 2 x ,···,k M x ) and calculate the discretization interval Δk when the speaker array is configured by the plurality of speaker units at equal intervals. x a wave number calculation unit in the x-axis direction that calculates the predetermined sampling frequency f s and using the DFT point number F and frequency index l=0 to F-1 (l is an integer), the formula: ω=2πf s a wave number calculation unit that calculates the angular frequency ω from l / F, and calculates the wave number k from the equation: k=ω / c using the angular frequency ω and the sound speed c; x , the wave number k calculated by the wave number calculation unit, and the preset reproduction boundary y ref and using the value y0, the formula: TIFF0007792263000013.tif13170, a program for causing a reproduction sound field angular spectrum calculation unit to calculate a reproduction sound field angular spectrum, a desired sound field angular spectrum calculation unit to calculate a desired sound field angular spectrum, and a division unit to obtain 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 reproduction sound field angular spectrum calculated by the reproduction sound field angular spectrum calculation unit, s ,y s ) and the sound source signal s(n), the wave number k x At each wall position y w Sound pressure angular spectrum P^ d (k x ,y w ) is calculated, and the wave number k calculated by the wave number calculation unit and the wave number k calculated by the x-axis direction wave number calculation unit are calculated. x , and the sound pressure angular spectrum P^ d (k x ,y w ) to obtain the formula: TIFF0007792263000014.tif48170, calculate the vector P^, and x , the previously set y-axis particle velocity transmittance distribution angular spectrum H(k x ), the desired boundary y des and the wall position y w Using the formula: TIFF0007792263000015.tif23170(k n x and k m x where n, m = 1 to M (n, m are integers) and the angular spectrum Ĥ(k n x ,k m x ) and calculate the formula: TIFF0007792263000016.tif30170, the matrix H is calculated, and the discretization interval Δk calculated by the x-axis direction wave number calculation unit is calculated. x , and multiplying the vector P^ and the matrix H to obtain the desired sound field angular spectrum. [Effects of the Invention]

[0029] As described above, according to the present invention, when generating a drive signal for a speaker array using SDM, if a wall with an opening is placed between the sound source and the desired boundary of the desired sound field, it is possible to reproduce a wavefront that takes into account the influence of the wall. As a result, it is possible to reproduce the occlusion or change in tone due to the presence of the wall. [Brief explanation of the drawings]

[0030] [Figure 1] 1A and 1B are diagrams illustrating an example of the configuration of a sound field assumed in an embodiment of the present invention. [Figure 2] 1 is a block diagram illustrating an example of the configuration of a sound field reproduction device according to an embodiment of the present invention. [Figure 3] FIG. 3 is a block diagram showing an example of the configuration of a drive signal calculation unit. [Figure 4] 10 is a flowchart illustrating an example of processing by a drive signal calculation unit. [Figure 5] 10 is a block diagram showing an example of the configuration of a reproduced sound field angular spectrum calculation unit. FIG. [Figure 6] 10 is a block diagram showing an example of the configuration of a desired sound field angular spectrum calculation unit. FIG. [Figure 7] FIG. 1 is a diagram showing an example of speaker array arrangement in wave field synthesis (WFS). [Figure 8] FIG. 1 is a diagram showing an example of speaker array placement in boundary sound field control (BoSC). [Figure 9] FIG. 1 is a diagram showing an example of speaker array arrangement in Higher Order Ambisonics (HOA). [Figure 10] FIG. 1 is a diagram showing an example of speaker array arrangement in the spectral division method (SDM). [Figure 11] FIG. 1 is a diagram illustrating a general SDM. [Figure 12] FIG. 2 is a diagram illustrating a reproduction area and a reproduction boundary. DETAILED DESCRIPTION OF THE INVENTION

[0031] The present invention is characterized in that, when generating a driving signal for a speaker array using SDM, a sound field is formulated taking into account the influence of walls.

[0032] [Speaker array drive signal] First, a speaker array driving signal generated by a sound field reproduction device according to an embodiment of the present invention will be described.

[0033] FIG. 1 is a diagram illustrating an example of the configuration of a sound field assumed in an embodiment of the present invention. The coordinates of a point sound source in the xyz coordinate space are (x s ,y s , z0), and z0 and ω are omitted in the equations below.

[0034] y=y w >y s A thin shielding plate (wall 50 with an opening) is installed at the point y <yw Angular sound pressure spectrum P^ of sound waves arriving from the direction d (k x ,y w ) is assumed to be known.

[0035] At this time, y <y w Angular particle velocity spectrum W(k) in the y-axis direction (when no wall 50 exists) x ,y w ) is given by:

number

[0036] In addition, y=y w The equivalent rate distribution of particle velocity in the y-axis direction is defined as h(x), and its angular spectrum expression is the y-axis particle velocity transmittance distribution angular spectrum H(k x ) Then, y>y on the straight wall 50 w Particle velocity angular spectrum W^(k) in the y-axis direction (when a wall 50 exists) x ,y w ) can be expressed as a convolution in the wavenumber domain (spatial frequency domain) as follows:

number

[0037] Furthermore, the wave field in the wavenumber domain is extrapolated to obtain the particle velocity angular spectrum W^(k x ,y w ) to the sound pressure angular spectrum P~(k x ,y des ) to get y=y des >y w The sound pressure angular spectrum P(k) of the desired sound field in x ,y des ) is expressed as follows:

number

[0038] To actually calculate the above integral on a computer, k x It is necessary to perform discretization and finite length truncation for P^ d (k x ,y w ) and P~(k x ,y des ), for |k x In |≦k, it represents a plane wave component that propagates far away, and |k x It represents the evanescent wave component that decays rapidly in |>k.

[0039] Therefore, P^ d (k x ,y w )=0,P~(k x ,y des )=0(|k x In practice, there is no problem even if the calculation is performed under the condition |>k). In other words, the integral range of the above formula (11) can be truncated to a finite value, and further, by performing discrete approximation, the above formula (11) can be expressed as the following formula using the matrix H and the vector P^.

number

[0040] where Δk x is the wave number k in the x-axis direction x The discretization interval of k x is |k n x Let |≦k (n is a positive integer and indicates the index of the discretized wave number. The same applies to m, which will be described later. m, n = 1 to M), and define the following equation, where M is the number of speaker units.

number

number

number

[0041] where H￣(k n x ,k m x ) is as follows:

number

[0042] Therefore, the above equation (4) is k x and P^ in the above equation (5) d (k x ,y des , z0,ω), P~ obtained in the above equation (12), that is, P~(k x ,y des ) (z0, ω are omitted) to obtain the drive signal.

[0043] That is, the above equation (4) is expressed as k x The angular spectrum of the driving signal D^(k x , y0, z0, ω) is the sound pressure angular spectrum of the desired sound field, P~(k x ,y des ) into the reproduced sound field angular spectrum G^(k x ,y ref -y0,z0,ω)

[0044] In this way, the sound field reproduction device according to the embodiment of the present invention calculates the sound pressure angular spectrum P~(k x ,y des ) is calculated. Then, the sound field reproduction device calculates the equation (4) as k x and P^ in the above equation (5) d (k x ,y des , z0,ω), P~(k x ,y des ) to obtain the sound pressure angular spectrum P~(k x ,ydes ) and the reproduced sound field angular spectrum G^(k x ,y ref -y0,z0,ω) to the angular spectrum D^(k x , y0, z0, ω) are calculated.

[0045] [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.

[0046] The excitation signal calculation unit 10 calculates virtual sound source coordinates (x s ,y s ) and the sound source signal s(n) are input, and the wall position y w y-axis particle velocity transmittance distribution angular spectrum H(k x ), wall position y w , desired boundary y des , the reproduction boundary y ref , number of speaker units M, speaker unit spacing Δx, sampling frequency f s , and the number of DFT points F, which is the number of DFT (or FFT) taps. n is the time sample.

[0047] The virtual sound source (point sound source) has fixed virtual sound source coordinates (x s ,y s ) or a directional, moving virtual sound source located at coordinates (x s ,y s ) may be a sound source that moves and rotates along an arbitrary trajectory.

[0048] The drive signal calculation unit 10 calculates the virtual sound source coordinates (x s ,y s ), sound source signal s(n), wall position y w Using the above, the wall position y wSound pressure angular spectrum P^ d (k x ,y w ) is calculated. Then, the drive signal calculation unit 10 calculates the wall position y w Sound pressure angular spectrum P^ d (k x ,y w ), y-axis particle velocity transmittance distribution angular spectrum H(k x ) and so on, the desired sound field angular spectrum P~(k x ,y des The drive signal calculation unit 10 calculates the division formula of the SDM (the above formula (4)) by k x and the wave number k in the x-axis direction is calculated by substituting the equation (12) in place of the equation (5) (corresponding to the equation (4)). x Angular spectrum D^(k x , y0, z0, ω) and outputs them to the spatial frequency domain discrete inverse Fourier transform unit 11. Details of the drive signal calculation unit 10 will be described later.

[0049] The spatial frequency domain discrete inverse Fourier transform unit 11 receives the wave number k in the x-axis direction from the drive signal calculation unit 10. x Angular spectrum D^(k x , y0, z0, ω), as well as the number of speaker units M and the speaker unit interval Δx that are set in advance.

[0050] The spatial frequency domain discrete inverse Fourier transform unit 11 calculates the angular spectrum D^(k x , y0, z0, ω), to obtain the driving signals D~(x,ω) in the time frequency domain 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-1)Δx from the number of speaker units M and the speaker unit spacing Δx, and performs the calculation of the above formula (8).

[0051] The spatial frequency domain discrete inverse Fourier transform unit 11 outputs the drive signal D~(x, ω) in the time frequency domain for each speaker unit to the time frequency domain discrete inverse Fourier transform unit 12.

[0052] The time-frequency domain discrete inverse Fourier transform unit 12 receives the time-frequency domain drive signals 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 signals D~(x,ω) to obtain time-domain speaker drive signals for each speaker unit that constitutes the speaker array 100. As a result, speaker drive signals for reproducing a sound field that takes into account the influence of the wall 50 are generated.

[0053] 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 constituting the speaker array 100 .

[0054] (Drive signal calculation unit 10) Next, a detailed description will be given of the drive signal calculation section 10 shown in Fig. 2. Fig. 3 is a block diagram showing an example of the configuration of the drive signal calculation section 10, and Fig. 4 is a flowchart showing an example of processing by the drive signal calculation section 10.

[0055] The drive signal calculation unit 10 includes an x-axis direction wave number calculation unit 20, a wave number calculation unit 21, a reproduced sound field angular spectrum calculation unit 22, a desired sound field angular spectrum calculation unit 23, and a division unit 24. These components will be described using the flowchart in FIG.

[0056] The drive signal calculation unit 10 calculates the number of speaker units M, the speaker unit interval Δx, the sampling frequency f s , number of DFT points F, reproduction boundary y ref , the particle velocity transmittance distribution angular spectrum H(k x ), wall position y w and the desired boundary y des is input (step S401).

[0057] (x-axis direction wave number calculation unit 20) The x-axis direction wave number calculation unit 20 receives the preset number of speaker units M and the speaker unit interval Δx. Then, the x-axis direction wave number calculation unit 20 calculates the wave number k in the x-axis direction using the number of speaker units M, the speaker unit interval Δx, and a predetermined index m. x (=2πm / (ΔxM)) and the wave number k in the x-axis direction x The discretization interval Δk x is calculated (step S402).

[0058] Here, the x-axis direction wave number calculation unit 20 calculates M wave numbers k corresponding to m=0 to M−1. x (=k 1 x ,k 2 x ,···,k M x ) is calculated. Also, it is assumed that adjacent speaker units constituting the speaker array 100 are equally spaced, and two wave numbers k in the x-axis direction corresponding to any adjacent speaker units are calculated as follows: x Using the wave number k in the x-axis direction x The discretization interval Δk x is calculated.

[0059] The x-axis wave number calculation unit 20 calculates the wave number k in the x-axis direction. x (=k 1 x ,k 2 x ,···,k M x ) to the reproduced sound field angular spectrum calculation unit 22 and the desired sound field angular spectrum calculation unit 23, and the wave number k in the x-axis direction is calculated. x The discretization interval Δk x is output to the desired sound field angular spectrum calculation unit 23.

[0060] (Wave number calculation unit 21) The wave number calculation unit 21 calculates the wave number at a preset sampling frequency f sand the number of DFT points F. Then, the wave number calculation unit 21 inputs the sampling frequency f s and the number of DFT points F, and the predetermined frequency index l, the wave number k (= ω / c) and the angular frequency ω (= 2πf s l / F) is calculated (step S403), where c is the speed of sound, and the frequency index l is l=0 to F-1 (l is an integer).

[0061] The wave number calculation unit 21 outputs the wave number k and the angular frequency ω to the reproduced sound field angular spectrum calculation unit 22 and the desired sound field angular spectrum calculation unit 23.

[0062] (Reproduced sound field angular spectrum calculation unit 22) The reproduced sound field angular spectrum calculation unit 22 calculates the reproduced sound field angular spectrum by using the preset reproduction boundary y ref is input, and the wave number k in the x-axis direction is calculated from the x-axis wave number calculation unit 20. x is input, and the wave number k and angular frequency ω are input from the wave number calculation unit 21.

[0063] The reproduced sound field angular spectrum calculation unit 22 calculates the reproduction boundary y ref , wave number k in the x-axis direction x , wave number k, angular frequency ω, and a preset value y0, the reproduced sound field angular spectrum G^(k x ,y ref -y0, z0, ω) is calculated (step S404).

[0064] As described above, the speaker array 100 is arranged on y=y0 of z=z0, which is set in advance in the xyz space, and the value y0 indicates the y value in the xyz space where the speaker array 100 is arranged. Then, the reproduced sound field angular spectrum calculation unit 22 calculates the reproduced sound field angular spectrum G^(k x ,y ref −y0, z0, ω) is output to the division unit 24.

[0065] 5 is a block diagram showing an example of the configuration of the reproduced sound field angular spectrum calculation unit 22. The reproduced sound field angular spectrum calculation unit 22 includes a calculation unit 30 and a memory 31.

[0066] The calculation unit 30 calculates the preset reproduction boundary y ref is input, and the wave number k in the x-axis direction (for the number M of speaker units) is calculated from the x-axis wave number calculation unit 20. x The wave number k and the angular frequency ω are input from the wave number calculation unit 21.

[0067] The calculation unit 30 calculates the reproduced sound field angular spectrum using the following formula, and stores the reproduced sound field angular spectrum in the memory 31.

number

[0068] Specifically, the calculation unit 30 calculates the wave number k in the x-axis direction from the square of the wave number k. x Subtract the squared value of , find the square root of the subtraction result, and calculate the reproduction boundary y ref The calculation unit 30 subtracts the value y0 from the result of the subtraction, calculates the absolute value of the subtraction result, and multiplies the square root by the absolute value. Then, the calculation unit 30 calculates the zero-order Hankel function H0 (2) The reproduced sound field angular spectrum is calculated by multiplying the calculation result by (-i / 4), where i is the imaginary unit. This reproduced sound field angular spectrum is calculated by multiplying the calculation result by (-i / 4), where i is the imaginary unit. ref Angular spectrum G^(k x ,y ref -y0,z0,ω).

[0069] In this way, the reproduced sound field angular spectrum is ref Therefore, the reproduced sound field angular spectrum is calculated using only the angular spectrum D^(k x , y0, z0, ω) are calculated in advance and stored in the memory 31.

[0070] (Desired sound field angular spectrum calculation unit 23) 3 and 4, the desired sound field angular spectrum calculation unit 23 calculates the virtual sound source coordinates (x s ,y s ) and the sound source signal s(n) (step S405). In addition, the desired sound field angular spectrum calculation unit 23 inputs the preset reproduction boundary y ref , the particle velocity transmittance distribution angular spectrum H(k x ), wall position y w and the desired boundary y des Furthermore, the desired sound field angular spectrum calculation unit 23 receives the wave number k in the x-axis direction from the x-axis direction wave number calculation unit 20. x and the wave number k in the x-axis direction x The discretization interval Δk x is input, and the wave number k and angular frequency ω are input from the wave number calculation unit 21.

[0071] The desired sound field angular spectrum calculation unit 23 calculates the virtual sound source coordinates (x s ,y s ) and the source signal s(n), and the wall position y w For example, if the point sound source is an omnidirectional stationary sound source, the wall position y can be calculated by the above equation (5) (or equation (18) to be described later). w Sound pressure angular spectrum P^ d (k x ,y w , z0,ω) is calculated (step S406).

[0072] The desired sound field angular spectrum calculation unit 23 calculates the desired sound field angular spectrum at the wall position y w Sound pressure angular spectrum P^ d (k x ,y w , z0,ω), the particle velocity transmittance distribution angular spectrum H(k x ), etc., the desired sound field angular spectrum P~(k x ,y des , z0, ω) (step S407). Then, the desired sound field angular spectrum calculation unit 23 calculates the desired sound field angular spectrum P~(k x ,ydes , z0,ω) is output to the division unit 24.

[0073] 6 is a block diagram showing an example of the configuration of the desired sound field angular spectrum calculation unit 23. This desired sound field angular spectrum calculation unit 23 includes a sound pressure angular spectrum calculation unit 40, a matrix H calculation unit 41, a vector P^ calculation unit 42, and a multiplication unit 43.

[0074] The sound pressure angular spectrum calculation unit 40 calculates the virtual sound source coordinates (x s ,y s ) and the sound source signal s(n) are input, and the preset wall position y w The sound pressure angular spectrum calculation unit 40 also receives the wave number k in the x-axis direction from the x-axis direction wave number calculation unit 20. x is input, and the wave number k and angular frequency ω are input from the wave number calculation unit 21.

[0075] For example, when the point sound source is an omnidirectional stationary sound source, the sound pressure angular spectrum calculation unit 40 calculates the desired boundary y des Instead of wall position y w The virtual sound source coordinates (x s ,y s ) and the source signal s(n), as well as the wave number k and the wave number k in the x-axis direction x and wall position y w Using the wall position y w Sound pressure angular spectrum P^ d (k x ,y w , z0,ω) is calculated.

number

[0076] This wall position y w Sound pressure angular spectrum P^ d (k x ,y w ,z0,ω) is the wave number k in the x-axis direction x(=k 1 x ,k 2 x ,···,k M x ) and P^ in the above formula (15) d (k 1 x ,y w ),P^ d (k 2 x ,y w ),···,P^ d (k M x ,y w ) is equivalent to

[0077] The sound pressure angular spectrum calculation unit 40 calculates the wave number k in the x-axis direction. x Wall position y w Sound pressure angular spectrum P^ d (k x ,y w , z0,ω) is output to the vector P^ calculation unit 42.

[0078] In the above example, the sound pressure angular spectrum calculation unit 40 calculates the sound pressure angular spectrum P^ using the above equation (18) which indicates the case where the point sound source is an omnidirectional stationary sound source. d (k x ,y w , z0, ω). On the other hand, when the point sound source has directionality and moves and rotates along an arbitrary trajectory, the sound pressure angular spectrum calculation unit 40 calculates the sound pressure angular spectrum P^ d (k x ,y w , z0,ω) may be calculated.

[0079] In this case, the sound pressure angular spectrum P^ d (k x ,y w For the calculation process of (z0,ω), please refer to Japanese Patent Application No. 2021-134942, which was filed by the same applicant and inventor as the present patent application and was not published at the time of filing the present patent application. d (k x ,yw , z0,ω) varies depending on various conditions such as the stationary or moving state of the point sound source and its directivity. d (k x ,y w , z0,ω) is not limited to the calculation process.

[0080] The matrix H calculation unit 41 calculates the particle velocity transmittance distribution angular spectrum H(k x ), wall position y w and the desired boundary y des is input, and the wave number k in the x-axis direction is calculated from the x-axis wave number calculation unit 20. x is input, and the wave number k and angular frequency ω are input from the wave number calculation unit 21.

[0081] The matrix H calculation unit 41 calculates the wave number k in the x-axis direction. x Two combinations (k n x ,k m x ) for each wave number k in the x-axis direction n x ,k m x , wave number k, wall position y w and the desired boundary y des Using the above equation (16), the angular spectrum H̃(k n x ,k m x ) and calculate the matrix H shown in the above equation (14) which is configured using these as elements. n and m are integers from 1 to M.

[0082] Specifically, the matrix H calculation unit 41 calculates the wave number k in the x-axis direction as shown in the above equation (16). x Two combinations (k n x ,k m x ) from the square of the wave number k to the wave number k in the x-axis direction n x Subtract the squared value of , find the square root of the subtraction result, and calculate the desired boundary y des From wall position y wis subtracted to obtain the subtraction result, (-i) is multiplied by the square root and the subtraction result to obtain the multiplication result, and an exponential function is obtained with the natural number e as the base and the multiplication result as the exponent.

[0083] The matrix H calculation unit 41 divides the exponential function by the square root to obtain the division result, and calculates the wave number k in the x-axis direction. n x wave number k in the x-axis direction m x The result obtained by subtracting k x " the particle velocity transmittance distribution angular spectrum H(k x Then, the matrix H calculation unit 41 calculates the y-axis direction particle velocity transmittance distribution angular spectrum H(k x ”), the angular spectrum H̃(k n x ,k m x ) is found.

[0084] In this way, the matrix H calculation unit 41 calculates the combination of n=1, . . . , M and m=1, . . . , M (k 1 x ,k 1 x ),···,(k 1 x ,k M x ),(k 2 x ,k 1 x ),···,(k 2 x ,k M x ),···,(k M x ,k 1 x ),···,(k M x ,k M x ), the angular spectrum H̃(k 1 x ,k 1 x ),···,H￣(k1 x ,k M x ),H￣(k 2 x ,k 1 x ),···,H￣(k 2 x ,k M x ),···,H￣(k M x ,k 1 x ),···,H￣(k M x ,k M x ) is calculated. As a result, the matrix H shown in the above equation (14) is calculated.

[0085] The matrix H calculation unit 41 outputs the matrix H shown in the above equation (14) to the multiplication unit 43.

[0086] The vector P^ calculation unit 42 calculates the wave number k in the x-axis direction from the sound pressure angular spectrum calculation unit 40. x Sound pressure angle spectrum P^ d (k x ,y w , z0,ω) and the preset wall position y w The vector P^ calculation unit 42 also receives the wave number k in the x-axis direction from the x-axis direction wave number calculation unit 20. x is input, and the wave number k and angular frequency ω are input from the wave number calculation unit 21.

[0087] The vector P^ calculation unit 42 calculates the wave number k in the x-axis direction. x For each wave number k, the square of the wave number k is calculated. x Then, the squared value of the vector P^ shown in the above formula (15) is subtracted, and the square root of the subtraction result is found. 2 -k 1 x 2 ),√(k 2 -k 2 x 2 ),···,√(k 2 -k Mx 2 ) is calculated.

[0088] The vector P^ calculation unit 42 calculates the wave number k in the x-axis direction. x For each square root, the sound pressure angle spectrum P^ d (k x ,y w ,z0,ω) where P^ d (k 1 x ,y w ),P^ d (k 2 x ,y w ),···,P^ d (k M x ,y w ) to obtain each element of the vector P^ shown in the above equation (15). In this way, the vector P^ shown in the above equation (15) is calculated.

[0089] The vector P^ calculation unit 42 outputs the vector P^ shown in the above equation (15) to the multiplication unit 43.

[0090] The multiplication unit 43 receives the matrix H from the matrix H calculation unit 41, the vector P^ from the vector P^ calculation unit 42, and the wave number k in the x-axis direction from the x-axis direction wave number calculation unit 20. x The discretization interval Δk x Enter.

[0091] The multiplication unit 43 calculates the wave number k in the x-axis direction. x The discretization interval Δk x , by multiplying the matrix H and the vector P̂, the desired sound field angular spectrum P̂(k x ,y des , z0, ω) is calculated. This gives the desired sound field angular spectrum P~(k x ,y des , z0,ω) is calculated.

[0092] The multiplication unit 43 multiplies the desired sound field angular spectrum P~(k x ,ydes , z0,ω) is output to the division unit 24.

[0093] (Division section 24) 3 and 4, the division unit 24 receives the reproduced sound field angular spectrum G^(k x ,y ref -y0, z0, ω) is input, and the desired sound field angle spectrum P~(k x ,y des ,z0,ω).

[0094] The division unit 24 calculates the desired sound field angular spectrum P~(k x ,y des ,z0,ω) reproduces the sound field angular spectrum G^(k x ,y ref -y0,z0,ω) to obtain the angular spectrum of the drive signal D^(k x , y0, z0, ω). Then, the division unit 24 calculates the angular spectrum D^(k x , y0, z0, ω) to the spatial frequency domain discrete inverse Fourier transform unit 11 (step S408).

[0095] Unless a predetermined termination condition is satisfied (step S409: N), the drive signal calculation unit 10 proceeds to step S405. Then, the drive signal calculation unit 10 calculates the angular spectrum D^(k x , y0, z0, ω) is calculated by repeating steps S405 to S408. On the other hand, if a predetermined termination condition is satisfied (step S409: Y), the drive signal calculation unit 10 terminates the process.

[0096] As described above, according to the sound field reproduction device 1 of the embodiment of the present invention, the reproduced sound field angular spectrum calculation unit 22 of the drive signal calculation unit 10 calculates the reproduced sound field angular spectrum G^(k x ,y ref-y0,z0,ω) is calculated.

[0097] The desired sound field angular spectrum calculation unit 23 of the drive signal calculation unit 10 calculates the virtual sound source coordinates (x s ,y s ) and the source signal s(n), as well as the wave number k and the wave number k in the x-axis direction x and wall position y w Using the above equation (18), the wall position y w Sound pressure angular spectrum P^ d (k x ,y w , z0,ω) is calculated.

[0098] The desired sound field angular spectrum calculation unit 23 of the drive signal calculation unit 10 calculates the wave number k in the x-axis direction. x Two combinations (k n x ,k m x ) for each wave number k in the x-axis direction n x ,k m x , wave number k, wall position y w and the desired boundary y des Using the above equation (16), the angular spectrum H̃(k n x ,k m x ) to obtain the matrix H shown in the above equation (14).

[0099] The desired sound field angular spectrum calculation unit 23 of the drive signal calculation unit 10 calculates the wave number k, the wave number k x and sound pressure angular spectrum P^ d (k x ,y w ,z0,ω) where P^ d (k 1 x ,y w ),P^ d (k 2 x ,y w ),···,P^ d (k M x ,y w) is used to calculate the vector P̂ according to the above equation (15).

[0100] The desired sound field angular spectrum calculation unit 23 of the drive signal calculation unit 10 calculates the wave number k in the x-axis direction. x The discretization interval Δk x Using the matrix H and the vector P̂, the desired sound field angular spectrum P̂(k x ,y des , z0,ω) is calculated.

[0101] The division unit 24 of the drive signal calculation unit 10 calculates the desired sound field angular spectrum P~(k x ,y des , z0,ω) and the reproduced sound field angular spectrum G^(k x ,y ref -y0,z0,ω) and the wave number k in the x-axis direction is calculated using the SDM division formula. x Angular spectrum D^(k x , y0, z0, ω) are calculated.

[0102] The spatial frequency domain discrete inverse Fourier transform unit 11 calculates the wave number k in the x-axis direction. x Angular spectrum D^(k x , y0, z0, ω), a driving signal D~(x, ω) in the time frequency domain for each speaker unit is obtained.

[0103] 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 signals to the speaker units that make up the speaker array 100.

[0104] As described above, in the embodiment of the present invention, when a drive signal for speaker array 100 using SDM is generated, if a wall with an opening is placed between the sound source of the desired sound field and the desired boundary, it is possible to reproduce a wavefront that takes into account the influence of the wall. As a result, it is also possible to reproduce occlusion or changes in tone due to the presence of the wall.

[0105] Although the present invention has been described above with reference to the embodiments, the present invention is not limited to the above-described embodiments and can be modified in various ways without departing from the technical concept thereof.

[0106] 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.

[0107] 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.

[0108] These programs are stored in the storage medium and are read and executed by the CPU. These programs can also be stored in a storage medium such as a magnetic disk (e.g., a floppy disk, a hard disk), an optical disk (e.g., a CD-ROM, a DVD), or a semiconductor memory and distributed, or can be transmitted and received via a network. [Explanation of symbols]

[0109] 1. Sound field reproduction device 10 Drive signal calculation unit 11 Spatial frequency domain discrete inverse Fourier transform section 12 Time-frequency domain discrete inverse Fourier transform section 20 x-axis wave number calculation unit 21 Wave number calculation section 22 Reproduced sound field angle spectrum calculation unit 23 Desired sound field angle spectrum calculation unit 24 Division part 30 Calculation Unit 31 memory 40 Sound pressure angle spectrum calculation section 41 Matrix H calculation section 42 Vector P^ calculation section 43 Multiplication section 50 Wall 100, 100-1, 100-2, 100-3, 100-4 speaker array 101 Infinite Linear Sound Source (x s ,y s ) Virtual sound source coordinates s(n) sound source signal y w wall position y des desired boundary y ref Reproduced boundary H(k x ) y-axis particle velocity transmittance distribution angular spectrum M Number of speaker units Δx Speaker unit spacing f s Sampling Frequency F DFT score k wavenumber k x Wave number in the x-axis direction Δk x k x The discretization interval of ω angular frequency D^(k x ,y0,z0,ω) Angular spectrum of the driving signal P~(k x ,y des ),P~(k x ,y des ,z0,ω) Desired sound field angular spectrum G^(k x ,y ref -y0,z0,ω) Reproduced sound field angular spectrum

Claims

1. A sound field reproduction device that generates a drive signal for reproducing a sound field formed by a desired sound source using a speaker array consisting of a plurality of speaker units, The speaker array is located at z=z in the xyz space. 0 y = y 0 where k is the wave number and k x is the wave number in the x-axis direction, Δk x The wave number k x The discretization interval of ω is the angular frequency, c is the speed of sound, (x s , y s ) is the virtual sound source coordinates which are the coordinates of the desired sound source, s(n) is the sound source signal of the desired sound source for each time sample n, y w is the wall position, H(k x ) at the wall position y w y-axis particle velocity transmittance distribution angular spectrum, y des is the desired boundary, y ref is the reproduction boundary, M is the number of speaker units, Δx is the speaker unit spacing, f s is the sampling frequency, F is the number of DFT points, H 0 (2) is the Hankel function of the second kind, and i is the imaginary unit. Using the preset number of speaker units M, the speaker unit interval Δx, and an index m=0 to M−1 (m is an integer), the equation: k x = 2πm / (ΔxM), the wave number k x (= k 1 x , k 2 x , ..., k M x ) and calculate the discretization interval Δk when the speaker array is configured by the plurality of speaker units at equal intervals. x an x-axis direction wave number calculation unit that calculates The preset sampling frequency f s and the number of DFT points F and frequency index l = 0 to F-1 (l is an integer), the formula: ω = 2πf s a wave number calculation unit that calculates the angular frequency ω by l / F, and calculates the wave number k by the equation: k = ω / c using the angular frequency ω and the sound speed c; The wave number k calculated by the x-axis direction wave number calculation unit x , the wave number k calculated by the wave number calculation unit, and the preset reproduction boundary y ref and the value y 0 Using the formula: a reproduced sound field angular spectrum calculation unit that calculates a reproduced sound field angular spectrum by a desired sound field angular spectrum calculation unit that calculates a desired sound field angular spectrum; a division unit that obtains 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, The desired sound field angular spectrum calculation unit The virtual sound source coordinates (x s , y s ) and the sound source signal s(n), the wave number k x The wall position y w Sound pressure angular spectrum P^ d (k x , y w ) is calculated, The wave number k calculated by the wave number calculation unit, the wave number k calculated by the x-axis direction wave number calculation unit x , and the sound pressure angular spectrum P^ d (k x , y w ) to obtain the formula: The vector P^ is calculated by The wave number k, the wave number k x , the previously set y-axis direction particle velocity transmittance distribution angular spectrum H(k x ), the desired boundary y des and the wall position y w Using the formula: (k n x and k m x where n, m = 1 to M (n, m are integers) Thus, the angular spectrum H(k n x , k m x ) is calculated, formula: The matrix H is calculated by The discretization interval Δk calculated by the x-axis direction wave number calculation unit x and multiplying the vector P^ and the matrix H to obtain the desired sound field angular spectrum.

2. A computer constituting a sound field reproduction device that generates a drive signal for reproducing a sound field formed by a desired sound source using a speaker array consisting of a plurality of speaker units, The speaker array is located at z=z in the xyz space. 0 y = y 0 where k is the wave number and k x is the wave number in the x-axis direction, Δk x The wave number k x The discretization interval of ω is the angular frequency, c is the speed of sound, (x s , y s ) is the virtual sound source coordinates which are the coordinates of the desired sound source, s(n) is the sound source signal of the desired sound source for each time sample n, y w is the wall position, H(k x ) at the wall position y w y-axis particle velocity transmittance distribution angular spectrum, y des is the desired boundary, y ref is the reproduction boundary, M is the number of speaker units, Δx is the speaker unit spacing, f s is the sampling frequency, F is the number of DFT points, H 0 (2) is the Hankel function of the second kind, and i is the imaginary unit. Using the preset number of speaker units M, the speaker unit interval Δx, and an index m=0 to M−1 (m is an integer), the equation: k x = 2πm / (ΔxM), the wave number k x (= k 1 x , k 2 x , ..., k M x ) and calculate the discretization interval Δk when the speaker array is configured by the plurality of speaker units at equal intervals. x An x-axis direction wave number calculation unit that calculates The preset sampling frequency f s and the number of DFT points F and frequency index l = 0 to F-1 (l is an integer), the formula: ω = 2πf s a wave number calculation unit that calculates the angular frequency ω by l / F, and calculates the wave number k by the equation: k = ω / c using the angular frequency ω and the sound speed c; The wave number k calculated by the x-axis direction wave number calculation unit x , the wave number k calculated by the wave number calculation unit, and the preset reproduction boundary y ref and the value y 0 Using the formula: a reproduced sound field angular spectrum calculation unit that calculates a reproduced sound field angular spectrum by a desired sound field angular spectrum calculation unit that calculates a desired sound field angular spectrum; and a program for causing a division unit to function as a division unit that obtains 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, The desired sound field angular spectrum calculation unit The virtual sound source coordinates (x s , y s ) and the sound source signal s(n), the wave number k x The wall position y w Sound pressure angular spectrum P^ d (k x , y w ) is calculated, The wave number k calculated by the wave number calculation unit, the wave number k calculated by the x-axis direction wave number calculation unit x , and the sound pressure angular spectrum P^ d (k x , y w ) to obtain the formula: The vector P^ is calculated by The wave number k, the wave number k x , the previously set y-axis direction particle velocity transmittance distribution angular spectrum H(k x ), the desired boundary y des and the wall position y w Using the formula: (k n x and k m x where n, m = 1 to M (n, m are integers) Thus, the angular spectrum H(k n x , k m x ) was calculated using the formula: The matrix H is calculated by The discretization interval Δk calculated by the x-axis direction wave number calculation unit x and multiplying the vector P^ and the matrix H to obtain the desired sound field angular spectrum.

Images