1 324891a5 2006-06-25 devnull #include <u.h>
  2 324891a5 2006-06-25 devnull #include <libc.h>
  3 324891a5 2006-06-25 devnull #include <draw.h>
  4 324891a5 2006-06-25 devnull #include <memdraw.h>
  5 324891a5 2006-06-25 devnull 
  6 324891a5 2006-06-25 devnull /*
  7 324891a5 2006-06-25 devnull  * ellipse(dst, c, a, b, t, src, sp)
  8 324891a5 2006-06-25 devnull  *   draws an ellipse centered at c with semiaxes a,b>=0
  9 324891a5 2006-06-25 devnull  *   and semithickness t>=0, or filled if t<0.  point sp
 10 324891a5 2006-06-25 devnull  *   in src maps to c in dst
 11 324891a5 2006-06-25 devnull  *
 12 324891a5 2006-06-25 devnull  *   very thick skinny ellipses are brushed with circles (slow)
 13 324891a5 2006-06-25 devnull  *   others are approximated by filling between 2 ellipses
 14 324891a5 2006-06-25 devnull  *   criterion for very thick when b<a: t/b > 0.5*x/(1-x)
 15 324891a5 2006-06-25 devnull  *   where x = b/a
 16 324891a5 2006-06-25 devnull  */
 17 324891a5 2006-06-25 devnull 
 18 324891a5 2006-06-25 devnull typedef struct Param	Param;
 19 324891a5 2006-06-25 devnull typedef struct State	State;
 20 324891a5 2006-06-25 devnull 
 21 324891a5 2006-06-25 devnull static	void	bellipse(int, State*, Param*);
 22 324891a5 2006-06-25 devnull static	void	erect(int, int, int, int, Param*);
 23 324891a5 2006-06-25 devnull static	void	eline(int, int, int, int, Param*);
 24 324891a5 2006-06-25 devnull 
 25 324891a5 2006-06-25 devnull struct Param {
 26 324891a5 2006-06-25 devnull 	Memimage	*dst;
 27 324891a5 2006-06-25 devnull 	Memimage	*src;
 28 324891a5 2006-06-25 devnull 	Point			c;
 29 324891a5 2006-06-25 devnull 	int			t;
 30 324891a5 2006-06-25 devnull 	Point			sp;
 31 324891a5 2006-06-25 devnull 	Memimage	*disc;
 32 324891a5 2006-06-25 devnull 	int			op;
 33 324891a5 2006-06-25 devnull };
 34 324891a5 2006-06-25 devnull 
 35 324891a5 2006-06-25 devnull /*
 36 324891a5 2006-06-25 devnull  * denote residual error by e(x,y) = b^2*x^2 + a^2*y^2 - a^2*b^2
 37 324891a5 2006-06-25 devnull  * e(x,y) = 0 on ellipse, e(x,y) < 0 inside, e(x,y) > 0 outside
 38 324891a5 2006-06-25 devnull  */
 39 324891a5 2006-06-25 devnull 
 40 324891a5 2006-06-25 devnull struct State {
 41 324891a5 2006-06-25 devnull 	int	a;
 42 324891a5 2006-06-25 devnull 	int	x;
 43 324891a5 2006-06-25 devnull 	vlong	a2;	/* a^2 */
 44 324891a5 2006-06-25 devnull 	vlong	b2;	/* b^2 */
 45 324891a5 2006-06-25 devnull 	vlong	b2x;	/* b^2 * x */
 46 324891a5 2006-06-25 devnull 	vlong	a2y;	/* a^2 * y */
 47 324891a5 2006-06-25 devnull 	vlong	c1;
 48 324891a5 2006-06-25 devnull 	vlong	c2;	/* test criteria */
 49 324891a5 2006-06-25 devnull 	vlong	ee;	/* ee = e(x+1/2,y-1/2) - (a^2+b^2)/4 */
 50 324891a5 2006-06-25 devnull 	vlong	dxe;
 51 324891a5 2006-06-25 devnull 	vlong	dye;
 52 324891a5 2006-06-25 devnull 	vlong	d2xe;
 53 324891a5 2006-06-25 devnull 	vlong	d2ye;
 54 324891a5 2006-06-25 devnull };
 55 324891a5 2006-06-25 devnull 
 56 324891a5 2006-06-25 devnull static
 57 324891a5 2006-06-25 devnull State*
 58 324891a5 2006-06-25 devnull newstate(State *s, int a, int b)
 59 324891a5 2006-06-25 devnull {
 60 324891a5 2006-06-25 devnull 	s->x = 0;
 61 324891a5 2006-06-25 devnull 	s->a = a;
 62 324891a5 2006-06-25 devnull 	s->a2 = (vlong)(a*a);
 63 324891a5 2006-06-25 devnull 	s->b2 = (vlong)(b*b);
 64 324891a5 2006-06-25 devnull 	s->b2x = (vlong)0;
 65 324891a5 2006-06-25 devnull 	s->a2y = s->a2*(vlong)b;
 66 324891a5 2006-06-25 devnull 	s->c1 = -((s->a2>>2) + (vlong)(a&1) + s->b2);
 67 324891a5 2006-06-25 devnull 	s->c2 = -((s->b2>>2) + (vlong)(b&1));
 68 324891a5 2006-06-25 devnull 	s->ee = -s->a2y;
 69 324891a5 2006-06-25 devnull 	s->dxe = (vlong)0;
 70 324891a5 2006-06-25 devnull 	s->dye = s->ee<<1;
 71 324891a5 2006-06-25 devnull 	s->d2xe = s->b2<<1;
 72 324891a5 2006-06-25 devnull 	s->d2ye = s->a2<<1;
 73 324891a5 2006-06-25 devnull 	return s;
 74 324891a5 2006-06-25 devnull }
 75 324891a5 2006-06-25 devnull 
 76 324891a5 2006-06-25 devnull /*
 77 324891a5 2006-06-25 devnull  * return x coord of rightmost pixel on next scan line
 78 324891a5 2006-06-25 devnull  */
 79 324891a5 2006-06-25 devnull static
 80 324891a5 2006-06-25 devnull int
 81 324891a5 2006-06-25 devnull step(State *s)
 82 324891a5 2006-06-25 devnull {
 83 324891a5 2006-06-25 devnull 	while(s->x < s->a) {
 84 324891a5 2006-06-25 devnull 		if(s->ee+s->b2x <= s->c1 ||	/* e(x+1,y-1/2) <= 0 */
 85 324891a5 2006-06-25 devnull 		   s->ee+s->a2y <= s->c2) {	/* e(x+1/2,y) <= 0 (rare) */
 86 fa325e9b 2020-01-10 cross 			s->dxe += s->d2xe;
 87 fa325e9b 2020-01-10 cross 			s->ee += s->dxe;
 88 324891a5 2006-06-25 devnull 			s->b2x += s->b2;
 89 fa325e9b 2020-01-10 cross 			s->x++;
 90 324891a5 2006-06-25 devnull 			continue;
 91 324891a5 2006-06-25 devnull 		}
 92 fa325e9b 2020-01-10 cross 		s->dye += s->d2ye;
 93 fa325e9b 2020-01-10 cross 		s->ee += s->dye;
 94 324891a5 2006-06-25 devnull 		s->a2y -= s->a2;
 95 324891a5 2006-06-25 devnull 		if(s->ee-s->a2y <= s->c2) {	/* e(x+1/2,y-1) <= 0 */
 96 fa325e9b 2020-01-10 cross 			s->dxe += s->d2xe;
 97 fa325e9b 2020-01-10 cross 			s->ee += s->dxe;
 98 324891a5 2006-06-25 devnull 			s->b2x += s->b2;
 99 324891a5 2006-06-25 devnull 			return s->x++;
100 324891a5 2006-06-25 devnull 		}
101 324891a5 2006-06-25 devnull 		break;
102 324891a5 2006-06-25 devnull 	}
103 fa325e9b 2020-01-10 cross 	return s->x;
104 324891a5 2006-06-25 devnull }
105 324891a5 2006-06-25 devnull 
106 324891a5 2006-06-25 devnull void
107 324891a5 2006-06-25 devnull memellipse(Memimage *dst, Point c, int a, int b, int t, Memimage *src, Point sp, int op)
108 324891a5 2006-06-25 devnull {
109 324891a5 2006-06-25 devnull 	State in, out;
110 324891a5 2006-06-25 devnull 	int y, inb, inx, outx, u;
111 324891a5 2006-06-25 devnull 	Param p;
112 324891a5 2006-06-25 devnull 
113 324891a5 2006-06-25 devnull 	if(a < 0)
114 324891a5 2006-06-25 devnull 		a = -a;
115 324891a5 2006-06-25 devnull 	if(b < 0)
116 324891a5 2006-06-25 devnull 		b = -b;
117 324891a5 2006-06-25 devnull 	p.dst = dst;
118 324891a5 2006-06-25 devnull 	p.src = src;
119 324891a5 2006-06-25 devnull 	p.c = c;
120 324891a5 2006-06-25 devnull 	p.t = t;
121 324891a5 2006-06-25 devnull 	p.sp = subpt(sp, c);
122 324891a5 2006-06-25 devnull 	p.disc = nil;
123 324891a5 2006-06-25 devnull 	p.op = op;
124 324891a5 2006-06-25 devnull 
125 324891a5 2006-06-25 devnull 	u = (t<<1)*(a-b);
126 324891a5 2006-06-25 devnull 	if(b<a && u>b*b || a<b && -u>a*a) {
127 324891a5 2006-06-25 devnull /*	if(b<a&&(t<<1)>b*b/a || a<b&&(t<<1)>a*a/b)	# very thick */
128 324891a5 2006-06-25 devnull 		bellipse(b, newstate(&in, a, b), &p);
129 324891a5 2006-06-25 devnull 		return;
130 324891a5 2006-06-25 devnull 	}
131 324891a5 2006-06-25 devnull 
132 324891a5 2006-06-25 devnull 	if(t < 0) {
133 324891a5 2006-06-25 devnull 		inb = -1;
134 324891a5 2006-06-25 devnull 		newstate(&out, a, y = b);
135 fa325e9b 2020-01-10 cross 	} else {
136 324891a5 2006-06-25 devnull 		inb = b - t;
137 324891a5 2006-06-25 devnull 		newstate(&out, a+t, y = b+t);
138 324891a5 2006-06-25 devnull 	}
139 324891a5 2006-06-25 devnull 	if(t > 0)
140 324891a5 2006-06-25 devnull 		newstate(&in, a-t, inb);
141 324891a5 2006-06-25 devnull 	inx = 0;
142 324891a5 2006-06-25 devnull 	for( ; y>=0; y--) {
143 324891a5 2006-06-25 devnull 		outx = step(&out);
144 324891a5 2006-06-25 devnull 		if(y > inb) {
145 324891a5 2006-06-25 devnull 			erect(-outx, y, outx, y, &p);
146 324891a5 2006-06-25 devnull 			if(y != 0)
147 324891a5 2006-06-25 devnull 				erect(-outx, -y, outx, -y, &p);
148 324891a5 2006-06-25 devnull 			continue;
149 324891a5 2006-06-25 devnull 		}
150 324891a5 2006-06-25 devnull 		if(t > 0) {
151 324891a5 2006-06-25 devnull 			inx = step(&in);
152 324891a5 2006-06-25 devnull 			if(y == inb)
153 324891a5 2006-06-25 devnull 				inx = 0;
154 324891a5 2006-06-25 devnull 		} else if(inx > outx)
155 324891a5 2006-06-25 devnull 			inx = outx;
156 324891a5 2006-06-25 devnull 		erect(inx, y, outx, y, &p);
157 324891a5 2006-06-25 devnull 		if(y != 0)
158 324891a5 2006-06-25 devnull 			erect(inx, -y, outx, -y, &p);
159 324891a5 2006-06-25 devnull 		erect(-outx, y, -inx, y, &p);
160 324891a5 2006-06-25 devnull 		if(y != 0)
161 324891a5 2006-06-25 devnull 			erect(-outx, -y, -inx, -y, &p);
162 324891a5 2006-06-25 devnull 		inx = outx + 1;
163 324891a5 2006-06-25 devnull 	}
164 324891a5 2006-06-25 devnull }
165 324891a5 2006-06-25 devnull 
166 324891a5 2006-06-25 devnull static Point p00 = {0, 0};
167 324891a5 2006-06-25 devnull 
168 324891a5 2006-06-25 devnull /*
169 324891a5 2006-06-25 devnull  * a brushed ellipse
170 324891a5 2006-06-25 devnull  */
171 324891a5 2006-06-25 devnull static
172 324891a5 2006-06-25 devnull void
173 324891a5 2006-06-25 devnull bellipse(int y, State *s, Param *p)
174 324891a5 2006-06-25 devnull {
175 324891a5 2006-06-25 devnull 	int t, ox, oy, x, nx;
176 324891a5 2006-06-25 devnull 
177 324891a5 2006-06-25 devnull 	t = p->t;
178 324891a5 2006-06-25 devnull 	p->disc = allocmemimage(Rect(-t,-t,t+1,t+1), GREY1);
179 324891a5 2006-06-25 devnull 	if(p->disc == nil)
180 324891a5 2006-06-25 devnull 		return;
181 324891a5 2006-06-25 devnull 	memfillcolor(p->disc, DTransparent);
182 324891a5 2006-06-25 devnull 	memellipse(p->disc, p00, t, t, -1, memopaque, p00, p->op);
183 324891a5 2006-06-25 devnull 	oy = y;
184 324891a5 2006-06-25 devnull 	ox = 0;
185 324891a5 2006-06-25 devnull 	nx = x = step(s);
186 324891a5 2006-06-25 devnull 	do {
187 324891a5 2006-06-25 devnull 		while(nx==x && y-->0)
188 324891a5 2006-06-25 devnull 			nx = step(s);
189 324891a5 2006-06-25 devnull 		y++;
190 324891a5 2006-06-25 devnull 		eline(-x,-oy,-ox, -y, p);
191 324891a5 2006-06-25 devnull 		eline(ox,-oy,  x, -y, p);
192 324891a5 2006-06-25 devnull 		eline(-x,  y,-ox, oy, p);
193 324891a5 2006-06-25 devnull 		eline(ox,  y,  x, oy, p);
194 324891a5 2006-06-25 devnull 		ox = x+1;
195 324891a5 2006-06-25 devnull 		x = nx;
196 324891a5 2006-06-25 devnull 		y--;
197 324891a5 2006-06-25 devnull 		oy = y;
198 324891a5 2006-06-25 devnull 	} while(oy > 0);
199 324891a5 2006-06-25 devnull }
200 324891a5 2006-06-25 devnull 
201 324891a5 2006-06-25 devnull /*
202 324891a5 2006-06-25 devnull  * a rectangle with closed (not half-open) coordinates expressed
203 324891a5 2006-06-25 devnull  * relative to the center of the ellipse
204 324891a5 2006-06-25 devnull  */
205 324891a5 2006-06-25 devnull static
206 324891a5 2006-06-25 devnull void
207 324891a5 2006-06-25 devnull erect(int x0, int y0, int x1, int y1, Param *p)
208 324891a5 2006-06-25 devnull {
209 324891a5 2006-06-25 devnull 	Rectangle r;
210 324891a5 2006-06-25 devnull 
211 324891a5 2006-06-25 devnull /*	print("R %d,%d %d,%d\n", x0, y0, x1, y1); */
212 324891a5 2006-06-25 devnull 	r = Rect(p->c.x+x0, p->c.y+y0, p->c.x+x1+1, p->c.y+y1+1);
213 324891a5 2006-06-25 devnull 	memdraw(p->dst, r, p->src, addpt(p->sp, r.min), memopaque, p00, p->op);
214 324891a5 2006-06-25 devnull }
215 324891a5 2006-06-25 devnull 
216 324891a5 2006-06-25 devnull /*
217 324891a5 2006-06-25 devnull  * a brushed point similarly specified
218 324891a5 2006-06-25 devnull  */
219 324891a5 2006-06-25 devnull static
220 324891a5 2006-06-25 devnull void
221 324891a5 2006-06-25 devnull epoint(int x, int y, Param *p)
222 324891a5 2006-06-25 devnull {
223 324891a5 2006-06-25 devnull 	Point p0;
224 324891a5 2006-06-25 devnull 	Rectangle r;
225 324891a5 2006-06-25 devnull 
226 324891a5 2006-06-25 devnull /*	print("P%d %d,%d\n", p->t, x, y);	*/
227 324891a5 2006-06-25 devnull 	p0 = Pt(p->c.x+x, p->c.y+y);
228 324891a5 2006-06-25 devnull 	r = Rpt(addpt(p0, p->disc->r.min), addpt(p0, p->disc->r.max));
229 324891a5 2006-06-25 devnull 	memdraw(p->dst, r, p->src, addpt(p->sp, r.min), p->disc, p->disc->r.min, p->op);
230 324891a5 2006-06-25 devnull }
231 324891a5 2006-06-25 devnull 
232 fa325e9b 2020-01-10 cross /*
233 324891a5 2006-06-25 devnull  * a brushed horizontal or vertical line similarly specified
234 324891a5 2006-06-25 devnull  */
235 324891a5 2006-06-25 devnull static
236 324891a5 2006-06-25 devnull void
237 324891a5 2006-06-25 devnull eline(int x0, int y0, int x1, int y1, Param *p)
238 324891a5 2006-06-25 devnull {
239 324891a5 2006-06-25 devnull /*	print("L%d %d,%d %d,%d\n", p->t, x0, y0, x1, y1); */
240 324891a5 2006-06-25 devnull 	if(x1 > x0+1)
241 324891a5 2006-06-25 devnull 		erect(x0+1, y0-p->t, x1-1, y1+p->t, p);
242 324891a5 2006-06-25 devnull 	else if(y1 > y0+1)
243 324891a5 2006-06-25 devnull 		erect(x0-p->t, y0+1, x1+p->t, y1-1, p);
244 324891a5 2006-06-25 devnull 	epoint(x0, y0, p);
245 324891a5 2006-06-25 devnull 	if(x1-x0 || y1-y0)
246 324891a5 2006-06-25 devnull 		epoint(x1, y1, p);
247 324891a5 2006-06-25 devnull }