comparison src/matlab/ifgram.m @ 0:c52bc3e8d3ad tip

user: boblsturm branch 'default' added README.md added assets/.DS_Store added assets/playButton.jpg added assets/stopButton.png added assets/swapButton.jpg added data/.DS_Store added data/fiveoctaves.mp3 added data/glock2.wav added data/sinScale.mp3 added data/speech_female.mp3 added data/sweep.wav added nimfks.m.lnk added src/.DS_Store added src/matlab/.DS_Store added src/matlab/AnalysisCache.m added src/matlab/CSS.m added src/matlab/DataHash.m added src/matlab/ExistsInCache.m added src/matlab/KLDivCost.m added src/matlab/LoadFromCache.m added src/matlab/SA_B_NMF.m added src/matlab/SaveInCache.m added src/matlab/Sound.m added src/matlab/SynthesisCache.m added src/matlab/chromagram_E.m added src/matlab/chromagram_IF.m added src/matlab/chromagram_P.m added src/matlab/chromsynth.m added src/matlab/computeSTFTFeat.m added src/matlab/controller.m added src/matlab/decibelSliderReleaseCallback.m added src/matlab/drawClickCallBack.m added src/matlab/fft2chromamx.m added src/matlab/hz2octs.m added src/matlab/ifgram.m added src/matlab/ifptrack.m added src/matlab/istft.m added src/matlab/nimfks.fig added src/matlab/nimfks.m added src/matlab/nmfFn.m added src/matlab/nmf_beta.m added src/matlab/nmf_divergence.m added src/matlab/nmf_euclidean.m added src/matlab/prune_corpus.m added src/matlab/rot_kernel.m added src/matlab/templateAdditionResynth.m added src/matlab/templateDelCb.m added src/matlab/templateScrollCb.m
author boblsturm
date Sun, 18 Jun 2017 06:26:13 -0400
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:c52bc3e8d3ad
1 function [F,D] = ifgram(X, N, W, H, SR)
2 % [F,D] = ifgram(X, N, W, H, SR) Instantaneous frequency by phase deriv.
3 % X is a 1-D signal. Process with N-point FFTs applying a W-point
4 % window, stepping by H points; return (N/2)+1 channels with the
5 % instantaneous frequency (as a proportion of the sampling rate)
6 % obtained as the time-derivative of the phase of the complex spectrum
7 % as described by Toshihiro Abe et al in ICASSP'95, Eurospeech'97
8 % Same arguments and some common code as dpwebox/stft.m.
9 % Calculates regular STFT as side effect - returned in D.
10 % after 1998may02 dpwe@icsi.berkeley.edu
11 % 2001-03-05 dpwe@ee.columbia.edu revised version
12 % 2001-12-13 dpwe@ee.columbia.edu Fixed to work when N != W
13 % $Header: $
14
15 % Copyright (c) 2006 Columbia University.
16 %
17 % This file is part of LabROSA-coversongID
18 %
19 % LabROSA-coversongID is free software; you can redistribute it and/or modify
20 % it under the terms of the GNU General Public License version 2 as
21 % published by the Free Software Foundation.
22 %
23 % LabROSA-coversongID is distributed in the hope that it will be useful, but
24 % WITHOUT ANY WARRANTY; without even the implied warranty of
25 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
26 % General Public License for more details.
27 %
28 % You should have received a copy of the GNU General Public License
29 % along with LabROSA-coversongID; if not, write to the Free Software
30 % Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA
31 % 02110-1301 USA
32 %
33 % See the file "COPYING" for the text of the license.
34
35 if nargin < 2; N = 256; end
36 if nargin < 3; W = N; end
37 if nargin < 4; H = W/2; end
38 if nargin < 5; SR = 1; end
39
40 s = length(X);
41 % Make sure it's a single row
42 if size(X,1) > 1
43 X = X';
44 end
45
46 %win = [0,hanning(W-1)'];
47 win = 0.5*(1-cos([0:(W-1)]/W*2*pi));
48
49 % Window for discrete differentiation
50 T = W/SR;
51 dwin = -pi / T * sin([0:(W-1)]/W*2*pi);
52
53 % sum(win) takes out integration due to window, 2 compensates for neg frq
54 norm = 2/sum(win);
55
56 % How many complete windows?
57 nhops = 1 + floor((s - W)/H);
58
59 F = zeros(1 + N/2, nhops);
60 D = zeros(1 + N/2, nhops);
61
62 nmw1 = floor( (N-W)/2 );
63 nmw2 = N-W - nmw1;
64
65 ww = 2*pi*[0:(N-1)]*SR/N;
66
67 for h = 1:nhops
68 u = X((h-1)*H + [1:W]);
69 % if(h==0)
70 % plot(u)
71 % end
72 % Apply windows now, while the length is right
73 wu = win.*u;
74 du = dwin.*u;
75
76 % Pad or truncate samples if N != W
77 if N > W
78 wu = [zeros(1,nmw1),wu,zeros(1,nmw2)];
79 du = [zeros(1,nmw1),du,zeros(1,nmw2)];
80 end
81 if N < W
82 wu = wu(-nmw1+[1:N]);
83 du = du(-nmw1+[1:N]);
84 end
85 % FFTs of straight samples plus differential-weighted ones
86 t1 = fft(fftshift(du));
87 t2 = fft(fftshift(wu));
88 % Scale down to factor out length & window effects
89 D(:,h) = t2(1:(1 + N/2))'*norm;
90
91 % Calculate instantaneous frequency from phase of differential spectrum
92 t = t1 + j*(ww.*t2);
93 a = real(t2);
94 b = imag(t2);
95 da = real(t);
96 db = imag(t);
97 instf = (1/(2*pi))*(a.*db - b.*da)./((a.*a + b.*b)+(abs(t2)==0));
98 % 1/2pi converts rad/s into cycles/s
99 % sampling rate already factored in as constant in dwin & ww
100 % so result is in Hz
101
102 F(:,h) = instf(1:(1 + N/2))';
103
104 end;
105