/* Richard A. DeVenezia
 * www.devenezia.com
 * 5/16/2004
 */

ods listing ;

%let pixels=600;
%let pixels=143;

goptions reset=all target=png device=png xpixels=&pixels ypixels=&pixels;
goptions gsfname=gout;

filename gout "%sysfunc(pathname(WORK))\bezier-curve.png";
filename gout "\\extreme\samples\bezier-curve.png";

data _null_;

 /* Adapted from http://www.codetoad.com/vb_bezier.asp
  * A cubic Bezier curve is defined by four points.
  *
  * (x0,y0) & (x3,y3) are endpoints and
  * (x1,y1) & (x2,y2) are control points.
  *
  * The following equations define the points
  * on the curve.  Both are evaluated for an
  * arbitrary number of values of t between 0 and 1.
  *
  * X(t) = ax * t ^ 3 + bx * t ^ 2 + cx * t + x0
  *
  * X1 = x0 + cx / 3
  * X2 = X1 + (cx + bx) / 3
  * x3 = x0 + cx + bx + ax
  *
  * Y(t) = ay * t ^ 3 + by * t ^ 2 + cy * t + y0
  *
  * Y1 = y0 + cy / 3
  * Y2 = Y1 + (cy + by) / 3
  * y3 = y0 + cy + by + ay
  */

  libname = 'WORK';
  memname = 'BEZIER';
  entname = 'CURVE';

  rc = gset ('CATALOG', libname, memname);

  rc = ginit();
  rc = graph('CLEAR', entname);

  f = 1;

  rc = gset ('WINDOW', 1, -f,-f, f,f);
  rc = gset ('TRANSNO', 1);

  rc = gset ('COLREP', 1, 'CXAAAAAA');
  rc = gset ('COLREP', 2, 'BLACK');

  rc = gset ('COLREP', 1, 'CXCCCCFF');
  rc = gset ('COLREP', 2, 'CX0000FF');

  r = .95 * f;

  rc = gset ('LINCOLOR', 2);
* rc = gset ('LINWIDTH', 1);
* rc = gset ('LINTYPE', 1);
  rc = gdraw('ARC', 0,0,r, 0,360);

  rc = gset ('LINCOLOR', 1);
* rc = gset ('LINWIDTH', 1);
* rc = gset ('LINTYPE', 1);

  rc = gset ('ASF', 'LINWIDTH', 'INDIVIDUAL');

  nVerts=12;

  do i = 1 to nVerts;

    theta1 = 2 * constant('PI') * (i-1) / nVerts;
    theta2 = 2 * constant('PI') * (i  ) / nVerts;

    x0 = r * cos(theta1) ;
    y0 = r * sin(theta1) ;
    x3 = r * cos(theta2) ;
    y3 = r * sin(theta2) ;

    * aff - arc flattening factor none=0..1=straight;
    do aff = 0 to 1-1e-6 by .25;
      link drawBez;
    end;
  end;

  rc = graph('UPDATE');*,'NOSHOW');
  rc = gterm();

  stop;

drawBez:
  array x[0:3] x0-x3;
  array y[0:3] y0-y3;

  x1 = aff * (x0 + x3) / 2;
  y1 = aff * (y0 + y3) / 2;
  x2 = x1;
  y2 = y1;

  cx = 3 * (x(1) - x(0));
  bx = 3 * (x(2) - x(1)) - cx;
  ax = x(3) - x(0) - cx - bx;

  cy = 3 * (y(1) - y(0));
  by = 3 * (y(2) - y(1)) - cy;
  ay = y(3) - y(0) - cy - by;


  xtp = x(0);
  ytp = y(0);

  iterations = 20;

  do t = 0 to 1 by 1/iterations;

    xt = ax * t ** 3 + bx * t ** 2 + cx * t + x(0);
    yt = ay * t ** 3 + by * t ** 2 + cy * t + y(0);

*   rc = gdraw('LINE', 2, xtp,xtp, ytp,ytp);
    rc = gdraw('LINE', 2, xtp,xt , ytp,yt );

    xtp = xt;
    ytp = yt;
  end;

return;

run;

options noxwait;
x "start %sysfunc(pathname(gout))";