Description of frq2bark

0001 function [b,c] = frq2bark(f,m)
0002 %FRQ2BARK  Convert Hertz to BARK frequency scale BARK=(FRQ)
0003 %       bark = frq2bark(frq) converts a vector of frequencies (in Hz)
0004 %       to the corresponding values on the BARK scale.
0005 % Inputs: f  matrix of frequencies in Hz
0006 %         m  mode options
0007 %            'h'   use high frequency correction from [1]
0008 %            'l'   use low frequency correction from [1]
0009 %            'H'   do not apply any high frequency correction
0010 %            'L'   do not apply any low frequency correction
0011 %            'z'   use the expressions from Zwicker et al. (1980) for b and c
0012 %            's'   use the expression from Schroeder et al. (1979)
0013 %            'u'   unipolar version: do not force b to be an odd function
0014 %                  This has no effect on the default function which is odd anyway
0015 %            'g'   plot a graph
0016 %
0017 % Outputs: b  bark values
0018 %          c  Critical bandwidth: d(freq)/d(bark)
0019 
0020 %   The Bark scale is named in honour of Barkhausen, the creator
0021 %   of the unit of loudness level [2]. Criitical band k extends
0022 %   from bark2frq(k-1) to bark2frq(k).
0023 %
0024 %   There are many published formulae approximating the Bark scale.
0025 %   The default is the one from [1] but with a correction at high and
0026 %   low frequencies to give a better fit to [2] with a continuous derivative
0027 %   and ensure that 0 Hz = 0 Bark.
0028 %   The h and l mode options apply the corrections from [1] which are
0029 %   not as good and do not give a continuous derivative. The H and L
0030 %   mode options suppress the correction entirely to give a simple formula.
0031 %   The 's' option uses the less accurate formulae from [3] which have been
0032 %   widely used in the lterature.
0033 %   The 'z' option uses the formulae from [4] in which the c output
0034 %   is not exactly the reciprocal of the derivative of the bark function.
0035 %
0036 %   [1] H. Traunmuller, Analytical Expressions for the
0037 %       Tonotopic Sensory Scale”, J. Acoust. Soc. Am. 88,
0038 %       1990, pp. 97-100.
0039 %   [2] E. Zwicker, Subdivision of the audible frequency range into
0040 %       critical bands, J Accoust Soc Am 33, 1961, p248.
0041 %   [3] M. R. Schroeder, B. S. Atal, and J. L. Hall. Optimizing digital
0042 %       speech coders by exploiting masking properties of the human ear.
0043 %       J. Acoust Soc Amer, 66 (6): 1647–1652, 1979. doi: 10.1121/1.383662.
0044 %   [4] E. Zwicker and E. Terhardt.  Analytical expressions for
0045 %       critical-band rate and critical bandwidth as a function of frequency.
0046 %       J. Acoust Soc Amer, 68 (5): 1523–1525, Nov. 1980.
0047 
0048 %   The following code reproduces the graphs 3(c) and 3(d) from [1].
0049 %       b0=(0:0.5:24)';
0050 %       f0=[[2 5 10 15 20 25 30 35 40 45 51 57 63 70 77 ...
0051 %           84 92 100 108 117 127 137 148 160 172 185 200 ...
0052 %           215 232 250 270 290 315]*10 [34 37 40 44 48 53 ...
0053 %           58 64 70 77 85 95 105 120 135 155]*100]';
0054 %       b1=frq2bark(f0);      b2=frq2bark(f0,'lh');
0055 %       b3=frq2bark(f0,'LH'); b4=frq2bark(f0,'z');
0056 %       plot(b0,[b0 b1 b2 b3 b4]-repmat(b0,1,5));
0057 %       xlabel('Frequency (Bark)'); ylabel('Error (Bark)');
0058 %       legend('Exact','voicebox','Traunmuller1990', ...
0059 %              'Traunmuller1983','Zwicker1980','Location','South');
0060 
0061 %      Copyright (C) Mike Brookes 2006-2010
0062 %      Version: $Id: frq2bark.m 713 2011-10-16 14:45:43Z dmb $
0063 %
0064 %   VOICEBOX is a MATLAB toolbox for speech processing.
0065 %   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
0066 %
0067 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0068 %   This program is free software; you can redistribute it and/or modify
0069 %   it under the terms of the GNU General Public License as published by
0070 %   the Free Software Foundation; either version 2 of the License, or
0071 %   (at your option) any later version.
0072 %
0073 %   This program is distributed in the hope that it will be useful,
0074 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0075 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0076 %   GNU General Public License for more details.
0077 %
0078 %   You can obtain a copy of the GNU General Public License from
0079 %   http://www.gnu.org/copyleft/gpl.html or by writing to
0080 %   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
0081 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0082 persistent A B C D P Q R S T U
0083 if isempty(P)
0084     A=26.81;
0085     B=1960;
0086     C=-0.53;
0087     D=A*B;
0088     P=(0.53/(3.53)^2);
0089     Q=0.25;
0090     R=20.4;
0091     xy=2;
0092     S=0.5*Q/xy;
0093     T=R+0.5*xy;
0094     U=T-xy;
0095 end
0096 if nargin<2
0097     m=' ';
0098 end
0099 if any(m=='u')
0100     g=f;
0101 else
0102     g=abs(f);
0103 end
0104 if any(m=='z')
0105     b=13*atan(0.00076*g)+3.5*atan((f/7500).^2);
0106     c=25+75*(1+1.4e-6*f.^2).^0.69;
0107 elseif any(m=='s')
0108     b=7*log(g/650+sqrt(1+(g/650).^2));
0109     c=cosh(b/7)*650/7;
0110 else
0111     b=A*g./(B+g)+C;
0112     d=D*(B+g).^(-2);
0113     if any(m=='l')
0114         m1=(b<2);
0115         d(m1)=d(m1)*0.85;
0116         b(m1)=0.3+0.85*b(m1);
0117     elseif ~any(m=='L')
0118         m1=(b<3);
0119         b(m1)=b(m1)+P*(3-b(m1)).^2;
0120         d(m1)=d(m1).*(1-2*P*(3-b(m1)));
0121     end
0122     if any(m=='h')
0123         m1=(b>20.1);
0124         d(m1)=d(m1)*1.22;
0125         b(m1)=1.22*b(m1)-4.422;
0126     elseif ~any(m=='H')
0127         m2=(b>T);
0128         m1=(b>U) & ~m2;
0129         b(m1)=b(m1)+S*(b(m1)-U).^2;
0130         b(m2)=(1+Q)*b(m2)-Q*R;
0131         d(m2)=d(m2).*(1+Q);
0132         d(m1)=d(m1).*(1+2*S*(b(m1)-U));
0133     end
0134     c=d.^(-1);
0135 end
0136 if ~any(m=='u')
0137     b=b.*sign(f);          % force to be odd
0138 end
0139 
0140 if ~nargout || any(m=='g')
0141     subplot(212)
0142     semilogy(f,c,'-r');
0143     ha=gca;
0144     ylabel(['Critical BW (' yticksi 'Hz)']);
0145     xlabel(['Frequency (' xticksi 'Hz)']);
0146     subplot(211)
0147     plot(f,b,'x-b');
0148     hb=gca;
0149     ylabel('Bark');
0150     xlabel(['Frequency (' xticksi 'Hz)']);
0151     linkaxes([ha hb],'x');
0152 end
frq2bark

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SOURCE CODE
