  

%wsurtexact
  n=25;
  n1=15;
  s=sum(combnk(1:n,n1),2);
  b = unique(s');
for i = 1:length(b)
fb(i) = sum(ismember(s,b(i)));
end
pro=fb/nchoosek(n, n1);
cp=cumsum(pro);
%plots
close all force
figure(1)
plot(b, cp, 'k-')
hold on
plot(b, normcdf(b, n1*(n+1)/2, sqrt(n1*(n-n1)*(n+1)/12)),'r-')
plot(b, normcdf(b+1/2, n1*(n+1)/2, sqrt(n1*(n-n1)*(n+1)/12)),'g-')

figure(2)
b1=b(60:70);
cp1=cp(60:70);
plot(b1, cp1, 'k-')
hold on
plot(b1, normcdf(b1, n1*(n+1)/2, sqrt(n1*(n-n1)*(n+1)/12)),'r-')
plot(b1, normcdf(b1+1/2, n1*(n+1)/2, sqrt(n1*(n-n1)*(n+1)/12)),'g-')
legend('1','2','3','location','NW')
figure(3)
stem(b, pro,'fill','-');
hold on
plot(b, normpdf(b, n1*(n+1)/2, sqrt(n1*(n-n1)*(n+1)/12)),'r-')

