view dsp/wigner.m @ 61:eff6bddf82e3 tip

Finally implemented perceptual brightness thing.

function wig = wigner(sig,m,l,k,win)
% Wigner -- Wigner-Ville Distribution. This version
%           operates on an entire signal in a vector. 
%  Usage
%    wig = wigner(sig,m,l,k,win)
%  Inputs
%    sig    1-d signal
%    m      half window length (must be power of 2)
%    l      number of frequencies (ie rows in output)
%    k      window step (defaults to 1)
%    win    half window of size m+1
%  Outputs
%   wig     complex-valued matrix representing the Wigner-Ville
%           distribution of zero-extended signal with rows corresponding
%           to frequencies and columns corresponding to times.
%
%  Side Effects
%   Image Plot of the Wigner-Ville distribution
%

sig = sig(:);
n   = length(sig);

if nargin<4 k   = 1; end
if nargin<3 l   = m+1; end
if nargin<5 win = ones(m+1,1); end
   
f   = [zeros(m,1); sig; zeros(m,1)];
wig = zeros(l,floor(n/k));
ix  = 0:m;

col=1;
for t=1:k:n,
    tplus    = m + t + ix;
    tminus   = m + t - ix;
    R = win .* ( f(tplus) .* f(tminus) );
    x = [ R(1:m); R(m+1:-1:2) ];
    y = fft(x);
    wig(:,col) = real(y(1:l));
    col = col+1;
end