|
@@ -3,13 +3,19 @@
|
|
|
%% sphear_fir_render_directions: render impulse responses throught the FIR filter A2B matrix
|
|
%% sphear_fir_render_directions: render impulse responses throught the FIR filter A2B matrix
|
|
|
%%
|
|
%%
|
|
|
%% usage:
|
|
%% usage:
|
|
|
-%% [] = sphear_fir_render_directions(FIR, fmin, fmax, N, Npre, Fs, smooth)
|
|
|
|
|
|
|
+%% [] = sphear_fir_render_directions(FIR, fmin, fmax, N, Npre, Nfft, plot_on, smooth)
|
|
|
%%
|
|
%%
|
|
|
%% FIR: A2B matrix of FIR filters
|
|
%% FIR: A2B matrix of FIR filters
|
|
|
%% fmin, fmax: range of frequencies to use
|
|
%% fmin, fmax: range of frequencies to use
|
|
|
%% N: number of samples to filter
|
|
%% N: number of samples to filter
|
|
|
%% Npre: number of samples before the main impulse peak
|
|
%% Npre: number of samples before the main impulse peak
|
|
|
-%% Fs: sampling rate
|
|
|
|
|
|
|
+%% Nfft: size of FFT
|
|
|
|
|
+%%
|
|
|
|
|
+%% optional parameters:
|
|
|
|
|
+%%
|
|
|
|
|
+%% plot_on: plot frequency response of components, default is "epsc"
|
|
|
|
|
+%% can be one of 0, "eps", "epsc" (color), "png", "pngc" (color)
|
|
|
|
|
+%% smooth: number of samples to average, default is 0
|
|
|
%%
|
|
%%
|
|
|
%% Copyright 2016, Fernando Lopez-Lezcano, All Rights Reserved
|
|
%% Copyright 2016, Fernando Lopez-Lezcano, All Rights Reserved
|
|
|
%% nando@ccrma.stanford.edu
|
|
%% nando@ccrma.stanford.edu
|
|
@@ -27,7 +33,19 @@
|
|
|
%% You should have received a copy of the GNU General Public License
|
|
%% You should have received a copy of the GNU General Public License
|
|
|
%% along with this program. If not, see <http://www.gnu.org/licenses/>.
|
|
%% along with this program. If not, see <http://www.gnu.org/licenses/>.
|
|
|
|
|
|
|
|
-function [BFF] = sphear_fir_render_directions(FIR, fmin, fmax, N, Npre, Fs, smooth)
|
|
|
|
|
|
|
+function [BFF] = sphear_fir_render_directions(FIR, fmin, fmax, N, Npre, Nfft, plot_on, smooth)
|
|
|
|
|
+
|
|
|
|
|
+ %% no smoothing or plotting by default
|
|
|
|
|
+ if (!exist("smooth"))
|
|
|
|
|
+ smooth = 0;
|
|
|
|
|
+ end
|
|
|
|
|
+ if (!exist("plot_on"))
|
|
|
|
|
+ plot_on = "epsc";
|
|
|
|
|
+ else
|
|
|
|
|
+ if (plot_on == 1)
|
|
|
|
|
+ plot_on = "epsc"
|
|
|
|
|
+ end
|
|
|
|
|
+ end
|
|
|
|
|
|
|
|
[ref, Fs] = wavread('../cal/2016.04.20/export/reference-001-1.wav');
|
|
[ref, Fs] = wavread('../cal/2016.04.20/export/reference-001-1.wav');
|
|
|
%% find highest peak location and define sample range
|
|
%% find highest peak location and define sample range
|
|
@@ -38,72 +56,69 @@ function [BFF] = sphear_fir_render_directions(FIR, fmin, fmax, N, Npre, Fs, smoo
|
|
|
ref_fft = fft(tukeywin(size(range)(2), 0.02) .* ref(range), N);
|
|
ref_fft = fft(tukeywin(size(range)(2), 0.02) .* ref(range), N);
|
|
|
printf("reference: peak %f at %d, sample range [%d %d]\n", value, index, r_start, r_end);
|
|
printf("reference: peak %f at %d, sample range [%d %d]\n", value, index, r_start, r_end);
|
|
|
|
|
|
|
|
- capnames = cellstr(['LF'; 'RF'; 'LB'; 'RB']);
|
|
|
|
|
-
|
|
|
|
|
- minindex = N;
|
|
|
|
|
- maxpeak = 0;
|
|
|
|
|
-
|
|
|
|
|
- Nfft = 4096;
|
|
|
|
|
-
|
|
|
|
|
- f = figure;
|
|
|
|
|
- hold on;
|
|
|
|
|
- grid on;
|
|
|
|
|
- axis([fmin, fmax, -30, 0]);
|
|
|
|
|
-
|
|
|
|
|
- %% title(sprintf("azimuth = %05.1f degrees", angle));
|
|
|
|
|
- tics('y', [6, 0, -3, -6, -9, -12, -15, -20, -25, -30]);
|
|
|
|
|
- tics('x', [300, 500, 1000, 2000, 3000, 5000, 7000, 10000, 15000, 20000]);
|
|
|
|
|
-
|
|
|
|
|
- idx = 1;
|
|
|
|
|
-
|
|
|
|
|
- for m = [1, 5, 9, 13, 3, 7, 11, 15]
|
|
|
|
|
- %% render b format components
|
|
|
|
|
- bf(idx,:,:) = sphear_a_to_b(m, FIR, range, N);
|
|
|
|
|
-
|
|
|
|
|
- %% calculate B format frequency responses
|
|
|
|
|
- for b = 1:4
|
|
|
|
|
- bfft(idx,b,:) = fastsmooth(20*log10(abs(fft(bf(idx,b,:), Nfft))), smooth, 2, 1);
|
|
|
|
|
- end
|
|
|
|
|
- fr = (1:Nfft/2) * Fs/Nfft;
|
|
|
|
|
- semilogx(fr, bfft(idx,1,1:end/2), 'k');
|
|
|
|
|
- if (idx > 4)
|
|
|
|
|
- col = 'b';
|
|
|
|
|
- else
|
|
|
|
|
- col = 'r';
|
|
|
|
|
- end
|
|
|
|
|
- if (rem(idx, 2) == 1)
|
|
|
|
|
- semilogx(fr, bfft(idx,2,1:end/2), col);
|
|
|
|
|
|
|
+ %% plot if requested
|
|
|
|
|
+ if (plot_on != 0)
|
|
|
|
|
+ if (strcmp(plot_on, "epsc") || strcmp(plot_on, "pngc"))
|
|
|
|
|
+ %% color plot
|
|
|
|
|
+ wopts = 'k;W;';
|
|
|
|
|
+ popts = 'r;X;';
|
|
|
|
|
+ dopts = 'b;Y;';
|
|
|
|
|
+ copt = '-color';
|
|
|
|
|
+ cext = "_color";
|
|
|
else
|
|
else
|
|
|
- semilogx(fr, bfft(idx,3,1:end/2), col);
|
|
|
|
|
|
|
+ %% mono plot
|
|
|
|
|
+ wopts = '-;W;';
|
|
|
|
|
+ popts = '--;X;';
|
|
|
|
|
+ dopts = '-.;Y;';
|
|
|
|
|
+ copt = '-mono';
|
|
|
|
|
+ cext = "_mono";
|
|
|
end
|
|
end
|
|
|
|
|
+ f = figure;
|
|
|
|
|
+ hold on;
|
|
|
|
|
+ grid on;
|
|
|
|
|
+ axis([fmin, fmax, -30, 0]);
|
|
|
|
|
+ tics('y', [6, 0, -3, -6, -9, -12, -15, -20, -25, -30]);
|
|
|
|
|
+ tics('x', [300, 500, 1000, 2000, 3000, 5000, 7000, 10000, 15000, 20000]);
|
|
|
|
|
|
|
|
- idx = idx + 1;
|
|
|
|
|
- end
|
|
|
|
|
|
|
+ idx = 1;
|
|
|
|
|
+ for m = [1, 5, 9, 13, 3, 7, 11, 15]
|
|
|
|
|
+ %% render b format components
|
|
|
|
|
+ bf(idx,:,:) = sphear_a_to_b(m, FIR, range, N);
|
|
|
|
|
|
|
|
- print(f, '-deps','-color',["figures/sphear_directions.eps"]);
|
|
|
|
|
-end
|
|
|
|
|
|
|
+ %% calculate B format frequency responses
|
|
|
|
|
+ for b = 1:4
|
|
|
|
|
+ if (smooth > 0)
|
|
|
|
|
+ bfft(idx,b,:) = fastsmooth(20*log10(abs(fft(bf(idx,b,:), Nfft))), smooth, 2, 1);
|
|
|
|
|
+ else
|
|
|
|
|
+ bfft(idx,b,:) = 20*log10(abs(fft(bf(idx,b,:), Nfft)));
|
|
|
|
|
+ end
|
|
|
|
|
+ end
|
|
|
|
|
+ fr = (1:Nfft/2) * Fs/Nfft;
|
|
|
|
|
+ %% W
|
|
|
|
|
+ semilogx(fr, bfft(idx,1,1:end/2) * 10^(6/20), wopts);
|
|
|
|
|
+ %% print principal and diagonal directions with different colors
|
|
|
|
|
+ if (idx > 4)
|
|
|
|
|
+ col = popts;
|
|
|
|
|
+ else
|
|
|
|
|
+ col = dopts;
|
|
|
|
|
+ end
|
|
|
|
|
+ if (rem(idx, 2) == 1)
|
|
|
|
|
+ semilogx(fr, bfft(idx,2,1:end/2), col);
|
|
|
|
|
+ else
|
|
|
|
|
+ semilogx(fr, bfft(idx,3,1:end/2), col);
|
|
|
|
|
+ end
|
|
|
|
|
|
|
|
|
|
+ idx = idx + 1;
|
|
|
|
|
+ end
|
|
|
|
|
|
|
|
-function [bf] = sphear_a_to_b(m, FIR, range, N);
|
|
|
|
|
- capnames = cellstr(['LF'; 'RF'; 'LB'; 'RB']);
|
|
|
|
|
- %% read in capsule signals and apply FIR filter
|
|
|
|
|
- for c = 1:4
|
|
|
|
|
- cs(c,:) = wavread(["../cal/2016.04.20/export/a-format-0" sprintf("%02.0f", m) "-" char(capnames(c)) ".wav"]);
|
|
|
|
|
- cr(c,:) = tukeywin(size(range)(2), 0.0) .* cs(c,:)'(range);
|
|
|
|
|
- for b = 1:4
|
|
|
|
|
- bs(c,b,:) = fftfilt(reshape(FIR(b,c,:), 1, [])', cr(c,:), N);
|
|
|
|
|
|
|
+ l = legend('W (+6dB)', 'XY principal', 'XY diagonal');
|
|
|
|
|
+ legend(l, 'location', 'west');
|
|
|
|
|
+ if (strcmp(plot_on, 'eps') || strcmp(plot_on, 'epsc') || plot_on == 1)
|
|
|
|
|
+ print(f, '-deps', '-tight', copt,["figures/sphear_directions.eps"]);
|
|
|
|
|
+ else
|
|
|
|
|
+ print(f, '-deps', copt, ["figures/sphear_directions.eps"]);
|
|
|
end
|
|
end
|
|
|
- end
|
|
|
|
|
- %% add filtered capsule output to b format components
|
|
|
|
|
- for b = 1:4
|
|
|
|
|
- bf(b,:) = bs(1,b,:) + bs(2,b,:) + bs(3,b,:) + bs(4,b,:);
|
|
|
|
|
|
|
+
|
|
|
end
|
|
end
|
|
|
end
|
|
end
|
|
|
|
|
|
|
|
-
|
|
|
|
|
-function [i, v] = findmaxpeak(signal);
|
|
|
|
|
- [peaks, locs] = findpeaks(signal, "DoubleSided");
|
|
|
|
|
- [~, index] = max(peaks);
|
|
|
|
|
- i = locs(index);
|
|
|
|
|
- v = peaks(index);
|
|
|
|
|
-end
|
|
|