Реализации алгоритмов/Алгоритм Коэна — Сазерленда

Алгоритм Коэна — Сазерленда (англ. Cohen–Sutherland) — алгоритм отсечения отрезков, то есть алгоритм, позволяющий определить часть отрезка, которая пересекает прямоугольник.

C

Примечание для тех, кто не знаком с языком Си: переменные a, b, c — не сами точки, а указатели на точки. То есть после присвоения «с = a», c и a указывают на одну и ту же ячейку памяти, так что «c->x» становится «псевдонимом» «a->x», и «c->y» становится «псевдонимом» «a->y».

#define LEFT  1  /* двоичное 0001 */
#define RIGHT 2  /* двоичное 0010 */
#define BOT   4  /* двоичное 0100 */
#define TOP   8  /* двоичное 1000 */

/* точка */
struct point {
	double x, y;
};

/* прямоугольник */
struct rect {
	double x_min, y_min, x_max, y_max;
};

/* вычисление кода точки
   r : указатель на struct rect; p : указатель на struct point */
#define vcode(r, p) \
	((((p)->x < (r)->x_min) ? LEFT : 0)  +  /* +1 если точка левее прямоугольника */ \
	 (((p)->x > (r)->x_max) ? RIGHT : 0) +  /* +2 если точка правее прямоугольника */\
	 (((p)->y < (r)->y_min) ? BOT : 0)   +  /* +4 если точка ниже прямоугольника */  \
	 (((p)->y > (r)->y_max) ? TOP : 0))     /* +8 если точка выше прямоугольника */

/* если отрезок ab не пересекает прямоугольник r, функция возвращает -1;
   если отрезок ab пересекает прямоугольник r, функция возвращает 0 и отсекает
   те части отрезка, которые находятся вне прямоугольника */
int cohen_sutherland (const struct rect *r, struct point *a, struct point *b)
{
	int code_a, code_b, code; /* код концов отрезка */
	struct point *c; /* одна из точек */

	code_a = vcode(r, a);
	code_b = vcode(r, b);
	
	/* пока одна из точек отрезка вне прямоугольника */
	while (code_a | code_b) {
		/* если обе точки с одной стороны прямоугольника, то отрезок не пересекает прямоугольник */
		if (code_a & code_b)
			return -1;

		/* выбираем точку c с ненулевым кодом */
		if (code_a) {
			code = code_a;
			c = a;
		} else {
			code = code_b;
			c = b;
		}

		/* если c левее r, то передвигаем c на прямую x = r->x_min
		   если c правее r, то передвигаем c на прямую x = r->x_max */
		if (code & LEFT) {
			c->y += (a->y - b->y) * (r->x_min - c->x) / (a->x - b->x);
			c->x = r->x_min;
		} else if (code & RIGHT) {
			c->y += (a->y - b->y) * (r->x_max - c->x) / (a->x - b->x);
			c->x = r->x_max;
		}/* если c ниже r, то передвигаем c на прямую y = r->y_min
		    если c выше r, то передвигаем c на прямую y = r->y_max */
		else if (code & BOT) {
			c->x += (a->x - b->x) * (r->y_min - c->y) / (a->y - b->y);
			c->y = r->y_min;
		} else if (code & TOP) {
			c->x += (a->x - b->x) * (r->y_max - c->y) / (a->y - b->y);
			c->y = r->y_max;
		}

		/* обновляем код */
		if (code == code_a) {
                        /* a = c */
			code_a = vcode(r,a);
		} else {
                        /* b = c */
			code_b = vcode(r,b);
                }
	}

	/* оба кода равны 0, следовательно обе точки в прямоугольнике */
	return 0;
}

Реализация алгоритма для трёхмерной модели идентична двумерной реализации, за исключением того, что вместо четырёхразрядного кода применяется шестиразрядный (дополнительные два бита глубины).

#define BOTTOM 1   // 00 000001
#define LEFT   2   // 00 000010
#define TOP    4   // 00 000100
#define RIGHT  8   // 00 001000
#define BACK  16   // 00 010000
#define FRONT 32   // 00 100000

typedef struct Area {
	float xlt,ylt,zlt; // Координаты левого верхнего угла ближней грани
	float xrb,yrb,zrb; // Координаты правого нижнего угла дальней грани
} Area;

int GetCode ( const float point[1][4], const Area &anArea ) {
	int code = 0;
	if ( point[0][1] > anArea.ylt )     // Точка выше области отсечения
		code |= TOP;
	else if ( point[0][1] < anArea.yrb )// Точка ниже области отсечения
		code |= BOTTOM;
	if ( point[0][0] > anArea.xrb )     // Точка правее области отсечения
		code |= RIGHT;
	else if ( point[0][0] < anArea.xlt )// Точка левее области отсечения
		code |= LEFT;
	if ( point[0][2] < anArea.zlt )     // Точка перед областью отсечения
		code |= FRONT;
	else if ( point[0][2] > anArea.zrb )// Точка заобластью отсечения	
		code |= BACK;
	return code;
}

// Трёхмерное отсечение методом Коэна-Сазерленда.
// beg - координаты начала отрезка [xn,yn,zn,1]
// end - координаты конца отрезка  [xk,yk,zk,1]
// anArea - координаты отсекающей области
void CS_Clip( const float beg[1][4], const float end[1][4], const Area &anArea )
{
	// Коды концов отрезка
	int outcode0 = 0, outcode1 = 0, outcodeOut = 0;
		
	char codes[2][7] = {"000000","000000"};
	float temp[2][1][4];

	// Флаг видимости и флаг завершения отсечения
	bool accept = false, done = false;

	outcode0 = GetCode( beg, codes[0],anArea ); // Вычисляем начальные коды концов отрезка
	outcode1 = GetCode( end, codes[1],anArea );
	do
	{
		float x, y, z;
		if (!( outcode0 | outcode1 )) { // Отрезок полностью видимый
			accept = true;
			done = true;
		}
		else if ( outcode0 & outcode1 ) // Отрезок полностью невидимый
			done = true;
		else {                          // Отрезок частично видимый. Вычисление точек пересечения отрезка и области отсечения
			outcodeOut = outcode0 ? outcode0: outcode1;
			if( outcodeOut & TOP ) {
				x = beg[0][0]+(end[0][0]-beg[0][0])*(anArea.ylt-beg[0][1])/(end[0][1]-beg[0][1]);
				z = beg[0][2]+(end[0][2]-beg[0][2])*(anArea.ylt-beg[0][1])/(end[0][1]-beg[0][1]);
				y = anArea.ylt;
			}
			else if( outcodeOut & BOTTOM ) {
				x = beg[0][0]+(end[0][0]-beg[0][0])*(anArea.yrb-beg[0][1])/(end[0][1]-beg[0][1]);
				z = beg[0][2]+(end[0][2]-beg[0][2])*(anArea.yrb-beg[0][1])/(end[0][1]-beg[0][1]);
				y = anArea.yrb;
			}
			else if( outcodeOut & RIGHT ) {
				y = beg[0][1]+(end[0][1]-beg[0][1])*(anArea.xrb-beg[0][0])/(end[0][0]-beg[0][0]);
				z = beg[0][2]+(end[0][2]-beg[0][2])*(anArea.xrb-beg[0][0])/(end[0][0]-beg[0][0]);
				x = anArea.xrb;
			}
			else if( outcodeOut & LEFT ) {
				y = beg[0][1]+(end[0][1]-beg[0][1])*(anArea.xlt-beg[0][0])/(end[0][0]-beg[0][0]);
				z = beg[0][2]+(end[0][2]-beg[0][2])*(anArea.xlt-beg[0][0])/(end[0][0]-beg[0][0]);
				x = anArea.xlt;
			}
			else if( outcodeOut & FRONT ) {
				x = beg[0][0]+(end[0][0]-beg[0][0])*(anArea.zlt-beg[0][2])/(end[0][2]-beg[0][2]);
				y = beg[0][1]+(end[0][1]-beg[0][1])*(anArea.zlt-beg[0][2])/(end[0][2]-beg[0][2]);
				z = anArea.zrb;
			}
			else if( outcodeOut & BACK ) {
				x = beg[0][0]+(end[0][0]-beg[0][0])*(anArea.zrb-beg[0][2])/(end[0][2]-beg[0][2]);
				y = beg[0][1]+(end[0][1]-beg[0][1])*(anArea.zrb-beg[0][2])/(end[0][2]-beg[0][2]);
				z = anArea.zlt;
			}
			if ( outcodeOut == outcode0 ) {
				temp[0][0][0] = x;
				temp[0][0][1] = y;
				temp[0][0][2] = z;
				outcode0 = GetCode( temp[0],codes[0],anArea );
			}
			else {
				temp[1][0][0] = x;
				temp[1][0][1] = y;
				temp[1][0][2] = z;
				outcode1 = GetCode( temp[1], codes[1],anArea );
			}
		}
	}while ( !done );	
	if( accept ) { // Отрезок полностью видимый. Вывод отрезка на экран.
		LINE(temp[0],temp[1]);
	}
}

Python

Свободная реализация на Python