DM 'LOG;CLEAR;OUT;CLEAR;';
OPTIONS NODATE NONUMBER LS=80 PS=60;

***********************************************************;
***  Operating Characteristic (OC) Curve for a p Chart  ***;
***********************************************************;

 data circ_one;
   input batch fail @@;
   cards;
  1   7   2   6   3   6   4   9    5   2
  6  11   7   8   8   8   9   6   10  19
 11   7  12   5  13   7  14   5   15   8
 16  13  17   7  18  14  19  19   20   5
 21   7  22   5  23   7  24   5   25  11
 26   4  27   6  28   3  29  11   30   3
 ;

 proc shewhart data=circ_one;
    where batch^=10 and batch^=19;
    pchart fail*batch / subgroupn = 500
                        nochart
                        outindex  = 'Revised Limits'
                        outlimits = faillim2;
 run;

 data ocpchart;
    set faillim2;
    keep beta fraction;
    nucl=_limitn_*_uclp_;
    nlcl=_limitn_*_lclp_;
    do p=0 to 500;
       fraction=p/1000;
       if nucl=floor(nucl) then
          adjust=probbnml(fraction,_limitn_,nucl) -
                 probbnml(fraction,_limitn_,nucl-1);
       else adjust=0;
       if nlcl=0 then
          beta=1 - probbeta(fraction,nucl,_limitn_-nucl+1)+adjust;
       else beta=probbeta(fraction,nlcl,_limitn_-nlcl+1) -
                 probbeta(fraction,nucl,_limitn_-nucl+1) +
                 adjust;
       if beta >= 0.001 then output;
       end;
    call symput('lcl', put(_lclp_,5.3));
    call symput('mean',put(_p_,   5.3));
    call symput('ucl', put(_uclp_,5.3));
 run;

 title "OC Curve for p Chart With LCL=&LCL, p0=&MEAN and UCL=&UCL";
 symbol i=j w=2 v=none c=white;
 proc gplot data=ocpchart;
    plot beta*fraction /
       vaxis  = axis1
       haxis  = axis2
       cframe = gray
       frame
       autovref
       autohref
       lvref  = 2
       lhref  = 2
       vzero
       hzero;
    label fraction='Fraction Nonconforming'
          beta ='Beta';
 axis1 offset=(0,.5) minor=none order=0 to 1.0 by 0.1;
 axis2 offset=(0,0)  minor=none order=0 to 0.06 by 0.005;
 run;