#include "trop.h"

//=======================================================================

class HalfCord

{

	int *m_array;

	int m_radius;

public:

	HalfCord(int radius)

	{

		assert(radius >= 0);

		m_radius = radius;

		m_array = new int[m_radius + 1];

		memset(m_array, 0, (m_radius + 1) * sizeof(int));

		float dCircle = 1.25f - m_radius; //  inizializza decision variable

		int y = m_radius;				  //  inizializzazione indice scanline

		int x = 0;						  //  inizializzazione indice colonna

		do {

			m_array[y] = tmax(x, m_array[y]);

			m_array[x] = y;

			if (dCircle <= 0) {

				dCircle = dCircle + 2 * x + 3;

			} else {

				y--;

				dCircle = dCircle + 2 * (x - y) + 5;

			}

			x++;

		} while (y >= x);

	}

	~HalfCord()

	{

		delete[] m_array;

	}

	inline int getCord(int x)

	{

		assert(0 <= x && x <= m_radius);

		return m_array[x];

	};

private:

	// not implemented

	HalfCord(const HalfCord &);

	HalfCord &operator=(const HalfCord &);

};

//=======================================================================

void TRop::brush(

	TRaster32P ras,

	const TPoint &aa,

	const TPoint &bb,

	int radius,

	const TPixel32 &col)

{

	TPoint a = aa;

	TPoint b = bb;

	if (a.y > b.y)

		tswap(a, b); //  a e' piu' in basso di b

	int lx = ras->getLx();

	int ly = ras->getLy();

	ras->lock();

	// ----- radius = 0

	if (radius == 0) {

		//  k = +1/-1 se il rettangolo e' inclinato positivamente (0<=m)/negativamente (m<0)

		//  (se k<0 viene fatta una riflessione sulle ascisse prima di tornare alle

		//  coordinate "di schermo")

		int k = 1;

		int dy = b.y - a.y;

		int dx = b.x - a.x;

		if (dx < 0) {

			dx = -dx;

			k = -1;

		}

		assert(dx >= 0);

		assert(dy >= 0);

		double m; //  m sara' definita solo per dx!=0)

		if (dx > 0) {

			m = dy / (double)dx;

		}

		//double length = sqrt(dx*dx + dy*dy);

		const int alpha = dy, beta = -dx;

		const int incE = alpha;

		const int incNE = alpha + beta;

		const int incN = beta;

		//  N.B. le coordinate sono relative ad un sist. di rif. con l'origine in a

		//  l'eq. della retta e' alpha * x + beta * y = 0

		int yMin = tmax(a.y, 0) - a.y;		//  clipping y + cambio  riferimento

		int yMax = tmin(b.y, ly - 1) - a.y; //  (trasporto dell'origine in a)

		if (dx > 0 && m <= 1) {

			//  midpoint algorithm

			TPoint segm;

			if (dy == 0) //  segmento orizzontale: inizializza segm

			{

				segm.x = 0;

				segm.y = yMin;

			} else //  0

			{

				segm.x = tceil((yMin - 0.5) / m);

				segm.y = yMin;

			}

			int dSegm = tfloor(alpha * (segm.x + 1) + beta * (segm.y + 0.5));

			while (segm.y <= yMax) {

				int count = 0;					  //  i trati orizzontali di segm vengono disegnati in "blocco"

				while (dSegm < 0 && segm.x <= dx) //  Est:  segm.x<=dx evita il ciclo

				{								  //  infinito quando m=0 (incE=0)

					dSegm = dSegm + incE;

					segm.x++;

					count++;

				}

				//  NordEst

				int xMin, xMax;

				if (k > 0) {

					xMin = tmax(a.x + segm.x - count, a.x, 0); //  clipping x + ritorno alle

					xMax = tmin(a.x + segm.x, b.x, lx - 1);	//  coordinate "di schermo"

				} else {

					xMin = tmax(a.x - segm.x, a.x - dx, 0);			//  clipping x + riflessione + ritorno

					xMax = tmin(a.x - segm.x + count, a.x, lx - 1); //  alle  coordinate "di schermo"

				}

				TPixel32 *p = ras->pixels(segm.y + a.y) + xMin;

				TPixel32 *q = p + (xMax - xMin);

				while (p <= q)

					*p++ = col;

				dSegm = dSegm + incNE;

				segm.x++;

				segm.y++;

			}

		} else //  m>1 oppure segmento verticale

		{

			//  midpoint algorithm

			TPoint segm;

			if (dx == 0) //  segmento verticale: inizializza segm

			{

				segm.x = 0;

				segm.y = yMin;

			} else //  m>1 :  inizializza segm

			{

				segm.x = tround(yMin / m);

				segm.y = yMin;

			}

			int dSegm = tfloor(alpha * (segm.x + 0.5) + beta * (segm.y + 1));

			while (segm.y <= yMax) {

				int xMin, xMax;

				if (k > 0) {

					xMin = tmax(a.x + segm.x, 0);	  //  clipping x + ritorno alle

					xMax = tmin(a.x + segm.x, lx - 1); //  coordinate "di schermo"

				} else {

					xMin = tmax(a.x - segm.x, 0);	  //  clipping x + riflessione + ritorno

					xMax = tmin(a.x - segm.x, lx - 1); //  alle  coordinate "di schermo"

				}

				TPixel32 *p = ras->pixels(segm.y + a.y) + xMin;

				TPixel32 *q = p + (xMax - xMin);

				while (p <= q)

					*p++ = col;

				if (dSegm <= 0) //  NordEst

				{

					dSegm = dSegm + incNE;

					segm.x++;

				} else //  Nord

				{

					dSegm = dSegm + incN;

				}

				segm.y++;

			}

		}

		ras->unlock();

		return;

	}

	HalfCord halfCord(radius);

	int x, y;

	// ----- punti iniziali coincidenti: disegna un cerchio

	if (a == b) {

		int yMin = tmax(a.y - radius, 0);	  //  clipping y

		int yMax = tmin(a.y + radius, ly - 1); //  clipping y

		for (y = yMin; y <= yMax; y++) {

			int deltay = abs(y - a.y);

			int xMin = tmax(a.x - halfCord.getCord(deltay), 0);		 //  clipping x

			int xMax = tmin(a.x + halfCord.getCord(deltay), lx - 1); //  clipping x

			TPixel32 *p = ras->pixels(y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

		}

		ras->unlock();

		return;

	}

	// -----  rettangolo orizzontale (a.y = b.y, a.x != b.x)

	if (a.y == b.y) {

		int yMin = tmax((a.y - radius), 0);		 //  clipping y

		int yMax = tmin((a.y + radius), ly - 1); //  clipping y

		int xLeft = tmin(a.x, b.x);

		int xRight = tmax(a.x, b.x);

		for (y = yMin; y <= yMax; y++) {

			int deltay = abs(y - a.y);

			int xMin = tmax(xLeft - halfCord.getCord(deltay), 0);		//  clipping x

			int xMax = tmin(xRight + halfCord.getCord(deltay), lx - 1); //  clipping x

			TPixel32 *p = ras->pixels(y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

		}

		ras->unlock();

		return;

	}

	// -----  rettangolo verticale (a.x = b.x, a.y != b.y)

	if (a.x == b.x) {

		int xMin = tmax(a.x - radius, 0);	  //  clipping x

		int xMax = tmin(a.x + radius, lx - 1); //  clipping x

		for (x = xMin; x <= xMax; x++) {

			int deltax = abs(x - a.x);

			int yMin = tmax(a.y - halfCord.getCord(deltax), 0);		 //  clipping y

			int yMax = tmin(b.y + halfCord.getCord(deltax), ly - 1); //  clipping y

			if (yMin <= yMax) {

				TPixel32 *p = ras->pixels(yMin) + x;

				TPixel32 *q = ras->pixels(yMax) + x;

				int wrap = ras->getWrap();

				while (p <= q) {

					*p = col;

					p += wrap;

				}

			}

		}

		ras->unlock();

		return;

	}

	// -----  rettangolo inclinato

	//	k = +1/-1 se il rettangolo e' inclinato positivamente/negativamente

	int k = 1;

	int dx = b.x - a.x;

	if (dx < 0) {

		dx = -dx;

		k = -1;

	}

	int dy = b.y - a.y;

	assert(dx > 0);

	assert(dy > 0);

	double length = sqrt((double)(dx * dx + dy * dy));

	const double m = dy / (double)dx;

	//punto di tangenza superiore nel sistema di riferimento del cerchio

	TPointD up(-radius * dy / length, radius * dx / length);

	//semi-ampiezza orizzontale delle "calotte" circolari

	int halfAmplCap = tfloor(-up.x);

	//  A meno di intersezioni relative tra le diverse zone:

	//  le scanline della "calotta" circolare superiore sono (b.y+cutExt,b.y+radius]

	//  le scanline del trapezoide circolare superiore sono [b.y-cutIn,b.y+cutExt]

	//  le scanline del parallelogramma sono (a.y+cutIn,b.y-cutIn)

	//  le scanline del trapezoide circolare inferiore sono [a.y-cutExt,a.y+cutIn]

	//  le scanline della "calotta" circolare inferiore sono [a.y-radius,a.y-cutExt)

	int cutExt, cutIn;

	// vertici del parallelogramma

	TPointD rightUp;

	TPointD rightDown;

	TPointD leftUp;

	TPointD leftDown;

	double mParall; //coeff. angolare parallelogramma

	//  NOTA BENE:  halfAmplCap=0 <=> (radius=0 (caso a parte) , 1)

	if (radius > 1) {

		for (cutExt = radius; cutExt >= 0 && halfCord.getCord(cutExt) <= halfAmplCap; cutExt--)

			;

		cutIn = cutExt; //  vedi else successivo

		rightUp.x = dx + halfCord.getCord(cutIn);

		rightUp.y = dy - cutIn;

		rightDown.x = halfCord.getCord(cutIn);

		rightDown.y = -cutIn;

		leftUp.x = dx - halfCord.getCord(cutIn);

		leftUp.y = dy + cutIn;

		leftDown.x = -halfCord.getCord(cutIn);

		leftDown.y = cutIn;

		mParall = dy / (double)dx;

	} else //  N.B. cutExt != cutIn solo quando radius=1

	{

		cutExt = radius; //  radius=1 => halfAmplCap=0 (non ci sono mai le "calotte" circolari)

		cutIn = 0;		 //  anche per radius=1 il limite "interno" dei trapezoidi circolari e' < radius

		rightUp.x = dx - up.x;

		rightUp.y = dy - up.y;

		rightDown.x = -up.x;

		rightDown.y = -up.y;

		leftUp.x = dx + up.x;

		leftUp.y = dy + up.y;

		leftDown.x = up.x;

		leftDown.y = up.y;

		mParall = m;

	}

	// -----  riempie "calotte" circolari

	// -----  riempie "calotta" circolare inferiore

	int yMin = tmax(a.y - radius, 0);		   //  clipping y

	int yMax = tmin(a.y - cutExt - 1, ly - 1); //  clipping y

	for (y = yMin; y <= yMax; y++) {

		int r = halfCord.getCord(a.y - y);

		int xMin = tmax(a.x - r, 0);	  //  clipping x

		int xMax = tmin(a.x + r, lx - 1); //  clipping x

		TPixel32 *p = ras->pixels(y) + xMin;

		TPixel32 *q = p + (xMax - xMin);

		while (p <= q)

			*p++ = col;

	}

	// -----  riempie "calotta" circolare superiore

	yMin = tmax(b.y + cutExt + 1, 0);  //  clipping y

	yMax = tmin(b.y + radius, ly - 1); //  clipping y

	for (y = yMin; y <= yMax; y++) {

		int r = halfCord.getCord(y - b.y);

		int xMin = tmax(b.x - r, 0);	  //  clipping x

		int xMax = tmin(b.x + r, lx - 1); //  clipping x

		TPixel32 *p = ras->pixels(y) + xMin;

		TPixel32 *q = p + (xMax - xMin);

		while (p <= q)

			*p++ = col;

	}

	// -----  riempie trapezoidi

	// (se k<0 viene fatta una riflessione sulle ascisse prima di tornare alle

	// coordinate "di schermo")

	//  limite destro assoluto delle scanline trapezoide:

	int xSegmMax = tround(dx - up.x); //  coordinata x del punto di tangenza inferiore sul cerchio superiore

	//  limite sinistro assoluto delle scanline:

	int xSegmMin = tround(up.x); //  coordinata x del punto di tangenza superiore sul cerchio inferiore

	// -----  riempie trapezoide inferiore

	// N.B. le coordinate sono relative ad un sist. di rif. con l'origine sul centro

	// del cerchio inferiore

	yMin = tmax(a.y - cutExt, 0) - a.y;						 //  clipping y

	yMax = tmin(a.y + cutIn, b.y - cutIn - 1, ly - 1) - a.y; //  clipping y

	// l'eq. della retta e' alpha * x + beta * y + gammaRight = 0

	const int alpha = dy, beta = -dx;

	const double gammaRight = rightDown.y * dx - rightDown.x * dy;

	const int incE = alpha;

	const int incNE = alpha + beta;

	const int incN = beta;

	if (m <= 1) {

		//  midpoint algorithm; le scanline vengono disegnate solo

		//  sul NordEst. L'ultima scanline non viene disegnata

		TPoint segmRight(tceil((yMin + 0.5 - rightDown.y) / mParall + rightDown.x) - 1, yMin);

		int dSegmRight = tfloor(alpha * (segmRight.x + 1) + beta * (segmRight.y + 0.5) + gammaRight);

		while (segmRight.y <= yMax) {

			if (dSegmRight < 0) //  Est

			{

				dSegmRight = dSegmRight + incE;

				segmRight.x++;

			} else //  NordEst

			{

				int xMin, xMax;

				if (k > 0) {

					xMin = tmax(a.x - halfCord.getCord(abs(segmRight.y)), 0); //  clipping x

					xMax = tmin(a.x + tmin(segmRight.x, xSegmMax), lx - 1);   //  clipping x

				} else {

					xMin = tmax(a.x - tmin(segmRight.x, xSegmMax), 0);			   //  clipping x + ritorno alle

					xMax = tmin(a.x + halfCord.getCord(abs(segmRight.y)), lx - 1); //   coordinate "di schermo"

				}

				TPixel32 *p = ras->pixels(segmRight.y + a.y) + xMin;

				TPixel32 *q = p + (xMax - xMin);

				while (p <= q)

					*p++ = col;

				dSegmRight = dSegmRight + incNE;

				segmRight.x++;

				segmRight.y++;

			}

		}

	} else //  m>1

	{

		//  midpoint algorithm; le scanline vengono disegnate sempre

		TPoint segmRight(tround((yMin - rightDown.y) / mParall + rightDown.x), yMin);

		int dSegmRight = tfloor(alpha * (segmRight.x + 0.5) + beta * (segmRight.y + 1) + gammaRight);

		while (segmRight.y <= yMax) {

			int xMin, xMax;

			if (k > 0) {

				xMin = tmax(a.x - halfCord.getCord(abs(segmRight.y)), 0); //  clipping x

				xMax = tmin(a.x + segmRight.x, lx - 1);					  //  clipping x

			} else {

				xMin = tmax(a.x - segmRight.x, 0);							   //  clipping x + ritorno alle coordinate

				xMax = tmin(a.x + halfCord.getCord(abs(segmRight.y)), lx - 1); //  "di schermo"

			}

			TPixel32 *p = ras->pixels(segmRight.y + a.y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

			if (dSegmRight <= 0) //  NordEst

			{

				dSegmRight = dSegmRight + incNE;

				segmRight.x++;

			} else //  Nord

			{

				dSegmRight = dSegmRight + incN;

			}

			segmRight.y++;

		}

	}

	// -----  riempie trapezoide superiore

	//  N.B. le coordinate sono relative ad un sist. di rif. con l'origine sul centro

	//  del cerchio superiore

	yMin = tmax(b.y - cutIn, a.y + cutIn + 1, 0) - b.y; //  clipping y

	yMax = tmin(b.y + cutExt, ly - 1) - b.y;			//  clipping y

	//   l'eq. della retta e' alpha * x + beta * y + gammaLeft = 0

	const double gammaLeft = leftDown.y * dx - leftDown.x * dy;

	if (m <= 1) {

		//  midpoint algorithm; le scanline vengono disegnate solo

		//  sul NordEst. L'ultima scanline non viene disegnata

		TPoint segmLeft(tceil((yMin - 0.5 - leftDown.y) / mParall + leftDown.x), yMin);

		int dSegmLeft = tfloor(alpha * (segmLeft.x + 1) + beta * (segmLeft.y + 0.5) + gammaLeft);

		while (segmLeft.y <= yMax) {

			int xMin, xMax;

			if (k > 0) {

				xMin = tmax(b.x + tmax(segmLeft.x, xSegmMin - dx), 0);		  //  clipping x

				xMax = tmin(b.x + halfCord.getCord(abs(segmLeft.y)), lx - 1); //  clipping x

			} else {

				xMin = tmax(b.x - halfCord.getCord(abs(segmLeft.y)), 0);	//  clipping x + ritorno alle

				xMax = tmin(b.x - tmax(segmLeft.x, xSegmMin - dx), lx - 1); //   coordinate "di schermo"

			}

			TPixel32 *p = ras->pixels(segmLeft.y + b.y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

			while (dSegmLeft < 0) {

				dSegmLeft = dSegmLeft + incE;

				segmLeft.x++;

			}

			dSegmLeft = dSegmLeft + incNE;

			segmLeft.x++;

			segmLeft.y++;

		}

	} else //  m>1

	{

		//  midpoint algorithm; le scanline vengono disegnate sempre

		TPoint segmLeft(tround((yMin - leftDown.y) / mParall + leftDown.x), yMin);

		int dSegmLeft = tfloor(alpha * (segmLeft.x + 0.5) + beta * (segmLeft.y + 1) + gammaLeft);

		while (segmLeft.y <= yMax) {

			int xMin, xMax;

			if (k > 0) {

				xMin = tmax(b.x + segmLeft.x, 0);							  //  clipping x

				xMax = tmin(b.x + halfCord.getCord(abs(segmLeft.y)), lx - 1); //  clipping x

			} else {

				xMin = tmax(b.x - halfCord.getCord(abs(segmLeft.y)), 0); //  clipping x + ritorno alle

				xMax = tmin(b.x - segmLeft.x, lx - 1);					 //   coordinate "di schermo"

			}

			TPixel32 *p = ras->pixels(segmLeft.y + b.y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

			if (dSegmLeft <= 0) //  NordEst

			{

				dSegmLeft = dSegmLeft + incNE;

				segmLeft.x++;

			} else //  Nord

			{

				dSegmLeft = dSegmLeft + incN;

			}

			segmLeft.y++;

		}

	}

	// -----  parallelogramma (in alternativa a "parallelogrammoide circolare")

	// N.B. le coordinate sono relative ad un sist. di rif. con l'origine sul centro

	// del cerchio inferiore

	//  retta destra di equaz.   alpha * x + beta * y + gammaRight = 0

	//  retta sinistra di equaz. alpha * x + beta * y + gammaLeft = 0

	yMin = tmax(a.y + cutIn + 1, 0) - a.y;		//clipping y

	yMax = tmin(b.y - cutIn - 1, ly - 1) - a.y; //clipping y

	if (m <= 1) {

		//  midpoint algorithm; le scanline vengono disegnate solo

		//  sul NordEst. L'ultima scanline non viene disegnata

		TPoint segmRight(tceil((yMin + 0.5 - rightDown.y) / mParall + rightDown.x) - 1, yMin);

		TPoint segmLeft = TPoint(tceil((yMin - 0.5 - leftDown.y) / mParall + leftDown.x), yMin);

		int dSegmRight = tfloor(alpha * (segmRight.x + 1) + beta * (segmRight.y + 0.5) + gammaRight);

		int dSegmLeft = tfloor(alpha * (segmLeft.x + 1) + beta * (segmLeft.y + 0.5) + gammaLeft);

		while (segmRight.y <= yMax) {

			if (dSegmRight < 0) //  segmRight a Est

			{

				dSegmRight = dSegmRight + incE;

				segmRight.x++;

			} else //  segmRight a NordEst

			{

				int xMin, xMax;

				if (k > 0) {

					xMin = tmax(a.x + tmax(segmLeft.x, xSegmMin), 0);		//  clipping x

					xMax = tmin(a.x + tmin(segmRight.x, xSegmMax), lx - 1); //  clipping x

				} else {

					xMin = tmax(a.x - tmin(segmRight.x, xSegmMax), 0);	 //  clipping x + ritorno alle

					xMax = tmin(a.x - tmax(segmLeft.x, xSegmMin), lx - 1); //   coordinate "di schermo"

				}

				TPixel32 *p = ras->pixels(segmRight.y + a.y) + xMin;

				TPixel32 *q = p + (xMax - xMin);

				while (p <= q)

					*p++ = col;

				dSegmRight = dSegmRight + incNE;

				segmRight.x++;

				segmRight.y++;

				while (dSegmLeft < 0) //  segmLeft a Est

				{

					dSegmLeft = dSegmLeft + incE;

					segmLeft.x++;

				}

				//  segmLeft a NordEst

				dSegmLeft = dSegmLeft + incNE;

				segmLeft.x++;

				segmLeft.y++;

			}

		}

	} else //  m>1

	{

		//  midpoint algorithm; le scanline vengono disegnate sempre

		TPoint segmRight(tround((yMin - rightDown.y) / mParall + rightDown.x), yMin);

		TPoint segmLeft(tround((yMin - leftDown.y) / mParall + leftDown.x), yMin);

		int dSegmRight = tfloor(alpha * (segmRight.x + 0.5) + beta * (segmRight.y + 1) + gammaRight);

		int dSegmLeft = tfloor(alpha * (segmLeft.x + 0.5) + beta * (segmLeft.y + 1) + gammaLeft);

		while (segmRight.y <= yMax) {

			int xMin, xMax;

			if (k > 0) {

				xMin = tmax(a.x + segmLeft.x, 0);		//  clipping x

				xMax = tmin(a.x + segmRight.x, lx - 1); //  clipping x

			} else {

				xMin = tmax(a.x - segmRight.x, 0);	 //  clipping x + ritorno alle

				xMax = tmin(a.x - segmLeft.x, lx - 1); //   coordinate "di schermo"

			}

			TPixel32 *p = ras->pixels(segmRight.y + a.y) + xMin;

			TPixel32 *q = p + (xMax - xMin);

			while (p <= q)

				*p++ = col;

			if (dSegmRight <= 0) //  segmRight a NordEst

			{

				dSegmRight = dSegmRight + incNE;

				segmRight.x++;

			} else //  segmRight a Nord

			{

				dSegmRight = dSegmRight + incN;

			}

			segmRight.y++;

			if (dSegmLeft <= 0) //  segmLeft a NordEst

			{

				dSegmLeft = dSegmLeft + incNE;

				segmLeft.x++;

			} else //  segmLeft a Nord

			{

				dSegmLeft = dSegmLeft + incN;

			}

		}

	}

	// ----  parallelogrammoide circolare (in alternativa a parallelogramma)

	// N.B. coordinate di schermo (riflessione per k<0 )

	yMin = tmax(b.y - cutIn, 0);

	yMax = tmin(a.y + cutIn, ly - 1);

	for (y = yMin; y <= yMax; y++) {

		int xMin, xMax;

		if (k > 0) {

			xMin = tmax(a.x - halfCord.getCord(abs(y - a.y)), 0);	  //  clipping x

			xMax = tmin(b.x + halfCord.getCord(abs(b.y - y)), lx - 1); //  clipping x

		} else {

			xMin = tmax(b.x - halfCord.getCord(abs(b.y - y)), 0);	  //  clipping x + ritorno alle

			xMax = tmin(a.x + halfCord.getCord(abs(y - a.y)), lx - 1); //   coordinate "di schermo"

		}

		TPixel32 *p = ras->pixels(y) + xMin;

		TPixel32 *q = p + (xMax - xMin);

		while (p <= q)

			*p++ = col;

	}

	ras->unlock();

}