# include <kapa_internal.h>

static int bWeight;
static int bType;
static bDrawColor bColor;
static bDrawBuffer *bBuffer;
void bDrawCircleSingle (double xc, double yc, double radius);

bDrawBuffer *bDrawBufferCreate (int Nx, int Ny) {

  int i, j;
  bDrawColor white;
  bDrawBuffer *buffer;

  white = KapaColorByName ("white");

  ALLOCATE (buffer, bDrawBuffer, 1);
  buffer[0].Nx = Nx;
  buffer[0].Ny = Ny;

  ALLOCATE (buffer[0].pixels, bDrawColor *, Ny);
  for (i = 0; i < Ny; i++) {
    ALLOCATE (buffer[0].pixels[i], bDrawColor, Nx);
    for (j = 0; j < Nx; j++) {
      buffer[0].pixels[i][j] = white;
    }
  }
  return (buffer);
}

void bDrawBufferFree (bDrawBuffer *buffer) {

  int i;

  for (i = 0; i < buffer[0].Ny; i++) {
    free (buffer[0].pixels[i]);
  }
  free (buffer[0].pixels);
  free (buffer);
  return;
}

void bDrawSetBuffer (bDrawBuffer *buffer) {
  bBuffer = buffer;
  return;
}

void bDrawSetStyle (bDrawColor color, int lw, int lt) {
  bColor = color;
  bWeight = lw;
  bType = lt;
  return;
}

// draw a point in the current color 
void bDrawPoint (int x, int y) {

  if (x < 0) return;
  if (y < 0) return;
  if (x >= bBuffer[0].Nx) return;
  if (y >= bBuffer[0].Ny) return;
  bBuffer[0].pixels[y][x] = bColor;
  return;
}

// draw a point in the current color 
void bDrawPointf (float x, float y) {

  bDrawPoint (ROUND(x), ROUND(y));
  return;
}

void bDrawTriOpen (double x1, double y1, double x2, double y2, double x3, double y3) {

  bDrawLine (x1, y1, x2, y2);
  bDrawLine (x2, y2, x3, y3);
  bDrawLine (x3, y3, x1, y1);

  return;
}

void bDrawRectOpen (double x1, double y1, double x2, double y2) {

  int X1, Y1, X2, Y2;

  if (x1 > x2) SWAP (x1, x2);
  if (y1 > y2) SWAP (y1, y2);

  X1 = MIN (MAX (ROUND (x1), 0), bBuffer[0].Nx - 1);
  X2 = MIN (MAX (ROUND (x2), 1), bBuffer[0].Nx);

  Y1 = MIN (MAX (ROUND (y1), 0), bBuffer[0].Ny - 1);
  Y2 = MIN (MAX (ROUND (y2), 1), bBuffer[0].Ny);

  bDrawLineHorizontal (X1, X2, Y1);
  bDrawLineHorizontal (X1, X2, Y2);
  bDrawLineVertical (X1, Y1, Y2);
  bDrawLineVertical (X2, Y1, Y2);
  return;
}

void bDrawRectFill (double x1, double y1, double x2, double y2) {

  int i;
  int X1, Y1, X2, Y2;

  if (x1 > x2) SWAP (x1, x2);
  if (y1 > y2) SWAP (y1, y2);

  X1 = MIN (MAX (ROUND (x1), 0), bBuffer[0].Nx - 1);
  X2 = MIN (MAX (ROUND (x2), 1), bBuffer[0].Nx);

  Y1 = MIN (MAX (ROUND (y1), 0), bBuffer[0].Ny - 1);
  Y2 = MIN (MAX (ROUND (y2), 0), bBuffer[0].Ny);

  for (i = Y1; i < Y2; i++) {
    bDrawLineHorizontal (X1, X2, i);
  } 
  return;
}

// identify the quadrant and draw the correct line
void bDrawLine (double x1, double y1, double x2, double y2) {

  int FlipDirect, FlipCoords;
  int X1, Y1, X2, Y2, dX, dY;

  /* rather than draw the line from float positions, we find the closest
     integer end-points and draw the line between those pixels */ 

  X1 = ROUND(x1);
  Y1 = ROUND(y1);
  X2 = ROUND(x2);
  Y2 = ROUND(y2);

  dX = X2 - X1;
  dY = Y2 - Y1;

  FlipCoords = (abs(dX) < abs(dY));
  FlipDirect = FlipCoords ? (y1 > y2) : (x1 > x2);

  if (!FlipDirect && !FlipCoords) bDrawLineWeight (X1, Y1, X2, Y2, FALSE);
  if ( FlipDirect && !FlipCoords) bDrawLineWeight (X2, Y2, X1, Y1, FALSE);
  if (!FlipDirect &&  FlipCoords) bDrawLineWeight (Y1, X1, Y2, X2, TRUE);
  if ( FlipDirect &&  FlipCoords) bDrawLineWeight (Y2, X2, Y1, X1, TRUE);

  return;
}

// draw a series of lines to give the line weight
void bDrawLineWeight (int X1, int Y1, int X2, int Y2, int swapcoords) {

  int dN, dNs, dNe;

  dNs = -0.5*(bWeight - 1); 
  /* 0, 0, 0, -1, -1, -2, -2 */

  dNe = +0.5*bWeight + 1; 
  /* 1, 1, 2, 2, 2, 3, 3 */

  for (dN = dNs; dN < dNe; dN++) {
    bDrawLineBresen (X1, Y1 + dN, X2, Y2 + dN, swapcoords);
  }
  return;
}

// use the Bresenham line drawing technique
// integer-only Bresenham line-draw version which is fast
void bDrawLineBresen (int X1, int Y1, int X2, int Y2, int swapcoords) {

  int X, Y, dX, dY;
  int e, e2;
  int N, DashOn;

  dX = X2 - X1;
  dY = Y2 - Y1;

  DashOn = TRUE;

  Y = Y1;
  e = 0;
  for (X = X1, N = 0; X <= X2; X++, N++) {
    if (bType == 1) { DashOn = (N % 10) < 5; }
    if (bType == 2) { DashOn = (N % 6) < 3; }
    if (swapcoords) {
      if (DashOn) bDrawPoint (Y,X);
    } else {
      if (DashOn) bDrawPoint (X,Y);
    }
    e += dY;
    e2 = 2 * e;
    if (e2 > dX) {
      Y++;
      e -= dX;
    } 
    if (e2 < -dX) {
      Y--;
      e += dX;
    }
  }
  return;
}

void bDrawLineHorizontal (int X1, int X2, int Y) {
  
  int i;

  for (i = X1; i < X2; i++) {
    bBuffer[0].pixels[Y][i] = bColor;
  }
  return;
}

void bDrawLineVertical (int X, int Y1, int Y2) {
  
  int i;

  for (i = Y1; i < Y2; i++) {
    bBuffer[0].pixels[i][X] = bColor;
  }
  return;
}

void bDrawTriFill (double x1, double y1, double x2, double y2, double x3, double y3) {

  bDrawTriOpen (x1, y1, x2, y2, x3, y3);
  return;
}

void bDrawArc (double Xc, double Yc, double Xr, double Yr, double Ts, double Te) {

  float t, dt;
  int x, y;

  /* drawing a complete circle */
//  if ((fabs(Te - Ts) > 360.0) && (Xr == Yr)) {
//    bDrawCircle (Xc, Yc, Xr);
//    return;
//  }

  /* only draw a single loop */
  if (fabs(Te - Ts) > 360.0) {
    Te = 360.0;
    Ts = 0.0;
  }

  /* smallest angle is 1/Rmax */
  dt = MAX (fabs(Xr * sin(Ts*RAD_DEG)), fabs(Yr * cos(Ts*RAD_DEG)));
  dt = 1.0 / dt;

  for (t = Ts*RAD_DEG; t <= Te*RAD_DEG; t += dt) {
    x = Xr*cos(t) + Xc;
    y = Yr*sin(t) + Yc;

    /* we could use the value of MAX(dy/dt,dx/dt) to set dt */
    bDrawPoint (x,y);

    dt = MAX (fabs(Xr * sin(t)), fabs(Yr * cos(t)));
    dt = 1.0 / dt;
  }
  return;
}

// draw a series of circles to give line weight
void bDrawCircle (double xc, double yc, double radius) {

  int dN, dNs, dNe;

  dNs = -0.5*(bWeight - 1); 
  /* 0, 0, 0, -1, -1, -2, -2 */

  dNe = +0.5*bWeight + 1; 
  /* 1, 1, 2, 2, 2, 3, 3 */

  for (dN = dNs; dN < dNe; dN++) {
    bDrawCircleSingle (xc, yc, radius + dN);
  }
  return;
}

// draw a pure circle  
void bDrawCircleSingle (double xc, double yc, double radius) {

  int Xc, Yc, Radius;
  int x, y, d;

  Xc = ROUND(xc);
  Yc = ROUND(yc);
  Radius = ROUND(radius);

  x = 0;
  y = Radius;

  // d = 3 - 2*Radius;
  d = 5 - 4*radius;

  while (x <= y) {
    bDrawPoint (Xc+x, Yc+y);
    bDrawPoint (Xc+x, Yc-y);
    bDrawPoint (Xc-x, Yc+y);
    bDrawPoint (Xc-x, Yc-y);
    bDrawPoint (Xc+y, Yc+x);
    bDrawPoint (Xc+y, Yc-x);
    bDrawPoint (Xc-y, Yc+x);
    bDrawPoint (Xc-y, Yc-x);

    if (d < 0) {
      // d = d + 4*x + 6;
      d = d + 8*x + 4;
    } else {
      // d = d + 4*(x-y) + 10;
      d = d + 8*(x-y) + 8;
      y--;
    }
    x++;
  }
}

// draw a pure circle  
void bDrawCircleFill (double xc, double yc, double radius) {

  int Xc, Yc, Radius;
  int x, y, d;

  Xc = ROUND(xc);
  Yc = ROUND(yc);
  Radius = ROUND(radius);

  x = 0;
  y = Radius;

  // d = 3 - 2*Radius;
  d = 5 - 4*radius;

  while (x <= y) {
    bDrawLineHorizontal (Xc-x, Xc+x, Yc+y);
    bDrawLineHorizontal (Xc-x, Xc+x, Yc-y);
    bDrawLineHorizontal (Xc-y, Xc+y, Yc+x);
    bDrawLineHorizontal (Xc-y, Xc+y, Yc-x);

    if (d < 0) {
      // d = d + 4*x + 6;
      d = d + 8*x + 4;
    } else {
      // d = d + 4*(x-y) + 10;
      d = d + 8*(x-y) + 8;
      y--;
    }
    x++;
  }
}

/* 
the discriminant of inside or outside the circle is:

f(x,y) = x^2 + y^2 - r^2

- negative: (x,y) inside circle
- positive: (x,y) outside circle

given d(0) = f(x,y); find d(1) = f(x+1,y):
d(1) = (x+1)^2 + y^2 - r^2
d(1) = d(0) + 2x + 1
(use d(1) if d(0) < 0, ie inside circle)

given d(0) = f(x,y); find d(1) = f(x+1,y-1):
d(1) = (x+1)^2 + (y-1)^2 - r^2
d(1) = d(0) + 2x + 1 - 2y + 1
d(1) = d(0) + 2(x-y) + 2

* init d to f(1,r-1/2) instead of r, keeping the effective boundary 
  between two pixels (also, an inside point)

f(1,r-1/2) = 1 + (r-1/2)^2 - r^2 = 5/4 - r
f(1,r-1/4) = 1 + (r-1/4)^2 - r^2 = 1 + r/2 + 1/16 = 17/16 + r/2

f(1,r-x)   = 1 + (r-x)^2 - r^2 = 1 - 2xr + x^2 = A (3 - 2r)

1+x^2 = 3A
-2x = -2A

A = x

1 - 3x + x^2 = 0

(3 +/- sqrt(9 - 4))/2 = (3 - sqrt(5))/2

multiply all d values by 4 to get integer tests:

d'(0) = 5 - 4r
d'(in) = d' + 8x + 4
d'(out) = d' + 8(x-y) + 8

*/

/* This is the Bresenham line-drawing algorithm for 1st and 4th quandrant
   vectors with positive and negative slopes < 1. this is the sequence if we use
   float errors and tests it is easy to understand, but slower than it could be

  { 
    float e, m;
    m = dY / dX;
    Y = Y1;
    e = 0;
    for (X = X1; X <= X2; X++) {
      plot (X,Y);
      e += m;
      if (e > 0.5) {
	Y++;
	e -= 1.0;
      }
      if (e < -0.5) {
	Y--;
	e += 1.0;
      }
    }
  }
*/
