function filter = hohmann2002_filter(arg1,arg2,arg3,arg4,arg5)
%hohmann2002_filter Construct a Gammatone filter within the HOHMANN2002 framework
% Usage: hohmann2002_filter(a_tilde,gamma_order);
% hohmann2002_filter(fs, fc, gamma_order);
% hohmann2002_filter(fs, fc, gamma_order, bandwidth_factor);
% hohmann2002_filter(fs, fc, bw, attenuation_db, gamma_order);
%
% Input parameters:
% a_tilde : Complex-valued filter coefficients
% fs : Sampling rate (in Hz)
% fc : Centre frequency (in Hz)
% bw : Filter bandwidth (in Hz)
% gamma_order : Gammatone filter order
% bandwidth_factor : Bandwidth factor with respect to 1 ERB. Default: 1.0
%
% HOHMANN2002_FILTER(a_tilde, gamma_order) specifies the complex-valued
% filter coefficients directly as at the beginning of Sec. 2.3 in Hohmann
% (2002).
%
% HOHMANN2002_FILTER(fs, fc, gamma_order) computes
% filter coefficients for sampling rate fs, center frequency fc, and order of
% the gammatone filter gamma_order. The filter bandwidth will be 1 ERB.
% Filter coefficients are computed according to Eq. 13 and 14 from Hohmann (2002).
%
% HOHMANN2002_FILTER(fs, fc, gamma_order, bandwidth_factor) uses filter
% bandwidth of bandwidth_factor ERB (instead of 1 ERB).
%
% HOHMANN2002_FILTER(fs, fc, bw, attenuation_db, gamma_order) designs
% the filter such that for the bandwidth bw, the attenuation attenuation_db
% at the corner frequencies of the filter is reached. Filter coefficient
% are computed as in Eq. 11 and 12 from Hohmann (2002)
%
% References:
% V. Hohmann. Frequency analysis and synthesis using a gammatone
% filterbank. Acta Acustica united with Acoustica, 88(3):433--442, 2002.
%
%
% Url: http://amtoolbox.sourceforge.net/amt-0.10.0/doc/modelstages/hohmann2002_filter.php
% Copyright (C) 2009-2020 Piotr Majdak and the AMT team.
% This file is part of Auditory Modeling Toolbox (AMT) version 1.0.0
%
% 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 <http://www.gnu.org/licenses/>.
% AUTHOR: Universitaet Oldenburg, tp
% Adapted to AMT (PM, Jan 2016) from function gfb_filter_new
filter.type = 'gfb_Filter';
switch nargin
case 2 % a_tilde, gamma_order
filter.coefficient = arg1;
filter.gamma_order = arg2;
case 3 % fs, fc, gamma_order
fs = arg1;
fc = arg2;
filter.gamma_order = arg3;
filter.bandwidth_factor = 1.0;
filter.L=24.7; % Eq. 17 in [Hohmann 2002]
filter.Q=9.265; % Eq. 17 in [Hohmann 2002]
audiological_erb = (filter.L + fc / filter.Q) * filter.bandwidth_factor; % Eq. 13 in [Hohmann 2002]
go = filter.gamma_order;
a_gamma = (pi * factorial(2*go - 2) * 2^-(2*go - 2) / factorial(go - 1)^2); % Eq. 14, line 3 [Hohmann 2002]
b = audiological_erb / a_gamma; % Eq. 14, line 2 [Hohmann 2002]
lambda = exp(-2 * pi * b / fs); % Eq. 14, line 1 [Hohmann 2002]
beta = 2 * pi * fc / fs; % Eq. 10 [Hohmann 2002]
filter.coefficient = lambda * exp(1i * beta); % Eq 1, line 2 [Hohmann 2002]
case 4 % fs, fc, gamma_order, bandwidth_factor
fs = arg1;
fc = arg2;
filter.gamma_order = arg3;
filter.bandwidth_factor = arg4;
filter.L=24.7; % see equation (17) in [Hohmann 2002]
filter.Q=9.265; % see equation (17) in [Hohmann 2002]
audiological_erb = (filter.L + fc / filter.Q) * filter.bandwidth_factor; % Eq. 13 in [Hohmann 2002]
go = filter.gamma_order;
a_gamma = (pi * factorial(2*go-2) * 2^-(2*go-2) / factorial(go-1)^2); % Eq. 14, line 3 [Hohmann 2002]
b = audiological_erb / a_gamma; % Eq. 14, line 2 [Hohmann 2002]
lambda = exp(-2 * pi * b / fs); % Eq. 14, line 1 [Hohmann 2002]
beta = 2 * pi * fc / fs; % Eq. 10 [Hohmann 2002]
filter.coefficient = lambda * exp(1i * beta); % Eq. 1, line 2 [Hohmann 2002]
case 5 % fs, fc, bw, attenuation_db, gamma_order
fs = arg1;
fc = arg2;
bw = arg3;
attenuation_db = arg4;
filter.gamma_order = arg5;
phi = pi * bw / fs; % Eq. 12, line 4 [Hohmann 2002]
u = -attenuation_db/filter.gamma_order; % Eq. 12, line 3 [Hohmann 2002]
p = (-2 + 2 * 10^(u/10) * cos(phi)) / (1 - 10^(u/10)); % Eq. 12, line 2 [Hohmann 2002]
lambda = -p/2 - sqrt(p*p/4 - 1); % Eq. 12, line 1 [Hohmann 2002]
beta = 2*pi*fc/fs; % Eq. 10, [Hohmann 2002]
filter.coefficient = lambda * exp(1i*beta); % Eq. 1, line 2 [Hohmann 2002]
otherwise
error ('Number of required arguments: 2, 3, 4 or 5');
end
% normalization factor from section 2.2 (text) [Hohmann 2002]:
filter.normalization_factor = 2 * (1-abs(filter.coefficient)) ^ filter.gamma_order;
filter.state = zeros(1, filter.gamma_order);