2123c29f0fe6b9b6a63d141b3cb52efaa685e057
[ctsim.git] / libctsupport / clip.cpp
1 /*****************************************************************************
2 ** FILE IDENTIFICATION
3 **
4 **   Name:         clip.c
5 **   Purpose:      Line clipping routines
6 **   Programmer:   Kevin Rosenberg
7 **   Date Started: 1984
8 **
9 ** OVERVIEW
10 **      Routines to clip lines against objects
11 **      All routines get the endpoints of the line, and
12 **      the SNARK size of the object (u,v)
13 **
14 **  This is part of the CTSim program
15 **  Copyright (C) 1983-2000 Kevin Rosenberg
16 **
17 **  $Id: clip.cpp,v 1.5 2000/12/06 01:46:43 kevin Exp $
18 **
19 **  This program is free software; you can redistribute it and/or modify
20 **  it under the terms of the GNU General Public License (version 2) as
21 **  published by the Free Software Foundation.
22 **
23 **  This program is distributed in the hope that it will be useful,
24 **  but WITHOUT ANY WARRANTY; without even the implied warranty of
25 **  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
26 **  GNU General Public License for more details.
27 **
28 **  You should have received a copy of the GNU General Public License
29 **  along with this program; if not, write to the Free Software
30 **  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
31 ******************************************************************************/
32
33 #include "ctsupport.h"
34
35
36 /* NAME
37  *    clip_segment              Clip against a segment of a circle
38  *
39  * SYNOPSIS
40  *    clip_segment (x1, y1, x2, y2, u, v)
41  *    double& x1,*y1,*x2,*y2    Endpoints of line
42  *    double u,v                Dimensions of segment
43  */
44
45 int 
46 clip_segment (double& x1, double& y1, double& x2, double& y2, const double u, const double v)
47 {
48   double xc1 = x1 * u;
49   double yc1 = y1 * v;
50   double xc2 = x2 * u;
51   double yc2 = y2 * v;
52
53   if (yc1 > 0. && yc2 > 0.)     // reject lines above y-axis 
54     return false;
55
56   double radius = sqrt (u * u + v * v);
57
58   if (clip_circle (xc1, yc1, xc2, yc2, 0.0, v, radius, 0.0, 0.0) == false)
59     return false;
60
61   if (yc1 > 0. && yc2 > 0.)     // trivial reject above y-axis 
62     return false;
63
64   // clip above x-axis 
65   if (yc1 > 0.) {
66     xc1 = xc1 + (xc2-xc1)*(0.0-yc1)/(yc2-yc1);
67     yc1 = 0.0;
68   } else if (yc2 > 0.) {
69     xc2 = xc1 + (xc2-xc1)*(0.0-yc1)/(yc2-yc1);
70     yc2 = 0.0;
71   }
72
73   x1 = xc1 / u;
74   y1 = yc1 / v;
75   x2 = xc2 / u;
76   y2 = yc2 / v;
77
78   return true;
79 }
80
81
82 /* NAME
83  *    clip_sector               Clip a line against a sector of a circle
84  *
85  * SYNOPSIS
86  *    clip_sector (x1, y1, x2, y2, u, v)
87  *    double& x1,*y1,*x2,*y2    Endpoints of line
88  *    double u,v                Size of sector
89  */
90
91 int 
92 clip_sector (double& x1, double& y1, double& x2, double& y2, const double u, const double v)
93 {
94   double xc1 = x1 * u;
95   double yc1 = y1 * v;
96   double xc2 = x2 * u;
97   double yc2 = y2 * v;
98   
99   double radius = sqrt (u * u + v * v);
100
101   if (clip_circle (xc1, yc1, xc2, yc2, 0.0, v, radius, 0.0, 0.0) == false)
102     return false;
103
104   if (clip_triangle (xc1, yc1, xc2, yc2, u, v, false) == false)
105     return false;
106   
107   x1 = xc1 / u;
108   y1 = yc1 / v;
109   x2 = xc2 / u;
110   y2 = yc2 / v;
111   
112   return true;
113 }
114
115
116 /* NAME
117  *    clip_circle               Clip a line against a circle
118  *
119  * SYNOPSIS
120  *    clip_circle (x1,y1,x2,y2,cx,cy,radius,t1,t2)
121  *    double& x1,*y1,*x2,*y2    Endpoints of line to be clipped
122  *    double cx,cy              Center of circle
123  *    double radius             Radius of circle
124  *    double t1,t2              Starting & stopping angles of clipping
125  */
126
127 int 
128 clip_circle (double& x1, double& y1, double& x2, double& y2, const double cx, const double cy, const double radius, double t1, double t2)
129 {
130   double xc1 = x1;
131   double yc1 = y1;
132   double xc2 = x2;
133   double yc2 = y2;
134   double ccx = cx;
135   double ccy = cy;
136
137   double xtrans = -xc1;                 // move (x1,y1) to origin 
138   double ytrans = -yc1;
139
140   xc1 += xtrans;
141   yc1 += ytrans;
142   xc2 += xtrans;
143   yc2 += ytrans;
144   ccx += xtrans;
145   ccy += ytrans;
146
147   double theta = -atan2 (yc2, xc2);     // rotate line to lie on x-axis
148   GRFMTX_2D rotmtx;
149   rot_mtx2 (rotmtx, theta);     // xc1,yc1 is at origin, no need to rot 
150   xform_mtx2 (rotmtx, xc2, yc2);
151   xform_mtx2 (rotmtx, ccx, ccy);
152   t1 += theta;                  // rotate start and stop angles 
153   t2 += theta;
154   t1 = normalizeAngle (t1);
155   t2 = normalizeAngle (t2);
156
157   if (xc2 < -D_EPSILON || fabs(yc2) > F_EPSILON) {
158     sys_error (ERR_SEVERE, "Internal error in clip_circle\n x1=%6.2f, y1=%6.2f, x2=%6.2f, y2=%6.2f, xc2=%6.2f, yc2=%6.2f, theta=%6.2f", x1, y1, x2, y2, xc2, yc2, theta);
159     return false;
160   }
161
162   if (fabs(ccy) > radius)               /* check if can reject */
163     return false;
164   
165   double temp = sqrt (radius * radius - ccy * ccy);
166   double xcmin = ccx - temp;
167   double xcmax = ccx + temp;
168   
169   if (fabs(t2 - t1) < D_EPSILON) {
170     if (xc1 < xcmin)
171       xc1 = xcmin;
172     if (xc2 > xcmax)
173       xc2 = xcmax;
174   } else if (t1 < t2) {
175     if (t1 < PI && t2 > PI)
176       if (xc1 < xcmin)
177         xc1 = xcmin;
178   } else if (t1 > t2) {
179     if (t1 < PI)
180       if (xc1 < xcmin)
181         xc1 = xcmin;
182     if (xc2 > xcmax)
183       xc2 = xcmax;
184   }
185
186   rot_mtx2 (rotmtx, -theta);
187   xform_mtx2 (rotmtx, xc1, yc1);
188   xform_mtx2 (rotmtx, xc2, yc2);
189   xc1 += -xtrans;
190   yc1 += -ytrans;
191   xc2 += -xtrans;
192   yc2 += -ytrans;
193
194   x1 = xc1;
195   y1 = yc1;
196   x2 = xc2;
197   y2 = yc2;
198
199   return true;
200 }
201
202
203 /* NAME
204  *    clip_triangle             Clip a line against a triangle
205  *
206  * SYNOPSIS
207  *    clip_triangle (x1, y1, x2, y2, u, v, clip_xaxis)
208  *    double& x1, *y1, *x2, *y2 Endpoints of line
209  *    double u, v               Size of 1/2 base len & height
210  *    int clip_xaxis            Boolean flag whether to clip against x axis
211  *                              (Use true for all triangles)
212  *                              (false if used internally by sector clipping routine)
213  *
214  * DESCRIPTION
215  *              x
216  *             /|\              Note that vertices of triangle are
217  *            / | \                 (-u, 0)
218  *           /  |  \                (u, 0)
219  *          /   |   \               (0, v)
220  *         /    | v  \
221  *        /     |     \
222  *       +------+------+
223  *            (0,0)  u
224  *
225  * NOTES
226  *   1) Inside of this routine, values of (u,v) are assumed to be (1,1)
227  *
228  *   2) Derivation of clipping equations:
229  *      Using parametric equations for the line
230  *          xv = x1 + t * (x2 - x1)
231  *          yv = y1 + t * (y2 - y1)
232  *      so,
233  *          t  = (xv - x1) / (x2 - x1)
234  *          yv = y1 + (xv - x1) * (y2 - y1) / (x2 - x1)
235  *          yv = y1 + (xv - x1) * dy / dx
236  *
237  *      Now, find the intersections with the following clipping boundries:
238  *          yv = v - (v/u) * xv         (yv = mx + b)
239  *          yv = v + (v/u) * xv         (m = v/u, b = v);
240  */
241
242 static int tcode (const double x, const double y, const double m, const double b, const int clip_xaxis);
243
244 int 
245 clip_triangle (double& x1, double& y1, double& x2, double& y2, const double u, const double v, const int clip_xaxis)
246 {
247   double m = v / u;      // slope of triangle lines
248   double b = v;          // y-intercept of triangle lines 
249
250   int c1 = tcode (x1, y1, m, b, clip_xaxis);
251   int c2 = tcode (x2, y2, m, b, clip_xaxis);
252
253 #ifdef DEBUG
254   printf ("x1:%6.2f  y1:%6.2f  code1:%2d  x2:%6.2f  y2:%6.2f code2:%2d\n",
255            x1, y1, c1, x2, y2, c2);
256 #endif
257   while ( c1 || c2 ) {
258     if ( c1 & c2 ) {
259       return false;                     // trivial reject 
260     }
261     int c = c1;
262     if (c1 == 0)
263       c = c2;
264
265     double x = 0, y = 0;
266     if (c & 1) {                        // below 
267       x = x1 + (x2-x1)*(0.0-y1)/(y2-y1);
268       y = 0.0;
269     } else if (c & 2) {                 // right 
270       double dx, dy;
271       dx = x2 - x1;
272       dy = y2 - y1;
273       if (fabs(dx) > D_EPSILON)
274         x = (-y1 + b + x1 * dy/dx) / (m + dy/dx);
275       else
276         x = x1;
277       y = -m * x + b;
278     } else if (c & 4) {                 /* left */
279       double dx, dy;
280       dx = x2 - x1;
281       dy = y2 - y1;
282       if (fabs(dx) > D_EPSILON) {
283         x = (y1 - b - x1 * dy/dx);
284         x /= (m - dy/dx);
285       } else
286         x = x1;
287       y = m * x + b;
288     }
289     
290     if (c == c1) {
291       x1=x; y1=y; c1=tcode (x1,y1,m,b,clip_xaxis);
292     } else {
293       x2=x; y2=y; c2=tcode (x2,y2,m,b,clip_xaxis);
294     }
295 #ifdef DEBUG
296     printf ("x1:%6.2f  y1:%6.2f  code1:%2d  x2:%6.2f  y2:%6.2f code2:%2d\n", x1, y1, c1, x2, y2, c2);
297 #endif
298   }
299
300   return true;          /* we have clipped the line, and it is good */
301 }
302
303
304 /* compute region code */
305 static int 
306 tcode (const double x, const double y, const double m, const double b, const int clip_xaxis)
307 {
308   int c = 0;
309
310   if (clip_xaxis && y < 0.)     // below triange 
311     c = 1;
312
313   if (y > -m * x + b + D_EPSILON)               // right of triangle 
314     c += 2;
315   if (y > m * x + b + D_EPSILON)                // left of triangle 
316     c += 4;
317   
318   return (c);
319 }  
320
321
322 /* NAME
323  *    clip_rect                 Clip a line against a rectangle
324  *
325  * SYNOPSIS
326  *    clip_rect (x1, y1, x2, y2, rect)
327  *    double& x1, *y1, *x2, *y2 Endpoints of line
328  *    double rect[4]            Rectangle to clip against
329  *                              ordered xmin, ymin, xmax, ymax
330  */
331
332 static int rectcode (double x, double y, const double rect[4]);
333
334 bool\r
335 clip_rect (double& x1, double& y1, double& x2, double& y2, const double rect[4])
336 {
337   double x = 0, y = 0;
338
339   int c1 = rectcode (x1, y1, rect);
340   int c2 = rectcode (x2, y2, rect);
341
342   while (c1 || c2) {
343     if (c1 & c2)
344       return false;                     // trivial reject 
345     int c = c1;
346     if (c1 == 0)
347       c = c2;
348     if (c & 1) {                        // left 
349       y = y1 + (y2-y1)*(rect[0]-x1)/(x2-x1);
350       x = rect[0];
351     } else if (c & 2) {                 // right 
352       y = y1 + (y2-y1)*(rect[2]-x1)/(x2-x1);
353       x = rect[2];
354     } else if (c & 4) {                 // bottom 
355       x = x1 + (x2-x1)*(rect[1]-y1)/(y2-y1);
356       y = rect[1];
357     } else if (c & 8) {                 // top 
358       x = x1 + (x2-x1)*(rect[3]-y1)/(y2-y1);
359       y = rect[3];
360     }
361     
362     if (c == c1) {
363       x1=x; y1=y; c1=rectcode(x1,y1,rect);
364     } else {
365       x2=x; y2=y; c2=rectcode(x2,y2,rect);
366     }
367   }
368   return true;          // we have clipped the line, and it is good 
369 }
370
371
372 /* NAME
373  *   rectcode                   INTERNAL routine to return position of
374  *                              point relative to a rectangle
375  *
376  * SYNOPSIS
377  *    c = rectcode (x, y, rect)
378  *    int c                     Position of point relative to the window
379  *    double x, y               Point to test against window
380  *    double rect[4]            Coordinates of rectangle extent
381  *                              Ordered [xmin, ymin, xmax, ymax]
382  */
383
384 static int 
385 rectcode (double x, double y, const double rect[4]) 
386 {
387   int c = 0;
388
389   if (x < rect[0])
390     c = 1;
391   else if (x > rect[2])
392     c = 2;
393   if (y < rect[1])
394     c += 4;
395   else if (y > rect[3])
396     c += 8;
397   return (c);
398 }