Impulse Response Plot and Receiver Dataset
The impulse response (IR) in a room acoustics simulation can be analyzed by collecting the ray information using the Receiver dataset and then postprocess using the Impulse Response plot. Room acoustic metrics like reverberation times, early energy metrics, and speech intelligibility can also be postprocessed using the Energy Decay subfeature. The IR can be determined for sources and boundary data (absorption coefficients, scattering parameters, source directivity, volume attenuation etc.) given in octave, 1/3 octave, or 1/6 octave bands.
Preparing a Room Acoustics Simulation
When setting up a Ray Acoustics simulation and preparing it for a room acoustics simulation where the impulse response (IR) is an important output, certain considerations and steps need to be taken. They are as follows:
•
Decide on the interpretation of the frequency in the ray acoustics simulation. To be able to reconstruct the IR with the dedicated plot, the frequency should represent the center frequency of octaves, 1/3 octave, or 1/6 octave bands. The center frequency is also used as input for boundary conditions (absorption coefficients and scattering parameter), for sources (source strength and directivity), and for the volume absorption in the medium, typically air.
•
Create a parameter under Global Definitions > Parameters for the band center frequency, for example, f0.
•
Set up interpolation functions for the properties that depend on the frequency. The data can easily be stored in a Spreadsheet formatted file and interpolation can be set to Nearest neighbor. An example, of the format for a (.txt) file, that specifies absorption coefficients for some boundaries in a simulation could be (data is invented):
%f0, wall, entrance, window, floor, diffusers, seats
125 0.1 0.4 0.01 0.02 0.1 0.05
250 0.07 0.2 0.01 0.03 0.1 0.06
500 0.05 0.12 0.02 0.05 0.1 0.07
1000 0.04 0.07 0.02 0.1 0.1 0.15
2000 0.04 0.05 0.02 0.3 0.1 0.35
4000 0.05 0.05 0.03 0.5 0.1 0.4
8000 0.08 0.05 0.04 0.5 0.1 0.6
16000 0.1 0.05 0.05 0.5 0.1 0.6
•
Make sure that either the Compute intensity and power or the Compute power option is selected in the Intensity Computation section. The Compute power option requires fewer degrees of freedom and is more robust. If it is selected the ray intensity (rac.I and rac.logI) cannot be visualized in ray trajectories plots, but both the Sound Pressure Level Calculation on boundaries and the impulse response can be computed.
•
Make sure that Count reflections is selected in the Additional Variables section.
•
Also make sure that the parameter that you created for the band center frequency (say f0) is used as the Ray frequency under the Ray Properties node.
•
Set up the ray acoustics model with appropriate boundary conditions and sources.
•
To be able to postprocess the IR, a Parametric Sweep needs to be used in the study around the Ray Tracing study step. The sweep parameter should be the frequency defined in the parameters. To easily create a Parameter value list representing the center frequencies of the bands use the ISO preferred frequencies entry method.
•
For the Times, specified under the Ray Tracing study step, only enter 0 and the final simulation time (the final time should be around the estimated reverberation time). COMSOL uses small internal steps that accurately account for all reflections. The so-called Extra Time Steps are used when reconstructing the IR and other data (see below).
•
Run the simulation.
Receiver Dataset
Use a Receiver 2D () or Receiver 3D dataset (), selected from the More 2D Datasets submenu and the More 3D Datasets submenu, respectively, to collect the data necessary to visualize the impulse response using an Impulse Response Plot. The Impulse Response plot uses the data from a Receiver dataset as input
Data
Select the appropriate ray dataset in the Dataset list. Then select the frequency parameters that should be used by the receiver; select one frequency for a single band analysis or all for the broad band analysis. This is the typical behavior if a parametric sweep has been set up, as described in Preparing a Room Acoustics Simulation.
Receiver
The receiver behaves as a (transparent) sphere with a given radius and located at a center position (the microphone needs to have a certain finite size such that there is a reasonable probability of rays interacting with it). By changing the center position the receiver can be moved around without the need to run the simulation again.
•
Under Center, specify the x, y, and (3D only) z coordinates for the center of the receiver.
•
Under Radius, specify the radius of the receiver. From the Radius input list, choose Expression to determine the radius using an expression (see below) or choose User defined to enter a value for the radius in the Radius field (SI unit: m). Different theories exist for the appropriate size of the receiver.
The following built-in expression determines the radius R of the receiver:
where you enter the values for the Number of rays, N; Room volume, V (SI unit: m2); and Source-receiver distance, dSR (SI unit: m), to determine the radius. See Ref. 19 for details on the expression.
Directivity
From the Directivity type list, choose Omnidirectional (the default), User defined (dB), or User defined (linear) and enter a directivity expression in the Expression field.
•
The User defined (dB) option applies the given gain to the received signal (positive or negative gain).
•
The User defined (linear) option applies a linear gain. If a negative value is used the phase of the arriving ray is inverted.
The expression can depend on the ray direction vector components (rac.nix, rac.niy, rac.niz) in order to set up a complex directivity pattern.
Extra Time Steps
Use the default All option. The impulse response plot has been tuned to use this option in order to get accurate and fast evaluation of arrival time and ray power.
Advanced
For the Interpolation between time steps list always use the default Linear interpolation option. The Cubic interpolation is not applicable for the impulse response computation and can lead to erroneous results.
Under Normal variables and Other variables, if desired, you can change the default names of the created variables:
•
The normal direction of the receiver surface (nx, ny, and, in 3D, nz).
•
The Distance traveled, that is the distance traveled by a ray inside the receiver.
•
The receiver Volume.
•
The Directivity, that can be used to visualize the expression for the directivity.
•
The First ray arrival time.
Evaluation and Export
The main purpose for the Receiver dataset is to be used as input for the Impulse Response Plot. However, the data generated by the dataset can also be exported and used in an external analysis tool or software (ray power and arrival times).
Right-click the dataset and select Add Data to Export, then enter the desired data. For example, export the arrival time of the rays t, the frequency of the rays rac.f, the power of the rays rac.Q, the distance traveled by the rays in receiver re1dist, and so forth. The exported data can be stored and sorted in several formats.
When working with phase alignment of loudspeakers, it can be useful to evaluate the arrival time if the first ray. This is easily done in a Derived Values > Global Evaluation feature that points to the receiver dataset. Simply evaluate the special variable re1first.
See also the Receiver 2D and Receiver 3D section in the COMSOL Multiphysics Reference Manual.
Impulse Response Plot
In room acoustics applications, the impulse response (IR) of a source-receiver configuration represents one of the most important postprocessing results. To create the plot based on the data collected by the Receiver dataset, add an Impulse Response subnode () to a 1D Plot Group to create the impulse response plot. After defining the characteristics of the impulse response, click Plot () to create the plot. The Energy Decay () subnode can be added to perform an analysis of the impulse response to get room acoustic metrics and the energy/level decay curves. The first rendering of the IR plot can take some time. As the data is cached any subsequent changes to the plot will be fast.
The Small Concert Hall Acoustics model is a tutorial on how to compute the impulse response using the Ray Acoustics interface. Application Library path: Acoustics_Module/Building_and_Room_Acoustics/ small_concert_hall
In the ray tracing methodology it is, as discussed above, typical to characterize sources and boundaries in frequency bands. This also means that some of the frequency content is lost. This frequency content is reconstructed when determining or reconstructing the IR. For each ray that interacts with the receiver, a small piece of signal with appropriate amplitude, arrival time, phase, and frequency content is added to the IR (the filter bank kernel). The sum of all the contributions reconstruct the full IR.
Data
From the Dataset list, choose the appropriate Receiver Dataset. Then from the Frequency interpretation list choose an interpretation of the frequency: Octave (the default), 1/3 octave, or 1/6 octave. This selections should coincide with the frequency interpretation used in boundary conditions and sources in the underlying ray acoustics simulation, as described in Preparing a Room Acoustics Simulation. The frequency content and resolution of the impulse response is based on this selection.
Expressions
Set up the variables that are necessary for reconstructing the impulse response (typically use the defaults). These variables are used to define and add the contribution of each ray, that intersects the receiver dataset, to the impulse response.
•
Define the frequency variable for the rays (SI unit: Hz) in the Frequency field. The default is rac.f.
•
Define the power variable for the rays (SI unit: kg·m2/s2) in the Power field. The default is rac.Q.
•
Define the density (SI unit: kg/m3) in the Density field. The default is rac.rho.
•
Define the speed of sound (SI unit: m/s) in the Speed of Sound field. The default is rac.c.
•
Define the number of reflections undergone by the rays (SI unit: 1) in the Number of reflections field. The default is rac.Nrefl.
x-Axis Data
From the Transformation list, choose a transformation of the data on the x-axis: None (the default), for no transformation, to show the time domain signal, or Frequency spectrum. The latter will perform an FFT of the IR and depict the power spectra as function of frequency.
Advanced
The following advanced settings are available to control some aspects of the IR reconstruction as well as the underlying filter kernel used.
•
Define a sampling frequency fs (SI unit: Hz). The default is 44100.
•
Define the Zero passing length. The default is 11025.
The next two options are legacy options (to achieve behavior for COMSOL Multiphysics releases older than version 5.6).
•
Select Remove noncausal signal to remove the signal (force it to zero) prior to the arrival of the first ray. Using this option will alter the energy content of the early energy component of the IR.
•
Select Use fully randomized phase to use the legacy method for reconstructing the IR. The fully randomized phase method (see Ref. 20) can lead to unwanted cancellations of the direct sound and early energy contributions.
The following options control the filter kernel used to reconstruct the impulse response. Advanced options allow to modify the kernel, enter a user-defined analytical kernel, or a user-defined kernel based on an interpolation function. The default is to use a Brick-wall with Kaiser window filter kernel (see Ref. 21).
•
Select Show the filters to show and inspect the filter bank in the Graphics window, click Plot () to create the plot. To switch between time domain and frequency domain representation change the Transformation option in the x-Axis Data section.
•
Select the Filter kernel to use Brick-wall with Kaiser window (the default) or User defined. When User defined is selected enter the expression for the filter kernel hfc[m]. Use and set up parameters in the Parameters table (as you would use variables). The parameter fc is reserved for the band center frequency; fs is reserved for the sampling frequency; and m is the zero centered index. Note that per default the User defined option sets up the rick-wall with Kaiser window filter kernel for octaves.
•
If the default Brick-wall with Kaiser window is selected. You can control the ripple factor δ (the default is 0.05) and the passband slope factor E (the default is 1.03).
Coloring and Style, and Legends
You can also make change in the Coloring and Style section or add and format a legend in the Legends section.
See also the Impulse Response section tn the COMSOL Multiphysics Reference Manual.
Energy Decay and Objective Room Acoustic Metrics
For further analysis of the impulse response, add the Energy Decay () subnode to the Impulse Response plot. Once selections have been made, click Plot () to create the plot.
Display
Select the Band type as Broadband (the default) or Individual bands. For the Individual bands option select the Band frequency, either an individual band or All. For All, all the level or energy decay curves are plotted.
Select the Plot as Energy decay (the default), Level decay, or Modulation transfer function, to define the type of graphical analysis of the impulse response to be plotted in the Graphics window.
Table
Select the room acoustic metrics to be included in the results table, per default all objective quality metrics are selected. The default name of the table is Objective Quality Metrics. The table displays the results in columns with the band center frequency fc in the first column. The results table can also be located under the Results > Tables node. A plot of the results located in the table can be created using the Table Graph plot.
The speech transmission index (STI) metric is defined as a single value that combines information from 14 modulation transfer function bands and 7 octave bands. To get the value specified by the standard solve the model for octaves from 125 Hz to 8 kHz. When evaluating use the Band type set to the default Broadband option.
Custom Plots
The energy response depicted as a reflectogram or echogram can be plotted using a Ray (Plot) and depict the expression re1dist*rac.Q/re1vol. This will give the intensity received by the “microphone” for each ray.