function y = var(x,w)
%VAR    Variance. 
%   For vectors, VAR(X) returns the variance of X.  
%   For matrices, VAR(X) is a row vector containing the
%   variance of each column of X. 
%
%   VAR(X) normalizes by N-1 where N is the sequence length.  This
%   makes VAR(X) the best unbiased estimate of the variance if X
%   is a sample from a normal distribution.
%
%   VAR(X,1) normalizes by N and produces the second moment of the
%   sample about its mean.
%
%   VAR(X,W) computes the variance using the weight vector, W.
%   The number of elements in W must equal the number of rows 
%   in X unless W = 1 which is treated as a short-cut for a 
%   vector of ones.
%
%   The elements of W must be positive. VAR normalizes W by
%   dividing each element in W by the sum of all its elements.
%
%   The variance is the square of the standard deviation
%   (STD). 
%
%   See also MEAN, STD, COV and CORRCOEF.

%   VAR supports both common definitions of variance simply.
%   Let SS be the sum of the squared deviations elements of a 
%   vector X from their mean. The default behavior (no weight vector) is:
%   VAR(X) = SS/(N-1) where N is the number of elements in X. But,
%   VAR(X,1) = SS/N

%   B.A. Jones 9-30-92
%   Revised 2-13-95 BJ
%   Copyright 1984-2000 The MathWorks, Inc. 
%   $Revision: 1.5 $  $Date: 2000/06/01 02:10:46 $

if min(size(x)) == 1
    x = x(:);
end

[m,n] = size(x);

if isempty(x)
   y = NaN;
   return
end

if m == 1
    y = 0;
    return
end

if nargin == 2
    if min(size(w)) ~= 1
        error('The second argument must be a vector.');
    else
        w = w(:);
    end
    k = find(w < 0);
    if any(k)
        error('The second argument must be positive.');
    end

   if w == 1
      w = w(ones(m,1),1);
   elseif max(size(w)) ~= m
       error('The sizes of the arguments are not compatible.');
    end

    w = w/sum(w);
else
    w = 1/(m - 1);
end 

if nargin == 2
   w = w(:,ones(n,1));
   avg = sum(w.*x);
else
   avg = mean(x);
end

if n == 1
   cx = x - avg; 
   y = sum((conj(cx).*cx).*w);
else
   centerx = x - avg(ones(m,1),:);
   y = sum((conj(centerx) .* centerx).*w);
end
