背景有点儿
我希望能够播放尽可能接近我已校准的录音的声音。据我了解的理论，我需要找到扬声器的传递函数，并解卷积我想随其发出的信号。像这样（在频域中）：
X -> H -> XH
其中
X是发射信号H是扬声器传递函数，而XH是X乘以H。现在，除法（./）应该给我H。现在，要发出校准信号，应将其除以
H：X/H -> H -> X
完成的操作
在三脚架上相距1 m放置了扬声器和已校准的麦克风。
记录了30+线性扫描150KHz-20KHz，长20ms，并以500记录了KS / s。
使用下面的Matlab / Octave脚本对齐和平均后的信号，在该脚本下可以看到结果信号。
files = dir('Mandag*');
rng = [1.5e6, 1.52e6];
[X, fs] = wavread(files(1).name, rng);
X = X(:,1);
for i=2:length(files)
[Y, fs] = wavread(files(i).name, rng);
sig = Y(:,1);
[x, off] = max(xcorr(X', sig'));
off = length(X) - off;
if(off < 0)
sig = [zeros(1, -off), sig(1:end+off)'];
elseif (off > 0)
sig = [sig(off:end)', zeros(1, off-1)];
end
X = X + sig';
end
X = X/length(files);
XH并进行了上述计算，结果似乎是合理的。下面是H（紫色）和X/H（绿色）的归一化图。
该图已被截断为相关频率。
请让我知道我是否以错误的方式进行操作。域，我假设这将是一个简单的
X/H和ifft(X./H)，但是到目前为止，我的所有尝试都未能获得任何合理的答案。可以在此处以mat7-binary格式找到频率向量wavwrite，Hf和H。也许我只是累了，这里有一个简单的解决方案，但此刻我看不到它。任何帮助/建议都非常感谢。
#1 楼
在查看吉姆·克莱（Jim Clay）在评论中提到的参考文献后找到了答案，谢谢吉姆。我犯了一个错误，仅考虑幅度会导致零相位信号并且无法合理使用对于发射，至少在此设置中没有。
我最终使用的代码可以在下面看到。
该脚本遵循保留时域信号的命名约定
% Align and sum all files called Mandag*
files = dir('Mandag*');
% Where in the recordings the signal is
rng = [1.5e6, 1.52e6];
% Initialize the xh vector
[xh, fs] = wavread(files(1).name, rng);
xh = xh(:,1);
for i=2:length(files)
y = wavread(files(i).name, rng);
y = y(:,1);
% Determine offset between xh and y
[~, off] = max(xcorr(xh', y'));
off = length(xh) - off;
% Shift signal appropriately
if(off < 0)
y = [zeros(1, -off), y(1:end+off)'];
elseif (off > 0)
y = [y(off:end)', zeros(1, off-1)];
end
xh = xh + y';
end
% Average
xh = xh/length(files);
% Location of the 20ms signal
xh = xh(2306:12306-1);
% Normalize
xh = xh / max(xh);
% Apply a moving average filter on xh to reduce noise. Window size of 4 was
% experimentally determined to give the best results
n = 4;
B = zeros(n, 1);
for i=1:n
B(i) = 1/n;
end
xh = filter(B, 1, xh);
xh = xh / max(xh);
x = wavread('sweep.wav');
x = x(1:2:end); % Sweep generated @ 1MHz, decimate
% to have same length as xh
% Transform x into frequency domain and determine H
X = fft(x);
H = fft(xh) ./ X;
% Vector indices to choose only frequencies of interest
starti = 20e3 / 50;
endi = 100e3 / 50;
rng = starti:endi;
irng = (length(x) - endi) : (length(x) - starti);
% Zero out unwanted frequencies
X = [zeros(1, starti - 1 ), X( rng)', zeros(1, length(X)/2 - endi) ...
zeros(1, length(X)/2 - endi), X(irng)', zeros(1, starti - 1 )]';
% Deconvolve x with h
X_deconv_H = X ./ H;
% Transform X/H to time domain and normalise
x_deconv_h = real(ifft(X_deconv_H));
x_deconv_h = x_deconv_h / max(x_deconv_h);
% Save the deconvolved sweep
wavwrite(x_deconv_h, fs, 'deconvolved_sweep.wav');
% Generate spectrograms of xh and x_deconv_h
winsize = 512;
overlap = round(.99 * winsize);
figure(1)
specgram(xh, winsize, fs, hann(winsize), overlap)
colorbar
figure(2)
specgram(x_deconv_h, winsize, fs, hann(winsize), overlap)
colorbar
x conv h的频谱图如下所示：
尽管去卷积信号中有一些噪声，但这些对我来说似乎很合理。
下一个测试将看是否发射
x deconv h除非扬声器无法发出那些频率，否则类似x_deconv_y。用测试结果更新
我们使用对数下扫来重述上述测量。这些结果似乎表明该方法可行。
验证测试包括发射
x并期望取回X / H，即在所有频率下能量相等。由于最差的输出频率约为比最好的要弱20dB，最高的输出电平会低得多。
发出的信号：
所记录信号的时间序列和频谱图如下所示：
评论
$ \ begingroup $
对此有任何更新吗？您发出的信号看起来如何？
$ \ endgroup $
–兰斯
2012年8月1日14:44
$ \ begingroup $
@Busk：感谢您的关注。由于该设备正在其他地方使用，因此我还没有机会对其进行测试。完成测试后，我将发布结果。
$ \ endgroup $
–雷神
2012年8月1日在16:04
$ \ begingroup $
@Busk：我们现在已经对其进行了测试，它似乎可以工作:-)。
$ \ endgroup $
–雷神
2012年10月31日8:33
评论
这个线程-dsp.stackexchange.com/questions/953/…-和这个dsp.stackexchange.com/questions/2705/…-对您可能有用。是的，发现我的错误，谢谢吉姆。我只考虑信号的幅度，这导致了零相位时间信号。看来我现在得到了正确的结果，我将其添加为答案。