Skip to content

Commit da6de06

Browse files
feliximmohrhagenw
feliximmohr
authored and
hagenw
committed
Add voronoi interpolation to IR processing (#171)
1 parent ac82079 commit da6de06

File tree

5 files changed

+811
-82
lines changed

5 files changed

+811
-82
lines changed

‎SFS_config.m

Lines changed: 20 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -384,18 +384,28 @@
384384
% Settings regarding all the stuff with impulse responses from the SFS_ir and
385385
% SFS_binaural_synthesis folders
386386
%
387-
% Use interpolation to get the desired HRTF or BRIR for binaural simulation. If this
388-
% is disabled, the HRTF/BRIR returned by a nearest neighbour search is used instead.
387+
% Use interpolation to get the desired HRTF or BRIR for binaural simulation.
388+
% If this is disabled, the HRTF/BRIR returned by a nearest neighbour search is
389+
% used instead.
389390
conf.ir.useinterpolation = true; % boolean
390-
% You can choose the way the points for interpolation are selected. Depending on the
391-
% geometry of the measured HRTF/BRIR data set, the interpolation will be done between
392-
% two or three HRTFs. Available methods:
393-
% 'nearestneighbour' - Interpolation between nearest neighbours. This only works
394-
% for interpolation points on a circle in the horizontal plane.
391+
% You can choose the way the points for interpolation are selected. Depending on
392+
% the geometry of the measured HRTF/BRIR data set and the interpolation method,
393+
% a different number of HRTFs are selected for interpolation. See validation
394+
% script test_interpolation_point_selection.m for examples. Available methods:
395+
% 'nearestneighbour' - Interpolation between nearest neighbours. This only
396+
% works for interpolation points on a circle in the
397+
% horizontal plane.
395398
% 'delaunay' - Interpolation between surrounding points according to
396-
% Delaunay triangulation. This only works for interpolation
397-
% points on a sphere. See validation script
398-
% test_interpolation_point_selection.m for examples.
399+
% Delaunay triangulation. This only works for
400+
% interpolation points on a sphere.
401+
% 'voronoi' - For spherical Voronoi interpolation, all coordinates
402+
% are projected on the unit sphere. Voronoi regions
403+
% (including their area) on the sphere are calculated
404+
% for all measured positions included in the HRTF/BRIR
405+
% data set. The query position is added to the point
406+
% cloud and the regions are calculated, again.
407+
% The weights result from the area "stolen" by the query
408+
% position from each of the other coordinates
399409
conf.ir.interpolationpointselection = 'nearestneighbour';
400410
% You can choose between the following interpolation methods:
401411
% 'simple' - Interpolation in the time domain performed samplewise. This

‎SFS_general/check_dimensionality.m

Lines changed: 104 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,104 @@
1+
function [x0, xs, dim, eq_idx] = check_dimensionality(x0, xs, tol, gamma)
2+
%CHECK_DIMENSIONALITY checks dimensionality and rotates the grid to its
3+
%principal axes.
4+
%
5+
% Usage: [x0,xs,dim,eq_idx] = check_dimensionality(x0,xs,tol)
6+
%
7+
% Input parameters:
8+
% x0 - point cloud on a sphere around the origin / m [nx3]
9+
% xs - desired direction as point in R^3 / m [1x3]
10+
% tol - tolerance for equality check (default: 1e-6)
11+
% gamma - scalar in 0 < gamma << 1 (default: 0.1)
12+
%
13+
% Output parameters:
14+
% x0 - original point cloud rotated to principal axes / m [nx3]
15+
% xs - original query point rotated to principal axes / m [1x3]
16+
% dim - dimensionality, either 1, 2, 2.5 or 3
17+
% eq_idx - indice of point from point cloud equal to query point
18+
% within tolerance specified in tol or -1 otherwise
19+
%
20+
% 1D: xs is colinear with or equal to one x0
21+
% 2D: all x0 and xs are in one plane
22+
% 2.5D: all x0 are in one plane, but xs is not
23+
% 3D: otherwise
24+
%
25+
% See also: findvoronoi, test_interpolation_point_selection
26+
27+
%*****************************************************************************
28+
% The MIT License (MIT) *
29+
% *
30+
% Copyright (c) 2010-2018 SFS Toolbox Developers *
31+
% *
32+
% Permission is hereby granted, free of charge, to any person obtaining a *
33+
% copy of this software and associated documentation files (the "Software"), *
34+
% to deal in the Software without restriction, including without limitation *
35+
% the rights to use, copy, modify, merge, publish, distribute, sublicense, *
36+
% and/or sell copies of the Software, and to permit persons to whom the *
37+
% Software is furnished to do so, subject to the following conditions: *
38+
% *
39+
% The above copyright notice and this permission notice shall be included in *
40+
% all copies or substantial portions of the Software. *
41+
% *
42+
% THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
43+
% IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, *
44+
% FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL *
45+
% THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
46+
% LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING *
47+
% FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER *
48+
% DEALINGS IN THE SOFTWARE. *
49+
% *
50+
% The SFS Toolbox allows to simulate and investigate sound field synthesis *
51+
% methods like wave field synthesis or higher order ambisonics. *
52+
% *
53+
% http://sfstoolbox.org sfstoolbox@gmail.com *
54+
%*****************************************************************************
55+
56+
57+
%% ===== Checking of input parameters ==================================
58+
nargmin = 2;
59+
nargmax = 4;
60+
narginchk(nargmin,nargmax);
61+
if nargin < nargmax
62+
gamma = 0.1; % inverse of aspect ratio of principal axes
63+
end
64+
if nargin < nargmax-1
65+
tol = 1e-6; % in case no tolerance is provided as input arg
66+
end
67+
68+
69+
%% ===== Computation ====================================================
70+
71+
% In case no 1D or 2D case is detected
72+
dim = 3;
73+
74+
% Check for 1D case (equality\collinearity within tolerance) and return idx
75+
eq = vector_norm(bsxfun(@minus,x0,xs),2);
76+
eq_idx = find(eq<=tol);
77+
if ~isempty(eq_idx)
78+
dim = 1;
79+
warning('SFS:check_dimensionality', ...
80+
'%s: Query point is apparently colinear with or equal to one grid point.', ...
81+
upper(mfilename))
82+
return
83+
else
84+
eq_idx = -1;
85+
end
86+
87+
% Check for 2D or 2.5D case and rotate to principal axes
88+
x0mean = mean(x0,1);
89+
[~,S,V] = svd(bsxfun(@minus,x0,x0mean));
90+
x0mean = x0mean*V;
91+
x0 = x0*V;
92+
xs = xs*V;
93+
S = diag(S);
94+
if S(end)/S(end-1) < gamma
95+
dim = 2.5;
96+
warning('SFS:check_dimensionality',...
97+
'%s: Grid is apparently two-dimensional. ',upper(mfilename));
98+
% Check for 2D case
99+
if abs(xs(3) - x0mean(3)) < tol
100+
dim = 2;
101+
x0(:,3) = x0(:,3) - x0mean(3);
102+
xs(3) = xs(3) - x0mean(3);
103+
end
104+
end

0 commit comments

Comments
 (0)