diff --git a/README.md b/README.md index 73f8ca2..235a34e 100644 --- a/README.md +++ b/README.md @@ -49,10 +49,10 @@ pip3 install . Usage -------- -See `examples` folder (they require `libsndfile` module installed). You need to `cd examples` and run `python3 mesh_sim.py` (we recommend starting with this one). This script demonstrates two equivalent ways to define the environment for sound propagation, and save the impulse response as an audio file. You can use a `.obj` file with an optional `.mtl` file with the same name to define the room geometry and materials. In this case, the `.mtl` file has two extra rows compared with conventional `.mtl` file used for visual rendering: +See `examples` folder (extra modules may be required). You need to `cd examples` and run `python3 mesh_sim.py` (we recommend starting with this one). This script demonstrates two equivalent ways to define the environment for sound propagation, and save the impulse response as an audio file. You can use a `.obj` file with an optional `.mtl` file with the same name to define the room geometry and materials. In this case, the `.mtl` file has two extra rows compared with conventional `.mtl` file used for visual rendering: ``` sound_a 0.5 0.6 0.6 0.7 0.75 0.8 0.9 0.9 # sound absorption coefficients, for 8 octave bands [62.5, 125, 250, 500, 1000, 2000, 4000, 8000]Hz -sound_s 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 # sound scattering coefficients, if you don't know the details of diffuse/specular reflections, keep it 0.5 +sound_s 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 # sound scattering coefficients, if you don't know the details of diffuse/specular reflections, keep it low ``` or directly create a shoebox shaped room using our API: ``` diff --git a/examples/compare_image.py b/examples/compare_image.py index b6d1169..58a3823 100644 --- a/examples/compare_image.py +++ b/examples/compare_image.py @@ -1,13 +1,14 @@ import numpy as np import pygsound as ps import matplotlib.pyplot as plt -import rirgenerator as rg +import rir_generator as rir +# using https://github.com/audiolabs/rir-generator/ def main(): l = 10 w = 6 - h = 2 + h = 4 absorb = 0.3 reflec = np.sqrt(1.0 - absorb) @@ -17,7 +18,7 @@ def main(): ctx.channel_type = ps.ChannelLayoutType.mono scene = ps.Scene() - mesh = ps.createbox(l, w, h, absorb, 0.5) + mesh = ps.createbox(l, w, h, absorb, 0.1) scene.setMesh(mesh) src_coord = [1, 1, 0.5] @@ -29,39 +30,34 @@ def main(): lis = ps.Listener(lis_coord) lis.radius = 0.01 - rir_gs = scene.computeIR(src, lis, ctx) - rir_img = rg.rir_generator(343, 16000, lis_coord, src_coord, [l, w, h], beta=[reflec]*6) - rir_img = rir_img / max(abs(rir_img[0])) - - max_cnt = min(len(rir_gs['samples']), len(rir_img[0])) * 2 - - fig, axs = plt.subplots(4, 1, figsize=(10, 20)) + rir_gs = scene.computeIR([src], [lis], ctx) + rir_img = rir.generate( + c=343, # Sound velocity (m/s) + fs=ctx.sample_rate, # Sample frequency (samples/s) + r=[lis_coord], + s=src_coord, # Source position [x y z] (m) + L=[l, w, h], # Room dimensions [x y z] (m) + beta=[reflec] * 6, # Reflection coefficients + ) + rir_img = rir_img[:, 0] / max(abs(rir_img)) + rir_gs = rir_gs['samples'][0][0][0] + max_cnt = min(len(rir_gs), len(rir_img)) + + fig, axs = plt.subplots(2, 1, figsize=(16, 10)) axs[0].set_title('Image Method (waveform)') - axs[0].plot(rir_img[0], linewidth=0.5) + axs[0].plot(rir_img, linewidth=0.5) axs[0].set_xlabel('Sample') axs[0].set_xlim(0, max_cnt) axs[0].set_ylim(-1, 1) axs[0].set_ylabel('Amplitude') axs[1].set_title('Geometric Sound Propagation (waveform)') - axs[1].plot(rir_gs['samples'], linewidth=0.5) + axs[1].plot(rir_gs, linewidth=0.5) axs[1].set_xlabel('Sample') axs[1].set_xlim(0, max_cnt) axs[1].set_ylim(-1, 1) axs[1].set_ylabel('Amplitude') - axs[2].set_title('Image Method (spectrogram)') - axs[2].specgram(rir_img[0], mode='magnitude', NFFT=1024, Fs=16000, noverlap=512) - axs[2].set_xlim(0, max_cnt / 16000) - axs[2].set_xlabel('Times (s)') - axs[2].set_ylabel('Frequency (Hz)') - - axs[3].set_title('Geometric Sound Propagation (spectrogram)') - axs[3].specgram(rir_gs['samples'], mode='magnitude', NFFT=1024, Fs=16000, noverlap=512) - axs[3].set_xlim(0, max_cnt / 16000) - axs[3].set_xlabel('Times (s)') - axs[3].set_ylabel('Frequency (Hz)') - fig.tight_layout() plt.savefig('img_comparison.png') plt.show() diff --git a/examples/cube.mtl b/examples/cube.mtl index f77b71c..d20823f 100644 --- a/examples/cube.mtl +++ b/examples/cube.mtl @@ -11,4 +11,4 @@ Ni 1.000000 d 1.000000 illum 2 sound_a 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 -sound_s 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 +sound_s 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 diff --git a/examples/custom_array.py b/examples/custom_array.py index 302d8b2..f759043 100644 --- a/examples/custom_array.py +++ b/examples/custom_array.py @@ -1,11 +1,10 @@ -import os import numpy as np import pygsound as ps from wavefile import WaveWriter, Format def compute_array(meshpath, src_coord, lis_coord, r, s, micarray): - mesh = ps.loadobj(meshpath, os.path.join(os.path.dirname(meshpath), ''), r, s) + mesh = ps.loadobj(meshpath, r, s) ctx = ps.Context() ctx.diffuse_count = 20000 ctx.specular_count = 2000 @@ -19,14 +18,13 @@ def compute_array(meshpath, src_coord, lis_coord, r, s, micarray): res = {} res_buffer = [] rate = 0 - abandon_flag = False for offset in micarray: lis = ps.Listener((offset + lis_coord).tolist()) lis.radius = 0.01 - res_ch = scene.computeIR(src, lis, ctx) + res_ch = scene.computeIR([src], [lis], ctx) rate = res_ch['rate'] - sa = res_ch['samples'] + sa = res_ch['samples'][0][0][0] res['rate'] = rate res_buffer.append(sa) res['samples'] = np.zeros((len(res_buffer), np.max([len(ps) for ps in res_buffer]))) @@ -43,7 +41,7 @@ def main(): src_coord = [1, 1, 1] lis_coord = [0.1, 0.1, 0.1] - res = compute_array("cube.obj", src_coord, lis_coord, 0.5, 0.5, custom_array) + res = compute_array("cube.obj", src_coord, lis_coord, 0.5, 0.1, custom_array) with WaveWriter('custom_array.wav', channels=np.shape(res['samples'])[0], samplerate=int(res['rate'])) as w: w.write(np.array(res['samples'])) diff --git a/examples/mesh_sim.py b/examples/mesh_sim.py index 969ecc4..4fd898e 100644 --- a/examples/mesh_sim.py +++ b/examples/mesh_sim.py @@ -10,7 +10,7 @@ def main(): ctx.specular_count = 2000 ctx.channel_type = ps.ChannelLayoutType.stereo - mesh1 = ps.loadobj("cube.obj", "") # if the second argument is empty, the code will infer the .mtl name using .obj name + mesh1 = ps.loadobj("cube.obj") scene = ps.Scene() scene.setMesh(mesh1) @@ -19,29 +19,26 @@ def main(): src = ps.Source(src_coord) src.radius = 0.01 + src.power = 1.0 lis = ps.Listener(lis_coord) lis.radius = 0.01 - res = scene.computeMultichannelIR(src, lis, ctx) - - with WaveWriter('test1.wav', channels=np.shape(res['samples'])[0], samplerate=int(res['rate'])) as w1: - w1.write(np.array(res['samples'])) + res = scene.computeIR([src], [lis], ctx) # you may pass lists of sources and listeners to get N_src x N_lis IRs + audio_data = np.array(res['samples'][0][0]) # the IRs are indexed by [i_src, i_lis, i_channel] + with WaveWriter('test1.wav', channels=audio_data.shape[0], samplerate=int(res['rate'])) as w1: + w1.write(audio_data) print("IR using .obj input written to test1.wav.") # Simulation using a shoebox definition - ctx = ps.Context() - ctx.diffuse_count = 20000 - ctx.specular_count = 2000 - ctx.channel_type = ps.ChannelLayoutType.stereo - - mesh2 = ps.createbox(10, 6, 2, 0.5, 0.5) + mesh2 = ps.createbox(10, 6, 2, 0.5, 0.1) scene = ps.Scene() scene.setMesh(mesh2) - res = scene.computeMultichannelIR(src, lis, ctx) - with WaveWriter('test2.wav', channels=np.shape(res['samples'])[0], samplerate=int(res['rate'])) as w: - w2.write(np.array(res['samples'])) + res = scene.computeIR([src_coord], [lis_coord], ctx) # use default source and listener settings if you only pass coordinates + audio_data = np.array(res['samples'][0][0]) + with WaveWriter('test2.wav', channels=audio_data.shape[0], samplerate=int(res['rate'])) as w2: + w2.write(audio_data) print("IR using shoebox input written to test2.wav.") diff --git a/setup.py b/setup.py index 3416c19..386ee38 100644 --- a/setup.py +++ b/setup.py @@ -72,7 +72,7 @@ def build_extension(self, ext): setup( name='pygsound', - version='0.2', + version='0.3', author='Zhenyu Tang, Carl Schissler, Dinesh Manocha', author_email='zhy@umd.edu', description='A room impulse response simulator using for geometric sound propagation', @@ -92,6 +92,5 @@ def build_extension(self, ext): install_requires=[ 'Cython', 'numpy', - 'wavefile', ], ) diff --git a/src/GSound/gsound/gsMeshRequest.cpp b/src/GSound/gsound/gsMeshRequest.cpp index b7ae87f..ddb97d2 100755 --- a/src/GSound/gsound/gsMeshRequest.cpp +++ b/src/GSound/gsound/gsMeshRequest.cpp @@ -74,7 +74,7 @@ MeshRequest:: MeshRequest() weldTolerance( 0.00001f ), simplifyTolerance( 0.0001f ), minDiffractionEdgeAngle( 0.0f ), - minDiffractionEdgeLength( 0.0f ), + minDiffractionEdgeLength( 0.5f ), diffuseResolution( 0.5f ), edgeResolution( 0.5f ), minRaysPerEdge( 1 ), diff --git a/src/GSound/gsound/gsMeshRequest.h b/src/GSound/gsound/gsMeshRequest.h index 82a540a..c9e083b 100755 --- a/src/GSound/gsound/gsMeshRequest.h +++ b/src/GSound/gsound/gsMeshRequest.h @@ -125,7 +125,7 @@ class MeshRequest * If the angle between normals of two neighboring triangles * is less than this value, then the edge that they share is not diffracting. * Thus, a larger diffraction angle will result in less diffracting edges while - * lower thresold will result in more edges. + * lower threshold will result in more edges. */ Real minDiffractionEdgeAngle; diff --git a/src/pygsound/src/Context.cpp b/src/pygsound/src/Context.cpp index e05fe06..cc0b5ab 100644 --- a/src/pygsound/src/Context.cpp +++ b/src/pygsound/src/Context.cpp @@ -12,7 +12,7 @@ Context::Context() prop_request.flags.set(gs::PropagationFlags::DIFFUSE, true); prop_request.flags.set(gs::PropagationFlags::DIFFRACTION, true); prop_request.flags.set(gs::PropagationFlags::SOURCE_DIRECTIVITY, false); - prop_request.flags.set(gs::PropagationFlags::DOPPLER_SORTING, true); + prop_request.flags.set(gs::PropagationFlags::DOPPLER_SORTING, false); prop_request.flags.set(gs::PropagationFlags::ADAPTIVE_QUALITY, false); prop_request.flags.set(gs::PropagationFlags::AIR_ABSORPTION, true); prop_request.flags.set(gs::PropagationFlags::ADAPTIVE_IR_LENGTH, true); diff --git a/src/pygsound/src/Listener.hpp b/src/pygsound/src/Listener.hpp index 1074713..bc8cfe7 100644 --- a/src/pygsound/src/Listener.hpp +++ b/src/pygsound/src/Listener.hpp @@ -36,7 +36,7 @@ class Listener friend std::ostream &operator<<( std::ostream &_strm, const Listener &_list ); - gs::SoundListener m_listener; + gs::SoundListener m_listener; }; #endif // INC_LISTENER_HPP diff --git a/src/pygsound/src/Scene.cpp b/src/pygsound/src/Scene.cpp index 13973e4..a5fe06c 100644 --- a/src/pygsound/src/Scene.cpp +++ b/src/pygsound/src/Scene.cpp @@ -21,65 +21,109 @@ Scene::setMesh( SoundMesh &_mesh ) } py::dict -Scene::computeIR( SoundSource &_source, Listener &_listener, Context &_context ) +Scene::computeIR( std::vector &_sources, std::vector &_listeners, Context &_context ) { - m_scene.addSource( &_source.m_source); - m_scene.addListener( &_listener.m_listener ); - if (m_scene.getObjectCount() == 0){ - std::cout << "object count is zero, cannot propagate sound!" << std::endl; - } + int n_src = _sources.size(); + int n_lis = _listeners.size(); + + for (SoundSource& p : _sources){ + m_scene.addSource(&p.m_source); + } + for (Listener& p : _listeners){ + m_scene.addListener(&p.m_listener); + } + + if (m_scene.getObjectCount() == 0){ + std::cerr << "object count is zero, cannot propagate sound!" << std::endl; + } + propagator.propagateSound(m_scene, _context.internalPropReq(), sceneIR); - const gs::SoundSourceIR& sourceIR = sceneIR.getListenerIR(0).getSourceIR(0); - gs::ImpulseResponse result; - result.setIR(sourceIR, _listener.m_listener, _context.internalIRReq()); - auto rate = result.getSampleRate(); + py::list IRPairs(n_src); + auto rate = _context.getSampleRate(); + for (int i_src = 0; i_src < n_src; ++i_src){ + py::list srcSamples(n_src); + for (int i_lis = 0; i_lis < n_lis; ++i_lis){ + const gs::SoundSourceIR& sourceIR = sceneIR.getListenerIR(i_lis).getSourceIR(i_src); + gs::ImpulseResponse result; + result.setIR(sourceIR, *m_scene.getListener(i_lis), _context.internalIRReq()); + auto numOfChannels = int(result.getChannelCount()); - auto *sample = result.getChannel( 0 ); - std::vector samples(sample, sample+result.getLengthInSamples()); + py::list samples; + for (int ch = 0; ch < numOfChannels; ch++) + { + auto *sample_ch = result.getChannel(ch); + std::vector samples_ch(sample_ch, sample_ch+result.getLengthInSamples()); + samples.append(samples_ch); + } + srcSamples[i_lis] = samples; + } + IRPairs[i_src] = srcSamples; + } - m_scene.clearSources(); - m_scene.clearListeners(); + m_scene.clearSources(); + m_scene.clearListeners(); - py::dict ret; - ret["rate"] = rate; - ret["samples"] = samples; - return ret; + py::dict ret; + ret["rate"] = rate; + ret["samples"] = IRPairs; // index by [i_src, i_lis, i_channel] -// return result; + return ret; } py::dict -Scene::computeMultichannelIR( SoundSource &_source, Listener &_listener, Context &_context ) +Scene::computeIR( std::vector> &_sources, std::vector> &_listeners, Context &_context, + float src_radius, float src_power, float lis_radius) { - m_scene.addSource( &_source.m_source); - m_scene.addListener( &_listener.m_listener ); - if (m_scene.getObjectCount() == 0){ - std::cout << "object count is zero, cannot propagate sound!" << std::endl; - } - propagator.propagateSound(m_scene, _context.internalPropReq(), sceneIR); - const gs::SoundSourceIR& sourceIR = sceneIR.getListenerIR(0).getSourceIR(0); -// gs::ImpulseResponse result; -// result.setIR(sourceIR, _listener.m_listener, _context.internalIRReq()); - auto result = std::make_shared(); - result->setIR(sourceIR, _listener.m_listener, _context.internalIRReq()); - auto rate = result->getSampleRate(); - - auto numOfChannels = int(result->getChannelCount()); - assert(numOfChannels > 0); - - py::list samples; - for (int ch = 0; ch < numOfChannels; ch++) - { - auto *sample_ch = result->getChannel(ch); - std::vector samples_ch(sample_ch, sample_ch+result->getLengthInSamples()); - samples.append(samples_ch); + int n_src = _sources.size(); + int n_lis = _listeners.size(); + + // listener propagation is most expensive, so swap them for computation if there are more listeners + bool swapBuffer = false; + auto src_pos = std::ref(_sources); + auto lis_pos = std::ref(_listeners); + if (n_src < n_lis){ + swapBuffer = true; + src_pos = std::ref(_listeners); + lis_pos = std::ref(_sources); + std::swap(n_src, n_lis); } + + std::vector sources; + std::vector listeners; + + for (auto p : src_pos.get()){ + auto source = new SoundSource(p); + source->setRadius(src_radius); + source->setPower(src_power); + sources.push_back( *source ); + } + for (auto p : lis_pos.get()){ + auto listener = new Listener(p); + listener->setRadius(lis_radius); + listeners.push_back( *listener ); + } + + py::dict _ret = computeIR(sources, listeners, _context); + py::list IRPairs = _ret["samples"]; + py::dict ret; - ret["rate"] = rate; - ret["samples"] = samples; - m_scene.clearSources(); - m_scene.clearListeners(); + ret["rate"] = _ret["rate"]; - return ret; + // swap the IR buffer if sources and listeners have been swapped + if (swapBuffer){ + py::list swapIRPairs; + for (int i_lis = 0; i_lis < n_lis; ++i_lis) { + py::list srcSamples; + for (int i_src = 0; i_src < n_src; ++i_src){ + srcSamples.append(IRPairs[i_src].cast()[i_lis]); + } + swapIRPairs.append(srcSamples); + } + ret["samples"] = swapIRPairs; + }else{ + ret["samples"] = IRPairs; // index by [i_src, i_lis, i_channel] + } + + return ret; } diff --git a/src/pygsound/src/Scene.hpp b/src/pygsound/src/Scene.hpp index 0638763..fd17db6 100644 --- a/src/pygsound/src/Scene.hpp +++ b/src/pygsound/src/Scene.hpp @@ -24,8 +24,9 @@ class Scene void setMesh( SoundMesh &_mesh ); - py::dict computeIR( SoundSource &_source, Listener &_listener, Context &_context ); - py::dict computeMultichannelIR( SoundSource &_source, Listener &_listener, Context &_context ); + py::dict computeIR( std::vector &_sources, std::vector &_listeners, Context &_context ); + py::dict computeIR( std::vector> &_sources, std::vector> &_listeners, Context &_context, + float src_radius = 0.01, float src_power = 1.0, float lis_radius = 0.01); public: diff --git a/src/pygsound/src/SoundMesh.cpp b/src/pygsound/src/SoundMesh.cpp index 8127df1..4a9b6ae 100644 --- a/src/pygsound/src/SoundMesh.cpp +++ b/src/pygsound/src/SoundMesh.cpp @@ -20,26 +20,26 @@ namespace oms = om::sound; namespace omm = om::math; std::shared_ptr< SoundMesh > -SoundMesh::loadObj( const std::string &_path, std::string _basepath, float _forceabsorp, float _forcescatter ) +SoundMesh::loadObj( const std::string &_path, float _forceabsorp, float _forcescatter ) { tinyobj::attrib_t attrib; std::vector< tinyobj::shape_t > shapes; std::vector< tinyobj::material_t > materials; std::string err; - if (_basepath.empty()) - { - auto loc = _path.rfind('/'); - _basepath = _path; - _basepath.erase(loc + 1); - //std::cout << "mesh path: " << _path << std::endl << "inferred base path: " << _basepath << std::endl; - } + std::string basepath; + auto loc = _path.rfind('/'); + basepath = _path; + basepath.erase(loc + 1); - bool errv = tinyobj::LoadObj( &attrib, &shapes, &materials, &err, _path.c_str(), _basepath.c_str() ); + + bool errv = tinyobj::LoadObj( &attrib, &shapes, &materials, &err, _path.c_str(), basepath.c_str() ); if ( !errv ) throw std::runtime_error( err ); else if ( !err.empty() ) std::cerr << "WARNING: " << err << "\n"; + if ( materials.empty() ) + std::cerr << "WARNING: no material loaded for " << _path << "\n"; std::vector< gs::SoundVertex > verts; std::vector< gs::SoundTriangle > tris; @@ -164,10 +164,12 @@ SoundMesh::loadObj( const std::string &_path, std::string _basepath, float _forc auto ret = std::make_shared< SoundMesh >(); gs::SoundMeshPreprocessor preprocessor; - + gs::MeshRequest meshRequest; + meshRequest.minDiffractionEdgeAngle = 30; + meshRequest.minDiffractionEdgeLength = 0.5; if ( !preprocessor.processMesh( &verts[0], verts.size(), &tris[0], tris.size(), - &mats[0], mats.size(), gs::MeshRequest(), ret->m_mesh ) ) + &mats[0], mats.size(), meshRequest, ret->m_mesh ) ) throw std::runtime_error( "Cannot preprocess sound mesh!" ); diff --git a/src/pygsound/src/SoundMesh.hpp b/src/pygsound/src/SoundMesh.hpp index 67794cd..4ea21ac 100644 --- a/src/pygsound/src/SoundMesh.hpp +++ b/src/pygsound/src/SoundMesh.hpp @@ -11,9 +11,9 @@ class SoundMesh { public: - static std::shared_ptr< SoundMesh > loadObj( const std::string &_path, std::string _basepath, float _forceabsorp = -1.0, float _forcescatter = -1.0 ); - static std::shared_ptr< SoundMesh > createBox( float _width, float _length, float _height, float _absorp = 0.5, float _scatter = 0.2 ); - static std::shared_ptr< SoundMesh > createBox( float _width, float _length, float _height, std::vector _absorp, float _scatter = 0.2 ); + static std::shared_ptr< SoundMesh > loadObj( const std::string &_path, float _forceabsorp = -1.0, float _forcescatter = -1.0 ); + static std::shared_ptr< SoundMesh > createBox( float _width, float _length, float _height, float _absorp = 0.5, float _scatter = 0.1 ); + static std::shared_ptr< SoundMesh > createBox( float _width, float _length, float _height, std::vector _absorp, float _scatter = 0.1 ); gsound::SoundMesh &mesh() { return m_mesh; } diff --git a/src/pygsound/src/module.cpp b/src/pygsound/src/module.cpp index 345e9ca..0e399a7 100644 --- a/src/pygsound/src/module.cpp +++ b/src/pygsound/src/module.cpp @@ -42,20 +42,21 @@ PYBIND11_MODULE(pygsound, ps) .def(py::init<>()); ps.def( "loadobj", &SoundMesh::loadObj, "A function to load mesh and materials", - py::arg("_path"), py::arg("_basepath"), py::arg("_forceabsorp") = -1.0, py::arg("_forcescatter") = -1.0 ); + py::arg("_path"), py::arg("_forceabsorp") = -1.0, py::arg("_forcescatter") = -1.0 ); ps.def( "createbox", py::overload_cast(&SoundMesh::createBox), "A function to create a simple shoebox mesh", py::arg("_width"), py::arg("_length"), py::arg("_height"), - py::arg("_absorp") = 0.5, py::arg("_scatter") = 0.5 ); + py::arg("_absorp") = 0.5, py::arg("_scatter") = 0.1 ); ps.def( "createbox", py::overload_cast, float>(&SoundMesh::createBox), "A function to create a simple shoebox mesh", py::arg("_width"), py::arg("_length"), py::arg("_height"), - py::arg("_absorp"), py::arg("_scatter") = 0.5 ); + py::arg("_absorp"), py::arg("_scatter") = 0.1 ); - -py::class_< Scene, std::shared_ptr< Scene > >( ps, "Scene" ) + py::class_< Scene, std::shared_ptr< Scene > >( ps, "Scene" ) .def(py::init<>()) .def( "setMesh", &Scene::setMesh ) - .def( "computeIR", &Scene::computeIR ) - .def( "computeMultichannelIR", &Scene::computeMultichannelIR ); + .def( "computeIR", py::overload_cast&, std::vector&, Context&>(&Scene::computeIR), + "A function to calculate IRs based on pre-defined sources and listeners", py::arg("_sources"), py::arg("_listeners"), py::arg("_context")) + .def( "computeIR", py::overload_cast>&, std::vector>&, Context&, float, float, float>(&Scene::computeIR), + "A function to calculate IRs based on source and listener locations", py::arg("_sources"), py::arg("_listeners"), py::arg("_context"), py::arg("src_radius") = 0.01, py::arg("src_power") = 1.0, py::arg("lis_radius") = 0.01 ); py::class_< SoundSource, std::shared_ptr< SoundSource > >( ps, "Source" ) .def( py::init>() ) diff --git a/tests/rir_test.py b/tests/rir_test.py index 99adb0b..441fcab 100644 --- a/tests/rir_test.py +++ b/tests/rir_test.py @@ -4,40 +4,87 @@ import numpy as np +def check_ir(samples): + assert np.argmax(np.fabs(samples)) > 0, "Max val is at beginning" + assert np.max(np.fabs(samples)) > 0, "IR max is zero" + assert ~np.isnan(samples).any(), "IR contains NAN" + + +def same_ir(ir1, ir2): + start_id = np.argmax(np.abs(ir1)) + assert start_id == np.argmax(np.abs(ir2)), "IRs have different start sample" + minlen = min([len(ir1), len(ir2), start_id + 1000]) + corrcoef_thresh = 0.9 + corrcoef = np.corrcoef(ir1[start_id:minlen], ir2[start_id:minlen])[0, 1] + assert corrcoef > corrcoef_thresh, "IRs are not similar" + + class MainTest(unittest.TestCase): - def test_rir(self): + def test_exception(self): + self.assertRaises(AssertionError, check_ir, [1, 0, 0]) + self.assertRaises(AssertionError, check_ir, [0, 0, 0]) + self.assertRaises(AssertionError, check_ir, [np.nan, 0, 0]) + self.assertRaises(AssertionError, same_ir, [1, 0, 0], [0, 1, 0]) + self.assertRaises(AssertionError, same_ir, [1, 0, 0], [1, 0.5, -0.5]) + + @staticmethod + def test_rir(): seed = 0 np.random.seed(seed) low = 0.5 high = 0.99 - N = 20 + N = 10 cnt = 0 - bad_cnt = 0 margin = 0.1 - bad_configs = [] roomdims = np.random.uniform(0.5, 10.0, (N, 3)).tolist() lis_locs = [[np.random.uniform(margin, x - margin) for x in roomdim] for roomdim in roomdims] src_locs = [[np.random.uniform(margin, x - margin) for x in roomdim] for roomdim in roomdims] alphas = 1 - np.random.uniform(low, high, (N, )) - while (cnt < N): + while cnt < N: alpha = alphas[cnt] - scatter = 0.5 - tasks = [[src_locs[cnt], lis_locs[cnt], 'seed{}_{}'.format(seed, cnt)]] - status = compute_scene_ir_absorb(roomdims[cnt], tasks, alpha) - self.assertEqual(status, 0) - if status: - bad_cnt += 1 - bad_configs.append([cnt] + roomdims[cnt] + src_locs[cnt] + lis_locs[cnt] + [alpha]) + tasks = [[src_locs[cnt], lis_locs[cnt]]] + compute_scene_ir_absorb(roomdims[cnt], tasks, alpha) cnt += 1 - print('{}/{} bad IR'.format(bad_cnt, cnt)) + + @staticmethod + def test_rir_pairs(): + roomdim = [10, 10, 10] + src_locs = [[0.5, 0.5, 0.5], [9.5, 9.5, 9.5]] + lis_locs = [[2.5, 0.5, 0.5], [5.0, 5.0, 5.0], [9.5, 0.5, 0.5]] + alpha = 0.5 + + mesh = ps.createbox(roomdim[0], roomdim[1], roomdim[2], alpha, 0.5) + + ctx = ps.Context() + ctx.diffuse_count = 20000 + ctx.specular_count = 2000 + ctx.threads_count = min(multiprocessing.cpu_count(), 8) + + scene = ps.Scene() + scene.setMesh(mesh) + + channel = ps.ChannelLayoutType.mono + ctx.channel_type = channel + ctx.sample_rate = 16000 + + src_lis_res = scene.computeIR(src_locs, lis_locs, ctx) + lis_src_res = scene.computeIR(lis_locs, src_locs, ctx) + + for i_src in range(len(src_locs)): + for i_lis in range(len(lis_locs)): + ir1 = src_lis_res['samples'][i_src][i_lis][0] + ir2 = lis_src_res['samples'][i_lis][i_src][0] + check_ir(ir1) + check_ir(ir2) + same_ir(ir1, ir2) # IRs should be similar by reciprocity def compute_scene_ir_absorb(roomdim, tasks, r): # Initialize scene mesh try: - mesh = ps.createbox(roomdim[0], roomdim[1], roomdim[2], r, 0.5) + mesh = ps.createbox(roomdim[0], roomdim[1], roomdim[2], r, 0.1) except Exception as e: print(str(e)) @@ -53,38 +100,18 @@ def compute_scene_ir_absorb(roomdim, tasks, r): ctx.channel_type = channel ctx.sample_rate = 16000 - status = 0 - for task in tasks: src_coord = task[0] lis_coord = task[1] - wavname = task[2] - cnt = 0 src = ps.Source(src_coord) src.radius = 0.01 + src.power = 1 lis = ps.Listener(lis_coord) lis.radius = 0.01 - # lis.channel_layout_type = ps.ChannelLayoutType.mono - - a = np.array(src_coord) - b = np.array(lis_coord) - dist = np.linalg.norm(a - b) - direct_idx = int(ctx.sample_rate * dist / 343) - res = scene.computeIR(src, lis, ctx) - res['samples'] = np.atleast_2d(res['samples']) - if (np.argmax(np.fabs(res['samples'])) == 0): - status = 1 - wavname += '_startzero.wav' - elif (np.max(np.fabs(res['samples'])) == 0): - status = 2 - wavname += '_zeromax.wav' - elif (np.isnan(res['samples']).any()): - status = 3 - wavname += '_nan.wav' - - return status + res = scene.computeIR([src], [lis], ctx) + check_ir(res['samples'][0][0][0]) if __name__ == "__main__":