function ospf=ostat(pf, cf, n, j) 
% ostat(pf, cf, n, j)  generates the probability function for the jth of n 
%                      order statistic of the probability and cumulative 
%                      functions 
% 
%   If you'd like nice results, the functions had better be matched 
if j>n 
   error('j cannot be >n') 
end 
if j<1 
   error('j cannot be <1') 
end 
if ispf(pf)==0 
   error('invalid probability function') 
end 
if iscf(cf)==0 
   error('invalid cumulative function') 
end 
x=xvals(pf); 
p=pvals(pf); 
P=pvals(cf); 
nc=[n-j j-1 1]; 
if isdiscrt(pf)==0 
   newp=((P.^(j-1)).*((1-P).^(n-j)).*p).*multinom(n, nc); 
end 
ospf=[pf(1,1) 1; x' newp']; 
if ispf(ospf)==0 
   newp=newp./ttlprob(ospf); 
   ospf=[pf(1,1) 1; x' newp']; 
end