Zoom H2で音源方向認識

Zoom H2というのは4chのPCM録音ができるポータブル・レコーダーなのですが、せっかく4つのマイク位置が固定になっているので、それを使って遊んでみようと思い、Matlabで非常に単純な方法で音源方向の認識を行うMatlabルーチンを作ってみました。H2には前方は90度に、後方は120度に開いたマイクが搭載されています。便宜的に、全てのマイクが同じ位置に置かれた理想的な単一指向マイクだとして、音源到来方向を4chに録音された音の音圧比から求めます。

途中にマイクの指向性を描いたり(デバッグのため)、窓関数をかけたり(メインループ部分をSTFTルーチンから持ってきた名残)、いろいろと無駄な部分がありますが、マラカス単一音源で試してみたら、ちゃんと方向を示してくれました。せっかくだから、これをMax/MSPで拡張して音源方向を向くマイク※とか話者にカメラを向ける三脚※とかを作ろうかと思ったのに、Zoom H2をコンピュータにUSB接続すると2chマイクとしてしか機能してくれませんでした。残念。

※こういうアイディアは、何年も前から真剣に研究している人々がいます。話者にピントを合わせたり、その話者の声以外を消したり、話者が複数人でも大丈夫だったり・・・日本音響学会にもいろいろと面白いデモが出ています。

% sound direction analyzer for Zoom H2 four channel recording
% (good for only one source)
%
% Four cardioid microphones are placed inside a housing enclosure
% coincidently, with angles between front and rear two microphones are 90
% and 120 degrees.  Let center be 0 degrees, and then FL=45, FR=285,
% RL=120, RR=240 degrees (counter-clock-wise from the center-front).

% 2008-09-11 by JikanBae



%% read sound files
[xF,fs] = wavread('SR000F.wav');
[xR,fs] = wavread('SR000R.wav');
x = [xF xR];
fc = 1000/(fs/2);
[N, Wn] = buttord([fc/(2^(1/2)) fc*(2^(1/2))], [fc/(2^(1/1)) fc*(2^(1/1))], 3, 60);
[B, A] = butter(N, Wn);
[B, A] = butter(6, [707 1414]/(fs/2));
x = filter(B, A, [xF xR]);



%% make a look-up-table of angles
% a=2 for hypocardioid,
% 1 for cardioid
% 0.5 for supercardioid,
% 0.33 for hypercardioid,
% 0 for bi-directional
a = 1.0;   b = 1.0;
th = [0:360]' / 180 * pi;
th = th(1:end-1);
rFL = (a + cos(th -  45/180*pi).^b) / (1+a);
rFR = (a + cos(th +  45/180*pi).^b) / (1+a);
rRL = (a + cos(th - 120/180*pi).^b) / (1+a);
rRR = (a + cos(th + 120/180*pi).^b) / (1+a);
figure(1);
subplot(2,2,1); polar(th, rFL);
subplot(2,2,2); polar(th, rFR);
subplot(2,2,3); polar(th, rRL);
subplot(2,2,4); polar(th, rRR);
r = [rFL rFR rRL rRR];






%% prepare variables for windowing
winTime = 300; % in ms (CHANGE HERE)
winSize = 2 ^ (nextpow2(winTime/1000 * fs));
totalSize = ceil(size(x,1) / winSize) * winSize;

%% set sample length to be nice
x = [x; zeros(totalSize - size(x,1), size(x,2))];
y = zeros(size(x));

%% separate into frame and do process
w = hann(winSize, 'periodic');
tic
k = 1;
amp = zeros(totalSize / winSize, 4);
amp2 = zeros(totalSize / winSize, 4);
ths = zeros(totalSize / winSize, 1);
for m = 1:winSize/2:(totalSize-winSize/2)
  zin = x(m:m+winSize-1,:);
  amp(k,:) = rms(zin);
  % find the best match with the look-up-table
  [val, ind] = max(abs(r / amp(k,:)));
  ths(k) = th(ind);
  
  k = k+1;
end
toc

% unwrap angles and display as a graph
thsa = unwrap(ths) / pi*180;
thsa = mod(thsa+180, 360)-180;
figure(2);
plot(0:(winSize/2)/fs:(totalSize-winSize/2-1)/fs, thsa, 'k.');
set(gca, 'XLim', [0 length(x)-1]/fs);
set(gca, 'YLim', [-180 180]);
set(gca, 'YTick', [-180:30:180]);
xlabel('Time (sec)');
ylabel('Incident Angle (deg)');
grid on;

%% polar plot of sound movement
figure(3);
ma = mean(amp,2);
comet(ma.*cos(unwrap(ths)), ma.*sin(unwrap(ths)));