-
Notifications
You must be signed in to change notification settings - Fork 17
/
vbip.m
148 lines (127 loc) · 4.13 KB
/
vbip.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
function GainMtx = vbip(src_dirs, ls_groups, ls_invMtx, spread, num_spread_src, num_spread_rings3d)
%VBIP Computes vector-base intensity panning gains for a set of directions
%
% INPUTS:
%
% src_dirs: panning direction in degrees, vector for 2D, [Nsrc x 2] matrix
% for 3D, in [azi elev] convention
% ls_groups: valid pairs (for 2D) or triplets (for 3D triplets) returned
% by findLsPairs() or findLsTriplets()
% ls_InvMtx: matrix of loudspeaker inversions returned by invertLsMtx()
% spread: value of spread in degrees of the panning gains for MDAP
% Additional spreading parameters like number of spread sources and
% rings can be enabled easily, see this code and getSpreadSrcDirs()
% num_spread_src: see getSpreadSrcDirs()
% num_spread_rings3d: see getSpreadSrcDirs()
%
% OUTPUTS:
%
% GainMtx: (Nsrc x Nspeaker) VBIP gain matrix.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Archontis Politis, 1/11/2015
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if nargin<4
spread = 0;
num_spread_src = 0;
num_spread_rings3d = 0;
elseif nargin<5
num_spread_src = 8;
num_spread_rings3d = 1;
elseif nargin<6
num_spread_rings3d = 1;
end
src_num = size(src_dirs,1);
dim = size(ls_groups,2);
ls_num = max(ls_groups(:));
GainMtx = zeros(src_num, ls_num);
% 3D case
if dim == 3
for ns=1:src_num
gains = max(0, vbip3(src_dirs(ns,1), src_dirs(ns,2), ls_groups, ls_invMtx, spread, num_spread_src, num_spread_rings3d));
GainMtx(ns,:) = gains;
end
% 2D case
elseif dim == 2
for ns=1:src_num
gains = max(0, vbip2(src_dirs(ns), ls_groups, ls_invMtx, spread, num_spread_src));
GainMtx(ns,:) = gains;
end
end
end
%%%%%% COMPUTE GAIN FACTORS CORRESPONDING TO INPUT ANGLE, 2D
function gains = vbip2(azi, ls_groups, ls_invMtx, spread, num_spread_src)
if nargin<4, spread = 0; end
ls_num = max(ls_groups(:));
energies = zeros(1,ls_num);
if spread
U_spread = getSpreadSrcDirs(azi, spread, num_spread_src);
for ns = 1:size(U_spread, 1);
u_ns = U_spread(ns,:);
e_ns = zeros(1,ls_num);
for i=1:size(ls_groups,1);
e_tmp(1) = ls_invMtx(i,1:2) * u_ns';
e_tmp(2) = ls_invMtx(i,3:4) * u_ns';
if min(e_tmp) > -0.001
e_ns(ls_groups(i,:)) = e_tmp/sum(e_tmp);
end
end
energies = energies + e_ns;
end
else
azi_rad = azi*pi/180;
u = [cos(azi_rad) sin(azi_rad)];
energies = zeros(1,ls_num);
for i=1:size(ls_groups,1);
e_tmp(1) = ls_invMtx(i,1:2) * u';
e_tmp(2) = ls_invMtx(i,3:4) * u';
if min(e_tmp) > -0.001
energies(ls_groups(i,:)) = e_tmp/sum(e_tmp);
end
end
end
energies=energies./sum(energies);
gains = sqrt(energies);
end
function gains = vbip3(azi, elev, ls_groups, ls_invMtx, spread, num_spread_src, num_spread_rings3d)
if nargin<5, spread = 0; end
ls_num = max(ls_groups(:));
energies = zeros(1,ls_num);
if spread
U_spread = getSpreadSrcDirs([azi elev], spread, num_spread_src, num_spread_rings3d);
for ns = 1:size(U_spread, 1);
u_ns = U_spread(ns,:);
e_ns = zeros(1,ls_num);
for i=1:size(ls_groups,1);
e_tmp(1) = ls_invMtx(i,1:3) * u_ns';
e_tmp(2) = ls_invMtx(i,4:6) * u_ns';
e_tmp(3) = ls_invMtx(i,7:9) * u_ns';
if min(e_tmp) > -0.001
e_ns(ls_groups(i,:)) = e_tmp/sum(e_tmp);
energies = energies + e_ns;
break
end
end
end
else
azi_rad = azi*pi/180;
elev_rad = elev*pi/180;
u = [cos(azi_rad)*cos(elev_rad) sin(azi_rad)*cos(elev_rad) sin(elev_rad)];
energies = zeros(1,ls_num);
for i=1:size(ls_groups,1);
e_tmp(1) = ls_invMtx(i,1:3) * u';
e_tmp(2) = ls_invMtx(i,4:6) * u';
e_tmp(3) = ls_invMtx(i,7:9) * u';
if min(e_tmp) > -0.001
energies(ls_groups(i,:)) = e_tmp/sum(e_tmp);
break
end
end
end
energies = energies/sum(energies);
energies(energies<0) = 0;
gains = sqrt(energies);
end