总共有四个函数分别如下,导入采集到的音频文件,设置好代码读取路径就可以实现了。
第一个函数:采集系统
% In this program, the signals are divided in pieces of lengData samples
% and process separately. This program calls two functions:
% f_adap
% get_angle
% At the end it returns the direction in degrees.
% CONSTANT VALUES
i=0;
fs=44100; % Sampling frequency (Hz)
d_micro=0.1; % Distance between microphones (m)
c=340; % Speed of sound (m/s)
muestrasMAX=ceil(d_microfs/c); % Maximum number of samples Nmax
DESP=ceil(muestrasMAX1.5); % Delay we insert in the micro 2
% We leave 50% of margin of error.
lengData=441000.5; % Number of samples in which the voice
% signals are divided to be processed.
% All values depend on fs and d_micro in
% case it was necessary to change them.
% Importing the file to process.
signal1=audioread(‘60.wav’)’;% LOOP WHERE THE DIFFERENT PARTS OF THE IMPORTED FILE ARE PROCESSED
for k=lengData:lengData:length(signal1)
tic; % Measure of time tic;toc;
signal=signal1(1,k-(lengData-1):k); % Signal in MIC B
d=signal1(2,k-(lengData-1):k); % Signal in MIC A
% NORMALIZATION PROCESS
M1=max(abs(signal)); % Maximum of channel 1
M2=max(abs(d)); % Maximum of channel 2
M3=max(M1,M2); % Normalization value
signal=signal/M32; % Normalizing
d=d/M32;
% LMS ALGORITHM
hDESP=[zeros(1,DESP) 1]; % Filter to delay the signal DESP samples.
d1=conv(hDESP,d);
P=50; % Parameters of the algorithm
mu=0.0117;
h0=zeros(1,P); h0(1)=0; % Initialazing the adaptative filter
[h y e]=f_adap(signal,d1,h0,mu); % Recursive function calculating the
% coefficients of the filter h(n)
% PROCESSING THE FILTER BEFORE THE FREQUENCY ANALYSIS.
h1=[zeros(1,DESP-muestrasMAX-3),h(DESP-muestrasMAX-2:length(h))];
h1(DESP+muestrasMAX+2:length(h1))=0;
h1(DESP+1)=h1(DESP+1)/2;
[B,I]=sort(h1,‘descend’);
H1=[zeros(1,I(1)-3),h(I(1)-2:I(1)+2),zeros(1,length(h)-(I(1)+2))];
% FREQUENCY ANALYSIS TO OBTAIN THE DELAY (IN SAMPLES)
% 1-FFT
lh=128; % Length of the FFT
H=fft(h1,lh); % FFT of the filter h(n)
% 2-ANGLE(+UNWRAP)
alpha=angle(fftshift(H)); % Obtaining the phase
q=unwrap(angle(fftshift(H)));
% 3-SLOPE
M=diff(q); % Obtaining the slope of the phase
% 4-SLOPE’S AVERAGE
lM=length(M)+2; % The slope M1 is not a unique value,
p1=floor(lM/2-4); % it’s an array. So we calculate the
p2=ceil(lM/2+4); % average of the values, K.
K=mean(M(p1:p2));
Nprime=(-Klh/(2*pi)); % Number of samples before
% substracting DESP.
% 5-SAMPLES
if Nprime<0 % Two possible cases: negative or positive
N=Nprime+lh;
N=N-DESP;
else
N=Nprime;
N=N-DESP;
end
% CALLING THE FUNCTION WHICH RETURNS THE ANGLE
angleGRAD1=get_angle(N,fs,d_micro);
if isreal(angleGRAD1)==1 % Security measures in case
angleGRAD1 % the number is complex
i=i+1;
else
angleGRAD1=real(angleGRAD1)
i=i+1;
end
timeElapsed=toc; % Time is kept in variable timeElapsed
timeElapsed;
i
end