Mercurial > hg > camir-aes2014
comparison toolboxes/MIRtoolbox1.3.2/MIRToolbox/mirfilterbank.m @ 0:e9a9cd732c1e tip
first hg version after svn
author | wolffd |
---|---|
date | Tue, 10 Feb 2015 15:05:51 +0000 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:e9a9cd732c1e |
---|---|
1 function varargout = mirfilterbank(orig,varargin) | |
2 % b = mirfilterbank(x) performs a filterbank decomposition of an audio | |
3 % waveform. | |
4 % Optional arguments: | |
5 % mirfilterbank(...,t) selects a type of filterbank. | |
6 % Possible values: | |
7 % t = 'Gammatone' (default for audio files). | |
8 % Requires the Auditory Toolbox. | |
9 % mirfilterbank(...,'Lowest',f): lowest frequency in Hz | |
10 % (default: 50) | |
11 % t = '2Channels' proposed in (Tolonen & Karjalainen, 2000) | |
12 % mirfilterbank(...,'NbChannels',N), or simply filterbank(x,N): | |
13 % specifies the number of channels in the bank. | |
14 % (default: N = 10) | |
15 % mirfilterbank(...,'Channel',c) only output the channels whose | |
16 % ranks are indicated in the array c. | |
17 % (default: c = (1:N)) | |
18 % mirfilterbank(...,'Manual',f) specifies a set of low-pass, band- | |
19 % pass and high-pass eliptic filters (Scheirer, 1998). | |
20 % The series of cut-off frequencies as to be specified | |
21 % as next parameter f. | |
22 % If this series of frequencies f begins with -Inf, | |
23 % the first filter is low-pass. | |
24 % If this series of frequencies f ends with Inf, | |
25 % the last filter is high-pass. | |
26 % mirfilterbank(...,'Order',o) specifies the order of the filter. | |
27 % (Default: o = 4) (Scheirer, 1998) | |
28 % mirfilterbank(...,'Hop',h) specifies the degree of spectral | |
29 % overlapping between successive channels. | |
30 % If h = 1 (default value), the filters are non-overlapping. | |
31 % If h = 2, the filters are half-overlapping. | |
32 % If h = 3, the spectral hop factor between successive | |
33 % filters is a third of the whole frequency region, etc. | |
34 % mirfilterbank(...,p) specifies predefined filterbanks, all | |
35 % implemented using elliptic filters, by default of order 4. | |
36 % Possible values: | |
37 % p = 'Mel' (mel-scale) | |
38 % p = 'Bark' (bark-scale) | |
39 % p = 'Scheirer' proposed in (Scheirer, 1998) corresponds to | |
40 % 'Manual',[-Inf 200 400 800 1600 3200 Inf] | |
41 % p = 'Klapuri' proposed in (Klapuri, 1999) corresponds to | |
42 % 'Manual',44*[2.^ ([ 0:2, ( 9+(0:17) )/3 ]) ] | |
43 | |
44 nCh.key = 'NbChannels'; | |
45 nCh.type = 'Integer'; | |
46 nCh.default = 10; | |
47 option.nCh = nCh; | |
48 | |
49 Ch.key = {'Channel','Channels'}; | |
50 Ch.type = 'Integer'; | |
51 Ch.default = 0; | |
52 option.Ch = Ch; | |
53 | |
54 lowF.key = 'Lowest'; | |
55 lowF.type = 'Integer'; | |
56 lowF.default = 50; | |
57 option.lowF = lowF; | |
58 | |
59 freq.key = 'Manual'; | |
60 freq.type = 'Integer'; | |
61 freq.default = NaN; | |
62 option.freq = freq; | |
63 | |
64 overlap.key = 'Hop'; | |
65 overlap.type = 'Boolean'; | |
66 overlap.default = 1; | |
67 option.overlap = overlap; | |
68 | |
69 filtertype.type = 'String'; | |
70 filtertype.choice = {'Gammatone','2Channels',0}; | |
71 filtertype.default = 'Gammatone'; | |
72 option.filtertype = filtertype; | |
73 | |
74 filterorder.key = 'Order'; | |
75 filterorder.type = 'Integer'; | |
76 filterorder.default = 4; | |
77 option.filterorder = filterorder; | |
78 | |
79 presel.type = 'String'; | |
80 presel.choice = {'Scheirer','Klapuri','Mel','Bark'}; | |
81 presel.default = ''; | |
82 option.presel = presel; | |
83 | |
84 specif.option = option; | |
85 | |
86 specif.eachchunk = @eachchunk; | |
87 specif.combinechunk = @combinechunk; | |
88 | |
89 varargout = mirfunction(@mirfilterbank,orig,varargin,nargout,specif,@init,@main); | |
90 | |
91 | |
92 function [x type] = init(x,option) | |
93 type = 'miraudio'; | |
94 | |
95 | |
96 function b = main(x,option,postoption) | |
97 if iscell(x) | |
98 x = x{1}; | |
99 end | |
100 f = get(x,'Sampling'); | |
101 | |
102 if strcmpi(option.presel,'Scheirer') | |
103 option.freq = [-Inf 200 400 800 1600 3200 Inf]; | |
104 elseif strcmpi(option.presel,'Klapuri') | |
105 option.freq = 44*[2.^([0:2,(9+(0:17))/3])]; | |
106 elseif strcmpi(option.presel,'Mel') | |
107 lowestFrequency = 133.3333; | |
108 linearFilters = 13; | |
109 linearSpacing = 66.66666666; | |
110 logFilters = 27; | |
111 logSpacing = 1.0711703; | |
112 totalFilters = linearFilters + logFilters; | |
113 cepstralCoefficients = 13; | |
114 option.freq = lowestFrequency + (0:linearFilters-1)*linearSpacing; | |
115 option.freq(linearFilters+1:totalFilters+2) = ... | |
116 option.freq(linearFilters) * logSpacing.^(1:logFilters+2); | |
117 | |
118 option.overlap = 2; | |
119 elseif strcmpi(option.presel,'Bark') | |
120 option.freq = [10 20 30 40 51 63 77 92 108 127 148 172 200 232 ... | |
121 270 315 370 440 530 640 770 950 1200 1550]*10; %% Hz | |
122 end | |
123 if not(isnan(option.freq)) | |
124 option.filtertype = 'Manual'; | |
125 end | |
126 | |
127 d = get(x,'Data'); | |
128 if size(d{1}{1},3) > 1 | |
129 warning('WARNING IN MIRFILTERBANK: The input data is already decomposed into channels. No more channel decomposition.'); | |
130 if option.Ch | |
131 cx = get(x,'Channels'); | |
132 db = cell(1,length(d)); | |
133 for k = 1:length(d) | |
134 db{k} = cell(1,length(d{k})); | |
135 for l = 1:length(d{k}) | |
136 for i = 1:length(option.Ch) | |
137 db{k}{l}(:,:,i) = d{k}{l}(:,:,find(cx{k} == option.Ch(i))); | |
138 end | |
139 end | |
140 end | |
141 b = set(x,'Data',db,'Channels',option.Ch); | |
142 else | |
143 b = x; | |
144 end | |
145 else | |
146 i = 1; | |
147 while i <= length(d) | |
148 if size(d{i}{1},2) > 1 | |
149 warning('WARNING IN MIRFILTERBANK: The frame decomposition should be performed after the filterbank decomposition.'); | |
150 disp(['Suggestion: Use the ' char(34) 'Frame' char(34) ' option available to some of the MIRtoolbox functions.']) | |
151 break | |
152 end | |
153 i = i+1; | |
154 end | |
155 nCh = option.nCh; | |
156 Ch = option.Ch; | |
157 if strcmpi(option.filtertype, 'Gammatone'); | |
158 [tmp x] = gettmp(x); %get(x,'Tmp'); | |
159 if not(Ch) | |
160 Ch = 1:nCh; | |
161 end | |
162 erb = cell(1,length(d)); | |
163 for i = 1:length(d) | |
164 erb{i} = cell(1,length(d{i})); | |
165 nch{i} = Ch; | |
166 for j = 1:length(d{i}) | |
167 try | |
168 coef = MakeERBFilters(f{i},nCh,option.lowF); | |
169 catch | |
170 error(['ERROR IN MIRFILTERBANK: Auditory Toolbox needs to be installed.']); | |
171 end | |
172 [erb{i}{j} tmp] = ERBfilterbank(d{i}{j},coef,Ch,tmp); | |
173 end | |
174 end | |
175 b = set(x,'Data',erb,'Channels',nch);%,'Tmp',tmp); | |
176 clear erb | |
177 b = settmp(b,tmp); | |
178 elseif strcmpi(option.filtertype, '2Channels'); | |
179 if not(Ch) | |
180 Ch = 1:2; | |
181 end | |
182 output = cell(1,length(d)); | |
183 try | |
184 for i = 1:length(d) | |
185 output{i} = cell(1,length(d{i})); | |
186 [bl,al] = butter(4,[70 1000]/f{i}*2); | |
187 k = 1; | |
188 if ismember(1,Ch) | |
189 Hdl = dfilt.df2t(bl,al); | |
190 %fvtool(Hdl) | |
191 Hdl.PersistentMemory = true; | |
192 for j = 1:length(d{i}) | |
193 low = filter(Hdl,d{i}{j}); | |
194 output{i}{j}(:,:,k) = low; | |
195 end | |
196 k = k+1; | |
197 end | |
198 if ismember(2,Ch) | |
199 if f{i} < 20000 | |
200 [bh,ah] = butter(2,1000/f{i}*2,'high'); | |
201 else | |
202 [bh,ah] = butter(2,[1000 10000]/f{i}*2); | |
203 end | |
204 Hdlh = dfilt.df2t(bl,al); | |
205 %fvtool(Hdlh) | |
206 Hdh = dfilt.df2t(bh,ah); | |
207 %fvtool(Hdh) | |
208 Hdlh.PersistentMemory = true; | |
209 Hdh.PersistentMemory = true; | |
210 for j = 1:length(d{i}) | |
211 high = filter(Hdh,d{i}{j}); | |
212 high = max(high,0); | |
213 high = filter(Hdlh,high); | |
214 output{i}{j}(:,:,k) = high; | |
215 end | |
216 end | |
217 nch{i} = Ch; | |
218 end | |
219 b = set(x,'Data',output,'Channels',nch); | |
220 catch | |
221 warning(['WARNING IN MIRFILTERBANK: Signal Processing Toolbox (version 6.2.1, or higher) not installed: no filterbank decomposition.']); | |
222 b = x; | |
223 end | |
224 elseif strcmpi(option.filtertype,'Manual'); | |
225 output = cell(1,length(d)); | |
226 for i = 1:length(d) | |
227 freqi = option.freq; | |
228 j = 1; | |
229 while j <= length(freqi) | |
230 if not(isinf(freqi(j))) && freqi(j)>f{i}/2 | |
231 if j == length(freqi) | |
232 freqi(j) = Inf; | |
233 else | |
234 freqi(j) = []; | |
235 j = j-1; | |
236 end | |
237 end | |
238 j = j+1; | |
239 end | |
240 output{i} = cell(1,length(d{i})); | |
241 step = option.overlap; | |
242 for j = 1:length(freqi)-step | |
243 if isinf(freqi(j)) | |
244 [z{j},p{j},k{j}] = ellip(option.filterorder,3,40,... | |
245 freqi(j+step)/f{i}*2); | |
246 elseif isinf(freqi(j+step)) | |
247 [z{j},p{j},k{j}] = ellip(option.filterorder,3,40,... | |
248 freqi(j)/f{i}*2,'high'); | |
249 else | |
250 [z{j},p{j},k{j}] = ellip(option.filterorder,3,40,... | |
251 freqi([j j+step])/f{i}*2); | |
252 end | |
253 end | |
254 for j = 1:length(z) | |
255 [sos,g] = zp2sos(z{j},p{j},k{j}); | |
256 Hd(j) = dfilt.df2tsos(sos,g); | |
257 Hd(j).PersistentMemory = true; | |
258 for h = 1:length(d{i}) | |
259 output{i}{h}(:,:,j) = filter(Hd(j),d{i}{h}); | |
260 end | |
261 end | |
262 %fvtool(Hd) | |
263 nch{i} = 1:length(freqi)-step; | |
264 end | |
265 b = set(x,'Data',output,'Channels',nch); | |
266 end | |
267 end | |
268 | |
269 | |
270 function [output tmp] = ERBfilterbank(x, fcoefs, chans, tmp) | |
271 % The following is based on the source code from Auditory Toolbox | |
272 % (A part that I could not call directly from MIRtoolbox, and that I | |
273 % generalize to chunk decomposition) | |
274 | |
275 % (Malcolm Slaney, August 1993, (c) 1998 Interval Research Corporation) | |
276 | |
277 % Process an input waveform with a gammatone filter bank. | |
278 % The fcoefs parameter, which completely specifies the Gammatone filterbank, | |
279 % should be designed with the MakeERBFilters function. | |
280 | |
281 A0 = fcoefs(:,1); | |
282 A11 = fcoefs(:,2); | |
283 A12 = fcoefs(:,3); | |
284 A13 = fcoefs(:,4); | |
285 A14 = fcoefs(:,5); | |
286 A2 = fcoefs(:,6); | |
287 B0 = fcoefs(:,7); | |
288 B1 = fcoefs(:,8); | |
289 B2 = fcoefs(:,9); | |
290 gain= fcoefs(:,10); | |
291 lc = length(chans); | |
292 output = zeros(size(x,1),size(x,2),lc); | |
293 if isempty(tmp) | |
294 emptytmp = 1; | |
295 else | |
296 emptytmp = 0; | |
297 end | |
298 for i = 1:lc | |
299 chan = length(gain)-chans(i)+1; % Revert the channels order | |
300 aa1 = [A0(chan)/gain(chan) A11(chan)/gain(chan) A2(chan)/gain(chan)]; | |
301 bb1 = [B0(chan) B1(chan) B2(chan)]; | |
302 aa2 = [A0(chan) A12(chan) A2(chan)]; | |
303 bb2 = [B0(chan) B1(chan) B2(chan)]; | |
304 aa3 = [A0(chan) A12(chan) A2(chan)]; | |
305 bb3 = [B0(chan) B1(chan) B2(chan)]; | |
306 aa4 = [A0(chan) A13(chan) A2(chan)]; | |
307 bb4 = [B0(chan) B1(chan) B2(chan)]; | |
308 if emptytmp | |
309 [y1 tmp(:,:,i,1)] = filter(aa1,bb1,x); | |
310 [y2 tmp(:,:,i,2)] = filter(aa2,bb2,y1); | |
311 [y3 tmp(:,:,i,3)] = filter(aa3,bb3,y2); | |
312 [y4 tmp(:,:,i,4)] = filter(aa4,bb4,y3); | |
313 else | |
314 [y1 tmp(:,:,i,1)] = filter(aa1,bb1,x,tmp(:,:,i,1)); | |
315 [y2 tmp(:,:,i,2)] = filter(aa2,bb2,y1,tmp(:,:,i,2)); | |
316 [y3 tmp(:,:,i,3)] = filter(aa3,bb3,y2,tmp(:,:,i,3)); | |
317 [y4 tmp(:,:,i,4)] = filter(aa4,bb4,y3,tmp(:,:,i,4)); | |
318 end | |
319 output(:,:,i) = y4; | |
320 end | |
321 | |
322 | |
323 function [y orig] = eachchunk(orig,option,missing) | |
324 y = mirfilterbank(orig,option); | |
325 | |
326 | |
327 function y = combinechunk(old,new) | |
328 do = get(old,'Data'); | |
329 to = get(old,'Time'); | |
330 dn = get(new,'Data'); | |
331 tn = get(new,'Time'); | |
332 y = set(old,'Data',{{[do{1}{1};dn{1}{1}]}},... | |
333 'Time',{{[to{1}{1};tn{1}{1}]}}); |