Simulation of Room Acoustics using 2-D Digital Waveguide Meshes

Copyright© 2006 IEEE. Reprinted from Proceedings of 31st IEEE International Conference on Acoustics, Speech, and Signal Processing, Vol. 5, pp. 313-316, Toulouse, France, May 2006. This material is posted here with permission of the IEEE. Such permission of the IEEE does not in any way imply IEEE endorsement of your university's web sites products or services. Internal or personal use of this material is permitted. However, permission to reprint/republish this material for advertising or promotional purposes or for creating new collective works for resale or redistribution must be obtained from the IEEE by writing to [email protected]. By choosing to view this document, you agree to all provisions of the copyright laws protecting it. SIMULATION OF ROOM ACOUSTICS USING 2-D DIGITAL WAVEGUIDE MESHES Antti Kelloniemi1) , Vesa Välimäki2) , and Lauri Savioja1) 1) Helsinki University of Technology Telecommunications Software and Multimedia Laboratory P.O. Box 5400, FI-02015 TKK, Finland, emails: [email protected], [email protected] 2) Helsinki University of Technology Laboratory of Acoustics and Audio Signal Processing P.O. Box 3000, FI-02015 TKK, Finland, email: [email protected] ABSTRACT A novel method for simulation of acoustic spaces, such as concert halls or listening rooms, using several 2-D digital waveguide mesh simulations is discussed. The advantages of this approach include reduced computational load, reduced memory usage, and simplified model structure in comparison to a 3-D waveguide mesh simulation. In approximating the modal frequencies of rooms, all the most important lowest modes get modeled, but some higher modes are missing. The proposed method is useful for finding low-frequency modes and for detecting changes in modal distribution when a sound source is moved, for example. As an acoustic visualization tool, the method is superior over a 3-D simulation in that it can isolate a certain layer of the acoustic wave field and it automatically hides waves that propagate in other directions and thus confuse the visualization. 1. INTRODUCTION Modeling of room acoustics is commonly based on geometrical methods, such as the ray tracing or the image-source method [1, 2]. These methods are computationally very intensive and their modeling error becomes large at low frequencies when the wavelength is large compared to the room dimensions. At the low frequencies, the modal frequencies and their spatial distribution are the most interesting features of room acoustics. For shoebox rooms these can be solved analytically, but for more complicated geometries numerical methods are needed [3]. One method for the accurate modeling of the low frequency modal behavior is the digital waveguide mesh [4, 5, 6]. A reason for this method not being popular in room acoustics modeling is its large computational load in the 3-D cases. The computation is significantly less demanding in 2-D as discussed in Section 2. Thus we propose the use of The first author has received funding from the Academy of Finland (project no. 201050) and the Nokia Foundation. The authors would like to thank Tapio Lokki for the reference measurement data. several 2-D meshes to approximately simulate the modes of a room. In Section 3, a case study is presented. The benefits of using 2-D meshes for visualization are investigated in Section 4. 2. MODELING ROOM ACOUSTICS WITH THE DIGITAL WAVEGUIDE MESH 2.1. Digital waveguide mesh algorithm The digital waveguide mesh is a multi-dimensional extension of 1-D digital waveguides, which are widely used in musical instrument modeling [7]. The mesh is constructed of scattering junctions at regular grid locations. The signal values are passed between neighboring junctions by bi-directional unit delays. Of multiple possible mesh topologies, we are using the rectilinear grid, as it is relatively easy to handle and suits well with the rectangular shapes often encountered in rooms. Numerical dispersion errors have been minimized by use of interpolated mesh structure and frequency warping [8]. Signal value p at a certain junction at time step n is calculated as P 2 l Yl hl pl (n − 1) P pc (n) = − pc (n − 2), (1) l Yl hl where Y and h are the admittance and interpolation coefficients of an interconnection, p represents the signal value at a junction, subscript c denotes the junction to be calculated and index l runs through all its axial and diagonal neigboring junctions. In homogenous media the admittances Yl of all connections are equal, so the denominator in (1) reduces to a constant [5]. At the boundaries, reflection coefficients are simulated using the admittance coefficients and the mesh is truncated with an absorbing boundary condition [9]. Because of the rectilinear mesh structure, the simulation is limited to frequencies below one quarter of the sampling frequency √ c N , (2) fs = ∆x where N represents the dimensionality of the mesh, c is the wave propagation speed and ∆x is the spatial interval between two neighboring junctions in the mesh [10]. x×y 20 0 −20 20 2.2. Reduction from a 3-D mesh to multiple 2-D meshes where Lx , Ly , and Lz are the dimensions of the space along its main axis and nx , ny , and nz are the index numbers of the modes [3]. For the axial 1-D standing waves occurring between two parallel walls, two out of the three n-indices in (3) are set to zero. For the 2-D modes, one of the n-indices equals zero, and for the 3-D modes all the indices have non-zero values. As a substitution of the full 3-D model, three 2-D meshes are set to the sizes of the walls of the room. The reflection coefficients at the mesh boundaries are set to approximate the characteristics of the respective wall. The source and receiver locations are projected on the meshes, the simulation is run for a set number of time steps and finally the received responses are added together. As there are three dimensions in a room but the model comprising of three 2-D meshes has six dimensions altogether, each room dimension gets modeled twice. Thus all the 1-D modes are included twice, the 2-D modes get simulated once and the 3-D modes are not included in the results. If a more accurate model is wanted, more wave propagation paths can be simulated using additional meshes. One way to construct such a mesh is to set one of its edges to be of length Lx and setting the other edge length to be a combination of the two other lengths Ly and Lz as ,s 1 1 + 2 . (4) Lyz = 1 L2y Lz This can be seen also as constricting two of the indices in (3) to have equal values ny = nz . The same equation is used for projecting the source and receiver locations on the new mesh. To cover more index combinations, more 2-D meshes of different sizes could be used. However, a great number of modes would then be modeled two or more times. 60 80 100 120 140 160 180 200 y×z 20 Magnitude / dB A three-dimensional digital waveguide mesh is computationally heavy and its memory consumption is large. For example simulation of a lecture hall of 15 × 7 × 2.5 m3 with ∆x = 0.02 m corresponding to fs ≈ 29 kHz would require 33 × 106 3-D junctions excluding the boundaries. This would take 920 × 106 additions and 130 × 106 multiplications per time step and a total of 66 × 106 memory locations. Thus we decided to approximate the room shape with multiple 2-D meshes. The room shape has to be close to rectangular for this approximation to be valid. In a rectangular space, standing waves occur at frequencies s   2  2 2 ny nx c nz + + , (3) f= 2 Lx Ly Lz 40 0 −20 20 40 60 80 100 120 140 160 180 200 x×z 20 0 −20 20 40 60 80 100 120 140 160 180 200 z × xy 20 0 −20 20 40 60 80 100 120 Frequency / Hz 140 160 180 200 Fig. 1. The impulse responses of each of the four 2-D meshes when simulating the response in receiving point A. As the lowest, most important modes get modeled already with the three basic meshes, we decided to use only one additional mesh. One dimension of it was combined from the two largest room dimensions by the rule in (4), and the second mesh dimension was set to equal the smallest room dimension. A more detailed discussion of the model is found in Section 3. This kind of model would require only 1.0×106 2-D junctions for modeling the same space with the same sampling frequency as in the previous 3-D example. This would imply 2.0 × 106 memory locations, updating of which would require 9.8 × 106 additions and 5.9 × 106 multiplications per time step. Compared to the full 3-D model, the memory consumption is thus diminished by 97% and the computational cost by 99%, assuming that additions and multiplications are equally expensive. 3. APPROXIMATION OF THE MODAL FREQUENCIES As a case study we investigated the acoustic modeling of lecture hall T3 at Helsinki University of Technology. The modeling results were compared with measured responses. The hall was empty of furniture during the measurements. The floor dimensions are 12 m × 7.3 m and the roof dimensions are 11.3 m × 6.6 m, because the left side wall and the back wall are sloping. The height of the room is 2.8 meters but a soft absorbent plate is hung at 0.40 m below the ceiling, covering the full area from the back wall to 3.0 m from the front wall. The sound source was placed 2.7 m away from the left side wall and from the front wall at the height of 1.2 m. The (a) 20 1 0 −20 Magnitude / dB 20 40 60 80 100 120 140 160 180 200 (b) 20 0 −20 20 40 60 80 100 120 140 160 180 0.5 200 (c) 20 0 −20 20 40 60 80 100 120 140 160 180 200 2 (d) 20 0 0 −20 20 40 60 80 100 120 Frequency / Hz 140 160 180 200 Fig. 2. The simulated responses of the lecture hall at receiving points A (a) and B (b), and the measured responses at points A (c) and B (d). impulse response was recorded at two locations, point A at 10 m and point B at 5.0 m away from the front wall. Both points were 2.0 m from the left side wall at the height of 1.7 m. The hall was roughly simulated using rectangular meshes with junction spacing ∆x = 0.10 m corresponding to fs ≈ 4.8 kHz. The mesh dimensions are listed in Table 1. The reflection coefficient of the ceiling was set to 0.90, and that of all the other surfaces was set to 0.97. The difference of distance attenuation of 1/r in 2-D compared to that of 1/r2 in 3-D had to be taken into account by setting the reflection coefficients to slightly lower values than those of the surface materials in the hall. The impulse response of a −3 dB/octave lowpass filter was used as excitation signal. The simulation was run for 214 time steps, which corresponds to about 3.5 s. Figure 1 shows the responses of each 2-D mesh separately for the first receiving point. The responses simulated using four 2-D meshes at two receiving points are shown in Fig. 2 along with the measured responses. The 30 lowest theoretical modal frequencies are marked with dashed vertical lines in all frequency plots. It can be seen that the low frequency modes of the simulated response match well with the analytical and measured mode frequencies. However, the error Table 1. Dimensions of the 2-D meshes in the room mode approximation. Mesh Dimensions x×y 116 × 69 y×z 69 × 25 x×z 116 × 25 z × xy 25 × 60 Fig. 3. Instantaneous absolute sound pressure distribution at three time steps on a 11.6 × 6.9 m2 plane excited with a 0.1 s burst of a 43.2 Hz sine signal at the upper right hand corner. In the upper series the absorbing boundary with reflection coefficient value r = 0.70 is on the vertical boundary marked with dashed line, in the lower series it is on the horizontal boundary. At all other boundaries r = 0.98. in heights and Q-values of the modal peaks are significant because of the roughness of the model. Some modes are missing because they are related to 3-D wave propagation that is not modeled with the 2-D layers chosen for this simulation. Also the reflection coefficient values and room dimensions were only roughly estimated. Still the differences seen between the modeling results in the two receiving points are similar with the differences in measurement results. See, for example, the level of the lowest mode that decreases 12 dB in the simulation results and 16 dB in the measurement results when moving from the receiving point A to B. 4. 2-D DIGITAL WAVEGUIDE MESH AS A VISUALIZATION TOOL As signal values passing through a digital waveguide mesh present actual sound pressure values, wave propagation can be followed by visualizing the values at each junction point. The visualization task is demanding if the full state of a three dimensional mesh is drawn, as already after few reflections the wavefront has typically split into multiple waves moving in many directions. One possible solution is to choose only one layer of junctions of the 3-D model to be visualized. This way the picture becomes more understandable, but still the waves passing through the chosen layer seem to come from nowhere and to vanish unexpectedly, which makes following of the visualization quite difficult [5]. If we choose to model the space again with 2-D meshes 6m 9m 9m 6m 1 0.5 the 2-D method outperforms the 3-D construction when used for visualization: besides being much faster, modeling in two dimensions leaves out waves traveling through the selected layer and makes it easier to follow the sound propagation on the plane of interest. 4m 4m 6. REFERENCES 2 0 3m 3m Fig. 4. Sound pressure distribution in a room depends much of the location of the sound source marked with a black dot. A plane is excited with a 0.1 s burst of 140 Hz sine signal. Absolute sound pressure distribution is plotted at time instant 0.2 s. In the left image the source is located at (1.3, 1.0) m from the upper left hand corner of the room, in the right image it is located at (1.8, 1.0) m. The reflection coefficient of the walls is r = 0.98. as proposed in Section 2, only the wavefronts moving at the level of the chosen layer are to be seen [11]. By investigating multiple such layers independently of each other, the visualizations become much more usable compared to the use of 3-D models. Some reflections will not appear in the visualizations, but the most important wavefronts, standing waves and flutter echoes are seen more clearly now. For example, the waves propagating on different layers can be seen in Figs. 3 and 4. The absolute sound pressure distribution is plotted. A single Hanning-windowed sine frequency burst of 0.1 s was used as input signal. In Fig. 3 it is shown how the location of absorbing wall material affects the amplitude of the standing wave. As the input signal has the frequency of a standing wave born vertically in the figures, absoption on the vertical boundary is less effective compared to the absorption on the horizontal boundary. Figure 4 illustrates how moving of the sound source 0.5 m away from the sloping wall affects the form of the sound field. For example the amplitude of the standing wave close to the vertical wall on the left hand side is significantly different in the two shown cases. 5. CONCLUSIONS The idea of simulating the sound field in a room using twodimensional digital waveguide meshes was introduced. When one three-dimensional simulation is replaced with at least three 2-D meshes, the computational load and memory consumption are diminished dramatically. A semi-rectangular hall was studied, where the new method was shown to produce modes at correct frequencies. As a result of roughness of the model, error is caused in amplitudes and Q-values of the modal peaks. Some modes not occurring at the selected 2-D planes are left out of the simulation results. On the other hand, [1] A. Krokstad, S. Strom, and S. Sorsdal, “Calculating the acoustical room response by the use of a ray tracing technique,” J. Sound Vib., vol. 8, no. 1, pp. 118–125, 1968. [2] J. B. Allen and D. A. Berkley, “Image method for efficiently simulating small-room acoustics,” J. Acoust. Soc. Am., vol. 65, no. 4, pp. 943–950, 1979. [3] H. Kuttruff, Room acoustics, Spon Press, London, 4th edition, 2000. [4] J. O. Smith, “Physical modeling using digital waveguides,” Comp. Music J., vol. 16, no. 4, pp. 74–91, 1992. [5] L. Savioja, T. Rinne, and T. Takala, “Simulation of room acoustics with a 3-D finite difference mesh,” in Proc. Int. Computer Music Conf. (ICMC), Aarhus, Denmark, September 1994, pp. 463–466. [6] L. Savioja, J. Backman, A. Järvinen, and T. Takala, “Waveguide mesh method for low-frequency simulation of room acoustics,” in Proc. 15th Int. Congr. Acoust. (ICA’95), Trondheim, Norway, June 1995, vol. 2, pp. 637–640. [7] S. Van Duyne and J. O. Smith, “The 2-D digital waveguide mesh,” in Proc. IEEE Workshop on Applications of Signal Processing to Audio and Acoustics (WASPAA’93), New Paltz, NY, October 1993. [8] L. Savioja and V. Välimäki, “Reducing the dispersion error in the digital waveguide mesh using interpolation and frequency warping techniques,” IEEE Trans. Speech and Audio Processing, vol. 8, pp. 184 – 194, March 2000. [9] A. Kelloniemi, “Improved adjustable boundary condition for the 2-D digital waveguide mesh,” in Proc. Int. Conf. Digital Audio Effects (DAFx’05), Madrid, Spain, September 2005, pp. 237–242. [10] L. Savioja, J. Huopaniemi, T. Lokki, and R. Väänänen, “Creating interactive virtual acoustic environments,” J. Audio Eng. Soc., vol. 47, no. 9, pp. 675–705, September 1999. [11] T. Yokota, S. Sakamoto, and H. Tachibana, “Visualization of sound propagation and scattering in rooms,” Acoust. Sci. & Tech., vol. 23, no. 1, pp. 40–46, 2002.