Modelling Diffusion at the Boundary of a Digital Waveguide MeshSkip other details (including permanent urls, DOI, citation information)
This work is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 3.0 License. Please contact firstname.lastname@example.org to use this work in a way not covered by the license. :
For more information, read Michigan Publishing's access and usage policy.
Page 00000492 Modelling Diffusion at the Boundary of a Digital Waveguide Mesh Joel Laird, Paul Masri, Nishan Canagarajah Digital Music Research Group, University of Bristol 5.01 Merchant Venturers Building, Woodland Road, Bristol BS8 1UB, United Kingdom Tel: +44 (0)117 954-5192; Fax: +44 (0)117 954-5206 email: Joel.Laird@bristol.ac.uk; URL: http://www.fen.bris.ac.uk/elec/dmr/ Abstract Digital waveguide meshes are well established for modelling 2-D and 3-D spaces. In a previous publication  an extension to the mesh perimeter was proposed, to enable the modelling of curved boundaries. This paper presents a method to incorporate the effects of diffusion at these boundaries, which enables the modelling of rough surfaces, such as room walls or musical instrument bodies. The technique works by randomly varying the incident angle of travelling waves, just before they are reflected. This variation in angle is achieved using circulant matrices, to convolve the signals arriving at the perimeter of the waveguide mesh. The resulting diffusion model is lossless and controllable, enabling the modelling of floor and wall surfaces inside a room, or abstract effects such as time-varying diffusion. 1 Introduction Digital waveguide meshes model the propagation of waves through 2-D surfaces and 3-D spaces. Being ideally suited to DSP implementation, they have useful applications in musical instrument synthesis and room acoustics modelling. Alternative methods exist for modelling room acoustics, such as virtual image-source  and raytracing . However, only the digital waveguide mesh implicitly models diffraction at low frequencies . Nevertheless, accurate modelling at the mesh boundaries is not trivial. In a previous publication a method was presented for modelling curved surfaces, such that the uneven edges of the mesh structure are compensated for . This is briefly reviewed in section 2. Although this enables meshes to be constructed with arbitrary shapes, all boundaries are assumed to be smooth. In this paper, a method is presented that models diffuse reflection, so that surfaces of arbitrary roughness can be specified. Two-dimensional meshes have been considered throughout, though the techniques discussed are applicable to three dimensions. Similarly, a triangular lattice structure is assumed, which is desirable because of its almost direction independent dispersion error . 2 Modelling Smooth Boundaries The lattice structure of a digital waveguide mesh causes its edge to be jagged. To overcome this problem, the solution presented in  involves forming the mesh so that it just fits within the desired boundary. 'Rimguides' are then extended from the perimeter nodes of this mesh to the closest point on the boundary, where they are terminated with a reflection. Rimguides are non-integer length waveguides, each consisting of a pure (integer) delay and a first order all-pass filter to model the fractional part . The rimguide impedances are adjusted, so that they cancel out the impedance mismatch between the perimeter nodes and nodes within the mesh. The resultant mesh behaves as though the boundary were smooth; see figure 1. When applied to a circular membrane, the technique was shown to enable accurate modelling at low sampling rates . 3 Modelling Diffusion At smooth boundaries the reflections are specular, or 'mirror-like', such that the angle of reflection is equal to the angle of incidence. Conversely, for diffuse reflections which occur at rough boundaries, energy is scattered in almost every direction regardless of the angle of incidence . Rather than modifying the reflection mechanism of the digital waveguide mesh, a solution was found by preWaveguide Mesh \\ 4 perimeter waveguides 3 perimeter waveguides 5 perimeter " waveguides Rimguides Boudary Figure 1. Modelling a Smooth, Curved, Boundary using Rimguides - 492 - ICMC Proceedings 1999
Page 00000493 altering the angle of incidence. By randomly varying this angle over time, the effect of diffusion can be simulated. When altering the direction of wave propagation it is important to ensure that signal strength and power is conserved. These requirements can be met by using circulant matrices, which have also been used to model diffusion in a class of artificial reverberators known as feedback delay networks . The stability of any matrix can be ensured by distributing its eigenvalues along the unit circle. For a circulant matrix in particular, calculating the coefficients is easy - the discrete inverse Fourier transform of the eigenvalues computes the first row of the matrix. 3.1 Using Circulant Matrices to Rotate Wavefronts Within a waveguide mesh, signals travel into a node from directions evenly distributed about a circle. These may be treated as vector components of a planar wavefront, which collectively describe its direction and amplitude. Convolving them with a circulant matrix can cause a rotation in the propagating wavefront's direction. This may be demonstrated mathematically. A planar wavefront of amplitude S travelling in the direction 6', can be~ represented by vectors, of magnitude si, whose directions, Pi, evenly span a circle: mapping s P P S S2 2 S (a) (b) s "2 Figure 2. (a) Ideal and (b) Real representations of a planar wave at a 3-waveguide perimeter node N=5, X = [ e " ej21, e-ej2J, e- ] (4) etc... Performing the discrete inverse Fourier transform on these eigenvalues gives the coefficients, x, of the circulant matrix A: X0 X1... XN-j A x-l xo x. XN-2 (5) When this matrix is applied to the signal amplitudes, si, the result is the convolved signal amplitudes s,: soi sF 1 A -= s' LN- S (6) i r SP os, 1) Si ~-5Lsin 0J (1) Recombining these signals results in a new signal of amplitude S', travelling in the direction 0': COS Io-- where PN= j.2 sin iN i = 01,...,N-I, N is the number of vectors used to represent the planar wave of amplitude S. See figure 2(a) for an example where N = 3. The values si are effectively the amplitudes of signals that are travelling through waveguides into a node within the mesh. The unit vectors, P,, correspond to the directions of these waveguides at that node. In order to rotate the planar wave through angle (p, the following eigenvalues, X, are used: For N = 3, N-1 1=0 i=o (7) (8) It can be shown that S' = S and 9$ = e. + (,. 3.2 Rotating Wavefronts at the Digital Waveguide Mesh Boundary To simulate diffusion, the wavefront rotation is applied at the perimeter nodes, just prior to reflection. Therefore, along the waveguides that connect perimeter nodes to the rest of the mesh. it is the amplitudes of the signals travelling toward the perimeter that will be convolved. N = 4, x=[ e1 e-i(P X= e' -" e-'jj S For a triangular mesh structure, the number of (2) F signals, N, is usually 3, 4 or 5 (not including the rimguide). However, as can be seen in figure 1. they are not evenly spaced about each node. This results (3) in errors upon rotation. ICMC Proceedings 1999 -493 -
Page 00000494 Perimeter Node Structure: Uncorrected Rotation Errors: < It. 60. C.501. *too,.... |" 0.3....... (Pr 0 30' (degrees) 100, 501 I" - S.-120 0 -. 120.ins ga 'o (degrees) (a) so 5 6 4P0 '46D,.120 60 3 o0 12 (b).30 0 60 (p. 120 e 0so 180 (c) Corrected Rotation Errors:.. 1000 O -12 e0 re) -0 in.(d.e s"60o ( r 30"....0" on 0 (degrees) O so (degrees) (d) lo. (P 30 so 10 m 00 50.1.90 sPr 60 120 9 o 0 ie (e) (f) Figure 3. Plots of the Rotation Error for Different Perimeter Node Structures. (a), (b) and (c) are uncorrected, (d), (e) and (f) have been corrected for small values of (p The difference in vector geometry between the ideal case (discussed in the previous section) and the real case of a perimeter node has two effects upon the calculation of wavefront rotation. Firstly, there is a unique (but nonlinear) mapping between the domains, so that an angle in the ideal case can be translated to an angle in the real case. Secondly, as a result of this nonlinearity, a rotation in the real case is dependent not only upon the amount of rotation in the ideal case, but also the initial angle. Mapping of Wavefront Angles between Domains For a wavefront, S, travelling in direction., in the ideal domain, the individual signal amplitudes, si, can be calculated using equation 1. These signal amplitudes, applied in the real domain to a perimeter node, result in a wavefront travelling in the direction: =, =ZXSi, mi (9) i=o where Pi are the unit vectors for the incoming waveguides to the perimeter node. This achieves a mapping from a wavefront angle in the ideal domain to its angle in the real domain, as shown in figure 2. Analysis of Rotation Errors After convolving the signal amplitudes with the circulant matrix, the rotated angle of the wavefront in the real domain can be found: (N-1 Si mi i=O (10) The resulting error in the rotation angle, (T -9, )- qpr, for wavefronts of different incident angles, is plotted in figures 3(a), 3(b) and 3(c). It can be seen that the greatest errors occur when N = 3, and the rotation is the most accurate when N = 5. 3.3 Improving the Accuracy of the Rotation By applying a function to the rotation angle, (p,, it is possible to improve the accuracy for small values of rotation. This was done by finding the gradient, g=dG~eldp,r, when 9, =0 and Pr =0. The value Pr /g could then be substituted in place of (Pr in the calculation of the matrix eigenvalues (equation (2), (3) or (4)). The improved response is shown in figures 3(d), 3(e) and 3(f). The resulting errors are now only a reduction in rotation angle. This is due to the fact that for positive angles of rotation the errors are negative, and vice-versa. 3.4 Randomising the Rotation to Model Diffusion To apply the rotation method to the modelling of diffusion, the rotation angles can be randomly altered for each sample step. By specifying the maximum - 494 - ICMC Proceedings 1999
Page 00000495 angle of rotation, the amount of diffusion can be contr-olled. A uniform probability distribution, used to randomise the rotation angle, models diffusion where the energy is scattered evenly in all directions, within ihc desired range. Al ternati ve probability distributions could be used to give different diffusion behaviour. For example, a normal distribution would make small scattering angles more likely than large ones. 4 Results A lossless circular mesh was constructed, to model a disc of radiius 30cm with wave speed 300ms'', at a sampling! ralte of 22050Hz. This model was excited * with an off'-centre impulse, and the diffusion range was increased from 00 to ~ 900 over time. A diff'usion mlethod of uniform probability distribution and corrected rotation angle was used. The spectrogram of the output (pressures summed over the entire area of thie disc, to simulate an external listening point) is shown in figure 4. Whilst the reflections are specular the energy in each r~esonant mode remains constant. It is interesting to see that as the diffusion starts to take effect, the energy ente~rs the high-frequency modes that weren't initially excited. As diffusion increases further, significant (random) amplitude modulations are caused at each mode. This results in a widening in freqluency of the modes and the presence of additional noise. By the end of the simulation, noise is the dominant sound. 5 Conclusions The method presented in this paper, to model diffusion, conserved energy yet enabled a range of no diffusion to almost total diffusion. This adds an extension to the digital waveguide mesh, so that the effect of rough surface boundaries can be included. The diffusion was achieved by randomly varying angles of reflection, at each sample step. The variation in reflection angles was done by convolving, with circulant matrices, the signals in waveguides, that travel towards the perimeter nodes. It is possible to simultaneously model different amounts of diffusion, and alter this over time. Frequency dependency was ignored, however, it may be included in the model by replacing the coefficients;, xi, within the circulant matrix A, with filters. Care would have to be taken to ensure that all filters within a row of the matri x were complementary, so that it would be lossless. It may be possible to make the diffusion dependent on the angle of incidence, for specialised surfaces. 'However in this case, non-circulant matrices would Diffusion Range (~cp, ) Figure 4l. Spectrogram of Sound from a Lossless Circular Mesh, with increasing Diffusion have to be used. Ensuring that these matrices conserve signal strength and power makes their design much more complicated. 6 Acknowledgements The authors would like to thank Premier Percussion Ltd. for the funding of this work. 7 References [I] J. Laird; P. Masri; C.N. Canagarajah. 1998. "Efficient and Accurate Synthesis of Circular Membranes using Digital Waveguides" in Proc. lEE Colloquium on A udio and Music Technology: The Challenge of Creative DSP, 12/1-12/6  J.B. Allen; D.A. Berkley. 1979. "Image Method for Efficiently Simulating Small-Room Acoustics" in Journal of the Acoustical Society of America (J.Acous.Soc.Am.): 65(4), pp. 943 -950  A. Kulowski. 1985. "Algorithmic Representation of the Ray Tracing Technique" in Applied Acoustics: 1 8(6), pp. 449-469  L. Savioja; T.J. Rinne; T. Takala. 1994. "Simulation of Room Acoustics with a 3-D Finite Difference Mesh" in Proc. International Computer Music Conference (ICMC), pp. 463 -466  S.A. Van Duyne; J.O. Smith. 1996. "The 3-D Tetrahedral Digital Waveguide Mesh with Musical Applications" in Proc. international Computer Music Conference (ICMC), pp. 9-16  T.I. Laakso; V. Viilim~ki; M. Karjalainen; U.K. Lamne. 1996. "Splitting the Unit Delay" in IEEE Signal Processing Magazine, Jan., pp 30-60 [7) D.E. Hall. 1987. "Room Acoustics" (ch.12) in Basic Acoustics. p.1 86. PubI. John Wiley & Sons, inc.  D. Rocchesso; J.O. Smith. 1994. "Circulant Feedback Delay Networks for Sound Synthesis and Processing" in Proc. international Computer Music Conference (IC'MC), pp 378 -381 ICMC Proceedings 1999-49 (c~95 -