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