Nicolas Tsingos

Ingrid Carlbom

Gary Elko

Thomas Funkhouser

Robert Kubli

Validation of Acoustical Simulations in the "Bell Labs Box"

Introduction

Computer simulated sound propagation through 3D environments is important in many applications, including computer-aided design, training, and virtual reality. In many cases, the accuracy of the acoustical simulation is critical to the success of the application. For example, in concert hall and factory design (where OSHA sound limits must be met), the accuracy of the simulation may save costly re-engineering after construction. In virtual environments, experiments have shown that more accurate acoustic modeling provides a stronger sense of presence [START_REF] Durlach | Virtual reality scientific and technological challenges[END_REF]. Furthermore, auditory cues aid in formation of spatial impressions, separation of simultaneous sound signals, and localization of objects [START_REF] Blauert | Spatial Hearing : The Psychophysics of Human Sound Localization[END_REF], such as when a soldier locates an enemy in a training exercise or a fire-fighter locates a person stranded in a burning building. In contrast, incorrect auditory cues can lead to negative training [START_REF] Welch | Perceptual Modification: Adapting to Altered Sensory Environments[END_REF][START_REF] Boff | Handbook of Perception and Human Performance[END_REF].

Although several systems (e.g., [START_REF]CATT-Acoustic[END_REF][START_REF] Naylor | Odeon -another hybrid room acoustical model[END_REF][START_REF]Bose Modeler[END_REF]) are available for computing sound propagation in 3D environments, there has not been detailed evaluation of their accuracy. The primary reason is that the acoustics of most real world environments is very complex, and thus detailed quantitative comparison of measured with computed impulse responses is difficult. As a result, acousticians have resorted to comparing gross statistics of impulse responses (e.g., reverberation time, clarity) [START_REF] Kuttruff | Room Acoustics[END_REF] and/or using human listening tests for validation of computer simulations [START_REF] Ando | Concert Hall Acoustics[END_REF]. This is similar to the situation in computer graphics during the mid-1980s. At that time, there were several global illumination algorithms that simulated the propagation of light through a 3D environment and produced an image for a given camera (e.g., ray £ The first author is now with the REVES research project at INRIA Sophia-Antipolis, France. Nicolas.Tsingos@sophia.inria.fr carlbom@research.bell-labs.com Ý Summit, NJ, USA. gwe@ieee.org Þ funk@cs.princeton.edu tracing [START_REF] Whitted | An improved illumination model for shaded display[END_REF] and/or radiosity [START_REF] Goral | Modeling the interaction of light between diffuse surfaces[END_REF]). Yet, there were very few results regarding the accuracy of the simulations. Computer-generated pictures and corresponding photographs rarely matched, and little insight into the causes of the mismatch could be derived from the comparisons due to the complexity of light transport in most scenes.

In response to this situation, researchers at Cornell University [START_REF] Goral | Modeling the interaction of light between diffuse surfaces[END_REF][START_REF] Michael | The hemi-cube, a radiosity solution for complex environments[END_REF][START_REF] Meyer | An experimental evaluation of computer graphics imagery[END_REF] constructed a simple real-world scene, the "Cornell Box", in which the lighting, geometry, and reflectance properties of every surface were carefully measured and duplicated in a 3D computer graphics model. Using this simple and controlled experimental setup, they were able to make meaningful comparisons between global illumination simulations and photographs of the scene and provide explanations for the differences.

Motivated by the success of the Cornell Box, we have built a simple experimental setup for validating sound propagation simulations. We call it the "Bell Labs Box." The setup comprises a simple configuration of planar surfaces with a speaker and a microphone (see Figure 1). While the basic configuration of the room is a simple six-sided box, it is constructed with reconfigurable panels that can be inserted or removed to create a variety of interesting geometries, including ones with diffracting panels (see Figure 2).

The key ideas behind this experimental setup are simplicity and control. First, the room has just a few planar surfaces, and thus possible sequences of reflections and diffractions are easily computed in a simulation and recognized in an impulse response. Second, the speaker radiation pattern, microphone directivity pattern, and reflectance properties of every surface are independently measured in an anechoic chamber, providing simulation parameters that closely match the real world environment. Finally, the room is configurable, providing a mechanism by which acoustical simulations can be validated with a variety of geometric effects. These features enable us to compare the results of simulations with measurement data in a simple and controlled setting, thereby allowing more detailed quantitative validation than has previously been possible.

Background and Previous Work

The problem addressed by this work is the validation of a solution to an integral equation expressing the wavefield (acoustic pressure and particle velocity) at some point in space in terms of the wavefield at other points (or equivalently on surrounding surfaces). For light simulations, the wave equation is described by Kajiya's rendering equation [START_REF] Kajiya | The rendering equation[END_REF]. For sound simulations, it is described by the Helmoltz-Kirchoff integral theorem [START_REF] Born | Principles of Optics[END_REF], which incorporates time and phase dependencies.

Since sound and light are both wave phenomena, the methods for their simulation and validation share many features. Yet, sound has characteristics different from light which introduce new and interesting problems:

¯Wavelength: the wavelengths of audible sound are five to seven orders of magnitude longer than visible light, ranging between 0.02 and 17 meters (for 20KHz and 20Hz, respectively). Diffraction of sound occurs around obstacles of the same size as the wavelength (such as tables) and reflections are primarily specular for large, flat surfaces (such as walls), while small objects (like coffee mugs) have little effect on the sound field (for all but the highest wavelengths). As a result, when compared to computer graphics, acoustics simulations tend to use 3D models with far less geometric detail. But, they must find propagation paths with specular reflections and diffractions efficiently, and they must consider the effects for different obstacles at a range of wavelengths.

¯Speed: at 343 meters per second, the speed of sound in air is six orders of magnitude less than light, and sound propagation delays are perceptible to humans. Thus, acoustic models must compute the exact time/frequency distribution of the propagation paths, and the source sound must be auralized by convolution with the corresponding impulse response. This digital filter represents the delay and amplitude of the sound arriving along different propagation paths. In contrast, the propagation delay of light can be ignored and only the energy steady-state response must be computed.

¯Coherence: sound is a coherent wave phenomenon, and interference between out-of-phase waves can be significant. Accordingly, acoustical simulations must consider phase when summing the cumulative contribution of many propagation paths to a receiver. More specifically, since the phase of the wave traveling along each propagation path is determined by the path length, acoustical models must compute accurate path lengths (up to a small percentage of the wavelength). In contrast, most light sources (except lasers) emit largely incoherent waves, and thus lighting simulations simply sum the power of different propagation paths.

¯Dynamic range: the human ear is sensitive to five orders of magnitude difference in sound amplitude and arrival time differences allow some high-order reflections to be audible [START_REF] Zwicker | Psychoacoustics: Facts and Models[END_REF][START_REF] Blauert | Spatial Hearing : The Psychophysics of Human Sound Localization[END_REF]. Therefore, as compared to computer graphics, acoustical simulations usually aim to compute several times more reflections, and the statistical time/frequency effects of late sound reverberation are much more significant than for global illumination. Despite these differences, the problems of simulation and validation are similar for both sound and light. For simulation, the main difficulty arises from the wavefield discontinuities caused by occlusions, caustics, and specular highlights, resulting in large variations over small portions of the integration domain (i.e., surfaces and/or directions). Due to these discontinuities, no general-purpose, analytic formula can describe the wavefield at a given point, and solutions must rely upon sampling or subdivision of the integration domain into components that can be solved efficiently and accurately. Traditionally, four approaches have been used to address this problem: finite or boundary element methods [START_REF] Goral | Modeling the interaction of light between diffuse surfaces[END_REF][START_REF] Moore | An Approach to the Analysis of Sound in Auditoria[END_REF][START_REF]Boundary Element Methods in Acoustics[END_REF], recursive ray tracing [START_REF] Krockstadt | Calculating the acoustical room response by the use of a ray tracing technique[END_REF][START_REF] Whitted | An improved illumination model for shaded display[END_REF], Monte Carlo path tracing [START_REF] Kajiya | The rendering equation[END_REF][START_REF] Naylor | Odeon -another hybrid room acoustical model[END_REF], and beam tracing [START_REF]CATT-Acoustic[END_REF][START_REF] Martin | Binaural simulation of concert halls : A new approach for the binaural reverberation process[END_REF][START_REF] Funkhouser | A beam tracing approach to acoustic modeling for interactive virtual environments[END_REF][START_REF] Heckbert | Beam tracing polygonal objects[END_REF]. All four methods have been used for both sound and light.

For validation, the common problems are: (1) defining a quantitative measure of accuracy, and (2) understanding the causes of simulation errors. The first problem is generally addressed in acoustics with statistical measures of impulse responses and with human listening tests. The second problem has largely been explained by conjecture. For example, recent "round robin" studies applied several computer simulation tools to the same concert hall model and compared their output to in-situ measurements [START_REF] Vorlander | International round robin on room acoustical computer simulations[END_REF][START_REF] Bork | A comparison of room simulation software -The 2nd round robin on room acoustical computer simulation[END_REF]. Although, these studies provide a nice reference for comparing simulations, they differ from our work in several respects. First, they focus on complicated concert halls, rather than simple rooms like ours, and thus it is difficult to evaluate the correctness of any of the simulations (i.e., the acoustics are so complex that none of the tools get the right answer). Second, they compare only gross statistical measures of impulse responses (e.g., T30, EDT, D50, C80, TS, G, LF, LFC, IACC) (ISO 3328), and thus it is difficult to determine the causes for simulation errors. Is it because propagation due to edge diffraction is ignored or because angle-dependent reflection functions were not modeled correctly? Without detailed examination of simple impulse responses, these questions are difficult to answer.

Our work is different than most previous validation studies in that we focus only on simple and controlled acoustical environments. In contrast to previous simulation systems, we not only use a simple 3D model of the environment, but the environment itself is very simple (e.g., a box with a few configurable panels). In this respect, our work is similar to recent detailed comparisons of geometrical simulations and measurements [START_REF] Suh | Measurment of transcient responses of rooms and comparison with geometrical acoustic models[END_REF][START_REF] Torres | Computation of edge diffraction for more accurate room acoustics auralization[END_REF]. These studies show that geometrical room acoustics can lead to satisfying simulations. But Suh and Nelson [START_REF] Suh | Measurment of transcient responses of rooms and comparison with geometrical acoustic models[END_REF] do not take into account sound diffraction, while Torres et al [START_REF] Torres | Computation of edge diffraction for more accurate room acoustics auralization[END_REF] do not use reverberant environments. In contrast, we present a simulation system and measurement environment which treats the combined effects of specular reflection and diffraction in a reverberant environment. As a result, we are generally able to simulate the sound propagation quite accurately, and we are able to evaluate computed impulse responses with great accuracy. For instance, we can often identify the sequence of surface reflections and edge diffractions that cause a particular peak to appear in an impulse response. Our goal is to obtain a detailed understanding of the limitations of the simulations in these simple environments, and then apply the acquired knowledge to improve our simulations with increased model complexity.

The details of our methods and validation studies appear in the following three sections. First, in Section 3, we explain the algorithms used in our computer simulation system. Then, in Section 4, we describe the construction of the Bell Labs Box and explain how we measure impulse responses. Third, in Section 5, we compare simulations with measurements and discuss our results. While the details of each step are non-trivial, the main message of the paper is that it is only possible to truly understand complicated processes in simple and controlled tests. Only when we have mastered the simple cases (e.g., a box with a diffracting panel) can we consider understanding more complex ones (e.g., a full concert hall).

Simulation System

Our sound simulation system is based on the beam tracing method described in [START_REF] Funkhouser | A beam tracing approach to acoustic modeling for interactive virtual environments[END_REF][START_REF] Funkhouser | Real-time acoustic modeling for distributed virtual environments[END_REF][START_REF] Tsingos | Modeling acoustics in virtual environments using the uniform theory of diffraction[END_REF]. One major advantage of this method over alternative approaches is that it finds all propagation paths with arbitrary combinations of transmission, specular reflection, and diffraction without spatial aliasing. This feature is particularly important for our studies with diffracting panels inside the Bell Labs Box.

The system takes as input: 1) a 3D environment described as a set of polygons with frequency-dependent impedances or bidirectional scattering filters, 2) a speaker described by its location and angular radiation pattern, 3) a microphone described by its location and optional angular directivity pattern. The system outputs a simulated impulse response.

The system executes in four steps, as shown in Figure 3. In the first step, we build a spatial subdivision data structure representing a binary space partition of 3D space into convex polyhedral cells. The purpose of this step is to decompose space into cells whose boundaries are aligned with polygons of the 3D input model and whose adjacencies are stored explicitly in a graph structure in order to enable efficient traversals of 3D space during beam tracing.

In the second step, we trace the convex polyhedral beams representing different propagation sequences through cells of the spatial subdivision. The traversal starts in the cell containing the speaker location with a beam representing the entire cell. Next, it visits adjacent cells iteratively, considering different permutations of transmissions, specular reflections and diffractions due to the faces and edges on the boundary of the "current" cell. As the algorithm traverses a cell boundary into a new cell, a copy of the current convex pyramidal beam is "clipped" to include the region of space passing through the convex polygonal boundary to model transmissions. At each reflecting cell boundary, a copy of the transmission beam is mirrored across the plane supporting the cell boundary to model specular reflections. At each diffracting wedge, a new beam is spawned whose source is the edge and whose extent includes all rays predicted by the Geometrical Theory of Diffraction [START_REF] Keller | Geometrical theory of diffraction[END_REF]. The traversal along any sequence terminates when either the length of the shortest path within the beam or the cumulative attenuation exceed some user-specified thresholds. The traversal may also be terminated when the total number of beams traced or the elapsed time exceed other thresholds.

In the third step, for each beam containing the microphone location, we compute the shortest propagation path from the speaker to the microphone along the sequence of transmissions, diffractions, and specular reflections represented by the beam (see Figure 4). The intersections with specularly reflecting faces are uniquely determined by the locations of the speaker, microphone, and the intersections with diffracting edges (diffraction points). The diffrac-tion points are found by solving a non-linear system of equations expressing equal angle constraints at diffracting edges. Once the diffraction points are found, we construct a piecewise-linear polyline representing the path along which sound travels from source to receiver along the propagation sequence, from which a length-, angle-, and frequency-dependent filter is computed. In the final step, for each valid propagation path from the speaker to the microphone, we add its contribution to the simulated impulse response. Our implementation includes source and material filtering effects derived from either measurements (see Section 4) or analytical models. We compute diffraction coefficients using the Uniform Theory of Diffraction [START_REF] Keller | Geometrical theory of diffraction[END_REF][START_REF] Kouyoumjian | A uniform geometrical theory of diffraction for an edge in a perfectly conducting surface[END_REF][START_REF] Kawai | Sound diffraction by a many sided barrier or pillar[END_REF] (see [START_REF] Tsingos | Modeling acoustics in virtual environments using the uniform theory of diffraction[END_REF] for details). Atmospheric scattering is also taken into account following the ISO 9013-1 specifications. All calculations are performed in complex Fourier domain at the sampling rate resolution. Our current system uses a sampling rate of 51200 Hz (the sampling rate of our source and material measurements). Thus, for a one second long response we are computing 25600 complex coefficients per path. The simulated transfer function is the sum of the coefficients for all paths, and the final impulse response is the inverse Fourier transform of the sum. After all beams have been traced up to a user-specified termination criterion and the contribution of all propagation paths have been summed, the resulting impulse response is output for comparison to measurements.

The Bell Labs Box is a small 2.19¢3.03¢2.42 meter enclosure (a volume of 16.058 meters ¿ ). It comprises a set of configurable panels, a sound source, a microphone and is constructed using standard residential housing techniques (see Figure 5). Sound measurements proceed in two phases. First, in an anechoic chamber, we measure the responses of the sound source, the microphone, and the surface panels independently. Second, we configure the panels into a room and make measurements of the sound propagation from the source to the microphone at different locations within the room that we later compare to simulation results (see Section 5).

Our sound source is a Brüel&Kjaer artificial mouth type 4227 speaker [4] (see Figure 6). During the first phase, we measure its directional responses by placing it on a rotator in the Bell Labs anechoic chamber [2], collecting responses at every 5 degrees in the azimuthal plane (this source has a revolution symmetry by design). Our listening device is a Brüel&Kjaer 4134 1/2-inch microphone, calibrated for free-field recording [4]. Based on the manufacturer's specifications, we assume its frequency response to be flat below 10KHz. We also assume the device to be perfectly omnidirectional below 10KHz. This defines the upper bound on the frequency for which measurement and simulations can be compared. For directional listening devices our system would measure their angle-dependent responses using a setup similar to the one used for sound sources.

Measurements of both the sound source and panel responses are acquired with a Siglab measurement system [START_REF] Siglab | spectral dynamics[END_REF] using repeated chirp stimuli. We use a 51200 Hz sampling rate for all measurements.

The surface panels of the Bell Labs Box are made of 3/4-inch melamine, smooth surface facing in, fastened to 3/4-inch FRP plywood panels. All seams are lapped and caulked. Before assembly, we measure their reflection responses using a bidirectional measurement rig in the anechoic chamber (see Figure 7). The rig has six degrees of freedom (three for both source and microphone). Thus, it allows for collecting bidirectional impulse responses of the scattering created by the panel. This setup can be used to measure the acoustic characteristics of any material and derive useful parameters, such as the complex acoustic impedance, which can be considered as an intrinsic material property for hard surfaces. Using the acoustic impedance, complex frequency dependent specular reflection coefficient are derived using the classic model for plane wave reflection which depends only on the impedance and incident angle [START_REF] Kuttruff | Room Acoustics[END_REF].

Panel Source Microphone

Figure 7: Measurement setup for bidirectional material properties in Bell Labs' anechoic chamber. In particular, our six degree of freedom rig allows measurement of the baffle reflection characteristics at different incident directions for source and microphone.

Validation Results

In order to evaluate our methods, we ran simulations and measurements for two different geometric panel configurations (see Figure 8). The first is a simple box-shaped enclosure comprising six rectangular panels. The second is the same box-shaped enclosure with a single rectangular panel spanning from floor to ceiling along one half the interior of the box with the speaker and microphone on opposite sides of the panel, as shown in Figure 2. The first configuration is very simple yielding only specular reflections, and we use it as a baseline for validation and comparison. The second configuration is a reverberant environment with diffraction and thus incorporates propagation paths combining both edge diffraction and specular reflection. This is a more difficult case for which detail validation results have not previously appeared in the literature.

For both configurations, we measured an impulse response using the Brüel& Kjaer artificial mouth and 1/2-inch microphone connected to the audio outputs and inputs of a MOTU 828 multichannel firewire audio interface [START_REF] Of | [END_REF]. The MOTU interface was connected to an off-the-shelf laptop running Windows. The source signal was a repeated chirp stimuli, and the sampling rate was 48KHz. The output signal used to feed the speaker was also fed back into the interface as a reference. We low-pass filtered the resulting response to get an actual bandwidth of 10KHz. We simulated all possible propagation paths combining up to 10th order specular reflection in the first configuration (empty box), and up to 4th order specular reflection and 2nd order edge diffraction in the second configuration (with diffracting baffle). These simulation parameters were chosen as a reasonable compromise between simulation accuracy and computational expense.

The simulated and measured impulse responses for each configuration are shown in Figures 9 and10. Note that fine band simulation with measured source and material characteristics allows us to compare the waveform of the simulated and measured early responses. As can be seen, our simulation is able to capture the temporal structure of the impulse response and is also able to capture the effects of the diffraction by the edge of the panel (e.g., first peak in Figure 10), an effect commonly ignored by most acoustic simulation systems. For comparison we also plotted a simulated impulse response without diffraction effects in Figure 10. Omitting diffraction prevents accurate modeling of the early part of the re-sponse, but has little influence on the late part which is dominated by specular reflection. A significant amount of energy is absent from the early response when reflections alone are simulated. This can introduce errors when evaluating perceptual criteria based on early to late energy ratios (e.g., clarity) [START_REF] Leo | Concert and Opera Halls: How They Sound[END_REF]. We also present in Fig- ure 11 a comparison between spectrograms for the simulated and measured impulse responses of Figures 9 and10. From these plots, we can conclude that our simulation performs well for the different frequencies. However, differences appear at low frequencies for which geometrical acoustics is no longer a valid model of sound propagation. Table 1 shows detailed simulation results for the test configuration with the diffracting baffle. Specifically, the columns list the path id, length, time delay, and sequence of scattering events, for the 60 shortest propagation paths (out of 1358 simulated). Every face and edge in the 3D model is assigned a unique identifier (see Figure 12) to allow detailed analysis of each simulated propagation sequence. For instance, looking in the fourth column of the second row, (s 4) denotes a specular reflection off surface 4 and (d 23) denotes a diffraction off edge 23. The ability to relate the propagation delay to the sequence of events along each path is a powerful tool for analyzing the impulse response since it makes it possible to derive correspondences between features in measured and simulated responses. For instance, we labeled some key features of the simulated response in Figure 10 with the corresponding path identifiers in Table 1. The directly diffracted contribution (the shortest possible path) is clearly identified (#1), followed by contributions of the back and side wall occluded by the panel (#2 and #3). Then, a strong specular reflection out of the other side-wall reaches the listener (#5). The strong reflected and diffracted contribution out of the back wall (path #2) is attenuated in our simulation. This might be due to slight inaccuracies in the measured source and microphone positions which results in interference of this path with the diffracted contribution reflected off the ceiling (path #3). This type of detailed analysis makes it possible to explain discrepancies between our simulation results and our measurements, and allows us to improve further our simulation models.

Conclusion

In this paper, we have described our initial experiences in building and using the Bell Labs Box, a simple enclosure for validating sound propagation simulations. We find that it is possible to achieve remarkably good matches between simulated and measured impulse responses in a carefully controlled experimental environment. Moreover, since the environment is so simple, we are not only able to validate the simulations at a gross level (as has been done before), but gain insight into the causes for simulation errors from detailed analysis of the differences between simulated and measured responses. We believe such detailed validation is necessary before the results of acoustics simulations can be understood for more complex environments.

Future Work

The study described in this paper is just the beginning of a long path toward understanding the validity of computer simulations of acoustical environments.

In the near future, we plan to extend our experiments to more complex environments. Specifically, we intend to consider a series of incremental steps of added complexity (e.g., inserting more surfaces, putting boxes in the enclosure, changing surface materials), measuring and validating the setup after each step. Our general approach is to validate simulations with gradually more complex For each path we assign a unique identifier and give the corresponding path length, propagation delay and sequence of events along the path.

environments so that we can understand and quantify the limitations of our simulations in detail. Eventually, we hope to validate simulations of concert halls and other more complex real world environments. Also, in the near-term, we plan to use our experimental setup to evaluate the accuracy of different models for reflection and diffraction in the Bell Labs Box. For instance, we could compare the Geometrical Theory of Diffraction [START_REF] Keller | Geometrical theory of diffraction[END_REF] with the Biot-Medwin-Tolstoy edge diffraction formulation [START_REF] Torres | Computation of edge diffraction for more accurate room acoustics auralization[END_REF], which may lead to more accurate simulations at greater computational expense. Or, we could evaluate more accurate reflection models, such as Thomasson's exact spherical wave reflection model off planar surfaces [START_REF] Thomasson | Reflection of waves from a point source by an impedance boundary[END_REF]. The results of these experiments may lead to an understanding of the trade-offs between computational expense and accuracy of different models.

In the longer-term, it will be important to consider psychoacoustical evaluation of acoustic models. We plan to build a system that will use validated acoustic simulations to investigate the psychoacoustic effects of varying acoustic modeling parameters. Our system will allow a user to interactively change acoustics parameters with real-time auralization and visualization feedback. With such a system, it may be possible to address psychoacoustic questions such as: "how many reflections are psychoacoustically important to model?," "which surface reflection model provides a psychoacoustically better approximation?," and "which among conflicting aural and visual cues are dominant in an interactive virtual environment?" We believe that the answers to such questions are of critical importance to future designers of 3D virtual environments.

Figure 1 :

 1 Figure 1: Left: a schematic view of the Bell Labs Box with measurement apparatus. Right: a computer graphics rendering of the Bell Labs Box lit by two point light sources.

Figure 2 :

 2 Figure 2: Additional panels can be mounted inside the Bell Labs box to study the effects of sound diffraction.

  (a) spatial subdivision (b) polyhedral beam tracing (c) path construction (d) impulse response generation

Figure 3 :

 3 Figure 3: The four phases of our simulation system.

Figure 4 :

 4 Figure 4: All possible 45 propagation paths combining two specular reflections and one diffraction around the edge of the panel in a model of the Bell Labs Box.

Figure 5 :

 5 Figure 5: (a) Outside view of the Bell Labs Box. (b) An inside view of the six-sided Bell Labs Box (without diffraction panel).

Figure 6 :

 6 Figure 6: Our reference source, the Brüel&Kjaer artificial mouth.

Figure 9 :Figure 10 :

 910 Figure 9: Comparison of a simulated early impulse response (top) including the first ten orders of specular reflection and a measured response (bottom) in the Bell Labs Box without the baffle. Our simulator computed the contribution of 1524 propagation paths in 738 seconds.

Figure 8 :

 8 Figure 8: Our two measurement configurations.

Figure 12 :

 12 Figure 12: In order to track the effect of different propagation sequences on the resulting impulse response, we assign a unique identifier to every edge (cyan) and face (yellow) in our model.

Figure 11 :

 11 Figure 11: Comparison between spectrograms of simulated and measured early responses. The two spectrograms on the left correspond to the (a) measured and (b) simulated impulse responses for the empty box (Figure 9), while those on the right correspond to the (c) measured and (d) simulated responses for the box with a baffle (Figure 10).

Table 1 :

 1 The 60 shortest propagation sequences (out of 1358) including four specular reflections and two diffractions off the edge of the panel in the Bell Labs Box with diffracting baffle.

	Id	Path length (m) Delay (ms) Propagation sequence
	1	1.6758	4.88	(d 23)
	2	2.7205	7.92	(s 4) (d 23)
	3	2.7416	7.98	(d 23) (s 7)
	4	3.0022	8.74	(s 5) (d 23)
	5	3.0154	8.78	(s 8)
	6	3.0603	8.91	(d 23) (s 8)
	7	3.1522	9.18	(d 23) (s 9)
	8	3.1530	9.18	(d 23) (s 10)
	9	3.3714	9.82	(s 3) (d 23)
	10	3.4798	10.13	(s 4) (s 6) (d 23)
	11	3.5356	10.30	(d 23) (s 3) (d 23)
	12	3.5356	10.30	(d 23) (s 8) (d 23)
	13	3.5645	10.38	(s 4) (s 5) (d 23)
	14	3.6860	10.73	(d 23) (s 2)
	15	3.7042	10.79	(s 5) (s 6) (d 23)
	16	3.7150	10.82	(s 7) (s 8)
	17	3.7515	10.92	(d 23) (s 7) (s 8)
	18	3.7515	10.92	(s 11) (s 4) (d 23)
	19	3.8117	11.10	(s 4) (s 1) (d 23)
	20	3.8274	11.15	(d 23) (s 7) (s 10)
	21	3.8637	11.25	(s 4) (s 3) (d 23)
	22	3.9156	11.40	(d 23) (s 4) (d 23)
	23	4.0093	11.68	(s 3) (s 6) (d 23)
	24	4.0176	11.70	(s 5) (s 1) (d 23)
	25	4.0275	11.73	(s 9) (s 8)
	26	4.0612	11.83	(d 23) (s 9) (s 8)
	27	4.1051	11.95	(s 4) (d 23) (s 8)
	28	4.1315	12.03	(d 23) (s 9) (s 10)
	29	4.1483	12.08	(d 23) (s 3) (s 6) (d 23)
	30	4.1483	12.08	(d 23) (s 8) (s 7) (d 23)
	31	4.1730	12.15	(s 4) (s 5) (s 6) (d 23)
	32	4.1956	12.22	(d 23) (s 10) (d 23)
	33	4.1956	12.22	(d 23) (s 5) (d 23)
	34	4.1978	12.22	(s 4) (d 23) (s 10)
	35	4.2772	12.46	(d 23) (s 7) (s 2)
	36	4.3005	12.52	(s 3) (s 1) (d 23)
	37	4.3338	12.62	(s 11) (s 4) (s 6) (d 23)
	38	4.3497	12.67	(s 11) (s 5) (s 4) (d 23)
	39	4.3709	12.73	(d 23) (s 8) (s 2)
	40	4.3868	12.77	(s 5) (d 23) (s 8)
	41	4.4304	12.90	(d 23) (s 3) (s 1) (d 23)
	42	4.4308	12.90	(d 23) (s 10) (s 2)
	43	4.4313	12.90	(s 4) (s 3) (s 6) (d 23)
	44	4.4535	12.97	(s 4) (s 5) (s 1) (d 23)
	45	4.4766	13.04	(d 23) (s 4) (s 6) (d 23)
	46	4.4795	13.04	(s 5) (d 23) (s 10)
	47	4.5513	13.25	(d 23) (s 9) (s 2)
	48	4.5804	13.34	(s 4) (d 23) (s 3) (d 23)
	49	4.5804	13.34	(s 4) (d 23) (s 8) (d 23)
	50	4.5851	13.35	(s 11) (s 4) (s 3) (d 23)
	51	4.5872	13.36	(d 23) (s 3) (s 4) (d 23)
	52	4.5872	13.36	(d 23) (s 4) (s 3) (d 23)
	53	4.6045	13.41	(s 11) (s 4) (s 1) (d 23)
	54	4.6433	13.52	(s 4) (d 23) (s 7) (s 8)
	55	4.6964	13.68	(s 4) (s 3) (s 1) (d 23)
	56	4.7235	13.75	(d 23) (s 5) (s 6) (d 23)
	57	4.7254	13.76	(s 4) (d 23) (s 7) (s 10)
	58	4.7308	13.78	(s 4) (d 23) (s 2)
	59	4.7392	13.80	(d 23) (s 4) (s 1) (d 23)
	60	4.7561	13.85	(s 3) (d 23) (s 8)

Acknowledgments

Thomas Funkhouser is partially funded by a National Science Foundation CAREER grant (CCR-0093343) and an Alfred P. Sloan Fellowship. The authors would like to thank the anonymous reviewers for their helpful comments.


Oliver Deussen

Carsten Colditz

Marc Stamminger

George Drettakis

Interactive Visualization of Complex Plant Ecosystems

Keywords: I.3.3 [Picture/Image Generation]: Display algorithms-[I.3.7]: Three-Dimensional Graphics and Realism-Animation Synthetic Plants, Ecosystems, Point-based rendering, Level-of-detail Algorithms

We present a method for interactive rendering of large outdoor scenes. Complex polygonal plant models and whole plant populations are represented by relatively small sets of point and line primitives. This enables us to show landscapes faithfully using only a limited percentage of primitives. In addition, a hierarchical data structure allows us to smoothly reduce the geometrical representation to any desired number of primitives. The scene is hierarchically divided into local portions of geometry to achieve large reduction factors for distant regions. Additionally, the data reduction is adapted to the visual importance of geometric objects. This allows us to maintain the visual fidelity of the representation while reducing most of the geometry drastically. With our system, we are able to interactively render very complex landscapes with good visual quality.

INTRODUCTION

Interactive rendering of large outdoor scenes is an important task for applications such as landscaping, architecture or computer games as well as virtual reality and simulation. The display of rich natural vegetation enhances comprehension of spatial relations in a planned landscape, to judge the visual impact of new buildings or to imagine the shape of a landscape in the future.

In recent years, a number of plant generation systems have been developed that allow the user to generate the geometry of single plants with high quality. Various aspects in designing plant shapes, growth of plants and their interaction have been investigated by Prusinkiewicz et al. [START_REF] Prusinkiewicz | The Algorithmic Beauty of Plants[END_REF] based on L-Systems. Other approaches use parameterised algorithms [START_REF] De Reffye | Plant models faithful to botanical structure and development[END_REF] or a set of components like the xfrog system [START_REF] Lintermann | Interactive modeling of plants[END_REF].

Based on these models, it is possible to generate large landscapes with many different plants. An open system architecture for such a system as well as algorithms for distributing and rendering complex scenes was presented by Deussen et al. [START_REF] Deussen | Realistic modeling and rendering of plant ecosystems[END_REF], while commercial products like Bryce (http://www.metacreations.com), Animatek World Builder (http://www.digi-elements.com) use plants to populate whole planets. The limiting factor for interactive rendering of scenes is the huge amount of geometry that has to be processed for each frame.

Vegetation geometry differs from other objects in a significant manner. While man-made objects such as houses, streets or cars are usually represented by a limited set of relatively large connected surfaces, vegetation usually consists of many small and isolated surfaces. Thus, most traditional level-of-detail schemes cannot be applied effectively to such scenes.

On the other hand, in recent years new rendering techniques have been proposed that allow efficient rendering of complex objects using the most simple primitives one can imagine: points. These methods seem to be appropriate for plant scenes because collections of discrete points can convey the shape of distant plants naturally.

The importance of interactive yet realistic rendering of these very complex ecosystem models cannot be overstated. For landscape design, or any kind of architectural or public works assessment, interactive viewing opens enormous potential for interactive experimentation and trial-and-error evaluation, which has previously been very hard to achieve for models of this scale and complexity. Such an interactive renderer is an invaluable tool in all architecture/design applications, as well as for virtual reality or video games.

In this paper, we present a system for interactively rendering lage outdoor scenes. Our contributions are the introduction of an appropriate method for pre-processing the model data, an efficient rendering algorithm using a point-and line-based level-of-detail approach and optimization of the scene using hierarchies.

For pre-processing we introduce a combined polygon/point and line representation; in contrast to previous methods we use lines in addition to points. Long and thin parts of plants such as leaves or branches are represented by lines, while more compact objects are approximated by points. This enables us to better convey the shape for both object categories in distant views.

Using our representation, the number of displayed vertices per second can be minimized keeping image quality high. Additionally, a new method of importance reduction allows us to represent visual important parts of a scene in a higher qualilty than others. Doing so, we can drastically reduce unimportant parts of the scene while keeping the visual impression stable. Although we focus on representing plants, the method can be combined with other levelof-detail approaches for different kinds of geometry to form a full purpose system for interactive display of complex scenes.

For rendering, we introduce a novel approach for blending through representations, which is both efficient and results in high quality visual results. Blending is done in a way that more and more polygons of a plant are represented by points and lines. Level-of-detail is obtained by reducing the displayed number of points and lines in correspondence to the viewing distance. Both polygonal data, points, as well as line approximations are stored in vertex arrays that allow us to efficiently display any desired portion of data.

We optimize the performance of our system using two types of hierarchies: Single plants such as trees that have high geometric complexity are spatially subdivided using an octree data structure. On the other hand, collections of many small plants are grouped locally and each group is entirely represented by a point and line representation. This enables us to display scenes of nearly arbitrary size efficiently.

After discussing related work we proceed by giving a system overview and describe how the data has to be pre-processed. The two forms of hierarchy are then introduced and rendering issues are presented. We conclude with results and some ideas of further work.

RELATED WORK

Previous approaches that are related to our work can be divided into two domains: level-of-detail representations and control for complex data and point-based rendering methods.

Level-of-detail (LOD) control: Among the many existing LOD methods one can distinguish between static methods that store an object by a small set of discrete representations that are blended into each other and dynamic methods that allow finer control of the geometry generated. In combination with outdoor scenes, the most simple way to form a discrete LOD is to use billboards that contain images of single plants [START_REF]Iris performer programmer's guide[END_REF] and replace the plant geometry for distant views. If the viewpoint changes, no parallax can be seen which makes the scenes flat and unrealistic. Due to its simplicity, this method is currently used in games and VR applications. Improvements are sets of impostors [START_REF] Schaufler | A three dimensional image cache for virtual reality[END_REF] and sprites with depth [START_REF] Shade | Layered depth images[END_REF] that offer a limited amount of parallax.

Max et al. [START_REF] Max | Rendering trees from precomputed Z-buffer views[END_REF][START_REF] Max | Hierarchical imagebased rendering using texture mapping hardware[END_REF] hierarchically replace parts of trees by precomputed views and is able to smoothly blend between these. Unfortunately, the proposed technique is too inefficient for interactive display. Perbert and Cani [START_REF] Perbert | Animating prairies in real-time[END_REF] explore an efficient rendering and animation method for meadows based on impostors.

Dynamic LOD representations for smooth objects were presented by Hoppe et al. [START_REF] Hoppe | Progressive meshes[END_REF]. Objects are represented by a base geometry and a set of vertex split operations that refine the triangle mesh. Methods applied to terrain rendering were presented by Duchaineau et al. [START_REF] Duchaineau | ROAMing terrain: Real-time optimally adapting meshes[END_REF] and by Hoppe [START_REF] Hoppe | View-dependent refinement of progressive meshes[END_REF]. Unfortunately, foliage cannot be processed this way. The discrete structure of the leaves allows edge collapses only to a very limited extend. Nevertheless, Rossignac and Borrel reduced geometry by collapsing nearby vertices regardless of the connection structure [START_REF] Rossignac | Multi-resolution 3D approximations for rendering complex scenes[END_REF]. This can be applied to the foliage but changes the visual appearance drastically.

Point rendering:

The first general idea of using points as display primitives was presented by Levoy and Whitted [START_REF] Levoy | The use of points as display primitives[END_REF]. In the same year, Reeves and Blau [START_REF] Reeves | Approximate and probabilistic algorithms for shading and rendering structured particle systems[END_REF] rendered complex trees using a set of small disks representing the foliage. Weber and Penn presented a level-of-detail representation of trees using sets of points for the leaves and lines for the tree skeleton [START_REF] Weber | Creation and rendering of realistic trees[END_REF]. As mentioned above, we extend this idea to model whole plant populations. Shade et al. [START_REF] Shade | Layered depth images[END_REF] sampled geometric plant models using a modified ray tracing procedure. Based on the resulting point sets they were able to render images from various view points. This works well for distant views but results in holes and a grainy look for close-ups.

An LOD method based on points is presented by Rusinkiewicz and Levoy [START_REF] Rusinkiewicz | Qsplat: A mulitresolution point rendering system for large meshes[END_REF]. A hierarchy of point representations allows smooth reduction of the amount of geometry for smooth surfaces. Pfister et al. [START_REF] Pfister | Surfels: Surface elements as rendering primitives[END_REF] present sophisticated rendering methods for point represen-tations. Their method works in software and therefore is too slow for interactive rendering highly complex scenes.

Cohen et al. [START_REF] Cohen | Hybrid simplification: Combining multi-resolution polygon and point rendering[END_REF] combine a LOD representation for polygonal objects with point rendering using a multi-resolution graph which describes simplification operations that lead to a point representation. By displaying the prefix of this graph, geometry reduction is performed. In this manner, they are able to smoothly blend between geometry and point representations. A similar transition method is presented by Chen and Nguyen [START_REF] Chen | Pop: A hybrid point and polygon rendering system for large data[END_REF]. However, for our highly complex objects an explicit representation of refinement operations seems to be too inefficient. Given a single tree with half a million of leaves (see figure 7(b)) the resulting graph would need a large amount of memory. Instead, we use vertex arrays that store the geometry of plants in an appropriate manner that allows us to reduce geometry by showing only the first part of such an array. [START_REF] Stamminger | Interactive sampling and rendering for complex and procedural geometry[END_REF] present a point based method for plants using random sample positions on the surface. Each plant is described by a set of points obtained from these positions. A similar approach is described in [START_REF] Wand | The randomized z-buffer algorithm: Interactive rendering of highly complex scenes[END_REF], where a plant is rendered by a random set of points, which is just dense enough not to exhibit holes. As the authors determine the point positions individually for each frame, aliasing is introduced.

Stamminger and Drettakis

SYSTEM OVERVIEW

Plants for our system are generated from surface oriented plant descriptions produced by the plant modellers described previously. In a pre-processing step these models are converted into point and line sets. Both the polygonal description and its point or line approximation are stored in a data file for each plant. A separate editor, similar to the one described in [START_REF] Deussen | Realistic modeling and rendering of plant ecosystems[END_REF], is used to combine the plants to generate complex scenes. This involves model quantization to reduce the number of different plant geometries. For a single plant population in most cases a set of five to ten different plant models is sufficient, the whole landscape may require several dozens.

Our scenes are described by a set of plant models and so called eco files that store positions of plants and instancing information. Instances of the plant models are spread over the terrain to form the plant population. Eco files can be used hierarchically, i.e., one eco file may by a combination of reference to other eco files.

The rendering system reads in eco files and plant geometries. For the entire content of each eco file a point approximation is determined. In distant views only this very sparse point representation is used. When zooming in, all the plants of the eco file are displayed. Visibility culling on the basis of bounding boxes is used to determine the visible plants. The system decides whether a polygonal representation will be used or its approximation as a function of the projected size of a plant. Figure 1 displays the system outline.

In the following, level-of-detail control is done with respect to the vertex count. This count describes the number of vertices that max-imally can be processed by the graphics processor per second. In most cases, our scenes consist of many small triangles. The pixel fill rate of the graphics processor does not play a significant role in this case and the performance is thus fully determined by the vertex count.

PRE-PROCESSING OF MODELS

During modelling of a plant, each plant is divided into elements that are approximated by points and others that are approximated by lines in distant views. The user decides which parts of a plant model are represented by which primitive. This is only done during modelling the plants in an early stage. If several instances of a plant type with slightly varying geometry are used in a scene, modelling, with its corresponding effort, is only performed once and therefore is negligible. Another attribute that has to be set by the user is the visual importance of the plant parts that is later used for showing visually important plant parts with higher precision than others.

Point Representation

Point representations are typically used for leaves of a tree, petals etc. By default, for an object of n triangles, 2n points are generated, reflecting the fact that a separate triangle requires about as much time for rendering as three points and therefore an approximation by two points needs less rendering time. The 2n points are distributed randomly over the object, in no particular order, but in a way that equal surfaces receive equal numbers of points. This means that small triangles might not be attributed any points, whereas big triangles can be represented by many. All the random points are stored in a single list. During rendering, only a prefix of this list with length depending on the projected size of the plant is drawn [START_REF] Stamminger | Interactive sampling and rendering for complex and procedural geometry[END_REF] (see Fig. 2). A different sampling technique is applied if all triangles of the object have about the same area. In this case, we generate two sample points on each of the triangles. This corresponds to stratified sampling and leads to more even sampling. Furthermore, the triangles and points are randomly reordered in the same order, so we still know the point-triangle correspondence. This allow us to render part of the object by triangles, and the rest by points. In this manner, we can blend between polygonal and point representations smoothly, largely avoiding popping artefacts, as we will describe later.

Line Representations

Thin and long structures like twigs or several types of leaves are better represented by a line than by a set of points [START_REF] Weber | Creation and rendering of realistic trees[END_REF]. We generate a line set for these object parts already in the modelling phase, when all information about the object structure is available. The line set is also randomly reordered, and stored in a list. In the same manner as for the points, for each view only a prefix of this list is rendered, according to the current viewing distance. We will discuss this issue further in the rendering section.

Importance Reduction

Best visual results can be achieved if the plants are not uniformly reduced. Visual important parts like petals or other plant elements with special colors should be reduced more slowly than than the rest of the plant. In Figure 4 an example is shown. The daisies in the meadow are visually important even though their geometry is small in comparison to the rest. If all geometry is reduced uniformly, the reduction is more noticeable, and the perceived quality much worse than if these elements are handled separately.

Currently, an importance factor is set interactively by the user. In the future this might be automated by color analysis of the geometric primitives: in most cases objects seem to be important if their color is clearly distinguishable from the majority of element colors in a plant.

Memory Requirements

The additional point and line data roughly doubles the memory requirements for storing our models. This is why we will try to fully automate the entire pre-processing step. In this case the plant models can be stored as procedural plant descriptions in the form of L-Systems or xfrog-files with very low memory needs. The extraction and approximation process then has to be done for each run in a pre-processing step. However, the memory needs for our plant models should are not very high, due to the approximate instancing mechanism in which a complex scene needs about several dozens of different models requiring usually less than a hundred megabytes of memory. For the scenes shown on the accompanying video, the entire geometric data was compact enough to be stored in AGP memory, which increases the rendering throughput by a factor of almost two.

RENDERING PLANT MODELS

After eco files and plant models are input into the rendering system, the geometry (polygons, points and lines) is stored in vertex arrays. During display, for each part of a plant its representation-triangles or points and lines-and the number of displayed primitives is determined dynamically.

Before rendering an object, visibility culling is applied. First, the culling is performed on each of the eco files, then on the plant model instances. If eco files are organized in a hierarchy, the culling is also applied hierarchically.

In the following, we assume that the original plant model consists of n triangles T1, .., Tn with area A in world space. Some of these triangles are approximated by n l lines L1, .., Ln l with overall length l, others are approximated by np points P1, .., Pn p . We now split the set of triangles T1, .., Tn of the model description into two sets, one set TP = T1, .., T k that is approximated by points and the set TL = T k+1 , .., Tn approximated by lines. The set TL has an area A l , the set TP has an area Ap in world space. Clearly we have:

A = Ap + A l .

Rendering Points

In the rendering stage, point and line representations are handled separately: For the point approximation the number of points is calculated that is required to render the current object faithfully, i.e. without holes and with correct coverage. Here each point of the approximation is represented by a small subimage, a so called splat, usually one pixel in diameter.

The user controls the speed versus quality tradeoff by setting a sample distance parameter d , that defines the average distance between two neighboring point samples on the image plane. A large value of d results in a larger distance and thus fewer samples and vice versa. Consequently, the required point splatting area is A sp = d 2 (in practice, we will use slightly bigger splat areas to avoid holes; see below).

The number of points that is needed for representing the plant faithfully depends on the surface Ap, as well as its average distance to the virtual camera r. In order to avoid scaling factors, we define the image plane to be at distance one from the camera. The approximate projected area A p of the triangle set TP on the image plane and the required point number p can be calculated by (see also [START_REF] Stamminger | Interactive sampling and rendering for complex and procedural geometry[END_REF]):

A p = 1 2 Ap r 2 p = cp A p A sp np (1)
In this formula it is assumed that the triangles are double sided, i.e., a leaf is modelled by a single surface and not two separate surfaces for front and back. The factor 1/2 accounts for leaf orientation under the assumption that all leaves are randomly distributed. 1/r 2 corresponds to perspective foreshortening. The splat factor cp allows us to increase the splat size in order to avoid gaps that appear due to the random sampling. Values from 1.2 to 1.5 are usually sufficient. If the number of required points p > np the model part is rendered polygonally, otherwise by its point approximation using p points.

Rendering Lines

For the plant part that is to be approximated by lines the approach is slightly different. Here we have two requirements for the approximation: the projected area needed to represent the object faithfully and the geometric approximation error of the lines. This is due to the fact that two line sets with equal projected area can approximate the object differently.

Let the line set L1, .., Ln l approximate the triangle set TL with a maximum distance l in world space. In image space the value is compared with the sampling distance. If l /r ≥ d we render the set by triangles, otherwise -as in the point case -we compute the image area if all lines were drawn. Doing so we need the projected area A l of the polygonal representation of corresponding triangle set T L. This is determined by projecting it similarly to Eq. ( 1). If l is the length of all lines in world space, the lines will have length l/2r in the image plane, where again 1/2 accounts for the (random) orientation, and 1/r for the perspective foreshortening (1/r instead of 1/r 2 , because we project length and not area!). The image plane area covered by the lines drawn with line width d is thus

A sl = l d = ld 2r , (2)
We are now able to determine the ratio between projected area of T L and its line approximation q l = A l /A sl . If q l ≤ 1 we draw q l n l lines, otherwise the model part is drawn using polygons. Note that this assumes, that all line segments have about the same length, which is the case for the lines generated by our modelling tool.

Blending Representations

In the above section we described how and when to switch between a polygonal description and its approximation by points or lines. Depending on the selected sample distance d , this switch in representations may become visible (popping). However, this effect can easily be avoided by not switching between point and triangle representations immediately, but instead continuously replacing more and more triangles by points or vice versa.

As described above, for certain objects, we determined two random points on each triangle for the point representation. The point list is scrambled in the same way as the triangles, so that the first 2n points always correspond to the first n triangles. The selected level of detail depends on p, the number of requested points (see Eq. ( 1) and also Fig. 6).

Figure 6: Rendering a plant in mixed polygon/line/point representations

If p is larger than three times the number of triangles k, that form the part of the plant that is to be approximated by points, the full triangle set is rendered. On the other hand, if p < 2k, we render the first p points of the point array. For 2k ≤ p < 3k, we blend between both representations by rendering the first p-2k triangles, and the points corresponding to the remaining triangles. Due to the ordering enforced during the pre-process, these points can be rendered efficiently as vertex arrays.

If we count the number of rendered vertices (number of rendered points plus three times the number of triangles), we see that p vertices are always rendered, but never more than 3k, so the rendering time is roughly linear in the number of pixels covered by the object, and essentially independent of object complexity.

The blending between polygons and their approximation works particulary well for the isolated surfaces which constitute plants. In Figure 5(d)-(f) the representation change is presented. While subfigure (d) shows the polygonal model, in (e) a mixed representation is given and in (f) a point approximation.

Implementation Issues

Image display can be altered by two parameters: the sampling distance d and the splat factor cp in Eq. 1. Fig. 5(a)-(c) demonstrates the influence of both parameters. If the sampling distance is increased, elements are converted into points earlier, resulting in fewer, but larger points. If the splat size is enlarged, resolution is decreased but the larger area of the splats allows a stronger reduction.

A noteworthy side effect of calculating the representations in the described way is that the sampling distance d can be varied dynamically to ensure a desired vertex complexity. This enables animations with nearly constant frame rates.

Another nice aspect of our rendering method is that due to the fixed size of the display elements a point representation often causes fewer aliasing artefacts than a polygonal representation with very small polygons. In this sense, our point representation causes aliasing due to its discrete structure, but due to the stable display of points in animations the point approximations of plants are visually more stable than their polygonal representation.

OPTIMIZATION

As mentioned above, in our system vegetation is represented by plant models and eco files with instancing information. In order to deal with the huge complexity of nature, we introduced hierarchical data structures in two ways: for the efficient rendering of complex individual plants and for the representation of very large numbers of small plants.

Since our rendering system is tightly coupled to the modeller, we directly work on its output, and need to get fast visual feedback.

Resorting the objects, building completely new scene hierarchies, or prefiltering of the point sets is thus not possible. We essentially have to use the structure of the scene as provided by the modeller, and enhance it with lightweight hierarchical information. In the following we describe our solution, which adds hierarchical information even to complex scenes in a few seconds.

Representation of Complex Plants

Eqs. 1 and 2 provide a global mechanism for representing a plant geometry polygonally or by points and lines. This is in contrast to the local schemes presented in [START_REF] Cohen | Hybrid simplification: Combining multi-resolution polygon and point rendering[END_REF][START_REF] Chen | Pop: A hybrid point and polygon rendering system for large data[END_REF]. The advantage of the global decision is that it is very efficient. Nonetheless, this approach has a drawback for large and complex plant models: If the viewer is close to a single leaf of a complex tree, due to Eq. ( 1) all the leaves are displayed polygonally because of the large projected size of the nearby leaves. This is not optimal for the distant leaves of the same plant that might be better represented by points or lines.

To reduce this effect, large and complex objects are spatially subdivided by an octree of small depth. All the triangles in a cell are represented by a polygon set and a point and line set. If the camera comes closer, the content of nearby cells is shown polygonally while distant cells can be approximated. The maple in Figure 7(b) is a typical example of a tree that benefits from such a structure.

Grouping of Plants

Rendering large scenes requires the storage and display of millions of plants. In this case, the iteration over all plants tends to be the limiting factor even if each plant is displayed by only a single point or line.

To deal with plant models of such size, vast plant populations are stored in a set of quadratic tiles that fit together. This results in a local grouping of plants. If a single tile has too many instances, it is further subdivided into sub-tiles, until a given number of instances per tile is achieved.

Each tile is represented by an eco file. For each eco file a rough point representation for the geometry of the entire tile is computed by randomly collecting points from the geometry approximations of the plant instances. Additionally, eco files can be organized hierarchically i.e., an eco file contains several eco files with positional information. If this is the case, the rough point approximations of the children are merged to form a point approximation of the overall content of the eco file. Doing so, we are able to display large areas.

In contrast to a quad-tree data structure this model-oriented spatial data structure gives us more freedom to organize different plant populations by several eco files that, for instance, obtain different importance parameters. This is why we prefer this data organization. Nevertheless, it is relatively easy to use a quad-tree data structure for the same purpose.

Caching

For performance reasons it is very important to load as much data as possible into memory which is quickly accessible by the graphics processor. Free AGP memory is usually quite limited, e.g., about 18MB on our GeForce3. Therefore, we determine how often a plant model is used in the current scene. The geometric data of frequently used models is primarily loaded into AGP memory. The previously described instancing scheme that limits the number of different models in our scenes helps in this respect by minimizing the amount of geometry needed. We are thus able to achieve a vertex rate of about eleven million vertices per second for rendering our scenes on a GeForce3.

Shadow Computation

Shadows are very important for correct light interaction in outdoor scenes. In combination with our interactive system we implemented two shadow map mechanisms. The first is to obtain a standard shadow map in a pre-processing step. This can be used efficiently during rendering but does not allow any dynamic objects and results in pixelisation artifacts for the large scenes we want to show.

The second possibility is to use the perspective shadow map proposed by Stamminger and Drettakis recently [START_REF] Stamminger | Perspective shadow maps[END_REF]. In this case, the map is obtained for each image in a way that the pixel resolution is concentrated to areas where it is needed. This is done by computing the map after the perspective projection of the scene. The approach allows us to use a relatively small map without visible artefacts. On the other hand, the map has to be obtained for each frame, which reduces performance. When computing a conventional high resolution shadow map the scene has to be used without any level-of-detail approximation -a procedure that might last for some seconds-the perspective shadow map can be generated in level-of-detail mode seen from the camera. The frame rate is thus about 50% of the frame rate without shadow computation.

RESULTS AND FUTURE WORK

We have applied the method to a number of plant models and eco systems. With modern graphics hardware based on PC's -we used a 1,2 GHz Pentium with NVidia Geforce3 graphics processor-over eleven million vertices per second can be drawn. The vertex rate is enough for interactively displaying a number of fairly complex plant scenes as can be seen in the accompanying video.

In Figure 7(a) an outdoor scene is shown, that consists of over 120 million polygons. The original scene generated by Deussen et al. [START_REF] Deussen | Realistic modeling and rendering of plant ecosystems[END_REF] is a part of what we use in our system. Their conventional renderer needed about 75 minutes for image generation of 30 million triangles computing 9 samples per pixel. In our system we achieve eight to ten frames per second for the much larger scene albeit without anti-aliasing. Future work will include the automatic texture synthesis for distant eco files and a simple approach for simulating wind. Using the ingredients described, a variety of users from fields like such architecture, landscaping, simulators, and game design can use scenes with rich vegetation for a real virtual reality.

Figure 1 :

 1 Figure1: System structure, following the outline given in[START_REF] Deussen | Realistic modeling and rendering of plant ecosystems[END_REF]

Figure 2 :

 2 Figure 2: Point representation of a pine. Polygonal model and representation by 13,000, 6,500, 3,250, and 1,625 points. Low sampling density has been selected to visualize the point models.

Figure 3 :

 3 Figure 3: Representing a small palm by lines. Upper row: polygonal model and representation with 106 lines. Lower row: 1,500 plants distributed on a plane, original scene size is 12M triangles.

Figure 4 :

 4 Figure 4: Importance reduction: a) uniform reduction; b) importance reduction, the visually important daisies are represented polygonally.

Figure 5 :

 5 Figure 5: Left column: Influence of sampling distance d and splat factor cp: a) Sampling distance of about one pixel. Plants in the front are drawn by polygons, those in the back by points. b) With a larger sampling distance, more plants are rendered by points. An insufficient splat factor leads to holes in the rendering. c) Same sample distance as in b), but a larger splat factor leads to bigger splats, such that holes are reduced. Right column: Mixed representation of plants. d) a tree rendered by triangles, e) a mix of triangles and points, f) points only.

Figure 7 (Figure 7 :

 77 Figure 7(b) shows another scene with a large tree and 13 thousand sunflowers, 70 million triangles in total. This scene is displayed typically at three to four frames per second. The relatively low performance stems from the models that allow no strict reduction

Table 1 :

 1 Table 1 some values are given. Complexity of different scenes and count of elements in average during interactive rendering. Typical frame rates; if camera is very close smaller rates can result.

		Size	LOD	LOD	LOD	Frame
	Model (triang) Triang. Lines Points	rate
			(av.)	(av.)	(av.)	(av.)
	640 trees
	as in Fig. 2	10M	500K	100K 400K	9-16 Hz

ACKNOWLEDGEMENTS

Marc Stamminger has been supported by a Marie-Curie postdoctoral fellowship while doing this work.


Thomas Funkhouser

; Tsingos

; Carlbom

; Elko

; Sondhi

; West

Modeling Sound Reflection and Diffraction in Architectural Environments with Beam Tracing

Keywords: Reference: 43, 55, Ka

This paper describes and analyzes a beam tracing method for computing sound propagation paths from sound sources to receivers in architectural environments. The algorithm traces polyhedral beams from the location of each sound source through a precomputed data structure encoding spatial adjacencies. The main features of the method are: 1) it enumerates all potential sequences of specular reflection, diffraction, and transmission up to user-specified termination criteria, 2) it enables evaluation of early propagation paths at interactive rates, and 3) it scales to support large, "densely-occluded" architectural environments. It is being used to support real-time auralization in immersive virtual environment applications.

INTRODUCTION

Geometric acoustic modeling tools are commonly used for design and simulation of 3D architectural environments. For example, acousticians use CAD tools to evaluate the acoustic properties of proposed auditorium designs, factory planners predict the sound levels at different positions on factory floors, and video game designers provide spatialized sound effects in interactive virtual environment systems.

One major challenge in geometric acoustic modeling is accurate and efficient computation of early sound propagation paths. As sound travels from source to receiver via a multitude of paths containing reflections, transmissions, and diffractions, accurate simulation is extremely compute intensive. Most prior systems for computing early reflections have been based on image source methods [START_REF] Allen | Image Method for Efficiently Simulating Small-Room Acoustics[END_REF][START_REF] Borish | Extension of the Image Model to Arbitrary Polyhedra[END_REF] and/or ray tracing [START_REF] Krockstadt | Calculating the Acoustical Room Response by the Use of a Ray Tracing Technique[END_REF], and therefore they do not generally scale well to support large 3D environments, and/or they fail to find many significant propagation paths (e.g., ones due to diffractions).

In this paper, we describe a beam tracing method that computes early propagation paths in 3D models of architectural environments. An important feature of this method is that it finds all paths combining transmission, specular reflection, and edge diffraction up to user-specified termination criteria (e.g., all paths arriving within a certain time delay) without sampling or exhaustive enumeration. As a result, it ensures that no propagation paths are missed due to undersampling, and it executes efficiently for large 3D polygonal models. We currently use this method to simulate simple architectural spaces and to spatialize sounds in interactive virtual environment applications.

The remainder of the paper is organized as follows. Section 2 briefly reviews related work in geometric acoustic modeling. Section 3 contains a description of our beam tracing method. Sections 4 and 5 present experimental results. Finally Section 6 contains a brief conclusion, along with a discussion of applications, limitations, and topics for future work.

BACKGROUND AND RELATED WORK

The problem of modeling sound propagation can be solved using the Helmoltz-Kirchoff integral equation, which expresses the wave field at every point in space in terms of the wave field at other points (or equivalently on surrounding surfaces). The main challenge in solving this equation is that the wave field has discontinuities due to occlusions, diffractions, and specular highlights, which induce large variations over small portions of the integration domain (i.e. surfaces and/or directions). Due to these discontinuities, no analytic formula can be applied to compute the wave field at a given point in the general case, and solution methods must rely upon sampling or subdivision of the integration domain into separate parts that can be treated efficiently and accurately. Traditionally, four approaches have been used to address this problem: boundary element methods, Monte Carlo path tracing, ray tracing, and beam tracing [START_REF] Kleiner | Auralization -An Overview[END_REF][START_REF] Kuttruff | Room Acoustics, 3rd Edition[END_REF].

Boundary element methods find an approximate solution to the wave equation by solving the system of equations resulting from discretizing the surfaces into patches [START_REF] Moore | An Approach to the Analysis of Sound in Auditoria[END_REF]. The major problem with this method is that the surface mesh must be fine enough to account for phase differences. Also, it is difficult to adapt the mesh of each surface to capture the irregularities and discontinuities of the wave field, and it is expensive to compute extended form factors, solve threepoint transport problems, or store directional basis functions for scenes with non-diffuse surfaces. These problems make boundary element methods impractical for sound simulations at all but the lowest audible frequencies.

Monte Carlo path tracing methods [START_REF] Kajiya | The Rendering Equation[END_REF] solve the integral equation by sampling randomly generated paths from the source to the receiver. In order to avoid large errors (variance) in the computed solutions, the sampling density must fit the shape of the integrand, which is difficult to achieve due to discontinuities in the wave field. Although methods have been proposed to achieve robust and efficient integration, it is still difficult to find the most significant paths, and a large number of samples are needed to render solutions with sufficient fidelity. Since each sample requires a separate set of ray-scene intersection tests, and many of the samples are rejected due to occlusions, these methods generally require long compute times, and thus they are not suitable for interactive applications.

Ray tracing methods [START_REF] Krockstadt | Calculating the Acoustical Room Response by the Use of a Ray Tracing Technique[END_REF] find propagation paths between a source and receiver by generating rays emanating from the source (or receiver) position and following them individually as they propagate through the environment. Although this method is very general and simple to implement, it is subject to aliasing artifacts as the space of rays is sampled discretely [START_REF] Lehnert | Systematic errors of the ray-tracing algorithm[END_REF]. For instance, receiver position and diffracting edges are often approximated by volumes of space (in order to admit intersections with infinitely thin rays), which can lead to false hits and paths counted multiple times. More often, important propagation paths may be missed by all samples. Moreover, ray tracing is very compute intensive, usually taking minutes to hours to compute a receiver-dependent solution.

Image source methods [START_REF] Allen | Image Method for Efficiently Simulating Small-Room Acoustics[END_REF][START_REF] Borish | Extension of the Image Model to Arbitrary Polyhedra[END_REF] compute specular reflection paths by considering virtual sources generated by mirroring the location of the audio source over each polygonal surface of the environment. The key idea is that a direct path from each virtual source has the same directionality and length as a specular reflection path. Thus, specular reflection paths can be modeled up to any order by recursive generation of virtual sources. This method is simple for rectangular rooms [START_REF] Allen | Image Method for Efficiently Simulating Small-Room Acoustics[END_REF]. However, in general environments, N R virtual sources must be generated for R reflections over N polygons. Moreover, for every new receiver location, each of the O(N R ) virtual sources must be checked to see if it is visible to the receiver, since the specular reflection path might be blocked by a polygon or intersect a mirroring plane outside a polygon [START_REF] Borish | Extension of the Image Model to Arbitrary Polyhedra[END_REF]. As a result, this method is practical only for computing very few specular reflections from stationary sources in simple environments.

Beam tracing methods [START_REF] Heckbert | Beam Tracing Polygonal Objects[END_REF][START_REF] Dadoun | The Geometry of Beam Tracing[END_REF] find propagation paths from a source by tracing beams (i.e., bundles of rays) through a 3D polyhedral environment. In general, a set of beams is constructed that completely covers the space of rays from the source. For each beam, polygons are considered for intersection in order from front to back. As each intersecting polygon P is detected, the original beam is clipped to remove the shadow region. Additionally, new beams may be created representing various types of surface scattering. For instance, a transmission beam may be constructed matching the shadow region, a specular reflection beam may be constructed by mirroring the transmission beam over the plane supporting P, and a diffraction beam may be constructed to represent scattering from a wedge. These new beams are then traced through the environment until some application-specific termination criteria are met (e.g., a maximum propagation delay). The main advantage of this approach is that it finds all propagation paths up to the termination criteria. The disadvantage is that the geometric operations required for beam tracing are more complex than for the individual paths.

BEAM TRACING APPROACH

Our approach is a hybrid between beam tracing and path tracing. The general strategy is to trace beams that partition the space of rays into topologically distinct bundles corresponding to different sequences of scattering events at surfaces of the 3D scene (propagation sequences), and then use them to guide sampling in a path tracing algorithm. Conceptually, the approach proceeds in two phases. In the first phase, beams are traced to enumerate potential propagation sequences comprising different permutations of diffuse reflections, specular reflections, diffractions, and transmissions. The result of this phase is a set of conservative beams, each containing the bundle of rays potentially traveling along a different propagation sequence. In the second phase, the wave field is evaluated at a point by solution methods that use traced beams to identify subdomains for analytic integration and for structured sampling in a path tracing algorithm. The motivation for this approach is that attributes of a relatively small number of beams traced during the first phase can provide useful information about the wave field that can be used to guide and accelerate evaluation of samples during the second phase. Some advantages of this approach are as follows:

• Fast Enumeration of propagation sequences: beam tracing provides an efficient method for finding potential sequences of surface scattering events, without sampling as in ray tracing or path tracing, and without exhaustive search of all surface permutations, as in image source methods. In particular, it provides a deterministic algorithm to enumerate every propagation path up to a specified termination criteria without risk of undersampling errors (unlike ray tracing), and it scales to support large architectural environments with many surfaces (unlike image source methods).

• Spatial coherence in ray intersections: a priori beam tracing can improve the efficiency of multiple ray intersection tests. Due to coherence in ray space, each beam-surface intersection represents an infinite number of ray-surface intersections. Thus, once a beam has been traced along a certain sequence of surface intersections, generating ray paths following the same sequence becomes a sampling problem rather than a combinatorial problem, and the expensive computation of casting rays through a scene can be avoided for individual samples.

• Progressive refinement: characteristics of beams can be used to guide progressive refinement and importance sampling strategies. For instance, bounded estimates of the path lengths or the energy carried by different beams can be used to prioritize beam tracing steps [START_REF] Funkhouser | Real-time Acoustic Modeling for Distributed Virtual Environments[END_REF][START_REF] Min | Priority-Driven Acoustic Modeling for Virtual Environments[END_REF] and/or to determine the number of samples to evaluate for each propagation sequence.

• Interactive path generation: the result of a beam tracing computation can be computed offline and stored in a data structure (e.g., a beam tree) for later evaluation of propagation paths at interactive rates. This feature enables visualization and/or auralization of sound traveling from a stationary source to an arbitrarily moving receiver location (or vice versa) in an interactive application [START_REF] Funkhouser | A Beam Tracing Approach to Acoustic Modeling for Interactive Virtual Environments[END_REF].

The challenge of implementing a beam tracing system is to perform polygon sorting, polygon intersections, and beam transformations efficiently and robustly. Although several methods have been proposed to accelerate these geometric operations, including ones based on BSP-trees [START_REF] Dadoun | The Geometry of Beam Tracing[END_REF], cell adjacency graphs [START_REF] Fortune | Topological Beam Tracing[END_REF][START_REF] Funkhouser | A Beam Tracing Approach to Acoustic Modeling for Interactive Virtual Environments[END_REF][START_REF] Teller | Visibility Computations in Densely Occuded Polyhedral Environments[END_REF], and layers of 2D triangulations, there is still work to be done for general scattering phenomena (e.g., diffraction and diffuse reflection) and general 3D scenes (e.g., environments with fine detail).

This paper outlines our recent work on using beam tracing for modeling early propagation paths combining specular reflection and edge diffraction in interior architectural environments. The next section briefly describes our methods, while the following two sections present results of experiments evaluating its efficiency and accuracy. Please see [START_REF] Funkhouser | A Beam Tracing Approach to Acoustic Modeling for Interactive Virtual Environments[END_REF][START_REF] Funkhouser | Real-time Acoustic Modeling for Distributed Virtual Environments[END_REF][START_REF] Tsingos | Modeling Acoustics in Virtual Environments Using the Uniform Theory of Diffraction[END_REF] for more details).

SIMULATION SYSTEM

Our simulation system takes as input: 1) a model of a 3D environment described as a set of polygons with frequency-dependent impedances or bidirectional scattering filters, 2) a model of a speaker described by a location and angular radiation pattern, and 3) a model of a microphone described by a location and optional angular directivity pattern. The result of the simulation is an estimated impulse response.

(c) Shortest propagation path generated for the corresponding sequence of reflecting surfaces (green), transmitting portals (purple), and diffracting wedges (magenta). (d)

The procedure repeated for multiple beams reaching the receiver.

The system executes in four steps, as shown in Figure 1. In the first step, we build a spatial subdivision data structure representing a binary space partition of 3D space into convex polyhedral cells. The purpose of this step is to decompose space into cells whose boundaries are aligned with polygons of the 3D input model and whose adjacencies are stored explicitly in a graph structure in order to enable efficient traversals of 3D space during beam tracing.

In the second step, we trace the convex polyhedral beams representing different propagation sequences through cells of the spatial subdivision. The traversal starts in the cell containing the speaker location with a beam representing the entire cell. Next, it visits adjacent cells iteratively, considering different permutations of transmissions, specular reflections and diffractions due to the faces and edges on the boundary of the "current" cell. As the algorithm traverses a cell boundary into a new cell, a copy of the current convex pyramidal beam is "clipped" to include the region of space passing through the convex polygonal boundary to model transmissions. At each reflecting cell boundary, a copy of the transmission beam is mirrored across the plane supporting the cell boundary to model specular reflections. At each diffracting wedge, a new beam is spawned whose source is the edge and whose extent includes all rays predicted by the Geometrical Theory of Diffraction [START_REF] Keller | Geometrical Theory of Diffraction[END_REF]. The traversal along any sequence terminates when either the length of the shortest path within the beam or the cumulative attenuation exceed some user-specified thresholds.

The traversal may also be terminated when the total number of beams traced or the elapsed time exceed other thresholds.

In the third step, for each beam containing the microphone location, we compute the shortest propagation path from the speaker to the microphone along the sequence of transmissions, diffractions, and specular reflections represented by the beam. The intersections with specularly reflecting faces are uniquely determined by the locations of the speaker, microphone, and the intersections with diffracting edges (diffraction points). The diffraction points are found by solving a non-linear system of equations expressing equal angle constraints at diffracting edges. Once the diffraction points are found, we construct a piecewise-linear polyline representing the path along which sound travels from source to receiver along the propagation sequence, from which a length-, angle-, and frequency-dependent filter is computed.

In the final step, for each valid propagation path from the speaker to the microphone, we add its contribution to the simulated impulse response. Our implementation includes source and material filtering effects derived from either measurements or analytical models. We compute diffraction coefficients using the Uniform Theory of Diffraction [START_REF] Keller | Geometrical Theory of Diffraction[END_REF][START_REF] Kouyoumjian | A Uniform Geometrical Theory of Diffraction for an Edge in a Perfectly Conducting Surface[END_REF]. Air absorption due to scattering is also taken into account following the ISO 9013-1 specifications. All calculations are performed in complex Fourier domain at the sampling rate resolution. Our current system uses a sampling rate of 51200 Hz (the sampling rate of our source and material measurements). Thus, for a one second long response we are computing 25600 complex coefficients per path. The simulated transfer function is the sum of the coefficients for all paths, and the final impulse response is the inverse Fourier transform of the sum.

After all beams have been traced up to a user-specified termination criterion and the contribution of all propagation paths have been summed, the resulting impulse response is output for auralization or comparison to measurements.

COMPUTATIONAL RESULTS

In order to evaluate whether our beam tracing framework is computationally practical for large architectural environments, we tested it for modeling specular reflection and diffraction in several 3D models (see Figure 2). For each model, we traced beams from a stationary sound source (located at the red/white dots) and stored them in a beam tree data structure. Then, as the sound receiver was moved interactively (located at the yellow dots), we updated propagation paths comprising specular refection and diffraction at interactive rates. All the tests were run on a Silicon Graphics Onyx2 workstation with a 195MHz R10000 processor. Figure 3 shows a visualization of the simulated results for one 3D model representing a real architectural design (Soda Hall at UC Berkeley). In this case, the model contained 1,762 polygons describing the walls, ceilings, and floors of the building interior. The spatial subdivision constructed 752 cells in 7.1 seconds. Then, the beam tracing algorithm computed 50,000 beams (Figure 3a) comprising specular reflections and diffractions from a source point (white dot) in 13.2 seconds. Next, 75-95 propagation paths (Figure 3b) were updated 20 times per second to a receiver location (pink dot) moving under interactive user control. Finally, echograms (Figure 3c) were computed for specific receiver locations. The system is able to provide interactive updates while visualizing and auralizing sound propagation in this architectural environment.

VALIDATION RESULTS

In order to provide an evaluation of the accuracy of our simulation methods, we compared impulse responses computed with our system versus measurements in two simple environments (see Figure 4). The first is a simple box-shaped enclosure comprising six rectangular panels. The second is the same box-shaped enclosure with a single rectangular panel spanning from floor to ceiling along one half the interior of the box with the speaker and microphone on opposite sides of the panel. The first configuration is very simple, yielding only specular reflections --we use it as a baseline for validation and comparison. The second configuration is a reverberant environment with significant propagation paths combining both edge diffraction and specular reflection. This is a more difficult case, which has not been widely studied. We focus on these simple environments in order to be able to compare the early parts of simulated and measured impulse responses directly, which would be difficult in a complex environment. For both configurations, we measured the impulse response of the environment using a Bruel & Kjaer artificial mouth and a 1/2-inch microphone connected to the audio outputs and inputs of a MOTU 828 multi-channel firewire audio interface [MOTU]. The MOTU interface was connected to an off-the-shelf laptop running Windows. The source signal was a repeated chirp stimuli, and the sampling rate was 48KHz. The output signal used to feed the speaker was also fed back into the interface as a reference. We low-pass filtered the resulting response to get an actual bandwidth of 10KHz.

We simulated all possible 1,524 propagation paths with up to 10th order specular reflection in the first configuration (empty box), and 1,358 paths combining up to 4th order specular reflection and 2nd order edge diffraction in the second configuration (with diffracting baffle). These simulation parameters were chosen as a reasonable compromise between simulation accuracy and computational expense (738 and 631 seconds, respectively). Models of both the sound source and panel responses were acquired through measurements in an anechoic chamber.

The simulated and measured impulse responses for each configuration are shown in Figure 5. Note that the fine band simulation with measured source and material characteristics allows us to compare the waveform of the simulated and measured early responses (numbers on peaks of the impulse response on the top-right correspond to identifiers for specific simulated propagation paths). As can be seen, our simulation is able to capture the temporal structure of the impulse response and is also able to capture the effects of the diffraction by the edge of the panel (e.g., first peak of the top curve in Figure 5b). For comparison, we also plotted a simulated impulse response without diffraction effects (the bottom curve in Figure 5b). We also present in Figure 6 a comparison between spectrograms for the simulated and measured impulse responses. Looking at these plots, we expect that the simulations match the measurements well enough to be useful for prediction or immersive simulation. Please refer to [START_REF] Tsingos | Validation of acoustical simulations in the Bell Labs Box[END_REF] for more details.

DISCUSSION & CONCLUSION

In this paper, we have described a beam tracing approach to simulating sound propagation in architectural environments. Its advantageous features are: 1) it enumerates all potential sequences of specular reflection, diffraction, and transmission up to user-specified termination criteria, 2) it enables evaluation of early propagation paths at interactive rates, and 3) it scales to support large, "densely-occluded" architectural environments (e.g., buildings and cities).

The approach is not without limitations. In particular, the current implementation is practical only for coarsely detailed 3D models, since beams get fragmented in scenes with many free-space cell boundaries [START_REF] Fortune | Topological Beam Tracing[END_REF][START_REF] Funkhouser | A Beam Tracing Approach to Acoustic Modeling for Interactive Virtual Environments[END_REF][START_REF] Teller | Visibility Computations in Densely Occuded Polyhedral Environments[END_REF]. Also, it models the building structure at only one scale, which may be appropriate for some sound wavelengths, but not others. Finally, our auralization methods are tuned for interactive applications and thus employ approximations that are quick to compute (e.g., the Geometric Theory of Diffraction) but less accurate than alternative models. In future work, we plan to investigate how our beam tracing methods can be used for enumerating propagation sequences in off-line simulations with more accurate models to improve the efficiency and accuracy of existing acoustics prediction systems.

Figure 1 :

 1 Figure 1: Overview of simulation process: a) Virtual environment (office cubicles) with source S, receiver R, and spatial subdivision marked in pink. (b) Sample reflected and diffracted sound beam (cyan) reaching the receiver. (c) Shortest propagation path generated for the corresponding sequence of reflecting surfaces (green), transmitting portals (purple), and diffracting wedges (magenta). (d) The procedure repeated for multiple beams reaching the receiver.

Figure 2 :

 2 Figure 2: Visualizations of sound simulations in a city (left) and a concert hall (right).

Figure 3 :

 3 Figure 3: Visualization of a) propagation paths (yellow lines) traced from source (white dot) to receiver (orange dot), b) power plot, and c) echogram simulated with 50,000 beams in model of a building with 1,762 polygons.

Figure 4 :

 4 Figure 4: Two simple environments used for validation experiments. One is a simple six-sided box (left), while the other contains an interior baffle (middle). An interior view of the baffle is shown in the image on the right.

Figure 5 :

 5 Figure 5: Comparison of simulated (above) and measured (below) impulse responses for the empty box (left) and box with a diffracting baffle (right). The plot on the right includes curves for simulations both with (blue) and without (yellow) diffraction.

Figure 6 :

 6 Figure 6: Spectrograms of measured (left) and simulated (right) early responses in box with a baffle.

ACKNOWLEDGEMENTS

The authors thank Sid Ahuja for his support of this work. Thomas Funkhouser is also supported by a National Science Foundation CAREER grant (CCR-0093343).