dicklyon@534
|
1 % Copyright 2012, Google, Inc.
|
dicklyon@534
|
2 % Author: Richard F. Lyon
|
dicklyon@534
|
3 %
|
dicklyon@534
|
4 % This Matlab file is part of an implementation of Lyon's cochlear model:
|
dicklyon@534
|
5 % "Cascade of Asymmetric Resonators with Fast-Acting Compression"
|
dicklyon@534
|
6 % to supplement Lyon's upcoming book "Human and Machine Hearing"
|
dicklyon@534
|
7 %
|
dicklyon@534
|
8 % Licensed under the Apache License, Version 2.0 (the "License");
|
dicklyon@534
|
9 % you may not use this file except in compliance with the License.
|
dicklyon@534
|
10 % You may obtain a copy of the License at
|
dicklyon@534
|
11 %
|
dicklyon@534
|
12 % http://www.apache.org/licenses/LICENSE-2.0
|
dicklyon@534
|
13 %
|
dicklyon@534
|
14 % Unless required by applicable law or agreed to in writing, software
|
dicklyon@534
|
15 % distributed under the License is distributed on an "AS IS" BASIS,
|
dicklyon@534
|
16 % WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
|
dicklyon@534
|
17 % See the License for the specific language governing permissions and
|
dicklyon@534
|
18 % limitations under the License.
|
dicklyon@534
|
19
|
dicklyon@534
|
20 function [zY, state] = CARFAC_CAR_Step(x_in, CAR_coeffs, state)
|
dicklyon@534
|
21 % function [zY, state] = CARFAC_CAR_Step(x_in, CAR_coeffs, state)
|
dicklyon@534
|
22 %
|
dicklyon@534
|
23 % One sample-time update step for the filter part of the CARFAC.
|
dicklyon@534
|
24
|
dicklyon@534
|
25 % Most of the update is parallel; finally we ripple inputs at the end.
|
dicklyon@534
|
26
|
dicklyon@534
|
27 % Local nonlinearity zA and AGC feedback zB reduce pole radius:
|
dicklyon@534
|
28 zA = state.zA_memory;
|
dicklyon@534
|
29 zB = state.zB_memory + state.dzB_memory; % AGC interpolation
|
dicklyon@534
|
30 r1 = CAR_coeffs.r1_coeffs;
|
dicklyon@534
|
31 g = state.g_memory + state.dg_memory; % interp g
|
dicklyon@534
|
32 v_offset = CAR_coeffs.v_offset;
|
dicklyon@534
|
33 v2_corner = CAR_coeffs.v2_corner;
|
dicklyon@534
|
34 v_damp_max = CAR_coeffs.v_damp_max;
|
dicklyon@534
|
35
|
dicklyon@534
|
36 % zB and zA are "extra damping", and multiply zr (compressed theta):
|
dicklyon@534
|
37 r = r1 - CAR_coeffs.zr_coeffs .* (zA + zB);
|
dicklyon@534
|
38
|
dicklyon@534
|
39 % now reduce state by r and rotate with the fixed cos/sin coeffs:
|
dicklyon@534
|
40 z1 = r .* (CAR_coeffs.a0_coeffs .* state.z1_memory - ...
|
dicklyon@534
|
41 CAR_coeffs.c0_coeffs .* state.z2_memory);
|
dicklyon@534
|
42 % z1 = z1 + inputs;
|
dicklyon@534
|
43 z2 = r .* (CAR_coeffs.c0_coeffs .* state.z1_memory + ...
|
dicklyon@534
|
44 CAR_coeffs.a0_coeffs .* state.z2_memory);
|
dicklyon@534
|
45
|
dicklyon@534
|
46 % update the "velocity" for cubic nonlinearity, into zA:
|
dicklyon@534
|
47 zA = (((state.z2_memory - z2) .* CAR_coeffs.velocity_scale) + ...
|
dicklyon@534
|
48 v_offset) .^ 2;
|
dicklyon@534
|
49 % soft saturation to make it more like an "essential" nonlinearity:
|
dicklyon@534
|
50 zA = v_damp_max * zA ./ (v2_corner + zA);
|
dicklyon@534
|
51
|
dicklyon@534
|
52 zY = CAR_coeffs.h_coeffs .* z2; % partial output
|
dicklyon@534
|
53
|
dicklyon@534
|
54 % Ripple input-output path, instead of parallel, to avoid delay...
|
dicklyon@534
|
55 % this is the only part that doesn't get computed "in parallel":
|
dicklyon@534
|
56 in_out = x_in;
|
dicklyon@534
|
57 for ch = 1:length(zY)
|
dicklyon@534
|
58 % could do this here, or later in parallel:
|
dicklyon@534
|
59 z1(ch) = z1(ch) + in_out;
|
dicklyon@534
|
60 % ripple, saving final channel outputs in zY
|
dicklyon@534
|
61 in_out = g(ch) * (in_out + zY(ch));
|
dicklyon@534
|
62 zY(ch) = in_out;
|
dicklyon@534
|
63 end
|
dicklyon@534
|
64
|
dicklyon@534
|
65 % put new state back in place of old
|
dicklyon@534
|
66 % (z1 and z2 are genuine temps; the others can update by reference in C)
|
dicklyon@534
|
67 state.z1_memory = z1;
|
dicklyon@534
|
68 state.z2_memory = z2;
|
dicklyon@534
|
69 state.zA_memory = zA;
|
dicklyon@534
|
70 state.zB_memory = zB;
|
dicklyon@534
|
71 state.zY_memory = zY;
|
dicklyon@534
|
72 state.g_memory = g;
|
dicklyon@534
|
73
|