#ifndef lint
static char sccsid[] = "@(#)trappath.c 9.2 88/01/19 Copyright 1985 Sun Micro";
#endif

/*
 * Copyright (c) 1985 by Sun Microsystems, Inc. 
 */

/*-
	Module to make trapezons out of paths.

	trappath.c, Wed Jul 31 17:35:30 1985

 */

#include <stdio.h>
#ifdef REF
#include <sys/types.h>
#include <ref/config.h>
#endif
#include "fract.h"
#include "shape.h"
#include "qalloc.h"

#define sfrsq(a) frsq(frrsh(a, 4))


/*
 * Returns the trapezons that correspond to the shape described by the
 * circular list of arcs l.  The endpoints of arcs adjacent in the list
 * are expected to be coincident: l should describe a single closed curve. 
 */
struct shape_trapezon *
looptotrap(l, wmask)
    struct shape_arc *l;
{
    register struct shape_arc *p;
    fract       miny;
    short       globalmaxy;
    struct shape_arc *topmost;
    struct shape_dcurve *chains = 0;
    int         chaincnt = 0;
    struct shape_trapezon *finished_traps = 0;
    struct shape_trapezon *last_trap;

    /*
     * Run through all the arcs, converting curvature to sharpness and
     * subdividing arcs which contain a vertical or horizontal minima or
     * maxima
     */
    if ((p = l) == 0)
	return 0;
    {	fract gxmin, gxmax, gymin, gymax;
	gxmin = gxmax = p->A.x;
	gymin = gymax = p->A.y;
	do {
	    if (gxmin > p->A.x) gxmin = p->A.x;
	    if (gxmax < p->A.x) gxmax = p->A.x;
	    if (gymin > p->A.y) gymin = p->A.y;
	    if (gymax < p->A.y) gymax = p->A.y;
	    if (gxmin > p->B.x) gxmin = p->B.x;
	    if (gxmax < p->B.x) gxmax = p->B.x;
	    if (gymin > p->B.y) gymin = p->B.y;
	    if (gymax < p->B.y) gymax = p->B.y;
	    if (gxmin > p->C.x) gxmin = p->C.x;
	    if (gxmax < p->C.x) gxmax = p->C.x;
	    if (gymin > p->C.y) gymin = p->C.y;
	    if (gymax < p->C.y) gymax = p->C.y;
	    p = p->next;
	} while (p != l);
#define huge(fr) ((unsigned)cfloorfr(fr) + (1<<14) >= (1<<15))
	if (huge(gxmax) || huge(gxmin) || huge(gymin) || huge(gymax)) {
	    freearc(l);
	    return 0;
	}
    }
    do {
	if (p->sh < 0) {
	    register    denom = frlsh(sfrsq(p->C.x - p->B.x)
				      + sfrsq(p->C.y - p->B.y)
				      + sfrsq(p->A.x - p->B.x)
				      + sfrsq(p->A.y - p->B.y),
				      1);
	    if (denom)
		p->sh = frmul(-p->sh,
			      frsqrt(frdiv(sfrsq(p->A.x - p->C.x)
					   + sfrsq(p->A.y - p->C.y),
					   denom)));
	    else
		p->sh = 0 /* FRHUGE */;
	}
	if (p->sh) {
	    if (p->A.y != p->B.y && p->B.y != p->C.y &&
		    (p->A.y < p->B.y) != (p->B.y < p->C.y))
		shape_subdiv(p, 1);
	    if (p->A.x != p->B.x && p->B.x != p->C.x &&
		    (p->A.x < p->B.x) != (p->B.x < p->C.x))
		shape_subdiv(p, 0);
	}
	p = p->next;
    } while (p != l);
    topmost = 0;
    miny = l->A.y + 1;

#ifdef DEBUG
    check_arclist(l);
#endif

    do {			/* throw out horizontal arcs */
	register struct shape_arc *n = p->next;
	while (n->A.y == n->C.y) {
	    p->next = n->next;
	    if (n == p)	{	/* Degenerate */
		qtfree(n);
		return 0;
	    }
	    if (n == l)
		l = p;
	    qtfree(n);
	    n = p->next;
	}
	if (miny > n->A.y) {	/* Find the topmost arc */
	    miny = n->A.y;
	    topmost = n;
	}
	if (miny > n->C.y) {
	    miny = n->C.y;
	    topmost = n;
	}
	p = n;
    } while (p != l);
    p = topmost;
    if (p == 0)
	return 0;
    if (p->next == p) {
	qtfree(p);
	return 0;
    }
    globalmaxy = roundfr(p->A.y);

#ifdef DEBUG
    check_arclist(topmost);
#endif

    {
	do {
	    register    dir = p->A.y - p->C.y;	/* upwards or downwards */
	    register struct shape_arc *nchain;
	    register struct shape_dcurve *ch;
	    fract       miny = p->A.y;
	    fract       minyx = p->A.x;
	    fract       maxy = p->A.y;
	    fract       maxyx = p->A.x;
	    if (p->A.y < p->C.y) {
		miny = p->A.y;
		minyx = p->A.x;
		maxy = p->C.y;
		maxyx = p->C.x;
	    }
	    else {
		miny = p->C.y;
		minyx = p->C.x;
		maxy = p->A.y;
		maxyx = p->A.x;
	    }
	    nchain = p->next;
	    do {
		register    ndir = nchain->A.y - nchain->C.y;
		if ((ndir ^ dir) < 0)
		    break;
		if (ndir < 0) {
		    if (nchain->A.y < miny)
			if (nchain->C.y != miny)
			    break;
			else
			    miny = nchain->A.y, minyx = nchain->A.x;
		    else if (nchain->C.y > maxy)
			if (nchain->A.y != maxy)
			    break;
			else
			    maxy = nchain->C.y, maxyx = nchain->C.x;
		    else
			break;
		}
		else {
		    if (nchain->C.y < miny)
			if (nchain->A.y != miny)
			    break;
			else
			    miny = nchain->C.y, minyx = nchain->C.x;
		    else if (nchain->A.y > maxy)
			if (nchain->C.y != maxy)
			    break;
			else
			    maxy = nchain->A.y, maxyx = nchain->A.y;
		    else
			break;
		}
		nchain = nchain->next;
	    } while (nchain != topmost);
	    if (miny >= maxy) {
		p = nchain;
		continue;
	    }
	    {
		register    size = roundfr(maxy - miny) + 3;
		if (size > 4000 || size <= 0) {
		    p = nchain;
		    continue;
		}
		size = sizeof(struct shape_dcurve) + size * sizeof(short);
		QALLOC(ch, size, (struct shape_dcurve *));
	    }
	    ch->y0 = roundfr(miny);
	    ch->x0 = roundfr(minyx);
	    ch->y1 = roundfr(maxy);
	    ch->x1 = roundfr(maxyx);
	    ch->x = (short *) (((char *) ch) + sizeof(struct shape_dcurve));
	    if (ch->y1 > globalmaxy)
		globalmaxy = ch->y1;
	    ch->direction = dir < 0 ? -1 : 1;
	    ch->next = chains;
	    ch->left = ch->x0;
	    ch->right = ch->x0;
	    chains = ch;
	    chaincnt++;
	    do {
		extern      arc_subdivision_depth;
		arc_subdivision_depth = 0;
		if (dir > 0) {
		    struct fpoint T;
		    T = p->A;
		    p->A = p->C;
		    p->C = T;
		}
		if (cfloorfr(p->A.x) < ch->left)
		    ch->left = cfloorfr(p->A.x);
		if (cfloorfr(p->C.x) < ch->left)
		    ch->left = cfloorfr(p->C.x);
		if (cfloorfr(p->A.x) > ch->right)
		    ch->right = cfloorfr(p->A.x);
		if (cfloorfr(p->C.x) > ch->right)
		    ch->right = cfloorfr(p->C.x);
		arc_to_chain(p, ch);
		p = p->next;
	    } while (p != nchain);
	} while (p != topmost);
    }

#ifdef DEBUG
    check_arclist(topmost);
    check_arclist(l);
#endif

/*-    assert(chaincnt >= 2); */
    /*- printloop(topmost);
        printchain(chains); */
    if (chaincnt >= 2) {
	struct chainstate {
	    struct shape_dcurve *chain;
	    short      *x;
	    short       y0,
	                y1;
	};
	struct trappair {
	    struct shape_trapezon *inuse,
	               *waiting;
	}          *traps;
	struct chainstate *state;
	register struct chainstate *cst;
	int         statesize = sizeof(struct chainstate) * chaincnt;
	int         trapssize = (chaincnt /* / 2  */ + 1) * sizeof(struct trappair);
	QALLOC(state, statesize, (struct chainstate *));
	QALLOC(traps, trapssize, (struct trappair *));
	bzero(traps, trapssize);
	/* assert((chaincnt&1)==0); */
	{
	    register struct shape_dcurve *cp;
	    for (cp = chains, cst = state; cp; cp = cp->next, cst++) {
		cst->chain = cp;
		cst->x = cp->x + 1;
		cst->y0 = cp->y0;
		cst->y1 = cp->y1;
	    }
	    assert((int) cst <= ((int) state) + statesize);
	}
	{
	    register struct chainstate *p2 = &state[chaincnt] - 2;
	    while (p2 >= state) {
		for (cst = state; cst <= p2; cst++)
		    if (cst[0].y0 == cst[1].y0 ? cst[0].x[0] > cst[1].x[0] : cst[0].y0 > cst[1].y0) {
			struct chainstate T;
			T = cst[0];
			cst[0] = cst[1];
			cst[1] = T;
		    }
		p2--;
	    }
/*-	    assert(state[0].y0 == state[1].y0); */
	}
	{
	    struct chainstate *first = state;
	    register struct chainstate *rest = state;
	    struct chainstate *limit = &state[chaincnt];
	    register    y = first->y0;
	    while (rest < limit && rest->y0 <= y)
		rest++;
	    while (first < limit) {
		register struct chainstate *cs;
		int         ylimit;
		/*-		assert(rest - first >= 2); */
		cs = first;
		ylimit = cs++->y1;
		while (cs < rest) {
		    if (cs->y1 < ylimit)
			ylimit = cs->y1;
		    cs++;
		}
		if (cs < limit && cs->y0 < ylimit)
		    ylimit = cs->y0;
		if (rest == first + 2) {	/* Only two falls: fast case */
		    register struct shape_trapseg *tseg;
		    register short *x0p = first[0].x;
		    register short *x1p = first[1].x;
		    register struct shape_trapezon *tt;
		    if (first[0].x[0] > first[1].x[0]) {
			struct chainstate tchain;
			tchain = first[0];
			first[0] = first[1];
			first[1] = tchain;
		    }
		    if ((tt = traps[0].waiting) == 0) {
			register    size = sizeof(struct shape_trapezon) +
			sizeof(struct shape_trapseg) * (globalmaxy - y + 1);
			QALLOC(traps[0].waiting, size, (struct shape_trapezon *));
			tt = traps[0].waiting;
			if (finished_traps == 0)
			    finished_traps = tt;
			else
			    last_trap->next = tt;
			last_trap = tt;
			tt->next = 0;
			tt->refcnt = 1;
			tt->top = y;
			tt->left = tt->right = *x0p;
			tt->bottom = 0;
			tt->used = 0;
			tt->size = globalmaxy - y + 2;
			tseg = tt->segs;
		    }
		    else {
			tseg = &tt->segs[tt->used];
			if (tseg->y != y)
			    tseg++;
		    }
		    if (first[0].chain->left < tt->left)
			tt->left = first[0].chain->left;
		    if (first[0].chain->right > tt->right)
			tt->right = first[0].chain->right;
		    if (first[1].chain->left < tt->left)
			tt->left = first[1].chain->left;
		    if (first[1].chain->right > tt->right)
			tt->right = first[1].chain->right;
		    if (first[0].chain->right <= first[1].chain->left)
			while (y < ylimit) {	/* No crossing possible */
			    tseg->x1 = *x0p++;
			    tseg->x2 = *x1p++;
			    if (tseg[0].x1 != tseg[-1].x1 || tseg[0].x2 != tseg[-1].x2 || tseg == tt->segs) {
				tseg->y = y;
				tseg++;
			    }
			    y++;
			}
		    else
			while (y < ylimit) {
			    if (*x0p <= *x1p) {
				tseg->x1 = *x0p++;
				tseg->x2 = *x1p++;
			    }
			    else {
				tseg->x1 = *x1p++;
				tseg->x2 = *x0p++;
			    }
			    if (tseg[0].x1 != tseg[-1].x1 || tseg[0].x2 != tseg[-1].x2 || tseg == tt->segs) {
				tseg->y = y;
				tseg++;
			    }
			    y++;
			}
		    first[0].x = x0p;
		    first[1].x = x1p;
		    if (*x0p > *x1p) {
			struct chainstate tchain;
			tchain = first[0];
			first[0] = first[1];
			first[1] = tchain;
		    }
		    tt->used = tseg - tt->segs;
		    tt->bottom = y;
		    tseg->y = y;
		    tseg->x1 = 0;
		    tseg->x2 = 0;
		    traps[0].inuse = 0;
		    assert(tt->used < tt->size);
		}
		else if (first + 1 >= rest) {	/* only one fall, only happens
						 * in bizarre degenerate
						 * cases: flush */
		    first = rest;
		    if (rest < limit)
			y = rest->y0;
		}
		else
		 /* More than two crossings per scan line */ if (y < ylimit)
		    while (1) {
			int         winding_number = 0;
			struct chainstate *leftedge = 0;
			struct trappair *tp = traps;
			cs = first;
			while (cs < rest) {
			    winding_number += cs->chain->direction;
			    if ((winding_number & wmask) != 0) {	/* inside */
				if (leftedge == 0)
				    leftedge = cs;
			    }
			    else if (leftedge) {	/* outside */
				if (leftedge->x[-1] < cs->x[0]) {
				    register struct shape_trapezon *t;
				    assert((int) tp < ((int) traps) + trapssize - sizeof(*tp));
				    if ((t = tp->inuse) == 0)
					if ((t = tp->waiting) != 0) {
					    tp->inuse = t;
					    if (t->segs[t->used].y != y)
						t->used++;
					    {
						register    T;
						if ((T = leftedge->chain->left) < t->left)
						    t->left = T;
						if ((T = cs->chain->left) < t->left)
						    t->left = T;
						if ((T = leftedge->chain->right) > t->right)
						    t->right = T;
						if ((T = cs->chain->right) > t->right)
						    t->right = T;
					    }
					}
					else {
					    register    size = sizeof(struct shape_trapezon) +
					    sizeof(struct shape_trapseg) * (globalmaxy - y + 1);
					    QALLOC(t, size, (struct shape_trapezon *));
					    tp->waiting = t;
					    tp->inuse = t;
					    if (finished_traps == 0)
						finished_traps = t;
					    else
						last_trap->next = t;
					    last_trap = t;
					    t->next = 0;
					    t->refcnt = 1;
					    t->top = y;
					    t->bottom = 0;
					    t->used = 0;
					    t->size = globalmaxy - y + 2;
					    t->left = min(leftedge->chain->left, cs->chain->left);
					    t->right = max(leftedge->chain->right, cs->chain->right);
					}
				    {
					register struct shape_trapseg *seg;
					seg = &t->segs[t->used];
					seg->x1 = leftedge->x[-1];
					seg->x2 = cs->x[0];
					seg->y = y;
					if (t->used == 0 || seg[-1].x1 != seg[0].x1
						|| seg[-1].x2 != seg[0].x2)
					    t->used++;
					assert(t->used < t->size);
				    }
				    tp++;
				}
				leftedge = 0;
			    }
			    cs->x++;
			    cs++;
			}
			y++;
			while (tp->inuse) {
			    register struct shape_trapseg *seg = &tp->inuse->segs[tp->inuse->used];
			    seg->x1 = 0;
			    seg->x2 = 0;
			    seg->y = y;
			    tp->inuse->bottom = y;
			    tp->inuse = 0;
			}
			if (y >= ylimit)
			    break;
			for (cs = first + 1; cs < rest; cs++)
			    if (cs[0].x[0] < cs[-1].x[0]) {
				register struct chainstate *csp = cs - 1;
				struct chainstate T;
				while (csp > first && cs[0].x[0] < csp[-1].x[0])
				    csp--;
				T = *cs;
				bcopy(csp, csp + 1, ((int) cs) - ((int) csp));
				*csp = T;
			    }
		    }
		for (cs = first; cs < rest; cs++)
		    if (cs->y1 <= y) {
			if (cs != first)
			    bcopy(first, first + 1, ((int) cs) - (int) first);
			first++;
		    }
		while (rest < limit && rest->y0 <= y) {
		    if (rest->x[0] < rest[-1].x[0]) {
			struct chainstate T;
			T = rest[0];
			for (cs = rest - 1; cs >= first && cs->x[0] >= T.x[0]; cs--)
			    cs[1] = cs[0];
			cs[1] = T;
		    }
		    rest++;
		}
		{
		    register struct trappair *tp = traps;
		    while (tp->inuse) {
			register struct shape_trapseg *seg = &tp->waiting->segs[tp->inuse->used];
			seg->x1 = 0;
			seg->x2 = 0;
			seg->y = y;
			tp->inuse->bottom = y;
			tp->inuse = 0;
			tp++;
		    }
		}
	    }
	}
	QFREE(state);
	QFREE(traps);
    }
    {
	register struct shape_dcurve *ch,
	           *cn;
	for (ch = chains; ch; ch = cn) {
	    cn = ch->next;
	    QFREE(ch);
	}
    }

#ifdef DEBUG
    check_arclist(l);
#endif

    freearc(l);
    return finished_traps;
}

/*
 * Given an arc a1 subdivide it and form a second arc, linked to it, that
 * is broken at the horizontal or vertical minima or maxima, which is
 * presumed to exist 
 */
static
shape_subdiv(a1, hor)
    register struct shape_arc *a1;
    int         hor;
{
    register struct shape_arc *a2;
    fract       mu;
    int         swap = 0;
    {
	register fract a, c;
	if (hor) {		/* hor == 1 means find horizontal point */
	    a = a1->B.y - a1->A.y;
	    c = a1->C.y - a1->B.y;
	}
	else {
	    a = a1->B.x - a1->A.x;
	    c = a1->C.x - a1->B.x;
	}
	/*
	 * a and c have opposite signs (unless both are 0) Negate the
	 * negative one if any. 
	 */
	if (a < 0)
	    a = -a;
	if (c < 0)
	    c = -c;
	if (a > c) {
	    register fract t;
	    t = a, a = c, c = t, swap = 1;
	}
	if (c == 0)
	    return;
	/* now 0 <= a <= c, 0 < c */
	mu = frdiv(a, c);
	assert(mu<=fracti(1));
    }
    {
	register fract s, t;
	fract       lambda, rho, Ss, St, sqrtFz, Dx, Dy, Dz, Ex, Ey, Ez, Fx, Fy, Fz;
	if (mu < fraction(1, 256) || a1->sh == 0)
	    lambda = frmul(a1->sh, mu);
	else {
	    rho = frdiv(-(fracti(1) - mu), 2 * a1->sh);
	    lambda = rho + frsqrt(frmul(rho, rho) + mu);
	}
	t = frdiv(lambda, lambda + fracti(1));
	if (swap)
	    s = t, t = fracti(1) - s;
	else
	    s = fracti(1) - t;
	assert(t <= fracti(1));
	assert(s <= fracti(1));
	qtalloc(a2, (struct shape_arc *));
	a2->next = a1->next, a1->next = a2, a2->C = a1->C;
	St = frmul(a1->sh, t), Ss = a1->sh - St;
	Dx = frmul(s, a1->A.x) + frmul(St, a1->B.x);
	Dy = frmul(s, a1->A.y) + frmul(St, a1->B.y);
	Dz = s + St;
	Ex = frmul(Ss, a1->B.x) + frmul(t, a1->C.x);
	Ey = frmul(Ss, a1->B.y) + frmul(t, a1->C.y);
	Ez = Ss + t;
	Fx = frmul(s, Dx) + frmul(t, Ex);
	Fy = frmul(s, Dy) + frmul(t, Ey);
	Fz = frmul(s, Dz) + frmul(t, Ez);
	sqrtFz = frsqrt(Fz);
	a2->A.x = frdiv(Fx, Fz);
	a2->A.y = frdiv(Fy, Fz);
	a2->sh = frdiv(Ez, sqrtFz);
	a1->sh = frdiv(Dz, sqrtFz);
	if (hor) {
	    a2->B.x = a1->B.x == a1->C.x ? a2->C.x : frdiv(Ex, Ez);
	    a2->B.y = a2->A.y;
	    a1->B.x = a1->B.x == a1->A.x ? a1->A.x : frdiv(Dx, Dz);
	    a1->B.y = a2->A.y;
	}
	else {
	    a2->B.x = a2->A.x;
	    a2->B.y = a1->B.y == a1->C.y ? a2->C.y : frdiv(Ey, Ez);
	    a1->B.x = a2->A.x;
	    a1->B.y = a1->B.y == a1->A.y ? a1->A.y : frdiv(Dy, Dz);
	    /*-	a2->B.x = !hor? a2->A.x: a1->B.x == a1->C.x? a2->C.x: frdiv(Ex,Ez);
    	    	a2->B.y =  hor? a2->A.y: a1->B.y == a1->C.y? a2->C.y: frdiv(Ey,Ez);
    	    	a1->B.x = !hor? a2->A.x: a1->B.x == a1->A.x? a1->A.x: frdiv(Dx,Dz);
    	    	a1->B.y =  hor? a2->A.y: a1->B.y == a1->A.y? a1->A.y: frdiv(Dy,Dz); */
	}
	a1->C = a2->A;
    }
}

#ifdef DEBUG
check_arclist(l)
    register struct shape_arc *l;
{
#ifdef undef
    register struct shape_arc *p;
    int         cnt = 30;
    qalloc_check();
    p = l;
    do {
	p = p->next;
	assert(p->A.x >= 0 && p->A.y >= 0);
	assert(p->B.x >= 0 && p->B.y >= 0);
	assert(p->C.x >= 0 && p->C.y >= 0);
    } while (p != l && --cnt >= 0);
    assert(cnt >= 0);
#endif
}

#endif

compute_bbox_from_traps(sh)
    register struct shape *sh;
{
    register struct shape_trapezon *p;
    if ((p = sh->trapset) && !sh->is_rect) {
	sh->pos.x = p->left;
	sh->pos.y = p->top;
	sh->size.x = p->right - p->left;
	sh->size.y = p->bottom - p->top;
	while (p = p->next) {
	    register    T;
	    if (sh->pos.x > (T = p->left))
		sh->pos.x = T;
	    if (sh->pos.y > (T = p->top))
		sh->pos.y = T;
	    if (sh->size.x > (T = p->right - p->left))
		sh->size.x = T;
	    if (sh->size.y > (T = p->bottom - p->top))
		sh->size.y = T;
	}
    }
}