changeset 528:741187dc780f

More updates

--- a/trunk/matlab/bmm/carfac/CARFAC_Design.m	Wed Mar 07 23:18:32 2012 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_Design.m	Sat Mar 10 05:05:35 2012 +0000
@@ -97,10 +97,10 @@
     'v2_corner', 0.2, ...  % corner for essential nonlin
     'v_damp_max', 0.01, ... % damping delta damping from velocity nonlin
     'min_zeta', 0.12, ...
-    'first_pole_theta', 0.78*pi, ...
-    'zero_ratio', sqrt(2), ...
-    'ERB_per_step', 0.3333, ... % assume G&M's ERB formula
-    'min_pole_Hz', 40 );
+    'first_pole_theta', 0.85*pi, ...
+    'zero_ratio', sqrt(2), ... % how far zero is above pole
+    'ERB_per_step', 0.5, ... % assume G&M's ERB formula
+    'min_pole_Hz', 30 );
 end
 
 if nargin < 1
@@ -126,8 +126,11 @@
 end
 % now we have n_ch, the number of channels, and pole_freqs array
 
+max_channels_per_octave = log(2) / log(pole_freqs(1)/pole_freqs(2));
+
 CF = struct( ...
   'fs', fs, ...
+  'max_channels_per_octave', max_channels_per_octave, ...
   'filter_params', CF_filter_params, ...
   'AGC_params', CF_AGC_params, ...
   'IHC_params', CF_IHC_params, ...
@@ -175,7 +178,8 @@
 % different possible interpretations for min-damping r:
 % r = exp(-theta * CF_filter_params.min_zeta).
 % Using sin gives somewhat higher Q at highest thetas.
-r = (1 - sin(theta) * filter_params.min_zeta);
+ff = 5;  % fudge factor for theta distortion; at least 1.0
+r = (1 - ff*sin(theta/ff) * filter_params.min_zeta);
 filter_coeffs.r_coeffs = r;
 
 % undamped coupled-form coefficients:
@@ -186,10 +190,16 @@
 h = sin(theta) .* f;
 filter_coeffs.h_coeffs = h;
 
-% unity gain at min damping, radius r:
-filter_coeffs.g_coeffs = (1 - 2*r.*cos(theta) + r.^2) ./ ...
+% % unity gain at min damping, radius r:
+g = (1 - 2*r.*cos(theta) + r.^2) ./ ...
   (1 - 2*r .* cos(theta) + h .* r .* sin(theta) + r.^2);
+% or assume r is 1, for the zero-damping gain g0:
+g0 = (2 - 2*cos(theta)) ./ ...
+  (2 - 2 * cos(theta) + h .* sin(theta));
 
+filter_coeffs.g_coeffs = g0;
+% make coeffs that can correct g0 to make g based on (1 - r).^2:
+filter_coeffs.gr_coeffs = ((g ./ g0) - 1) ./ ((1 - r).^2);
 
 %% the AGC design coeffs:
 function AGC_coeffs = CARFAC_DesignAGC(AGC_params, fs)

--- a/trunk/matlab/bmm/carfac/CARFAC_FilterStep.m	Wed Mar 07 23:18:32 2012 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_FilterStep.m	Sat Mar 10 05:05:35 2012 +0000
@@ -33,7 +33,8 @@
 v2_corner = filter_coeffs.v2_corner;
 v_damp_max = filter_coeffs.v_damp_max;
 
-r = r0 - filter_coeffs.c_coeffs .* (zA + zB);
+% zB and zA are "extra damping", and multiply c or sin(theta):
+r = r0 - filter_coeffs.c_coeffs .* (zA + zB); 
 
 % now reduce state by r and rotate with the fixed cos/sin coeffs:
 z1 = r .* (filter_coeffs.a_coeffs .* state.z1_memory - ...
@@ -47,8 +48,12 @@
   v_offset) .^ 2;
 zA = v_damp_max * zA ./ (v2_corner + zA);  % make it more like an "essential" nonlinearity
 
+% Adjust gain for r variation:
+g = filter_coeffs.g_coeffs;
+g = g .* (1 + filter_coeffs.gr_coeffs .* (1 - r).^2);
+
 % Get outputs from inputs and new z2 values:
-zY = filter_coeffs.g_coeffs .* (inputs + filter_coeffs.h_coeffs .* z2);
+zY = g .* (inputs + filter_coeffs.h_coeffs .* z2);
 
 % put new state back in place of old
 state.z1_memory = z1;
@@ -57,7 +62,7 @@
 state.zB_memory = zB;
 state.zY_memory = zY;
 
-% accum the straight hwr version for the sake of AGC range:
-hwr_detect = max(0, zY);  % detect with HWR
-state.detect_accum = state.detect_accum + hwr_detect;
+% % accum the straight hwr version for the sake of AGC range:
+% hwr_detect = max(0, zY);  % detect with HWR
+% state.detect_accum = state.detect_accum + hwr_detect;
 

--- a/trunk/matlab/bmm/carfac/CARFAC_Run.m	Wed Mar 07 23:18:32 2012 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_Run.m	Sat Mar 10 05:05:35 2012 +0000
@@ -130,13 +130,12 @@
       plot(CF.AGC_state(mic).AGC_memory(:, 1), 'k-', 'LineWidth', 1)
       maxes(mic) = max(CF.AGC_state(mic).AGC_memory(:));
       hold on
-      stage1 = 4; % as opposed to stage
       for stage = 1:3;
-        plot(2^(stage1-1) * (CF.AGC_state(mic).AGC_memory(:, stage) - ...
+        plot(2^(stage-1) * (CF.AGC_state(mic).AGC_memory(:, stage) - ...
           2 * CF.AGC_state(mic).AGC_memory(:, stage+1)));
       end
       stage = 4;
-      plot(2^(stage1-1) * CF.AGC_state(mic).AGC_memory(:, stage));
+      plot(2^(stage-1) * CF.AGC_state(mic).AGC_memory(:, stage));
     end
     axis([0, CF.n_ch+1, -0.01, max(maxes) + 0.01]);
     drawnow

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_Run_Linear.m	Sat Mar 10 05:05:35 2012 +0000
@@ -0,0 +1,57 @@
+% Copyright 2012, Google, Inc.
+% Author Richard F. Lyon
+%
+% This Matlab file is part of an implementation of Lyon's cochlear model:
+% "Cascade of Asymmetric Resonators with Fast-Acting Compression"
+% to supplement Lyon's upcoming book "Human and Machine Hearing"
+%
+% Licensed under the Apache License, Version 2.0 (the "License");
+% you may not use this file except in compliance with the License.
+% You may obtain a copy of the License at
+%
+%     http://www.apache.org/licenses/LICENSE-2.0
+%
+% Unless required by applicable law or agreed to in writing, software
+% distributed under the License is distributed on an "AS IS" BASIS,
+% WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+% See the License for the specific language governing permissions and
+% limitations under the License.
+
+function [naps, CF] = CARFAC_Run_Linear(CF, input_waves, extra_damping)
+% function [naps, CF] = CARFAC_Run_Linear(CF, input_waves, extra_damping)
+%
+% This function runs the CARFAC; that is, filters a 1 or more channel
+% sound input to make one or more neural activity patterns (naps);
+% however, unlike CARFAC_Run, it forces it to be linear, and gives a
+% linear (not detected) output.
+
+saved_v_damp_max = CF.filter_coeffs.v_damp_max;
+CF.filter_coeffs.v_damp_max = 0.00;  % make it linear for now
+
+[n_samp, n_mics] = size(input_waves);
+n_ch = CF.n_ch;
+
+if n_mics ~= CF.n_mics
+  error('bad number of input_waves channels passed to CARFAC_Run')
+end
+
+for mic = 1:CF.n_mics
+  % for the state of the AGC interpolator:
+  CF.filter_state(mic).zB_memory(:) = extra_damping;  % interpolator state
+  CF.filter_state(mic).dzB_memory(:) = 0;  % interpolator slope
+end
+
+naps = zeros(n_samp, n_ch, n_mics);
+
+for k = 1:n_samp
+  % at each time step, possibly handle multiple channels
+  for mic = 1:n_mics
+    [filters_out, CF.filter_state(mic)] = CARFAC_FilterStep( ...
+      input_waves(k, mic), CF.filter_coeffs, CF.filter_state(mic));
+    naps(k, :, mic) = filters_out;  % linear
+  end
+  % skip IHC and AGC updates
+end
+
+CF.filter_coeffs.v_damp_max = saved_v_damp_max;
+

--- a/trunk/matlab/bmm/carfac/CARFAC_Transfer_Functions.m	Wed Mar 07 23:18:32 2012 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_Transfer_Functions.m	Sat Mar 10 05:05:35 2012 +0000
@@ -18,11 +18,11 @@
 % limitations under the License.
 
 function [complex_transfns_freqs, ...
-  stage_numerators, stage_denominators] = ...
-  CARFAC_Transfer_Functions(CF, freqs, to_channels, from_channels)
+  stage_numerators, stage_denominators] = CARFAC_Transfer_Functions( ...
+  CF, freqs, to_channels, from_channels)
 % function [complex_transfns_freqs, ...
-%   stage_z_numerators, stage_z_denominators] = ...
-%   CARFAC_Transfer_Functions(CF, freqs, to_channels, from_channels)
+%   stage_numerators, stage_denominators] = CARFAC_Transfer_Functions( ...
+%   CF, freqs, to_channels, from_channels)
 % Return transfer functions as polynomials in z (nums & denoms);
 % And evaluate them at freqs if it's given, to selected output,
 %   optionally from selected starting points (from 0, input, by default).
@@ -63,7 +63,7 @@
     from_channels = 0;  % tranfuns from input, called channel 0.
   end
   if length(from_channels) == 1
-    from_channels = from_channels * ones(length(to_channels));
+    from_channels = from_channels * ones(1,length(to_channels));
   end
   % Default to cum gain of 1 (log is 0), from input channel 0:
   from_cum = zeros(length(to_channels), length(z_row));
@@ -88,24 +88,35 @@
 
 
 function [stage_numerators, stage_denominators] = ...
-  CARFAC_Rational_Functions(CF, chans)
+  CARFAC_Rational_Functions(CF)
 % function [stage_z_numerators, stage_z_denominators] = ...
 %   CARFAC_Rational_Functions(CF, chans)
 % Return transfer functions of all stages as rational functions.
 
-if nargin < 2
-  n_ch = CF.n_ch;
-  chans = 1:n_ch;
+n_ch = CF.n_ch;
+coeffs = CF.filter_coeffs;
+min_zeta = CF.filter_params.min_zeta;
+
+a0 = coeffs.a_coeffs;
+c0 = coeffs.c_coeffs;
+
+% get r, adapted if we have state:
+r =  coeffs.r_coeffs;
+if isfield(CF, 'filter_state')
+  state = CF.filter_state;
+  zB = state.zB_memory; % current extra damping
+  r = r - c0 .* zB;
 else
-  n_ch = length(chans);
+  zB = 0;
 end
-coeffs = CF.filter_coeffs;
-r = coeffs.r_coeffs(chans);
-a = coeffs.a_coeffs(chans) .* r;
-c = coeffs.c_coeffs(chans) .* r;
+
+a = a0 .* r;
+c = c0 .* r;
 r2 = r .* r;
-h = coeffs.h_coeffs(chans);
-g = coeffs.g_coeffs(chans);
+h = coeffs.h_coeffs;
+g0 = coeffs.g_coeffs;
+g = g0 .* (1 + coeffs.gr_coeffs .* (1 - r).^2);
+
 stage_denominators = [ones(n_ch, 1), -2 * a, r2];
 stage_numerators = [g .* ones(n_ch, 1), g .* (-2 * a + h .* c), g .* r2];
 

--- a/trunk/matlab/bmm/carfac/CARFAC_hacking.m	Wed Mar 07 23:18:32 2012 +0000
+++ b/trunk/matlab/bmm/carfac/CARFAC_hacking.m	Sat Mar 10 05:05:35 2012 +0000
@@ -52,7 +52,7 @@
 %   nap = deskew(nap);  % deskew doesn't make much difference
 
   if n_mics == 1  % because this hack doesn't work for binarual yet
-    MultiScaleSmooth(nap_decim, 10);
+    MultiScaleSmooth(nap_decim, 4);
   end
 
 %   nap_decim = nap;