1function W=ref_wignervilledist(f,g) 2%-*- texinfo -*- 3%@deftypefn {Function} ref_wignervilledist 4%@verbatim 5%REF_WIGNERVILLEDIST Reference wigner-Ville distribution 6% Usage: W=ref_wignervilledist(f) 7% 8% REF_WIGNERVILLEDIST(f,g) computes the Wigner-Ville distribution of f and g. 9%@end verbatim 10%@strong{Url}: @url{http://ltfat.github.io/doc/reference/ref_wignervilledist.html} 11%@end deftypefn 12 13% Copyright (C) 2005-2016 Peter L. Soendergaard <peter@sonderport.dk>. 14% This file is part of LTFAT version 2.3.1 15% 16% This program is free software: you can redistribute it and/or modify 17% it under the terms of the GNU General Public License as published by 18% the Free Software Foundation, either version 3 of the License, or 19% (at your option) any later version. 20% 21% This program is distributed in the hope that it will be useful, 22% but WITHOUT ANY WARRANTY; without even the implied warranty of 23% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 24% GNU General Public License for more details. 25% 26% You should have received a copy of the GNU General Public License 27% along with this program. If not, see <http://www.gnu.org/licenses/>. 28 29% AUTHOR: Jordy van Velthoven 30 31if ~all(length(f)==length(g)) 32 error('%s: f and g must have the same length.', upper(mfilename)); 33end; 34 35L = length(f); 36H = floor(L/2); 37R = zeros(L,L); 38W = zeros(L,L); 39 40% Compute the analytic representation of f 41if (nargin == 1) 42 if isreal(f) 43 z = fft(f); 44 z(2:L-H) = 2*z(2:L-H); 45 z(H+2:L) = 0; 46 z1 = ifft(z); 47 z2 = z1; 48 else 49 z1 = f; 50 z2 = z1; 51 end 52elseif (nargin == 2) 53 if isreal(f) || isreal(g) 54 z1 = fft(f); 55 z1(2:L-H) = 2*z1(2:L-H); 56 z1(H+2:L) = 0; 57 z1 = ifft(z1); 58 z2 = fft(g); 59 z2(2:L-H) = 2*z2(2:L-H); 60 z2(H+2:L) = 0; 61 z2 = ifft(z2); 62 else 63 z1 = f; 64 z2 = g; 65 end 66end 67 68% Compute instantaneous autocorrelation matrix R 69for l = 0 : L-1; 70 for m = -min([L-l, l, round(L/2)-1]) : min([L-l, l, round(L/2)-1]); 71 R(mod(L+m,L)+1, l+1) = z1(mod(l+m, L)+1).*conj(z2(mod(l-m, L)+1)); 72 end 73end 74 75% Compute the Wigner-Ville distribution 76for hh=0:L-1 77 for ii=0:L-1 78 for jj = 0:L-1 79 W(hh+1, ii+1) = W(hh+1, ii+1) + R(jj+1, ii+1) .* exp(-2*pi*i*hh*jj/L); 80 end 81 end 82end 83 84