1
% This script generates a sqrt(j k) pre-equalization filter
2
% for Wave Field Synthesis
5
% Deutsche Telekom Laboratories
8
% ===== parameters ========================================================
10
% sampling frequency of filter
13
% aliasing frequency of system (= upper frequency limit of equalization filter)
16
% lower frequency limit of filter (= frequency when subwoofer is active)
19
% number of coefficients for filter
23
% ===== variables =========================================================
24
f = linspace(0,f_s/2,441*10);
26
idx_al=max(find(f<f_al));
27
idx_sw=max(find(f<f_sw));
29
H = ones(1,length(f));
32
% -------------------------------------------------------------------------
33
% computation of pre-equalization filter
34
% -------------------------------------------------------------------------
37
H(1:idx_al) = sqrt(2*pi*f(1:idx_al))./sqrt(2*pi*f_al);
38
H(1:idx_sw) = H(idx_sw)*ones(1,idx_sw);
41
h = firls(Nfilt,2*f/f_s,H);
43
% truncate length to power of 2
46
% -------------------------------------------------------------------------
47
% plot resulting filter characteristics
48
% -------------------------------------------------------------------------
50
Hfilt=fftshift(fft(h));
51
Hfilt=Hfilt(Nfilt/2+1:end);
52
f2 = linspace(0,f_s/2,length(Hfilt));
55
plot(f,20*log10(abs(H)));
57
plot(f2,20*log10(abs(Hfilt)),'ro-');
61
xlabel('frequency -> [Hz]');
62
ylabel('magnitude response -> [dB]');
63
legend('desired response','filter response','Location','SouthEast');
64
axis([0 2*f_al -20 2]);
70
% -------------------------------------------------------------------------
71
% save filter coefficients
72
% -------------------------------------------------------------------------
73
fname = sprintf('wfs_prefilter_%d_%d_%d.wav',f_sw,f_al,f_s);
74
disp(['Wrote pre-equalization filter into file: ' fname]);
75
wavwrite(h,f_s,fname);