📄 configurewhereisit.m
字号:
function configurewhereisit(fo)
%CONFIGUREWHEREISIT Configuration for WHEREISIT.
% CONFIGUREWHEREISIT will display a series of dialogs which will prompt
% you to place a speaker at a specified location. Once completed,
% the resulting data is displayed.
%
% configurewhereisit uses a single tone at 500 Hz
% configurewhereisit(fo) uses a single tone at frequency fo Hz
% configurewhereisit(wavfilename) plays the contents of wavfilename instead
%
% See also WHEREISIT, MAKEASOUND, TRACKER.
% Loren Dean (original version)
% Scott Hirsch (this version)
filename=''; %Name of wav file for sound source (default: single tone)
if nargin==0
fo = 500;
elseif isstr(fo) & exist(fo,'file'); %Specify file name
filename=fo;
%We still want to know the fundamental frequency. We'll guess from spectral
% analysis of the contents of the wav file
[y,fs]=wavread(filename);
[Pxx,F]=psd(y(:,1),min([length(y) 1024*4]),fs);
[P_max,Max_ind] = max(Pxx);
%watch out for big DC offset
if Max_ind<4 %Any of the first few frequencies
[P_max,Max_ind] = max(Pxx(4:end));
Max_ind=Max_ind+3;
end;
fo = F(Max_ind);
end;
%Grab a handle to the Tracker window to find out which products we have
set(0,'ShowHiddenHandles','on');
tracker_fig = findobj('Name','Tracker');
tracker_handles = guidata(tracker_fig);
has = tracker_handles.has;
set(0,'ShowHiddenHandles','off');
%Create the figure
fh = figure('Name','Configure Acoustic Tracker', ...
'NumberTitle','off');
h=msgbox(['We will now configure this experiment by placing the ' ...
'speaker at two different locations. You will be ' ...
'prompted with three more dialogs. Follow the ' ...
'instructions given with each dialog.'], ...
'Where is it? Configuration');
waitfor(h)
% Get the data when the speaker is placed next to the left microphone
[data(1,:),trace]=LRetrieveData( ...
['Please place the speaker next to the left microphone and ' ...
'press OK when it is there. A tone will sound for 1 second ' ...
'after which you will be prompted with another dialog.'], ...
'Left Microphone',fo,filename);
%Plot the data - top left of figure
time = linspace(0,1,length(trace));
subplot(4,3,1);
plot(time,trace(:,1));
ylabel('Left Mic');
title('Left Tone');
ax = axis;
subplot(4,3,4);
plot(time,trace(:,2));
ylabel('Right Mic');
xlabel('Time (s)');
axis(ax);
if data(1,1) < data(1,2) %Mics are switched
errordlg('Your mics are crossed. Please switch the left and right mics.');
%close(fh)
return
end;
% Get the data when the speaker is placed next to the right microphone
[data(2,:),trace]=LRetrieveData( ...
['Please place the speaker next to the right microphone and ' ...
'press OK when it is there. A tone will sound for 1 second ' ...
'after which you will be prompted with another dialog.'], ...
'Right Microphone',fo,filename);
%Plot the data - top right
subplot(4,3,3);
plot(time,trace(:,1));
title('Right Tone');
axis(ax)
subplot(4,3,6);
plot(time,trace(:,2));
xlabel('Time (s)');
axis(ax);
% Get the data when the speaker is placed in the middle. We'll use this
% for mic calibration and for placement guess
[data(3,:), trace]=LRetrieveData( ...
['Please place the speaker mid-way between the two microphones.' ...
'A tone will sound for 1 second ' ...
'after which you will be prompted with another dialog.'], ...
'Center',fo,filename);
%Plot the data - top center
subplot(4,3,2);
plot(time,trace(:,1));
title('Center Tone');
axis(ax)
subplot(4,3,5);
plot(time,trace(:,2));
xlabel('Time (s)');
axis(ax);
%Calibrate the microphones so that they have about the same respone
%There are a bunch of ways we could do this with the data collected
%The ideal would be to measure the response of both mics subjected
% to an identical sound.
%If the third cal is done with the speaker truly in the middle, this is
% a good way to go.
if 0
cal = data(3,:);
end;
%We could also use the first two tones, assuming that the speaker is held
% in about the same position relative to each mic for the two cals.
cal = data([1 5]); %First mic, first test; second mic, second test
%Apply calibration factor to mic in second channel.
data(:,2) = data(:,2)*cal(1)/cal(2);
%Fit 1/r curve to data
x=[0 .5 1]';
y = data([1 3 2],:);
if has.cfit %Has curve fitting toolbox
warning off %Hide warnings about not specifying the start point
LeftMic = fit(x,y(:,1),'rat01'); %rat01: rational, 0 order num, 1st order den
RightMic = fit(x,flipud(y(:,2)),'rat01'); %Use same 1/r, then reverse
warning on
%Plot results of curve fit
subplot(2,2,3);
plot(LeftMic,x,y(:,1));
legend('data','1/r fit')
ylabel('Left Mic');
subplot(2,2,4);
plot(RightMic,x,flipud(y(:,2)));
legend('data','1/r fit',2)
ylabel('Right Mic');
set(gca,'XDir','reverse');
else
%People who don't have the curve fitting toolbox make my life more difficult!
%Try to fit y = p1/(x+q1)
%To use MATLAB's polyfit:
% y = pl / (x + q1)
% y_inv = 1/y = (x + q1) / p1
% y_inv = (1/p1)*x^1 + (q1/p1)*x^0
% y_inv = c1*x^1 + c0*x^0 %Fit this
% then, p1 = 1/c1;
% q1 = co/c1;
C_Left = polyfit(x,1./y(:,1),1);
C_Right = polyfit(x,1./flipud(y(:,2)),1);
p1_Left = 1/C_Left(1);
q1_Left = C_Left(2) / C_Left(1);
p1_Right = 1/C_Right(1);
q1_Right = C_Right(2) / C_Right(1);
%Build LeftMic, RightMic structures
LeftMic.p1 = p1_Left;
LeftMic.q1 = q1_Left;
RightMic.p1 = p1_Right;
RightMic.q1 = q1_Right;
%Plot results of curve fit
x_fit = linspace(0,1);
LeftMic_fit = 1./polyval(C_Left,x_fit);
RightMic_fit = 1./polyval(C_Right,x_fit);
subplot(2,2,3);
plot(x,y(:,1),'r.',x_fit,LeftMic_fit,'b');
legend('data','1/r fit')
ylabel('Left Mic');
subplot(2,2,4);
plot(x,flipud(y(:,2)),'r.',x_fit,RightMic_fit,'b');
legend('data','1/r fit',2)
ylabel('Right Mic');
set(gca,'XDir','reverse');
end;
% Save the data
save whereisitconfiguration data fo LeftMic RightMic cal filename
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%% LRetrieveData %%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%
function [magnitude,data]=LRetrieveData(msgString,msgTitle,fo,filename)
% Function to make a sound, retrieve data for 1 second through two
% microphones and then determine the magnitude of the incoming data.
% Display a message box to prompt the user where to place the speaker
h=msgbox(msgString,msgTitle);
waitfor(h)
% Start playing a sound
if exist(filename,'file')
makeasound(filename);
else
makeasound(fo);
end;
% Create an input object with two channels and configure it to
% run for 1 second
AI=analoginput('winsound');
addchannel(AI,1:2);
AI.SampleRate=44100; % Hz
AI.SamplesPerTrigger=AI.SampleRate*2; % Run for 1 second
% Start the object and get the data out from its buffer
% The object will automatically stop when it has acquired the
% number of samples specified by SamplesPerTrigger above.
start(AI)
data=getdata(AI);
% Stop playing the sound and delete the input object
makeasound stop
delete(AI)
% Determine the magnitude of the sound but only look at the data in the
% middle of the dataset so that any transients are removed.
range=length(data)/4:3/4*length(data);
magnitude = localRMS(data(range,1:2));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function out=localRMS(in,dim);
%RMS Compute rms value of a steady signal
[nr,nc]=size(in);
if nargin<1
help(mfilename);
elseif nargin<2
[junk,dim]=max([nr nc]);
end;
out=sqrt(mean(in.^2,dim));
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -