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-2001 Kevin Rosenberg
  16 **
  17 **  $Id$
  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 bool
  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 bool
  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 bool
 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 bool
 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 #if 0
 254   printf ("x1:%6.2f  y1:%6.2f  code1:%2d  x2:%6.2f  y2:%6.2f code2:%2d\n", x1, y1, c1, x2, y2, c2);
 255 #endif
 256   while ( c1 || c2 ) {
 257     if ( c1 & c2 ) {
 258       return false;                     // trivial reject
 259     }
 260     int c = c1;
 261     if (c1 == 0)
 262       c = c2;
 263
 264     double x = 0, y = 0;
 265     if (c & 1) {                        // below
 266       x = x1 + (x2-x1)*(0.0-y1)/(y2-y1);
 267       y = 0.0;
 268     } else if (c & 2) {                 // right
 269       double dx, dy;
 270       dx = x2 - x1;
 271       dy = y2 - y1;
 272       if (fabs(dx) > D_EPSILON)
 273         x = (-y1 + b + x1 * dy/dx) / (m + dy/dx);
 274       else
 275         x = x1;
 276       y = -m * x + b;
 277     } else if (c & 4) {                 /* left */
 278       double dx, dy;
 279       dx = x2 - x1;
 280       dy = y2 - y1;
 281       if (fabs(dx) > D_EPSILON) {
 282         x = (y1 - b - x1 * dy/dx);
 283         x /= (m - dy/dx);
 284       } else
 285         x = x1;
 286       y = m * x + b;
 287     }
 288
 289     if (c == c1) {
 290       x1=x; y1=y; c1=tcode (x1,y1,m,b,clip_xaxis);
 291     } else {
 292       x2=x; y2=y; c2=tcode (x2,y2,m,b,clip_xaxis);
 293     }
 294 #if 0
 295     printf ("x1:%6.2f  y1:%6.2f  code1:%2d  x2:%6.2f  y2:%6.2f code2:%2d\n", x1, y1, c1, x2, y2, c2);
 296 #endif
 297   }
 298
 299   return true;          /* we have clipped the line, and it is good */
 300 }
 301
 302
 303 /* compute region code */
 304 static int
 305 tcode (const double x, const double y, const double m, const double b, const int clip_xaxis)
 306 {
 307   int c = 0;
 308
 309   if (clip_xaxis && y < 0.)     // below triange
 310     c = 1;
 311
 312   if (y > -m * x + b + D_EPSILON)               // right of triangle
 313     c += 2;
 314   if (y > m * x + b + D_EPSILON)                // left of triangle
 315     c += 4;
 316
 317   return (c);
 318 }
 319
 320
 321 /* NAME
 322  *    clip_rect                 Clip a line against a rectangle
 323  *
 324  * SYNOPSIS
 325  *    clip_rect (x1, y1, x2, y2, rect)
 326  *    double& x1, *y1, *x2, *y2 Endpoints of line
 327  *    double rect[4]            Rectangle to clip against
 328  *                              ordered xmin, ymin, xmax, ymax
 329  */
 330
 331 static int rectcode (double x, double y, const double rect[4]);
 332
 333 bool
 334 clip_rect (double& x1, double& y1, double& x2, double& y2, const double rect[4])
 335 {
 336   double x = 0, y = 0;
 337
 338   int c1 = rectcode (x1, y1, rect);
 339   int c2 = rectcode (x2, y2, rect);
 340
 341   while (c1 || c2) {
 342     if (c1 & c2)
 343       return false;                     // trivial reject
 344     int c = c1;
 345     if (c1 == 0)
 346       c = c2;
 347     if (c & 1) {                        // left
 348       y = y1 + (y2-y1)*(rect[0]-x1)/(x2-x1);
 349       x = rect[0];
 350     } else if (c & 2) {                 // right
 351       y = y1 + (y2-y1)*(rect[2]-x1)/(x2-x1);
 352       x = rect[2];
 353     } else if (c & 4) {                 // bottom
 354       x = x1 + (x2-x1)*(rect[1]-y1)/(y2-y1);
 355       y = rect[1];
 356     } else if (c & 8) {                 // top
 357       x = x1 + (x2-x1)*(rect[3]-y1)/(y2-y1);
 358       y = rect[3];
 359     }
 360
 361     if (c == c1) {
 362       x1=x; y1=y; c1=rectcode(x1,y1,rect);
 363     } else {
 364       x2=x; y2=y; c2=rectcode(x2,y2,rect);
 365     }
 366   }
 367   return true;          // we have clipped the line, and it is good
 368 }
 369
 370
 371 /* NAME
 372  *   rectcode                   INTERNAL routine to return position of
 373  *                              point relative to a rectangle
 374  *
 375  * SYNOPSIS
 376  *    c = rectcode (x, y, rect)
 377  *    int c                     Position of point relative to the window
 378  *    double x, y               Point to test against window
 379  *    double rect[4]            Coordinates of rectangle extent
 380  *                              Ordered [xmin, ymin, xmax, ymax]
 381  */
 382
 383 static int
 384 rectcode (double x, double y, const double rect[4])
 385 {
 386   int c = 0;
 387
 388   if (x < rect[0])
 389     c = 1;
 390   else if (x > rect[2])
 391     c = 2;
 392   if (y < rect[1])
 393     c += 4;
 394   else if (y > rect[3])
 395     c += 8;
 396   return (c);
 397 }