%pdfdat.m
% M-file to calculate to prob.density functon(pdf) of data
% Then does a Chi Square goodness of fit.
  mu=mean(data);
  stdd=std(data);
  z=(data-mu)/stdd;    %use the standard normal variable,z

  numbins=14;
  [N,x]=hist(z,numbins);    % do a histogram over 14 bins
  nlength=sum(N);
  l=length(data);
  cdiff=diff(x);
  c=cdiff(1);
  Norm=(1/c)*N/(l); 
%  ii=find(N>=5);      % do chisquare when at least 5 values in each bin
%  ii=find(N>=2);      % do chisquare when at least 2 values in each bin
 ii=find(N>=0.02*l*c);   % do chisquare when at least 2% in each bin
  Nrm=(1/c)*N(ii)/(l);  % change to normalized relative freq histogram
  fi = Nrm;           % normalize so sum(fi*c)==1

  x1=x(1)-c;
  xl=x(numbins)+c;
  xx=[x1 x xl];

  gaussx=normpdf(xx,0,1);

  %h = bar(x,normpdf(x,0,1));
  Fi = gaussian(x)';
  Fi = Fi(ii);

Fi=Fi*l*c;
fi=fi*l*c;


iii = find(Fi>=5);  
i2 = find(Fi<5);
%iii = 1:length(Fi);
%i2 = [];    
i3 = i2(find(i2<(length(Fi)/2)));
i4 = i2(find(i2>=(length(Fi)/2)));


first=iii(1);				% Add up small outliers in Fi and fi
last=iii(length(iii));
if (~isempty(i4)) & (sum(Fi(i4))<5)         
 Fi(last)=Fi(last)+sum(Fi(i4));
 Fi(i4)=[];
 fi(last)=fi(last)+sum(fi(i4));                        
 fi(i4)=[];
elseif (~isempty(i4))
 Fi(last+1)=sum(Fi(i4));
 Fi=Fi(1:last+1);
 fi(last+1)=sum(fi(i4));
 fi=fi(1:last+1);
end

if (~isempty(i3)) & (sum(Fi(i3))<5)         
 Fi(first)=Fi(first)+sum(Fi(i3));
 Fi(i3)=[];
 fi(first)=fi(first)+sum(fi(i3));                        
 fi(i3)=[];
elseif (~isempty(i3))
 Fi(first-1)=sum(Fi(i3));
 Fi=Fi((first-1):(length(Fi)));
 fi(first-1)=sum(fi(i3));
 fi=fi((first-1):(length(fi)));
end


  Xsquare=sum(((fi(:)-Fi(:)).^2)./Fi);   

fi=fi/l/c;
Fi=Fi/l/c;
  nbin = length(Fi);
  
%the areas should be 1 in a pdf - check this
  dx = c;
  Sfi = sum(fi*c);
  SFi = sum(Fi*c);
