comparison Code/Descriptors/Matlab/MPEG7/FromWeb/VoiceSauce/func_GetA1A2A3.m @ 4:92ca03a8fa99 tip

Update to ICASSP 2013 benchmark

+function [A1, A2, A3] = func_GetA1A2A3(y, Fs, F0, F1, F2, F3, variables, textgridfile)
+% [A1, A2, A3] = func_GetA1A2A3(y, Fs, F0, F1, F2, F3, variables, textgridfile)
+% Input:  y - input signal
+%         Fs - sampling frequency
+%         F0 - fundamental frequency (vector)
+%         F1-F3 - formant frequencies (vectors)
+%         variables - settings from VS
+%         textgridfile - textgrid filename if available (optional)
+% Output: A1, A2, A3 - uncorrected amplitudes of the spectra at the formant
+% frequencies.
+%
+% Author: Yen-Liang Shue, Speech Processing and Auditory Perception Laboratory, UCLA
+% Copyright UCLA SPAPL 2009
+
+N_periods = variables.Nperiods;
+sampleshift = (Fs / 1000 * variables.frameshift);
+
+A1 = zeros(length(F0), 1) * NaN;
+A2 = zeros(length(F0), 1) * NaN;
+A3 = zeros(length(F0), 1) * NaN;
+
+if (nargin == 7)
+
+    for k=1:length(F0)
+        ks = round(k * sampleshift);
+        
+        if (ks <= 0 || ks > length(y))
+            continue;
+        end
+        
+        F0_curr = F0(k);
+        
+        if (isnan(F0_curr) || F0_curr == 0) 
+            continue; 
+        end
+
+        N0_curr = 1 / F0_curr * Fs;
+               
+        ystart = round(ks - N_periods/2*N0_curr);
+        yend = round(ks + N_periods/2*N0_curr) - 1;
+        
+        if (ystart <= 0) 
+            continue; 
+        end;
+        
+        if (yend > length(y)) 
+            continue; 
+        end;
+        
+        yseg = y(ystart:yend);
+        
+        % proceed if no NaNs detected
+        if (~isnan(F1(k)) && ~isnan(F2(k)) && ~isnan(F3(k)) && ~isnan(N0_curr))
+            [A1_, fmax] = ana_GetMagnitudeMax(yseg, F1(k), Fs, 8192);
+            [A2_, fmax] = ana_GetMagnitudeMax(yseg, F2(k), Fs, 8192);
+            [A3_, fmax] = ana_GetMagnitudeMax(yseg, F3(k), Fs, 8192);
+
+            A1(k) = A1_;
+            A2(k) = A2_;
+            A3(k) = A3_;
+        end    
+    end
+else
+    % use textgrids - get ignore string list
+    tbuffer = variables.tbuffer;
+    ignorelabels = textscan(variables.TextgridIgnoreList, '%s', 'delimiter', ',');
+    ignorelabels = ignorelabels{1};   
+    
+    [labels, start, stop] = func_readTextgrid(textgridfile);
+    labels_tmp = [];
+    start_tmp = [];
+    stop_tmp = [];
+    
+    for k=1:length(variables.TextgridTierNumber)
+        inx = variables.TextgridTierNumber(k);
+        if (inx <= length(labels))
+            labels_tmp = [labels_tmp; labels{inx}];
+            start_tmp = [start_tmp; start{inx}];
+            stop_tmp = [stop_tmp; stop{inx}];
+        end
+    end
+    
+    labels = labels_tmp;
+    start = start_tmp;
+    stop = stop_tmp; 
+    
+    
+    for m=1:length(start)
+        switch(labels{m})
+            case ignorelabels
+                continue;  % skip anything that is within the ignore list
+        end
+        
+        kstart = round((start(m) * 1000 - tbuffer) / variables.frameshift);
+        kstop = round((stop(m) * 1000 + tbuffer) / variables.frameshift);
+        kstart(kstart <= 0) = 1;
+        kstop(kstop > length(F0)) = length(F0);
+        
+        ystart = round(kstart * sampleshift);
+        ystop = round(kstop * sampleshift);
+        ystart(ystart <= 0) = 1;
+        ystop(ystop > length(y)) = length(y);
+        
+        ysegment = y(ystart:ystop);
+        F0seg = F0(kstart:kstop);
+        F1seg = F1(kstart:kstop);
+        F2seg = F2(kstart:kstop);
+        F3seg = F3(kstart:kstop);
+               
+        for k=1:length(F0seg)
+            ks = round(k*sampleshift);
+            
+            if (ks <= 0 || ks > length(ysegment))
+                continue;
+            end
+            
+            F0_curr = F0seg(k);
+            if (isnan(F0_curr) || F0_curr == 0)
+                continue;
+            end
+            
+            N0_curr = Fs/F0_curr;
+            ysstart = round(ks - N_periods/2 * N0_curr);
+            ysend = round(ks + N_periods/2 * N0_curr) - 1;
+            
+            if (isnan(F1seg(k)) || isnan(F2seg(k)) || isnan(F3seg(k)))
+                continue;
+            end
+            
+            if (ysstart <= 0) 
+                continue; 
+            end;
+            if (ysend > length(ysegment)) 
+                continue; 
+            end;
+            
+            yseg = ysegment(ysstart:ysend);
+            
+            [A1_, fmax] = ana_GetMagnitudeMax(yseg, F1seg(k), Fs, 8192);
+            [A2_, fmax] = ana_GetMagnitudeMax(yseg, F2seg(k), Fs, 8192);
+            [A3_, fmax] = ana_GetMagnitudeMax(yseg, F3seg(k), Fs, 8192);
+            
+            inx = kstart + k - 1;
+            if (inx <= length(A1))
+                A1(inx) = A1_;
+                A2(inx) = A2_;
+                A3(inx) = A3_;
+            end                
+        
+        end
+    
+    end
+end
+
+
+%------------------------------------------------------------------------
+function [M,fmax] = ana_GetMagnitudeMax(x, Fx, Fs, fftlen)
+% Get maximal spectral magnitude of signal x around position Fx Hz in dB
+% Fx can be a vector of frequency points
+% Note that the bigger fftlen the better the resolution
+
+if (isnan(Fx))
+  M = NaN;
+  fmax = NaN;
+else
+
+  len = length(x);
+  hamlen = min(fftlen, len);
+  %X = fft(hamming(hamlen).*x(1:hamlen), fftlen);
+  factor = 1;  %/length(x);  % not exactly Kosher but compatible to dfs_magn()
+  X = fft(x,fftlen);
+  X(X==0) = 0.000000001;  % guard against log(0) warnings
+  X = 20*log10(factor*abs(X(1:fftlen/2, :)));
+  fstep = Fs / fftlen;
+  lowf = Fx-0.1*Fx;
+  if (lowf<0) lowf=0; end
+  highf = Fx+0.1*Fx;
+  if (highf>Fs/2-fstep) highf=Fs/2-fstep; end
+
+  for cnt = 1:length(Fx)
+      [M(cnt),pos] = max(X(1+round(lowf(cnt)/fstep):1+round(highf(cnt)/fstep), :));
+      fmax(cnt) = (pos-1+round(lowf(cnt)/fstep))*fstep;
+  end
+  
+end
+
+% f=0:fstep:(Fs/2-fstep);
+% plot(f,X); hold;
+% plot(pos, M, 'rx'); hold off;
+% fprintf('End get_magnitude_max()\n');