function p = parzen_pdf(X,Samples,sigma,W)
% p = parzen_pdf(X,Samples,sigma)
% p = parzen_pdf(X,Samples,sigma,Weights)
%
% Parzen estimator of 1-D distribution F, from which Samples have
% been drawn, over the domain X.
%
% X       = 1-d vector of pdf domain values.
% Samples = [Nsamples,Npdfs] vector
% sigma   = scalar width of (gaussian) parzen window.
% Weights = [Nsamples,1] vector, weighting to be applied to each sample.
% returns p = [d,Npdfs] vector of pdf over X.

% Liam Pedersen, CMU, 1998
% revised Liam Pedersen, NASA ARC, August 1999 to include weights.

[N, NP] = size(Samples);

if (nargin < 4) W = ones(N,1); end

p = zeros(length(X),NP);
dx = diff(X);
dx = [dx 0]; 
for i = 1:N
   for j = 1:NP
      ParzenWindow = gaussian(X,Samples(i,j),sigma);
   	I = sum(dx.*ParzenWindow);
      p(:,j) = p(:,j) + W(i)*ParzenWindow'/I;
   end	
end
for j = 1:NP
   S = sum(dx.*p(:,j)');
   p(:,j) = p(:,j)/S;
end
