annotate Code/Descriptors/Matlab/MPEG7/FromWeb/VoiceSauce/func_GetCPP.m @ 4:92ca03a8fa99 tip

Update to ICASSP 2013 benchmark
author Dawn Black
date Wed, 13 Feb 2013 11:02:39 +0000
parents
children
rev   line source
Dawn@4 1 function CPP = func_GetCPP(y, Fs, F0, variables)
Dawn@4 2 % CPP = func_GetCPP(y, Fs, F0, variables)
Dawn@4 3 % Input: y, Fs - from wavread
Dawn@4 4 % F0 - vector of fundamental frequencies
Dawn@4 5 % variables - global settings
Dawn@4 6 % Output: CPP vector
Dawn@4 7 % Notes: Calculates the CPP according to the formula from Hillenbrand, 1994.
Dawn@4 8 % Currently does not use textgrid information. For very very long files this
Dawn@4 9 % may be a process-bottleneck.
Dawn@4 10 %
Dawn@4 11 % Author: Yen-Liang Shue, Speech Processing and Auditory Perception Laboratory, UCLA
Dawn@4 12 % Copyright UCLA SPAPL 2009
Dawn@4 13
Dawn@4 14 N_periods = variables.Nperiods_EC;
Dawn@4 15 sampleshift = (Fs / 1000 * variables.frameshift);
Dawn@4 16
Dawn@4 17 CPP = zeros(length(F0), 1) * NaN;
Dawn@4 18
Dawn@4 19 N_ms = round(Fs / 1000); % Quefrency below 1ms are ignored as per Hillenbrand
Dawn@4 20
Dawn@4 21 for k=1:length(F0)
Dawn@4 22 ks = round(k * sampleshift);
Dawn@4 23
Dawn@4 24 if (ks <= 0 || ks > length(y))
Dawn@4 25 continue;
Dawn@4 26 end
Dawn@4 27
Dawn@4 28 F0_curr = F0(k);
Dawn@4 29
Dawn@4 30 if (isnan(F0_curr))
Dawn@4 31 continue;
Dawn@4 32 end
Dawn@4 33
Dawn@4 34 N0_curr = Fs / F0_curr;
Dawn@4 35
Dawn@4 36 ystart = round(ks - N_periods/2*N0_curr);
Dawn@4 37 yend = round(ks + N_periods/2*N0_curr) - 1;
Dawn@4 38
Dawn@4 39 if (ystart <= 0)
Dawn@4 40 continue;
Dawn@4 41 end;
Dawn@4 42
Dawn@4 43 if (yend > length(y))
Dawn@4 44 continue;
Dawn@4 45 end;
Dawn@4 46
Dawn@4 47 yseg = y(ystart:yend);
Dawn@4 48 yseg = yseg .* hamming(length(yseg)); %windowing
Dawn@4 49 YSEG = fft(yseg);
Dawn@4 50 yseg_c = ifft(log(abs(YSEG)));
Dawn@4 51 yseg_c_db = 10*log10(yseg_c .^ 2);
Dawn@4 52
Dawn@4 53 yseg_c_db = yseg_c_db(1:floor(length(yseg)/2));
Dawn@4 54 [vals, inx] = func_pickpeaks(yseg_c_db(N_ms:end)', 2*N0_curr);
Dawn@4 55 [pos, pinx] = min(abs(inx - N0_curr));
Dawn@4 56
Dawn@4 57 if (~isempty(pinx))
Dawn@4 58 inx = inx(pinx(1));
Dawn@4 59 vals = yseg_c_db(inx+N_ms-1);
Dawn@4 60 p = polyfit([N_ms:length(yseg_c_db)], yseg_c_db(N_ms:end)', 1);
Dawn@4 61 base_val = p(1) * (inx+N_ms-1) + p(2);
Dawn@4 62 CPP(k) = vals - base_val;
Dawn@4 63
Dawn@4 64 % figure;
Dawn@4 65 % plot(yseg_c_db); hold on;
Dawn@4 66 % plot(inx+N_ms-1, vals, 'rx');
Dawn@4 67 % plot([N_ms:length(yseg_c_db)], polyval(p,[N_ms:length(yseg_c_db)]), 'g');
Dawn@4 68
Dawn@4 69 end
Dawn@4 70
Dawn@4 71 end
Dawn@4 72