总共有四个函数分别如下，导入采集到的音频文件，设置好代码读取路径就可以实现了。

第一个函数：采集系统

% 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