1 /* 2 ** License Applicability. Except to the extent portions of this file are 3 ** made subject to an alternative license as permitted in the SGI Free 4 ** Software License B, Version 1.1 (the "License"), the contents of this 5 ** file are subject only to the provisions of the License. You may not use 6 ** this file except in compliance with the License. You may obtain a copy 7 ** of the License at Silicon Graphics, Inc., attn: Legal Services, 1600 8 ** Amphitheatre Parkway, Mountain View, CA 94043-1351, or at: 9 ** 10 ** http://oss.sgi.com/projects/FreeB 11 ** 12 ** Note that, as provided in the License, the Software is distributed on an 13 ** "AS IS" basis, with ALL EXPRESS AND IMPLIED WARRANTIES AND CONDITIONS 14 ** DISCLAIMED, INCLUDING, WITHOUT LIMITATION, ANY IMPLIED WARRANTIES AND 15 ** CONDITIONS OF MERCHANTABILITY, SATISFACTORY QUALITY, FITNESS FOR A 16 ** PARTICULAR PURPOSE, AND NON-INFRINGEMENT. 17 ** 18 ** Original Code. The Original Code is: OpenGL Sample Implementation, 19 ** Version 1.2.1, released January 26, 2000, developed by Silicon Graphics, 20 ** Inc. The Original Code is Copyright (c) 1991-2000 Silicon Graphics, Inc. 21 ** Copyright in any portions created by third parties is as indicated 22 ** elsewhere herein. All Rights Reserved. 23 ** 24 ** Additional Notice Provisions: The application programming interfaces 25 ** established by SGI in conjunction with the Original Code are The 26 ** OpenGL(R) Graphics System: A Specification (Version 1.2.1), released 27 ** April 1, 1999; The OpenGL(R) Graphics System Utility Library (Version 28 ** 1.3), released November 4, 1998; and OpenGL(R) Graphics with the X 29 ** Window System(R) (Version 1.3), released October 19, 1998. This software 30 ** was created using the OpenGL(R) version 1.2.1 Sample Implementation 31 ** published by SGI, but has not been independently verified as being 32 ** compliant with the OpenGL(R) version 1.2.1 Specification. 33 */ 34 35 /* 36 * mapdesc.c++ 37 * 38 */ 39 40 //#include <stdio.h> 41 //#include "glimports.h" 42 //#include "mystdio.h" 43 //#include "myassert.h" 44 //#include "mystring.h" 45 #include "mymath.h" 46 #include "backend.h" 47 //#include "nurbsconsts.h" 48 #include "mapdesc.h" 49 50 Mapdesc::Mapdesc( long _type, int _israt, int _ncoords, Backend& b ) 51 : backend( b ) 52 { 53 type = _type; 54 isrational = _israt; 55 ncoords = _ncoords; 56 hcoords = _ncoords + (_israt ? 0 : 1 ); 57 inhcoords = _ncoords - (_israt ? 1 : 0 ); 58 mask = ((1<<(inhcoords*2))-1); 59 next = 0; 60 61 assert( hcoords <= MAXCOORDS ); 62 assert( inhcoords >= 1 ); 63 64 pixel_tolerance = 1.0; 65 error_tolerance = 1.0; 66 bbox_subdividing = N_NOBBOXSUBDIVISION; 67 culling_method = N_NOCULLING; 68 sampling_method = N_NOSAMPLING; 69 clampfactor = N_NOCLAMPING; 70 minsavings = N_NOSAVINGSSUBDIVISION; 71 s_steps = 0.0; 72 t_steps = 0.0; 73 maxrate = ( s_steps < 0.0 ) ? 0.0 : s_steps; 74 maxsrate = ( s_steps < 0.0 ) ? 0.0 : s_steps; 75 maxtrate = ( t_steps < 0.0 ) ? 0.0 : t_steps; 76 identify( bmat ); 77 identify( cmat ); 78 identify( smat ); 79 for( int i = 0; i != inhcoords; i++ ) 80 bboxsize[i] = 1.0; 81 } 82 83 void 84 Mapdesc::setBboxsize( INREAL *mat ) 85 { 86 for( int i = 0; i != inhcoords; i++ ) 87 bboxsize[i] = (REAL) mat[i]; 88 } 89 90 void 91 Mapdesc::identify( REAL dest[MAXCOORDS][MAXCOORDS] ) 92 { 93 memset( dest, 0, sizeof( REAL ) * MAXCOORDS * MAXCOORDS ); 94 for( int i=0; i != hcoords; i++ ) 95 dest[i][i] = 1.0; 96 } 97 98 void 99 Mapdesc::surfbbox( REAL bb[2][MAXCOORDS] ) 100 { 101 backend.surfbbox( type, bb[0], bb[1] ); 102 } 103 104 void 105 Mapdesc::copy( REAL dest[MAXCOORDS][MAXCOORDS], long n, INREAL *src, 106 long rstride, long cstride ) 107 { 108 assert( n >= 0 ); 109 for( int i=0; i != n; i++ ) 110 for( int j=0; j != n; j++ ) 111 dest[i][j] = src[i*rstride + j*cstride]; 112 } 113 114 /*-------------------------------------------------------------------------- 115 * copyPt - copy a homogeneous point 116 *-------------------------------------------------------------------------- 117 */ 118 void 119 Mapdesc::copyPt( REAL *d, REAL *s ) 120 { 121 assert( hcoords > 0 ); 122 switch( hcoords ) { 123 case 4: 124 d[3] = s[3]; 125 d[2] = s[2]; 126 d[1] = s[1]; 127 d[0] = s[0]; 128 break; 129 case 3: 130 d[2] = s[2]; 131 d[1] = s[1]; 132 d[0] = s[0]; 133 break; 134 case 2: 135 d[1] = s[1]; 136 d[0] = s[0]; 137 break; 138 case 1: 139 d[0] = s[0]; 140 break; 141 case 5: 142 d[4] = s[4]; 143 d[3] = s[3]; 144 d[2] = s[2]; 145 d[1] = s[1]; 146 d[0] = s[0]; 147 break; 148 default: 149 memcpy( d, s, hcoords * sizeof( REAL ) ); 150 break; 151 } 152 } 153 154 /*-------------------------------------------------------------------------- 155 * sumPt - compute affine combination of two homogeneous points 156 *-------------------------------------------------------------------------- 157 */ 158 void 159 Mapdesc::sumPt( REAL *dst, REAL *src1, REAL *src2, REAL alpha, REAL beta ) 160 { 161 assert( hcoords > 0 ); 162 switch( hcoords ) { 163 case 4: 164 dst[3] = src1[3] * alpha + src2[3] * beta; 165 dst[2] = src1[2] * alpha + src2[2] * beta; 166 dst[1] = src1[1] * alpha + src2[1] * beta; 167 dst[0] = src1[0] * alpha + src2[0] * beta; 168 break; 169 case 3: 170 dst[2] = src1[2] * alpha + src2[2] * beta; 171 dst[1] = src1[1] * alpha + src2[1] * beta; 172 dst[0] = src1[0] * alpha + src2[0] * beta; 173 break; 174 case 2: 175 dst[1] = src1[1] * alpha + src2[1] * beta; 176 dst[0] = src1[0] * alpha + src2[0] * beta; 177 break; 178 case 1: 179 dst[0] = src1[0] * alpha + src2[0] * beta; 180 break; 181 case 5: 182 dst[4] = src1[4] * alpha + src2[4] * beta; 183 dst[3] = src1[3] * alpha + src2[3] * beta; 184 dst[2] = src1[2] * alpha + src2[2] * beta; 185 dst[1] = src1[1] * alpha + src2[1] * beta; 186 dst[0] = src1[0] * alpha + src2[0] * beta; 187 break; 188 default: { 189 for( int i = 0; i != hcoords; i++ ) 190 dst[i] = src1[i] * alpha + src2[i] * beta; 191 } 192 break; 193 } 194 } 195 196 /*-------------------------------------------------------------------------- 197 * clipbits - compute bit-vector indicating point/window position 198 * of a (transformed) homogeneous point 199 *-------------------------------------------------------------------------- 200 */ 201 unsigned int 202 Mapdesc::clipbits( REAL *p ) 203 { 204 assert( inhcoords >= 0 ); 205 assert( inhcoords <= 3 ); 206 207 int nc = inhcoords; 208 REAL pw = p[nc]; 209 REAL nw = -pw; 210 unsigned int bits = 0; 211 212 if( pw == 0.0 ) return mask; 213 214 if( pw > 0.0 ) { 215 switch( nc ) { 216 case 3: 217 if( p[2] <= pw ) bits |= (1<<5); 218 if( p[2] >= nw ) bits |= (1<<4); 219 if( p[1] <= pw ) bits |= (1<<3); 220 if( p[1] >= nw ) bits |= (1<<2); 221 if( p[0] <= pw ) bits |= (1<<1); 222 if( p[0] >= nw ) bits |= (1<<0); 223 return bits; 224 case 2: 225 if( p[1] <= pw ) bits |= (1<<3); 226 if( p[1] >= nw ) bits |= (1<<2); 227 if( p[0] <= pw ) bits |= (1<<1); 228 if( p[0] >= nw ) bits |= (1<<0); 229 return bits; 230 case 1: 231 if( p[0] <= pw ) bits |= (1<<1); 232 if( p[0] >= nw ) bits |= (1<<0); 233 return bits; 234 default: { 235 int bit = 1; 236 for( int i=0; i<nc; i++ ) { 237 if( p[i] >= nw ) bits |= bit; 238 bit <<= 1; 239 if( p[i] <= pw ) bits |= bit; 240 bit <<= 1; 241 } 242 abort(); 243 break; 244 } 245 } 246 } else { 247 switch( nc ) { 248 case 3: 249 if( p[2] <= nw ) bits |= (1<<5); 250 if( p[2] >= pw ) bits |= (1<<4); 251 if( p[1] <= nw ) bits |= (1<<3); 252 if( p[1] >= pw ) bits |= (1<<2); 253 if( p[0] <= nw ) bits |= (1<<1); 254 if( p[0] >= pw ) bits |= (1<<0); 255 return bits; 256 case 2: 257 if( p[1] <= nw ) bits |= (1<<3); 258 if( p[1] >= pw ) bits |= (1<<2); 259 if( p[0] <= nw ) bits |= (1<<1); 260 if( p[0] >= pw ) bits |= (1<<0); 261 return bits; 262 case 1: 263 if( p[0] <= nw ) bits |= (1<<1); 264 if( p[0] >= pw ) bits |= (1<<0); 265 return bits; 266 default: { 267 int bit = 1; 268 for( int i=0; i<nc; i++ ) { 269 if( p[i] >= pw ) bits |= bit; 270 bit <<= 1; 271 if( p[i] <= nw ) bits |= bit; 272 bit <<= 1; 273 } 274 abort(); 275 break; 276 } 277 } 278 } 279 return bits; 280 } 281 282 /*-------------------------------------------------------------------------- 283 * xformRational - transform a homogeneous point 284 *-------------------------------------------------------------------------- 285 */ 286 void 287 Mapdesc::xformRational( Maxmatrix mat, REAL *d, REAL *s ) 288 { 289 assert( hcoords >= 0 ); 290 291 if( hcoords == 3 ) { 292 REAL x = s[0]; 293 REAL y = s[1]; 294 REAL z = s[2]; 295 d[0] = x*mat[0][0]+y*mat[1][0]+z*mat[2][0]; 296 d[1] = x*mat[0][1]+y*mat[1][1]+z*mat[2][1]; 297 d[2] = x*mat[0][2]+y*mat[1][2]+z*mat[2][2]; 298 } else if( hcoords == 4 ) { 299 REAL x = s[0]; 300 REAL y = s[1]; 301 REAL z = s[2]; 302 REAL w = s[3]; 303 d[0] = x*mat[0][0]+y*mat[1][0]+z*mat[2][0]+w*mat[3][0]; 304 d[1] = x*mat[0][1]+y*mat[1][1]+z*mat[2][1]+w*mat[3][1]; 305 d[2] = x*mat[0][2]+y*mat[1][2]+z*mat[2][2]+w*mat[3][2]; 306 d[3] = x*mat[0][3]+y*mat[1][3]+z*mat[2][3]+w*mat[3][3]; 307 } else { 308 for( int i=0; i != hcoords; i++ ) { 309 d[i] = 0; 310 for( int j = 0; j != hcoords; j++ ) 311 d[i] += s[j] * mat[j][i]; 312 } 313 } 314 } 315 316 /*-------------------------------------------------------------------------- 317 * xformNonrational - transform a inhomogeneous point to a homogeneous point 318 *-------------------------------------------------------------------------- 319 */ 320 void 321 Mapdesc::xformNonrational( Maxmatrix mat, REAL *d, REAL *s ) 322 { 323 if( inhcoords == 2 ) { 324 REAL x = s[0]; 325 REAL y = s[1]; 326 d[0] = x*mat[0][0]+y*mat[1][0]+mat[2][0]; 327 d[1] = x*mat[0][1]+y*mat[1][1]+mat[2][1]; 328 d[2] = x*mat[0][2]+y*mat[1][2]+mat[2][2]; 329 } else if( inhcoords == 3 ) { 330 REAL x = s[0]; 331 REAL y = s[1]; 332 REAL z = s[2]; 333 d[0] = x*mat[0][0]+y*mat[1][0]+z*mat[2][0]+mat[3][0]; 334 d[1] = x*mat[0][1]+y*mat[1][1]+z*mat[2][1]+mat[3][1]; 335 d[2] = x*mat[0][2]+y*mat[1][2]+z*mat[2][2]+mat[3][2]; 336 d[3] = x*mat[0][3]+y*mat[1][3]+z*mat[2][3]+mat[3][3]; 337 } else { 338 assert( inhcoords >= 0 ); 339 for( int i=0; i != hcoords; i++ ) { 340 d[i] = mat[inhcoords][i]; 341 for( int j = 0; j < inhcoords; j++ ) 342 d[i] += s[j] * mat[j][i]; 343 } 344 } 345 } 346 347 /*-------------------------------------------------------------------------- 348 * xformAndCullCheck - transform a set of points that may be EITHER 349 * homogeneous or inhomogeneous depending on the map description and 350 * check if they are either completely inside, completely outside, 351 * or intersecting the viewing frustrum. 352 *-------------------------------------------------------------------------- 353 */ 354 int 355 Mapdesc::xformAndCullCheck( 356 REAL *pts, int uorder, int ustride, int vorder, int vstride ) 357 { 358 assert( uorder > 0 ); 359 assert( vorder > 0 ); 360 361 unsigned int inbits = mask; 362 unsigned int outbits = 0; 363 364 REAL *p = pts; 365 for( REAL *pend = p + uorder * ustride; p != pend; p += ustride ) { 366 REAL *q = p; 367 for( REAL *qend = q + vorder * vstride; q != qend; q += vstride ) { 368 REAL cpts[MAXCOORDS]; 369 xformCulling( cpts, q ); 370 unsigned int bits = clipbits( cpts ); 371 outbits |= bits; 372 inbits &= bits; 373 if( ( outbits == (unsigned int)mask ) && ( inbits != (unsigned int)mask ) ) return CULL_ACCEPT; 374 } 375 } 376 377 if( outbits != (unsigned int)mask ) { 378 return CULL_TRIVIAL_REJECT; 379 } else if( inbits == (unsigned int)mask ) { 380 return CULL_TRIVIAL_ACCEPT; 381 } else { 382 return CULL_ACCEPT; 383 } 384 } 385 386 /*-------------------------------------------------------------------------- 387 * cullCheck - check if a set of homogeneous transformed points are 388 * either completely inside, completely outside, 389 * or intersecting the viewing frustrum. 390 *-------------------------------------------------------------------------- 391 */ 392 int 393 Mapdesc::cullCheck( REAL *pts, int uorder, int ustride, int vorder, int vstride ) 394 { 395 unsigned int inbits = mask; 396 unsigned int outbits = 0; 397 398 REAL *p = pts; 399 for( REAL *pend = p + uorder * ustride; p != pend; p += ustride ) { 400 REAL *q = p; 401 for( REAL *qend = q + vorder * vstride; q != qend; q += vstride ) { 402 unsigned int bits = clipbits( q ); 403 outbits |= bits; 404 inbits &= bits; 405 if( ( outbits == (unsigned int)mask ) && ( inbits != (unsigned int)mask ) ) return CULL_ACCEPT; 406 } 407 } 408 409 if( outbits != (unsigned int)mask ) { 410 return CULL_TRIVIAL_REJECT; 411 } else if( inbits == (unsigned int)mask ) { 412 return CULL_TRIVIAL_ACCEPT; 413 } else { 414 return CULL_ACCEPT; 415 } 416 } 417 418 /*-------------------------------------------------------------------------- 419 * cullCheck - check if a set of homogeneous transformed points are 420 * either completely inside, completely outside, 421 * or intersecting the viewing frustrum. 422 *-------------------------------------------------------------------------- 423 */ 424 int 425 Mapdesc::cullCheck( REAL *pts, int order, int stride ) 426 { 427 unsigned int inbits = mask; 428 unsigned int outbits = 0; 429 430 REAL *p = pts; 431 for( REAL *pend = p + order * stride; p != pend; p += stride ) { 432 unsigned int bits = clipbits( p ); 433 outbits |= bits; 434 inbits &= bits; 435 if( ( outbits == (unsigned int)mask ) && ( inbits != (unsigned int)mask ) ) return CULL_ACCEPT; 436 } 437 438 if( outbits != (unsigned int)mask ) { 439 return CULL_TRIVIAL_REJECT; 440 } else if( inbits == (unsigned int)mask ) { 441 return CULL_TRIVIAL_ACCEPT; 442 } else { 443 return CULL_ACCEPT; 444 } 445 } 446 447 /*-------------------------------------------------------------------------- 448 * xformSampling - transform a set of points that may be EITHER 449 * homogeneous or inhomogeneous depending on the map description 450 * into sampling space 451 *-------------------------------------------------------------------------- 452 */ 453 void 454 Mapdesc::xformSampling( REAL *pts, int order, int stride, REAL *sp, int outstride ) 455 { 456 xformMat( smat, pts, order, stride, sp, outstride ); 457 } 458 459 void 460 Mapdesc::xformBounding( REAL *pts, int order, int stride, REAL *sp, int outstride ) 461 { 462 xformMat( bmat, pts, order, stride, sp, outstride ); 463 } 464 465 /*-------------------------------------------------------------------------- 466 * xformCulling - transform a set of points that may be EITHER 467 * homogeneous or inhomogeneous depending on the map description 468 * into culling space 469 *-------------------------------------------------------------------------- 470 */ 471 void 472 Mapdesc::xformCulling( REAL *pts, int order, int stride, REAL *cp, int outstride ) 473 { 474 xformMat( cmat, pts, order, stride, cp, outstride ); 475 } 476 477 /*-------------------------------------------------------------------------- 478 * xformCulling - transform a set of points that may be EITHER 479 * homogeneous or inhomogeneous depending on the map description 480 * into culling space 481 *-------------------------------------------------------------------------- 482 */ 483 void 484 Mapdesc::xformCulling( REAL *pts, 485 int uorder, int ustride, 486 int vorder, int vstride, 487 REAL *cp, int outustride, int outvstride ) 488 { 489 xformMat( cmat, pts, uorder, ustride, vorder, vstride, cp, outustride, outvstride ); 490 } 491 492 /*-------------------------------------------------------------------------- 493 * xformSampling - transform a set of points that may be EITHER 494 * homogeneous or inhomogeneous depending on the map description 495 * into sampling space 496 *-------------------------------------------------------------------------- 497 */ 498 void 499 Mapdesc::xformSampling( REAL *pts, 500 int uorder, int ustride, 501 int vorder, int vstride, 502 REAL *sp, int outustride, int outvstride ) 503 { 504 xformMat( smat, pts, uorder, ustride, vorder, vstride, sp, outustride, outvstride ); 505 } 506 507 void 508 Mapdesc::xformBounding( REAL *pts, 509 int uorder, int ustride, 510 int vorder, int vstride, 511 REAL *sp, int outustride, int outvstride ) 512 { 513 xformMat( bmat, pts, uorder, ustride, vorder, vstride, sp, outustride, outvstride ); 514 } 515 516 void 517 Mapdesc::xformMat( 518 Maxmatrix mat, 519 REAL * pts, 520 int order, 521 int stride, 522 REAL * cp, 523 int outstride ) 524 { 525 if( isrational ) { 526 REAL *pend = pts + order * stride; 527 for( REAL *p = pts ; p != pend; p += stride ) { 528 xformRational( mat, cp, p ); 529 cp += outstride; 530 } 531 } else { 532 REAL *pend = pts + order * stride; 533 for( REAL *p = pts ; p != pend; p += stride ) { 534 xformNonrational( mat, cp, p ); 535 cp += outstride; 536 } 537 } 538 } 539 540 void 541 Mapdesc::xformMat( Maxmatrix mat, REAL *pts, 542 int uorder, int ustride, 543 int vorder, int vstride, 544 REAL *cp, int outustride, int outvstride ) 545 { 546 if( isrational ) { 547 REAL *pend = pts + uorder * ustride; 548 for( REAL *p = pts ; p != pend; p += ustride ) { 549 REAL *cpts2 = cp; 550 REAL *qend = p + vorder * vstride; 551 for( REAL *q = p; q != qend; q += vstride ) { 552 xformRational( mat, cpts2, q ); 553 cpts2 += outvstride; 554 } 555 cp += outustride; 556 } 557 } else { 558 REAL *pend = pts + uorder * ustride; 559 for( REAL *p = pts ; p != pend; p += ustride ) { 560 REAL *cpts2 = cp; 561 REAL *qend = p + vorder * vstride; 562 for( REAL *q = p; q != qend; q += vstride ) { 563 xformNonrational( mat, cpts2, q ); 564 cpts2 += outvstride; 565 } 566 cp += outustride; 567 } 568 } 569 } 570 571 /*-------------------------------------------------------------------------- 572 * subdivide - subdivide a curve along an isoparametric line 573 *-------------------------------------------------------------------------- 574 */ 575 576 void 577 Mapdesc::subdivide( REAL *src, REAL *dst, REAL v, int stride, int order ) 578 { 579 REAL mv = 1.0 - v; 580 581 for( REAL *send=src+stride*order; src!=send; send-=stride, dst+=stride ) { 582 copyPt( dst, src ); 583 REAL *qpnt = src + stride; 584 for( REAL *qp=src; qpnt!=send; qp=qpnt, qpnt+=stride ) 585 sumPt( qp, qp, qpnt, mv, v ); 586 } 587 } 588 589 /*-------------------------------------------------------------------------- 590 * subdivide - subdivide a patch along an isoparametric line 591 *-------------------------------------------------------------------------- 592 */ 593 594 void 595 Mapdesc::subdivide( REAL *src, REAL *dst, REAL v, 596 int so, int ss, int to, int ts ) 597 { 598 REAL mv = 1.0 - v; 599 600 for( REAL *slast = src+ss*so; src != slast; src += ss, dst += ss ) { 601 REAL *sp = src; 602 REAL *dp = dst; 603 for( REAL *send = src+ts*to; sp != send; send -= ts, dp += ts ) { 604 copyPt( dp, sp ); 605 REAL *qp = sp; 606 for( REAL *qpnt = sp+ts; qpnt != send; qp = qpnt, qpnt += ts ) 607 sumPt( qp, qp, qpnt, mv, v ); 608 } 609 } 610 } 611 612 613 #define sign(x) ((x > 0) ? 1 : ((x < 0.0) ? -1 : 0)) 614 615 /*-------------------------------------------------------------------------- 616 * project - project a set of homogeneous coordinates into inhomogeneous ones 617 *-------------------------------------------------------------------------- 618 */ 619 int 620 Mapdesc::project( REAL *src, int rstride, int cstride, 621 REAL *dest, int trstride, int tcstride, 622 int nrows, int ncols ) 623 { 624 int s = sign( src[inhcoords] ); 625 REAL *rlast = src + nrows * rstride; 626 REAL *trptr = dest; 627 for( REAL *rptr=src; rptr != rlast; rptr+=rstride, trptr+=trstride ) { 628 REAL *clast = rptr + ncols * cstride; 629 REAL *tcptr = trptr; 630 for( REAL *cptr = rptr; cptr != clast; cptr+=cstride, tcptr+=tcstride ) { 631 REAL *coordlast = cptr + inhcoords; 632 if( sign( *coordlast ) != s ) return 0; 633 REAL *tcoord = tcptr; 634 for( REAL *coord = cptr; coord != coordlast; coord++, tcoord++ ) { 635 *tcoord = *coord / *coordlast; 636 } 637 } 638 } 639 return 1; 640 } 641 642 /*-------------------------------------------------------------------------- 643 * project - project a set of homogeneous coordinates into inhomogeneous ones 644 *-------------------------------------------------------------------------- 645 */ 646 int 647 Mapdesc::project( REAL *src, int stride, REAL *dest, int tstride, int ncols ) 648 { 649 int s = sign( src[inhcoords] ); 650 651 REAL *clast = src + ncols * stride; 652 for( REAL *cptr = src, *tcptr = dest; cptr != clast; cptr+=stride, tcptr+=tstride ) { 653 REAL *coordlast = cptr + inhcoords; 654 if( sign( *coordlast ) != s ) return 0; 655 for( REAL *coord = cptr, *tcoord = tcptr; coord != coordlast; coord++, tcoord++ ) 656 *tcoord = *coord / *coordlast; 657 } 658 659 return 1; 660 } 661 662 int 663 Mapdesc::bboxTooBig( 664 REAL *p, 665 int rstride, 666 int cstride, 667 int nrows, 668 int ncols, 669 REAL bb[2][MAXCOORDS] ) 670 { 671 REAL bbpts[MAXORDER][MAXORDER][MAXCOORDS]; 672 const int trstride = sizeof(bbpts[0]) / sizeof(REAL); 673 const int tcstride = sizeof(bbpts[0][0]) / sizeof(REAL); 674 675 // points have been transformed, therefore they are homogeneous 676 // project points 677 int val = project( p, rstride, cstride, 678 &bbpts[0][0][0], trstride, tcstride, nrows, ncols ); 679 if( val == 0 ) return -1; 680 681 // compute bounding box 682 bbox( bb, &bbpts[0][0][0], trstride, tcstride, nrows, ncols ); 683 684 // find out if bounding box can't fit in unit cube 685 if( bbox_subdividing == N_BBOXROUND ) { 686 for( int k=0; k != inhcoords; k++ ) 687 if( ceilf(bb[1][k]) - floorf(bb[0][k]) > bboxsize[k] ) return 1; 688 } else { 689 for( int k=0; k != inhcoords; k++ ) 690 if( bb[1][k] - bb[0][k] > bboxsize[k] ) return 1; 691 } 692 return 0; 693 } 694 695 void 696 Mapdesc::bbox( 697 REAL bb[2][MAXCOORDS], 698 REAL *p, 699 int rstride, 700 int cstride, 701 int nrows, 702 int ncols ) 703 { 704 int k; 705 for( k=0; k != inhcoords; k++ ) 706 bb[0][k] = bb[1][k] = p[k]; 707 708 for( int i=0; i != nrows; i++ ) 709 for( int j=0; j != ncols; j++ ) 710 for( k=0; k != inhcoords; k++ ) { 711 REAL x = p[i*rstride + j*cstride + k]; 712 if( x < bb[0][k] ) bb[0][k] = x; 713 else if( x > bb[1][k] ) bb[1][k] = x; 714 } 715 } 716 717 /*-------------------------------------------------------------------------- 718 * calcVelocityRational - calculate upper bound on first partial derivative 719 * of a homogeneous set of points and bounds on each row of points. 720 *-------------------------------------------------------------------------- 721 */ 722 REAL 723 Mapdesc::calcVelocityRational( REAL *p, int stride, int ncols ) 724 { 725 REAL tmp[MAXORDER][MAXCOORDS]; 726 727 assert( ncols <= MAXORDER ); 728 729 const int tstride = sizeof(tmp[0]) / sizeof(REAL); 730 731 if( project( p, stride, &tmp[0][0], tstride, ncols ) ) { 732 return calcPartialVelocity( &tmp[0][0], tstride, ncols, 1, 1.0 ); 733 } else { /* XXX */ 734 return calcPartialVelocity( &tmp[0][0], tstride, ncols, 1, 1.0 ); 735 } 736 } 737 738 /*-------------------------------------------------------------------------- 739 * calcVelocityNonrational - calculate upper bound on first partial 740 * derivative of a inhomogeneous set of points. 741 *-------------------------------------------------------------------------- 742 */ 743 REAL 744 Mapdesc::calcVelocityNonrational( REAL *pts, int stride, int ncols ) 745 { 746 return calcPartialVelocity( pts, stride, ncols, 1, 1.0 ); 747 } 748 749 int 750 Mapdesc::isProperty( long property ) 751 { 752 switch ( property ) { 753 case N_PIXEL_TOLERANCE: 754 case N_ERROR_TOLERANCE: 755 case N_CULLING: 756 case N_BBOX_SUBDIVIDING: 757 case N_S_STEPS: 758 case N_T_STEPS: 759 case N_SAMPLINGMETHOD: 760 case N_CLAMPFACTOR: 761 case N_MINSAVINGS: 762 return 1; 763 default: 764 return 0; 765 } 766 } 767 768 REAL 769 Mapdesc::getProperty( long property ) 770 { 771 switch ( property ) { 772 case N_PIXEL_TOLERANCE: 773 return pixel_tolerance; 774 case N_ERROR_TOLERANCE: 775 return error_tolerance; 776 case N_CULLING: 777 return culling_method; 778 case N_BBOX_SUBDIVIDING: 779 return bbox_subdividing; 780 case N_S_STEPS: 781 return s_steps; 782 case N_T_STEPS: 783 return t_steps; 784 case N_SAMPLINGMETHOD: 785 return sampling_method; 786 case N_CLAMPFACTOR: 787 return clampfactor; 788 case N_MINSAVINGS: 789 return minsavings; 790 default: 791 abort(); 792 return -1; //not necessary, needed to shut up compiler 793 } 794 } 795 796 void 797 Mapdesc::setProperty( long property, REAL value ) 798 { 799 800 switch ( property ) { 801 case N_PIXEL_TOLERANCE: 802 pixel_tolerance = value; 803 break; 804 case N_ERROR_TOLERANCE: 805 error_tolerance = value; 806 break; 807 case N_CULLING: 808 culling_method = value; 809 break; 810 case N_BBOX_SUBDIVIDING: 811 if( value <= 0.0 ) value = N_NOBBOXSUBDIVISION; 812 bbox_subdividing = value; 813 break; 814 case N_S_STEPS: 815 if( value < 0.0 ) value = 0.0; 816 s_steps = value; 817 maxrate = ( value < 0.0 ) ? 0.0 : value; 818 maxsrate = ( value < 0.0 ) ? 0.0 : value; 819 break; 820 case N_T_STEPS: 821 if( value < 0.0 ) value = 0.0; 822 t_steps = value; 823 maxtrate = ( value < 0.0 ) ? 0.0 : value; 824 break; 825 case N_SAMPLINGMETHOD: 826 sampling_method = value; 827 break; 828 case N_CLAMPFACTOR: 829 if( value <= 0.0 ) value = N_NOCLAMPING; 830 clampfactor = value; 831 break; 832 case N_MINSAVINGS: 833 if( value <= 0.0 ) value = N_NOSAVINGSSUBDIVISION; 834 minsavings = value; 835 break; 836 default: 837 abort(); 838 break; 839 } 840 } 841 842