Explorar o código

add more generic sphear_signal_plot_frequency_response, modify design example to use it

Fernando Pablo Lopez-Lezcano %!s(int64=9) %!d(string=hai) anos
pai
achega
3945776c56

+ 12 - 8
calibration/src/sphear_design_example.m

@@ -70,7 +70,7 @@ sphear_af_plot_frequency_response(AF, 300, 20000, -60, -10, 48000, 8192);
 %% calculate signal power in reference frequency interval for all
 %% measurements and capsules
 
-[R] = sphear_signal_power(AF, FS, 600, 1000, 8192);
+[R] = sphear_signal_power(AF, FS, 600, 1200, 8192);
 
 [R] = sphear_signal_power(AF, FS, 800, 2400, 8192);
 
@@ -88,7 +88,7 @@ sphear_plot_r(R, 0);
 				  
 %% plot frequency response of B-format components for each measurement direction (just a test)
 
-sphear_bf_plot_frequency_response(BF, 600, 1200, 300, 20000, -60, 0, 48000, 8192);
+sphear_bf_plot_frequency_response(BF, 600, 1200, 300, 20000, -50, 9, 48000, 8192);
 
 %% calculate power of B-format signals in logarithmically spaced frequency bands
 
@@ -156,7 +156,7 @@ sphear_bf_plot_filters(FC, W, X, Y, Z, U, V);
 %% plot the frequency response of the B-format components for each measurement
 %% direction, they should all be flat
 
-sphear_bf_plot_frequency_response(BFF, 600, 1200, 300, 20000, -35, -3, 48000, 8192);
+sphear_bf_plot_frequency_response(BFF, 600, 1200, 300, 20000, -30, 6, 48000, 8192);
 
 %% a revealing plot that highlights the cardinal direction response versus
 %% the diagonal direction response (incomplete code for the 8 capsule case)
@@ -182,13 +182,17 @@ sphear_write_simple_a2b_encoder(A2B, TFIR);
 
 %% plot B format components separately for a set of measurements
 
-sphear_capsule_plot_frequency_response(BFF, [1:16], 300, 20000, -50, 0, 48000, 8192);
+sphear_signal_plot_frequency_response(BFF, [1:16], 300, 20000, -50, 0, 48000, 8192);
 
-sphear_capsule_plot_frequency_response(BFF, [1, 5, 9, 13], 300, 20000, -50, 0, 48000, 8192);
+%% tetra cardinal directions
+sphear_signal_plot_frequency_response(BFF, [1, 5, 9, 13], 300, 20000, -50, 0, 48000, 8192);
+%% tetra diagonal directions
+sphear_signal_plot_frequency_response(BFF, [1, 5, 9, 13] + 2, 300, 20000, -50, 0, 48000, 8192);
 
-sphear_capsule_plot_frequency_response(BFF, [1, 3, 5, 7, 9, 11, 13, 15], 300, 20000, -50, 0, 48000, 8192);
-
-sphear_capsule_plot_frequency_response(BFF, [1, 3, 5, 7, 9, 11, 13, 15] + 1, 300, 20000, -50, 0, 48000, 8192);
+%% octa cardinal directions
+sphear_signal_plot_frequency_response(BFF, [1, 3, 5, 7, 9, 11, 13, 15], 300, 20000, -50, 0, 48000, 8192);
+% octa diagonal directions
+sphear_signal_plot_frequency_response(BFF, [1, 3, 5, 7, 9, 11, 13, 15] + 1, 300, 20000, -50, 0, 48000, 8192);
 
 %%%%%%%%% testing code to determine z axis wobble during calibration measurements
 

+ 91 - 0
calibration/src/sphear_signal_plot_frequency_response.m

@@ -0,0 +1,91 @@
+%% -*- mode: octave -*-
+%%
+%% *SpHEAR Project
+%%
+%% sphear_signal_plot_frequency_response: plot frequency response of a signal
+%%
+%% usage:
+%%   [] = sphear_signal_plot_frequency_response(SIG, fmin, fmax, ymin, ymax, FS, fft_size, title, filename, smooth)
+%%
+%% SIG: signals for each measurement
+%% fmin, fmax: range of frequencies to use
+%% ymin, ymax: y range to plot
+%% FS: sampling rate
+%% fft_size: size of fft used to calculate frequency response
+%%
+%% optional parameters:
+%%
+%% title: title to add to figure
+%% filename: prefix to filename of figure (figure is not written unless this is specified)
+%% smooth: number of samples to average, default is 0
+%%
+%% Copyright 2016-2017, Fernando Lopez-Lezcano, All Rights Reserved
+%% nando@ccrma.stanford.edu
+%%
+%%  This program is free software: you can redistribute it and/or modify
+%%  it under the terms of the GNU General Public License as published by
+%%  the Free Software Foundation, either version 3 of the License, or
+%%  (at your option) any later version.
+%%
+%%  This program is distributed in the hope that it will be useful,
+%%  but WITHOUT ANY WARRANTY; without even the implied warranty of
+%%  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+%%  GNU General Public License for more details.
+%%
+%%  You should have received a copy of the GNU General Public License
+%%  along with this program.  If not, see <htp://www.gnu.org/licenses/>.
+
+function [] = sphear_signal_plot_frequency_response(SIG, measurements, fmin, fmax, ymin, ymax, FS, fft_size, figure_title, filename, smooth)
+
+  %% no smoothing by default
+  if (!exist("smooth"))
+    smooth = 0;
+  end
+  for s = 1:size(SIG,2)
+    %% calculate frequency responses
+    for m = measurements
+      if (smooth > 0)
+	sfft(m,s,:) = fastsmooth(20*log10(abs(fft(SIG(m,s,:), fft_size))), smooth, 2, 1);
+      else
+	sfft(m,s,:) = 20*log10(abs(fft(SIG(m,s,:), fft_size)));
+      end
+    end
+    fig = figure;
+    hold on;
+    fr = (1:fft_size/2) * FS/fft_size;
+    grid on;
+    if (exist("figure_title"))
+      title(figure_title);
+    end
+    axis([fmin, fmax, ymin, ymax]);
+    tics('y', [9, 6, 3, 0, -3, -6, -9, -12, -15, -20, -25, -30, -40, -50]);
+    tics('x', [300, 500, 1000, 2000, 3000, 5000, 7000, 10000, 15000, 20000]);
+    ca=get (gcf, "currentaxes");
+    set(ca, "fontweight", "bold", "linewidth", 1.5)
+    pbaspect([16, 9, 1])
+    for m = measurements
+      if (m == 1)
+	col = "k";
+      else
+	if (m < 8)
+	  col = "r";
+	else
+	  if (m == 9)
+	    col = "g";
+	  else
+	    col = "b";
+	  end
+	end
+      end
+      semilogx(fr, sfft(m,s,1:end/2), col, "linewidth", 1.5);
+    end
+    hold off;
+    if (exist("filename"))
+      %% write file with eps plot
+      print(fig, '-deps', '-tight', '-color', '-depsc', '-S1920,1080', ["../data/current/figures/" filename ".eps"]);
+      %% write png file
+      print(fig, '-deps', '-tight', '-color', '-dpng', '-S1920,1080', ["../data/current/figures/" filename ".png"]);
+    end
+  end
+end
+