-
Notifications
You must be signed in to change notification settings - Fork 3
/
mlescands.m
127 lines (113 loc) · 5.88 KB
/
mlescands.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
function [mleim] = mlescands(sigar, envpar, procpar)
% This function processes signals from a microphone array relative to a 3-D
% measurment space containing the microphones and the sound sources.
% The mic positions field of view grid point with any uint as long as
% the units for parameters listed below are consistent (I like meters
% and seconds).
% The output is the liklihood of the sound source as
% as function of the grid. If only a 2-D grid is input for the Mic positions
% and FOV, then only a 2-D grid at the output is computed.
%
% [mleim, fovxa, fovya] = mlescan(sigar, envpar, procpar)
%
% Inputs
% "sigar" matrix where each column is the rf signal from each mic.
% "envpar.mpos" 3-row matrix indicating the x,y,z cartesion positions of each
% mic, where the first row contains the x positions, second y
% position, and third z positions. If an axis is missing
% it will be assume to be the z and an array in the x-y
% plane will be assumed.
%
% "envpar.xg" array of grid points for x-axis (from meshgird)
% "envpar.yg" array of grid points for y-axis (from meshgird)
% "envpar.zg" array of grid points for z-axis (from meshgird), if this
% variable is missing, the function will only do a 2-D scan.
% "procpar.fs" scalar is the sampling frequency of the rf mic signals
% "envpar.c" scalar is the speed of sound in the room.
% "procpar.win" scalar the correlation window in time units. This should
% be less than the segment lenght of "sigar" or else zeros
% will be padded to allow for relative shifts between signals.
% "procpar.mxdelay" Maximimum delay to consider in forming beamfields
% "procpar.mprs" Indeces of all mic pairs to be combined in beamfield creation
% "procpar.combine" Sting indicating how pairwise mic correlation should be combined
% mean => average, min = Minimum, med => Median, max =>
% Maximum, absmed = Median of Absolute values
dimfov = size(envpar.xg); % Get dimensions for FOV
ndims = length(dimfov); % Check to see how many dimensions were provided
mleim = zeros(dimfov); % Initialize output image with Zeros
% Obtain mic array information
[mr, mc] = size(envpar.mpos); % Determine # number of mics = mc
[sr,sc] = size(sigar); % Determine the mic signal array size [# of samples, # of mics]
winsamp = procpar.win*procpar.fs; % Convert processing window length to samples
temp = zeros(winsamp,sc); % Initialze correlation matrix
dlysamp = procpar.mxdelay*procpar.fs; % Convert maximum delay length to samples
totwin = ceil(winsamp+dlysamp); % Compute maximum data window need to
% accomodate maximum delay for
% beamforming
% If maxdelay requires data not in array, pad with zeros
if totwin > sr;
zpadd = totwin-sr+1;
sigar = [sigar; zeros(zpadd, sc)];
end
% If only 2 dimension space is specified use double loop
if ndims == 2
for kx=1:dimfov(1);
for ky = 1:dimfov(2);
% Distance of position from microphones
ds = envpar.mpos - [envpar.xg(kx,ky); envpar.yg(kx,ky)]*ones(1,mc); % Spatial
rst = sqrt(ds(1,:).^2 + ds(2,:).^2)/envpar.c; % Corresponding Time differences
[mt, nmic] = min(rst); % find closest mic to point being considered
% Make all delays relative to closest mic, setting closest mic to
% 0 delay
rstref = rst-mt;
at = rst.^2;
at = (min(at) ./ (abs(at)));
% Convert initial and points to sample indices for each mic signal
indbeg = round(rstref*procpar.fs) + 1;
indend = winsamp + round(rstref*procpar.fs);
% Load up matrix with aligned mic signal corresponding to mic
% position
for ksig = 1:sc
temp(:,ksig) = at(ksig)*sigar(indbeg(ksig):indend(ksig),ksig); %.*hamming(winsamp);
%temp(:,ksig) = temp(:,ksig); %/(sqrt(ksig)*std(temp(:,ksig)));
end
% figure(1)
% imagesc(temp); colorbar
% pause
dumr =sum(temp');
mleim(kx,ky) = dumr*dumr';
end
end
%
elseif ndims == 3
% Check to see if mics have 3 indeces for position
mcs = size(envpar.mpos);
% If not assume mics are in a reference zero plane
if mcs < 3;
envpar.mpos = [envpar.mpos; zeros(1,msc(2))];
end
for kx=1:dimfov(1)
for ky = 1:dimfov(2)
for kz = 1:dimfov(3)
% Distance of position from microphones
ds = envpar.mpos - [envpar.xg(ky,kx,kz); envpar.yg(ky,kx,kz); envpar.zg(ky,kx,kz);]*ones(1,mc); % Spatial
rst = sqrt(ds(1,:).^2 + ds(2,:).^2 + ds(3,:).^2)/envpar.c; % Corresponding Time differences
[mt, nmic] = min(rst); % find closest mic to point being considered
% Make all delays relative to closest mic, setting closest mic to
% 0 delay
rstref = rst-mt;
at = rst.^2;
at = (min(at) ./ (abs(at)));
% Convert initial and points to sample indices for each mic signal
indbeg = round(rstref*fs) + 1;
indend = winsamp + round(rstref*fs) + 1;
% Load up matrix with aligned mic signal corresponding to mic
% position
for ksig = 1:sc
temp(:,ksig) = at(ksig)*sigar(indbeg(ksig):indend(ksig),ksig);
end
mleim(ky,kx,kz) = (mean(mean(temp(:,procpar.mprs(kmps,1)))));
end
end
end
end