view ddm_example_1.m @ 4:72c011ed1977 tip

more elaborate example with non-stat. estimate explanation

input_file = 'path_to_file';
%[sig,fs] = wavread(input_file);

%%%%% DEMO %%%%%
fs = 44100;
sig = rand(100000,1);
%%%%% DEMO %%%%%

% frame length
N = 2^9;
% DDM internal FFT size, for zero-padding, must be => N
N_fft = N;
% quantization, can be anything, but N_q = N makes the most sense
N_q = N;
% degree 2 DDM (linear FM, quadratic AM)
ddm_dgr = 2;
% DDM bandwidth (in bins)
nr_bins = 5;
% AM/FM model functions - 2nd degree polynomial
dsc_idxs = [-(N-1)/2:(N-1)/2];
t = dsc_idxs / fs;
poly_model = poly_mat(t,ddm_dgr,N);

% hann and hann' window
w = 0.5+0.5*cos(2*pi*dsc_idxs/(N-1));
wd = -pi/fs*(N-1)*sin(2*pi*dsc_idxs/(N-1));
% example
offset = 306; % some sample offset in the audio file
bffr = sig(offset:offset+N-1); % get the buffer of lentgh N 


%%% execute DDM on the whole frequency range (FFT) %%%
% create linear systems
[A, b, sig_fft] = ddm_lin_sys_fft(nr_bins, ddm_dgr, w', wd', poly_model', bffr, N, N_fft, fs);
% fast pseudo-inverse
ddm_estimates = flipud(pinv_2_fast(A,b,nr_bins,N-nr_bins+1).');

% phase function estimates
ddm_frequency_estimates_for_each_fft_bin = imag(ddm_estimates(2,:));
ddm_linear_fm_estimates_for_each_fft_bin = imag(ddm_estimates(1,:));
% amplitude funciton estimates
ddm_linear_am_estimates_for_each_fft_bin = real(ddm_estimates(2,:));
ddm_2nd_dgr_am_estimates_for_each_fft_bin = real(ddm_estimates(1,:));

% re-assign the spectral magnitude 
ddm_spectrogram = reallocate(abs(sig_fft), ddm_estimates, N_q, fs);

% find the peak
[peak_val peak_bin] = max(ddm_spectrogram);

% frequency trajectory of the peak bin
frequency_trajectory = ddm_linear_fm_estimates_for_each_fft_bin(peak_bin) * t + ddm_frequency_estimates_for_each_fft_bin(peak_bin);